Fast Auxiliary Space Preconditioning  1.8.4 Feb/15/2016
itsolver_bcsr.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_block.h"
11 #include "fasp_functs.h"
12 #include "itsolver_util.inl"
13 
14 /*---------------------------------*/
15 /*-- Public Functions --*/
16 /*---------------------------------*/
17 
37  dvector *b,
38  dvector *x,
39  precond *pc,
40  itsolver_param *itparam)
41 {
42  const SHORT prtlvl = itparam->print_level;
43  const SHORT itsolver_type = itparam->itsolver_type;
44  const SHORT stop_type = itparam->stop_type;
45  const SHORT restart = itparam->restart;
46  const INT MaxIt = itparam->maxit;
47  const REAL tol = itparam->tol;
48 
49  REAL solver_start, solver_end, solver_duration;
50  INT iter = ERROR_SOLVER_TYPE;
51 
52  fasp_gettime(&solver_start);
53 
54 #if DEBUG_MODE > 0
55  printf("### DEBUG: %s ...... [Start]\n", __FUNCTION__);
56  printf("### DEBUG: rhs/sol size: %d %d\n", b->row, x->row);
57 #endif
58 
59  /* Safe-guard checks on parameters */
60  ITS_CHECK ( MaxIt, tol );
61 
62  switch (itsolver_type) {
63 
64  case SOLVER_BiCGstab:
65  if ( prtlvl > PRINT_NONE ) printf("\nCalling BiCGstab solver (Block CSR) ...\n");
66  iter=fasp_solver_bdcsr_pbcgs(A, b, x, pc, tol, MaxIt, stop_type, prtlvl);
67  break;
68 
69  case SOLVER_MinRes:
70  if ( prtlvl > PRINT_NONE ) printf("\nCalling MinRes solver (Block CSR) ...\n");
71  iter=fasp_solver_bdcsr_pminres(A, b, x, pc, tol, MaxIt, stop_type, prtlvl);
72  break;
73 
74  case SOLVER_GMRES:
75  if ( prtlvl > PRINT_NONE ) printf("\nCalling GMRES solver (Block CSR) ...\n");
76  iter=fasp_solver_bdcsr_pgmres(A, b, x, pc, tol, MaxIt, restart, stop_type, prtlvl);
77  break;
78 
79  case SOLVER_VGMRES:
80  if ( prtlvl > PRINT_NONE ) printf("Calling vGMRES solver (Block CSR) ...\n");
81  iter=fasp_solver_bdcsr_pvgmres(A, b, x, pc, tol, MaxIt, restart, stop_type, prtlvl);
82  break;
83 
84  case SOLVER_VFGMRES:
85  if ( prtlvl > PRINT_NONE ) printf("Calling FGMRES solver (Block CSR) ...\n");
86  iter=fasp_solver_bdcsr_pvfgmres(A, b, x, pc, tol, MaxIt, restart, stop_type, prtlvl);
87  break;
88 
89  default:
90  printf("### ERROR: Unknown itertive solver type %d!\n", itsolver_type);
91 
92  }
93 
94  if ( (prtlvl >= PRINT_MIN) && (iter >= 0) ) {
95  fasp_gettime(&solver_end);
96  solver_duration = solver_end - solver_start;
97  print_cputime("Iterative method", solver_duration);
98  }
99 
100 #if DEBUG_MODE > 0
101  printf("### DEBUG: %s ...... [Finish]\n", __FUNCTION__);
102 #endif
103 
104  return iter;
105 }
106 
124  dvector *b,
125  dvector *x,
126  itsolver_param *itparam)
127 {
128  const SHORT prtlvl = itparam->print_level;
129 
130  INT status = FASP_SUCCESS;
131  REAL solver_start, solver_end, solver_duration;
132 
133 #if DEBUG_MODE > 0
134  printf("### DEBUG: %s ...... [Start]\n", __FUNCTION__);
135 #endif
136 
137  // solver part
138  fasp_gettime(&solver_start);
139 
140  status = fasp_solver_bdcsr_itsolver(A,b,x,NULL,itparam);
141 
142  fasp_gettime(&solver_end);
143 
144  solver_duration = solver_end - solver_start;
145 
146  if ( prtlvl >= PRINT_MIN )
147  print_cputime("Krylov method totally", solver_duration);
148 
149 #if DEBUG_MODE > 0
150  printf("### DEBUG: %s ...... [Finish]\n", __FUNCTION__);
151 #endif
152 
153  return status;
154 }
155 
178  dvector *b,
179  dvector *x,
180  itsolver_param *itparam,
181  AMG_param *amgparam,
182  dCSRmat *A_diag)
183 {
184  const SHORT prtlvl = itparam->print_level;
185  const SHORT precond_type = itparam->precond_type;
186 
187  INT status = FASP_SUCCESS;
188  REAL setup_start, setup_end, setup_duration;
189  REAL solver_start, solver_end, solver_duration;
190 
191  const SHORT max_levels = amgparam->max_levels;
192  INT m, n, nnz, i;
193 
194  AMG_data **mgl = NULL;
195 
196 #if WITH_UMFPACK
197  void **LU_diag = (void **)fasp_mem_calloc(3, sizeof(void *));
198 #endif
199 
200 #if DEBUG_MODE > 0
201  printf("### DEBUG: %s ...... [Start]\n", __FUNCTION__);
202 #endif
203 
204  /* setup preconditioner */
205  fasp_gettime(&setup_start);
206 
207  /* diagonal blocks are solved exactly */
208  if ( precond_type > 20 && precond_type < 30 ) {
209 
210 #if WITH_UMFPACK
211  // Need to sort the diagonal blocks for UMFPACK format
212  dCSRmat A_tran;
213 
214  for (i=0; i<3; i++){
215 
216  A_tran = fasp_dcsr_create(A_diag[i].row, A_diag[i].col, A_diag[i].nnz);
217  fasp_dcsr_transz(&A_diag[i], NULL, &A_tran);
218  fasp_dcsr_cp(&A_tran, &A_diag[i]);
219 
220  printf("Factorization for %d-th diagnol: \n", i);
221  LU_diag[i] = fasp_umfpack_factorize(&A_diag[i], prtlvl);
222 
223  }
224 
225  fasp_dcsr_free(&A_tran);
226 #endif
227 
228  }
229 
230  /* diagonal blocks are solved by AMG */
231  else if ( precond_type > 30 && precond_type < 40 ) {
232 
233  mgl = (AMG_data **)fasp_mem_calloc(3, sizeof(AMG_data *));
234 
235  for (i=0; i<3; i++){
236 
237  mgl[i] = fasp_amg_data_create(max_levels);
238  m = A_diag[i].row; n = A_diag[i].col; nnz = A_diag[i].nnz;
239  mgl[i][0].A=fasp_dcsr_create(m,n,nnz); fasp_dcsr_cp(&A_diag[i],&mgl[i][0].A);
240  mgl[i][0].b=fasp_dvec_create(n); mgl[i][0].x=fasp_dvec_create(n);
241 
242  switch (amgparam->AMG_type) {
243  case SA_AMG: // Smoothed Aggregation AMG
244  status = fasp_amg_setup_sa(mgl[i], amgparam); break;
245  case UA_AMG: // Unsmoothed Aggregation AMG
246  status = fasp_amg_setup_ua(mgl[i], amgparam); break;
247  default: // Classical AMG
248  status = fasp_amg_setup_rs(mgl[i], amgparam); break;
249  }
250 
251  }
252 
253  }
254 
255  else {
256  fasp_chkerr(ERROR_SOLVER_PRECTYPE, __FUNCTION__);
257  }
258 
259  precond_block_data precdata;
260  precdata.Abcsr = A;
261  precdata.A_diag = A_diag;
262  precdata.r = fasp_dvec_create(b->row);
263 
264  /* diagonal blocks are solved exactly */
265  if ( precond_type > 20 && precond_type < 30 ) {
266 #if WITH_UMFPACK
267  precdata.LU_diag = LU_diag;
268 #endif
269  }
270  /* diagonal blocks are solved by AMG */
271  else if ( precond_type > 30 && precond_type < 40 ) {
272  precdata.amgparam = amgparam;
273  precdata.mgl = mgl;
274  }
275  else {
276  fasp_chkerr(ERROR_SOLVER_PRECTYPE, __FUNCTION__);
277  }
278 
279  precond prec; prec.data = &precdata;
280 
281  switch (precond_type)
282  {
283  case 21:
285  break;
286 
287  case 22:
289  break;
290 
291  case 23:
293  break;
294 
295  case 24:
297  break;
298 
299  case 31:
301  break;
302 
303  case 32:
305  break;
306 
307  case 33:
309  break;
310 
311  case 34:
313  break;
314 
315  default:
316  fasp_chkerr(ERROR_SOLVER_PRECTYPE, __FUNCTION__);
317  break;
318  }
319 
320  if ( prtlvl >= PRINT_MIN ) {
321  fasp_gettime(&setup_end);
322  setup_duration = setup_end - setup_start;
323  print_cputime("Setup totally", setup_duration);
324  }
325 
326 
327  // solver part
328  fasp_gettime(&solver_start);
329 
330  status=fasp_solver_bdcsr_itsolver(A,b,x, &prec,itparam);
331 
332  fasp_gettime(&solver_end);
333 
334  solver_duration = solver_end - solver_start;
335 
336  if ( prtlvl >= PRINT_MIN )
337  print_cputime("Krylov method totally", solver_duration);
338 
339  // clean
340  /* diagonal blocks are solved exactly */
341  if ( precond_type > 20 && precond_type < 30 ) {
342 #if WITH_UMFPACK
343  for (i=0; i<3; i++) fasp_umfpack_free_numeric(LU_diag[i]);
344 #endif
345  }
346  /* diagonal blocks are solved by AMG */
347  else if ( precond_type > 30 && precond_type < 40 ) {
348  for (i=0; i<3; i++) fasp_amg_data_free(mgl[i], amgparam);
349  if (mgl) fasp_mem_free(mgl);
350  }
351  else {
352  fasp_chkerr(ERROR_SOLVER_PRECTYPE, __FUNCTION__);
353  }
354 
355 #if DEBUG_MODE > 0
356  printf("### DEBUG: %s ...... [Finish]\n", __FUNCTION__);
357 #endif
358 
359  return status;
360 }
361 
384  dvector *b,
385  dvector *x,
386  itsolver_param *itparam,
387  AMG_param *amgparam,
388  dCSRmat *A_diag)
389 {
390  const SHORT prtlvl = itparam->print_level;
391  const SHORT precond_type = itparam->precond_type;
392 
393  INT status = FASP_SUCCESS;
394  REAL setup_start, setup_end, setup_duration;
395  REAL solver_start, solver_end, solver_duration;
396 
397 #if WITH_UMFPACK
398  void **LU_diag = (void **)fasp_mem_calloc(4, sizeof(void *));
399  INT i;
400 #endif
401 
402 #if DEBUG_MODE > 0
403  printf("### DEBUG: %s ...... [Start]\n", __FUNCTION__);
404 #endif
405 
406  /* setup preconditioner */
407  fasp_gettime(&setup_start);
408 
409  /* diagonal blocks are solved exactly */
410  if ( precond_type > 20 && precond_type < 30 ) {
411 
412 #if WITH_UMFPACK
413  // Need to sort the matrices local_A for UMFPACK format
414  dCSRmat A_tran;
415 
416  for (i=0; i<4; i++){
417 
418  A_tran = fasp_dcsr_create(A_diag[i].row, A_diag[i].col, A_diag[i].nnz);
419  fasp_dcsr_transz(&A_diag[i], NULL, &A_tran);
420  fasp_dcsr_cp(&A_tran, &A_diag[i]);
421 
422  printf("Factorization for %d-th diagnol: \n", i);
423  LU_diag[i] = fasp_umfpack_factorize(&A_diag[i], prtlvl);
424 
425  }
426 
427  fasp_dcsr_free(&A_tran);
428 #endif
429 
430  }
431  else {
432  fasp_chkerr(ERROR_SOLVER_PRECTYPE, __FUNCTION__);
433  }
434 
435  precond_block_data precdata;
436 
437  precdata.Abcsr = A;
438  precdata.A_diag = A_diag;
439 #if WITH_UMFPACK
440  precdata.LU_diag = LU_diag;
441 #endif
442  precdata.r = fasp_dvec_create(b->row);
443 
444  precond prec; prec.data = &precdata;
445 
446  switch (precond_type)
447  {
448  case 21:
450  break;
451 
452  case 22:
454  break;
455  }
456 
457  if ( prtlvl >= PRINT_MIN ) {
458  fasp_gettime(&setup_end);
459  setup_duration = setup_end - setup_start;
460  print_cputime("Setup totally", setup_duration);
461  }
462 
463  // solver part
464  fasp_gettime(&solver_start);
465 
466  status=fasp_solver_bdcsr_itsolver(A,b,x, &prec,itparam);
467 
468  fasp_gettime(&solver_end);
469 
470  solver_duration = solver_end - solver_start;
471 
472  if ( prtlvl >= PRINT_MIN )
473  print_cputime("Krylov method totally", solver_duration);
474 
475  // clean
476 #if WITH_UMFPACK
477  for (i=0; i<4; i++) fasp_umfpack_free_numeric(LU_diag[i]);
478 #endif
479 
480 #if DEBUG_MODE > 0
481  printf("### DEBUG: %s ...... [Finish]\n", __FUNCTION__);
482 #endif
483 
484  return status;
485 }
486 
510  dvector *b,
511  dvector *x,
512  itsolver_param *itparam,
513  INT NumLayers,
514  block_dCSRmat *Ai,
515  dCSRmat *local_A,
516  ivector *local_index)
517 {
518  const SHORT prtlvl = itparam->print_level;
519 
520  INT status = FASP_SUCCESS;
521  REAL setup_start, setup_end, setup_duration;
522  REAL solve_start, solve_end, solve_duration;
523 
524  /* setup preconditioner */
525  fasp_gettime(&setup_start);
526 
527  void **local_LU = NULL;
528 
529 #if DEBUG_MODE > 0
530  printf("### DEBUG: %s ...... [Start]\n", __FUNCTION__);
531 #endif
532 
533 #if WITH_UMFPACK
534  // Need to sort the matrices local_A for UMFPACK format
535  INT l;
536  dCSRmat A_tran;
537  local_LU = (void **)fasp_mem_calloc(NumLayers, sizeof(void *));
538 
539  for ( l=0; l<NumLayers; l++ ) {
540 
541  A_tran = fasp_dcsr_create(local_A[l].row, local_A[l].col, local_A[l].nnz);
542  fasp_dcsr_transz(&local_A[l], NULL, &A_tran);
543  fasp_dcsr_cp(&A_tran, &local_A[l]);
544 
545  printf("Factorization for layer %d: \n", l);
546  local_LU[l] = fasp_umfpack_factorize(&local_A[l], prtlvl);
547 
548  }
549 
550  fasp_dcsr_free(&A_tran);
551 #endif
552 
553  precond_sweeping_data precdata;
554  precdata.NumLayers = NumLayers;
555  precdata.A = A;
556  precdata.Ai = Ai;
557  precdata.local_A = local_A;
558  precdata.local_LU = local_LU;
559  precdata.local_index = local_index;
560  precdata.r = fasp_dvec_create(b->row);
561  precdata.w = (REAL *)fasp_mem_calloc(10*b->row,sizeof(REAL));
562 
563  precond prec; prec.data = &precdata;
564  prec.fct = fasp_precond_sweeping;
565 
566  if ( prtlvl >= PRINT_MIN ) {
567  fasp_gettime(&setup_end);
568  setup_duration = setup_end - setup_start;
569  print_cputime("Setup totally", setup_duration);
570  }
571 
572  /* solver part */
573  fasp_gettime(&solve_start);
574 
575  status = fasp_solver_bdcsr_itsolver(A,b,x, &prec,itparam);
576 
577  fasp_gettime(&solve_end);
578 
579  solve_duration = solve_end - solve_start;
580 
581  if ( prtlvl >= PRINT_MIN )
582  print_cputime("Krylov method totally", setup_duration+solve_duration);
583 
584  // clean
585 #if WITH_UMFPACK
586  for (l=0; l<NumLayers; l++) fasp_umfpack_free_numeric(local_LU[l]);
587 #endif
588 
589 #if DEBUG_MODE > 0
590  printf("### DEBUG: %s ...... [Finish]\n", __FUNCTION__);
591 #endif
592 
593  return status;
594 }
595 
596 /*---------------------------------*/
597 /*-- End of File --*/
598 /*---------------------------------*/
dvector fasp_dvec_create(const INT m)
Create dvector data space of REAL type.
Definition: vec.c:56
INT fasp_solver_bdcsr_pvgmres(block_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:393
block_dCSRmat * A
Definition: fasp_block.h:652
INT fasp_solver_bdcsr_pminres(block_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:499
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_chkerr(const SHORT status, const char *fctname)
Check error status and print out error messages before quit.
Definition: message.c:199
SHORT AMG_type
type of AMG method
Definition: fasp.h:586
void fasp_precond_block_upper_3(REAL *r, REAL *z, void *data)
block upper triangular preconditioning (3x3 block matrix, each diagonal block is solved exactly) ...
Definition: precond_bcsr.c:525
Parameters for AMG solver.
Definition: fasp.h:583
#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
Vector with n entries of INT type.
Definition: fasp.h:356
dCSRmat * A_diag
Definition: fasp_block.h:509
dvector b
pointer to the right-hand side at level level_num
Definition: fasp.h:744
AMG_data ** mgl
Definition: fasp_block.h:520
void fasp_precond_block_upper_3_amg(REAL *r, REAL *z, void *data)
block upper triangular preconditioning (3x3 block matrix, each diagonal block is solved AMG) ...
Definition: precond_bcsr.c:607
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
void * fasp_mem_calloc(LONGLONG size, INT type)
1M = 1024*1024
Definition: memory.c:60
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
#define INT
Definition: fasp.h:64
block_dCSRmat * Ai
Definition: fasp_block.h:653
#define FASP_SUCCESS
Definition of return status and error messages.
Definition: fasp_const.h:27
Preconditioner data and action.
Definition: fasp.h:999
#define ERROR_SOLVER_TYPE
Definition: fasp_const.h:47
INT restart
Definition: fasp.h:1112
void fasp_precond_block_SGS_3(REAL *r, REAL *z, void *data)
block symmetric GS preconditioning (3x3 block matrix, each diagonal block is solved exactly) ...
Definition: precond_bcsr.c:688
void fasp_mem_free(void *mem)
Free up previous allocated memory body.
Definition: memory.c:150
SHORT precond_type
Definition: fasp.h:1108
SHORT stop_type
Definition: fasp.h:1109
void * data
data for preconditioner, void pointer
Definition: fasp.h:1002
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_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
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
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 ERROR_SOLVER_PRECTYPE
Definition: fasp_const.h:48
void fasp_gettime(REAL *time)
Get system time.
Definition: timing.c:28
INT fasp_solver_bdcsr_krylov_block_3(block_dCSRmat *A, dvector *b, dvector *x, itsolver_param *itparam, AMG_param *amgparam, dCSRmat *A_diag)
Solve Ax = b by standard Krylov methods.
INT fasp_solver_bdcsr_krylov_sweeping(block_dCSRmat *A, dvector *b, dvector *x, itsolver_param *itparam, INT NumLayers, block_dCSRmat *Ai, dCSRmat *local_A, ivector *local_index)
Solve Ax = b by standard Krylov methods.
#define SOLVER_GMRES
Definition: fasp_const.h:106
#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
void fasp_precond_block_lower_3_amg(REAL *r, REAL *z, void *data)
block lower triangular preconditioning (3x3 block matrix, each diagonal block is solved by AMG) ...
Definition: precond_bcsr.c:353
Sparse matrix of REAL type in CSR format.
Definition: fasp.h:148
void fasp_precond_block_diag_3(REAL *r, REAL *z, void *data)
block diagonal preconditioning (3x3 block matrix, each diagonal block is solved exactly) ...
Definition: precond_bcsr.c:26
void fasp_precond_block_diag_3_amg(REAL *r, REAL *z, void *data)
block diagonal preconditioning (3x3 block matrix, each diagonal block is solved by AMG) ...
Definition: precond_bcsr.c:110
INT row
number of rows
Definition: fasp.h:345
REAL tol
Definition: fasp.h:1111
SHORT print_level
Definition: fasp.h:1113
Block REAL CSR matrix format.
Definition: fasp_block.h:84
#define SHORT
FASP integer and floating point numbers.
Definition: fasp.h:63
Data passed to the preconditioner for sweeping preconditioning.
Definition: fasp_block.h:648
INT fasp_solver_bdcsr_itsolver(block_dCSRmat *A, dvector *b, dvector *x, precond *pc, itsolver_param *itparam)
Solve Ax = b by standard Krylov methods.
Definition: itsolver_bcsr.c:36
dCSRmat A
pointer to the matrix at level level_num
Definition: fasp.h:735
INT fasp_solver_bdcsr_krylov(block_dCSRmat *A, dvector *b, dvector *x, itsolver_param *itparam)
Solve Ax = b by standard Krylov methods.
void print_cputime(const char *message, const REAL cputime)
Print CPU walltime.
Definition: message.c:165
INT fasp_solver_bdcsr_pgmres(block_dCSRmat *A, dvector *b, dvector *x, precond *pc, const REAL tol, const INT MaxIt, const SHORT restart, const SHORT stop_type, const SHORT prtlvl)
Preconditioned GMRES method for solving Au=b.
Definition: pgmres.c:356
INT fasp_solver_bdcsr_pvfgmres(block_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:712
Data for AMG solvers.
Definition: fasp.h:722
void fasp_precond_block_SGS_3_amg(REAL *r, REAL *z, void *data)
block symmetric GS preconditioning (3x3 block matrix, each diagonal block is solved exactly) ...
Definition: precond_bcsr.c:806
INT fasp_solver_bdcsr_pbcgs(block_dCSRmat *A, dvector *b, dvector *u, precond *pc, const REAL tol, const INT MaxIt, const SHORT stop_type, const SHORT prtlvl)
A preconditioned BiCGstab method for solving Au=b.
Definition: pbcgs.c:774
void fasp_amg_data_free(AMG_data *mgl, AMG_param *param)
Free AMG_data data memeory space.
Definition: init.c:185
INT fasp_solver_bdcsr_krylov_block_4(block_dCSRmat *A, dvector *b, dvector *x, itsolver_param *itparam, AMG_param *amgparam, dCSRmat *A_diag)
Solve Ax = b by standard Krylov methods.
Header file for FASP block matrices.
void fasp_precond_block_diag_4(REAL *r, REAL *z, void *data)
block diagonal preconditioning (4x4 block matrix, each diagonal block is solved exactly) ...
Definition: precond_bcsr.c:175
AMG_param * amgparam
Definition: fasp_block.h:521
#define SOLVER_VGMRES
Definition: fasp_const.h:107
block_dCSRmat * Abcsr
Definition: fasp_block.h:507
void fasp_precond_block_lower_3(REAL *r, REAL *z, void *data)
block lower triangular preconditioning (3x3 block matrix, each diagonal block is solved exactly) ...
Definition: precond_bcsr.c:271
#define SOLVER_MinRes
Definition: fasp_const.h:105
Data passed to the preconditioner for block preconditioning for block_dCSRmat format.
Definition: fasp_block.h:502
#define SA_AMG
Definition: fasp_const.h:163
void fasp_precond_block_lower_4(REAL *r, REAL *z, void *data)
block lower triangular preconditioning (4x4 block matrix, each diagonal block is solved exactly) ...
Definition: precond_bcsr.c:427
void fasp_dcsr_free(dCSRmat *A)
Free CSR sparse matrix data memory space.
Definition: sparse_csr.c:166
#define UA_AMG
Definition: fasp_const.h:164
void fasp_precond_sweeping(REAL *r, REAL *z, void *data)
sweeping preconditioner for Maxwell equations
Definition: precond_bcsr.c:916
#define PRINT_MIN
Definition: fasp_const.h:80
SHORT itsolver_type
Definition: fasp.h:1107