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