Fast Auxiliary Space Preconditioning 2.7.7 Aug/28/2022
SolWrapper.c
Go to the documentation of this file.
1
19#include "fasp.h"
20#include "fasp_block.h"
21#include "fasp_functs.h"
22
23/*---------------------------------*/
24/*-- Public Functions --*/
25/*---------------------------------*/
26
45void fasp_fwrapper_dcsr_pardiso_(INT* n, INT* nnz, INT* ia, INT* ja, REAL* a, REAL* b,
46 REAL* u, INT* ptrlvl)
47{
48 dCSRmat mat; // coefficient matrix
49 dvector rhs, sol; // right-hand-side, solution
50
51 // set up coefficient matrix
52 mat.row = *n;
53 mat.col = *n;
54 mat.nnz = *nnz;
55 mat.IA = ia;
56 mat.JA = ja;
57 mat.val = a;
58
59 rhs.row = *n;
60 rhs.val = b;
61 sol.row = *n;
62 sol.val = u;
63
64 fasp_dcsr_sort(&mat);
65
66 fasp_solver_pardiso(&mat, &rhs, &sol, *ptrlvl);
67}
68
90void fasp_fwrapper_dcsr_amg_(INT* n, INT* nnz, INT* ia, INT* ja, REAL* a, REAL* b,
91 REAL* u, REAL* tol, INT* maxit, INT* ptrlvl)
92{
93 dCSRmat mat; // coefficient matrix
94 dvector rhs, sol; // right-hand-side, solution
95 AMG_param amgparam; // parameters for AMG
96
97 // setup AMG parameters
98 fasp_param_amg_init(&amgparam);
99
100 amgparam.tol = *tol;
101 amgparam.print_level = *ptrlvl;
102 amgparam.maxit = *maxit;
103
104 // set up coefficient matrix
105 mat.row = *n;
106 mat.col = *n;
107 mat.nnz = *nnz;
108 mat.IA = ia;
109 mat.JA = ja;
110 mat.val = a;
111
112 rhs.row = *n;
113 rhs.val = b;
114 sol.row = *n;
115 sol.val = u;
116
117 fasp_solver_amg(&mat, &rhs, &sol, &amgparam);
118}
119
141void fasp_fwrapper_dcsr_krylov_ilu_(INT* n, INT* nnz, INT* ia, INT* ja, REAL* a,
142 REAL* b, REAL* u, REAL* tol, INT* maxit,
143 INT* ptrlvl)
144{
145 dCSRmat mat; // coefficient matrix
146 dvector rhs, sol; // right-hand-side, solution
147 ILU_param iluparam; // parameters for ILU
148 ITS_param itsparam; // parameters for itsolver
149
150 // setup ILU parameters
151 fasp_param_ilu_init(&iluparam);
152
153 iluparam.print_level = *ptrlvl;
154
155 // setup Krylov method parameters
156 fasp_param_solver_init(&itsparam);
157
158 itsparam.itsolver_type = SOLVER_VFGMRES;
159 itsparam.tol = *tol;
160 itsparam.maxit = *maxit;
161 itsparam.print_level = *ptrlvl;
162
163 // set up coefficient matrix
164 mat.row = *n;
165 mat.col = *n;
166 mat.nnz = *nnz;
167 mat.IA = ia;
168 mat.JA = ja;
169 mat.val = a;
170
171 rhs.row = *n;
172 rhs.val = b;
173 sol.row = *n;
174 sol.val = u;
175
176 fasp_solver_dcsr_krylov_ilu(&mat, &rhs, &sol, &itsparam, &iluparam);
177}
178
200void fasp_fwrapper_dcsr_krylov_amg_(INT* n, INT* nnz, INT* ia, INT* ja, REAL* a,
201 REAL* b, REAL* u, REAL* tol, INT* maxit,
202 INT* ptrlvl)
203{
204 dCSRmat mat; // coefficient matrix
205 dvector rhs, sol; // right-hand-side, solution
206 input_param inparam; // parameters from input files
207 AMG_param amgparam; // parameters for AMG
208 ITS_param itsparam; // parameters for itsolver
209 ILU_param iluparam; // parameters for ILU
210
212 char* inputfile = "ini/amg.dat"; // Added for fasp4ns 2022.04.08 --zcs
213 fasp_param_input(inputfile, &inparam);
214 fasp_param_init(&inparam, &itsparam, &amgparam, &iluparam, NULL);
215
216 itsparam.tol = *tol;
217 itsparam.maxit = *maxit;
218 itsparam.print_level = *ptrlvl;
219
220 // set up coefficient matrix
221 mat.row = *n;
222 mat.col = *n;
223 mat.nnz = *nnz;
224 mat.IA = ia;
225 mat.JA = ja;
226 mat.val = a;
227
228 rhs.row = *n;
229 rhs.val = b;
230 sol.row = *n;
231 sol.val = u;
232
233 fasp_solver_dcsr_krylov_amg(&mat, &rhs, &sol, &itsparam, &amgparam);
234}
235
258void fasp_fwrapper_dbsr_krylov_ilu_(INT* n, INT* nnz, INT* nb, INT* ia, INT* ja,
259 REAL* a, REAL* b, REAL* u, REAL* tol, INT* maxit,
260 INT* ptrlvl)
261{
262 dBSRmat mat; // coefficient matrix in BSR format
263 dvector rhs, sol; // right-hand-side, solution
264
265 ILU_param iluparam; // parameters for ILU
266 ITS_param itsparam; // parameters for itsolver
267
268 // setup ILU parameters
269 fasp_param_ilu_init(&iluparam);
270 iluparam.ILU_lfil = 0;
271 iluparam.print_level = *ptrlvl;
272
273 // setup Krylov method parameters
274 fasp_param_solver_init(&itsparam);
275
276 itsparam.itsolver_type = SOLVER_VFGMRES;
277 itsparam.tol = *tol;
278 itsparam.maxit = *maxit;
279 itsparam.print_level = *ptrlvl;
280
281 // set up coefficient matrix
282 mat.ROW = *n;
283 mat.COL = *n;
284 mat.NNZ = *nnz;
285 mat.nb = *nb;
286 mat.IA = ia;
287 mat.JA = ja;
288 mat.val = a;
289
290 rhs.row = *n * *nb;
291 rhs.val = b;
292 sol.row = *n * *nb;
293 sol.val = u;
294
295 // solve
296 fasp_solver_dbsr_krylov_ilu(&mat, &rhs, &sol, &itsparam, &iluparam);
297}
298
321void fasp_fwrapper_dbsr_krylov_amg_(INT* n, INT* nnz, INT* nb, INT* ia, INT* ja,
322 REAL* a, REAL* b, REAL* u, REAL* tol, INT* maxit,
323 INT* ptrlvl)
324{
325 dBSRmat mat; // coefficient matrix in CSR format
326 dvector rhs, sol; // right-hand-side, solution
327
328 AMG_param amgparam; // parameters for AMG
329 ITS_param itsparam; // parameters for itsolver
330
331 // setup AMG parameters
332 fasp_param_amg_init(&amgparam);
333 amgparam.AMG_type = UA_AMG;
334 amgparam.print_level = *ptrlvl;
335
336 // setup Krylov method parameters
337 fasp_param_solver_init(&itsparam);
338 itsparam.tol = *tol;
339 itsparam.print_level = *ptrlvl;
340 itsparam.maxit = *maxit;
341 itsparam.itsolver_type = SOLVER_VFGMRES;
342
343 // set up coefficient matrix
344 mat.ROW = *n;
345 mat.COL = *n;
346 mat.NNZ = *nnz;
347 mat.nb = *nb;
348 mat.IA = ia;
349 mat.JA = ja;
350 mat.val = a;
351
352 rhs.row = *n * *nb;
353 rhs.val = b;
354 sol.row = *n * *nb;
355 sol.val = u;
356
357 // solve
358 fasp_solver_dbsr_krylov_amg(&mat, &rhs, &sol, &itsparam, &amgparam);
359}
360
361/*---------------------------------*/
362/*-- End of File --*/
363/*---------------------------------*/
void fasp_param_input(const char *fname, input_param *inparam)
Read input parameters from disk file.
Definition: AuxInput.c:112
void fasp_param_ilu_init(ILU_param *iluparam)
Initialize ILU parameters.
Definition: AuxParam.c:495
void fasp_param_init(const input_param *iniparam, ITS_param *itsparam, AMG_param *amgparam, ILU_param *iluparam, SWZ_param *swzparam)
Initialize parameters, global variables, etc.
Definition: AuxParam.c:283
void fasp_param_solver_init(ITS_param *itsparam)
Initialize ITS_param.
Definition: AuxParam.c:473
void fasp_param_amg_init(AMG_param *amgparam)
Initialize AMG parameters.
Definition: AuxParam.c:407
void fasp_dcsr_sort(dCSRmat *A)
Sort each row of A in ascending order w.r.t. column indices.
Definition: BlaSparseCSR.c:385
INT fasp_solver_amg(dCSRmat *A, dvector *b, dvector *x, AMG_param *param)
Solve Ax = b by algebraic multigrid methods.
Definition: SolAMG.c:49
INT fasp_solver_dbsr_krylov_ilu(dBSRmat *A, dvector *b, dvector *x, ITS_param *itparam, ILU_param *iluparam)
Solve Ax=b by ILUs preconditioned Krylov methods.
Definition: SolBSR.c:289
INT fasp_solver_dbsr_krylov_amg(dBSRmat *A, dvector *b, dvector *x, ITS_param *itparam, AMG_param *amgparam)
Solve Ax=b by AMG preconditioned Krylov methods.
Definition: SolBSR.c:354
INT fasp_solver_dcsr_krylov_ilu(dCSRmat *A, dvector *b, dvector *x, ITS_param *itparam, ILU_param *iluparam)
Solve Ax=b by ILUs preconditioned Krylov methods.
Definition: SolCSR.c:593
INT fasp_solver_dcsr_krylov_amg(dCSRmat *A, dvector *b, dvector *x, ITS_param *itparam, AMG_param *amgparam)
Solve Ax=b by AMG preconditioned Krylov methods.
Definition: SolCSR.c:483
void fasp_fwrapper_dcsr_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: SolWrapper.c:90
void fasp_fwrapper_dbsr_krylov_ilu_(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 block ILU in BSR format.
Definition: SolWrapper.c:258
void fasp_fwrapper_dcsr_krylov_ilu_(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 ILUk.
Definition: SolWrapper.c:141
void fasp_fwrapper_dcsr_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: SolWrapper.c:200
void fasp_fwrapper_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 block AMG in BSR format.
Definition: SolWrapper.c:321
void fasp_fwrapper_dcsr_pardiso_(INT *n, INT *nnz, INT *ia, INT *ja, REAL *a, REAL *b, REAL *u, INT *ptrlvl)
Solve Ax=b by the Pardiso direct solver.
Definition: SolWrapper.c:45
INT fasp_solver_pardiso(dCSRmat *ptrA, dvector *b, dvector *u, const SHORT prtlvl)
Solve Ax=b by PARDISO directly.
Definition: XtrPardiso.c:45
Main header file for the FASP project.
#define REAL
Definition: fasp.h:67
#define INT
Definition: fasp.h:64
Header file for FASP block matrices.
#define SOLVER_VFGMRES
Definition: fasp_const.h:108
#define UA_AMG
Definition: fasp_const.h:164
Parameters for AMG methods.
Definition: fasp.h:447
SHORT print_level
print level for AMG
Definition: fasp.h:453
SHORT AMG_type
type of AMG method
Definition: fasp.h:450
REAL tol
stopping tolerance for AMG solver
Definition: fasp.h:459
INT maxit
max number of iterations of AMG
Definition: fasp.h:456
Parameters for ILU.
Definition: fasp.h:396
INT ILU_lfil
level of fill-in for ILUk
Definition: fasp.h:405
SHORT print_level
print level
Definition: fasp.h:399
Parameters for iterative solvers.
Definition: fasp.h:379
SHORT itsolver_type
Definition: fasp.h:382
SHORT print_level
Definition: fasp.h:381
REAL tol
Definition: fasp.h:388
INT maxit
Definition: fasp.h:387
Block sparse row storage matrix of REAL type.
Definition: fasp_block.h:34
INT COL
number of cols of sub-blocks in matrix A, N
Definition: fasp_block.h:40
INT NNZ
number of nonzero sub-blocks in matrix A, NNZ
Definition: fasp_block.h:43
REAL * val
Definition: fasp_block.h:57
INT nb
dimension of each sub-block
Definition: fasp_block.h:46
INT * IA
integer array of row pointers, the size is ROW+1
Definition: fasp_block.h:60
INT ROW
number of rows of sub-blocks in matrix A, M
Definition: fasp_block.h:37
INT * JA
Definition: fasp_block.h:64
Sparse matrix of REAL type in CSR format.
Definition: fasp.h:143
INT col
column of matrix A, n
Definition: fasp.h:149
REAL * val
nonzero entries of A
Definition: fasp.h:161
INT row
row number of matrix A, m
Definition: fasp.h:146
INT * IA
integer array of row pointers, the size is m+1
Definition: fasp.h:155
INT nnz
number of nonzero entries
Definition: fasp.h:152
INT * JA
integer array of column indexes, the size is nnz
Definition: fasp.h:158
Vector with n entries of REAL type.
Definition: fasp.h:346
REAL * val
actual vector entries
Definition: fasp.h:352
INT row
number of rows
Definition: fasp.h:349
Input parameters.
Definition: fasp.h:1111