Fast Auxiliary Space Preconditioning  1.8.4 Feb/15/2016
init.c
Go to the documentation of this file.
1 
8 #include "fasp.h"
9 #include "fasp_functs.h"
10 
11 /*---------------------------------*/
12 /*-- Public Functions --*/
13 /*---------------------------------*/
14 
26 {
27  pcdata->AMG_type = CLASSIC_AMG;
28  pcdata->print_level = PRINT_NONE;
29  pcdata->maxit = 500;
30  pcdata->max_levels = 20;
31  pcdata->tol = 1e-8;
32  pcdata->cycle_type = V_CYCLE;
33  pcdata->smoother = SMOOTHER_GS;
34  pcdata->smooth_order = CF_ORDER;
35  pcdata->presmooth_iter = 1;
36  pcdata->postsmooth_iter = 1;
37  pcdata->relaxation = 1.1;
38  pcdata->coarsening_type = 1;
39  pcdata->coarse_scaling = ON;
40  pcdata->amli_degree = 1;
42 }
43 
57 {
58  max_levels = MAX(1, max_levels); // at least allocate one level
59 
60  AMG_data *mgl = (AMG_data *)fasp_mem_calloc(max_levels,sizeof(AMG_data));
61 
62  INT i;
63  for ( i=0; i<max_levels; ++i ) {
64  mgl[i].max_levels = max_levels;
65  mgl[i].num_levels = 0;
66  mgl[i].near_kernel_dim = 0;
67  mgl[i].near_kernel_basis = NULL;
68  mgl[i].cycle_type = 0;
69  }
70 
71  return(mgl);
72 }
73 
87 {
88  max_levels = MAX(1, max_levels); // at least allocate one level
89 
90  AMG_data_bsr *mgl = (AMG_data_bsr *)fasp_mem_calloc(max_levels,sizeof(AMG_data_bsr));
91 
92  INT i;
93  for (i=0; i<max_levels; ++i) {
94  mgl[i].max_levels = max_levels;
95  mgl[i].num_levels = 0;
96  mgl[i].near_kernel_dim = 0;
97  mgl[i].near_kernel_basis = NULL;
98  mgl[i].A_nk = NULL;
99  mgl[i].P_nk = NULL;
100  mgl[i].R_nk = NULL;
101  }
102 
103  return(mgl);
104 }
105 
118 void fasp_ilu_data_alloc (const INT iwk,
119  const INT nwork,
120  ILU_data *iludata)
121 {
122 #if DEBUG_MODE > 0
123  printf("### DEBUG: %s ...... [Start]\n", __FUNCTION__);
124  printf("### DEBUG: iwk=%d, nwork=%d \n",iwk,nwork);
125 #endif
126  iludata->ijlu=(INT*)fasp_mem_calloc(iwk, sizeof(INT));
127 
128  iludata->luval=(REAL*)fasp_mem_calloc(iwk, sizeof(REAL));
129 
130  iludata->work=(REAL*)fasp_mem_calloc(nwork, sizeof(REAL));
131 #if DEBUG_MODE > 0
132  printf("### DEBUG: %s ...... %d [End]\n", __FUNCTION__,__LINE__);
133 #endif
134 
135  return;
136 }
137 
148 {
149  INT i;
150  fasp_dcsr_free(&Schwarz->A);
151 
152  for (i=0; i<Schwarz->nblk; ++i) fasp_dcsr_free (&((Schwarz->blk_data)[i]));
153 
154  Schwarz->nblk = 0;
155  fasp_mem_free (Schwarz->iblock);
156  fasp_mem_free (Schwarz->jblock);
157  fasp_dvec_free (&Schwarz->rhsloc1);
158  fasp_dvec_free (&Schwarz->xloc1);
159 
160  Schwarz->memt = 0;
161  fasp_mem_free (Schwarz->mask);
162  fasp_mem_free (Schwarz->maxa);
163 
164 #if WITH_MUMPS
165  if (Schwarz->mumps == NULL) return;
166 
167  for (i=0; i<Schwarz->nblk; ++i) fasp_mumps_free (&((Schwarz->mumps)[i]));
168 #endif
169 }
170 
186  AMG_param *param)
187 {
188  const INT max_levels = MAX(1,mgl[0].num_levels);
189 
190  INT i;
191 
192  for (i=0; i<max_levels; ++i) {
193  fasp_dcsr_free(&mgl[i].A);
194  fasp_dcsr_free(&mgl[i].P);
195  fasp_dcsr_free(&mgl[i].R);
196  fasp_dvec_free(&mgl[i].b);
197  fasp_dvec_free(&mgl[i].x);
198  fasp_dvec_free(&mgl[i].w);
199  fasp_ivec_free(&mgl[i].cfmark);
200  fasp_ilu_data_free(&mgl[i].LU);
201  fasp_Schwarz_data_free(&mgl[i].Schwarz);
202  }
203 
204  for (i=0; i<mgl->near_kernel_dim; ++i) {
205  if (mgl->near_kernel_basis[i]) fasp_mem_free(mgl->near_kernel_basis[i]);
206  mgl->near_kernel_basis[i] = NULL;
207  }
208 
209  // Clean direct solver data in necessary
210  switch (param->coarse_solver) {
211 
212 #if WITH_MUMPS
213  /* Destroy MUMPS direct solver on the coarsest level */
214  case SOLVER_MUMPS: {
215  mgl[max_levels-1].mumps.job = 3;
216  fasp_solver_mumps_steps(&mgl[max_levels-1].A, &mgl[max_levels-1].b,
217  &mgl[max_levels-1].x, &mgl[max_levels-1].mumps);
218  break;
219  }
220 #endif
221 #if WITH_UMFPACK
222  case SOLVER_UMFPACK: {
223  fasp_mem_free(mgl[max_levels-1].Numeric);
224  break;
225  }
226 #endif
227 
228 #if WITH_PARDISO
229  case SOLVER_PARDISO: {
230  fasp_pardiso_free_internal_mem(&mgl[max_levels-1].pdata);
231  }
232 #endif
233  default: // Do nothing!
234  break;
235  }
236 
238 
239  fasp_mem_free(mgl); mgl = NULL;
240 
241  if (param != NULL) {
242  if ( param->cycle_type == AMLI_CYCLE ) fasp_mem_free(param->amli_coef);
243  }
244 }
245 
257 {
258  const INT max_levels = mgl[0].max_levels;
259  INT i;
260 
261  for (i=0; i<max_levels; ++i) {
262  fasp_dbsr_free(&mgl[i].A);
263  fasp_dbsr_free(&mgl[i].P);
264  fasp_dbsr_free(&mgl[i].R);
265  fasp_dvec_free(&mgl[i].b);
266  fasp_dvec_free(&mgl[i].x);
267  fasp_dvec_free(&mgl[i].diaginv);
268  fasp_dvec_free(&mgl[i].diaginv_SS);
269  fasp_dcsr_free(&mgl[i].Ac);
270  fasp_dcsr_free(&mgl[i].PP);
271  fasp_mem_free(mgl[i].pw);
272  fasp_dbsr_free(&mgl[i].SS);
273  fasp_mem_free(mgl[i].sw);
274  fasp_dvec_free(&mgl[i].diaginv_SS);
275  fasp_ilu_data_free(&mgl[i].PP_LU);
276  fasp_dvec_free(&mgl[i].w);
277  fasp_ivec_free(&mgl[i].cfmark);
278  fasp_ilu_data_free(&mgl[i].LU);
279  }
280 
281  for (i=0; i<mgl->near_kernel_dim; ++i) {
282  if (mgl->near_kernel_basis[i]) fasp_mem_free(mgl->near_kernel_basis[i]);
283  mgl->near_kernel_basis[i] = NULL;
284  }
285 
287  fasp_mem_free(mgl); mgl = NULL;
288 }
289 
301 {
302  if (ILUdata==NULL) return;
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 
308  ILUdata->row = ILUdata->col = ILUdata->nzlu = ILUdata ->nwork = ILUdata->nb = 0;
309 }
310 
322 {
323  ILUdata->row = ILUdata->col = ILUdata->nzlu = 0;
324  ILUdata->ijlu = NULL; ILUdata->luval = NULL;
325 }
326 
338 {
339  pcdata->data = NULL;
340  pcdata->fct = NULL;
341 }
342 
343 /*---------------------------------*/
344 /*-- End of File --*/
345 /*---------------------------------*/
INT col
column of matrix LU, n
Definition: fasp.h:406
INT * ijlu
integer array of row pointers and column indexes, the size is nzlu
Definition: fasp.h:412
SHORT coarsening_type
switch of scaling of the coarse grid correction
Definition: fasp.h:834
AMG_data_bsr * fasp_amg_data_bsr_create(SHORT max_levels)
Create and initialize AMG_data data sturcture for AMG/SAMG (BSR format)
Definition: init.c:86
#define SMOOTHER_GS
Definition: fasp_const.h:184
SHORT coarse_scaling
switch of scaling of the coarse grid correction
Definition: fasp.h:840
Parameters for AMG solver.
Definition: fasp.h:583
void fasp_precond_data_null(precond_data *pcdata)
Initialize precond_data.
Definition: init.c:25
SHORT presmooth_iter
number of presmoothing
Definition: fasp.h:822
REAL relaxation
relaxation parameter for SOR smoother
Definition: fasp.h:828
void fasp_amg_data_bsr_free(AMG_data_bsr *mgl)
Free AMG_data_bsr data memeory space.
Definition: init.c:256
#define REAL
Definition: fasp.h:67
SHORT print_level
print level in AMG preconditioner
Definition: fasp.h:801
INT * iblock
row index of blocks
Definition: fasp.h:518
SHORT num_levels
number of levels in use <= max_levels
Definition: fasp.h:730
SHORT smoother
AMG smoother type.
Definition: fasp.h:816
INT maxit
max number of iterations of AMG preconditioner
Definition: fasp.h:804
Mumps_data * mumps
param for MUMPS
Definition: fasp.h:570
int fasp_solver_mumps_steps(dCSRmat *ptrA, dvector *b, dvector *u, Mumps_data *mumps)
Solve Ax=b by MUMPS in three steps.
#define SOLVER_GCG
Definition: fasp_const.h:109
INT nblk
number of blocks
Definition: fasp.h:515
void fasp_dvec_free(dvector *u)
Free vector data space of REAL type.
Definition: vec.c:139
INT * maxa
maxa
Definition: fasp.h:554
SHORT smooth_order
AMG smoother ordering.
Definition: fasp.h:819
void * fasp_mem_calloc(LONGLONG size, INT type)
1M = 1024*1024
Definition: memory.c:60
void(* fct)(REAL *, REAL *, void *)
action for preconditioner, void function pointer
Definition: fasp.h:1005
void fasp_Schwarz_data_free(Schwarz_data *Schwarz)
Free Schwarz_data data memeory space.
Definition: init.c:147
SHORT nl_amli_krylov_type
type of Krylov method used by Nonlinear AMLI cycle
Definition: fasp.h:846
void fasp_precond_null(precond *pcdata)
Initialize precond data.
Definition: init.c:337
INT near_kernel_dim
dimension of the near kernel for SAMG
Definition: fasp_block.h:261
SHORT AMG_type
type of AMG method
Definition: fasp.h:798
#define INT
Definition: fasp.h:64
Data for ILU setup.
Definition: fasp.h:400
Preconditioner data and action.
Definition: fasp.h:999
INT nwork
work space size
Definition: fasp.h:421
INT near_kernel_dim
dimension of the near kernel for SAMG
Definition: fasp.h:767
REAL ** near_kernel_basis
basis of near kernel space for SAMG
Definition: fasp.h:770
Data for multigrid levels. (BSR format)
Definition: fasp_block.h:198
void fasp_mem_free(void *mem)
Free up previous allocated memory body.
Definition: memory.c:150
#define V_CYCLE
Definition of cycle types.
Definition: fasp_const.h:175
void * data
data for preconditioner, void pointer
Definition: fasp.h:1002
void fasp_ivec_free(ivector *u)
Free vector data space of INT type.
Definition: vec.c:159
AMG_data * fasp_amg_data_create(SHORT max_levels)
Create and initialize AMG_data for classical and SA AMG.
Definition: init.c:56
REAL * luval
nonzero entries of LU
Definition: fasp.h:415
SHORT coarse_solver
coarse solver type
Definition: fasp.h:628
#define SOLVER_PARDISO
Definition: fasp_const.h:126
SHORT postsmooth_iter
number of postsmoothing
Definition: fasp.h:825
SHORT amli_degree
degree of the polynomial used by AMLI cycle
Definition: fasp.h:843
Main header file for FASP.
REAL * amli_coef
coefficients of the polynomial used by AMLI cycle
Definition: fasp.h:637
#define CLASSIC_AMG
Definition of AMG types.
Definition: fasp_const.h:162
#define SOLVER_MUMPS
Definition: fasp_const.h:125
INT row
row number of matrix LU, m
Definition: fasp.h:403
#define PRINT_NONE
Print level for all subroutines – not including DEBUG output.
Definition: fasp_const.h:79
dCSRmat * P_nk
Prolongation for near kernal.
Definition: fasp_block.h:273
INT * mask
mask
Definition: fasp.h:548
dvector rhsloc1
local right hand side
Definition: fasp.h:527
void fasp_ilu_data_alloc(const INT iwk, const INT nwork, ILU_data *iludata)
Allocate workspace for ILU factorization.
Definition: init.c:118
SHORT cycle_type
type of AMG cycle
Definition: fasp.h:604
SHORT max_levels
max number of AMG levels
Definition: fasp.h:807
dCSRmat A
pointer to the matrix
Definition: fasp.h:510
dCSRmat * R_nk
Resriction for near kernal.
Definition: fasp_block.h:276
INT cycle_type
cycle type
Definition: fasp.h:787
#define SHORT
FASP integer and floating point numbers.
Definition: fasp.h:63
INT max_levels
max number of levels
Definition: fasp_block.h:201
#define AMLI_CYCLE
Definition: fasp_const.h:177
Mumps_data mumps
data for MUMPS
Definition: fasp.h:784
REAL tol
tolerance for AMG preconditioner
Definition: fasp.h:810
INT memt
working space size
Definition: fasp.h:545
Data for Schwarz methods.
Definition: fasp.h:505
#define ON
Definition of switch.
Definition: fasp_const.h:73
INT num_levels
number of levels in use <= max_levels
Definition: fasp_block.h:204
Data for AMG solvers.
Definition: fasp.h:722
SHORT max_levels
max number of levels
Definition: fasp.h:727
dCSRmat * blk_data
matrix for each partition
Definition: fasp.h:557
#define SOLVER_UMFPACK
Definition: fasp_const.h:124
Data passed to the preconditioners.
Definition: fasp.h:795
void fasp_ilu_data_null(ILU_data *ILUdata)
Initialize ILU data.
Definition: init.c:321
void fasp_amg_data_free(AMG_data *mgl, AMG_param *param)
Free AMG_data data memeory space.
Definition: init.c:185
void fasp_ilu_data_free(ILU_data *ILUdata)
Create ILU_data sturcture.
Definition: init.c:300
INT nb
block size for BSR type only
Definition: fasp.h:418
SHORT cycle_type
AMG cycle type.
Definition: fasp.h:813
#define CF_ORDER
Definition: fasp_const.h:223
REAL * work
work space
Definition: fasp.h:424
REAL ** near_kernel_basis
basis of near kernel space for SAMG
Definition: fasp_block.h:264
dCSRmat * A_nk
Matrix data for near kernal.
Definition: fasp_block.h:270
INT nzlu
number of nonzero entries
Definition: fasp.h:409
INT job
work for MUMPS
Definition: fasp.h:467
void fasp_dbsr_free(dBSRmat *A)
Free memory space for BSR format sparse matrix.
Definition: sparse_bsr.c:133
void fasp_dcsr_free(dCSRmat *A)
Free CSR sparse matrix data memory space.
Definition: sparse_csr.c:166
#define MAX(a, b)
Definition of max, min, abs.
Definition: fasp.h:72
dvector xloc1
local solution
Definition: fasp.h:530
INT * jblock
column index of blocks
Definition: fasp.h:521