Fast Auxiliary Space Preconditioning  1.8.4 Feb/15/2016
itsolver_csr.c
Go to the documentation of this file.
1 
6 #include <time.h>
7 
8 #ifdef _OPENMP
9 #include <omp.h>
10 #endif
11 
12 #include "fasp.h"
13 #include "fasp_functs.h"
14 #include "itsolver_util.inl"
15 
16 /*---------------------------------*/
17 /*-- Public Functions --*/
18 /*---------------------------------*/
19 
40  dvector *b,
41  dvector *x,
42  precond *pc,
43  itsolver_param *itparam)
44 {
45  const SHORT prtlvl = itparam->print_level;
46  const SHORT itsolver_type = itparam->itsolver_type;
47  const SHORT stop_type = itparam->stop_type;
48  const SHORT restart = itparam->restart;
49  const INT MaxIt = itparam->maxit;
50  const REAL tol = itparam->tol;
51 
52  /* Local Variables */
53  REAL solver_start, solver_end, solver_duration;
54  INT iter;
55 
56  fasp_gettime(&solver_start);
57 
58 #if DEBUG_MODE > 0
59  printf("### DEBUG: %s ...... [Start]\n", __FUNCTION__);
60  printf("### DEBUG: rhs/sol size: %d %d\n", b->row, x->row);
61 #endif
62 
63  /* Safe-guard checks on parameters */
64  ITS_CHECK ( MaxIt, tol );
65 
66  /* Choose a desirable Krylov iterative solver */
67  switch ( itsolver_type ) {
68  case SOLVER_CG:
69  if ( prtlvl > PRINT_NONE ) printf("\nCalling PCG solver (CSR) ...\n");
70  iter = fasp_solver_dcsr_pcg(A, b, x, pc, tol, MaxIt, stop_type, prtlvl);
71  break;
72 
73  case SOLVER_BiCGstab:
74  if ( prtlvl > PRINT_NONE ) printf("\nCalling PBiCGstab solver (CSR) ...\n");
75  iter = fasp_solver_dcsr_pbcgs(A, b, x, pc, tol, MaxIt, stop_type, prtlvl);
76  break;
77 
78  case SOLVER_MinRes:
79  if ( prtlvl > PRINT_NONE ) printf("\nCalling PMinRes solver (CSR) ...\n");
80  iter = fasp_solver_dcsr_pminres(A, b, x, pc, tol, MaxIt, stop_type, prtlvl);
81  break;
82 
83  case SOLVER_GMRES:
84  if ( prtlvl > PRINT_NONE ) printf("\nCalling PGMRes solver (CSR) ...\n");
85  iter = fasp_solver_dcsr_pgmres(A, b, x, pc, tol, MaxIt, restart, stop_type, prtlvl);
86  break;
87 
88  case SOLVER_VGMRES:
89  if ( prtlvl > PRINT_NONE ) printf("\nCalling PVGMRes solver (CSR) ...\n");
90  iter = fasp_solver_dcsr_pvgmres(A, b, x, pc, tol, MaxIt, restart, stop_type, prtlvl);
91  break;
92 
93  case SOLVER_VFGMRES:
94  if ( prtlvl > PRINT_NONE ) printf("\nCalling PVFGMRes solver (CSR) ...\n");
95  iter = fasp_solver_dcsr_pvfgmres(A, b, x, pc, tol, MaxIt, restart, stop_type, prtlvl);
96  break;
97 
98  case SOLVER_GCG:
99  if ( prtlvl > PRINT_NONE ) printf("\nCalling PGCG solver (CSR) ...\n");
100  iter = fasp_solver_dcsr_pgcg(A, b, x, pc, tol, MaxIt, stop_type, prtlvl);
101  break;
102 
103  case SOLVER_GCR:
104  if ( prtlvl > PRINT_NONE ) printf("\nCalling PGCR solver (CSR) ...\n");
105  iter = fasp_solver_dcsr_pgcr(A, b, x, pc, tol, MaxIt, restart, stop_type, prtlvl);
106  break;
107 
108  default:
109  printf("### ERROR: Unknown itertive solver type %d!\n", itsolver_type);
110  return ERROR_SOLVER_TYPE;
111 
112  }
113 
114  if ( (prtlvl >= PRINT_SOME) && (iter >= 0) ) {
115  fasp_gettime(&solver_end);
116  solver_duration = solver_end - solver_start;
117  print_cputime("Iterative method", solver_duration);
118  }
119 
120 #if DEBUG_MODE > 0
121  printf("### DEBUG: %s ...... [Finish]\n", __FUNCTION__);
122 #endif
123 
124  return iter;
125 }
126 
144  dvector *b,
145  dvector *x,
146  itsolver_param *itparam)
147 {
148  const SHORT prtlvl = itparam->print_level;
149 
150  /* Local Variables */
151  INT status = FASP_SUCCESS;
152  REAL solver_start, solver_end, solver_duration;
153 
154 #if DEBUG_MODE > 0
155  printf("### DEBUG: %s ...... [Start]\n", __FUNCTION__);
156  printf("### DEBUG: matrix size: %d %d %d\n", A->row, A->col, A->nnz);
157  printf("### DEBUG: rhs/sol size: %d %d\n", b->row, x->row);
158 #endif
159 
160  fasp_gettime(&solver_start);
161 
162  status = fasp_solver_dcsr_itsolver(A,b,x,NULL,itparam);
163 
164  if ( prtlvl >= PRINT_MIN ) {
165  fasp_gettime(&solver_end);
166  solver_duration = solver_end - solver_start;
167  print_cputime("Krylov method totally", solver_duration);
168  }
169 
170 #if DEBUG_MODE > 0
171  printf("### DEBUG: %s ...... [Finish]\n", __FUNCTION__);
172 #endif
173 
174  return status;
175 }
176 
194  dvector *b,
195  dvector *x,
196  itsolver_param *itparam)
197 {
198  const SHORT prtlvl = itparam->print_level;
199 
200  /* Local Variables */
201  INT status = FASP_SUCCESS;
202  REAL solver_start, solver_end, solver_duration;
203 
204 #if DEBUG_MODE > 0
205  printf("### DEBUG: %s ...... [Start]\n", __FUNCTION__);
206  printf("### DEBUG: matrix size: %d %d %d\n", A->row, A->col, A->nnz);
207  printf("### DEBUG: rhs/sol size: %d %d\n", b->row, x->row);
208 #endif
209 
210  fasp_gettime(&solver_start);
211 
212  // setup preconditioner
213  dvector diag; fasp_dcsr_getdiag(0,A,&diag);
214 
215  precond pc;
216  pc.data = &diag;
217  pc.fct = fasp_precond_diag;
218 
219  // call iterative solver
220  status = fasp_solver_dcsr_itsolver(A,b,x,&pc,itparam);
221 
222  if ( prtlvl >= PRINT_MIN ) {
223  fasp_gettime(&solver_end);
224  solver_duration = solver_end - solver_start;
225  print_cputime("Diag_Krylov method totally", solver_duration);
226  }
227 
228  fasp_dvec_free(&diag);
229 
230 #if DEBUG_MODE > 0
231  printf("### DEBUG: %s ...... [Finish]\n", __FUNCTION__);
232 #endif
233 
234  return status;
235 }
236 
258  dvector *b,
259  dvector *x,
260  itsolver_param *itparam,
261  Schwarz_param *schparam)
262 {
263  Schwarz_param swzparam;
264  swzparam.Schwarz_mmsize = schparam->Schwarz_mmsize;
265  swzparam.Schwarz_maxlvl = schparam->Schwarz_maxlvl;
266  swzparam.Schwarz_type = schparam->Schwarz_type;
267  swzparam.Schwarz_blksolver = schparam->Schwarz_blksolver;
268 
269  const SHORT prtlvl = itparam->print_level;
270 
271  REAL setup_start, setup_end, setup_duration;
272  REAL solver_start, solver_end, solver_duration;
273  INT status = FASP_SUCCESS;
274 
275 #if DEBUG_MODE > 0
276  printf("### DEBUG: %s ...... [Start]\n", __FUNCTION__);
277  printf("### DEBUG: matrix size: %d %d %d\n", A->row, A->col, A->nnz);
278  printf("### DEBUG: rhs/sol size: %d %d\n", b->row, x->row);
279 #endif
280 
281  fasp_gettime(&setup_start);
282 
283  // setup preconditioner
285 
286  // symmetrize the matrix (get rid of this later)
287  Schwarz_data.A = fasp_dcsr_sympat(A);
288 
289  // construct Schwarz precondtioner
290  fasp_dcsr_shift (&Schwarz_data.A, 1);
291  fasp_Schwarz_setup(&Schwarz_data, &swzparam);
292 
293  fasp_gettime(&setup_end);
294  setup_duration = setup_end - setup_start;
295  printf("Schwarz_Krylov method setup costs %f seconds.\n", setup_duration);
296 
297  precond prec;
298  prec.data = &Schwarz_data;
299  prec.fct = fasp_precond_Schwarz;
300 
301  fasp_gettime(&solver_start);
302 
303  // solver part
304  status = fasp_solver_dcsr_itsolver(A,b,x,&prec,itparam);
305 
306  if ( prtlvl > PRINT_NONE ) {
307  fasp_gettime(&solver_end);
308  solver_duration = solver_end - solver_start;
309  printf("Schwarz_Krylov method totally costs %f seconds.\n", solver_duration);
310  }
311 
312 #if DEBUG_MODE > 0
313  printf("### DEBUG: %s ...... [Finish]\n", __FUNCTION__);
314 #endif
315 
316  fasp_Schwarz_data_free(&Schwarz_data);
317 
318  return status;
319 }
320 
339  dvector *b,
340  dvector *x,
341  itsolver_param *itparam,
342  AMG_param *amgparam)
343 {
344  const SHORT prtlvl = itparam->print_level;
345  const SHORT max_levels = amgparam->max_levels;
346  const INT nnz = A->nnz, m = A->row, n = A->col;
347 
348  /* Local Variables */
349  INT status = FASP_SUCCESS;
350  REAL solver_start, solver_end, solver_duration;
351 
352 #if DEBUG_MODE > 0
353  printf("### DEBUG: %s ...... [Start]\n", __FUNCTION__);
354  printf("### DEBUG: matrix size: %d %d %d\n", A->row, A->col, A->nnz);
355  printf("### DEBUG: rhs/sol size: %d %d\n", b->row, x->row);
356 #endif
357 
358  fasp_gettime(&solver_start);
359 
360  // initialize A, b, x for mgl[0]
361  AMG_data *mgl=fasp_amg_data_create(max_levels);
362  mgl[0].A=fasp_dcsr_create(m,n,nnz); fasp_dcsr_cp(A,&mgl[0].A);
363  mgl[0].b=fasp_dvec_create(n); mgl[0].x=fasp_dvec_create(n);
364 
365  // setup preconditioner
366  switch (amgparam->AMG_type) {
367 
368  case SA_AMG: // Smoothed Aggregation AMG
369  status = fasp_amg_setup_sa(mgl, amgparam); break;
370 
371  case UA_AMG: // Unsmoothed Aggregation AMG
372  status = fasp_amg_setup_ua(mgl, amgparam); break;
373 
374  default: // Classical AMG
375  status = fasp_amg_setup_rs(mgl, amgparam); break;
376 
377  }
378 
379 #if CHMEM_MODE
380  fasp_mem_usage();
381 #endif
382 
383  if (status < 0) goto FINISHED;
384 
385  // setup preconditioner
386  precond_data pcdata;
387  fasp_param_amg_to_prec(&pcdata,amgparam);
388  pcdata.max_levels = mgl[0].num_levels;
389  pcdata.mgl_data = mgl;
390 
391  precond pc; pc.data = &pcdata;
392 
393  if (itparam->precond_type == PREC_FMG) {
394  pc.fct = fasp_precond_famg; // Full AMG
395  }
396  else {
397  switch (amgparam->cycle_type) {
398  case AMLI_CYCLE: // AMLI cycle
399  pc.fct = fasp_precond_amli; break;
400  case NL_AMLI_CYCLE: // Nonlinear AMLI AMG
401  pc.fct = fasp_precond_nl_amli; break;
402  default: // V,W-Cycle AMG
403  pc.fct = fasp_precond_amg; break;
404  }
405  }
406 
407  // call iterative solver
408  status = fasp_solver_dcsr_itsolver(A, b, x, &pc, itparam);
409 
410  if ( prtlvl >= PRINT_MIN ) {
411  fasp_gettime(&solver_end);
412  solver_duration = solver_end - solver_start;
413  print_cputime("AMG_Krylov method totally", solver_duration);
414  }
415 
416 FINISHED:
417  fasp_amg_data_free(mgl, amgparam);
418 
419 #if DEBUG_MODE > 0
420  printf("### DEBUG: %s ...... [Finish]\n", __FUNCTION__);
421 #endif
422 
423  return status;
424 }
425 
444  dvector *b,
445  dvector *x,
446  itsolver_param *itparam,
447  ILU_param *iluparam)
448 {
449  const SHORT prtlvl = itparam->print_level;
450 
451  /* Local Variables */
452  INT status = FASP_SUCCESS;
453  REAL solver_start, solver_end, solver_duration;
454 
455 #if DEBUG_MODE > 0
456  printf("### DEBUG: %s ...... [Start]\n", __FUNCTION__);
457  printf("### DEBUG: matrix size: %d %d %d\n", A->row, A->col, A->nnz);
458  printf("### DEBUG: rhs/sol size: %d %d\n", b->row, x->row);
459 #endif
460 
461  fasp_gettime(&solver_start);
462 
463  // ILU setup for whole matrix
464  ILU_data LU;
465  if ( (status = fasp_ilu_dcsr_setup(A,&LU,iluparam)) < 0 ) goto FINISHED;
466 
467  // check iludata
468  if ( (status = fasp_mem_iludata_check(&LU)) < 0 ) goto FINISHED;
469 
470  // set preconditioner
471  precond pc;
472  pc.data = &LU;
473  pc.fct = fasp_precond_ilu;
474 
475  // call iterative solver
476  status = fasp_solver_dcsr_itsolver(A,b,x,&pc,itparam);
477 
478  if ( prtlvl >= PRINT_MIN ) {
479  fasp_gettime(&solver_end);
480  solver_duration = solver_end - solver_start;
481 
482  switch (iluparam->ILU_type) {
483  case ILUt:
484  print_cputime("ILUt_Krylov method totally", solver_duration);
485  break;
486  case ILUtp:
487  print_cputime("ILUtp_Krylov method totally", solver_duration);
488  break;
489  default: // ILUk
490  print_cputime("ILUk_Krylov method totally", solver_duration);
491  break;
492  }
493  }
494 
495 FINISHED:
496  fasp_ilu_data_free(&LU);
497 
498 #if DEBUG_MODE > 0
499  printf("### DEBUG: %s ...... [Finish]\n", __FUNCTION__);
500 #endif
501 
502  return status;
503 }
504 
528  dvector *b,
529  dvector *x,
530  itsolver_param *itparam,
531  ILU_param *iluparam,
532  dCSRmat *M)
533 {
534  const SHORT prtlvl = itparam->print_level;
535 
536  /* Local Variables */
537  REAL solver_start, solver_end, solver_duration;
538  INT status = FASP_SUCCESS;
539 
540 #if DEBUG_MODE > 0
541  printf("### DEBUG: %s ...... [Start]\n", __FUNCTION__);
542  printf("### DEBUG: matrix size: %d %d %d\n", A->row, A->col, A->nnz);
543  printf("### DEBUG: rhs/sol size: %d %d\n", b->row, x->row);
544 #endif
545 
546  fasp_gettime(&solver_start);
547 
548  // ILU setup for M
549  ILU_data LU;
550  if ( (status = fasp_ilu_dcsr_setup(M,&LU,iluparam))<0 ) goto FINISHED;
551 
552  // check iludata
553  if ( (status = fasp_mem_iludata_check(&LU))<0 ) goto FINISHED;
554 
555  // set precondtioner
556  precond pc;
557  pc.data = &LU;
558  pc.fct = fasp_precond_ilu;
559 
560  // call iterative solver
561  status = fasp_solver_dcsr_itsolver(A,b,x,&pc,itparam);
562 
563  if ( prtlvl >= PRINT_MIN ) {
564  fasp_gettime(&solver_end);
565  solver_duration = solver_end - solver_start;
566 
567  switch (iluparam->ILU_type) {
568  case ILUt:
569  print_cputime("ILUt_Krylov method totally", solver_duration);
570  break;
571  case ILUtp:
572  print_cputime("ILUtp_Krylov method totally", solver_duration);
573  break;
574  default: // ILUk
575  print_cputime("ILUk_Krylov method totally", solver_duration);
576  break;
577  }
578  }
579 
580 FINISHED:
581  fasp_ilu_data_free(&LU);
582 
583 #if DEBUG_MODE > 0
584  printf("### DEBUG: %s ...... [Finish]\n", __FUNCTION__);
585 #endif
586 
587  return status;
588 }
589 
612  dvector *b,
613  dvector *x,
614  itsolver_param *itparam,
615  AMG_param *amgparam,
616  dCSRmat *A_nk,
617  dCSRmat *P_nk,
618  dCSRmat *R_nk)
619 {
620  const SHORT prtlvl = itparam->print_level;
621  const SHORT max_levels = amgparam->max_levels;
622  const INT nnz = A->nnz, m = A->row, n = A->col;
623 
624  /* Local Variables */
625  INT status = FASP_SUCCESS;
626  REAL solver_start, solver_end, solver_duration;
627 
628 #if DEBUG_MODE > 0
629  printf("### DEBUG: %s ...... [Start]\n", __FUNCTION__);
630  printf("### DEBUG: matrix size: %d %d %d\n", A->row, A->col, A->nnz);
631  printf("### DEBUG: rhs/sol size: %d %d\n", b->row, x->row);
632 #endif
633 
634  fasp_gettime(&solver_start);
635 
636  // initialize A, b, x for mgl[0]
637  AMG_data *mgl=fasp_amg_data_create(max_levels);
638  mgl[0].A=fasp_dcsr_create(m,n,nnz); fasp_dcsr_cp(A,&mgl[0].A);
639  mgl[0].b=fasp_dvec_create(n); mgl[0].x=fasp_dvec_create(n);
640 
641  // setup preconditioner
642  switch (amgparam->AMG_type) {
643 
644  case SA_AMG: // Smoothed Aggregation AMG
645  status = fasp_amg_setup_sa(mgl, amgparam); break;
646 
647  case UA_AMG: // Unsmoothed Aggregation AMG
648  status = fasp_amg_setup_ua(mgl, amgparam); break;
649 
650  default: // Classical AMG
651  status = fasp_amg_setup_rs(mgl, amgparam); break;
652 
653  }
654 
655 #if CHMEM_MODE
656  fasp_mem_usage();
657 #endif
658 
659  if (status < 0) goto FINISHED;
660 
661  // setup preconditioner
662  precond_data pcdata;
663  fasp_param_amg_to_prec(&pcdata,amgparam);
664  pcdata.max_levels = mgl[0].num_levels;
665  pcdata.mgl_data = mgl;
666 
667  // near kernel space
668 #if WITH_UMFPACK // use UMFPACK directly
669  dCSRmat A_tran;
670  A_tran = fasp_dcsr_create(A_nk->row, A_nk->col, A_nk->nnz);
671  fasp_dcsr_transz(A_nk, NULL, &A_tran);
672  // It is equivalent to do transpose and then sort
673  // fasp_dcsr_trans(A_nk, &A_tran);
674  // fasp_dcsr_sort(&A_tran);
675  pcdata.A_nk = &A_tran;
676 #else
677  pcdata.A_nk = A_nk;
678 #endif
679 
680  pcdata.P_nk = P_nk;
681  pcdata.R_nk = R_nk;
682 
683  precond pc; pc.data = &pcdata;
685 
686  // call iterative solver
687  status = fasp_solver_dcsr_itsolver(A, b, x, &pc, itparam);
688 
689  if ( prtlvl >= PRINT_MIN ) {
690  fasp_gettime(&solver_end);
691  solver_duration = solver_end - solver_start;
692  print_cputime("AMG_NK_Krylov method totally", solver_duration);
693  }
694 
695 FINISHED:
696  fasp_amg_data_free(mgl, amgparam);
697 
698 #if WITH_UMFPACK // use UMFPACK directly
699  fasp_dcsr_free(&A_tran);
700 #endif
701 
702 #if DEBUG_MODE > 0
703  printf("### DEBUG: %s ...... [Finish]\n", __FUNCTION__);
704 #endif
705 
706  return status;
707 }
708 
709 /*---------------------------------*/
710 /*-- End of File --*/
711 /*---------------------------------*/
INT fasp_solver_dcsr_krylov_ilu(dCSRmat *A, dvector *b, dvector *x, itsolver_param *itparam, ILU_param *iluparam)
Solve Ax=b by ILUs preconditioned Krylov methods.
Definition: itsolver_csr.c:443
dvector fasp_dvec_create(const INT m)
Create dvector data space of REAL type.
Definition: vec.c:56
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
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
INT fasp_solver_dcsr_pgcg(dCSRmat *A, dvector *b, dvector *u, precond *pc, const REAL tol, const INT MaxIt, const SHORT stop_type, const SHORT prtlvl)
Preconditioned generilzed conjugate gradient (GCG) method for solving Au=b.
Definition: pgcg.c:44
SHORT num_levels
number of levels in use <= max_levels
Definition: fasp.h:730
INT fasp_solver_dcsr_pbcgs(dCSRmat *A, dvector *b, dvector *u, precond *pc, const REAL tol, const INT MaxIt, const SHORT stop_type, const SHORT prtlvl)
Preconditioned BiCGstab method for solving Au=b.
Definition: pbcgs.c:88
SHORT ILU_type
ILU type for decomposition.
Definition: fasp.h:380
INT fasp_Schwarz_setup(Schwarz_data *Schwarz, Schwarz_param *param)
Setup phase for the Schwarz methods.
INT fasp_solver_dcsr_pcg(dCSRmat *A, dvector *b, dvector *u, precond *pc, const REAL tol, const INT MaxIt, const SHORT stop_type, const SHORT prtlvl)
Preconditioned conjugate gradient method for solving Au=b.
Definition: pcg.c:84
void fasp_mem_usage()
Show total allocated memory currently.
Definition: memory.c:175
dvector b
pointer to the right-hand side at level level_num
Definition: fasp.h:744
#define SOLVER_GCG
Definition: fasp_const.h:109
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 Schwarz_maxlvl
maximal level for constructing the blocks
Definition: fasp.h:443
INT fasp_solver_dcsr_pminres(dCSRmat *A, dvector *b, dvector *u, precond *pc, const REAL tol, const INT MaxIt, const SHORT stop_type, const SHORT prtlvl)
A preconditioned minimal residual (Minres) method for solving Au=b.
Definition: pminres.c:92
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_Schwarz_data_free(Schwarz_data *Schwarz)
Free Schwarz_data data memeory space.
Definition: init.c:147
INT fasp_solver_dcsr_krylov(dCSRmat *A, dvector *b, dvector *x, itsolver_param *itparam)
Solve Ax=b by standard Krylov methods for CSR matrices.
Definition: itsolver_csr.c:143
void fasp_precond_nl_amli(REAL *r, REAL *z, void *data)
Nonlinear AMLI AMG preconditioner.
Definition: precond_csr.c:502
#define PRINT_SOME
Definition: fasp_const.h:81
#define INT
Definition: fasp.h:64
#define FASP_SUCCESS
Definition of return status and error messages.
Definition: fasp_const.h:27
Data for ILU setup.
Definition: fasp.h:400
Preconditioner data and action.
Definition: fasp.h:999
#define ERROR_SOLVER_TYPE
Definition: fasp_const.h:47
INT restart
Definition: fasp.h:1112
SHORT precond_type
Definition: fasp.h:1108
void fasp_precond_ilu(REAL *r, REAL *z, void *data)
ILU preconditioner.
Definition: precond_csr.c:185
SHORT stop_type
Definition: fasp.h:1109
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
void fasp_dcsr_transz(dCSRmat *A, INT *p, dCSRmat *AT)
Generalized transpose of A: (n x m) matrix given in dCSRmat format.
Definition: sparse_csr.c:1366
void fasp_dcsr_shift(dCSRmat *A, INT offset)
Re-index a REAL matrix in CSR format to make the index starting from 0 or 1.
Definition: sparse_csr.c:1085
AMG_data * fasp_amg_data_create(SHORT max_levels)
Create and initialize AMG_data for classical and SA AMG.
Definition: init.c:56
#define SOLVER_BiCGstab
Definition: fasp_const.h:104
#define SOLVER_VFGMRES
Definition: fasp_const.h:108
INT nnz
number of nonzero entries
Definition: fasp.h:157
SHORT Schwarz_type
type for Schwarz method
Definition: fasp.h:440
INT col
column of matrix A, n
Definition: fasp.h:154
Main header file for FASP.
SHORT max_levels
max number of levels of AMG
Definition: fasp.h:598
#define SOLVER_CG
Definition: fasp_const.h:103
dCSRmat * R_nk
Restriction for near kernel.
Definition: fasp.h:872
void fasp_gettime(REAL *time)
Get system time.
Definition: timing.c:28
#define SOLVER_GMRES
Definition: fasp_const.h:106
INT Schwarz_mmsize
maximal size of blocks
Definition: fasp.h:446
INT fasp_solver_dcsr_pgcr(dCSRmat *A, dvector *b, dvector *x, precond *pc, const REAL tol, const INT MaxIt, const SHORT restart, const SHORT stop_type, const SHORT prtlvl)
A preconditioned GCR method for solving Au=b.
Definition: pgcr.c:37
#define PRINT_NONE
Print level for all subroutines – not including DEBUG output.
Definition: fasp_const.h:79
INT row
row number of matrix A, m
Definition: fasp.h:151
Parameters passed to iterative solvers.
Definition: fasp.h:1105
INT Schwarz_blksolver
type of Schwarz block solver
Definition: fasp.h:449
void fasp_precond_diag(REAL *r, REAL *z, void *data)
Diagonal preconditioner z=inv(D)*r.
Definition: precond_csr.c:159
INT fasp_solver_dcsr_krylov_ilu_M(dCSRmat *A, dvector *b, dvector *x, itsolver_param *itparam, ILU_param *iluparam, dCSRmat *M)
Solve Ax=b by ILUs preconditioned Krylov methods: ILU of M as preconditioner.
Definition: itsolver_csr.c:527
Sparse matrix of REAL type in CSR format.
Definition: fasp.h:148
SHORT fasp_mem_iludata_check(ILU_data *iludata)
Check wether a ILU_data has enough work space.
Definition: memory.c:222
dCSRmat * A_nk
Matrix data for near kernel.
Definition: fasp.h:866
INT fasp_solver_dcsr_krylov_Schwarz(dCSRmat *A, dvector *b, dvector *x, itsolver_param *itparam, Schwarz_param *schparam)
Solve Ax=b by overlapping Schwarz Krylov methods.
Definition: itsolver_csr.c:257
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
INT row
number of rows
Definition: fasp.h:345
REAL tol
Definition: fasp.h:1111
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
INT fasp_solver_dcsr_pgmres(dCSRmat *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 for solving Au=b.
Definition: pgmres.c:53
SHORT print_level
Definition: fasp.h:1113
Parameters for Schwarz method.
Definition: fasp.h:434
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 print_cputime(const char *message, const REAL cputime)
Print CPU walltime.
Definition: message.c:165
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
INT fasp_solver_dcsr_krylov_amg(dCSRmat *A, dvector *b, dvector *x, itsolver_param *itparam, AMG_param *amgparam)
Solve Ax=b by AMG preconditioned Krylov methods.
Definition: itsolver_csr.c:338
#define PREC_FMG
Definition: fasp_const.h:141
INT fasp_solver_dcsr_krylov_diag(dCSRmat *A, dvector *b, dvector *x, itsolver_param *itparam)
Solve Ax=b by diagonal preconditioned Krylov methods.
Definition: itsolver_csr.c:193
Data for AMG solvers.
Definition: fasp.h:722
Data passed to the preconditioners.
Definition: fasp.h:795
INT fasp_solver_dcsr_krylov_amg_nk(dCSRmat *A, dvector *b, dvector *x, itsolver_param *itparam, AMG_param *amgparam, dCSRmat *A_nk, dCSRmat *P_nk, dCSRmat *R_nk)
Solve Ax=b by AMG preconditioned Krylov methods with an extra near kernel solve.
Definition: itsolver_csr.c:611
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
INT fasp_solver_dcsr_pvgmres(dCSRmat *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:51
#define NL_AMLI_CYCLE
Definition: fasp_const.h:178
dCSRmat * P_nk
Prolongation for near kernel.
Definition: fasp.h:869
#define SOLVER_VGMRES
Definition: fasp_const.h:107
INT fasp_solver_dcsr_itsolver(dCSRmat *A, dvector *b, dvector *x, precond *pc, itsolver_param *itparam)
Solve Ax=b by preconditioned Krylov methods for CSR matrices.
Definition: itsolver_csr.c:39
#define SOLVER_MinRes
Definition: fasp_const.h:105
#define SOLVER_GCR
Definition: fasp_const.h:110
#define SA_AMG
Definition: fasp_const.h:163
void fasp_dcsr_free(dCSRmat *A)
Free CSR sparse matrix data memory space.
Definition: sparse_csr.c:166
INT fasp_solver_dcsr_pvfgmres(dCSRmat *A, dvector *b, dvector *x, precond *pc, const REAL tol, const INT MaxIt, const SHORT restart, const SHORT stop_type, const SHORT prtlvl)
Solve "Ax=b" using PFGMRES(right preconditioned) iterative method in which the restart parameter can ...
Definition: pvfgmres.c:54
#define UA_AMG
Definition: fasp_const.h:164
#define ILUtp
Definition: fasp_const.h:150
#define PRINT_MIN
Definition: fasp_const.h:80
SHORT itsolver_type
Definition: fasp.h:1107
dCSRmat fasp_dcsr_sympat(dCSRmat *A)
Get symmetric part of a dCSRmat matrix.
Definition: sparse_csr.c:1232
#define ILUt
Definition: fasp_const.h:149