Fast Auxiliary Space Preconditioning  1.8.4 Feb/15/2016
amg.c
Go to the documentation of this file.
1 
6 #include <time.h>
7 
8 #include "fasp.h"
9 #include "fasp_functs.h"
10 
11 /*---------------------------------*/
12 /*-- Public Functions --*/
13 /*---------------------------------*/
14 
38  dvector *b,
39  dvector *x,
40  AMG_param *param)
41 {
42  const SHORT max_levels = param->max_levels;
43  const SHORT prtlvl = param->print_level;
44  const SHORT amg_type = param->AMG_type;
45  const SHORT cycle_type = param->cycle_type;
46  const INT nnz = A->nnz, m = A->row, n = A->col;
47 
48  // local variables
49  SHORT status;
50  AMG_data * mgl = fasp_amg_data_create(max_levels);
51  REAL AMG_start, AMG_end;
52 
53 #if DEBUG_MODE > 0
54  printf("### DEBUG: %s ...... [Start]\n", __FUNCTION__);
55 #endif
56 
57  if ( prtlvl > PRINT_NONE ) fasp_gettime(&AMG_start);
58 
59  // check matrix data
60  if ( m != n ) {
61  printf("### ERROR: A is not a square matrix!\n");
62  fasp_chkerr(ERROR_MAT_SIZE, __FUNCTION__);
63  }
64 
65  if ( nnz <= 0 ) {
66  printf("### ERROR: A has no nonzero entries!\n");
67  fasp_chkerr(ERROR_MAT_SIZE, __FUNCTION__);
68  }
69 
70  // Step 0: initialize mgl[0] with A, b and x
71  mgl[0].A = fasp_dcsr_create(m, n, nnz);
72  fasp_dcsr_cp(A, &mgl[0].A);
73 
74  mgl[0].b = fasp_dvec_create(n);
75  fasp_dvec_cp(b, &mgl[0].b);
76 
77  mgl[0].x = fasp_dvec_create(n);
78  fasp_dvec_cp(x, &mgl[0].x);
79 
80  // Step 1: AMG setup phase
81  switch (amg_type) {
82 
83  case SA_AMG: // Smoothed Aggregation AMG setup
84  if ( prtlvl > PRINT_NONE ) printf("\nCalling SA AMG ...\n");
85  status = fasp_amg_setup_sa(mgl, param); break;
86 
87  case UA_AMG: // Unsmoothed Aggregation AMG setup
88  if ( prtlvl > PRINT_NONE ) printf("\nCalling UA AMG ...\n");
89  status = fasp_amg_setup_ua(mgl, param); break;
90 
91  default: // Classical AMG setup
92  if ( prtlvl > PRINT_NONE ) printf("\nCalling classical AMG ...\n");
93  status = fasp_amg_setup_rs(mgl, param);
94 
95  }
96 
97  // Step 2: AMG solve phase
98  if ( status == FASP_SUCCESS ) { // call a multilevel cycle
99 
100  switch (cycle_type) {
101 
102  case AMLI_CYCLE: // AMLI-cycle
103  fasp_amg_solve_amli(mgl, param); break;
104 
105  case NL_AMLI_CYCLE: // Nonlinear AMLI-cycle
106  fasp_amg_solve_nl_amli(mgl, param); break;
107 
108  default: // V,W-cycles (determined by param)
109  fasp_amg_solve(mgl, param); break;
110 
111  }
112 
113  fasp_dvec_cp(&mgl[0].x, x);
114 
115  }
116 
117  else { // call a backup solver
118 
119  if ( prtlvl > PRINT_MIN ) {
120  printf("### WARNING: AMG setup failed!\n");
121  printf("### WARNING: Use a backup solver instead.\n");
122  }
123  fasp_solver_dcsr_spgmres (A, b, x, NULL, param->tol, param->maxit,
124  20, 1, prtlvl);
125 
126  }
127 
128  // clean-up memory
129  fasp_amg_data_free(mgl, param);
130 
131  // print out CPU time if needed
132  if ( prtlvl > PRINT_NONE ) {
133  fasp_gettime(&AMG_end);
134  print_cputime("AMG totally", AMG_end - AMG_start);
135  }
136 
137 #if DEBUG_MODE > 0
138  printf("### DEBUG: %s ...... [Finish]\n", __FUNCTION__);
139 #endif
140 
141  return;
142 }
143 
144 /*---------------------------------*/
145 /*-- End of File --*/
146 /*---------------------------------*/
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
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
Parameters for AMG solver.
Definition: fasp.h:583
INT maxit
max number of iterations of AMG
Definition: fasp.h:592
#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
INT fasp_amg_solve_amli(AMG_data *mgl, AMG_param *param)
AMLI – SOLVE phase.
Definition: amg_solve.c:125
dvector b
pointer to the right-hand side at level level_num
Definition: fasp.h:744
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_dvec_cp(dvector *x, dvector *y)
Copy dvector x to dvector y.
Definition: vec.c:345
void fasp_solver_amg(dCSRmat *A, dvector *b, dvector *x, AMG_param *param)
Solve Ax = b by algebraic multigrid methods.
Definition: amg.c:37
Vector with n entries of REAL type.
Definition: fasp.h:342
#define INT
Definition: fasp.h:64
#define FASP_SUCCESS
Definition of return status and error messages.
Definition: fasp_const.h:27
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
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
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
void fasp_gettime(REAL *time)
Get system time.
Definition: timing.c:28
#define PRINT_NONE
Print level for all subroutines – not including DEBUG output.
Definition: fasp_const.h:79
INT fasp_amg_solve_nl_amli(AMG_data *mgl, AMG_param *param)
Nonlinear AMLI – SOLVE phase.
Definition: amg_solve.c:209
INT row
row number of matrix A, m
Definition: fasp.h:151
Sparse matrix of REAL type in CSR format.
Definition: fasp.h:148
SHORT cycle_type
type of AMG cycle
Definition: fasp.h:604
#define SHORT
FASP integer and floating point numbers.
Definition: fasp.h:63
#define AMLI_CYCLE
Definition: fasp_const.h:177
INT fasp_amg_solve(AMG_data *mgl, AMG_param *param)
AMG – SOLVE phase.
Definition: amg_solve.c:36
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
SHORT print_level
print level for AMG
Definition: fasp.h:589
Data for AMG solvers.
Definition: fasp.h:722
#define ERROR_MAT_SIZE
Definition: fasp_const.h:33
REAL tol
stopping tolerance for AMG solver
Definition: fasp.h:595
void fasp_amg_data_free(AMG_data *mgl, AMG_param *param)
Free AMG_data data memeory space.
Definition: init.c:185
#define NL_AMLI_CYCLE
Definition: fasp_const.h:178
#define SA_AMG
Definition: fasp_const.h:163
INT fasp_solver_dcsr_spgmres(dCSRmat *A, dvector *b, dvector *x, precond *pc, const REAL tol, const INT MaxIt, SHORT restart, const SHORT stop_type, const SHORT prtlvl)
Preconditioned GMRES method for solving Au=b with safe-guard.
Definition: spgmres.c:46
#define UA_AMG
Definition: fasp_const.h:164
#define PRINT_MIN
Definition: fasp_const.h:80