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

Annotation of /trunk/paso/src/Smoother.c

Parent Directory Parent Directory | Revision Log Revision Log


Revision 3402 - (hide annotations)
Tue Dec 7 07:36:12 2010 UTC (8 years, 10 months ago) by gross
File MIME type: text/plain
File size: 15431 byte(s)
AMG support now block size > 3 (requires clapack or MKL).

additional diagnostics for AMG 


1 jfenwick 1851
2     /*******************************************************
3     *
4 jfenwick 2881 * Copyright (c) 2003-2010 by University of Queensland
5 jfenwick 1851 * 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 gross 3097 /* Paso: Gauss-Seidel */
18 jfenwick 1851
19     /**************************************************************/
20    
21     /* Author: artak@uq.edu.au */
22    
23     /**************************************************************/
24    
25     #include "Paso.h"
26 gross 3120 #include "Preconditioner.h"
27 jfenwick 1851 #include "PasoUtil.h"
28 gross 3097 #include "BlockOps.h"
29 jfenwick 1851
30     #include <stdio.h>
31    
32 gross 3094
33 jfenwick 1851 /**************************************************************/
34    
35 gross 3158 /* free all memory used by Smoother */
36 jfenwick 1851
37 gross 3158 void Paso_Preconditioner_Smoother_free(Paso_Preconditioner_Smoother * in) {
38 jfenwick 1851 if (in!=NULL) {
39 gross 3158 Paso_Preconditioner_LocalSmoother_free(in->localSmoother);
40 jfenwick 1851 MEMFREE(in);
41     }
42     }
43 gross 3158 void Paso_Preconditioner_LocalSmoother_free(Paso_Preconditioner_LocalSmoother * in) {
44 gross 3094 if (in!=NULL) {
45     MEMFREE(in->diag);
46     MEMFREE(in->pivot);
47 gross 3097 MEMFREE(in->buffer);
48 gross 3094 MEMFREE(in);
49     }
50     }
51 jfenwick 1851 /**************************************************************/
52    
53 gross 3094 /* constructs the symmetric Gauss-Seidel preconditioner
54 jfenwick 1851
55     */
56 gross 3158 Paso_Preconditioner_Smoother* Paso_Preconditioner_Smoother_alloc(Paso_SystemMatrix * A_p, const bool_t jacobi, const bool_t is_local, const bool_t verbose)
57 gross 3094 {
58    
59 jfenwick 1851 /* allocations: */
60 gross 3158 Paso_Preconditioner_Smoother* out=MEMALLOC(1,Paso_Preconditioner_Smoother);
61 jfenwick 3259 if (! Esys_checkPtr(out)) {
62 gross 3158 out->localSmoother=Paso_Preconditioner_LocalSmoother_alloc(A_p->mainBlock,jacobi,verbose);
63 gross 3094 out->is_local=is_local;
64 jfenwick 1851 }
65 jfenwick 3259 if (Esys_MPIInfo_noError(A_p->mpi_info)) {
66 jfenwick 1851 return out;
67 gross 3094 } else {
68 gross 3158 Paso_Preconditioner_Smoother_free(out);
69 jfenwick 1851 return NULL;
70     }
71     }
72 jfenwick 3259 Paso_Preconditioner_LocalSmoother* Paso_Preconditioner_LocalSmoother_alloc(Paso_SparseMatrix * A_p, const bool_t jacobi, bool_t verbose)
73 gross 3158 {
74 gross 3094
75     dim_t n=A_p->numRows;
76     dim_t n_block=A_p->row_block_size;
77     dim_t block_size=A_p->block_size;
78    
79 jfenwick 3259 double time0=Esys_timer();
80 gross 3094 /* allocations: */
81 gross 3158 Paso_Preconditioner_LocalSmoother* out=MEMALLOC(1,Paso_Preconditioner_LocalSmoother);
82 jfenwick 3259 if (! Esys_checkPtr(out)) {
83 gross 3094
84     out->diag=MEMALLOC( ((size_t) n) * ((size_t) block_size),double);
85     out->pivot=MEMALLOC( ((size_t) n) * ((size_t) n_block), index_t);
86 gross 3097 out->buffer=MEMALLOC( ((size_t) n) * ((size_t) n_block), double);
87 gross 3158 out->Jacobi=jacobi;
88 gross 3094
89 jfenwick 3259 if ( ! ( Esys_checkPtr(out->diag) || Esys_checkPtr(out->pivot) ) ) {
90 gross 3094 Paso_SparseMatrix_invMain(A_p, out->diag, out->pivot );
91     }
92    
93     }
94 jfenwick 3259 time0=Esys_timer()-time0;
95 gross 3094
96 jfenwick 3259 if (Esys_noError()) {
97 gross 3094 return out;
98     } else {
99 gross 3158 Paso_Preconditioner_LocalSmoother_free(out);
100 gross 3094 return NULL;
101     }
102     }
103 jfenwick 1851
104 gross 3097 /*
105 jfenwick 1851
106 gross 3097 performs a few sweeps of the from
107 jfenwick 1851
108 gross 3097 S (x_{k} - x_{k-1}) = b - A x_{k-1}
109 jfenwick 1851
110 gross 3097 where x_{0}=0 and S provides some approximatioon of A.
111 jfenwick 1851
112 gross 3097 Under MPI S is build on using A_p->mainBlock only.
113 gross 3158 if Smoother is local the defect b - A x_{k-1} is calculated using A_p->mainBlock only.
114 gross 3094
115 jfenwick 1851 */
116    
117 gross 3159 void Paso_Preconditioner_Smoother_solve(Paso_SystemMatrix* A_p, Paso_Preconditioner_Smoother * smoother, double * x, const double * b,
118     const dim_t sweeps, const bool_t x_is_initial)
119 gross 3094 {
120 gross 3097 const dim_t n= (A_p->mainBlock->numRows) * (A_p->mainBlock->row_block_size);
121    
122 gross 3158 double *b_new = smoother->localSmoother->buffer;
123     dim_t nsweeps=sweeps;
124     if (smoother->is_local) {
125 gross 3159 Paso_Preconditioner_LocalSmoother_solve(A_p->mainBlock,smoother->localSmoother,x,b,sweeps,x_is_initial);
126 gross 3097 } else {
127 gross 3159 if (! x_is_initial) {
128     Paso_Copy(n, x, b);
129     Paso_Preconditioner_LocalSmoother_Sweep(A_p->mainBlock,smoother->localSmoother,x);
130     nsweeps--;
131     }
132     while (nsweeps > 0 ) {
133 gross 3158 Paso_Copy(n, b_new, b);
134 gross 3097 Paso_SystemMatrix_MatrixVector((-1.), A_p, x, 1., b_new); /* b_new = b - A*x */
135 gross 3159 Paso_Preconditioner_LocalSmoother_Sweep(A_p->mainBlock,smoother->localSmoother,b_new);
136 gross 3158 Paso_AXPY(n, x, 1., b_new);
137     nsweeps--;
138 gross 3097 }
139    
140 gross 3094 }
141 gross 3097 }
142 gross 3159 void Paso_Preconditioner_LocalSmoother_solve(Paso_SparseMatrix* A_p, Paso_Preconditioner_LocalSmoother * smoother, double * x, const double * b,
143     const dim_t sweeps, const bool_t x_is_initial)
144 gross 3097 {
145     const dim_t n= (A_p->numRows) * (A_p->row_block_size);
146 gross 3158 double *b_new = smoother->buffer;
147     dim_t nsweeps=sweeps;
148 gross 3094
149 gross 3159 if (! x_is_initial) {
150     Paso_Copy(n, x, b);
151     Paso_Preconditioner_LocalSmoother_Sweep(A_p,smoother,x);
152     nsweeps--;
153     }
154 gross 3105
155 gross 3159 while (nsweeps > 0 ) {
156 gross 3158 Paso_Copy(n, b_new, b);
157 gross 3283
158     Paso_SparseMatrix_MatrixVector_CSR_OFFSET0((-1.), A_p, x, 1., b_new); /* b_new = b - A*x */
159 gross 3158 Paso_Preconditioner_LocalSmoother_Sweep(A_p,smoother,b_new);
160 gross 3159 Paso_AXPY(n, x, 1., b_new);
161 gross 3158 nsweeps--;
162 gross 3094 }
163     }
164 jfenwick 1851
165 gross 3158 void Paso_Preconditioner_LocalSmoother_Sweep(Paso_SparseMatrix* A, Paso_Preconditioner_LocalSmoother * smoother, double * x)
166 gross 3094 {
167     const dim_t nt=omp_get_max_threads();
168 gross 3158 if (smoother->Jacobi) {
169 gross 3402 Paso_BlockOps_solveAll(A->row_block_size,A->numRows,smoother->diag,smoother->pivot,x);
170 gross 3097 } else {
171 gross 3158 if (nt < 2) {
172     Paso_Preconditioner_LocalSmoother_Sweep_sequential(A,smoother,x);
173     } else {
174     Paso_Preconditioner_LocalSmoother_Sweep_colored(A,smoother,x);
175     }
176 gross 3094 }
177 jfenwick 1851 }
178    
179 gross 3097 /* inplace Gauss-Seidel sweep in seqential mode: */
180 gross 3094
181 gross 3158 void Paso_Preconditioner_LocalSmoother_Sweep_sequential(Paso_SparseMatrix* A_p, Paso_Preconditioner_LocalSmoother * smoother, double * x)
182 gross 3097 {
183     const dim_t n=A_p->numRows;
184     const dim_t n_block=A_p->row_block_size;
185 gross 3158 const double *diag = smoother->diag;
186 gross 3402 const index_t* pivot = smoother->pivot;
187     const dim_t block_len=A_p->block_size;
188 gross 3094
189 gross 3097 register dim_t i,k;
190     register index_t iptr_ik, mm;
191 gross 3402 register double rtmp;
192     int failed = 0;
193 gross 3097
194     const index_t* ptr_main = Paso_SparseMatrix_borrowMainDiagonalPointer(A_p);
195     /* forward substitution */
196    
197     if (n_block==1) {
198 gross 3379 x[0]*=diag[0];
199 gross 3098 for (i = 1; i < n; ++i) {
200 gross 3097 mm=ptr_main[i];
201 gross 3105 /* x_i=x_i-a_ik*x_k (with k<i) */
202 gross 3402 rtmp=x[i];
203 gross 3097 for (iptr_ik=A_p->pattern->ptr[i];iptr_ik<mm; ++iptr_ik) {
204 gross 3379 k=A_p->pattern->index[iptr_ik];
205 gross 3402 rtmp-=A_p->val[iptr_ik]*x[k];
206 gross 3097 }
207 gross 3402 x[i]=rtmp*diag[i];
208 gross 3097 }
209     } else if (n_block==2) {
210 gross 3402 Paso_BlockOps_MViP_2(&diag[0], &x[0]);
211 gross 3098 for (i = 1; i < n; ++i) {
212 gross 3097 mm=ptr_main[i];
213     for (iptr_ik=A_p->pattern->ptr[i];iptr_ik<mm; ++iptr_ik) {
214     k=A_p->pattern->index[iptr_ik];
215     Paso_BlockOps_SMV_2(&x[2*i], &A_p->val[4*iptr_ik], &x[2*k]);
216     }
217 gross 3402 Paso_BlockOps_MViP_2(&diag[4*i], &x[2*i]);
218 gross 3097 }
219     } else if (n_block==3) {
220 gross 3402 Paso_BlockOps_MViP_3(&diag[0], &x[0]);
221 gross 3098 for (i = 1; i < n; ++i) {
222 gross 3097 mm=ptr_main[i];
223     /* x_i=x_i-a_ik*x_k */
224     for (iptr_ik=A_p->pattern->ptr[i];iptr_ik<mm; ++iptr_ik) {
225     k=A_p->pattern->index[iptr_ik];
226     Paso_BlockOps_SMV_3(&x[3*i], &A_p->val[9*iptr_ik], &x[3*k]);
227     }
228 gross 3402 Paso_BlockOps_MViP_3(&diag[9*i], &x[3*i]);
229 gross 3097 }
230 gross 3402 } else {
231     Paso_BlockOps_solve_N(n_block, &x[0], &diag[0], &pivot[0], &failed);
232     for (i = 1; i < n; ++i) {
233     mm=ptr_main[i];
234     /* x_i=x_i-a_ik*x_k */
235     for (iptr_ik=A_p->pattern->ptr[i];iptr_ik<mm; ++iptr_ik) {
236     k=A_p->pattern->index[iptr_ik];
237     Paso_BlockOps_SMV_N(n_block, &x[n_block*i], &A_p->val[block_len*iptr_ik], &x[n_block*k]);
238     }
239     Paso_BlockOps_solve_N(n_block, &x[n_block*i], &diag[block_len*i], &pivot[n_block*i], &failed)
240     }
241     }
242 gross 3097 /* backward substitution */
243    
244     if (n_block==1) {
245 gross 3105 for (i = n-2; i > -1; --i) {
246 gross 3098 mm=ptr_main[i];
247 gross 3402 rtmp=x[i]*A_p->val[mm];
248 gross 3098 for (iptr_ik=mm+1; iptr_ik < A_p->pattern->ptr[i+1]; ++iptr_ik) {
249 gross 3402 k=A_p->pattern->index[iptr_ik];
250     rtmp-=A_p->val[iptr_ik]*x[k];
251 gross 3098 }
252 gross 3402 x[i]=diag[i]*rtmp;
253 gross 3097 }
254 gross 3098
255 gross 3097 } else if (n_block==2) {
256 gross 3105 for (i = n-2; i > -1; --i) {
257 gross 3098 mm=ptr_main[i];
258 gross 3402 Paso_BlockOps_MViP_2(&A_p->val[4*mm], &x[2*i]);
259 gross 3098 for (iptr_ik=mm+1; iptr_ik < A_p->pattern->ptr[i+1]; ++iptr_ik) {
260     k=A_p->pattern->index[iptr_ik];
261     Paso_BlockOps_SMV_2(&x[2*i], &A_p->val[4*iptr_ik], &x[2*k]);
262     }
263 gross 3402 Paso_BlockOps_MViP_2(&diag[i*4], &x[2*i]);
264 gross 3283
265 gross 3097 }
266     } else if (n_block==3) {
267 gross 3105 for (i = n-2; i > -1; --i) {
268 gross 3098 mm=ptr_main[i];
269 gross 3402 Paso_BlockOps_MViP_3(&A_p->val[9*mm], &x[3*i]);
270 gross 3098 for (iptr_ik=mm+1; iptr_ik < A_p->pattern->ptr[i+1]; ++iptr_ik) {
271     k=A_p->pattern->index[iptr_ik];
272     Paso_BlockOps_SMV_3(&x[3*i], &A_p->val[9*iptr_ik], &x[3*k]);
273     }
274 gross 3402 Paso_BlockOps_MViP_3(&diag[i*9], &x[3*i]);
275 gross 3097 }
276    
277 gross 3402 } else {
278     double *y=TMPMEMALLOC(n_block, double);
279     for (i = n-2; i > -1; --i) {
280     mm=ptr_main[i];
281     Paso_BlockOps_MV_N(n_block, &y[0], &A_p->val[block_len*mm], &x[n_block*i]);
282     for (iptr_ik=mm+1; iptr_ik < A_p->pattern->ptr[i+1]; ++iptr_ik) {
283     k=A_p->pattern->index[iptr_ik];
284     Paso_BlockOps_SMV_N(n_block, &y[0], &A_p->val[block_len*iptr_ik], &x[n_block*k]);
285     }
286     Paso_BlockOps_Cpy_N(n_block ,&x[n_block*i], &y[0]);
287     Paso_BlockOps_solve_N(n_block, &x[n_block*i], &diag[i*block_len], &pivot[i*n_block], &failed);
288     }
289     TMPMEMFREE(y);
290     }
291 gross 3097
292 gross 3402 if (failed > 0) {
293     Esys_setError(ZERO_DIVISION_ERROR, "Paso_Preconditioner_LocalSmoother_Sweep_sequential: non-regular main diagonal block.");
294     }
295    
296 gross 3097 return;
297     }
298 gross 3094
299 gross 3158 void Paso_Preconditioner_LocalSmoother_Sweep_colored(Paso_SparseMatrix* A_p, Paso_Preconditioner_LocalSmoother * smoother, double * x)
300 gross 3094 {
301     const dim_t n=A_p->numRows;
302     const dim_t n_block=A_p->row_block_size;
303 gross 3158 const double *diag = smoother->diag;
304     index_t* pivot = smoother->pivot;
305 gross 3402 const dim_t block_len=A_p->block_size;
306     double *y;
307 gross 3094
308     register dim_t i,k;
309     register index_t color,iptr_ik, mm;
310 gross 3402 register double rtmp;
311     int failed = 0;
312 gross 3094
313     const index_t* coloring = Paso_Pattern_borrowColoringPointer(A_p->pattern);
314     const dim_t num_colors = Paso_Pattern_getNumColors(A_p->pattern);
315     const index_t* ptr_main = Paso_SparseMatrix_borrowMainDiagonalPointer(A_p);
316 gross 3105
317 gross 3402 #pragma omp parallel private(mm, i,iptr_ik,k,rtmp, color, y)
318     {
319     if (n_block>3) {
320     y=TMPMEMALLOC(n_block, double);
321     } else {
322     y=NULL;
323 gross 3097 }
324 gross 3402 /* forward substitution */
325    
326     /* color = 0 */
327     if (n_block==1) {
328     #pragma omp for schedule(static)
329 gross 3094 for (i = 0; i < n; ++i) {
330 gross 3402 if (coloring[i]== 0 ) x[i]*=diag[i];
331 gross 3094 }
332 gross 3402 } else if (n_block==2) {
333     #pragma omp for schedule(static)
334 gross 3094 for (i = 0; i < n; ++i) {
335 gross 3402 if (coloring[i]== 0 ) Paso_BlockOps_MViP_2(&diag[i*4], &x[2*i]);
336 gross 3094 }
337 gross 3402 } else if (n_block==3) {
338     #pragma omp for schedule(static)
339     for (i = 0; i < n; ++i) {
340     if (coloring[i]==0) Paso_BlockOps_MViP_3(&diag[i*9], &x[3*i]);
341     }
342     } else {
343     #pragma omp for schedule(static)
344     for (i = 0; i < n; ++i) {
345     if (coloring[i]==0) Paso_BlockOps_solve_N(n_block, &x[n_block*i], &diag[block_len*i], &pivot[n_block*i], &failed);
346     }
347 gross 3097 }
348 gross 3094
349 gross 3402
350     for (color=1;color<num_colors;++color) {
351 gross 3105
352 gross 3402 if (n_block==1) {
353     #pragma omp for schedule(static)
354     for (i = 0; i < n; ++i) {
355     if (coloring[i]==color) {
356     /* x_i=x_i-a_ik*x_k */
357     rtmp=x[i];
358     for (iptr_ik=A_p->pattern->ptr[i];iptr_ik<A_p->pattern->ptr[i+1]; ++iptr_ik) {
359     k=A_p->pattern->index[iptr_ik];
360     if (coloring[k]<color) rtmp-=A_p->val[iptr_ik]*x[k];
361     }
362     x[i]=diag[i]*rtmp;
363 gross 3094 }
364     }
365 gross 3402 } else if (n_block==2) {
366     #pragma omp for schedule(static)
367     for (i = 0; i < n; ++i) {
368     if (coloring[i]==color) {
369     for (iptr_ik=A_p->pattern->ptr[i];iptr_ik<A_p->pattern->ptr[i+1]; ++iptr_ik) {
370     k=A_p->pattern->index[iptr_ik];
371     if (coloring[k]<color) Paso_BlockOps_SMV_2(&x[2*i], &A_p->val[4*iptr_ik], &x[2*k]);
372     }
373     Paso_BlockOps_MViP_2(&diag[4*i], &x[2*i]);
374 gross 3094 }
375     }
376 gross 3402 } else if (n_block==3) {
377     #pragma omp for schedule(static)
378     for (i = 0; i < n; ++i) {
379     if (coloring[i]==color) {
380     for (iptr_ik=A_p->pattern->ptr[i];iptr_ik<A_p->pattern->ptr[i+1]; ++iptr_ik) {
381     k=A_p->pattern->index[iptr_ik];
382     if (coloring[k]<color) Paso_BlockOps_SMV_3(&x[3*i], &A_p->val[9*iptr_ik], &x[3*k]);
383     }
384     Paso_BlockOps_MViP_3(&diag[9*i], &x[3*i]);
385 gross 3094 }
386     }
387 gross 3402 } else {
388     #pragma omp for schedule(static)
389     for (i = 0; i < n; ++i) {
390     if (coloring[i] == color) {
391     for (iptr_ik=A_p->pattern->ptr[i];iptr_ik<A_p->pattern->ptr[i+1]; ++iptr_ik) {
392     k=A_p->pattern->index[iptr_ik];
393     if (coloring[k]<color) Paso_BlockOps_SMV_N(n_block, &x[n_block*i], &A_p->val[block_len*iptr_ik], &x[n_block*k]);
394     }
395     Paso_BlockOps_solve_N(n_block, &x[n_block*i], &diag[block_len*i], &pivot[n_block*i], &failed);
396 gross 3097 }
397     }
398     }
399 gross 3402 } /* end of coloring loop */
400    
401     /* backward substitution */
402     for (color=(num_colors)-2 ;color>-1;--color) { /* Note: color=(num_colors)-1 is not required */
403     if (n_block==1) {
404     #pragma omp for schedule(static)
405     for (i = 0; i < n; ++i) {
406     if (coloring[i]==color) {
407     mm=ptr_main[i];
408     rtmp=A_p->val[mm]*x[i];
409     for (iptr_ik=A_p->pattern->ptr[i];iptr_ik<A_p->pattern->ptr[i+1]; ++iptr_ik) {
410     k=A_p->pattern->index[iptr_ik];
411     if (coloring[k]>color) rtmp-=A_p->val[iptr_ik]*x[k];
412     }
413     x[i]= rtmp*diag[i];
414 gross 3094 }
415     }
416 gross 3402 } else if (n_block==2) {
417     #pragma omp for schedule(static)
418     for (i = 0; i < n; ++i) {
419     if (coloring[i]==color) {
420     mm=ptr_main[i];
421     Paso_BlockOps_MViP_2(&A_p->val[4*mm], &x[2*i]);
422     for (iptr_ik=A_p->pattern->ptr[i];iptr_ik<A_p->pattern->ptr[i+1]; ++iptr_ik) {
423     k=A_p->pattern->index[iptr_ik];
424     if (coloring[k]>color) Paso_BlockOps_SMV_2(&x[2*i], &A_p->val[4*iptr_ik], &x[2*k]);
425     }
426     Paso_BlockOps_MViP_2(&diag[4*i], &x[2*i]);
427 gross 3094 }
428     }
429 gross 3402 } else if (n_block==3) {
430     #pragma omp for schedule(static)
431     for (i = 0; i < n; ++i) {
432     if (coloring[i]==color) {
433     mm=ptr_main[i];
434     Paso_BlockOps_MViP_3(&A_p->val[9*mm], &x[3*i]);
435     for (iptr_ik=A_p->pattern->ptr[i];iptr_ik<A_p->pattern->ptr[i+1]; ++iptr_ik) {
436     k=A_p->pattern->index[iptr_ik];
437     if (coloring[k]>color) Paso_BlockOps_SMV_3(&x[3*i], &A_p->val[9*iptr_ik], &x[3*k]);
438     }
439     Paso_BlockOps_MViP_3(&diag[9*i], &x[3*i]);
440 gross 3094 }
441     }
442 gross 3402 } else {
443     #pragma omp for schedule(static)
444     for (i = 0; i < n; ++i) {
445     if (coloring[i]==color) {
446     mm=ptr_main[i];
447     Paso_BlockOps_MV_N(n_block, &y[0], &A_p->val[block_len*mm], &x[n_block*i]);
448     for (iptr_ik=A_p->pattern->ptr[i];iptr_ik<A_p->pattern->ptr[i+1]; ++iptr_ik) {
449     k=A_p->pattern->index[iptr_ik];
450     if (coloring[k]>color) Paso_BlockOps_SMV_N(n_block, &y[0], &A_p->val[block_len*iptr_ik], &x[n_block*k]);
451     }
452     Paso_BlockOps_Cpy_N(n_block ,&x[n_block*i], &y[0]);
453     Paso_BlockOps_solve_N(n_block, &x[n_block*i], &diag[i*block_len], &pivot[i*n_block], &failed);
454 gross 3097 }
455     }
456     }
457     }
458 gross 3402 TMPMEMFREE(y);
459 gross 3094 }
460 gross 3402 if (failed > 0) {
461     Esys_setError(ZERO_DIVISION_ERROR, "Paso_Preconditioner_LocalSmoother_Sweep_colored: non-regular main diagonal block.");
462     }
463 gross 3094 return;
464 jfenwick 3259 }

  ViewVC Help
Powered by ViewVC 1.1.26