Fast Auxiliary Space Preconditioning 2.7.7 Aug/28/2022
SolBLC.c
Go to the documentation of this file.
1
17#include <math.h>
18#include <time.h>
19
20#include "fasp.h"
21#include "fasp_block.h"
22#include "fasp_functs.h"
23
24/*---------------------------------*/
25/*-- Declare Private Functions --*/
26/*---------------------------------*/
27
28#include "KryUtil.inl"
29
30/*---------------------------------*/
31/*-- Public Functions --*/
32/*---------------------------------*/
33
55 dvector *b,
56 dvector *x,
57 precond *pc,
58 ITS_param *itparam)
59{
60 const SHORT prtlvl = itparam->print_level;
61 const SHORT itsolver_type = itparam->itsolver_type;
62 const SHORT stop_type = itparam->stop_type;
63 const SHORT restart = itparam->restart;
64 const INT MaxIt = itparam->maxit;
65 const REAL tol = itparam->tol;
66
67 REAL solve_start, solve_end;
69
70#if DEBUG_MODE > 0
71 printf("### DEBUG: [-Begin-] %s ...\n", __FUNCTION__);
72 printf("### DEBUG: rhs/sol size: %d %d\n", b->row, x->row);
73#endif
74
75 fasp_gettime(&solve_start);
76
77 /* Safe-guard checks on parameters */
78 ITS_CHECK ( MaxIt, tol );
79
80 switch (itsolver_type) {
81
82 case SOLVER_BiCGstab:
83 iter=fasp_solver_dblc_pbcgs(A, b, x, pc, tol, MaxIt, stop_type, prtlvl);
84 break;
85
86 case SOLVER_MinRes:
87 iter=fasp_solver_dblc_pminres(A, b, x, pc, tol, MaxIt, stop_type, prtlvl);
88 break;
89
90 case SOLVER_GMRES:
91 iter=fasp_solver_dblc_pgmres(A, b, x, pc, tol, MaxIt, restart, stop_type, prtlvl);
92 break;
93
94 case SOLVER_VGMRES:
95 iter=fasp_solver_dblc_pvgmres(A, b, x, pc, tol, MaxIt, restart, stop_type, prtlvl);
96 break;
97
98 case SOLVER_VFGMRES:
99 iter=fasp_solver_dblc_pvfgmres(A, b, x, pc, tol, MaxIt, restart, stop_type, prtlvl);
100 break;
101
102 default:
103 printf("### ERROR: Unknown iterative solver type %d! [%s]\n",
104 itsolver_type, __FUNCTION__);
105 return ERROR_SOLVER_TYPE;
106
107 }
108
109 if ( (prtlvl >= PRINT_MIN) && (iter >= 0) ) {
110 fasp_gettime(&solve_end);
111 fasp_cputime("Iterative method", solve_end - solve_start);
112 }
113
114#if DEBUG_MODE > 0
115 printf("### DEBUG: [--End--] %s ...\n", __FUNCTION__);
116#endif
117
118 return iter;
119}
120
138 dvector *b,
139 dvector *x,
140 ITS_param *itparam)
141{
142 const SHORT prtlvl = itparam->print_level;
143
144 INT status = FASP_SUCCESS;
145 REAL solve_start, solve_end;
146
147#if DEBUG_MODE > 0
148 printf("### DEBUG: [-Begin-] %s ...\n", __FUNCTION__);
149#endif
150
151 // solver part
152 fasp_gettime(&solve_start);
153
154 status = fasp_solver_dblc_itsolver(A,b,x,NULL,itparam);
155
156 fasp_gettime(&solve_end);
157
158 if ( prtlvl >= PRINT_MIN )
159 fasp_cputime("Krylov method totally", solve_end - solve_start);
160
161#if DEBUG_MODE > 0
162 printf("### DEBUG: [--End--] %s ...\n", __FUNCTION__);
163#endif
164
165 return status;
166}
167
190 dvector *b,
191 dvector *x,
192 ITS_param *itparam,
193 AMG_param *amgparam,
194 dCSRmat *A_diag)
195{
196 const SHORT prtlvl = itparam->print_level;
197 const SHORT precond_type = itparam->precond_type;
198
199 INT status = FASP_SUCCESS;
200 REAL setup_start, setup_end;
201 REAL solve_start, solve_end;
202
203 const SHORT max_levels = amgparam->max_levels;
204 INT m, n, nnz, i;
205
206 AMG_data **mgl = NULL;
207
208#if WITH_UMFPACK
209 void **LU_diag = (void **)fasp_mem_calloc(3, sizeof(void *));
210#endif
211
212#if DEBUG_MODE > 0
213 printf("### DEBUG: [-Begin-] %s ...\n", __FUNCTION__);
214#endif
215
216 /* setup preconditioner */
217 fasp_gettime(&solve_start);
218 fasp_gettime(&setup_start);
219
220 /* diagonal blocks are solved exactly */
221 if ( precond_type > 20 && precond_type < 30 ) {
222#if WITH_UMFPACK
223 // Need to sort the diagonal blocks for UMFPACK format
224 dCSRmat A_tran;
225
226 for (i=0; i<3; i++){
227
228 A_tran = fasp_dcsr_create(A_diag[i].row, A_diag[i].col, A_diag[i].nnz);
229 fasp_dcsr_transz(&A_diag[i], NULL, &A_tran);
230 fasp_dcsr_cp(&A_tran, &A_diag[i]);
231
232 printf("Factorization for %d-th diagnol: \n", i);
233 LU_diag[i] = fasp_umfpack_factorize(&A_diag[i], prtlvl);
234
235 }
236
237 fasp_dcsr_free(&A_tran);
238#endif
239 }
240
241 /* diagonal blocks are solved by AMG */
242 else if ( precond_type > 30 && precond_type < 40 ) {
243
244 mgl = (AMG_data **)fasp_mem_calloc(3, sizeof(AMG_data *));
245
246 for (i=0; i<3; i++){
247
248 mgl[i] = fasp_amg_data_create(max_levels);
249 m = A_diag[i].row; n = A_diag[i].col; nnz = A_diag[i].nnz;
250 mgl[i][0].A=fasp_dcsr_create(m,n,nnz); fasp_dcsr_cp(&A_diag[i],&mgl[i][0].A);
251 mgl[i][0].b=fasp_dvec_create(n); mgl[i][0].x=fasp_dvec_create(n);
252
253 switch (amgparam->AMG_type) {
254 case SA_AMG: // Smoothed Aggregation AMG
255 status = fasp_amg_setup_sa(mgl[i], amgparam); break;
256 case UA_AMG: // Unsmoothed Aggregation AMG
257 status = fasp_amg_setup_ua(mgl[i], amgparam); break;
258 default: // Classical AMG
259 status = fasp_amg_setup_rs(mgl[i], amgparam); break;
260 }
261
262 fasp_chkerr(status, __FUNCTION__);
263 }
264
265 }
266
267 else {
268 fasp_chkerr(ERROR_SOLVER_PRECTYPE, __FUNCTION__);
269 }
270
271 precond_data_blc precdata;
272 precdata.Ablc = A;
273 precdata.A_diag = A_diag;
274 precdata.r = fasp_dvec_create(b->row);
275
276 /* diagonal blocks are solved exactly */
277 if ( precond_type > 20 && precond_type < 30 ) {
278#if WITH_UMFPACK
279 precdata.LU_diag = LU_diag;
280#endif
281 }
282 /* diagonal blocks are solved by AMG */
283 else if ( precond_type > 30 && precond_type < 40 ) {
284 precdata.amgparam = amgparam;
285 precdata.mgl = mgl;
286 }
287 else {
288 fasp_chkerr(ERROR_SOLVER_PRECTYPE, __FUNCTION__);
289 }
290
291 precond prec; prec.data = &precdata;
292
293 switch (precond_type) {
294 case 21:
295 prec.fct = fasp_precond_dblc_diag_3; break;
296
297 case 22:
298 prec.fct = fasp_precond_dblc_lower_3; break;
299
300 case 23:
301 prec.fct = fasp_precond_dblc_upper_3; break;
302
303 case 24:
304 prec.fct = fasp_precond_dblc_SGS_3; break;
305
306 case 31:
307 prec.fct = fasp_precond_dblc_diag_3_amg; break;
308
309 case 32:
311
312 case 33:
314
315 case 34:
316 prec.fct = fasp_precond_dblc_SGS_3_amg; break;
317
318 default:
319 fasp_chkerr(ERROR_SOLVER_PRECTYPE, __FUNCTION__); break;
320 }
321
322 if ( prtlvl >= PRINT_MIN ) {
323 fasp_gettime(&setup_end);
324 fasp_cputime("Setup totally", setup_end - setup_start);
325 }
326
327 // solve part
328 status = fasp_solver_dblc_itsolver(A,b,x, &prec,itparam);
329
330 fasp_gettime(&solve_end);
331
332 if ( prtlvl >= PRINT_MIN )
333 fasp_cputime("Krylov method totally", solve_end - solve_start);
334
335 // clean up
336 /* diagonal blocks are solved exactly */
337 if ( precond_type > 20 && precond_type < 30 ) {
338#if WITH_UMFPACK
339 for (i=0; i<3; i++) fasp_umfpack_free_numeric(LU_diag[i]);
340#endif
341 }
342 /* diagonal blocks are solved by AMG */
343 else if ( precond_type > 30 && precond_type < 40 ) {
344 for (i=0; i<3; i++) fasp_amg_data_free(mgl[i], amgparam);
345 fasp_mem_free(mgl); mgl = NULL;
346 }
347 else {
348 fasp_chkerr(ERROR_SOLVER_PRECTYPE, __FUNCTION__);
349 }
350
351#if DEBUG_MODE > 0
352 printf("### DEBUG: [--End--] %s ...\n", __FUNCTION__);
353#endif
354
355 return status;
356}
357
380 dvector *b,
381 dvector *x,
382 ITS_param *itparam,
383 AMG_param *amgparam,
384 dCSRmat *A_diag)
385{
386 const SHORT prtlvl = itparam->print_level;
387 const SHORT precond_type = itparam->precond_type;
388
389 INT status = FASP_SUCCESS;
390 REAL setup_start, setup_end;
391 REAL solve_start, solve_end;
392
393#if DEBUG_MODE > 0
394 printf("### DEBUG: [-Begin-] %s ...\n", __FUNCTION__);
395#endif
396
397 /* setup preconditioner */
398 fasp_gettime(&solve_start);
399 fasp_gettime(&setup_start);
400
401#if WITH_UMFPACK
402 void **LU_diag = (void **)fasp_mem_calloc(4, sizeof(void *));
403 INT i;
404#endif
405
406 /* diagonal blocks are solved exactly */
407 if ( precond_type > 20 && precond_type < 30 ) {
408
409#if WITH_UMFPACK
410 // Need to sort the matrices local_A for UMFPACK format
411 dCSRmat A_tran;
412
413 for (i=0; i<4; i++){
414
415 A_tran = fasp_dcsr_create(A_diag[i].row, A_diag[i].col, A_diag[i].nnz);
416 fasp_dcsr_transz(&A_diag[i], NULL, &A_tran);
417 fasp_dcsr_cp(&A_tran, &A_diag[i]);
418
419 printf("Factorization for %d-th diagnol: \n", i);
420 LU_diag[i] = fasp_umfpack_factorize(&A_diag[i], prtlvl);
421
422 }
423
424 fasp_dcsr_free(&A_tran);
425#endif
426
427 }
428 else {
429 fasp_chkerr(ERROR_SOLVER_PRECTYPE, __FUNCTION__);
430 }
431
432 precond_data_blc precdata;
433
434 precdata.Ablc = A;
435 precdata.A_diag = A_diag;
436#if WITH_UMFPACK
437 precdata.LU_diag = LU_diag;
438#endif
439 precdata.r = fasp_dvec_create(b->row);
440
441 precond prec; prec.data = &precdata;
442
443 switch (precond_type)
444 {
445 case 21:
447 break;
448
449 case 22:
451 break;
452 }
453
454 if ( prtlvl >= PRINT_MIN ) {
455 fasp_gettime(&setup_end);
456 fasp_cputime("Setup totally", setup_end - setup_start);
457 }
458
459 // solver part
460 status=fasp_solver_dblc_itsolver(A,b,x, &prec,itparam);
461
462 fasp_gettime(&solve_end);
463
464 if ( prtlvl >= PRINT_MIN )
465 fasp_cputime("Krylov method totally", solve_end - solve_start);
466
467 // clean
468#if WITH_UMFPACK
469 for (i=0; i<4; i++) fasp_umfpack_free_numeric(LU_diag[i]);
470#endif
471
472#if DEBUG_MODE > 0
473 printf("### DEBUG: [--End--] %s ...\n", __FUNCTION__);
474#endif
475
476 return status;
477}
478
502 dvector *b,
503 dvector *x,
504 ITS_param *itparam,
505 INT NumLayers,
506 dBLCmat *Ai,
507 dCSRmat *local_A,
508 ivector *local_index)
509{
510 const SHORT prtlvl = itparam->print_level;
511
512 INT status = FASP_SUCCESS;
513 REAL setup_start, setup_end;
514 REAL solve_start, solve_end;
515
516 void **local_LU = NULL;
517
518#if DEBUG_MODE > 0
519 printf("### DEBUG: [-Begin-] %s ...\n", __FUNCTION__);
520#endif
521
522 /* setup preconditioner */
523 fasp_gettime(&solve_start);
524 fasp_gettime(&setup_start);
525
526#if WITH_UMFPACK
527 // Need to sort the matrices local_A for UMFPACK format
528 INT l;
529 dCSRmat A_tran;
530 local_LU = (void **)fasp_mem_calloc(NumLayers, sizeof(void *));
531
532 for ( l=0; l<NumLayers; l++ ) {
533
534 A_tran = fasp_dcsr_create(local_A[l].row, local_A[l].col, local_A[l].nnz);
535 fasp_dcsr_transz(&local_A[l], NULL, &A_tran);
536 fasp_dcsr_cp(&A_tran, &local_A[l]);
537
538 printf("Factorization for layer %d: \n", l);
539 local_LU[l] = fasp_umfpack_factorize(&local_A[l], prtlvl);
540
541 }
542
543 fasp_dcsr_free(&A_tran);
544#endif
545
546 precond_data_sweeping precdata;
547 precdata.NumLayers = NumLayers;
548 precdata.A = A;
549 precdata.Ai = Ai;
550 precdata.local_A = local_A;
551 precdata.local_LU = local_LU;
552 precdata.local_index = local_index;
553 precdata.r = fasp_dvec_create(b->row);
554 precdata.w = (REAL *)fasp_mem_calloc(10*b->row,sizeof(REAL));
555
556 precond prec; prec.data = &precdata;
558
559 if ( prtlvl >= PRINT_MIN ) {
560 fasp_gettime(&setup_end);
561 fasp_cputime("Setup totally", setup_end - setup_start);
562 }
563
564 /* solver part */
565 status = fasp_solver_dblc_itsolver(A,b,x, &prec,itparam);
566
567 fasp_gettime(&solve_end);
568
569 if ( prtlvl >= PRINT_MIN )
570 fasp_cputime("Krylov method totally", solve_end - solve_start);
571
572 // clean
573#if WITH_UMFPACK
574 for (l=0; l<NumLayers; l++) fasp_umfpack_free_numeric(local_LU[l]);
575#endif
576
577#if DEBUG_MODE > 0
578 printf("### DEBUG: [--End--] %s ...\n", __FUNCTION__);
579#endif
580
581 return status;
582}
583
584/*---------------------------------*/
585/*-- End of File --*/
586/*---------------------------------*/
void fasp_mem_free(void *mem)
Free up previous allocated memory body and set pointer to NULL.
Definition: AuxMemory.c:155
void * fasp_mem_calloc(const unsigned int size, const unsigned int type)
Allocate, initiate, and check memory.
Definition: AuxMemory.c:65
void fasp_cputime(const char *message, const REAL cputime)
Print CPU walltime.
Definition: AuxMessage.c:179
void fasp_chkerr(const SHORT status, const char *fctname)
Check error status and print out error messages before quit.
Definition: AuxMessage.c:213
void fasp_gettime(REAL *time)
Get system time.
Definition: AuxTiming.c:36
dvector fasp_dvec_create(const INT m)
Create dvector data space of REAL type.
Definition: AuxVector.c:62
dCSRmat fasp_dcsr_create(const INT m, const INT n, const INT nnz)
Create CSR sparse matrix data memory space.
Definition: BlaSparseCSR.c:47
void fasp_dcsr_free(dCSRmat *A)
Free CSR sparse matrix data memory space.
Definition: BlaSparseCSR.c:184
void fasp_dcsr_transz(dCSRmat *A, INT *p, dCSRmat *AT)
Generalized transpose of A: (n x m) matrix given in dCSRmat format.
void fasp_dcsr_cp(const dCSRmat *A, dCSRmat *B)
copy a dCSRmat to a new one B=A
Definition: BlaSparseCSR.c:851
INT fasp_solver_dblc_pbcgs(dBLCmat *A, dvector *b, dvector *u, precond *pc, const REAL tol, const INT MaxIt, const SHORT StopType, const SHORT PrtLvl)
Preconditioned BiCGstab method for solving Au=b for BLC matrix.
Definition: KryPbcgs.c:713
INT fasp_solver_dblc_pgmres(dBLCmat *A, dvector *b, dvector *x, precond *pc, const REAL tol, const INT MaxIt, const SHORT restart, const SHORT StopType, const SHORT PrtLvl)
Preconditioned GMRES method for solving Au=b.
Definition: KryPgmres.c:675
INT fasp_solver_dblc_pminres(dBLCmat *A, dvector *b, dvector *u, precond *pc, const REAL tol, const INT MaxIt, const SHORT StopType, const SHORT PrtLvl)
A preconditioned minimal residual (Minres) method for solving Au=b.
Definition: KryPminres.c:475
INT fasp_solver_dblc_pvfgmres(dBLCmat *A, dvector *b, dvector *x, precond *pc, const REAL tol, const INT MaxIt, const SHORT restart, const SHORT StopType, const SHORT PrtLvl)
Solve "Ax=b" using PFGMRES (right preconditioned) iterative method in which the restart parameter can...
Definition: KryPvfgmres.c:714
INT fasp_solver_dblc_pvgmres(dBLCmat *A, dvector *b, dvector *x, precond *pc, const REAL tol, const INT MaxIt, const SHORT restart, const SHORT StopType, const SHORT PrtLvl)
Right preconditioned GMRES method in which the restart parameter can be adaptively modified during it...
Definition: KryPvgmres.c:757
SHORT fasp_amg_setup_rs(AMG_data *mgl, AMG_param *param)
Setup phase of Ruge and Stuben's classic AMG.
Definition: PreAMGSetupRS.c:52
SHORT fasp_amg_setup_sa(AMG_data *mgl, AMG_param *param)
Set up phase of smoothed aggregation AMG.
Definition: PreAMGSetupSA.c:63
SHORT fasp_amg_setup_ua(AMG_data *mgl, AMG_param *param)
Set up phase of unsmoothed aggregation AMG.
Definition: PreAMGSetupUA.c:55
void fasp_precond_dblc_upper_3(REAL *r, REAL *z, void *data)
block upper triangular preconditioning (3x3 blocks)
Definition: PreBLC.c:557
void fasp_precond_dblc_SGS_3(REAL *r, REAL *z, void *data)
Block symmetric GS preconditioning (3x3 blocks)
Definition: PreBLC.c:725
void fasp_precond_dblc_diag_4(REAL *r, REAL *z, void *data)
Block diagonal preconditioning (4x4 blocks)
Definition: PreBLC.c:191
void fasp_precond_dblc_lower_3(REAL *r, REAL *z, void *data)
block lower triangular preconditioning (3x3 blocks)
Definition: PreBLC.c:291
void fasp_precond_dblc_SGS_3_amg(REAL *r, REAL *z, void *data)
Block symmetric GS preconditioning (3x3 blocks)
Definition: PreBLC.c:838
void fasp_precond_dblc_diag_3_amg(REAL *r, REAL *z, void *data)
Block diagonal preconditioning (3x3 blocks)
Definition: PreBLC.c:126
void fasp_precond_dblc_lower_4(REAL *r, REAL *z, void *data)
block lower triangular preconditioning (4x4 blocks)
Definition: PreBLC.c:453
void fasp_precond_dblc_upper_3_amg(REAL *r, REAL *z, void *data)
block upper triangular preconditioning (3x3 blocks)
Definition: PreBLC.c:645
void fasp_precond_dblc_lower_3_amg(REAL *r, REAL *z, void *data)
block lower triangular preconditioning (3x3 blocks)
Definition: PreBLC.c:379
void fasp_precond_dblc_diag_3(REAL *r, REAL *z, void *data)
Block diagonal preconditioner (3x3 blocks)
Definition: PreBLC.c:38
void fasp_precond_dblc_sweeping(REAL *r, REAL *z, void *data)
Sweeping preconditioner for Maxwell equations.
Definition: PreBLC.c:939
AMG_data * fasp_amg_data_create(SHORT max_levels)
Create and initialize AMG_data for classical and SA AMG.
Definition: PreDataInit.c:64
void fasp_amg_data_free(AMG_data *mgl, AMG_param *param)
Free AMG_data data memeory space.
Definition: PreDataInit.c:101
INT fasp_solver_dblc_krylov(dBLCmat *A, dvector *b, dvector *x, ITS_param *itparam)
Solve Ax = b by standard Krylov methods.
Definition: SolBLC.c:137
INT fasp_solver_dblc_itsolver(dBLCmat *A, dvector *b, dvector *x, precond *pc, ITS_param *itparam)
Solve Ax = b by standard Krylov methods.
Definition: SolBLC.c:54
INT fasp_solver_dblc_krylov_block4(dBLCmat *A, dvector *b, dvector *x, ITS_param *itparam, AMG_param *amgparam, dCSRmat *A_diag)
Solve Ax = b by standard Krylov methods.
Definition: SolBLC.c:379
INT fasp_solver_dblc_krylov_block3(dBLCmat *A, dvector *b, dvector *x, ITS_param *itparam, AMG_param *amgparam, dCSRmat *A_diag)
Solve Ax = b by standard Krylov methods.
Definition: SolBLC.c:189
INT fasp_solver_dblc_krylov_sweeping(dBLCmat *A, dvector *b, dvector *x, ITS_param *itparam, INT NumLayers, dBLCmat *Ai, dCSRmat *local_A, ivector *local_index)
Solve Ax = b by standard Krylov methods.
Definition: SolBLC.c:501
Main header file for the FASP project.
#define REAL
Definition: fasp.h:67
#define SHORT
FASP integer and floating point numbers.
Definition: fasp.h:63
#define INT
Definition: fasp.h:64
Header file for FASP block matrices.
#define SOLVER_GMRES
Definition: fasp_const.h:106
#define FASP_SUCCESS
Definition of return status and error messages.
Definition: fasp_const.h:19
#define ERROR_SOLVER_TYPE
Definition: fasp_const.h:41
#define SOLVER_VGMRES
Definition: fasp_const.h:107
#define ERROR_SOLVER_PRECTYPE
Definition: fasp_const.h:42
#define SOLVER_BiCGstab
Definition: fasp_const.h:104
#define SA_AMG
Definition: fasp_const.h:163
#define SOLVER_MinRes
Definition: fasp_const.h:105
#define SOLVER_VFGMRES
Definition: fasp_const.h:108
#define PRINT_MIN
Definition: fasp_const.h:74
#define UA_AMG
Definition: fasp_const.h:164
Data for AMG methods.
Definition: fasp.h:790
dCSRmat A
pointer to the matrix at level level_num
Definition: fasp.h:803
dvector b
pointer to the right-hand side at level level_num
Definition: fasp.h:812
dvector x
pointer to the iterative solution at level level_num
Definition: fasp.h:815
Parameters for AMG methods.
Definition: fasp.h:447
SHORT AMG_type
type of AMG method
Definition: fasp.h:450
SHORT max_levels
max number of levels of AMG
Definition: fasp.h:462
Parameters for iterative solvers.
Definition: fasp.h:379
SHORT itsolver_type
Definition: fasp.h:382
SHORT print_level
Definition: fasp.h:381
SHORT precond_type
Definition: fasp.h:384
REAL tol
Definition: fasp.h:388
INT restart
Definition: fasp.h:386
SHORT stop_type
Definition: fasp.h:385
INT maxit
Definition: fasp.h:387
Block REAL CSR matrix format.
Definition: fasp_block.h:74
Sparse matrix of REAL type in CSR format.
Definition: fasp.h:143
INT col
column of matrix A, n
Definition: fasp.h:149
INT row
row number of matrix A, m
Definition: fasp.h:146
INT nnz
number of nonzero entries
Definition: fasp.h:152
Vector with n entries of REAL type.
Definition: fasp.h:346
INT row
number of rows
Definition: fasp.h:349
Vector with n entries of INT type.
Definition: fasp.h:361
Data for block preconditioners in dBLCmat format.
Definition: fasp_block.h:349
dBLCmat * Ablc
Definition: fasp_block.h:354
dCSRmat * A_diag
Definition: fasp_block.h:356
AMG_data ** mgl
Definition: fasp_block.h:368
AMG_param * amgparam
Definition: fasp_block.h:370
Data for sweeping preconditioner.
Definition: fasp_block.h:384
Preconditioner data and action.
Definition: fasp.h:1081
void * data
data for preconditioner, void pointer
Definition: fasp.h:1084
void(* fct)(REAL *, REAL *, void *)
action for preconditioner, void function pointer
Definition: fasp.h:1087