Fast Auxiliary Space Preconditioning 2.7.7 Aug/28/2022
PreDataInit.c
Go to the documentation of this file.
1
16#include "fasp.h"
17#include "fasp_functs.h"
18
19/*---------------------------------*/
20/*-- Public Functions --*/
21/*---------------------------------*/
22
34{
35 pcdata->AMG_type = CLASSIC_AMG;
36 pcdata->print_level = PRINT_NONE;
37 pcdata->maxit = 500;
38 pcdata->max_levels = 20;
39 pcdata->tol = 1e-8;
40 pcdata->cycle_type = V_CYCLE;
41 pcdata->smoother = SMOOTHER_GS;
42 pcdata->smooth_order = CF_ORDER;
43 pcdata->presmooth_iter = 1;
44 pcdata->postsmooth_iter = 1;
45 pcdata->relaxation = 1.1;
46 pcdata->coarsening_type = 1;
47 pcdata->coarse_scaling = ON;
48 pcdata->amli_degree = 1;
50}
51
65{
66 max_levels = MAX(1, max_levels); // at least allocate one level
67
68 AMG_data *mgl = (AMG_data *)fasp_mem_calloc(max_levels,sizeof(AMG_data));
69
70 INT i;
71 for ( i=0; i<max_levels; ++i ) {
72 mgl[i].max_levels = max_levels;
73 mgl[i].num_levels = 0;
74 mgl[i].near_kernel_dim = 0;
75 mgl[i].near_kernel_basis = NULL;
76 mgl[i].cycle_type = 0;
77#if MULTI_COLOR_ORDER
78 mgl[i].GS_Theta = 0.0E-2; //0.0; //1.0E-2;
79#endif
80 }
81
82 return(mgl);
83}
84
102 AMG_param *param)
103{
104 const INT max_levels = MAX(1,mgl[0].num_levels);
105
106 INT i;
107
108 switch (param->coarse_solver) {
109
110#if WITH_MUMPS
111 /* Destroy MUMPS direct solver on the coarsest level */
112 case SOLVER_MUMPS: {
113 mgl[max_levels-1].mumps.job = 3;
114 fasp_solver_mumps_steps(&mgl[max_levels-1].A, &mgl[max_levels-1].b,
115 &mgl[max_levels-1].x, &mgl[max_levels-1].mumps);
116 break;
117 }
118#endif
119
120#if WITH_UMFPACK
121 /* Destroy UMFPACK direct solver on the coarsest level */
122 case SOLVER_UMFPACK: {
123 fasp_mem_free(mgl[max_levels-1].Numeric); mgl[max_levels-1].Numeric = NULL;
124 break;
125 }
126#endif
127
128#if WITH_PARDISO
129 /* Destroy PARDISO direct solver on the coarsest level */
130 case SOLVER_PARDISO: {
131 fasp_pardiso_free_internal_mem(&mgl[max_levels-1].pdata);
132 break;
133 }
134
135#endif
136 default: // Do nothing!
137 break;
138 }
139
140 for ( i=0; i<max_levels; ++i ) {
141 fasp_ilu_data_free(&mgl[i].LU);
142 fasp_dcsr_free(&mgl[i].A);
143 if ( max_levels > 1 ) {
144 fasp_dcsr_free(&mgl[i].P);
145 fasp_dcsr_free(&mgl[i].R);
146 }
147 fasp_dvec_free(&mgl[i].b);
148 fasp_dvec_free(&mgl[i].x);
149 fasp_dvec_free(&mgl[i].w);
150 fasp_ivec_free(&mgl[i].cfmark);
151 fasp_swz_data_free(&mgl[i].Schwarz);
152 }
153
154 for ( i=0; i<mgl->near_kernel_dim; ++i ) {
155 fasp_mem_free(mgl->near_kernel_basis[i]); mgl->near_kernel_basis[i] = NULL;
156 }
157
159 fasp_mem_free(mgl); mgl = NULL;
160
161 if ( param == NULL ) return; // exit if no param given
162
163 if ( param->cycle_type == AMLI_CYCLE ) {
164 fasp_mem_free(param->amli_coef); param->amli_coef = NULL;
165 }
166
167}
168
182{
183 max_levels = MAX(1, max_levels); // at least allocate one level
184
185 AMG_data_bsr *mgl = (AMG_data_bsr *)fasp_mem_calloc(max_levels,sizeof(AMG_data_bsr));
186
187 INT i;
188 for (i=0; i<max_levels; ++i) {
189 mgl[i].max_levels = max_levels;
190 mgl[i].num_levels = 0;
191 mgl[i].near_kernel_dim = 0;
192 mgl[i].near_kernel_basis = NULL;
193 mgl[i].A_nk = NULL;
194 mgl[i].P_nk = NULL;
195 mgl[i].R_nk = NULL;
196 }
197
198 return(mgl);
199}
200
214{
215 const INT max_levels = MAX(1,mgl[0].num_levels);
216
217 INT i;
218
219 for ( i = 0; i < max_levels; ++i ) {
220
221 fasp_ilu_data_free(&mgl[i].LU);
222 fasp_dbsr_free(&mgl[i].A);
223 if ( max_levels > 1 ) {
224 fasp_dbsr_free(&mgl[i].P);
225 fasp_dbsr_free(&mgl[i].R);
226 }
227 fasp_dvec_free(&mgl[i].b);
228 fasp_dvec_free(&mgl[i].x);
229 fasp_dvec_free(&mgl[i].diaginv);
230 fasp_dvec_free(&mgl[i].diaginv_SS);
231 fasp_dcsr_free(&mgl[i].Ac);
232
233 fasp_ilu_data_free(&mgl[i].PP_LU);
234 fasp_dcsr_free(&mgl[i].PP);
235 fasp_dbsr_free(&mgl[i].SS);
236 fasp_dvec_free(&mgl[i].diaginv_SS);
237 fasp_dvec_free(&mgl[i].w);
238 fasp_ivec_free(&mgl[i].cfmark);
239
240 fasp_mem_free(mgl[i].pw); mgl[i].pw = NULL;
241 fasp_mem_free(mgl[i].sw); mgl[i].sw = NULL;
242 }
243
244 for ( i = 0; i < mgl->near_kernel_dim; ++i ) {
245 fasp_mem_free(mgl->near_kernel_basis[i]); mgl->near_kernel_basis[i] = NULL;
246 }
248 fasp_mem_free(mgl); mgl = NULL;
249}
250
265void fasp_ilu_data_create (const INT iwk,
266 const INT nwork,
267 ILU_data *iludata)
268{
269#if DEBUG_MODE > 0
270 printf("### DEBUG: [-Begin-] %s ...\n", __FUNCTION__);
271 printf("### DEBUG: iwk=%d, nwork=%d \n", iwk, nwork);
272#endif
273
274 iludata->ijlu=(INT*)fasp_mem_calloc(iwk, sizeof(INT));
275
276 if (iludata->type == ILUtp) iludata->iperm=(INT*)fasp_mem_calloc(iludata->row*2, sizeof(INT));
277
278 iludata->luval=(REAL*)fasp_mem_calloc(iwk, sizeof(REAL));
279
280 iludata->work=(REAL*)fasp_mem_calloc(nwork, sizeof(REAL));
281#if DEBUG_MODE > 0
282 printf("### DEBUG: %s ...... %d [End]\n", __FUNCTION__,__LINE__);
283#endif
284
285 return;
286}
287
301{
302 if ( iludata == NULL ) return; // There is nothing to do!
303
304 fasp_mem_free(iludata->ijlu); iludata->ijlu = NULL;
305 fasp_mem_free(iludata->luval); iludata->luval = NULL;
306 fasp_mem_free(iludata->work); iludata->work = NULL;
307 fasp_mem_free(iludata->ilevL); iludata->ilevL = NULL;
308 fasp_mem_free(iludata->jlevL); iludata->jlevL = NULL;
309 fasp_mem_free(iludata->ilevU); iludata->ilevU = NULL;
310 fasp_mem_free(iludata->jlevU); iludata->jlevU = NULL;
311
312 if ( iludata->type == ILUtp ) {
313
314 if ( iludata->A != NULL ) {
315 // To permute the matrix back to its original state use the loop:
316 INT k;
317 const INT nnz = iludata->A->nnz;
318 const INT *iperm = iludata->iperm;
319 for ( k = 0; k < nnz; k++ ) {
320 // iperm is in Fortran array format
321 iludata->A->JA[k] = iperm[ iludata->A->JA[k] ] -1;
322 }
323 }
324
325 fasp_mem_free(iludata->iperm); iludata->iperm = NULL;
326 }
327
328 iludata->row = iludata->col = iludata->nzlu = iludata->nwork = \
329 iludata->nb = iludata->nlevL = iludata->nlevU = 0;
330}
331
342{
343 INT i;
344
345 if ( swzdata == NULL ) return; // There is nothing to do!
346
347 fasp_dcsr_free(&swzdata->A);
348
349 for ( i=0; i<swzdata->nblk; ++i ) fasp_dcsr_free (&((swzdata->blk_data)[i]));
350
351 swzdata->nblk = 0;
352
353 fasp_mem_free (swzdata->iblock); swzdata->iblock = NULL;
354 fasp_mem_free (swzdata->jblock); swzdata->jblock = NULL;
355
356 fasp_dvec_free (&swzdata->rhsloc1);
357 fasp_dvec_free (&swzdata->xloc1);
358
359 swzdata->memt = 0;
360 fasp_mem_free (swzdata->mask); swzdata->mask = NULL;
361 fasp_mem_free (swzdata->maxa); swzdata->maxa = NULL;
362
363#if WITH_MUMPS
364 if ( swzdata->mumps == NULL ) return;
365
366 for ( i=0; i<swzdata->nblk; ++i ) fasp_mumps_free (&((swzdata->mumps)[i]));
367#endif
368}
369
370/*---------------------------------*/
371/*-- End of File --*/
372/*---------------------------------*/
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_dvec_free(dvector *u)
Free vector data space of REAL type.
Definition: AuxVector.c:145
void fasp_ivec_free(ivector *u)
Free vector data space of INT type.
Definition: AuxVector.c:164
void fasp_dbsr_free(dBSRmat *A)
Free memory space for BSR format sparse matrix.
Definition: BlaSparseBSR.c:146
void fasp_dcsr_free(dCSRmat *A)
Free CSR sparse matrix data memory space.
Definition: BlaSparseCSR.c:184
void fasp_precond_data_init(precond_data *pcdata)
Initialize precond_data.
Definition: PreDataInit.c:33
void fasp_swz_data_free(SWZ_data *swzdata)
Free SWZ_data data memeory space.
Definition: PreDataInit.c:341
void fasp_ilu_data_create(const INT iwk, const INT nwork, ILU_data *iludata)
Allocate workspace for ILU factorization.
Definition: PreDataInit.c:265
AMG_data * fasp_amg_data_create(SHORT max_levels)
Create and initialize AMG_data for classical and SA AMG.
Definition: PreDataInit.c:64
AMG_data_bsr * fasp_amg_data_bsr_create(SHORT max_levels)
Create and initialize AMG_data data sturcture for AMG/SAMG (BSR format)
Definition: PreDataInit.c:181
void fasp_amg_data_free(AMG_data *mgl, AMG_param *param)
Free AMG_data data memeory space.
Definition: PreDataInit.c:101
void fasp_ilu_data_free(ILU_data *iludata)
Create ILU_data sturcture.
Definition: PreDataInit.c:300
void fasp_amg_data_bsr_free(AMG_data_bsr *mgl)
Free AMG_data_bsr data memeory space.
Definition: PreDataInit.c:213
int fasp_solver_mumps_steps(dCSRmat *ptrA, dvector *b, dvector *u, Mumps_data *mumps)
Solve Ax=b by MUMPS in three steps.
Definition: XtrMumps.c:188
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 MAX(a, b)
Definition of max, min, abs.
Definition: fasp.h:74
#define INT
Definition: fasp.h:64
#define AMLI_CYCLE
Definition: fasp_const.h:179
#define SOLVER_PARDISO
Definition: fasp_const.h:126
#define SOLVER_MUMPS
Definition: fasp_const.h:125
#define ILUtp
Definition: fasp_const.h:150
#define V_CYCLE
Definition of cycle types.
Definition: fasp_const.h:177
#define SMOOTHER_GS
Definition: fasp_const.h:188
#define SOLVER_UMFPACK
Definition: fasp_const.h:124
#define CF_ORDER
Definition: fasp_const.h:234
#define CLASSIC_AMG
Definition of AMG types.
Definition: fasp_const.h:162
#define ON
Definition of switch.
Definition: fasp_const.h:67
#define PRINT_NONE
Print level for all subroutines – not including DEBUG output.
Definition: fasp_const.h:73
#define SOLVER_GCG
Definition: fasp_const.h:109
Data for multigrid levels in dBSRmat format.
Definition: fasp_block.h:146
INT near_kernel_dim
dimension of the near kernel for SAMG
Definition: fasp_block.h:209
REAL * pw
pointer to the auxiliary vectors for pressure block
Definition: fasp_block.h:185
dCSRmat * A_nk
Matrix data for near kernal.
Definition: fasp_block.h:218
REAL ** near_kernel_basis
basis of near kernel space for SAMG
Definition: fasp_block.h:212
INT num_levels
number of levels in use <= max_levels
Definition: fasp_block.h:152
dCSRmat * R_nk
Resriction for near kernal.
Definition: fasp_block.h:224
INT max_levels
max number of levels
Definition: fasp_block.h:149
dCSRmat * P_nk
Prolongation for near kernal.
Definition: fasp_block.h:221
REAL * sw
pointer to the auxiliary vectors for saturation block
Definition: fasp_block.h:191
Data for AMG methods.
Definition: fasp.h:790
INT near_kernel_dim
dimension of the near kernel for SAMG
Definition: fasp.h:835
INT cycle_type
cycle type
Definition: fasp.h:855
Mumps_data mumps
data for MUMPS
Definition: fasp.h:852
REAL ** near_kernel_basis
basis of near kernel space for SAMG
Definition: fasp.h:838
void * Numeric
pointer to the numerical factorization from UMFPACK
Definition: fasp.h:820
SHORT max_levels
max number of levels
Definition: fasp.h:795
SHORT num_levels
number of levels in use <= max_levels
Definition: fasp.h:798
Parameters for AMG methods.
Definition: fasp.h:447
REAL * amli_coef
coefficients of the polynomial used by AMLI cycle
Definition: fasp.h:501
SHORT coarse_solver
coarse solver type
Definition: fasp.h:492
SHORT cycle_type
type of AMG cycle
Definition: fasp.h:468
Data for ILU setup.
Definition: fasp.h:637
INT * ijlu
integer array of row pointers and column indexes, the size is nzlu
Definition: fasp.h:655
REAL * luval
nonzero entries of LU
Definition: fasp.h:658
INT col
column of matrix LU, n
Definition: fasp.h:649
INT * jlevU
mapping from row to color for upper triangle
Definition: fasp.h:702
INT nlevU
number of colors for upper triangle
Definition: fasp.h:690
INT nwork
work space size
Definition: fasp.h:664
INT row
row number of matrix LU, m
Definition: fasp.h:646
INT nlevL
number of colors for lower triangle
Definition: fasp.h:687
INT type
type of ILUdata
Definition: fasp.h:643
INT * iperm
permutation arrays for ILUtp
Definition: fasp.h:670
INT nzlu
number of nonzero entries
Definition: fasp.h:652
INT * jlevL
mapping from row to color for lower triangle
Definition: fasp.h:699
REAL * work
work space
Definition: fasp.h:667
INT * ilevL
number of vertices in each color for lower triangle
Definition: fasp.h:693
dCSRmat * A
pointer to the original coefficient matrix
Definition: fasp.h:640
INT * ilevU
number of vertices in each color for upper triangle
Definition: fasp.h:696
INT job
work for MUMPS
Definition: fasp.h:601
Data for Schwarz methods.
Definition: fasp.h:712
INT memt
working space size
Definition: fasp.h:752
INT * mask
mask
Definition: fasp.h:755
dCSRmat A
pointer to the original coefficient matrix
Definition: fasp.h:717
dvector rhsloc1
local right hand side
Definition: fasp.h:734
dvector xloc1
local solution
Definition: fasp.h:737
INT nblk
number of blocks
Definition: fasp.h:722
dCSRmat * blk_data
matrix for each partition
Definition: fasp.h:764
INT * jblock
column index of blocks
Definition: fasp.h:728
Mumps_data * mumps
param for MUMPS
Definition: fasp.h:777
INT * maxa
maxa
Definition: fasp.h:761
INT * iblock
row index of blocks
Definition: fasp.h:725
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
Data for preconditioners.
Definition: fasp.h:880
SHORT print_level
print level in AMG preconditioner
Definition: fasp.h:886
SHORT coarsening_type
switch of scaling of the coarse grid correction
Definition: fasp.h:919
SHORT nl_amli_krylov_type
type of Krylov method used by Nonlinear AMLI cycle
Definition: fasp.h:931
SHORT coarse_scaling
switch of scaling of the coarse grid correction
Definition: fasp.h:925
SHORT AMG_type
type of AMG method
Definition: fasp.h:883
REAL tol
tolerance for AMG preconditioner
Definition: fasp.h:895
REAL relaxation
relaxation parameter for SOR smoother
Definition: fasp.h:913
SHORT smoother
AMG smoother type.
Definition: fasp.h:901
SHORT cycle_type
AMG cycle type.
Definition: fasp.h:898
SHORT amli_degree
degree of the polynomial used by AMLI cycle
Definition: fasp.h:928
SHORT postsmooth_iter
number of postsmoothing
Definition: fasp.h:910
SHORT max_levels
max number of AMG levels
Definition: fasp.h:892
SHORT presmooth_iter
number of presmoothing
Definition: fasp.h:907
INT maxit
max number of iterations of AMG preconditioner
Definition: fasp.h:889
SHORT smooth_order
AMG smoother ordering.
Definition: fasp.h:904