/[escript]/trunk/paso/src/Solver_ILU.c
ViewVC logotype

Contents of /trunk/paso/src/Solver_ILU.c

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1811 - (show annotations)
Thu Sep 25 23:11:13 2008 UTC (11 years, 5 months ago) by ksteube
File MIME type: text/plain
File size: 25731 byte(s)
Copyright updated in all files

1
2 /*******************************************************
3 *
4 * Copyright (c) 2003-2008 by University of Queensland
5 * Earth Systems Science Computational Center (ESSCC)
6 * http://www.uq.edu.au/esscc
7 *
8 * Primary Business: Queensland, Australia
9 * Licensed under the Open Software License version 3.0
10 * http://www.opensource.org/licenses/osl-3.0.php
11 *
12 *******************************************************/
13
14
15 /**************************************************************/
16
17 /* Paso: ILU preconditioner with reordering */
18
19 /**************************************************************/
20
21 /* Copyrights by ACcESS Australia 2003,2004,2005 */
22 /* Author: gross@access.edu.au */
23
24 /**************************************************************/
25
26 #include "Paso.h"
27 #include "Solver.h"
28 #include "PasoUtil.h"
29
30 /**************************************************************/
31
32 /* free all memory used by ILU */
33
34 void Paso_Solver_ILU_free(Paso_Solver_ILU * in) {
35 if (in!=NULL) {
36 MEMFREE(in->colorOf);
37 MEMFREE(in->factors);
38 MEMFREE(in->main_iptr);
39 Paso_Pattern_free(in->pattern);
40 MEMFREE(in);
41 }
42 }
43
44 /**************************************************************/
45
46 /* constructs the incomplete block factorization of
47
48 */
49 Paso_Solver_ILU* Paso_Solver_getILU(Paso_SparseMatrix * A,bool_t verbose) {
50 dim_t n=A->numRows;
51 dim_t n_block=A->row_block_size;
52 index_t num_colors=0, *mis_marker=NULL;
53 register double A11,A12,A13,A21,A22,A23,A31,A32,A33,D;
54 register double S11,S12,S13,S21,S22,S23,S31,S32,S33;
55 register index_t i,iptr_main,iptr,iptr_ik,k,iptr_kj,j,iptr_ij,color,color2;
56 double time0,time_color,time_fac;
57 /* allocations: */
58 Paso_Solver_ILU* out=MEMALLOC(1,Paso_Solver_ILU);
59 if (Paso_checkPtr(out)) return NULL;
60 out->colorOf=MEMALLOC(n,index_t);
61 out->factors=MEMALLOC(A->len,double);
62 out->main_iptr=MEMALLOC(n,index_t);
63 out->pattern=Paso_Pattern_getReference(A->pattern);
64 out->n_block=n_block;
65 out->n=n;
66 if ( !(Paso_checkPtr(out->colorOf) || Paso_checkPtr(out->main_iptr) || Paso_checkPtr(out->factors)) ) {
67 time0=Paso_timer();
68 Paso_Pattern_color(A->pattern,&out->num_colors,out->colorOf);
69 time_color=Paso_timer()-time0;
70
71 if (Paso_noError()) {
72 time0=Paso_timer();
73 /* find main diagonal and copy matrix values */
74 #pragma omp parallel for schedule(static) private(i,iptr,iptr_main,k)
75 for (i = 0; i < n; ++i) {
76 for (iptr=A->pattern->ptr[i];iptr<A->pattern->ptr[i+1]; ++iptr) {
77 iptr_main=A->pattern->ptr[0]-1;
78 for (iptr=A->pattern->ptr[i];iptr<A->pattern->ptr[i+1]; iptr++) {
79 if (A->pattern->index[iptr]==i) iptr_main=iptr;
80 for (k=0;k<n_block*n_block;++k) out->factors[n_block*n_block*iptr+k]=A->val[n_block*n_block*iptr+k];
81 }
82 out->main_iptr[i]=iptr_main;
83 if (iptr_main==A->pattern->ptr[0]-1)
84 Paso_setError(VALUE_ERROR, "Paso_Solver_getILU: no main diagonal");
85 }
86 }
87 /* start factorization */
88
89 #pragma omp barrier
90 for (color=0;color<out->num_colors && Paso_noError();++color) {
91 if (n_block==1) {
92 #pragma omp parallel for schedule(static) private(i,color2,iptr_ik,k,iptr_kj,S11,j,iptr_ij,A11,iptr_main,D)
93 for (i = 0; i < n; ++i) {
94 if (out->colorOf[i]==color) {
95 for (color2=0;color2<color;++color2) {
96 for (iptr_ik=A->pattern->ptr[i];iptr_ik<A->pattern->ptr[i+1]; ++iptr_ik) {
97 k=A->pattern->index[iptr_ik];
98 if (out->colorOf[k]==color2) {
99 A11=out->factors[iptr_ik];
100 /* a_ij=a_ij-a_ik*a_kj */
101 for (iptr_kj=A->pattern->ptr[k];iptr_kj<A->pattern->ptr[k+1]; iptr_kj++) {
102 j=A->pattern->index[iptr_kj];
103 if (out->colorOf[j]>color2) {
104 S11=out->factors[iptr_kj];
105 for (iptr_ij=A->pattern->ptr[i];iptr_ij<A->pattern->ptr[i+1]; iptr_ij++) {
106 if (j==A->pattern->index[iptr_ij]) {
107 out->factors[iptr_ij]-=A11*S11;
108 break;
109 }
110 }
111 }
112 }
113 }
114 }
115 }
116 iptr_main=out->main_iptr[i];
117 D=out->factors[iptr_main];
118 if (ABS(D)>0.) {
119 D=1./D;
120 out->factors[iptr_main]=D;
121 /* a_ik=a_ii^{-1}*a_ik */
122 for (iptr_ik=A->pattern->ptr[i];iptr_ik<A->pattern->ptr[i+1]; ++iptr_ik) {
123 k=A->pattern->index[iptr_ik];
124 if (out->colorOf[k]>color) {
125 A11=out->factors[iptr_ik];
126 out->factors[iptr_ik]=A11*D;
127 }
128 }
129 } else {
130 Paso_setError(ZERO_DIVISION_ERROR, "Paso_Solver_getILU: non-regular main diagonal block.");
131 }
132 }
133 }
134 } else if (n_block==2) {
135 #pragma omp parallel for schedule(static) private(i,color2,iptr_ik,k,iptr_kj,S11,S21,S12,S22,j,iptr_ij,A11,A21,A12,A22,iptr_main,D)
136 for (i = 0; i < n; ++i) {
137 if (out->colorOf[i]==color) {
138 for (color2=0;color2<color;++color2) {
139 for (iptr_ik=A->pattern->ptr[i];iptr_ik<A->pattern->ptr[i+1]; ++iptr_ik) {
140 k=A->pattern->index[iptr_ik];
141 if (out->colorOf[k]==color2) {
142 A11=out->factors[iptr_ik*4 ];
143 A21=out->factors[iptr_ik*4+1];
144 A12=out->factors[iptr_ik*4+2];
145 A22=out->factors[iptr_ik*4+3];
146 /* a_ij=a_ij-a_ik*a_kj */
147 for (iptr_kj=A->pattern->ptr[k];iptr_kj<A->pattern->ptr[k+1]; iptr_kj++) {
148 j=A->pattern->index[iptr_kj];
149 if (out->colorOf[j]>color2) {
150 S11=out->factors[iptr_kj*4];
151 S21=out->factors[iptr_kj*4+1];
152 S12=out->factors[iptr_kj*4+2];
153 S22=out->factors[iptr_kj*4+3];
154 for (iptr_ij=A->pattern->ptr[i];iptr_ij<A->pattern->ptr[i+1]; iptr_ij++) {
155 if (j==A->pattern->index[iptr_ij]) {
156 out->factors[4*iptr_ij ]-=A11*S11+A12*S21;
157 out->factors[4*iptr_ij+1]-=A21*S11+A22*S21;
158 out->factors[4*iptr_ij+2]-=A11*S12+A12*S22;
159 out->factors[4*iptr_ij+3]-=A21*S12+A22*S22;
160 break;
161 }
162 }
163 }
164 }
165 }
166 }
167 }
168 iptr_main=out->main_iptr[i];
169 A11=out->factors[iptr_main*4];
170 A21=out->factors[iptr_main*4+1];
171 A12=out->factors[iptr_main*4+2];
172 A22=out->factors[iptr_main*4+3];
173 D = A11*A22-A12*A21;
174 if (ABS(D)>0.) {
175 D=1./D;
176 S11= A22*D;
177 S21=-A21*D;
178 S12=-A12*D;
179 S22= A11*D;
180 out->factors[iptr_main*4] = S11;
181 out->factors[iptr_main*4+1]= S21;
182 out->factors[iptr_main*4+2]= S12;
183 out->factors[iptr_main*4+3]= S22;
184 /* a_ik=a_ii^{-1}*a_ik */
185 for (iptr_ik=A->pattern->ptr[i];iptr_ik<A->pattern->ptr[i+1]; ++iptr_ik) {
186 k=A->pattern->index[iptr_ik];
187 if (out->colorOf[k]>color) {
188 A11=out->factors[iptr_ik*4 ];
189 A21=out->factors[iptr_ik*4+1];
190 A12=out->factors[iptr_ik*4+2];
191 A22=out->factors[iptr_ik*4+3];
192 out->factors[4*iptr_ik ]=S11*A11+S12*A21;
193 out->factors[4*iptr_ik+1]=S21*A11+S22*A21;
194 out->factors[4*iptr_ik+2]=S11*A12+S12*A22;
195 out->factors[4*iptr_ik+3]=S21*A12+S22*A22;
196 }
197 }
198 } else {
199 Paso_setError(ZERO_DIVISION_ERROR, "Paso_Solver_getILU: non-regular main diagonal block.");
200 }
201 }
202 }
203 } else if (n_block==3) {
204 #pragma omp parallel for schedule(static) private(i,color2,iptr_ik,k,iptr_kj,S11,S21,S31,S12,S22,S32,S13,S23,S33,j,iptr_ij,A11,A21,A31,A12,A22,A32,A13,A23,A33,iptr_main,D)
205 for (i = 0; i < n; ++i) {
206 if (out->colorOf[i]==color) {
207 for (color2=0;color2<color;++color2) {
208 for (iptr_ik=A->pattern->ptr[i];iptr_ik<A->pattern->ptr[i+1]; ++iptr_ik) {
209 k=A->pattern->index[iptr_ik];
210 if (out->colorOf[k]==color2) {
211 A11=out->factors[iptr_ik*9 ];
212 A21=out->factors[iptr_ik*9+1];
213 A31=out->factors[iptr_ik*9+2];
214 A12=out->factors[iptr_ik*9+3];
215 A22=out->factors[iptr_ik*9+4];
216 A32=out->factors[iptr_ik*9+5];
217 A13=out->factors[iptr_ik*9+6];
218 A23=out->factors[iptr_ik*9+7];
219 A33=out->factors[iptr_ik*9+8];
220 /* a_ij=a_ij-a_ik*a_kj */
221 for (iptr_kj=A->pattern->ptr[k];iptr_kj<A->pattern->ptr[k+1]; iptr_kj++) {
222 j=A->pattern->index[iptr_kj];
223 if (out->colorOf[j]>color2) {
224 S11=out->factors[iptr_kj*9 ];
225 S21=out->factors[iptr_kj*9+1];
226 S31=out->factors[iptr_kj*9+2];
227 S12=out->factors[iptr_kj*9+3];
228 S22=out->factors[iptr_kj*9+4];
229 S32=out->factors[iptr_kj*9+5];
230 S13=out->factors[iptr_kj*9+6];
231 S23=out->factors[iptr_kj*9+7];
232 S33=out->factors[iptr_kj*9+8];
233 for (iptr_ij=A->pattern->ptr[i];iptr_ij<A->pattern->ptr[i+1]; iptr_ij++) {
234 if (j==A->pattern->index[iptr_ij]) {
235 out->factors[iptr_ij*9 ]-=A11*S11+A12*S21+A13*S31;
236 out->factors[iptr_ij*9+1]-=A21*S11+A22*S21+A23*S31;
237 out->factors[iptr_ij*9+2]-=A31*S11+A32*S21+A33*S31;
238 out->factors[iptr_ij*9+3]-=A11*S12+A12*S22+A13*S32;
239 out->factors[iptr_ij*9+4]-=A21*S12+A22*S22+A23*S32;
240 out->factors[iptr_ij*9+5]-=A31*S12+A32*S22+A33*S32;
241 out->factors[iptr_ij*9+6]-=A11*S13+A12*S23+A13*S33;
242 out->factors[iptr_ij*9+7]-=A21*S13+A22*S23+A23*S33;
243 out->factors[iptr_ij*9+8]-=A31*S13+A32*S23+A33*S33;
244 break;
245 }
246 }
247 }
248 }
249 }
250 }
251 }
252 iptr_main=out->main_iptr[i];
253 A11=out->factors[iptr_main*9 ];
254 A21=out->factors[iptr_main*9+1];
255 A31=out->factors[iptr_main*9+2];
256 A12=out->factors[iptr_main*9+3];
257 A22=out->factors[iptr_main*9+4];
258 A32=out->factors[iptr_main*9+5];
259 A13=out->factors[iptr_main*9+6];
260 A23=out->factors[iptr_main*9+7];
261 A33=out->factors[iptr_main*9+8];
262 D = A11*(A22*A33-A23*A32)+ A12*(A31*A23-A21*A33)+A13*(A21*A32-A31*A22);
263 if (ABS(D)>0.) {
264 D=1./D;
265 S11=(A22*A33-A23*A32)*D;
266 S21=(A31*A23-A21*A33)*D;
267 S31=(A21*A32-A31*A22)*D;
268 S12=(A13*A32-A12*A33)*D;
269 S22=(A11*A33-A31*A13)*D;
270 S32=(A12*A31-A11*A32)*D;
271 S13=(A12*A23-A13*A22)*D;
272 S23=(A13*A21-A11*A23)*D;
273 S33=(A11*A22-A12*A21)*D;
274
275 out->factors[iptr_main*9 ]=S11;
276 out->factors[iptr_main*9+1]=S21;
277 out->factors[iptr_main*9+2]=S31;
278 out->factors[iptr_main*9+3]=S12;
279 out->factors[iptr_main*9+4]=S22;
280 out->factors[iptr_main*9+5]=S32;
281 out->factors[iptr_main*9+6]=S13;
282 out->factors[iptr_main*9+7]=S23;
283 out->factors[iptr_main*9+8]=S33;
284
285 /* a_ik=a_ii^{-1}*a_ik */
286 for (iptr_ik=A->pattern->ptr[i];iptr_ik<A->pattern->ptr[i+1]; ++iptr_ik) {
287 k=A->pattern->index[iptr_ik];
288 if (out->colorOf[k]>color) {
289 A11=out->factors[iptr_ik*9 ];
290 A21=out->factors[iptr_ik*9+1];
291 A31=out->factors[iptr_ik*9+2];
292 A12=out->factors[iptr_ik*9+3];
293 A22=out->factors[iptr_ik*9+4];
294 A32=out->factors[iptr_ik*9+5];
295 A13=out->factors[iptr_ik*9+6];
296 A23=out->factors[iptr_ik*9+7];
297 A33=out->factors[iptr_ik*9+8];
298 out->factors[iptr_ik*9 ]=S11*A11+S12*A21+S13*A31;
299 out->factors[iptr_ik*9+1]=S21*A11+S22*A21+S23*A31;
300 out->factors[iptr_ik*9+2]=S31*A11+S32*A21+S33*A31;
301 out->factors[iptr_ik*9+3]=S11*A12+S12*A22+S13*A32;
302 out->factors[iptr_ik*9+4]=S21*A12+S22*A22+S23*A32;
303 out->factors[iptr_ik*9+5]=S31*A12+S32*A22+S33*A32;
304 out->factors[iptr_ik*9+6]=S11*A13+S12*A23+S13*A33;
305 out->factors[iptr_ik*9+7]=S21*A13+S22*A23+S23*A33;
306 out->factors[iptr_ik*9+8]=S31*A13+S32*A23+S33*A33;
307 }
308 }
309 } else {
310 Paso_setError(ZERO_DIVISION_ERROR, "Paso_Solver_getILU: non-regular main diagonal block.");
311 }
312 }
313 }
314 } else {
315 Paso_setError(VALUE_ERROR, "Paso_Solver_getILU: block size greater than 3 is not supported.");
316 }
317 #pragma omp barrier
318 }
319 time_fac=Paso_timer()-time0;
320 }
321 }
322 if (Paso_noError()) {
323 if (verbose) {
324 printf("ILU: %d color used \n",out->num_colors);
325 printf("timing: ILU: coloring/elemination : %e/%e\n",time_color,time_fac);
326 }
327 return out;
328 } else {
329 Paso_Solver_ILU_free(out);
330 return NULL;
331 }
332 }
333
334 /**************************************************************/
335
336 /* apply ILU precondition b-> x
337
338 in fact it solves LUx=b in the form x= U^{-1} L^{-1}b
339
340 should be called within a parallel region
341 barrier synconization should be performed to make sure that the input vector available
342
343 */
344
345 void Paso_Solver_solveILU(Paso_Solver_ILU * ilu, double * x, double * b) {
346 register dim_t i,k;
347 register index_t color,iptr_ik,iptr_main;
348 register double S1,S2,S3,R1,R2,R3;
349 dim_t n_block=ilu->n_block;
350 dim_t n=ilu->n;
351
352
353 /* copy x into b*/
354 #pragma omp parallel for private(i) schedule(static)
355 for (i=0;i<n*n_block;++i) x[i]=b[i];
356 /* forward substitution */
357 for (color=0;color<ilu->num_colors;++color) {
358 if (n_block==1) {
359 #pragma omp parallel for schedule(static) private(i,iptr_ik,k,S1,R1,iptr_main)
360 for (i = 0; i < n; ++i) {
361 if (ilu->colorOf[i]==color) {
362 /* x_i=x_i-a_ik*x_k */
363 S1=x[i];
364 for (iptr_ik=ilu->pattern->ptr[i];iptr_ik<ilu->pattern->ptr[i+1]; ++iptr_ik) {
365 k=ilu->pattern->index[iptr_ik];
366 if (ilu->colorOf[k]<color) {
367 R1=x[k];
368 S1-=ilu->factors[iptr_ik]*R1;
369 }
370 }
371 iptr_main=ilu->main_iptr[i];
372 x[i]=ilu->factors[iptr_main]*S1;
373 }
374 }
375 } else if (n_block==2) {
376 #pragma omp parallel for schedule(static) private(i,iptr_ik,k,iptr_main,S1,S2,R1,R2)
377 for (i = 0; i < n; ++i) {
378 if (ilu->colorOf[i]==color) {
379 /* x_i=x_i-a_ik*x_k */
380 S1=x[2*i];
381 S2=x[2*i+1];
382 for (iptr_ik=ilu->pattern->ptr[i];iptr_ik<ilu->pattern->ptr[i+1]; ++iptr_ik) {
383 k=ilu->pattern->index[iptr_ik];
384 if (ilu->colorOf[k]<color) {
385 R1=x[2*k];
386 R2=x[2*k+1];
387 S1-=ilu->factors[4*iptr_ik ]*R1+ilu->factors[4*iptr_ik+2]*R2;
388 S2-=ilu->factors[4*iptr_ik+1]*R1+ilu->factors[4*iptr_ik+3]*R2;
389 }
390 }
391 iptr_main=ilu->main_iptr[i];
392 x[2*i ]=ilu->factors[4*iptr_main ]*S1+ilu->factors[4*iptr_main+2]*S2;
393 x[2*i+1]=ilu->factors[4*iptr_main+1]*S1+ilu->factors[4*iptr_main+3]*S2;
394 }
395
396 }
397 } else if (n_block==3) {
398 #pragma omp parallel for schedule(static) private(i,iptr_ik,iptr_main,k,S1,S2,S3,R1,R2,R3)
399 for (i = 0; i < n; ++i) {
400 if (ilu->colorOf[i]==color) {
401 /* x_i=x_i-a_ik*x_k */
402 S1=x[3*i];
403 S2=x[3*i+1];
404 S3=x[3*i+2];
405 for (iptr_ik=ilu->pattern->ptr[i];iptr_ik<ilu->pattern->ptr[i+1]; ++iptr_ik) {
406 k=ilu->pattern->index[iptr_ik];
407 if (ilu->colorOf[k]<color) {
408 R1=x[3*k];
409 R2=x[3*k+1];
410 R3=x[3*k+2];
411 S1-=ilu->factors[9*iptr_ik ]*R1+ilu->factors[9*iptr_ik+3]*R2+ilu->factors[9*iptr_ik+6]*R3;
412 S2-=ilu->factors[9*iptr_ik+1]*R1+ilu->factors[9*iptr_ik+4]*R2+ilu->factors[9*iptr_ik+7]*R3;
413 S3-=ilu->factors[9*iptr_ik+2]*R1+ilu->factors[9*iptr_ik+5]*R2+ilu->factors[9*iptr_ik+8]*R3;
414 }
415 }
416 iptr_main=ilu->main_iptr[i];
417 x[3*i ]=ilu->factors[9*iptr_main ]*S1+ilu->factors[9*iptr_main+3]*S2+ilu->factors[9*iptr_main+6]*S3;
418 x[3*i+1]=ilu->factors[9*iptr_main+1]*S1+ilu->factors[9*iptr_main+4]*S2+ilu->factors[9*iptr_main+7]*S3;
419 x[3*i+2]=ilu->factors[9*iptr_main+2]*S1+ilu->factors[9*iptr_main+5]*S2+ilu->factors[9*iptr_main+8]*S3;
420 }
421 }
422 }
423 }
424 /* backward substitution */
425 for (color=(ilu->num_colors)-1;color>-1;--color) {
426 if (n_block==1) {
427 #pragma omp parallel for schedule(static) private(i,iptr_ik,k,S1,R1)
428 for (i = 0; i < n; ++i) {
429 if (ilu->colorOf[i]==color) {
430 /* x_i=x_i-a_ik*x_k */
431 S1=x[i];
432 for (iptr_ik=ilu->pattern->ptr[i];iptr_ik<ilu->pattern->ptr[i+1]; ++iptr_ik) {
433 k=ilu->pattern->index[iptr_ik];
434 if (ilu->colorOf[k]>color) {
435 R1=x[k];
436 S1-=ilu->factors[iptr_ik]*R1;
437 }
438 }
439 x[i]=S1;
440 }
441 }
442 } else if (n_block==2) {
443 #pragma omp parallel for schedule(static) private(i,iptr_ik,k,S1,S2,R1,R2)
444 for (i = 0; i < n; ++i) {
445 if (ilu->colorOf[i]==color) {
446 /* x_i=x_i-a_ik*x_k */
447 S1=x[2*i];
448 S2=x[2*i+1];
449 for (iptr_ik=ilu->pattern->ptr[i];iptr_ik<ilu->pattern->ptr[i+1]; ++iptr_ik) {
450 k=ilu->pattern->index[iptr_ik];
451 if (ilu->colorOf[k]>color) {
452 R1=x[2*k];
453 R2=x[2*k+1];
454 S1-=ilu->factors[4*iptr_ik ]*R1+ilu->factors[4*iptr_ik+2]*R2;
455 S2-=ilu->factors[4*iptr_ik+1]*R1+ilu->factors[4*iptr_ik+3]*R2;
456 }
457 }
458 x[2*i]=S1;
459 x[2*i+1]=S2;
460 }
461 }
462 } else if (n_block==3) {
463 #pragma omp parallel for schedule(static) private(i,iptr_ik,k,S1,S2,S3,R1,R2,R3)
464 for (i = 0; i < n; ++i) {
465 if (ilu->colorOf[i]==color) {
466 /* x_i=x_i-a_ik*x_k */
467 S1=x[3*i ];
468 S2=x[3*i+1];
469 S3=x[3*i+2];
470 for (iptr_ik=ilu->pattern->ptr[i];iptr_ik<ilu->pattern->ptr[i+1]; ++iptr_ik) {
471 k=ilu->pattern->index[iptr_ik];
472 if (ilu->colorOf[k]>color) {
473 R1=x[3*k];
474 R2=x[3*k+1];
475 R3=x[3*k+2];
476 S1-=ilu->factors[9*iptr_ik ]*R1+ilu->factors[9*iptr_ik+3]*R2+ilu->factors[9*iptr_ik+6]*R3;
477 S2-=ilu->factors[9*iptr_ik+1]*R1+ilu->factors[9*iptr_ik+4]*R2+ilu->factors[9*iptr_ik+7]*R3;
478 S3-=ilu->factors[9*iptr_ik+2]*R1+ilu->factors[9*iptr_ik+5]*R2+ilu->factors[9*iptr_ik+8]*R3;
479 }
480 }
481 x[3*i]=S1;
482 x[3*i+1]=S2;
483 x[3*i+2]=S3;
484 }
485 }
486 }
487 #pragma omp barrier
488 }
489 return;
490 }
491 /*
492 * $Log$
493 * Revision 1.2 2005/09/15 03:44:40 jgs
494 * Merge of development branch dev-02 back to main trunk on 2005-09-15
495 *
496 * Revision 1.1.2.1 2005/09/05 06:29:50 gross
497 * These files have been extracted from finley to define a stand alone libray for iterative
498 * linear solvers on the ALTIX. main entry through Paso_solve. this version compiles but
499 * has not been tested yet.
500 *
501 *
502 */

Properties

Name Value
svn:eol-style native
svn:keywords Author Date Id Revision

  ViewVC Help
Powered by ViewVC 1.1.26