Fast Auxiliary Space Preconditioning 2.7.7 Aug/28/2022
PreMGCycle.c
Go to the documentation of this file.
1
17#include <math.h>
18#include <time.h>
19
20#include "fasp.h"
21#include "fasp_functs.h"
22
23/*---------------------------------*/
24/*-- Declare Private Functions --*/
25/*---------------------------------*/
26
27#include "PreMGUtil.inl"
28#include "PreMGSmoother.inl"
29
30/*---------------------------------*/
31/*-- Public Functions --*/
32/*---------------------------------*/
33
49 AMG_param *param)
50{
51 const SHORT prtlvl = param->print_level;
52 const SHORT amg_type = param->AMG_type;
53 const SHORT smoother = param->smoother;
54 const SHORT smooth_order = param->smooth_order;
55 const SHORT cycle_type = param->cycle_type;
56 const SHORT coarse_solver = param->coarse_solver;
57 const SHORT nl = mgl[0].num_levels;
58 const REAL relax = param->relaxation;
59 const REAL tol = param->tol * 1e-4;
60 const SHORT ndeg = param->polynomial_degree;
61
62 // Schwarz parameters
63 SWZ_param swzparam;
64 if ( param->SWZ_levels > 0 ) {
65 swzparam.SWZ_blksolver = param->SWZ_blksolver;
66 }
67
68 // local variables
69 REAL alpha = 1.0;
70 INT num_lvl[MAX_AMG_LVL] = {0}, l = 0;
71
72 // more general cycling types on each level --zcs 05/07/2020
73 INT ncycles[MAX_AMG_LVL] = {1};
74 SHORT i;
75 for ( i = 0; i < MAX_AMG_LVL; ++i ) ncycles[i] = 1; // initially V-cycle
76 switch(cycle_type) {
77 case 12:
78 for ( i = MAX_AMG_LVL-2; i > 0; i -= 2 ) ncycles[i] = 2;
79 break;
80 case 21:
81 for ( i = MAX_AMG_LVL-1; i > 0; i -= 2 ) ncycles[i] = 2;
82 break;
83 default:
84 for ( i = 0; i < MAX_AMG_LVL; i += 1 ) ncycles[i] = cycle_type;
85 }
86
87#if DEBUG_MODE > 0
88 printf("### DEBUG: [-Begin-] %s ...\n", __FUNCTION__);
89 printf("### DEBUG: n=%d, nnz=%d\n", mgl[0].A.row, mgl[0].A.nnz);
90#endif
91
92#if DEBUG_MODE > 1
93 printf("### DEBUG: AMG_level = %d, ILU_level = %d\n", nl, mgl->ILU_levels);
94#endif
95
96ForwardSweep:
97 while ( l < nl-1 ) {
98
99 num_lvl[l]++;
100
101 // pre-smoothing with ILU method
102 if ( l < mgl->ILU_levels ) {
103 fasp_smoother_dcsr_ilu(&mgl[l].A, &mgl[l].b, &mgl[l].x, &mgl[l].LU);
104 }
105
106 // or pre-smoothing with Schwarz method
107 else if ( l < mgl->SWZ_levels ) {
108 switch (mgl[l].Schwarz.SWZ_type) {
110 fasp_dcsr_swz_forward(&mgl[l].Schwarz, &swzparam, &mgl[l].x, &mgl[l].b);
111 fasp_dcsr_swz_backward(&mgl[l].Schwarz,&swzparam, &mgl[l].x, &mgl[l].b);
112 break;
113 default:
114 fasp_dcsr_swz_forward(&mgl[l].Schwarz, &swzparam, &mgl[l].x, &mgl[l].b);
115 break;
116 }
117 }
118
119 // or pre-smoothing with standard smoother
120 else {
121#if MULTI_COLOR_ORDER
122 // printf("fasp_smoother_dcsr_gs_multicolor, %s, %d\n", __FUNCTION__, __LINE__);
123 fasp_smoother_dcsr_gs_multicolor (&mgl[l].x, &mgl[l].A, &mgl[l].b, param->presmooth_iter,1);
124#else
125 fasp_dcsr_presmoothing(smoother, &mgl[l].A, &mgl[l].b, &mgl[l].x,
126 param->presmooth_iter, 0, mgl[l].A.row-1, 1,
127 relax, ndeg, smooth_order, mgl[l].cfmark.val);
128#endif
129 }
130
131 // form residual r = b - A x
132 fasp_darray_cp(mgl[l].A.row, mgl[l].b.val, mgl[l].w.val);
133 fasp_blas_dcsr_aAxpy(-1.0, &mgl[l].A, mgl[l].x.val, mgl[l].w.val);
134
135 // restriction r1 = R*r0
136 switch ( amg_type ) {
137 case UA_AMG:
138 fasp_blas_dcsr_mxv_agg(&mgl[l].R, mgl[l].w.val, mgl[l+1].b.val);
139 break;
140 default:
141 fasp_blas_dcsr_mxv(&mgl[l].R, mgl[l].w.val, mgl[l+1].b.val);
142 break;
143 }
144
145 // prepare for the next level
146 ++l; fasp_dvec_set(mgl[l].A.row, &mgl[l].x, 0.0);
147
148 }
149
150 // If AMG only has one level or we have arrived at the coarsest level,
151 // call the coarse space solver:
152 switch ( coarse_solver ) {
153
154#if WITH_PARDISO
155 case SOLVER_PARDISO: {
156 /* use Intel MKL PARDISO direct solver on the coarsest level */
157 fasp_pardiso_solve(&mgl[nl-1].A, &mgl[nl-1].b, &mgl[nl-1].x, &mgl[nl-1].pdata, 0);
158 break;
159 }
160#endif
161
162#if WITH_MUMPS
163 case SOLVER_MUMPS: {
164 // use MUMPS direct solver on the coarsest level
165 mgl[nl-1].mumps.job = 2;
166 fasp_solver_mumps_steps(&mgl[nl-1].A, &mgl[nl-1].b, &mgl[nl-1].x, &mgl[nl-1].mumps);
167 break;
168 }
169#endif
170
171#if WITH_UMFPACK
172 case SOLVER_UMFPACK: {
173 // use UMFPACK direct solver on the coarsest level
174 fasp_umfpack_solve(&mgl[nl-1].A, &mgl[nl-1].b, &mgl[nl-1].x, mgl[nl-1].Numeric, 0);
175 break;
176 }
177#endif
178
179#if WITH_SuperLU
180 case SOLVER_SUPERLU: {
181 // use SuperLU direct solver on the coarsest level
182 fasp_solver_superlu(&mgl[nl-1].A, &mgl[nl-1].b, &mgl[nl-1].x, 0);
183 break;
184 }
185#endif
186
187 default:
188 // use iterative solver on the coarsest level
189 fasp_coarse_itsolver(&mgl[nl-1].A, &mgl[nl-1].b, &mgl[nl-1].x, tol, prtlvl);
190
191 }
192
193 // BackwardSweep:
194 while ( l > 0 ) {
195
196 --l;
197
198 // find the optimal scaling factor alpha
199 if ( param->coarse_scaling == ON ) {
200 alpha = fasp_blas_darray_dotprod(mgl[l+1].A.row, mgl[l+1].x.val, mgl[l+1].b.val)
201 / fasp_blas_dcsr_vmv(&mgl[l+1].A, mgl[l+1].x.val, mgl[l+1].x.val);
202 alpha = MIN(alpha, 1.0); // Add this for safety! --Chensong on 10/04/2014
203 }
204
205 // prolongation u = u + alpha*P*e1
206 switch ( amg_type ) {
207 case UA_AMG:
208 fasp_blas_dcsr_aAxpy_agg(alpha, &mgl[l].P, mgl[l+1].x.val, mgl[l].x.val);
209 break;
210 default:
211 fasp_blas_dcsr_aAxpy(alpha, &mgl[l].P, mgl[l+1].x.val, mgl[l].x.val);
212 break;
213 }
214
215 // post-smoothing with ILU method
216 if ( l < mgl->ILU_levels ) {
217 fasp_smoother_dcsr_ilu(&mgl[l].A, &mgl[l].b, &mgl[l].x, &mgl[l].LU);
218 }
219
220 // post-smoothing with Schwarz method
221 else if ( l < mgl->SWZ_levels ) {
222 switch (mgl[l].Schwarz.SWZ_type) {
224 fasp_dcsr_swz_backward(&mgl[l].Schwarz,&swzparam, &mgl[l].x, &mgl[l].b);
225 fasp_dcsr_swz_forward(&mgl[l].Schwarz, &swzparam, &mgl[l].x, &mgl[l].b);
226 break;
227 default:
228 fasp_dcsr_swz_backward(&mgl[l].Schwarz,&swzparam, &mgl[l].x, &mgl[l].b);
229 break;
230 }
231 }
232
233 // post-smoothing with standard methods
234 else {
235#if MULTI_COLOR_ORDER
236 fasp_smoother_dcsr_gs_multicolor (&mgl[l].x, &mgl[l].A, &mgl[l].b, param->postsmooth_iter,-1);
237#else
238 fasp_dcsr_postsmoothing(smoother, &mgl[l].A, &mgl[l].b, &mgl[l].x,
239 param->postsmooth_iter, 0, mgl[l].A.row-1, -1,
240 relax, ndeg, smooth_order, mgl[l].cfmark.val);
241#endif
242 }
243
244 // General cycling on each level --zcs
245 if ( num_lvl[l] < ncycles[l] ) break;
246 else num_lvl[l] = 0;
247 }
248
249 if ( l > 0 ) goto ForwardSweep;
250
251#if DEBUG_MODE > 0
252 printf("### DEBUG: [--End--] %s ...\n", __FUNCTION__);
253#endif
254
255}
256
269 AMG_param *param)
270{
271 const SHORT prtlvl = param->print_level;
272 const SHORT nl = mgl[0].num_levels;
273 const SHORT smoother = param->smoother;
274 const SHORT cycle_type = param->cycle_type;
275 const SHORT coarse_solver = param->coarse_solver;
276 const REAL relax = param->relaxation;
277 INT steps = param->presmooth_iter;
278
279 // local variables
280 INT nu_l[MAX_AMG_LVL] = {0}, l = 0;
281 REAL alpha = 1.0;
282 INT i;
283
284 dvector r_nk, z_nk;
285
286 if ( mgl[0].A_nk != NULL ) {
287 fasp_dvec_alloc(mgl[0].A_nk->row, &r_nk);
288 fasp_dvec_alloc(mgl[0].A_nk->row, &z_nk);
289 }
290
291#if DEBUG_MODE > 0
292 printf("### DEBUG: [-Begin-] %s ...\n", __FUNCTION__);
293#endif
294
295#if DEBUG_MODE > 1
296 printf("### DEBUG: AMG_level = %d, ILU_level = %d\n", nl, mgl->ILU_levels);
297#endif
298
299ForwardSweep:
300 while ( l < nl-1 ) {
301 nu_l[l]++;
302 // pre smoothing
303 if ( l < mgl->ILU_levels ) {
304 fasp_smoother_dbsr_ilu(&mgl[l].A, &mgl[l].b, &mgl[l].x, &mgl[l].LU);
305 for ( i=0; i<steps; i++ )
306 fasp_smoother_dbsr_gs_ascend(&mgl[l].A, &mgl[l].b, &mgl[l].x, mgl[l].diaginv.val);
307 }
308 else {
309 if ( steps > 0 ) {
310 switch ( smoother ) {
311 case SMOOTHER_JACOBI:
312 for (i=0; i<steps; i++)
313 fasp_smoother_dbsr_jacobi1(&mgl[l].A, &mgl[l].b, &mgl[l].x,
314 mgl[l].diaginv.val);
315 break;
316 case SMOOTHER_GS:
317 for (i=0; i<steps; i++)
318 fasp_smoother_dbsr_gs_ascend(&mgl[l].A, &mgl[l].b, &mgl[l].x,
319 mgl[l].diaginv.val);
320 break;
321 case SMOOTHER_SGS:
322 for (i=0; i<steps; i++){
323 fasp_smoother_dbsr_gs_ascend(&mgl[l].A, &mgl[l].b, &mgl[l].x,
324 mgl[l].diaginv.val);
325 fasp_smoother_dbsr_gs_descend(&mgl[l].A, &mgl[l].b, &mgl[l].x,
326 mgl[l].diaginv.val);
327 }
328 break;
329 case SMOOTHER_SOR:
330 for (i=0; i<steps; i++)
331 fasp_smoother_dbsr_sor_ascend(&mgl[l].A, &mgl[l].b, &mgl[l].x,
332 mgl[l].diaginv.val, relax);
333 break;
334 case SMOOTHER_SSOR:
335 for (i=0; i<steps; i++) {
336 fasp_smoother_dbsr_sor_ascend(&mgl[l].A, &mgl[l].b, &mgl[l].x,
337 mgl[l].diaginv.val, relax);
338 }
339 fasp_smoother_dbsr_sor_descend(&mgl[l].A, &mgl[l].b, &mgl[l].x,
340 mgl[l].diaginv.val, relax);
341 break;
342 default:
343 printf("### ERROR: Unknown smoother type %d!\n", smoother);
344 fasp_chkerr(ERROR_SOLVER_TYPE, __FUNCTION__);
345 }
346 }
347 }
348
349 // extra kernel solve
350 if (mgl[l].A_nk != NULL) {
351
352 //--------------------------------------------
353 // extra kernel solve
354 //--------------------------------------------
355 // form residual r = b - A x
356 fasp_darray_cp(mgl[l].A.ROW*mgl[l].A.nb, mgl[l].b.val, mgl[l].w.val);
357 fasp_blas_dbsr_aAxpy(-1.0,&mgl[l].A, mgl[l].x.val, mgl[l].w.val);
358
359 // r_nk = R_nk*r
360 fasp_blas_dcsr_mxv(mgl[l].R_nk, mgl[l].w.val, r_nk.val);
361
362 // z_nk = A_nk^{-1}*r_nk
363#if WITH_UMFPACK // use UMFPACK directly
364 fasp_solver_umfpack(mgl[l].A_nk, &r_nk, &z_nk, 0);
365#else
366 fasp_coarse_itsolver(mgl[l].A_nk, &r_nk, &z_nk, 1e-12, 0);
367#endif
368
369 // z = z + P_nk*z_nk;
370 fasp_blas_dcsr_aAxpy(1.0, mgl[l].P_nk, z_nk.val, mgl[l].x.val);
371 }
372
373 // form residual r = b - A x
374 fasp_darray_cp(mgl[l].A.ROW*mgl[l].A.nb, mgl[l].b.val, mgl[l].w.val);
375 fasp_blas_dbsr_aAxpy(-1.0,&mgl[l].A, mgl[l].x.val, mgl[l].w.val);
376
377 // restriction r1 = R*r0
378 fasp_blas_dbsr_mxv(&mgl[l].R, mgl[l].w.val, mgl[l+1].b.val);
379
380 // prepare for the next level
381 ++l; fasp_dvec_set(mgl[l].A.ROW*mgl[l].A.nb, &mgl[l].x, 0.0);
382
383 }
384
385 // If AMG only has one level or we have arrived at the coarsest level,
386 // call the coarse space solver:
387 switch ( coarse_solver ) {
388
389#if WITH_PARDISO
390 case SOLVER_PARDISO: {
391 /* use Intel MKL PARDISO direct solver on the coarsest level */
392 fasp_pardiso_solve(&mgl[nl-1].Ac, &mgl[nl-1].b, &mgl[nl-1].x, &mgl[nl-1].pdata, 0);
393 break;
394 }
395#endif
396
397#if WITH_MUMPS
398 case SOLVER_MUMPS:
399 /* use MUMPS direct solver on the coarsest level */
400 mgl[nl-1].mumps.job = 2;
401 fasp_solver_mumps_steps(&mgl[nl-1].Ac, &mgl[nl-1].b, &mgl[nl-1].x, &mgl[nl-1].mumps);
402 break;
403#endif
404
405#if WITH_UMFPACK
406 case SOLVER_UMFPACK:
407 /* use UMFPACK direct solver on the coarsest level */
408 fasp_umfpack_solve(&mgl[nl-1].Ac, &mgl[nl-1].b, &mgl[nl-1].x, mgl[nl-1].Numeric, 0);
409 break;
410#endif
411
412#if WITH_SuperLU
413 case SOLVER_SUPERLU:
414 /* use SuperLU direct solver on the coarsest level */
415 fasp_solver_superlu(&mgl[nl-1].Ac, &mgl[nl-1].b, &mgl[nl-1].x, 0);
416 break;
417#endif
418
419 default: {
420 /* use iterative solver on the coarsest level */
421 const INT csize = mgl[nl-1].A.ROW*mgl[nl-1].A.nb;
422 const INT cmaxit = MIN(csize*csize, 200); // coarse level iteration number
423 const REAL ctol = param->tol; // coarse level tolerance
424 if ( fasp_solver_dbsr_pvgmres(&mgl[nl-1].A,&mgl[nl-1].b,&mgl[nl-1].x,
425 NULL,ctol,cmaxit,25,1,0) < 0 ) {
426 if ( prtlvl > PRINT_MIN ) {
427 printf("### WARNING: Coarse level solver did not converge!\n");
428 printf("### WARNING: Consider to increase maxit to %d!\n", 2*cmaxit);
429 }
430 }
431 }
432 }
433
434 // BackwardSweep:
435 while ( l > 0 ) {
436 --l;
437
438 // prolongation u = u + alpha*P*e1
439 if ( param->coarse_scaling == ON ) {
440 dvector PeH, Aeh;
441 PeH.row = Aeh.row = mgl[l].b.row;
442 PeH.val = mgl[l].w.val + mgl[l].b.row;
443 Aeh.val = PeH.val + mgl[l].b.row;
444
445 fasp_blas_dbsr_mxv (&mgl[l].P, mgl[l+1].x.val, PeH.val);
446 fasp_blas_dbsr_mxv (&mgl[l].A, PeH.val, Aeh.val);
447
448 alpha = (fasp_blas_darray_dotprod (mgl[l].b.row, Aeh.val, mgl[l].w.val))
449 / (fasp_blas_darray_dotprod (mgl[l].b.row, Aeh.val, Aeh.val));
450 alpha = MIN(alpha, 1.0); // Add this for safety! --Chensong on 10/04/2014
451 fasp_blas_darray_axpy (mgl[l].b.row, alpha, PeH.val, mgl[l].x.val);
452 }
453 else {
454 fasp_blas_dbsr_aAxpy(alpha, &mgl[l].P, mgl[l+1].x.val, mgl[l].x.val);
455 }
456
457 // extra kernel solve
458 if ( mgl[l].A_nk != NULL ) {
459 //--------------------------------------------
460 // extra kernel solve
461 //--------------------------------------------
462 // form residual r = b - A x
463 fasp_darray_cp(mgl[l].A.ROW*mgl[l].A.nb, mgl[l].b.val, mgl[l].w.val);
464 fasp_blas_dbsr_aAxpy(-1.0, &mgl[l].A, mgl[l].x.val, mgl[l].w.val);
465
466 // r_nk = R_nk*r
467 fasp_blas_dcsr_mxv(mgl[l].R_nk, mgl[l].w.val, r_nk.val);
468
469 // z_nk = A_nk^{-1}*r_nk
470#if WITH_UMFPACK // use UMFPACK directly
471 fasp_solver_umfpack(mgl[l].A_nk, &r_nk, &z_nk, 0);
472#else
473 fasp_coarse_itsolver(mgl[l].A_nk, &r_nk, &z_nk, 1e-12, 0);
474#endif
475
476 // z = z + P_nk*z_nk;
477 fasp_blas_dcsr_aAxpy(1.0, mgl[l].P_nk, z_nk.val, mgl[l].x.val);
478 }
479
480 // post-smoothing
481 if ( l < mgl->ILU_levels ) {
482 fasp_smoother_dbsr_ilu(&mgl[l].A, &mgl[l].b, &mgl[l].x, &mgl[l].LU);
483 for ( i=0; i<steps; i++ )
484 fasp_smoother_dbsr_gs_descend(&mgl[l].A, &mgl[l].b, &mgl[l].x,
485 mgl[l].diaginv.val);
486 }
487 else {
488 if ( steps > 0 ) {
489 switch ( smoother ) {
490 case SMOOTHER_JACOBI:
491 for (i=0; i<steps; i++)
492 fasp_smoother_dbsr_jacobi1(&mgl[l].A, &mgl[l].b, &mgl[l].x,
493 mgl[l].diaginv.val);
494 break;
495 case SMOOTHER_GS:
496 for (i=0; i<steps; i++)
497 fasp_smoother_dbsr_gs_descend(&mgl[l].A, &mgl[l].b, &mgl[l].x,
498 mgl[l].diaginv.val);
499 break;
500 case SMOOTHER_SGS:
501 for (i=0; i<steps; i++){
502 fasp_smoother_dbsr_gs_ascend(&mgl[l].A, &mgl[l].b, &mgl[l].x,
503 mgl[l].diaginv.val);
504 fasp_smoother_dbsr_gs_descend(&mgl[l].A, &mgl[l].b, &mgl[l].x,
505 mgl[l].diaginv.val);
506 }
507 break;
508 case SMOOTHER_SOR:
509 for (i=0; i<steps; i++)
510 fasp_smoother_dbsr_sor_descend(&mgl[l].A, &mgl[l].b, &mgl[l].x,
511 mgl[l].diaginv.val, relax);
512 break;
513 case SMOOTHER_SSOR:
514 for (i=0; i<steps; i++)
515 fasp_smoother_dbsr_sor_ascend(&mgl[l].A, &mgl[l].b, &mgl[l].x,
516 mgl[l].diaginv.val, relax);
517 fasp_smoother_dbsr_sor_descend(&mgl[l].A, &mgl[l].b, &mgl[l].x,
518 mgl[l].diaginv.val, relax);
519 break;
520 default:
521 printf("### ERROR: Unknown smoother type %d!\n", smoother);
522 fasp_chkerr(ERROR_SOLVER_TYPE, __FUNCTION__);
523 }
524 }
525 }
526
527 if ( nu_l[l] < cycle_type ) break;
528 else nu_l[l] = 0;
529 }
530
531 if ( l > 0 ) goto ForwardSweep;
532
533#if DEBUG_MODE > 0
534 printf("### DEBUG: [--End--] %s ...\n", __FUNCTION__);
535#endif
536
537}
538
539/*---------------------------------*/
540/*-- End of File --*/
541/*---------------------------------*/
void fasp_darray_cp(const INT n, const REAL *x, REAL *y)
Copy an array to the other y=x.
Definition: AuxArray.c:164
void fasp_chkerr(const SHORT status, const char *fctname)
Check error status and print out error messages before quit.
Definition: AuxMessage.c:213
void fasp_dvec_set(INT n, dvector *x, const REAL val)
Initialize dvector x[i]=val for i=0:n-1.
Definition: AuxVector.c:222
void fasp_dvec_alloc(const INT m, dvector *u)
Create dvector data space of REAL type.
Definition: AuxVector.c:105
REAL fasp_blas_darray_dotprod(const INT n, const REAL *x, const REAL *y)
Inner product of two arraies x and y.
Definition: BlaArray.c:741
void fasp_blas_darray_axpy(const INT n, const REAL a, const REAL *x, REAL *y)
y = a*x + y
Definition: BlaArray.c:93
void fasp_dcsr_swz_forward(SWZ_data *swzdata, SWZ_param *swzparam, dvector *x, dvector *b)
Schwarz smoother: forward sweep.
void fasp_dcsr_swz_backward(SWZ_data *swzdata, SWZ_param *swzparam, dvector *x, dvector *b)
Schwarz smoother: backward sweep.
void fasp_blas_dbsr_aAxpy(const REAL alpha, const dBSRmat *A, const REAL *x, REAL *y)
Compute y := alpha*A*x + y.
Definition: BlaSpmvBSR.c:348
void fasp_blas_dbsr_mxv(const dBSRmat *A, const REAL *x, REAL *y)
Compute y := A*x.
Definition: BlaSpmvBSR.c:910
void fasp_blas_dcsr_mxv_agg(const dCSRmat *A, const REAL *x, REAL *y)
Matrix-vector multiplication y = A*x (nonzeros of A = 1)
Definition: BlaSpmvCSR.c:437
void fasp_blas_dcsr_mxv(const dCSRmat *A, const REAL *x, REAL *y)
Matrix-vector multiplication y = A*x.
Definition: BlaSpmvCSR.c:241
void fasp_blas_dcsr_aAxpy_agg(const REAL alpha, const dCSRmat *A, const REAL *x, REAL *y)
Matrix-vector multiplication y = alpha*A*x + y (nonzeros of A = 1)
Definition: BlaSpmvCSR.c:724
REAL fasp_blas_dcsr_vmv(const dCSRmat *A, const REAL *x, const REAL *y)
vector-Matrix-vector multiplication alpha = y'*A*x
Definition: BlaSpmvCSR.c:834
void fasp_blas_dcsr_aAxpy(const REAL alpha, const dCSRmat *A, const REAL *x, REAL *y)
Matrix-vector multiplication y = alpha*A*x + y.
Definition: BlaSpmvCSR.c:493
void fasp_smoother_dbsr_gs_descend(dBSRmat *A, dvector *b, dvector *u, REAL *diaginv)
Gauss-Seidel relaxation in the descending order.
void fasp_smoother_dbsr_jacobi1(dBSRmat *A, dvector *b, dvector *u, REAL *diaginv)
Jacobi relaxation.
void fasp_smoother_dbsr_gs_ascend(dBSRmat *A, dvector *b, dvector *u, REAL *diaginv)
Gauss-Seidel relaxation in the ascending order.
void fasp_smoother_dbsr_sor_ascend(dBSRmat *A, dvector *b, dvector *u, REAL *diaginv, REAL weight)
SOR relaxation in the ascending order.
void fasp_smoother_dbsr_ilu(dBSRmat *A, dvector *b, dvector *x, void *data)
ILU method as the smoother in solving Au=b with multigrid method.
void fasp_smoother_dbsr_sor_descend(dBSRmat *A, dvector *b, dvector *u, REAL *diaginv, REAL weight)
SOR relaxation in the descending order.
void fasp_smoother_dcsr_ilu(dCSRmat *A, dvector *b, dvector *x, void *data)
ILU method as a smoother.
INT fasp_solver_dbsr_pvgmres(dBSRmat *A, dvector *b, dvector *x, precond *pc, const REAL tol, const INT MaxIt, const SHORT restart, const SHORT StopType, const SHORT PrtLvl)
Right preconditioned GMRES method in which the restart parameter can be adaptively modified during it...
Definition: KryPvgmres.c:413
void fasp_solver_mgcycle_bsr(AMG_data_bsr *mgl, AMG_param *param)
Solve Ax=b with non-recursive multigrid cycle.
Definition: PreMGCycle.c:268
void fasp_solver_mgcycle(AMG_data *mgl, AMG_param *param)
Solve Ax=b with non-recursive multigrid cycle.
Definition: PreMGCycle.c:48
int fasp_solver_mumps_steps(dCSRmat *ptrA, dvector *b, dvector *u, Mumps_data *mumps)
Solve Ax=b by MUMPS in three steps.
Definition: XtrMumps.c:188
int fasp_solver_superlu(dCSRmat *ptrA, dvector *b, dvector *u, const SHORT prtlvl)
Solve Au=b by SuperLU.
Definition: XtrSuperlu.c:47
INT fasp_solver_umfpack(dCSRmat *ptrA, dvector *b, dvector *u, const SHORT prtlvl)
Solve Au=b by UMFpack.
Definition: XtrUmfpack.c:44
Main header file for the FASP project.
#define MIN(a, b)
Definition: fasp.h:75
#define REAL
Definition: fasp.h:67
#define SHORT
FASP integer and floating point numbers.
Definition: fasp.h:63
#define INT
Definition: fasp.h:64
#define SMOOTHER_SOR
Definition: fasp_const.h:191
#define SOLVER_PARDISO
Definition: fasp_const.h:126
#define MAX_AMG_LVL
Definition: fasp_const.h:252
#define SMOOTHER_JACOBI
Definition of standard smoother types.
Definition: fasp_const.h:187
#define SOLVER_MUMPS
Definition: fasp_const.h:125
#define ERROR_SOLVER_TYPE
Definition: fasp_const.h:41
#define SOLVER_SUPERLU
Definition: fasp_const.h:123
#define SCHWARZ_SYMMETRIC
Definition: fasp_const.h:157
#define SMOOTHER_GS
Definition: fasp_const.h:188
#define SOLVER_UMFPACK
Definition: fasp_const.h:124
#define SMOOTHER_SGS
Definition: fasp_const.h:189
#define ON
Definition of switch.
Definition: fasp_const.h:67
#define SMOOTHER_SSOR
Definition: fasp_const.h:192
#define PRINT_MIN
Definition: fasp_const.h:74
#define UA_AMG
Definition: fasp_const.h:164
Data for multigrid levels in dBSRmat format.
Definition: fasp_block.h:146
dBSRmat A
pointer to the matrix at level level_num
Definition: fasp_block.h:155
Mumps_data mumps
data for MUMPS
Definition: fasp_block.h:231
dvector b
pointer to the right-hand side at level level_num
Definition: fasp_block.h:164
INT ILU_levels
number of levels use ILU smoother
Definition: fasp_block.h:203
INT num_levels
number of levels in use <= max_levels
Definition: fasp_block.h:152
dvector x
pointer to the iterative solution at level level_num
Definition: fasp_block.h:167
dvector w
temporary work space
Definition: fasp_block.h:228
Data for AMG methods.
Definition: fasp.h:790
dCSRmat A
pointer to the matrix at level level_num
Definition: fasp.h:803
Mumps_data mumps
data for MUMPS
Definition: fasp.h:852
dvector b
pointer to the right-hand side at level level_num
Definition: fasp.h:812
INT ILU_levels
number of levels use ILU smoother
Definition: fasp.h:829
dvector x
pointer to the iterative solution at level level_num
Definition: fasp.h:815
SHORT num_levels
number of levels in use <= max_levels
Definition: fasp.h:798
dvector w
temporary work space
Definition: fasp.h:849
ivector cfmark
pointer to the CF marker at level level_num
Definition: fasp.h:826
Parameters for AMG methods.
Definition: fasp.h:447
SHORT print_level
print level for AMG
Definition: fasp.h:453
SHORT polynomial_degree
degree of the polynomial smoother
Definition: fasp.h:489
SHORT coarse_scaling
switch of scaling of the coarse grid correction
Definition: fasp.h:495
SHORT AMG_type
type of AMG method
Definition: fasp.h:450
REAL tol
stopping tolerance for AMG solver
Definition: fasp.h:459
SHORT coarse_solver
coarse solver type
Definition: fasp.h:492
REAL relaxation
relaxation parameter for Jacobi and SOR smoother
Definition: fasp.h:486
SHORT smoother
smoother type
Definition: fasp.h:474
SHORT cycle_type
type of AMG cycle
Definition: fasp.h:468
INT SWZ_blksolver
type of Schwarz block solver
Definition: fasp.h:579
SHORT postsmooth_iter
number of postsmoothers
Definition: fasp.h:483
INT SWZ_levels
number of levels use Schwarz smoother
Definition: fasp.h:567
SHORT presmooth_iter
number of presmoothers
Definition: fasp.h:480
SHORT smooth_order
smoother order
Definition: fasp.h:477
INT job
work for MUMPS
Definition: fasp.h:601
Parameters for Schwarz method.
Definition: fasp.h:422
INT SWZ_blksolver
type of Schwarz block solver
Definition: fasp.h:437
INT nb
dimension of each sub-block
Definition: fasp_block.h:46
INT ROW
number of rows of sub-blocks in matrix A, M
Definition: fasp_block.h:37
INT row
row number of matrix A, m
Definition: fasp.h:146
INT nnz
number of nonzero entries
Definition: fasp.h:152
Vector with n entries of REAL type.
Definition: fasp.h:346
REAL * val
actual vector entries
Definition: fasp.h:352
INT row
number of rows
Definition: fasp.h:349
INT * val
actual vector entries
Definition: fasp.h:367