Fast Auxiliary Space Preconditioning  1.8.4 Feb/15/2016
wrapper.c
Go to the documentation of this file.
1 
6 #include "fasp.h"
7 #include "fasp_block.h"
8 #include "fasp_functs.h"
9 
10 /*---------------------------------*/
11 /*-- Public Functions --*/
12 /*---------------------------------*/
13 
36  INT *nnz,
37  INT *ia,
38  INT *ja,
39  REAL *a,
40  REAL *b,
41  REAL *u,
42  REAL *tol,
43  INT *maxit,
44  INT *ptrlvl)
45 {
46  dCSRmat mat; // coefficient matrix
47  dvector rhs, sol; // right-hand-side, solution
48  AMG_param amgparam; // parameters for AMG
49 
50  fasp_param_amg_init(&amgparam);
51  amgparam.tol = *tol;
52  amgparam.print_level = *ptrlvl;
53  amgparam.maxit = *maxit;
54 
55  mat.row = *n; mat.col = *n; mat.nnz = *nnz;
56  mat.IA = ia; mat.JA = ja; mat.val = a;
57 
58  rhs.row = *n; rhs.val = b;
59  sol.row = *n; sol.val = u;
60 
61  fasp_solver_amg(&mat, &rhs, &sol, &amgparam);
62 }
63 
86  INT *nnz,
87  INT *ia,
88  INT *ja,
89  REAL *a,
90  REAL *b,
91  REAL *u,
92  REAL *tol,
93  INT *maxit,
94  INT *ptrlvl)
95 {
96  dCSRmat mat; // coefficient matrix
97  dvector rhs, sol; // right-hand-side, solution
98  AMG_param amgparam; // parameters for AMG
99  itsolver_param itparam; // parameters for itsolver
100 
101  fasp_param_amg_init(&amgparam);
102  fasp_param_solver_init(&itparam);
103 
104  itparam.tol = *tol;
105  itparam.print_level = *ptrlvl;
106  itparam.maxit = *maxit;
107 
108  amgparam.print_level = *ptrlvl;
109 
110  mat.row = *n; mat.col = *n; mat.nnz = *nnz;
111  mat.IA = ia; mat.JA = ja; mat.val = a;
112 
113  rhs.row = *n; rhs.val = b;
114  sol.row = *n; sol.val = u;
115 
116  fasp_solver_dcsr_krylov_amg(&mat, &rhs, &sol, &itparam, &amgparam);
117 }
118 
144  INT nnz,
145  INT nb,
146  INT *ia,
147  INT *ja,
148  REAL *a,
149  REAL *b,
150  REAL *u,
151  REAL tol,
152  INT maxit,
153  INT ptrlvl)
154 {
155  dCSRmat mat; // coefficient matrix in CSR format
156  dBSRmat bsrmat; // coefficient matrix in BSR format
157  dvector rhs, sol; // right-hand-side, solution
158  AMG_param amgparam; // parameters for AMG
159  itsolver_param itparam; // parameters for itsolver
160  INT status = FASP_SUCCESS; // return parameter
161 
162  // setup AMG parameters
163  fasp_param_amg_init(&amgparam);
164  amgparam.AMG_type = UA_AMG;
165  amgparam.print_level = PRINT_NONE;
166 
167  amgparam.coarse_dof = 100;
168  amgparam.presmooth_iter = 1;
169  amgparam.postsmooth_iter = 1;
170 
171  amgparam.strong_coupled = 0.00;
172  amgparam.max_aggregation = 100;
173 
174  amgparam.ILU_type = ILUk;
175  amgparam.ILU_levels = 1;
176  amgparam.ILU_lfil = 1;
177 
178  // setup Krylov method parameters
179  fasp_param_solver_init(&itparam);
180  itparam.tol = tol;
181  itparam.print_level = ptrlvl;
182  itparam.maxit = maxit;
183 
184  itparam.itsolver_type = SOLVER_VFGMRES;
185  itparam.restart = 30;
186 
187  mat.row = n; mat.col = n; mat.nnz = nnz;
188  mat.IA = ia; mat.JA = ja; mat.val = a;
189 
190  // convert CSR to BSR format
191  bsrmat = fasp_format_dcsr_dbsr(&mat, nb);
192 
193  rhs.row = n; rhs.val = b;
194  sol.row = n; sol.val = u;
195 
196  // solve
197  status = fasp_solver_dbsr_krylov_amg(&bsrmat, &rhs, &sol, &itparam, &amgparam);
198 
199  // clean up
200  fasp_dbsr_free(&bsrmat);
201 
202  return status;
203 }
204 
230  INT nnz,
231  INT nb,
232  INT *ia,
233  INT *ja,
234  REAL *a,
235  REAL *b,
236  REAL *u,
237  REAL tol,
238  INT maxit,
239  INT ptrlvl)
240 {
241  dCOOmat coomat; // coefficient matrix in COO format
242  dCSRmat csrmat; // coefficient matrix in CSR format
243  dBSRmat bsrmat; // coefficient matrix in BSR format
244  dvector rhs, sol; // right-hand-side, solution
245  AMG_param amgparam; // parameters for AMG
246  itsolver_param itparam; // parameters for itsolver
247  INT status = FASP_SUCCESS; // return parameter
248 
249  // setup AMG parameters
250  fasp_param_amg_init(&amgparam);
251  amgparam.AMG_type = UA_AMG;
252  amgparam.print_level = ptrlvl;
253 
254  amgparam.coarse_dof = 100;
255  amgparam.presmooth_iter = 1;
256  amgparam.postsmooth_iter = 1;
257 
258  amgparam.strong_coupled = 0.00;
259  amgparam.max_aggregation = 100;
260 
261  amgparam.ILU_type = ILUk;
262  amgparam.ILU_levels = 1;
263  amgparam.ILU_lfil = 1;
264 
265  // setup Krylov method parameters
266  fasp_param_solver_init(&itparam);
267  itparam.tol = tol;
268  itparam.print_level = ptrlvl;
269  itparam.maxit = maxit;
270 
271  itparam.itsolver_type = SOLVER_VFGMRES;
272  itparam.restart = 30;
273 
274  // COO format
275  coomat.row = n; coomat.col = n; coomat.nnz = nnz;
276  coomat.rowind = ia; coomat.colind = ja; coomat.val = a;
277 
278  // convert COO to CSR format
279  fasp_format_dcoo_dcsr(&coomat,&csrmat);
280 
281  // convert CSR to BSR format
282  bsrmat = fasp_format_dcsr_dbsr(&csrmat, nb);
283 
284  // clean up
285  fasp_dcsr_free(&csrmat);
286 
287  rhs.row = n; rhs.val = b;
288  sol.row = n; sol.val = u;
289 
290  // solve
291  status = fasp_solver_dbsr_krylov_amg(&bsrmat, &rhs, &sol, &itparam, &amgparam);
292 
293  // clean up
294  fasp_dbsr_free(&bsrmat);
295 
296  return status;
297 }
298 
299 /*---------------------------------*/
300 /*-- End of File --*/
301 /*---------------------------------*/
INT max_aggregation
max size of each aggregate
Definition: fasp.h:673
INT * rowind
integer array of row indices, the size is nnz
Definition: fasp.h:221
SHORT AMG_type
type of AMG method
Definition: fasp.h:586
Parameters for AMG solver.
Definition: fasp.h:583
REAL strong_coupled
strong coupled threshold for aggregate
Definition: fasp.h:670
SHORT postsmooth_iter
number of postsmoothers
Definition: fasp.h:619
INT * JA
integer array of column indexes, the size is nnz
Definition: fasp.h:163
INT maxit
max number of iterations of AMG
Definition: fasp.h:592
INT coarse_dof
max number of coarsest level DOF
Definition: fasp.h:601
#define REAL
Definition: fasp.h:67
SHORT presmooth_iter
number of presmoothers
Definition: fasp.h:616
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
void fasp_fwrapper_krylov_amg_(INT *n, INT *nnz, INT *ia, INT *ja, REAL *a, REAL *b, REAL *u, REAL *tol, INT *maxit, INT *ptrlvl)
Solve Ax=b by Krylov method preconditioned by classic AMG.
Definition: wrapper.c:85
REAL * val
nonzero entries of A
Definition: fasp.h:166
REAL * val
nonzero entries of A
Definition: fasp.h:227
Block sparse row storage matrix of REAL type.
Definition: fasp_block.h:44
REAL * val
actual vector entries
Definition: fasp.h:348
INT col
column of matrix A, n
Definition: fasp.h:215
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
INT nnz
number of nonzero entries
Definition: fasp.h:218
#define INT
Definition: fasp.h:64
#define FASP_SUCCESS
Definition of return status and error messages.
Definition: fasp_const.h:27
INT restart
Definition: fasp.h:1112
INT row
row number of matrix A, m
Definition: fasp.h:212
Sparse matrix of REAL type in COO (or IJ) format.
Definition: fasp.h:209
INT fasp_wrapper_dcoo_dbsr_krylov_amg(INT n, INT nnz, INT nb, INT *ia, INT *ja, REAL *a, REAL *b, REAL *u, REAL tol, INT maxit, INT ptrlvl)
Solve Ax=b by Krylov method preconditioned by AMG (dcoo - > dbsr)
Definition: wrapper.c:229
#define ILUk
Type of ILU methods.
Definition: fasp_const.h:148
#define SOLVER_VFGMRES
Definition: fasp_const.h:108
INT nnz
number of nonzero entries
Definition: fasp.h:157
SHORT ILU_levels
number of levels use ILU smoother
Definition: fasp.h:682
INT col
column of matrix A, n
Definition: fasp.h:154
Main header file for FASP.
#define PRINT_NONE
Print level for all subroutines – not including DEBUG output.
Definition: fasp_const.h:79
SHORT fasp_format_dcoo_dcsr(dCOOmat *A, dCSRmat *B)
Transform a REAL matrix from its IJ format to its CSR format.
Definition: formats.c:27
INT row
row number of matrix A, m
Definition: fasp.h:151
Parameters passed to iterative solvers.
Definition: fasp.h:1105
Sparse matrix of REAL type in CSR format.
Definition: fasp.h:148
INT fasp_wrapper_dbsr_krylov_amg(INT n, INT nnz, INT nb, INT *ia, INT *ja, REAL *a, REAL *b, REAL *u, REAL tol, INT maxit, INT ptrlvl)
Solve Ax=b by Krylov method preconditioned by AMG (dcsr - > dbsr)
Definition: wrapper.c:143
void fasp_fwrapper_amg_(INT *n, INT *nnz, INT *ia, INT *ja, REAL *a, REAL *b, REAL *u, REAL *tol, INT *maxit, INT *ptrlvl)
Solve Ax=b by Ruge and Stuben's classic AMG.
Definition: wrapper.c:35
SHORT ILU_type
ILU type for smoothing.
Definition: fasp.h:685
INT row
number of rows
Definition: fasp.h:345
REAL tol
Definition: fasp.h:1111
void fasp_param_amg_init(AMG_param *amgparam)
Initialize AMG parameters.
Definition: parameters.c:390
SHORT print_level
Definition: fasp.h:1113
SHORT print_level
print level for AMG
Definition: fasp.h:589
INT * colind
integer array of column indices, the size is nnz
Definition: fasp.h:224
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
REAL tol
stopping tolerance for AMG solver
Definition: fasp.h:595
void fasp_param_solver_init(itsolver_param *itsparam)
Initialize itsolver_param.
Definition: parameters.c:455
Header file for FASP block matrices.
INT * IA
integer array of row pointers, the size is m+1
Definition: fasp.h:160
INT ILU_lfil
level of fill-in for ILUs and ILUk
Definition: fasp.h:688
void fasp_dbsr_free(dBSRmat *A)
Free memory space for BSR format sparse matrix.
Definition: sparse_bsr.c:133
dBSRmat fasp_format_dcsr_dbsr(dCSRmat *A, const INT nb)
Transfer a dCSRmat type matrix into a dBSRmat.
Definition: formats.c:721
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
SHORT itsolver_type
Definition: fasp.h:1107