Fast Auxiliary Space Preconditioning  1.8.4 Feb/15/2016
precond_csr.c
Go to the documentation of this file.
1 
6 #include "fasp.h"
7 #include "fasp_functs.h"
9 
10 #include "mg_util.inl"
11 
12 /*---------------------------------*/
13 /*-- Public Functions --*/
14 /*---------------------------------*/
15 
32 precond *fasp_precond_setup (const SHORT precond_type,
33  AMG_param *amgparam,
34  ILU_param *iluparam,
35  dCSRmat *A)
36 {
37  precond *pc = NULL;
38 
39  AMG_data *mgl = NULL;
40  precond_data *pcdata = NULL;
41  ILU_data *ILU = NULL;
42  dvector *diag = NULL;
43 
44  INT max_levels, nnz, m, n;
45 
46  switch (precond_type) {
47 
48  case PREC_AMG: // AMG preconditioner
49 
50  pc = (precond *)fasp_mem_calloc(1, sizeof(precond));
51  max_levels = amgparam->max_levels;
52  nnz=A->nnz, m=A->row, n=A->col;
53 
54  // initialize A, b, x for mgl[0]
55  mgl=fasp_amg_data_create(max_levels);
56  mgl[0].A=fasp_dcsr_create(m,n,nnz); fasp_dcsr_cp(A,&mgl[0].A);
57  mgl[0].b=fasp_dvec_create(n); mgl[0].x=fasp_dvec_create(n);
58 
59  // setup preconditioner
60  switch (amgparam->AMG_type) {
61  case SA_AMG: // Smoothed Aggregation AMG
62  fasp_amg_setup_sa(mgl, amgparam); break;
63  case UA_AMG: // Unsmoothed Aggregation AMG
64  fasp_amg_setup_ua(mgl, amgparam); break;
65  default: // Classical AMG
66  fasp_amg_setup_rs(mgl, amgparam); break;
67  }
68 
69  pcdata = (precond_data *)fasp_mem_calloc(1, sizeof(precond_data));
70  fasp_param_amg_to_prec(pcdata, amgparam);
71  pcdata->max_levels = mgl[0].num_levels;
72  pcdata->mgl_data = mgl;
73 
74  pc->data = pcdata;
75 
76  switch (amgparam->cycle_type) {
77  case AMLI_CYCLE: // AMLI cycle
78  pc->fct = fasp_precond_amli; break;
79  case NL_AMLI_CYCLE: // Nonlinear AMLI AMG
80  pc->fct = fasp_precond_nl_amli; break;
81  default: // V,W-Cycle AMG
82  pc->fct = fasp_precond_amg; break;
83  }
84 
85  break;
86 
87  case PREC_FMG: // FMG preconditioner
88 
89  pc = (precond *)fasp_mem_calloc(1, sizeof(precond));
90  max_levels = amgparam->max_levels;
91  nnz=A->nnz, m=A->row, n=A->col;
92 
93  // initialize A, b, x for mgl[0]
94  mgl=fasp_amg_data_create(max_levels);
95  mgl[0].A=fasp_dcsr_create(m,n,nnz); fasp_dcsr_cp(A,&mgl[0].A);
96  mgl[0].b=fasp_dvec_create(n); mgl[0].x=fasp_dvec_create(n);
97 
98  // setup preconditioner
99  switch (amgparam->AMG_type) {
100  case SA_AMG: // Smoothed Aggregation AMG
101  fasp_amg_setup_sa(mgl, amgparam); break;
102  case UA_AMG: // Unsmoothed Aggregation AMG
103  fasp_amg_setup_ua(mgl, amgparam); break;
104  default: // Classical AMG
105  fasp_amg_setup_rs(mgl, amgparam); break;
106  }
107 
108  pcdata = (precond_data *)fasp_mem_calloc(1, sizeof(precond_data));
109  fasp_param_amg_to_prec(pcdata, amgparam);
110  pcdata->max_levels = mgl[0].num_levels;
111  pcdata->mgl_data = mgl;
112 
113  pc->data = pcdata; pc->fct = fasp_precond_famg;
114 
115  break;
116 
117  case PREC_ILU: // ILU preconditioner
118 
119  pc = (precond *)fasp_mem_calloc(1, sizeof(precond));
120  ILU = (ILU_data *)fasp_mem_calloc(1, sizeof(ILU_data));
121  fasp_ilu_dcsr_setup(A, ILU, iluparam);
122  pc->data = ILU;
123  pc->fct = fasp_precond_ilu;
124 
125  break;
126 
127  case PREC_DIAG: // Diagonal preconditioner
128 
129  pc = (precond *)fasp_mem_calloc(1, sizeof(precond));
130  diag = (dvector *)fasp_mem_calloc(1, sizeof(dvector));
131  fasp_dcsr_getdiag(0, A, diag);
132 
133  pc->data = diag;
134  pc->fct = fasp_precond_diag;
135 
136  break;
137 
138  default: // No preconditioner
139 
140  break;
141 
142  }
143 
144  return pc;
145 }
146 
160  REAL *z,
161  void *data)
162 {
163  dvector *diag=(dvector *)data;
164  REAL *diagptr=diag->val;
165  INT i, m=diag->row;
166 
167  memcpy(z,r,m*sizeof(REAL));
168  for (i=0;i<m;++i) {
169  if (ABS(diag->val[i])>SMALLREAL) z[i]/=diagptr[i];
170  }
171 }
172 
186  REAL *z,
187  void *data)
188 {
189  ILU_data *iludata=(ILU_data *)data;
190  const INT m=iludata->row, mm1=m-1, memneed=2*m;
191  REAL *zz, *zr;
192 
193  if (iludata->nwork<memneed) goto MEMERR; // check this outside this subroutine!!
194 
195  zz = iludata->work;
196  zr = iludata->work+m;
197  fasp_array_cp(m, r, zr);
198 
199  {
200  INT i, j, jj, begin_row, end_row, mm2=m-2;
201  INT *ijlu=iludata->ijlu;
202  REAL *lu=iludata->luval;
203 
204  // forward sweep: solve unit lower matrix equation L*zz=zr
205  zz[0]=zr[0];
206 
207  for (i=1;i<=mm1;++i) {
208  begin_row=ijlu[i]; end_row=ijlu[i+1]-1;
209  for (j=begin_row;j<=end_row;++j) {
210  jj=ijlu[j];
211  if (jj<i) zr[i]-=lu[j]*zz[jj];
212  else break;
213  }
214  zz[i]=zr[i];
215  }
216 
217  // backward sweep: solve upper matrix equation U*z=zz
218  z[mm1]=zz[mm1]*lu[mm1];
219  for (i=mm2;i>=0;i--) {
220  begin_row=ijlu[i]; end_row=ijlu[i+1]-1;
221  for (j=end_row;j>=begin_row;j--) {
222  jj=ijlu[j];
223  if (jj>i) zz[i]-=lu[j]*z[jj];
224  else break;
225  }
226  z[i]=zz[i]*lu[i];
227  }
228  }
229 
230  return;
231 
232 MEMERR:
233  printf("### ERROR: Need %d memory, only %d available!\n", memneed, iludata->nwork);
234  exit(ERROR_ALLOC_MEM);
235 }
236 
250  REAL *z,
251  void *data)
252 {
253  ILU_data *iludata=(ILU_data *)data;
254  const INT m=iludata->row, mm1=m-1, memneed=2*m;
255  REAL *zz, *zr;
256 
257  if (iludata->nwork<memneed) goto MEMERR;
258 
259  zz = iludata->work;
260  zr = iludata->work+m;
261  fasp_array_cp(m, r, zr);
262 
263  {
264  INT i, j, jj, begin_row, end_row;
265  INT *ijlu=iludata->ijlu;
266  REAL *lu=iludata->luval;
267 
268  // forward sweep: solve unit lower matrix equation L*z=r
269  zz[0]=zr[0];
270  for (i=1;i<=mm1;++i) {
271  begin_row=ijlu[i]; end_row=ijlu[i+1]-1;
272  for (j=begin_row;j<=end_row;++j) {
273  jj=ijlu[j];
274  if (jj<i) zr[i]-=lu[j]*zz[jj];
275  else break;
276  }
277  zz[i]=zr[i];
278  }
279  }
280 
281  fasp_array_cp(m, zz, z);
282 
283  return;
284 
285 MEMERR:
286  printf("### ERROR: Need %d memory, only %d available!", memneed, iludata->nwork);
287  exit(ERROR_ALLOC_MEM);
288 }
289 
303  REAL *z,
304  void *data)
305 {
306  ILU_data *iludata=(ILU_data *)data;
307  const INT m=iludata->row, mm1=m-1, memneed=2*m;
308  REAL *zz;
309 
310  if (iludata->nwork<memneed) goto MEMERR;
311 
312  zz = iludata->work;
313  fasp_array_cp(m, r, zz);
314 
315  {
316  INT i, j, jj, begin_row, end_row, mm2=m-2;
317  INT *ijlu=iludata->ijlu;
318  REAL *lu=iludata->luval;
319 
320  // backward sweep: solve upper matrix equation U*z=zz
321  z[mm1]=zz[mm1]*lu[mm1];
322  for (i=mm2;i>=0;i--) {
323  begin_row=ijlu[i]; end_row=ijlu[i+1]-1;
324  for (j=end_row;j>=begin_row;j--) {
325  jj=ijlu[j];
326  if (jj>i) zz[i]-=lu[j]*z[jj];
327  else break;
328  }
329  z[i]=zz[i]*lu[i];
330  }
331 
332  }
333 
334  return;
335 
336 MEMERR:
337  printf("### ERROR: Need %d memory, only %d available!", memneed, iludata->nwork);
338  exit(ERROR_ALLOC_MEM);
339 }
340 
356  REAL *z,
357  void *data)
358 {
359  Schwarz_data * swzdata = (Schwarz_data *)data;
360  Schwarz_param * swzparam = swzdata->swzparam;
361  const INT swztype = swzdata->Schwarz_type;
362  const INT n = swzdata->A.row;
363 
364  dvector x, b;
365 
366  fasp_dvec_alloc(n, &x);
367  fasp_dvec_alloc(n, &b);
368  fasp_array_cp(n, r, b.val);
369 
370  fasp_dvec_set(n, &x, 0);
371 
372  switch (swztype) {
373  case 2:
374  fasp_dcsr_Schwarz_backward_smoother(swzdata, swzparam, &x, &b);
375  break;
376  case 3:
377  fasp_dcsr_Schwarz_forward_smoother(swzdata, swzparam, &x, &b);
378  fasp_dcsr_Schwarz_backward_smoother(swzdata, swzparam, &x, &b);
379  break;
380  default:
381  fasp_dcsr_Schwarz_forward_smoother(swzdata, swzparam, &x, &b);
382  break;
383  }
384 
385  fasp_array_cp(n, x.val, z);
386 }
387 
401  REAL *z,
402  void *data)
403 {
404  precond_data *pcdata=(precond_data *)data;
405  const INT m=pcdata->mgl_data[0].A.row;
406  const INT maxit=pcdata->maxit;
407  INT i;
408 
409  AMG_param amgparam; fasp_param_amg_init(&amgparam);
410  fasp_param_prec_to_amg(&amgparam,pcdata);
411 
412  AMG_data *mgl = pcdata->mgl_data;
413  mgl->b.row=m; fasp_array_cp(m,r,mgl->b.val); // residual is an input
414  mgl->x.row=m; fasp_dvec_set(m,&mgl->x,0.0);
415 
416  for (i=0;i<maxit;++i) fasp_solver_mgcycle(mgl,&amgparam);
417 
418  // We can also use a recursive version of MG:
419  // for (i=0;i<maxit;++i) fasp_solver_mgrecur(mgl,&amgparam,0);
420 
421  fasp_array_cp(m,mgl->x.val,z);
422 }
423 
437  REAL *z,
438  void *data)
439 {
440  precond_data *pcdata=(precond_data *)data;
441  const INT m=pcdata->mgl_data[0].A.row;
442  const INT maxit=pcdata->maxit;
443  INT i;
444 
445  AMG_param amgparam; fasp_param_amg_init(&amgparam);
446  fasp_param_prec_to_amg(&amgparam,pcdata);
447 
448  AMG_data *mgl = pcdata->mgl_data;
449  mgl->b.row=m; fasp_array_cp(m,r,mgl->b.val); // residual is an input
450  mgl->x.row=m; fasp_dvec_set(m,&mgl->x,0.0);
451 
452  for (i=0;i<maxit;++i) fasp_solver_fmgcycle(mgl,&amgparam);
453 
454  fasp_array_cp(m,mgl->x.val,z);
455 }
456 
470  REAL *z,
471  void *data)
472 {
473  precond_data *pcdata=(precond_data *)data;
474  const INT m=pcdata->mgl_data[0].A.row;
475  const INT maxit=pcdata->maxit;
476  INT i;
477 
478  AMG_param amgparam; fasp_param_amg_init(&amgparam);
479  fasp_param_prec_to_amg(&amgparam,pcdata);
480 
481  AMG_data *mgl = pcdata->mgl_data;
482  mgl->b.row=m; fasp_array_cp(m,r,mgl->b.val); // residual is an input
483  mgl->x.row=m; fasp_dvec_set(m,&mgl->x,0.0);
484 
485  for (i=0;i<maxit;++i) fasp_solver_amli(mgl,&amgparam,0);
486 
487  fasp_array_cp(m,mgl->x.val,z);
488 }
489 
503  REAL *z,
504  void *data)
505 {
506  precond_data *pcdata=(precond_data *)data;
507  const INT m=pcdata->mgl_data[0].A.row;
508  const INT maxit=pcdata->maxit;
509  const SHORT num_levels = pcdata->max_levels;
510  INT i;
511 
512  AMG_param amgparam; fasp_param_amg_init(&amgparam);
513  fasp_param_prec_to_amg(&amgparam,pcdata);
514 
515  AMG_data *mgl = pcdata->mgl_data;
516  mgl->b.row=m; fasp_array_cp(m,r,mgl->b.val); // residual is an input
517  mgl->x.row=m; fasp_dvec_set(m,&mgl->x,0.0);
518 
519  for (i=0;i<maxit;++i) fasp_solver_nl_amli(mgl, &amgparam, 0, num_levels);
520  fasp_array_cp(m,mgl->x.val,z);
521 }
522 
536  REAL *z,
537  void *data)
538 {
539  precond_data *pcdata=(precond_data *)data;
540  const INT m=pcdata->mgl_data[0].A.row;
541  const INT maxit=pcdata->maxit;
542  INT i;
543 
544  dCSRmat *A_nk = pcdata->A_nk;
545  dCSRmat *P_nk = pcdata->P_nk;
546  dCSRmat *R_nk = pcdata->R_nk;
547 
548  fasp_array_set(m, z, 0.0);
549 
550  // local variables
551  dvector r_nk, z_nk;
552  fasp_dvec_alloc(A_nk->row, &r_nk);
553  fasp_dvec_alloc(A_nk->row, &z_nk);
554 
555  //----------------------
556  // extra kernel solve
557  //----------------------
558  // r_nk = R_nk*r
559  fasp_blas_dcsr_mxv(R_nk, r, r_nk.val);
560 
561  // z_nk = A_nk^{-1}*r_nk
562 #if WITH_UMFPACK // use UMFPACK directly
563  fasp_solver_umfpack(A_nk, &r_nk, &z_nk, 0);
564 #else
565  fasp_coarse_itsolver(A_nk, &r_nk, &z_nk, 1e-12, 0);
566 #endif
567 
568  // z = z + P_nk*z_nk;
569  fasp_blas_dcsr_aAxpy(1.0, P_nk, z_nk.val, z);
570 
571  //----------------------
572  // AMG solve
573  //----------------------
574  AMG_param amgparam; fasp_param_amg_init(&amgparam);
575  fasp_param_prec_to_amg(&amgparam,pcdata);
576 
577  AMG_data *mgl = pcdata->mgl_data;
578  mgl->b.row=m; fasp_array_cp(m,r,mgl->b.val); // residual is an input
579  mgl->x.row=m; //fasp_dvec_set(m,&mgl->x,0.0);
580  fasp_array_cp(m, z, mgl->x.val);
581 
582  for ( i = 0; i < maxit; ++i ) fasp_solver_mgcycle(mgl,&amgparam);
583 
584  fasp_array_cp(m,mgl->x.val,z);
585 
586  //----------------------
587  // extra kernel solve
588  //----------------------
589  // r = r - A*z
590  fasp_blas_dcsr_aAxpy(-1.0, &(pcdata->mgl_data[0].A), z, mgl->b.val);
591 
592  // r_nk = R_nk*r
593  fasp_blas_dcsr_mxv(R_nk, mgl->b.val, r_nk.val);
594 
595  // z_nk = A_nk^{-1}*r_nk
596 #if WITH_UMFPACK // use UMFPACK directly
597  fasp_solver_umfpack(A_nk, &r_nk, &z_nk, 0);
598 #else
599  fasp_coarse_itsolver(A_nk, &r_nk, &z_nk, 1e-12, 0);
600 #endif
601 
602  // z = z + P_nk*z_nk;
603  fasp_blas_dcsr_aAxpy(1.0, P_nk, z_nk.val, z);
604 }
605 
619 void fasp_precond_free (const SHORT precond_type,
620  precond *pc)
621 {
622  switch (precond_type) {
623 
624  case PREC_AMG: // AMG preconditioner
625 
626  fasp_amg_data_free(((precond_data*)(pc->data))->mgl_data, NULL);
628  fasp_mem_free(pc);
629 
630  break;
631 
632  case PREC_FMG: // FMG preconditioner
633 
634  fasp_amg_data_free(((precond_data*)(pc->data))->mgl_data, NULL);
636  fasp_mem_free(pc);
637 
638  break;
639 
640  case PREC_ILU: // ILU preconditioner
641 
643  fasp_mem_free(pc);
644 
645  break;
646 
647  case PREC_DIAG: // Diagonal preconditioner
648 
649  fasp_dvec_free((dvector*)(pc->data));
650  fasp_mem_free(pc);
651 
652  break;
653 
654  default: // No preconditioner
655 
656  break;
657 
658  }
659 }
660 
661 /*---------------------------------*/
662 /*-- End of File --*/
663 /*---------------------------------*/
dvector fasp_dvec_create(const INT m)
Create dvector data space of REAL type.
Definition: vec.c:56
INT * ijlu
integer array of row pointers and column indexes, the size is nzlu
Definition: fasp.h:412
dCSRmat fasp_dcsr_create(const INT m, const INT n, const INT nnz)
Create CSR sparse matrix data memory space.
Definition: sparse_csr.c:34
void fasp_solver_fmgcycle(AMG_data *mgl, AMG_param *param)
#include "forts_ns.h"
Definition: fmgcycle.c:36
void fasp_dcsr_Schwarz_forward_smoother(Schwarz_data *Schwarz, Schwarz_param *param, dvector *x, dvector *b)
Schwarz smoother: forward sweep.
SHORT AMG_type
type of AMG method
Definition: fasp.h:586
Parameters for AMG solver.
Definition: fasp.h:583
void fasp_precond_amli(REAL *r, REAL *z, void *data)
AMLI AMG preconditioner.
Definition: precond_csr.c:469
#define REAL
Definition: fasp.h:67
SHORT fasp_amg_setup_ua(AMG_data *mgl, AMG_param *param)
Set up phase of unsmoothed aggregation AMG.
Definition: amg_setup_ua.c:38
void fasp_precond_amg(REAL *r, REAL *z, void *data)
AMG preconditioner.
Definition: precond_csr.c:400
SHORT num_levels
number of levels in use <= max_levels
Definition: fasp.h:730
INT maxit
max number of iterations of AMG preconditioner
Definition: fasp.h:804
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.
void fasp_dvec_free(dvector *u)
Free vector data space of REAL type.
Definition: vec.c:139
SHORT fasp_amg_setup_rs(AMG_data *mgl, AMG_param *param)
Setup phase of Ruge and Stuben's classic AMG.
Definition: amg_setup_rs.c:47
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_mem_calloc(LONGLONG size, INT type)
1M = 1024*1024
Definition: memory.c:60
void fasp_param_prec_to_amg(AMG_param *amgparam, precond_data *pcdata)
Set AMG_param with precond_data.
Definition: parameters.c:701
void(* fct)(REAL *, REAL *, void *)
action for preconditioner, void function pointer
Definition: fasp.h:1005
Vector with n entries of REAL type.
Definition: fasp.h:342
void fasp_precond_nl_amli(REAL *r, REAL *z, void *data)
Nonlinear AMLI AMG preconditioner.
Definition: precond_csr.c:502
#define INT
Definition: fasp.h:64
Data for ILU setup.
Definition: fasp.h:400
Preconditioner data and action.
Definition: fasp.h:999
INT Schwarz_type
Schwarz method type.
Definition: fasp.h:539
INT nwork
work space size
Definition: fasp.h:421
void fasp_mem_free(void *mem)
Free up previous allocated memory body.
Definition: memory.c:150
void fasp_precond_ilu(REAL *r, REAL *z, void *data)
ILU preconditioner.
Definition: precond_csr.c:185
void fasp_array_set(const INT n, REAL *x, const REAL val)
Set initial value for an array to be x=val.
Definition: array.c:48
void * data
data for preconditioner, void pointer
Definition: fasp.h:1002
void fasp_param_amg_to_prec(precond_data *pcdata, AMG_param *amgparam)
Set precond_data with AMG_param.
Definition: parameters.c:666
dvector x
pointer to the iterative solution at level level_num
Definition: fasp.h:747
SHORT fasp_amg_setup_sa(AMG_data *mgl, AMG_param *param)
Set up phase of smoothed aggregation AMG.
Definition: amg_setup_sa.c:48
void fasp_dcsr_cp(dCSRmat *A, dCSRmat *B)
copy a dCSRmat to a new one B=A
Definition: sparse_csr.c:723
void fasp_precond_amg_nk(REAL *r, REAL *z, void *data)
AMG with extra near kernel solve as preconditioner.
Definition: precond_csr.c:535
void fasp_precond_Schwarz(REAL *r, REAL *z, void *data)
get z from r by Schwarz
Definition: precond_csr.c:355
Schwarz_param * swzparam
param for Schwarz
Definition: fasp.h:573
AMG_data * fasp_amg_data_create(SHORT max_levels)
Create and initialize AMG_data for classical and SA AMG.
Definition: init.c:56
INT nnz
number of nonzero entries
Definition: fasp.h:157
REAL * luval
nonzero entries of LU
Definition: fasp.h:415
INT col
column of matrix A, n
Definition: fasp.h:154
precond * fasp_precond_setup(const SHORT precond_type, AMG_param *amgparam, ILU_param *iluparam, dCSRmat *A)
#include "forts_ns.h"
Definition: precond_csr.c:32
Main header file for FASP.
SHORT max_levels
max number of levels of AMG
Definition: fasp.h:598
dCSRmat * R_nk
Restriction for near kernel.
Definition: fasp.h:872
#define ERROR_ALLOC_MEM
Definition: fasp_const.h:37
INT row
row number of matrix LU, m
Definition: fasp.h:403
INT row
row number of matrix A, m
Definition: fasp.h:151
void fasp_precond_diag(REAL *r, REAL *z, void *data)
Diagonal preconditioner z=inv(D)*r.
Definition: precond_csr.c:159
Sparse matrix of REAL type in CSR format.
Definition: fasp.h:148
dCSRmat * A_nk
Matrix data for near kernel.
Definition: fasp.h:866
void fasp_dcsr_getdiag(INT n, dCSRmat *A, dvector *diag)
Get first n diagonal entries of a CSR matrix A.
Definition: sparse_csr.c:410
void fasp_precond_ilu_backward(REAL *r, REAL *z, void *data)
ILU preconditioner: only backward sweep.
Definition: precond_csr.c:302
#define PREC_AMG
Definition: fasp_const.h:140
INT row
number of rows
Definition: fasp.h:345
void fasp_param_amg_init(AMG_param *amgparam)
Initialize AMG parameters.
Definition: parameters.c:390
SHORT cycle_type
type of AMG cycle
Definition: fasp.h:604
SHORT max_levels
max number of AMG levels
Definition: fasp.h:807
dCSRmat A
pointer to the matrix
Definition: fasp.h:510
Parameters for Schwarz method.
Definition: fasp.h:434
#define ABS(a)
Definition: fasp.h:74
SHORT fasp_ilu_dcsr_setup(dCSRmat *A, ILU_data *iludata, ILU_param *iluparam)
Get ILU decomposition of a CSR matrix A.
Definition: ilu_setup_csr.c:50
AMG_data * mgl_data
AMG preconditioner data.
Definition: fasp.h:855
Parameters for ILU.
Definition: fasp.h:374
#define SHORT
FASP integer and floating point numbers.
Definition: fasp.h:63
#define AMLI_CYCLE
Definition: fasp_const.h:177
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
void fasp_precond_famg(REAL *r, REAL *z, void *data)
Full AMG preconditioner.
Definition: precond_csr.c:436
Data for Schwarz methods.
Definition: fasp.h:505
#define PREC_FMG
Definition: fasp_const.h:141
void fasp_precond_ilu_forward(REAL *r, REAL *z, void *data)
ILU preconditioner: only forward sweep.
Definition: precond_csr.c:249
void fasp_precond_free(const SHORT precond_type, precond *pc)
free preconditioner
Definition: precond_csr.c:619
Data for AMG solvers.
Definition: fasp.h:722
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
Data passed to the preconditioners.
Definition: fasp.h:795
void fasp_dvec_alloc(const INT m, dvector *u)
Create dvector data space of REAL type.
Definition: vec.c:99
void fasp_amg_data_free(AMG_data *mgl, AMG_param *param)
Free AMG_data data memeory space.
Definition: init.c:185
void fasp_ilu_data_free(ILU_data *ILUdata)
Create ILU_data sturcture.
Definition: init.c:300
void fasp_array_cp(const INT n, REAL *x, REAL *y)
Copy an array to the other y=x.
Definition: array.c:165
REAL * work
work space
Definition: fasp.h:424
void fasp_solver_nl_amli(AMG_data *mgl, AMG_param *param, INT level, INT num_levels)
Solve Ax=b with recursive nonlinear AMLI-cycle.
Definition: amlirecur.c:269
#define NL_AMLI_CYCLE
Definition: fasp_const.h:178
dCSRmat * P_nk
Prolongation for near kernel.
Definition: fasp.h:869
#define PREC_DIAG
Definition: fasp_const.h:139
void fasp_solver_amli(AMG_data *mgl, AMG_param *param, INT level)
Solve Ax=b with recursive AMLI-cycle.
Definition: amlirecur.c:45
#define SA_AMG
Definition: fasp_const.h:163
#define PREC_ILU
Definition: fasp_const.h:142
#define SMALLREAL
Definition: fasp_const.h:238
#define UA_AMG
Definition: fasp_const.h:164
void fasp_blas_dcsr_mxv(dCSRmat *A, REAL *x, REAL *y)
Matrix-vector multiplication y = A*x.
Definition: blas_csr.c:225