Fast Auxiliary Space Preconditioning  1.8.4 Feb/15/2016
itsolver_bsr.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 
38  dvector *b,
39  dvector *x,
40  precond *pc,
41  itsolver_param *itparam)
42 {
43  const SHORT prtlvl = itparam->print_level;
44  const SHORT itsolver_type = itparam->itsolver_type;
45  const SHORT stop_type = itparam->stop_type;
46  const SHORT restart = itparam->restart;
47  const INT MaxIt = itparam->maxit;
48  const REAL tol = itparam->tol;
49 
50  // Local variables
51  INT iter = ERROR_SOLVER_TYPE;
52  REAL solver_start, solver_end, solver_duration;
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  fasp_gettime(&solver_start);
60 
61  /* Safe-guard checks on parameters */
62  ITS_CHECK ( MaxIt, tol );
63 
64  switch (itsolver_type) {
65 
66  case SOLVER_CG:
67  if ( prtlvl > PRINT_NONE ) printf("\nCalling PCG solver (BSR) ...\n");
68  iter = fasp_solver_dbsr_pcg(A, b, x, pc, tol, MaxIt, stop_type, prtlvl);
69  break;
70 
71  case SOLVER_BiCGstab:
72  if ( prtlvl > PRINT_NONE ) printf("\nCalling BiCGstab solver (BSR) ...\n");
73  iter = fasp_solver_dbsr_pbcgs(A, b, x, pc, tol, MaxIt, stop_type, prtlvl);
74  break;
75 
76  case SOLVER_GMRES:
77  if ( prtlvl > PRINT_NONE ) printf("\nCalling GMRES solver (BSR) ...\n");
78  iter = fasp_solver_dbsr_pgmres(A, b, x, pc, tol, MaxIt, restart, stop_type, prtlvl);
79  break;
80 
81  case SOLVER_VGMRES:
82  if ( prtlvl > PRINT_NONE ) printf("\nCalling vGMRES solver (BSR) ...\n");
83  iter = fasp_solver_dbsr_pvgmres(A, b, x, pc, tol, MaxIt, restart, stop_type, prtlvl);
84  break;
85 
86  case SOLVER_VFGMRES:
87  if ( prtlvl > PRINT_NONE ) printf("\nCalling vFGMRes solver (BSR) ...\n");
88  iter = fasp_solver_dbsr_pvfgmres(A, b, x, pc, tol, MaxIt, restart, stop_type, prtlvl);
89  break;
90 
91  default:
92  printf("### ERROR: Unknown itertive solver type %d!\n", itsolver_type);
93 
94  }
95 
96  if ( (prtlvl > PRINT_MIN) && (iter >= 0) ) {
97  fasp_gettime(&solver_end);
98  solver_duration = solver_end - solver_start;
99  print_cputime("Iterative method", solver_duration);
100  }
101 
102 #if DEBUG_MODE > 0
103  printf("### DEBUG: %s ...... [Finish]\n", __FUNCTION__);
104 #endif
105 
106  return iter;
107 }
108 
126  dvector *b,
127  dvector *x,
128  itsolver_param *itparam)
129 {
130  const SHORT prtlvl = itparam->print_level;
131  INT status = FASP_SUCCESS;
132  REAL solver_start, solver_end;
133 
134 #if DEBUG_MODE > 0
135  printf("### DEBUG: %s ...... [Start]\n", __FUNCTION__);
136 #endif
137 
138  // solver part
139  fasp_gettime(&solver_start);
140 
141  status=fasp_solver_dbsr_itsolver(A,b,x,NULL,itparam);
142 
143  fasp_gettime(&solver_end);
144 
145  if ( prtlvl > PRINT_NONE ) {
146  REAL solver_duration = solver_end - solver_start;
147  print_cputime("Krylov method totally", solver_duration);
148  }
149 
150 #if DEBUG_MODE > 0
151  printf("### DEBUG: %s ...... [Finish]\n", __FUNCTION__);
152 #endif
153 
154  return status;
155 }
156 
177  dvector *b,
178  dvector *x,
179  itsolver_param *itparam)
180 {
181  const SHORT prtlvl = itparam->print_level;
182  INT status = FASP_SUCCESS;
183  REAL solver_start, solver_end;
184 
185  INT nb=A->nb,i,k;
186  INT nb2=nb*nb;
187  INT ROW=A->ROW;
188 
189 #ifdef _OPENMP
190  // variables for OpenMP
191  INT myid, mybegin, myend;
192  INT nthreads = FASP_GET_NUM_THREADS();
193 #endif
194  // setup preconditioner
195  precond_diagbsr diag;
196  fasp_dvec_alloc(ROW*nb2, &(diag.diag));
197 
198 #if DEBUG_MODE > 0
199  printf("### DEBUG: %s ...... [Start]\n", __FUNCTION__);
200 #endif
201 
202  // get all the diagonal sub-blocks
203 #ifdef _OPENMP
204  if (ROW > OPENMP_HOLDS) {
205 #pragma omp parallel for private(myid, mybegin, myend, i, k)
206  for (myid=0; myid<nthreads; ++myid) {
207  FASP_GET_START_END(myid, nthreads, ROW, &mybegin, &myend);
208  for (i = mybegin; i < myend; ++i) {
209  for (k = A->IA[i]; k < A->IA[i+1]; ++k) {
210  if (A->JA[k] == i)
211  memcpy(diag.diag.val+i*nb2, A->val+k*nb2, nb2*sizeof(REAL));
212  }
213  }
214  }
215  }
216  else {
217 #endif
218  for (i = 0; i < ROW; ++i) {
219  for (k = A->IA[i]; k < A->IA[i+1]; ++k) {
220  if (A->JA[k] == i)
221  memcpy(diag.diag.val+i*nb2, A->val+k*nb2, nb2*sizeof(REAL));
222  }
223  }
224 #ifdef _OPENMP
225  }
226 #endif
227 
228  diag.nb=nb;
229 
230 #ifdef _OPENMP
231 #pragma omp parallel for if(ROW>OPENMP_HOLDS)
232 #endif
233  for (i=0; i<ROW; ++i){
234  fasp_blas_smat_inv(&(diag.diag.val[i*nb2]), nb);
235  }
236 
237  precond *pc = (precond *)fasp_mem_calloc(1,sizeof(precond));
238  pc->data = &diag;
240 
241  // solver part
242  fasp_gettime(&solver_start);
243 
244  status=fasp_solver_dbsr_itsolver(A,b,x,pc,itparam);
245 
246  fasp_gettime(&solver_end);
247 
248  if ( prtlvl > PRINT_NONE ) {
249  REAL solver_duration = solver_end - solver_start;
250  print_cputime("Diag_Krylov method totally", solver_duration);
251  }
252 
253  // clean up
254  fasp_dvec_free(&(diag.diag));
255 
256 #if DEBUG_MODE > 0
257  printf("### DEBUG: %s ...... [Finish]\n", __FUNCTION__);
258 #endif
259 
260  return status;
261 }
262 
281  dvector *b,
282  dvector *x,
283  itsolver_param *itparam,
284  ILU_param *iluparam)
285 {
286  const SHORT prtlvl = itparam->print_level;
287  REAL solver_start, solver_end;
288  INT status = FASP_SUCCESS;
289 
290  ILU_data LU;
291  precond pc;
292 
293 #if DEBUG_MODE > 0
294  printf("### DEBUG: %s ...... [Start]\n", __FUNCTION__);
295  printf("### DEBUG: matrix size: %d %d %d\n", A->ROW, A->COL, A->NNZ);
296  printf("### DEBUG: rhs/sol size: %d %d\n", b->row, x->row);
297 #endif
298 
299  // ILU setup for whole matrix
300  if ( (status = fasp_ilu_dbsr_setup(A,&LU,iluparam)) < 0 ) goto FINISHED;
301 
302  // check iludata
303  if ( (status = fasp_mem_iludata_check(&LU)) < 0 ) goto FINISHED;
304 
305  // set preconditioner
306  pc.data = &LU; pc.fct = fasp_precond_dbsr_ilu;
307 
308  // solve
309  fasp_gettime(&solver_start);
310 
311  status=fasp_solver_dbsr_itsolver(A,b,x,&pc,itparam);
312 
313  fasp_gettime(&solver_end);
314 
315  if ( prtlvl > PRINT_NONE ) {
316  REAL solver_duration = solver_end - solver_start;
317  print_cputime("ILUk_Krylov method totally", solver_duration);
318  }
319 
320 FINISHED:
321  fasp_ilu_data_free(&LU);
322 
323 #if DEBUG_MODE > 0
324  printf("### DEBUG: %s ...... [Finish]\n", __FUNCTION__);
325 #endif
326 
327  return status;
328 }
329 
348  dvector *b,
349  dvector *x,
350  itsolver_param *itparam,
351  AMG_param *amgparam)
352 {
353  //--------------------------------------------------------------
354  // Part 1: prepare
355  // --------------------------------------------------------------
357  const SHORT prtlvl = itparam->print_level;
358  const SHORT max_levels = amgparam->max_levels;
359 
360  // return variable
361  INT status = FASP_SUCCESS;
362 
363  // data of AMG
364  AMG_data_bsr *mgl=fasp_amg_data_bsr_create(max_levels);
365 
366  // timing
367  REAL setup_start, setup_end, setup_duration;
368  REAL solver_start, solver_end, solver_duration;
369 
370 #if DEBUG_MODE > 0
371  printf("### DEBUG: %s ...... [Start]\n", __FUNCTION__);
372 #endif
373 
374  //--------------------------------------------------------------
375  //Part 2: set up the preconditioner
376  //--------------------------------------------------------------
377  fasp_gettime(&setup_start);
378 
379  // initialize A, b, x for mgl[0]
380  mgl[0].A = fasp_dbsr_create(A->ROW, A->COL, A->NNZ, A->nb, A->storage_manner);
381  mgl[0].b = fasp_dvec_create(mgl[0].A.ROW*mgl[0].A.nb);
382  mgl[0].x = fasp_dvec_create(mgl[0].A.COL*mgl[0].A.nb);
383 
384  fasp_dbsr_cp(A, &(mgl[0].A));
385 
386  switch (amgparam->AMG_type) {
387 
388  case SA_AMG: // Smoothed Aggregation AMG
389  status = fasp_amg_setup_sa_bsr(mgl, amgparam); break;
390 
391  default:
392  status = fasp_amg_setup_ua_bsr(mgl, amgparam); break;
393 
394  }
395 
396  if (status < 0) goto FINISHED;
397 
398  fasp_gettime(&setup_end);
399 
400  setup_duration = setup_end - setup_start;
401 
402  precond_data_bsr precdata;
403  precdata.print_level = amgparam->print_level;
404  precdata.maxit = amgparam->maxit;
405  precdata.tol = amgparam->tol;
406  precdata.cycle_type = amgparam->cycle_type;
407  precdata.smoother = amgparam->smoother;
408  precdata.presmooth_iter = amgparam->presmooth_iter;
409  precdata.postsmooth_iter = amgparam->postsmooth_iter;
410  precdata.coarsening_type = amgparam->coarsening_type;
411  precdata.relaxation = amgparam->relaxation;
412  precdata.coarse_scaling = amgparam->coarse_scaling;
413  precdata.amli_degree = amgparam->amli_degree;
414  precdata.amli_coef = amgparam->amli_coef;
415  precdata.tentative_smooth = amgparam->tentative_smooth;
416  precdata.max_levels = mgl[0].num_levels;
417  precdata.mgl_data = mgl;
418  precdata.A = A;
419 
420  if (status < 0) goto FINISHED;
421 
422  precond prec;
423  prec.data = &precdata;
424  switch (amgparam->cycle_type) {
425  case NL_AMLI_CYCLE: // Nonlinear AMLI AMG
427  break;
428  default: // V,W-Cycle AMG
429  prec.fct = fasp_precond_dbsr_amg;
430  break;
431  }
432 
433  if ( prtlvl >= PRINT_MIN ) {
434  print_cputime("BSR AMG setup", setup_duration);
435  }
436 
437  //--------------------------------------------------------------
438  // Part 3: solver
439  //--------------------------------------------------------------
440  fasp_gettime(&solver_start);
441 
442  status=fasp_solver_dbsr_itsolver(A,b,x,&prec,itparam);
443 
444  fasp_gettime(&solver_end);
445 
446  solver_duration = solver_end - solver_start;
447 
448  if ( prtlvl >= PRINT_MIN ) {
449  print_cputime("BSR AMG Krylov method totally", setup_duration+solver_duration);
450  }
451 
452 FINISHED:
454 
455 #if DEBUG_MODE > 0
456  printf("### DEBUG: %s ...... [Finish]\n", __FUNCTION__);
457 #endif
458 
459  if (status == ERROR_ALLOC_MEM) goto MEMORY_ERROR;
460  return status;
461 
462 MEMORY_ERROR:
463  printf("### ERROR: %s cannot allocate memory!\n", __FUNCTION__);
464  exit(status);
465 }
466 
489  dvector *b,
490  dvector *x,
491  itsolver_param *itparam,
492  AMG_param *amgparam,
493  dCSRmat *A_nk,
494  dCSRmat *P_nk,
495  dCSRmat *R_nk)
496 {
497  //--------------------------------------------------------------
498  // Part 1: prepare
499  // --------------------------------------------------------------
500  // parameters of iterative method
501  const SHORT prtlvl = itparam->print_level;
502  const SHORT max_levels = amgparam->max_levels;
503 
504  // return variable
505  INT status = FASP_SUCCESS;
506 
507  // data of AMG
508  AMG_data_bsr *mgl=fasp_amg_data_bsr_create(max_levels);
509 
510  // timing
511  REAL setup_start, setup_end, setup_duration;
512  REAL solver_start, solver_end, solver_duration;
513 
514 #if DEBUG_MODE > 0
515  printf("### DEBUG: %s ...... [Start]\n", __FUNCTION__);
516 #endif
517 
518  //--------------------------------------------------------------
519  //Part 2: set up the preconditioner
520  //--------------------------------------------------------------
521  fasp_gettime(&setup_start);
522 
523  // initialize A, b, x for mgl[0]
524  mgl[0].A = fasp_dbsr_create(A->ROW, A->COL, A->NNZ, A->nb, A->storage_manner);
525  fasp_dbsr_cp(A, &(mgl[0].A));
526  mgl[0].b = fasp_dvec_create(mgl[0].A.ROW*mgl[0].A.nb);
527  mgl[0].x = fasp_dvec_create(mgl[0].A.COL*mgl[0].A.nb);
528 
529  // near kernel space
530  mgl[0].A_nk = NULL;
531  mgl[0].P_nk = P_nk;
532  mgl[0].R_nk = R_nk;
533 
534  switch (amgparam->AMG_type) {
535 
536  case SA_AMG: // Smoothed Aggregation AMG
537  status = fasp_amg_setup_sa_bsr(mgl, amgparam); break;
538 
539  default:
540  status = fasp_amg_setup_ua_bsr(mgl, amgparam); break;
541 
542  }
543 
544  if (status < 0) goto FINISHED;
545 
546  fasp_gettime(&setup_end);
547 
548  setup_duration = setup_end - setup_start;
549 
550  precond_data_bsr precdata;
551  precdata.print_level = amgparam->print_level;
552  precdata.maxit = amgparam->maxit;
553  precdata.tol = amgparam->tol;
554  precdata.cycle_type = amgparam->cycle_type;
555  precdata.smoother = amgparam->smoother;
556  precdata.presmooth_iter = amgparam->presmooth_iter;
557  precdata.postsmooth_iter = amgparam->postsmooth_iter;
558  precdata.coarsening_type = amgparam->coarsening_type;
559  precdata.relaxation = amgparam->relaxation;
560  precdata.coarse_scaling = amgparam->coarse_scaling;
561  precdata.amli_degree = amgparam->amli_degree;
562  precdata.amli_coef = amgparam->amli_coef;
563  precdata.tentative_smooth = amgparam->tentative_smooth;
564  precdata.max_levels = mgl[0].num_levels;
565  precdata.mgl_data = mgl;
566  precdata.A = A;
567 
568 #if WITH_UMFPACK // use UMFPACK directly
569  dCSRmat A_tran;
570  A_tran = fasp_dcsr_create(A_nk->row, A_nk->col, A_nk->nnz);
571  fasp_dcsr_transz(A_nk, NULL, &A_tran);
572  // It is equivalent to do transpose and then sort
573  // fasp_dcsr_trans(A_nk, &A_tran);
574  // fasp_dcsr_sort(&A_tran);
575  precdata.A_nk = &A_tran;
576 #else
577  precdata.A_nk = A_nk;
578 #endif
579 
580  precdata.P_nk = P_nk;
581  precdata.R_nk = R_nk;
582 
583  if (status < 0) goto FINISHED;
584 
585  precond prec;
586  prec.data = &precdata;
587 
589 
590  if ( prtlvl >= PRINT_MIN ) {
591  print_cputime("BSR AMG setup", setup_duration);
592  }
593 
594  //--------------------------------------------------------------
595  // Part 3: solver
596  //--------------------------------------------------------------
597  fasp_gettime(&solver_start);
598 
599  status=fasp_solver_dbsr_itsolver(A,b,x,&prec,itparam);
600 
601  fasp_gettime(&solver_end);
602 
603  solver_duration = solver_end - solver_start;
604 
605  if ( prtlvl >= PRINT_MIN ) {
606  print_cputime("BSR AMG NK Krylov method totally", setup_duration+solver_duration);
607  }
608 
609 FINISHED:
611 
612 #if DEBUG_MODE > 0
613  printf("### DEBUG: %s ...... [Finish]\n", __FUNCTION__);
614 #endif
615 
616 #if WITH_UMFPACK // use UMFPACK directly
617  fasp_dcsr_free(&A_tran);
618 #endif
619  if (status == ERROR_ALLOC_MEM) goto MEMORY_ERROR;
620  return status;
621 
622 MEMORY_ERROR:
623  printf("### ERROR: %s cannot allocate memory!\n", __FUNCTION__);
624  exit(status);
625 }
626 
648  dvector *b,
649  dvector *x,
650  itsolver_param *itparam,
651  AMG_param *amgparam,
652  const INT nk_dim,
653  dvector *nk)
654 {
655  //--------------------------------------------------------------
656  // Part 1: prepare
657  // --------------------------------------------------------------
659  const SHORT prtlvl = itparam->print_level;
660  const SHORT max_levels = amgparam->max_levels;
661 
662  // local variable
663  INT i;
664 
665  // return variable
666  INT status = FASP_SUCCESS;
667 
668  // data of AMG
669  AMG_data_bsr *mgl=fasp_amg_data_bsr_create(max_levels);
670 
671  // timing
672  REAL setup_start, setup_end, setup_duration;
673  REAL solver_start, solver_end, solver_duration;
674 
675 #if DEBUG_MODE > 0
676  printf("### DEBUG: %s ...... [Start]\n", __FUNCTION__);
677 #endif
678 
679  //--------------------------------------------------------------
680  //Part 2: set up the preconditioner
681  //--------------------------------------------------------------
682  fasp_gettime(&setup_start);
683 
684  // initialize A, b, x for mgl[0]
685  mgl[0].A = fasp_dbsr_create(A->ROW, A->COL, A->NNZ, A->nb, A->storage_manner);
686  fasp_dbsr_cp(A, &(mgl[0].A));
687  mgl[0].b = fasp_dvec_create(mgl[0].A.ROW*mgl[0].A.nb);
688  mgl[0].x = fasp_dvec_create(mgl[0].A.COL*mgl[0].A.nb);
689 
690  /*-----------------------*/
691  /*-- setup null spaces --*/
692  /*-----------------------*/
693 
694  // null space for whole Jacobian
695  mgl[0].near_kernel_dim = nk_dim;
696  mgl[0].near_kernel_basis = (REAL **)fasp_mem_calloc(mgl->near_kernel_dim, sizeof(REAL*));
697 
698  for ( i=0; i < mgl->near_kernel_dim; ++i ) mgl[0].near_kernel_basis[i] = nk[i].val;
699 
700  switch (amgparam->AMG_type) {
701 
702  case SA_AMG: // Smoothed Aggregation AMG
703  status = fasp_amg_setup_sa_bsr(mgl, amgparam); break;
704 
705  default:
706  status = fasp_amg_setup_ua_bsr(mgl, amgparam); break;
707 
708  }
709 
710  if (status < 0) goto FINISHED;
711 
712  fasp_gettime(&setup_end);
713 
714  setup_duration = setup_end - setup_start;
715 
716  precond_data_bsr precdata;
717  precdata.print_level = amgparam->print_level;
718  precdata.maxit = amgparam->maxit;
719  precdata.tol = amgparam->tol;
720  precdata.cycle_type = amgparam->cycle_type;
721  precdata.smoother = amgparam->smoother;
722  precdata.presmooth_iter = amgparam->presmooth_iter;
723  precdata.postsmooth_iter = amgparam->postsmooth_iter;
724  precdata.coarsening_type = amgparam->coarsening_type;
725  precdata.relaxation = amgparam->relaxation;
726  precdata.coarse_scaling = amgparam->coarse_scaling;
727  precdata.amli_degree = amgparam->amli_degree;
728  precdata.amli_coef = amgparam->amli_coef;
729  precdata.tentative_smooth = amgparam->tentative_smooth;
730  precdata.max_levels = mgl[0].num_levels;
731  precdata.mgl_data = mgl;
732  precdata.A = A;
733 
734  if (status < 0) goto FINISHED;
735 
736  precond prec;
737  prec.data = &precdata;
738  switch (amgparam->cycle_type) {
739  case NL_AMLI_CYCLE: // Nonlinear AMLI AMG
741  break;
742  default: // V,W-Cycle AMG
743  prec.fct = fasp_precond_dbsr_amg;
744  break;
745  }
746 
747  if ( prtlvl >= PRINT_MIN ) {
748  print_cputime("BSR AMG setup", setup_duration);
749  }
750 
751  //--------------------------------------------------------------
752  // Part 3: solver
753  //--------------------------------------------------------------
754  fasp_gettime(&solver_start);
755 
756  status=fasp_solver_dbsr_itsolver(A,b,x,&prec,itparam);
757 
758  fasp_gettime(&solver_end);
759 
760  solver_duration = solver_end - solver_start;
761 
762  if ( prtlvl >= PRINT_MIN ) {
763  print_cputime("BSR AMG Krylov method totally", setup_duration+solver_duration);
764  }
765 
766 FINISHED:
768 
769 #if DEBUG_MODE > 0
770  printf("### DEBUG: %s ...... [Finish]\n", __FUNCTION__);
771 #endif
772 
773  if (status == ERROR_ALLOC_MEM) goto MEMORY_ERROR;
774  return status;
775 
776 MEMORY_ERROR:
777  printf("### ERROR: %s cannot allocate memory!\n", __FUNCTION__);
778  exit(status);
779 }
780 
781 /*---------------------------------*/
782 /*-- End of File --*/
783 /*---------------------------------*/
void fasp_precond_dbsr_amg_nk(REAL *r, REAL *z, void *data)
AMG with extra near kernel solve preconditioner.
Definition: precond_bsr.c:643
dvector fasp_dvec_create(const INT m)
Create dvector data space of REAL type.
Definition: vec.c:56
dCSRmat * A_nk
Matrix data for near kernal.
Definition: fasp_block.h:382
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
AMG_data_bsr * fasp_amg_data_bsr_create(SHORT max_levels)
Create and initialize AMG_data data sturcture for AMG/SAMG (BSR format)
Definition: init.c:86
SHORT AMG_type
type of AMG method
Definition: fasp.h:586
Parameters for AMG solver.
Definition: fasp.h:583
void fasp_dbsr_cp(dBSRmat *A, dBSRmat *B)
copy a dCSRmat to a new one B=A
Definition: sparse_bsr.c:181
SHORT postsmooth_iter
number of postsmoothers
Definition: fasp.h:619
REAL tentative_smooth
relaxation parameter for smoothing the tentative prolongation
Definition: fasp.h:676
SHORT smoother
smoother type
Definition: fasp.h:610
INT maxit
max number of iterations of AMG
Definition: fasp.h:592
REAL tol
tolerance for AMG preconditioner
Definition: fasp_block.h:326
SHORT coarsening_type
coarsening type
Definition: fasp.h:643
void fasp_amg_data_bsr_free(AMG_data_bsr *mgl)
Free AMG_data_bsr data memeory space.
Definition: init.c:256
SHORT amli_degree
degree of the polynomial used by AMLI cycle
Definition: fasp_block.h:356
#define REAL
Definition: fasp.h:67
SHORT amli_degree
degree of the polynomial used by AMLI cycle
Definition: fasp.h:634
dBSRmat fasp_dbsr_create(const INT ROW, const INT COL, const INT NNZ, const INT nb, const INT storage_manner)
Create BSR sparse matrix data memory space.
Definition: sparse_bsr.c:36
INT fasp_solver_dbsr_krylov_amg_nk(dBSRmat *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 with extra near kernel solve preconditioned Krylov methods.
Definition: itsolver_bsr.c:488
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
INT fasp_blas_smat_inv(REAL *a, const INT n)
Compute the inverse matrix of a small full matrix a.
Definition: smat.c:554
INT fasp_solver_dbsr_pgmres(dBSRmat *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:659
INT fasp_solver_dbsr_pbcgs(dBSRmat *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:431
INT nb
dimension of each sub-block
Definition: fasp_block.h:56
INT fasp_solver_dbsr_krylov_amg(dBSRmat *A, dvector *b, dvector *x, itsolver_param *itparam, AMG_param *amgparam)
Solve Ax=b by AMG preconditioned Krylov methods.
Definition: itsolver_bsr.c:347
REAL * amli_coef
coefficients of the polynomial used by AMLI cycle
Definition: fasp_block.h:359
void fasp_dvec_free(dvector *u)
Free vector data space of REAL type.
Definition: vec.c:139
Block sparse row storage matrix of REAL type.
Definition: fasp_block.h:44
REAL * val
actual vector entries
Definition: fasp.h:348
INT fasp_solver_dbsr_krylov_diag(dBSRmat *A, dvector *b, dvector *x, itsolver_param *itparam)
Solve Ax=b by diagonal preconditioned Krylov methods.
Definition: itsolver_bsr.c:176
dCSRmat * P_nk
Prolongation for near kernal.
Definition: fasp_block.h:385
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
SHORT cycle_type
AMG cycle type.
Definition: fasp_block.h:329
INT max_levels
max number of AMG levels
Definition: fasp_block.h:323
SHORT print_level
print level in AMG preconditioner
Definition: fasp_block.h:317
#define INT
Definition: fasp.h:64
dBSRmat * A
Matrix data.
Definition: fasp_block.h:377
#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
void FASP_GET_START_END(INT procid, INT nprocs, INT n, INT *start, INT *end)
Assign Load to each thread.
Definition: threads.c:83
Data for multigrid levels. (BSR format)
Definition: fasp_block.h:198
SHORT fasp_ilu_dbsr_setup(dBSRmat *A, ILU_data *iludata, ILU_param *iluparam)
Get ILU decoposition of a BSR matrix A.
Definition: ilu_setup_bsr.c:45
void fasp_precond_dbsr_nl_amli(REAL *r, REAL *z, void *data)
Nonlinear AMLI-cycle AMG preconditioner.
Definition: precond_bsr.c:607
SHORT stop_type
Definition: fasp.h:1109
REAL * val
Definition: fasp_block.h:67
SHORT fasp_amg_setup_sa_bsr(AMG_data_bsr *mgl, AMG_param *param)
Set up phase of smoothed aggregation AMG (BSR format)
Definition: amg_setup_sa.c:85
void * data
data for preconditioner, void pointer
Definition: fasp.h:1002
INT NNZ
number of nonzero sub-blocks in matrix A, NNZ
Definition: fasp_block.h:53
Data passed to the preconditioners.
Definition: fasp_block.h:311
INT fasp_solver_dbsr_krylov_ilu(dBSRmat *A, dvector *b, dvector *x, itsolver_param *itparam, ILU_param *iluparam)
Solve Ax=b by ILUs preconditioned Krylov methods.
Definition: itsolver_bsr.c:280
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
#define SOLVER_BiCGstab
Definition: fasp_const.h:104
#define SOLVER_VFGMRES
Definition: fasp_const.h:108
#define OPENMP_HOLDS
Definition: fasp_const.h:248
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
INT col
column of matrix A, n
Definition: fasp.h:154
INT storage_manner
storage manner for each sub-block
Definition: fasp_block.h:59
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
SHORT coarse_scaling
switch of scaling of the coarse grid correction
Definition: fasp_block.h:353
REAL * amli_coef
coefficients of the polynomial used by AMLI cycle
Definition: fasp.h:637
void fasp_gettime(REAL *time)
Get system time.
Definition: timing.c:28
#define ERROR_ALLOC_MEM
Definition: fasp_const.h:37
#define SOLVER_GMRES
Definition: fasp_const.h:106
SHORT smoother
AMG smoother type.
Definition: fasp_block.h:332
#define PRINT_NONE
Print level for all subroutines – not including DEBUG output.
Definition: fasp_const.h:79
Data passed to diagnal preconditioner for dBSRmat matrices.
Definition: fasp_block.h:293
REAL relaxation
relaxation parameter for SOR smoother
Definition: fasp.h:622
INT row
row number of matrix A, m
Definition: fasp.h:151
Parameters passed to iterative solvers.
Definition: fasp.h:1105
INT fasp_solver_dbsr_krylov_nk_amg(dBSRmat *A, dvector *b, dvector *x, itsolver_param *itparam, AMG_param *amgparam, const INT nk_dim, dvector *nk)
Solve Ax=b by AMG preconditioned Krylov methods with extra kernal space.
Definition: itsolver_bsr.c:647
void fasp_precond_dbsr_ilu(REAL *r, REAL *z, void *data)
ILU preconditioner.
Definition: precond_bsr.c:306
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
AMG_data_bsr * mgl_data
AMG preconditioner data.
Definition: fasp_block.h:368
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
INT COL
number of cols of sub-blocks in matrix A, N
Definition: fasp_block.h:50
SHORT print_level
Definition: fasp.h:1113
Parameters for ILU.
Definition: fasp.h:374
#define SHORT
FASP integer and floating point numbers.
Definition: fasp.h:63
SHORT coarsening_type
coarsening type
Definition: fasp_block.h:344
dvector diag
diagnal elements
Definition: fasp_block.h:299
INT fasp_solver_dbsr_pcg(dBSRmat *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:373
void print_cputime(const char *message, const REAL cputime)
Print CPU walltime.
Definition: message.c:165
SHORT print_level
print level for AMG
Definition: fasp.h:589
INT fasp_solver_dbsr_itsolver(dBSRmat *A, dvector *b, dvector *x, precond *pc, itsolver_param *itparam)
Solve Ax=b by preconditioned Krylov methods for BSR matrices.
Definition: itsolver_bsr.c:37
void fasp_precond_dbsr_amg(REAL *r, REAL *z, void *data)
AMG preconditioner.
Definition: precond_bsr.c:563
INT fasp_solver_dbsr_krylov(dBSRmat *A, dvector *b, dvector *x, itsolver_param *itparam)
Solve Ax=b by standard Krylov methods for BSR matrices.
Definition: itsolver_bsr.c:125
INT fasp_solver_dbsr_pvfgmres(dBSRmat *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:382
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
REAL tentative_smooth
smooth factor for smoothing the tentative prolongation
Definition: fasp_block.h:362
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
void fasp_ilu_data_free(ILU_data *ILUdata)
Create ILU_data sturcture.
Definition: init.c:300
INT maxit
max number of iterations of AMG preconditioner
Definition: fasp_block.h:320
INT ROW
number of rows of sub-blocks in matrix A, M
Definition: fasp_block.h:47
SHORT fasp_amg_setup_ua_bsr(AMG_data_bsr *mgl, AMG_param *param)
Set up phase of unsmoothed aggregation AMG (BSR format)
Definition: amg_setup_ua.c:69
INT * JA
Definition: fasp_block.h:74
SHORT presmooth_iter
number of presmoothing
Definition: fasp_block.h:338
void fasp_precond_dbsr_diag(REAL *r, REAL *z, void *data)
Diagonal preconditioner z=inv(D)*r.
Definition: precond_bsr.c:37
#define NL_AMLI_CYCLE
Definition: fasp_const.h:178
SHORT postsmooth_iter
number of postsmoothing
Definition: fasp_block.h:341
#define SOLVER_VGMRES
Definition: fasp_const.h:107
INT nb
dimension of each sub-block
Definition: fasp_block.h:296
dCSRmat * R_nk
Resriction for near kernal.
Definition: fasp_block.h:388
REAL relaxation
relaxation parameter for SOR smoother
Definition: fasp_block.h:347
INT * IA
integer array of row pointers, the size is ROW+1
Definition: fasp_block.h:70
#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
SHORT coarse_scaling
switch of scaling of the coarse grid correction
Definition: fasp.h:631
#define PRINT_MIN
Definition: fasp_const.h:80
SHORT itsolver_type
Definition: fasp.h:1107