Fast Auxiliary Space Preconditioning  1.8.4 Feb/15/2016
parameters.c
Go to the documentation of this file.
1 
6 #include <stdio.h>
7 
8 #include "fasp.h"
9 #include "fasp_functs.h"
10 
11 /*---------------------------------*/
12 /*-- Public Functions --*/
13 /*---------------------------------*/
14 
27 void fasp_param_set (int argc,
28  const char *argv[],
29  input_param *iniparam)
30 {
31  int arg_index = 1;
32  int print_usage = 0;
33  SHORT status = FASP_SUCCESS;
34 
35  // Option 1. set default input parameters
36  fasp_param_input_init(iniparam);
37 
38  while ( arg_index < argc ) {
39 
40  if ( strcmp(argv[arg_index], "-help") == 0 ) {
41  print_usage = 1; break;
42  }
43 
44  // Option 2. Get parameters from an ini file
45  else if ( strcmp(argv[arg_index], "-ini") == 0 ) {
46  arg_index++;
47  if ( arg_index >= argc ) {
48  printf("### ERROR: Missing ini file name!\n");
49  print_usage = 1; break;
50  }
51  strcpy(iniparam->inifile, argv[arg_index]);
52  fasp_param_input(iniparam->inifile,iniparam);
53  if ( ++arg_index >= argc ) break;
54  }
55 
56  // Option 3. Get parameters from command line input
57  else if ( strcmp(argv[arg_index], "-print") == 0 ) {
58  arg_index++;
59  if ( arg_index >= argc ) {
60  printf("### ERROR: Expecting print level (int between 0 and 10).\n");
61  print_usage = 1; break;
62  }
63  iniparam->print_level = atoi(argv[arg_index]);
64  if ( ++arg_index >= argc ) break;
65  }
66 
67  else if ( strcmp(argv[arg_index], "-output") == 0 ) {
68  arg_index++;
69  if ( arg_index >= argc ) {
70  printf("### ERROR: Expecting output type (0 or 1).\n");
71  print_usage = 1; break;
72  }
73  iniparam->output_type = atoi(argv[arg_index]);
74  if ( ++arg_index >= argc ) break;
75  }
76 
77  else if ( strcmp(argv[arg_index], "-solver") == 0 ) {
78  arg_index++;
79  if ( arg_index >= argc ) {
80  printf("### ERROR: Expecting solver type.\n");
81  print_usage = 1; break;
82  }
83  iniparam->solver_type = atoi(argv[arg_index]);
84  if ( ++arg_index >= argc ) break;
85  }
86 
87  else if ( strcmp(argv[arg_index], "-precond") == 0 ) {
88  arg_index++;
89  if ( arg_index >= argc ) {
90  printf("### ERROR: Expecting preconditioner type.\n");
91  print_usage = 1; break;
92  }
93  iniparam->precond_type = atoi(argv[arg_index]);
94  if ( ++arg_index >= argc ) break;
95  }
96 
97  else if ( strcmp(argv[arg_index], "-maxit") == 0 ) {
98  arg_index++;
99  if ( arg_index >= argc ) {
100  printf("### ERROR: Expecting max number of iterations.\n");
101  print_usage = 1; break;
102  }
103  iniparam->itsolver_maxit = atoi(argv[arg_index]);
104  if ( ++arg_index >= argc ) break;
105  }
106 
107  else if ( strcmp(argv[arg_index], "-tol") == 0 ) {
108  arg_index++;
109  if ( arg_index >= argc ) {
110  printf("### ERROR: Expecting tolerance for itsolver.\n");
111  print_usage = 1; break;
112  }
113  iniparam->itsolver_tol = atof(argv[arg_index]);
114  if ( ++arg_index >= argc ) break;
115  }
116 
117  else if ( strcmp(argv[arg_index], "-amgmaxit") == 0 ) {
118  arg_index++;
119  if ( arg_index >= argc ) {
120  printf("### ERROR: Expecting max number of iterations for AMG.\n");
121  print_usage = 1; break;
122  }
123  iniparam->AMG_maxit = atoi(argv[arg_index]);
124  if ( ++arg_index >= argc ) break;
125  }
126 
127  else if ( strcmp(argv[arg_index], "-amgtol") == 0 ) {
128  arg_index++;
129  if ( arg_index >= argc ) {
130  printf("### ERROR: Expecting tolerance for AMG.\n");
131  print_usage = 1; break;
132  }
133  iniparam->AMG_tol = atof(argv[arg_index]);
134  if ( ++arg_index >= argc ) break;
135  }
136 
137  else if ( strcmp(argv[arg_index], "-amgtype") == 0 ) {
138  arg_index++;
139  if ( arg_index >= argc ) {
140  printf("### ERROR: Expecting AMG type (1, 2, 3).\n");
141  print_usage = 1; break;
142  }
143  iniparam->AMG_type = atoi(argv[arg_index]);
144  if ( ++arg_index >= argc ) break;
145  }
146 
147  else if ( strcmp(argv[arg_index], "-amgcycle") == 0 ) {
148  arg_index++;
149  if ( arg_index >= argc ) {
150  printf("### ERROR: Expecting AMG cycle type (1, 2, 3).\n");
151  print_usage = 1; break;
152  }
153  iniparam->AMG_cycle_type = atoi(argv[arg_index]);
154  if ( ++arg_index >= argc ) break;
155  }
156 
157  else if ( strcmp(argv[arg_index], "-amgcoarsening") == 0 ) {
158  arg_index++;
159  if ( arg_index >= argc ) {
160  printf("### ERROR: Expecting AMG coarsening type.\n");
161  print_usage = 1; break;
162  }
163  iniparam->AMG_coarsening_type = atoi(argv[arg_index]);
164  if ( ++arg_index >= argc ) break;
165  }
166 
167  else if ( strcmp(argv[arg_index], "-amginterplation") == 0 ) {
168  arg_index++;
169  if ( arg_index >= argc ) {
170  printf("### ERROR: Expecting AMG interpolation type.\n");
171  print_usage = 1; break;
172  }
173  iniparam->AMG_interpolation_type = atoi(argv[arg_index]);
174  if ( ++arg_index >= argc ) break;
175  }
176 
177  else if ( strcmp(argv[arg_index], "-amgsmoother") == 0 ) {
178  arg_index++;
179  if ( arg_index >= argc ) {
180  printf("### ERROR: Expecting AMG smoother type.\n");
181  print_usage = 1; break;
182  }
183  iniparam->AMG_smoother = atoi(argv[arg_index]);
184  if ( ++arg_index >= argc ) break;
185  }
186 
187  else if ( strcmp(argv[arg_index], "-amgsthreshold") == 0 ) {
188  arg_index++;
189  if ( arg_index >= argc ) {
190  printf("### ERROR: Expecting AMG strong threshold.\n");
191  print_usage = 1; break;
192  }
193  iniparam->AMG_strong_threshold = atof(argv[arg_index]);
194  if ( ++arg_index >= argc ) break;
195  }
196 
197  else if ( strcmp(argv[arg_index], "-amgscouple") == 0 ) {
198  arg_index++;
199  if ( arg_index >= argc ) {
200  printf("### ERROR: Expecting AMG strong coupled threshold.\n");
201  print_usage = 1; break;
202  }
203  iniparam->AMG_strong_coupled = atof(argv[arg_index]);
204  if ( ++arg_index >= argc ) break;
205  }
206 
207  else {
208  print_usage = 1;
209  break;
210  }
211 
212  }
213 
214  if ( print_usage ) {
215 
216  printf("FASP command line options:\n");
217  printf("================================================================\n");
218  printf(" -ini [CharValue] : Ini file name\n");
219  printf(" -print [IntValue] : Print level\n");
220  printf(" -output [IntValue] : Output to screen or a log file\n");
221  printf(" -solver [IntValue] : Solver type\n");
222  printf(" -precond [IntValue] : Preconditioner type\n");
223  printf(" -maxit [IntValue] : Max number of iterations\n");
224  printf(" -tol [RealValue] : Tolerance for iterative solvers\n");
225  printf(" -amgmaxit [IntValue] : Max number of AMG iterations\n");
226  printf(" -amgtol [RealValue] : Tolerance for AMG methods\n");
227  printf(" -amgtype [IntValue] : AMG type\n");
228  printf(" -amgcycle [IntValue] : AMG cycle type\n");
229  printf(" -amgcoarsening [IntValue] : AMG coarsening type\n");
230  printf(" -amginterpolation [IntValue] : AMG interpolation type\n");
231  printf(" -amgsmoother [IntValue] : AMG smoother type\n");
232  printf(" -amgsthreshold [RealValue] : AMG strong threshold\n");
233  printf(" -amgscoupled [RealValue] : AMG strong coupled threshold\n");
234  printf(" -help : Brief help messages\n");
235 
236  exit(ERROR_INPUT_PAR);
237 
238  }
239 
240  // sanity checks
241  status = fasp_param_check(iniparam);
242 
243  // if meet unexpected input, stop the program
244  fasp_chkerr(status,__FUNCTION__);
245 
246 }
247 
270 void fasp_param_init (input_param *iniparam,
271  itsolver_param *itsparam,
272  AMG_param *amgparam,
273  ILU_param *iluparam,
274  Schwarz_param *schparam)
275 {
276 #if CHMEM_MODE
277  total_alloc_mem = 0; // initialize total memeory amount
278  total_alloc_count = 0; // initialize alloc count
279 #endif
280 
281  if (itsparam) fasp_param_solver_init(itsparam);
282 
283  if (amgparam) fasp_param_amg_init(amgparam);
284 
285  if (iluparam) fasp_param_ilu_init(iluparam);
286 
287  if (schparam) fasp_param_Schwarz_init(schparam);
288 
289  if (iniparam) {
290  if (itsparam) fasp_param_solver_set(itsparam,iniparam);
291  if (amgparam) fasp_param_amg_set(amgparam,iniparam);
292  if (iluparam) fasp_param_ilu_set(iluparam,iniparam);
293  if (schparam) fasp_param_Schwarz_set(schparam,iniparam);
294  }
295  else {
296  printf("### WARNING: No input specified. Use default values instead!\n");
297  }
298 }
299 
311 {
312  strcpy(iniparam->workdir,"../data/");
313 
314  // Input/output
315  iniparam->print_level = PRINT_SOME;
316  iniparam->output_type = 0;
317 
318  // Problem information
319  iniparam->problem_num = 10;
320  iniparam->solver_type = SOLVER_CG;
321  iniparam->precond_type = PREC_AMG;
322  iniparam->stop_type = STOP_REL_RES;
323 
324  // Solver parameters
325  iniparam->itsolver_tol = 1e-6;
326  iniparam->itsolver_maxit = 500;
327  iniparam->restart = 25;
328 
329  // ILU method parameters
330  iniparam->ILU_type = ILUk;
331  iniparam->ILU_lfil = 0;
332  iniparam->ILU_droptol = 0.001;
333  iniparam->ILU_relax = 0;
334  iniparam->ILU_permtol = 0.0;
335 
336  // Schwarz method parameters
337  iniparam->Schwarz_mmsize = 200;
338  iniparam->Schwarz_maxlvl = 2;
339  iniparam->Schwarz_type = 1;
340  iniparam->Schwarz_blksolver = SOLVER_DEFAULT;
341 
342  // AMG method parameters
343  iniparam->AMG_type = CLASSIC_AMG;
344  iniparam->AMG_levels = 20;
345  iniparam->AMG_cycle_type = V_CYCLE;
346  iniparam->AMG_smoother = SMOOTHER_GS;
347  iniparam->AMG_smooth_order = CF_ORDER;
348  iniparam->AMG_presmooth_iter = 1;
349  iniparam->AMG_postsmooth_iter = 1;
350  iniparam->AMG_relaxation = 1.0;
351  iniparam->AMG_coarse_dof = 500;
352  iniparam->AMG_coarse_solver = 0;
353  iniparam->AMG_tol = 1e-6;
354  iniparam->AMG_maxit = 1;
355  iniparam->AMG_ILU_levels = 0;
356  iniparam->AMG_Schwarz_levels = 0;
357  iniparam->AMG_coarse_scaling = OFF; // Require investigation --Chensong
358  iniparam->AMG_amli_degree = 1;
359  iniparam->AMG_nl_amli_krylov_type = 2;
360 
361  // Classical AMG specific
362  iniparam->AMG_coarsening_type = 1;
363  iniparam->AMG_interpolation_type = 1;
364  iniparam->AMG_max_row_sum = 0.9;
365  iniparam->AMG_strong_threshold = 0.3;
366  iniparam->AMG_truncation_threshold = 0.2;
367  iniparam->AMG_aggressive_level = 0;
368  iniparam->AMG_aggressive_path = 1;
369 
370  // Aggregation AMG specific
371  iniparam->AMG_aggregation_type = PAIRWISE;
372  iniparam->AMG_quality_bound = 8.0;
373  iniparam->AMG_pair_number = 2;
374  iniparam->AMG_strong_coupled = 0.25;
375  iniparam->AMG_max_aggregation = 9;
376  iniparam->AMG_tentative_smooth = 0.67;
377  iniparam->AMG_smooth_filter = ON;
378 }
379 
391 {
392  // General AMG parameters
393  amgparam->AMG_type = CLASSIC_AMG;
394  amgparam->print_level = PRINT_NONE;
395  amgparam->maxit = 1;
396  amgparam->tol = 1e-6;
397  amgparam->max_levels = 20;
398  amgparam->coarse_dof = 500;
399  amgparam->cycle_type = V_CYCLE;
400  amgparam->smoother = SMOOTHER_GS;
401  amgparam->smooth_order = CF_ORDER;
402  amgparam->presmooth_iter = 1;
403  amgparam->postsmooth_iter = 1;
404  amgparam->coarse_solver = SOLVER_DEFAULT;
405  amgparam->relaxation = 1.0;
406  amgparam->polynomial_degree = 3;
407  amgparam->coarse_scaling = OFF;
408  amgparam->amli_degree = 2;
409  amgparam->amli_coef = NULL;
410  amgparam->nl_amli_krylov_type = SOLVER_GCG;
411 
412  // Classical AMG specific
413  amgparam->coarsening_type = COARSE_RS;
414  amgparam->interpolation_type = INTERP_DIR;
415  amgparam->max_row_sum = 0.9;
416  amgparam->strong_threshold = 0.3;
417  amgparam->truncation_threshold = 0.2;
418  amgparam->aggressive_level = 0;
419  amgparam->aggressive_path = 1;
420 
421  // Aggregation AMG specific
422  amgparam->aggregation_type = PAIRWISE;
423  amgparam->quality_bound = 10.0;
424  amgparam->pair_number = 2;
425  amgparam->strong_coupled = 0.25;
426  amgparam->max_aggregation = 20;
427  amgparam->tentative_smooth = 0.67;
428  amgparam->smooth_filter = ON;
429 
430  // ILU smoother parameters
431  amgparam->ILU_type = ILUk;
432  amgparam->ILU_levels = 0;
433  amgparam->ILU_lfil = 0;
434  amgparam->ILU_droptol = 0.001;
435  amgparam->ILU_relax = 0;
436 
437  // Schwarz smoother parameters
438  amgparam->Schwarz_levels = 0; // how many levels will use Schwarz smoother
439  amgparam->Schwarz_mmsize = 200;
440  amgparam->Schwarz_maxlvl = 3; // blocksize -- vertices with smaller distance
441  amgparam->Schwarz_type = 1;
442  amgparam->Schwarz_blksolver = SOLVER_DEFAULT;
443 }
444 
456 {
457  itsparam->print_level = PRINT_NONE;
458  itsparam->itsolver_type = SOLVER_CG;
459  itsparam->precond_type = PREC_AMG;
460  itsparam->stop_type = STOP_REL_RES;
461  itsparam->maxit = 500;
462  itsparam->restart = 25;
463  itsparam->tol = 1e-6;
464 }
465 
477 {
478  iluparam->print_level = PRINT_NONE;
479  iluparam->ILU_type = ILUk;
480  iluparam->ILU_lfil = 2;
481  iluparam->ILU_droptol = 0.001;
482  iluparam->ILU_relax = 0;
483  iluparam->ILU_permtol = 0.01;
484 }
485 
499 {
500  schparam->print_level = PRINT_NONE;
501  schparam->Schwarz_type = 3;
502  schparam->Schwarz_maxlvl = 2;
503  schparam->Schwarz_mmsize = 200;
504  schparam->Schwarz_blksolver = 0;
505 }
506 
519  input_param *iniparam)
520 {
521  param->AMG_type = iniparam->AMG_type;
522  param->print_level = iniparam->print_level;
523 
524  if (iniparam->solver_type == SOLVER_AMG) {
525  param->maxit = iniparam->itsolver_maxit;
526  param->tol = iniparam->itsolver_tol;
527  }
528  else if (iniparam->solver_type == SOLVER_FMG) {
529  param->maxit = iniparam->itsolver_maxit;
530  param->tol = iniparam->itsolver_tol;
531  }
532  else {
533  param->maxit = iniparam->AMG_maxit;
534  param->tol = iniparam->AMG_tol;
535  }
536 
537  param->max_levels = iniparam->AMG_levels;
538  param->cycle_type = iniparam->AMG_cycle_type;
539  param->smoother = iniparam->AMG_smoother;
540  param->smooth_order = iniparam->AMG_smooth_order;
541  param->relaxation = iniparam->AMG_relaxation;
542  param->coarse_solver = iniparam->AMG_coarse_solver;
543  param->polynomial_degree = iniparam->AMG_polynomial_degree;
544  param->presmooth_iter = iniparam->AMG_presmooth_iter;
545  param->postsmooth_iter = iniparam->AMG_postsmooth_iter;
546  param->coarse_dof = iniparam->AMG_coarse_dof;
547  param->coarse_scaling = iniparam->AMG_coarse_scaling;
548  param->amli_degree = iniparam->AMG_amli_degree;
549  param->amli_coef = NULL;
550  param->nl_amli_krylov_type = iniparam->AMG_nl_amli_krylov_type;
551 
552  param->coarsening_type = iniparam->AMG_coarsening_type;
553  param->interpolation_type = iniparam->AMG_interpolation_type;
554  param->strong_threshold = iniparam->AMG_strong_threshold;
556  param->max_row_sum = iniparam->AMG_max_row_sum;
557  param->aggressive_level = iniparam->AMG_aggressive_level;
558  param->aggressive_path = iniparam->AMG_aggressive_path;
559 
560  param->aggregation_type = iniparam->AMG_aggregation_type;
561  param->pair_number = iniparam->AMG_pair_number;
562  param->quality_bound = iniparam->AMG_quality_bound;
563  param->strong_coupled = iniparam->AMG_strong_coupled;
564  param->max_aggregation = iniparam->AMG_max_aggregation;
565  param->tentative_smooth = iniparam->AMG_tentative_smooth;
566  param->smooth_filter = iniparam->AMG_smooth_filter;
567  param->smooth_filter = iniparam->AMG_smooth_filter;
568 
569  param->ILU_levels = iniparam->AMG_ILU_levels;
570  param->ILU_type = iniparam->ILU_type;
571  param->ILU_lfil = iniparam->ILU_lfil;
572  param->ILU_droptol = iniparam->ILU_droptol;
573  param->ILU_relax = iniparam->ILU_relax;
574  param->ILU_permtol = iniparam->ILU_permtol;
575 
576  param->Schwarz_levels = iniparam->AMG_Schwarz_levels;
577  param->Schwarz_mmsize = iniparam->Schwarz_mmsize;
578  param->Schwarz_maxlvl = iniparam->Schwarz_maxlvl;
579  param->Schwarz_type = iniparam->Schwarz_type;
580 }
581 
594  input_param *iniparam)
595 {
596  iluparam->print_level = iniparam->print_level;
597  iluparam->ILU_type = iniparam->ILU_type;
598  iluparam->ILU_lfil = iniparam->ILU_lfil;
599  iluparam->ILU_droptol = iniparam->ILU_droptol;
600  iluparam->ILU_relax = iniparam->ILU_relax;
601  iluparam->ILU_permtol = iniparam->ILU_permtol;
602 }
603 
616  input_param *iniparam)
617 {
618  schparam->print_level = iniparam->print_level;
619  schparam->Schwarz_type = iniparam->Schwarz_type;
620  schparam->Schwarz_maxlvl = iniparam->Schwarz_maxlvl;
621  schparam->Schwarz_mmsize = iniparam->Schwarz_mmsize;
622  schparam->Schwarz_blksolver = iniparam->Schwarz_blksolver;
623 }
624 
637  input_param *iniparam)
638 {
639  itsparam->print_level = iniparam->print_level;
640  itsparam->itsolver_type = iniparam->solver_type;
641  itsparam->precond_type = iniparam->precond_type;
642  itsparam->stop_type = iniparam->stop_type;
643  itsparam->restart = iniparam->restart;
644 
645  if ( itsparam->itsolver_type == SOLVER_AMG ) {
646  itsparam->tol = iniparam->AMG_tol;
647  itsparam->maxit = iniparam->AMG_maxit;
648  }
649  else {
650  itsparam->tol = iniparam->itsolver_tol;
651  itsparam->maxit = iniparam->itsolver_maxit;
652  }
653 }
654 
667  AMG_param *amgparam)
668 {
669  pcdata->AMG_type = amgparam->AMG_type;
670  pcdata->print_level = amgparam->print_level;
671  pcdata->maxit = amgparam->maxit;
672  pcdata->max_levels = amgparam->max_levels;
673  pcdata->tol = amgparam->tol;
674  pcdata->cycle_type = amgparam->cycle_type;
675  pcdata->smoother = amgparam->smoother;
676  pcdata->smooth_order = amgparam->smooth_order;
677  pcdata->presmooth_iter = amgparam->presmooth_iter;
678  pcdata->postsmooth_iter = amgparam->postsmooth_iter;
679  pcdata->coarsening_type = amgparam->coarsening_type;
680  pcdata->coarse_solver = amgparam->coarse_solver;
681  pcdata->relaxation = amgparam->relaxation;
682  pcdata->polynomial_degree = amgparam->polynomial_degree;
683  pcdata->coarse_scaling = amgparam->coarse_scaling;
684  pcdata->amli_degree = amgparam->amli_degree;
685  pcdata->amli_coef = amgparam->amli_coef;
686  pcdata->nl_amli_krylov_type = amgparam->nl_amli_krylov_type;
687  pcdata->tentative_smooth = amgparam->tentative_smooth;
688 }
689 
702  precond_data *pcdata)
703 {
704  amgparam->AMG_type = pcdata->AMG_type;
705  amgparam->print_level = pcdata->print_level;
706  amgparam->cycle_type = pcdata->cycle_type;
707  amgparam->smoother = pcdata->smoother;
708  amgparam->smooth_order = pcdata->smooth_order;
709  amgparam->presmooth_iter = pcdata->presmooth_iter;
710  amgparam->postsmooth_iter = pcdata->postsmooth_iter;
711  amgparam->relaxation = pcdata->relaxation;
712  amgparam->polynomial_degree = pcdata->polynomial_degree;
713  amgparam->coarse_solver = pcdata->coarse_solver;
714  amgparam->coarse_scaling = pcdata->coarse_scaling;
715  amgparam->amli_degree = pcdata->amli_degree;
716  amgparam->amli_coef = pcdata->amli_coef;
717  amgparam->nl_amli_krylov_type = pcdata->nl_amli_krylov_type;
718  amgparam->tentative_smooth = pcdata->tentative_smooth;
719  amgparam->ILU_levels = pcdata->mgl_data->ILU_levels;
720 }
721 
734  AMG_param *amgparam)
735 {
736  pcdata->AMG_type = amgparam->AMG_type;
737  pcdata->print_level = amgparam->print_level;
738  pcdata->maxit = amgparam->maxit;
739  pcdata->max_levels = amgparam->max_levels;
740  pcdata->tol = amgparam->tol;
741  pcdata->cycle_type = amgparam->cycle_type;
742  pcdata->smoother = amgparam->smoother;
743  pcdata->smooth_order = amgparam->smooth_order;
744  pcdata->presmooth_iter = amgparam->presmooth_iter;
745  pcdata->postsmooth_iter = amgparam->postsmooth_iter;
746  pcdata->coarse_solver = amgparam->coarse_solver;
747  pcdata->coarsening_type = amgparam->coarsening_type;
748  pcdata->relaxation = amgparam->relaxation;
749  pcdata->coarse_scaling = amgparam->coarse_scaling;
750  pcdata->amli_degree = amgparam->amli_degree;
751  pcdata->amli_coef = amgparam->amli_coef;
752  pcdata->nl_amli_krylov_type = amgparam->nl_amli_krylov_type;
753  pcdata->tentative_smooth = amgparam->tentative_smooth;
754 }
755 
768  precond_data_bsr *pcdata)
769 {
770  amgparam->AMG_type = pcdata->AMG_type;
771  amgparam->print_level = pcdata->print_level;
772  amgparam->cycle_type = pcdata->cycle_type;
773  amgparam->smoother = pcdata->smoother;
774  amgparam->smooth_order = pcdata->smooth_order;
775  amgparam->presmooth_iter = pcdata->presmooth_iter;
776  amgparam->postsmooth_iter = pcdata->postsmooth_iter;
777  amgparam->relaxation = pcdata->relaxation;
778  amgparam->coarse_solver = pcdata->coarse_solver;
779  amgparam->coarse_scaling = pcdata->coarse_scaling;
780  amgparam->amli_degree = pcdata->amli_degree;
781  amgparam->amli_coef = pcdata->amli_coef;
782  amgparam->nl_amli_krylov_type = pcdata->nl_amli_krylov_type;
783  amgparam->tentative_smooth = pcdata->tentative_smooth;
784  amgparam->ILU_levels = pcdata->mgl_data->ILU_levels;
785 }
786 
798 {
799 
800  if ( param ) {
801 
802  printf("\n Parameters in AMG_param\n");
803  printf("-----------------------------------------------\n");
804 
805  printf("AMG print level: %d\n", param->print_level);
806  printf("AMG max num of iter: %d\n", param->maxit);
807  printf("AMG type: %d\n", param->AMG_type);
808  printf("AMG tolerance: %.2e\n", param->tol);
809  printf("AMG max levels: %d\n", param->max_levels);
810  printf("AMG cycle type: %d\n", param->cycle_type);
811  printf("AMG coarse solver type: %d\n", param->coarse_solver);
812  printf("AMG scaling of coarse correction: %d\n", param->coarse_scaling);
813  printf("AMG smoother type: %d\n", param->smoother);
814  printf("AMG smoother order: %d\n", param->smooth_order);
815  printf("AMG num of presmoothing: %d\n", param->presmooth_iter);
816  printf("AMG num of postsmoothing: %d\n", param->postsmooth_iter);
817 
818  if ( param->smoother == SMOOTHER_SOR ||
819  param->smoother == SMOOTHER_SSOR ||
820  param->smoother == SMOOTHER_GSOR ||
821  param->smoother == SMOOTHER_SGSOR ) {
822  printf("AMG relax factor: %.4f\n", param->relaxation);
823  }
824 
825  if ( param->smoother == SMOOTHER_POLY ) {
826  printf("AMG polynomial smoother degree: %d\n", param->polynomial_degree);
827  }
828 
829  if ( param->cycle_type == AMLI_CYCLE ) {
830  printf("AMG AMLI degree of polynomial: %d\n", param->amli_degree);
831  }
832 
833  if ( param->cycle_type == NL_AMLI_CYCLE ) {
834  printf("AMG Nonlinear AMLI Krylov type: %d\n", param->nl_amli_krylov_type);
835  }
836 
837  switch (param->AMG_type) {
838  case CLASSIC_AMG:
839  printf("AMG coarsening type: %d\n", param->coarsening_type);
840  printf("AMG interpolation type: %d\n", param->interpolation_type);
841  printf("AMG dof on coarsest grid: %d\n", param->coarse_dof);
842  printf("AMG strong threshold: %.4f\n", param->strong_threshold);
843  printf("AMG truncation threshold: %.4f\n", param->truncation_threshold);
844  printf("AMG max row sum: %.4f\n", param->max_row_sum);
845  printf("AMG aggressive levels: %d\n", param->aggressive_level);
846  printf("AMG aggressive path: %d\n", param->aggressive_path);
847  break;
848 
849  default: // SA_AMG or UA_AMG
850  printf("Aggregation type: %d\n", param->aggregation_type);
851  if ( param->aggregation_type == PAIRWISE ) {
852  printf("Aggregation number of pairs: %d\n", param->pair_number);
853  printf("Aggregation quality bound: %.2f\n", param->quality_bound);
854  }
855  if ( param->aggregation_type == VMB ) {
856  printf("Aggregation AMG strong coupling: %.4f\n", param->strong_coupled);
857  printf("Aggregation AMG max aggregation: %d\n", param->max_aggregation);
858  printf("Aggregation AMG tentative smooth: %.4f\n", param->tentative_smooth);
859  printf("Aggregation AMG smooth filter: %d\n", param->smooth_filter);
860  }
861  break;
862  }
863 
864  if (param->ILU_levels>0) {
865  printf("AMG ILU smoother level: %d\n", param->ILU_levels);
866  printf("AMG ILU type: %d\n", param->ILU_type);
867  printf("AMG ILU level of fill-in: %d\n", param->ILU_lfil);
868  printf("AMG ILU drop tol: %e\n", param->ILU_droptol);
869  printf("AMG ILU relaxation: %f\n", param->ILU_relax);
870  }
871 
872  if (param->Schwarz_levels>0){
873  printf("AMG Schwarz smoother level: %d\n", param->Schwarz_levels);
874  printf("AMG Schwarz type: %d\n", param->Schwarz_type);
875  printf("AMG Schwarz forming block level: %d\n", param->Schwarz_maxlvl);
876  printf("AMG Schwarz maximal block size: %d\n", param->Schwarz_mmsize);
877  }
878 
879  printf("-----------------------------------------------\n\n");
880 
881  }
882  else {
883  printf("### WARNING: param has not been set!\n");
884  } // end if (param)
885 
886 }
887 
899 {
900  if ( param ) {
901 
902  printf("\n Parameters in ILU_param\n");
903  printf("-----------------------------------------------\n");
904  printf("ILU print level: %d\n", param->print_level);
905  printf("ILU type: %d\n", param->ILU_type);
906  printf("ILU level of fill-in: %d\n", param->ILU_lfil);
907  printf("ILU relaxation factor: %.4f\n", param->ILU_relax);
908  printf("ILU drop tolerance: %.2e\n", param->ILU_droptol);
909  printf("ILU permutation tolerance: %.2e\n", param->ILU_permtol);
910  printf("-----------------------------------------------\n\n");
911 
912  }
913  else {
914  printf("### WARNING: param has not been set!\n");
915  }
916 }
917 
929 {
930  if ( param ) {
931 
932  printf("\n Parameters in Schwarz_param\n");
933  printf("-----------------------------------------------\n");
934  printf("Schwarz print level: %d\n", param->print_level);
935  printf("Schwarz type: %d\n", param->Schwarz_type);
936  printf("Schwarz forming block level: %d\n", param->Schwarz_maxlvl);
937  printf("Schwarz maximal block size: %d\n", param->Schwarz_mmsize);
938  printf("Schwarz block solver type: %d\n", param->Schwarz_blksolver);
939  printf("-----------------------------------------------\n\n");
940 
941  }
942  else {
943  printf("### WARNING: param has not been set!\n");
944  }
945 }
946 
958 {
959  if ( param ) {
960 
961  printf("\n Parameters in itsolver_param\n");
962  printf("-----------------------------------------------\n");
963 
964  printf("Solver print level: %d\n", param->print_level);
965  printf("Solver type: %d\n", param->itsolver_type);
966  printf("Solver precond type: %d\n", param->precond_type);
967  printf("Solver max num of iter: %d\n", param->maxit);
968  printf("Solver tolerance: %.2e\n", param->tol);
969  printf("Solver stopping type: %d\n", param->stop_type);
970 
971  if (param->itsolver_type==SOLVER_GMRES ||
972  param->itsolver_type==SOLVER_VGMRES) {
973  printf("Solver restart number: %d\n", param->restart);
974  }
975 
976  printf("-----------------------------------------------\n\n");
977 
978  }
979  else {
980  printf("### WARNING: param has not been set!\n");
981  }
982 }
983 
984 /*---------------------------------*/
985 /*-- End of File --*/
986 /*---------------------------------*/
INT max_aggregation
max size of each aggregate
Definition: fasp.h:673
void fasp_param_input_init(input_param *iniparam)
Initialize input parameters.
Definition: parameters.c:310
SHORT coarsening_type
switch of scaling of the coarse grid correction
Definition: fasp.h:834
#define SMOOTHER_GS
Definition: fasp_const.h:184
void fasp_chkerr(const SHORT status, const char *fctname)
Check error status and print out error messages before quit.
Definition: message.c:199
SHORT AMG_type
type of AMG method
Definition: fasp.h:586
SHORT coarse_scaling
switch of scaling of the coarse grid correction
Definition: fasp.h:840
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
SHORT presmooth_iter
number of presmoothing
Definition: fasp.h:822
REAL tentative_smooth
relaxation parameter for smoothing the tentative prolongation
Definition: fasp.h:676
SHORT smoother
smoother type
Definition: fasp.h:610
INT maxit
max number of iterations of AMG
Definition: fasp.h:592
SHORT smooth_order
AMG smoother ordering.
Definition: fasp_block.h:335
REAL tol
tolerance for AMG preconditioner
Definition: fasp_block.h:326
#define ERROR_INPUT_PAR
Definition: fasp_const.h:31
SHORT coarsening_type
coarsening type
Definition: fasp.h:643
SHORT smooth_order
smoother order
Definition: fasp.h:613
REAL relaxation
relaxation parameter for SOR smoother
Definition: fasp.h:828
INT coarse_dof
max number of coarsest level DOF
Definition: fasp.h:601
SHORT amli_degree
degree of the polynomial used by AMLI cycle
Definition: fasp_block.h:356
SHORT amli_degree
degree of the polynomial used by AMLI cycle
Definition: fasp.h:634
REAL AMG_quality_bound
Definition: fasp.h:1091
SHORT presmooth_iter
number of presmoothers
Definition: fasp.h:616
SHORT print_level
print level in AMG preconditioner
Definition: fasp.h:801
SHORT AMG_postsmooth_iter
Definition: fasp.h:1070
REAL AMG_max_row_sum
Definition: fasp.h:1087
#define PAIRWISE
Definition of aggregation types.
Definition: fasp_const.h:169
INT itsolver_maxit
Definition: fasp.h:1045
SHORT ILU_type
ILU type for decomposition.
Definition: fasp.h:380
SHORT smoother
AMG smoother type.
Definition: fasp.h:816
SHORT AMG_smooth_filter
Definition: fasp.h:1097
SHORT polynomial_degree
degree of the polynomial smoother
Definition: fasp.h:625
INT maxit
max number of iterations of AMG preconditioner
Definition: fasp.h:804
INT AMG_coarse_dof
Definition: fasp.h:1071
SHORT AMG_presmooth_iter
Definition: fasp.h:1069
SHORT AMG_aggregation_type
Definition: fasp.h:1083
SHORT coarse_solver
coarse solver type for AMG
Definition: fasp_block.h:350
void fasp_param_init(input_param *iniparam, itsolver_param *itsparam, AMG_param *amgparam, ILU_param *iluparam, Schwarz_param *schparam)
Initialize parameters, global variables, etc.
Definition: parameters.c:270
#define SOLVER_GCG
Definition: fasp_const.h:109
#define SMOOTHER_GSOR
Definition: fasp_const.h:189
REAL ILU_relax
add the sum of dropped elements to diagonal element in proportion relax
Definition: fasp.h:389
SHORT solver_type
Definition: fasp.h:1041
REAL * amli_coef
coefficients of the polynomial used by AMLI cycle
Definition: fasp_block.h:359
REAL quality_bound
quality threshold for pairwise aggregation
Definition: fasp.h:607
SHORT smooth_order
AMG smoother ordering.
Definition: fasp.h:819
void fasp_param_ilu_set(ILU_param *iluparam, input_param *iniparam)
Set ILU_param with INPUT.
Definition: parameters.c:593
#define INTERP_DIR
Definition of interpolation types.
Definition: fasp_const.h:206
Input parameters.
Definition: fasp.h:1029
#define SMOOTHER_SSOR
Definition: fasp_const.h:188
INT Schwarz_maxlvl
maximal level for constructing the blocks
Definition: fasp.h:443
REAL AMG_relaxation
Definition: fasp.h:1067
REAL ILU_permtol
permuted if permtol*|a(i,j)| > |a(i,i)|
Definition: fasp.h:392
SHORT precond_type
Definition: fasp.h:1042
INT ILU_lfil
level of fill-in for ILUk
Definition: fasp.h:383
REAL ILU_permtol
Definition: fasp.h:1053
SHORT AMG_smooth_order
Definition: fasp.h:1066
SHORT AMG_smoother
Definition: fasp.h:1065
void fasp_param_prec_to_amg(AMG_param *amgparam, precond_data *pcdata)
Set AMG_param with precond_data.
Definition: parameters.c:701
INT aggressive_path
number of paths use to determine strongly coupled C points
Definition: fasp.h:664
INT aggressive_level
number of levels use aggressive coarsening
Definition: fasp.h:661
SHORT nl_amli_krylov_type
type of Krylov method used by Nonlinear AMLI cycle
Definition: fasp.h:846
SHORT ILU_type
Definition: fasp.h:1049
INT Schwarz_maxlvl
Definition: fasp.h:1057
#define SMOOTHER_SOR
Definition: fasp_const.h:187
char workdir[256]
Definition: fasp.h:1037
#define PRINT_SOME
Definition: fasp_const.h:81
SHORT cycle_type
AMG cycle type.
Definition: fasp_block.h:329
#define SMOOTHER_POLY
Definition: fasp_const.h:191
INT max_levels
max number of AMG levels
Definition: fasp_block.h:323
SHORT AMG_type
type of AMG method
Definition: fasp.h:798
SHORT print_level
print level in AMG preconditioner
Definition: fasp_block.h:317
#define FASP_SUCCESS
Definition of return status and error messages.
Definition: fasp_const.h:27
REAL strong_threshold
strong connection threshold for coarsening
Definition: fasp.h:652
REAL AMG_strong_threshold
Definition: fasp.h:1085
void fasp_param_ilu_init(ILU_param *iluparam)
Initialize ILU parameters.
Definition: parameters.c:476
REAL AMG_tentative_smooth
Definition: fasp.h:1096
INT restart
Definition: fasp.h:1112
SHORT AMG_nl_amli_krylov_type
Definition: fasp.h:1078
INT restart
Definition: fasp.h:1046
INT Schwarz_levels
number of levels use Schwarz smoother
Definition: fasp.h:700
SHORT precond_type
Definition: fasp.h:1108
unsigned INT total_alloc_count
Total allocated memory amount.
Definition: memory.c:35
#define V_CYCLE
Definition of cycle types.
Definition: fasp_const.h:175
char inifile[256]
Definition: fasp.h:1036
REAL * amli_coef
coefficients of the polynomial used by AMLI cycle
Definition: fasp.h:852
SHORT stop_type
Definition: fasp.h:1109
#define ILUk
Type of ILU methods.
Definition: fasp_const.h:148
void fasp_param_amg_to_prec(precond_data *pcdata, AMG_param *amgparam)
Set precond_data with AMG_param.
Definition: parameters.c:666
SHORT aggregation_type
aggregation type
Definition: fasp.h:646
Data passed to the preconditioners.
Definition: fasp_block.h:311
void fasp_param_set(int argc, const char *argv[], input_param *iniparam)
Read input from command-line arguments.
Definition: parameters.c:27
INT Schwarz_blksolver
Definition: fasp.h:1059
INT ILU_lfil
Definition: fasp.h:1050
INT AMG_maxit
Definition: fasp.h:1073
SHORT AMG_coarse_scaling
Definition: fasp.h:1076
SHORT AMG_interpolation_type
Definition: fasp.h:1084
SHORT coarse_solver
coarse solver type
Definition: fasp.h:628
REAL truncation_threshold
truncation threshold
Definition: fasp.h:658
void fasp_param_Schwarz_set(Schwarz_param *schparam, input_param *iniparam)
Set Schwarz_param with INPUT.
Definition: parameters.c:615
SHORT ILU_levels
number of levels use ILU smoother
Definition: fasp.h:682
SHORT Schwarz_type
type for Schwarz method
Definition: fasp.h:440
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
REAL ILU_permtol
permuted if permtol*|a(i,j)| > |a(i,i)|
Definition: fasp.h:697
Main header file for FASP.
SHORT max_levels
max number of levels of AMG
Definition: fasp.h:598
REAL ILU_relax
Definition: fasp.h:1052
#define SOLVER_CG
Definition: fasp_const.h:103
void fasp_param_ilu_print(ILU_param *param)
Print out ILU parameters.
Definition: parameters.c:898
SHORT coarse_scaling
switch of scaling of the coarse grid correction
Definition: fasp_block.h:353
REAL * amli_coef
coefficients of the polynomial used by AMLI cycle
Definition: fasp.h:637
INT AMG_aggressive_path
Definition: fasp.h:1089
#define CLASSIC_AMG
Definition of AMG types.
Definition: fasp_const.h:162
SHORT AMG_coarse_solver
Definition: fasp.h:1075
#define SOLVER_GMRES
Definition: fasp_const.h:106
INT Schwarz_mmsize
maximal size of blocks
Definition: fasp.h:446
REAL ILU_relax
relaxation for ILUs
Definition: fasp.h:694
SHORT nl_amli_krylov_type
type of Krylov method used by Nonlinear AMLI cycle
Definition: fasp.h:640
SHORT smoother
AMG smoother type.
Definition: fasp_block.h:332
#define PRINT_NONE
Print level for all subroutines – not including DEBUG output.
Definition: fasp_const.h:79
#define OFF
Definition: fasp_const.h:74
SHORT coarse_solver
coarse solver type for AMG
Definition: fasp.h:837
INT AMG_max_aggregation
Definition: fasp.h:1095
REAL relaxation
relaxation parameter for SOR smoother
Definition: fasp.h:622
SHORT fasp_param_check(input_param *inparam)
Simple check on input parameters.
Definition: input.c:25
REAL ILU_droptol
drop tolerance for ILUt
Definition: fasp.h:386
Parameters passed to iterative solvers.
Definition: fasp.h:1105
INT Schwarz_blksolver
type of Schwarz block solver
Definition: fasp.h:449
INT Schwarz_blksolver
type of Schwarz block solver
Definition: fasp.h:712
AMG_data_bsr * mgl_data
AMG preconditioner data.
Definition: fasp_block.h:368
void fasp_param_Schwarz_init(Schwarz_param *schparam)
Initialize Schwarz parameters.
Definition: parameters.c:498
REAL tentative_smooth
smooth factor for smoothing the tentative prolongation
Definition: fasp.h:849
#define STOP_REL_RES
Definition of iterative solver stopping criteria types.
Definition: fasp_const.h:131
INT ILU_levels
number of levels use ILU smoother
Definition: fasp.h:761
SHORT ILU_type
ILU type for smoothing.
Definition: fasp.h:685
#define PREC_AMG
Definition: fasp_const.h:140
REAL tol
Definition: fasp.h:1111
void fasp_param_amg_init(AMG_param *amgparam)
Initialize AMG parameters.
Definition: parameters.c:390
SHORT cycle_type
type of AMG cycle
Definition: fasp.h:604
SHORT max_levels
max number of AMG levels
Definition: fasp.h:807
INT AMG_pair_number
Definition: fasp.h:1090
SHORT smooth_filter
switch for filtered matrix used for smoothing the tentative prolongation
Definition: fasp.h:679
INT AMG_aggressive_level
Definition: fasp.h:1088
SHORT polynomial_degree
degree of the polynomial smoother
Definition: fasp.h:831
SHORT print_level
Definition: fasp.h:1113
Parameters for Schwarz method.
Definition: fasp.h:434
REAL max_row_sum
maximal row sum parameter
Definition: fasp.h:655
AMG_data * mgl_data
AMG preconditioner data.
Definition: fasp.h:855
Parameters for ILU.
Definition: fasp.h:374
#define SHORT
FASP integer and floating point numbers.
Definition: fasp.h:63
INT ILU_levels
number of levels use ILU smoother
Definition: fasp_block.h:255
#define AMLI_CYCLE
Definition: fasp_const.h:177
SHORT interpolation_type
interpolation type
Definition: fasp.h:649
void fasp_param_solver_print(itsolver_param *param)
Print out itsolver parameters.
Definition: parameters.c:957
SHORT coarsening_type
coarsening type
Definition: fasp_block.h:344
#define VMB
Definition: fasp_const.h:170
SHORT print_level
print leve
Definition: fasp.h:437
SHORT stop_type
Definition: fasp.h:1043
void fasp_param_amg_to_prec_bsr(precond_data_bsr *pcdata, AMG_param *amgparam)
Set precond_data_bsr with AMG_param.
Definition: parameters.c:733
INT Schwarz_mmsize
maximal block size
Definition: fasp.h:703
REAL tol
tolerance for AMG preconditioner
Definition: fasp.h:810
#define COARSE_RS
Definition of coarsening types.
Definition: fasp_const.h:197
INT pair_number
number of pairwise matchings
Definition: fasp.h:667
REAL itsolver_tol
Definition: fasp.h:1044
REAL AMG_tol
Definition: fasp.h:1072
#define SOLVER_AMG
Definition: fasp_const.h:120
INT problem_num
Definition: fasp.h:1038
REAL ILU_droptol
Definition: fasp.h:1051
SHORT print_level
print level for AMG
Definition: fasp.h:589
void fasp_param_input(const char *filenm, input_param *inparam)
Read input parameters from disk file.
Definition: input.c:102
#define ON
Definition of switch.
Definition: fasp_const.h:73
INT Schwarz_type
Definition: fasp.h:1058
unsigned INT total_alloc_mem
Definition: memory.c:34
void fasp_param_prec_to_amg_bsr(AMG_param *amgparam, precond_data_bsr *pcdata)
Set AMG_param with precond_data.
Definition: parameters.c:767
REAL AMG_strong_coupled
Definition: fasp.h:1094
SHORT AMG_type
Definition: fasp.h:1062
REAL tentative_smooth
smooth factor for smoothing the tentative prolongation
Definition: fasp_block.h:362
Data passed to the preconditioners.
Definition: fasp.h:795
SHORT nl_amli_krylov_type
type of krylov method used by Nonlinear AMLI cycle
Definition: fasp_block.h:365
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
INT maxit
max number of iterations of AMG preconditioner
Definition: fasp_block.h:320
INT AMG_Schwarz_levels
Definition: fasp.h:1079
SHORT cycle_type
AMG cycle type.
Definition: fasp.h:813
SHORT AMG_type
type of AMG method
Definition: fasp_block.h:314
#define CF_ORDER
Definition: fasp_const.h:223
INT Schwarz_maxlvl
maximal levels
Definition: fasp.h:706
#define SMOOTHER_SGSOR
Definition: fasp_const.h:190
SHORT presmooth_iter
number of presmoothing
Definition: fasp_block.h:338
#define NL_AMLI_CYCLE
Definition: fasp_const.h:178
SHORT AMG_cycle_type
Definition: fasp.h:1064
void fasp_param_amg_set(AMG_param *param, input_param *iniparam)
Set AMG_param from INPUT.
Definition: parameters.c:518
SHORT postsmooth_iter
number of postsmoothing
Definition: fasp_block.h:341
#define SOLVER_VGMRES
Definition: fasp_const.h:107
SHORT AMG_coarsening_type
Definition: fasp.h:1082
REAL ILU_droptol
drop tolerance for ILUt
Definition: fasp.h:691
REAL relaxation
relaxation parameter for SOR smoother
Definition: fasp_block.h:347
INT Schwarz_type
type of Schwarz method
Definition: fasp.h:709
void fasp_param_amg_print(AMG_param *param)
Print out AMG parameters.
Definition: parameters.c:797
SHORT print_level
Definition: fasp.h:1032
INT Schwarz_mmsize
Definition: fasp.h:1056
INT ILU_lfil
level of fill-in for ILUs and ILUk
Definition: fasp.h:688
#define SOLVER_FMG
Definition: fasp_const.h:121
SHORT AMG_amli_degree
Definition: fasp.h:1077
REAL AMG_truncation_threshold
Definition: fasp.h:1086
#define SOLVER_DEFAULT
Definition of solver types for iterative methods.
Definition: fasp_const.h:101
SHORT AMG_levels
Definition: fasp.h:1063
SHORT output_type
Definition: fasp.h:1033
SHORT coarse_scaling
switch of scaling of the coarse grid correction
Definition: fasp.h:631
SHORT AMG_ILU_levels
Definition: fasp.h:1074
void fasp_param_solver_set(itsolver_param *itsparam, input_param *iniparam)
Set itsolver_param with INPUT.
Definition: parameters.c:636
SHORT AMG_polynomial_degree
Definition: fasp.h:1068
SHORT itsolver_type
Definition: fasp.h:1107
void fasp_param_Schwarz_print(Schwarz_param *param)
Print out Schwarz parameters.
Definition: parameters.c:928
SHORT print_level
print level
Definition: fasp.h:377