OpenCAEPoro  v0.5.0
A simulator for multicomponent porous media flow
FaspSolver.cpp
Go to the documentation of this file.
1 
12 #include <math.h>
13 
14 #include "FaspSolver.hpp"
15 
16 void FaspSolver::SetupParam(const string& dir, const string& file)
17 {
18  solveDir = dir;
19  solveFile = file;
20 
21  string myfile = solveDir + solveFile;
22  InitParam(); // Set default solver parameters
23  ifstream ifs(myfile);
24  if (!ifs.is_open()) {
25  cout << "The input file " << myfile << " is missing!" << endl;
26  myfile = solveDir + "../conf/csr.fasp";
27  ifs.open(myfile);
28  if (!ifs.is_open()) {
29  cout << "The input file " << myfile << " is missing!" << endl;
30  cout << "Using the default parameters of FASP" << endl;
31  } else {
32  ifs.close();
33  cout << "Using the input file " << myfile << endl;
34  fasp_param_input(myfile.data(), &inParam);
35  }
36  } else {
37  ifs.close(); // if file has been opened, close it first
38  fasp_param_input(myfile.data(), &inParam);
39  }
40  fasp_param_init(&inParam, &itParam, &amgParam, &iluParam, &swzParam);
41 }
42 
43 void ScalarFaspSolver::Allocate(const vector<USI>& rowCapacity,
44  const OCP_USI& maxDim,
45  const USI& blockDim)
46 {
47  OCP_USI nnz = 0;
48  for (OCP_USI n = 0; n < maxDim; n++) {
49  nnz += rowCapacity[n];
50  }
51  A = fasp_dcsr_create(maxDim, maxDim, nnz);
52 }
53 
54 void ScalarFaspSolver::InitParam()
55 {
56  // Input/output
57  inParam.print_level = PRINT_MIN;
58  inParam.output_type = 0;
59 
60  // Problem information
61  inParam.solver_type = SOLVER_VFGMRES;
62  inParam.precond_type = PREC_AMG;
63  inParam.stop_type = STOP_REL_RES;
64 
65  // LinearSolver parameters
66  inParam.itsolver_tol = 1e-4;
67  inParam.itsolver_maxit = 100;
68  inParam.restart = 30;
69 
70  // ILU method parameters
71  inParam.ILU_type = ILUk;
72  inParam.ILU_lfil = 0;
73  inParam.ILU_droptol = 0.001;
74  inParam.ILU_relax = 0;
75  inParam.ILU_permtol = 0.0;
76 
77  // Schwarz method parameters
78  inParam.SWZ_mmsize = 200;
79  inParam.SWZ_maxlvl = 2;
80  inParam.SWZ_type = 1;
81  inParam.SWZ_blksolver = SOLVER_DEFAULT;
82 
83  // AMG method parameters
84  inParam.AMG_type = CLASSIC_AMG;
85  inParam.AMG_levels = 20;
86  inParam.AMG_cycle_type = V_CYCLE;
87  inParam.AMG_smoother = SMOOTHER_GS;
88  inParam.AMG_smooth_order = CF_ORDER;
89  inParam.AMG_presmooth_iter = 1;
90  inParam.AMG_postsmooth_iter = 1;
91  inParam.AMG_relaxation = 1.0;
92  inParam.AMG_coarse_dof = 500;
93  inParam.AMG_coarse_solver = 0;
94  inParam.AMG_tol = 1e-6;
95  inParam.AMG_maxit = 1;
96  inParam.AMG_ILU_levels = 0;
97  inParam.AMG_SWZ_levels = 0;
98  inParam.AMG_coarse_scaling = OFF;
99  inParam.AMG_amli_degree = 1;
100  inParam.AMG_nl_amli_krylov_type = 2;
101 
102  // Classical AMG specific
103  inParam.AMG_coarsening_type = 1;
104  inParam.AMG_interpolation_type = 1;
105  inParam.AMG_max_row_sum = 0.9;
106  inParam.AMG_strong_threshold = 0.3;
107  inParam.AMG_truncation_threshold = 0.2;
108  inParam.AMG_aggressive_level = 0;
109  inParam.AMG_aggressive_path = 1;
110 
111  // Aggregation AMG specific
112  inParam.AMG_aggregation_type = PAIRWISE;
113  inParam.AMG_quality_bound = 8.0;
114  inParam.AMG_pair_number = 2;
115  inParam.AMG_strong_coupled = 0.25;
116  inParam.AMG_max_aggregation = 9;
117  inParam.AMG_tentative_smooth = 0.67;
118  inParam.AMG_smooth_filter = ON;
119  inParam.AMG_smooth_restriction = ON;
120 }
121 
122 void ScalarFaspSolver::AssembleMat(const vector<vector<USI>>& colId,
123  const vector<vector<OCP_DBL>>& val,
124  const OCP_USI& dim,
125  const USI& blockDim,
126  vector<OCP_DBL>& rhs,
127  vector<OCP_DBL>& u)
128 {
129  // b & x
130  b.row = dim;
131  b.val = rhs.data();
132  x.row = dim;
133  x.val = u.data();
134  // A
135  OCP_USI nnz = 0;
136  for (OCP_USI i = 0; i < dim; i++) {
137  nnz += colId[i].size();
138  }
139 
140  A.row = dim;
141  A.col = dim;
142  A.nnz = nnz;
143 
144  // IA
145  OCP_USI count = 0;
146  A.IA[0] = 0;
147  for (OCP_USI i = 1; i < dim + 1; i++) {
148  USI nnz_Row = colId[i - 1].size();
149  A.IA[i] = A.IA[i - 1] + nnz_Row;
150 
151  for (USI j = 0; j < nnz_Row; j++) {
152  A.JA[count] = colId[i - 1][j];
153  A.val[count] = val[i - 1][j];
154  count++;
155  }
156  }
157 
158 #ifdef DEBUG
159  // check x and b ---- for test
160  for (int i = 0; i < dim; i++) {
161  if (!isfinite(b.val[i])) OCP_ABORT("sFasp b is infinite!");
162  if (!isfinite(x.val[i])) OCP_ABORT("sFasp x is infinite!");
163  }
164  // check A ---- for test
165  for (int i = 0; i < nnz; i++) {
166  if (!isfinite(A.val[i])) OCP_ABORT("sFasp A is infinite!");
167  }
168 #endif // DEBUG
169 }
170 
171 OCP_INT ScalarFaspSolver::Solve()
172 {
173  OCP_INT status = FASP_SUCCESS;
174 
175  const OCP_INT print_level = inParam.print_level;
176  const OCP_INT solver_type = inParam.solver_type;
177  const OCP_INT precond_type = inParam.precond_type;
178  const OCP_INT output_type = inParam.output_type;
179 
180  if (output_type) {
181  const char* outputfile = "out/Solver.out";
182  printf("Redirecting outputs to file: %s ...\n", outputfile);
183  freopen(outputfile, "w", stdout); // open a file for stdout
184  }
185 
186  // Preconditioned Krylov methods
187  if (solver_type >= 1 && solver_type <= 20) {
188 
189  // Using no preconditioner for Krylov iterative methods
190  if (precond_type == PREC_NULL) {
191  status = fasp_solver_dcsr_krylov(&A, &b, &x, &itParam);
192  }
193 
194  // Using diag(A) as preconditioner for Krylov iterative methods
195  else if (precond_type == PREC_DIAG) {
196  status = fasp_solver_dcsr_krylov_diag(&A, &b, &x, &itParam);
197  }
198 
199  // Using AMG as preconditioner for Krylov iterative methods
200  else if (precond_type == PREC_AMG || precond_type == PREC_FMG) {
201  if (print_level > PRINT_NONE) fasp_param_amg_print(&amgParam);
202  status = fasp_solver_dcsr_krylov_amg(&A, &b, &x, &itParam, &amgParam);
203  }
204 
205  // Using ILU as preconditioner for Krylov iterative methods
206  else if (precond_type == PREC_ILU) {
207  if (print_level > PRINT_NONE) fasp_param_ilu_print(&iluParam);
208  status = fasp_solver_dcsr_krylov_ilu(&A, &b, &x, &itParam, &iluParam);
209  }
210 
211  // Undefined iterative methods
212  else {
213  printf("### ERROR: Wrong preconditioner type %d!!!\n", precond_type);
214  status = ERROR_SOLVER_PRECTYPE;
215  }
216 
217  }
218 
219  // AMG as the iterative solver
220  else if (solver_type == SOLVER_AMG) {
221  if (print_level > PRINT_NONE) fasp_param_amg_print(&amgParam);
222  fasp_solver_amg(&A, &b, &x, &amgParam);
223  }
224 
225  // Full AMG as the iterative solver
226  else if (solver_type == SOLVER_FMG) {
227  if (print_level > PRINT_NONE) fasp_param_amg_print(&amgParam);
228  fasp_solver_famg(&A, &b, &x, &amgParam);
229  }
230 
231 #if WITH_MUMPS // use MUMPS directly
232  else if (solver_type == SOLVER_MUMPS) {
233  status = fasp_solver_mumps(&A, &b, &x, print_level);
234  if (status >= 0) status = 1; // Direct solver returns 1
235  }
236 #endif
237 
238 #if WITH_SuperLU // use SuperLU directly
239  else if (solver_type == SOLVER_SUPERLU) {
240  status = fasp_solver_superlu(&A, &b, &x, print_level);
241  if (status >= 0) status = 1; // Direct solver returns 1
242  }
243 #endif
244 
245 #if WITH_UMFPACK // use UMFPACK directly
246  else if (solver_type == SOLVER_UMFPACK) {
247  // Need to sort the matrix A for UMFPACK to work
248  dCSRmat A_trans = fasp_dcsr_create(A.row, A.col, A.nnz);
249  fasp_dcsr_transz(&A, NULL, &A_trans);
250  fasp_dcsr_sort(&A_trans);
251  status = fasp_solver_umfpack(&A_trans, &b, &x, print_level);
252  fasp_dcsr_free(&A_trans);
253  if (status >= 0) status = 1; // Direct solver returns 1
254  }
255 #endif
256 
257 #ifdef WITH_PARDISO // use PARDISO directly
258  else if (solver_type == SOLVER_PARDISO) {
259  fasp_dcsr_sort(&A);
260  status = fasp_solver_pardiso(&A, &b, &x, print_level);
261  if (status >= 0) status = 1; // Direct solver returns 1
262  }
263 #endif
264 
265  else {
266  printf("### ERROR: Wrong solver type %d!!!\n", solver_type);
267  status = ERROR_SOLVER_TYPE;
268  }
269 
270  if (status < 0) {
271  printf("### ERROR: Solver failed! Exit status = %d.\n\n", status);
272  }
273 
274  if (output_type) fclose(stdout);
275  return status;
276 }
277 
278 void VectorFaspSolver::Allocate(const vector<USI>& rowCapacity,
279  const OCP_USI& maxDim,
280  const USI& blockDim)
281 {
282  OCP_USI nnz = 0;
283  for (OCP_USI n = 0; n < maxDim; n++) {
284  nnz += rowCapacity[n];
285  }
286  A = fasp_dbsr_create(maxDim, maxDim, nnz, blockDim, 0);
287  Asc = fasp_dbsr_create(maxDim, maxDim, nnz, blockDim, 0);
288  fsc = fasp_dvec_create(maxDim * blockDim);
289  order = fasp_ivec_create(maxDim);
290  Dmat.resize(maxDim * blockDim * blockDim);
291 }
292 
293 void VectorFaspSolver::InitParam()
294 {
295  // Input/output
296  inParam.print_level = PRINT_MIN;
297  inParam.output_type = 0;
298 
299  // Problem information
300  inParam.solver_type = SOLVER_VFGMRES;
301  inParam.decoup_type = 1;
302  inParam.precond_type = 64;
303  inParam.stop_type = STOP_REL_RES;
304 
305  // Solver parameters
306  inParam.itsolver_tol = 1e-3;
307  inParam.itsolver_maxit = 100;
308  inParam.restart = 30;
309 
310  // ILU method parameters
311  inParam.ILU_type = ILUk;
312  inParam.ILU_lfil = 0;
313  inParam.ILU_droptol = 0.001;
314  inParam.ILU_relax = 0;
315  inParam.ILU_permtol = 0.0;
316 
317  // Schwarz method parameters
318  inParam.SWZ_mmsize = 200;
319  inParam.SWZ_maxlvl = 2;
320  inParam.SWZ_type = 1;
321  inParam.SWZ_blksolver = SOLVER_DEFAULT;
322 
323  // AMG method parameters
324  inParam.AMG_type = CLASSIC_AMG;
325  inParam.AMG_levels = 20;
326  inParam.AMG_cycle_type = V_CYCLE;
327  inParam.AMG_smoother = SMOOTHER_GS;
328  inParam.AMG_smooth_order = CF_ORDER;
329  inParam.AMG_presmooth_iter = 1;
330  inParam.AMG_postsmooth_iter = 1;
331  inParam.AMG_relaxation = 1.0;
332  inParam.AMG_coarse_dof = 500;
333  inParam.AMG_coarse_solver = 0;
334  inParam.AMG_tol = 1e-6;
335  inParam.AMG_maxit = 1;
336  inParam.AMG_ILU_levels = 0;
337  inParam.AMG_SWZ_levels = 0;
338  inParam.AMG_coarse_scaling = OFF; // Require investigation --Chensong
339  inParam.AMG_amli_degree = 1;
340  inParam.AMG_nl_amli_krylov_type = 2;
341 
342  // Classical AMG specific
343  inParam.AMG_coarsening_type = 1;
344  inParam.AMG_interpolation_type = 1;
345  inParam.AMG_max_row_sum = 0.9;
346  inParam.AMG_strong_threshold = 0.3;
347  inParam.AMG_truncation_threshold = 0.2;
348  inParam.AMG_aggressive_level = 0;
349  inParam.AMG_aggressive_path = 1;
350 
351  // Aggregation AMG specific
352  inParam.AMG_aggregation_type = PAIRWISE;
353  inParam.AMG_quality_bound = 8.0;
354  inParam.AMG_pair_number = 2;
355  inParam.AMG_strong_coupled = 0.25;
356  inParam.AMG_max_aggregation = 9;
357  inParam.AMG_tentative_smooth = 0.67;
358  inParam.AMG_smooth_filter = ON;
359  inParam.AMG_smooth_restriction = ON;
360 }
361 
362 void VectorFaspSolver::AssembleMat(const vector<vector<USI>>& colId,
363  const vector<vector<OCP_DBL>>& val,
364  const OCP_USI& dim,
365  const USI& blockDim,
366  vector<OCP_DBL>& rhs,
367  vector<OCP_DBL>& u)
368 {
369  const OCP_USI nrow = dim * blockDim;
370  // b & x
371  b.row = nrow;
372  b.val = rhs.data();
373  x.row = nrow;
374  x.val = u.data(); // x will be set to zero later
375 
376  // fsc & order
377  fsc.row = nrow;
378  order.row = nrow;
379 
380  // nnz
381  OCP_USI nnz = 0;
382  for (OCP_USI i = 0; i < dim; i++) {
383  nnz += colId[i].size();
384  }
385 
386  // Asc
387  Asc.ROW = dim;
388  Asc.COL = dim;
389  Asc.nb = blockDim;
390  Asc.NNZ = nnz;
391 
392  // A
393  A.ROW = dim;
394  A.COL = dim;
395  A.nb = blockDim;
396  A.NNZ = nnz;
397 
398  OCP_USI count1 = 0;
399  OCP_USI count2 = 0;
400  OCP_USI size_row;
401  A.IA[0] = 0;
402  for (OCP_USI i = 1; i < dim + 1; i++) {
403  USI nnb_Row = colId[i - 1].size();
404  A.IA[i] = A.IA[i - 1] + nnb_Row;
405 
406  for (USI j = 0; j < nnb_Row; j++) {
407  A.JA[count1] = colId[i - 1][j];
408  count1++;
409  }
410  size_row = nnb_Row * blockDim * blockDim;
411  for (USI k = 0; k < size_row; k++) {
412  A.val[count2] = val[i - 1][k];
413  count2++;
414  }
415  }
416 
417 #ifdef DEBUG
418  // check x and b ---- for test
419  for (int i = 0; i < nrow; i++) {
420  if (!isfinite(b.val[i])) OCP_ABORT("vFasp b is infinite!");
421  if (!isfinite(x.val[i])) OCP_ABORT("vFasp x is infinite!");
422  }
423  // check A ---- for test
424  const OCP_USI len = A.NNZ * blockDim * blockDim;
425  for (int i = 0; i < len; i++) {
426  if (!isfinite(A.val[i])) OCP_ABORT("vFasp A is infinite!");
427  }
428 #endif // DEBUG
429 }
430 
431 OCP_INT VectorFaspSolver::Solve()
432 {
433  OCP_INT status = FASP_SUCCESS;
434 
435  // Set local parameters
436  const OCP_INT print_level = inParam.print_level;
437  const OCP_INT solver_type = inParam.solver_type;
438  const OCP_INT precond_type = inParam.precond_type;
439  const OCP_INT output_type = inParam.output_type;
440 
441 #if WITH_FASP4BLKOIL || WITH_FASPCPR // Currently, only fasp4blkoil requires decoupling
442  const OCP_INT decoup_type = inParam.decoup_type;
443 #endif
444 
445  if (output_type) {
446  const char* outputfile = "../output/test.out";
447  printf("Redirecting outputs to file: %s ...\n", outputfile);
448  freopen(outputfile, "w", stdout); // open a file for stdout
449  }
450 
451  fasp_dvec_set(x.row, &x, 0);
452 
453  // Preconditioned Krylov methods
454  if (solver_type >= 1 && solver_type <= 10) {
455 
456  // Preconditioned Krylov methods in BSR format
457  switch (precond_type) {
458  case PC_NULL:
459  status = fasp_solver_dbsr_krylov(&A, &b, &x, &itParam);
460  break;
461  case PC_DIAG:
462  status = fasp_solver_dbsr_krylov_diag(&A, &b, &x, &itParam);
463  break;
464  case PC_BILU:
465  status = fasp_solver_dbsr_krylov_ilu(&A, &b, &x, &itParam, &iluParam);
466  break;
467 
468 #if WITH_FASPCPR
470  case PC_FASP1:
471  Decoupling(&A, &b, &Asc, &fsc, &order, Dmat.data(), decoup_type);
472  status = FASP_BSRSOL_ASCPR(&Asc, &fsc, &x, &itParam, &iluParam,
473  &amgParam, 0);
474  break;
475 
476  case PC_FASP1_SHARE:
477  Decoupling(&A, &b, &Asc, &fsc, &order, Dmat.data(), decoup_type);
478  status = FASP_BSRSOL_ASCPR(&Asc, &fsc, &x, &itParam, &iluParam,
480  break;
481 #endif
482 
483 #if WITH_FASP4BLKOIL
484  case PC_FASP1:
485  Decoupling(&A, &b, &Asc, &fsc, &order, Dmat.data(), decoup_type);
486 #if WITH_FASP4CUDA // zhaoli 2022.04.04
487  status = fasp_solver_dbsr_krylov_FASP1_cuda_interface(
488  &Asc, &fsc, &x, &itParam, &iluParam, &amgParam, NULL, &order);
489 #else
490  status = fasp_solver_dbsr_krylov_FASP1a(
491  &Asc, &fsc, &x, &itParam, &iluParam, &amgParam, NULL, &order);
492 #endif
493  break;
494  case PC_FASP1_SHARE: // zhaoli 2021.03.24
495  Decoupling(&A, &b, &Asc, &fsc, &order, Dmat.data(), decoup_type);
496 #if WITH_FASP4CUDA
497  status = fasp_solver_dbsr_krylov_FASP1_cuda_share_interface(
498  &Asc, &fsc, &x, &itParam, &iluParam, &amgParam, NULL, &order,
499  RESET_CONST);
500 #else
501  status = fasp_solver_dbsr_krylov_FASP1a_share_interface(
502  &Asc, &fsc, &x, &itParam, &iluParam, &amgParam, NULL, &order,
503  RESET_CONST);
504 #endif
505  break;
506  case PC_FASP2:
507  Decoupling(&A, &b, &Asc, &fsc, &order, Dmat.data(), decoup_type);
508  status = fasp_solver_dbsr_krylov_FASP2(
509  &Asc, &fsc, &x, &itParam, &iluParam, &amgParam, NULL, &order);
510  break;
511  case PC_FASP3:
512  Decoupling(&A, &b, &Asc, &fsc, &order, Dmat.data(), decoup_type);
513  status = fasp_solver_dbsr_krylov_FASP3(
514  &Asc, &fsc, &x, &itParam, &iluParam, &amgParam, NULL, &order);
515  break;
516  case PC_FASP4:
517  Decoupling(&A, &b, &Asc, &fsc, &order, Dmat.data(), decoup_type);
518 #if WITH_FASP4CUDA
519  status = fasp_solver_dbsr_krylov_FASP4_cuda(
520  &Asc, &fsc, &x, &itParam, &iluParam, &amgParam, NULL, &order);
521 #else
522  status = fasp_solver_dbsr_krylov_FASP4(
523  &Asc, &fsc, &x, &itParam, &iluParam, &amgParam, NULL, &order);
524 #endif
525  break;
526  case PC_FASP4_SHARE: // zhaoli 2021.04.24
527  Decoupling(&A, &b, &Asc, &fsc, &order, Dmat.data(), decoup_type);
528 #if WITH_FASP4CUDA // zhaoli 2022.08.03
529  status = fasp_solver_dbsr_krylov_FASP4_cuda_share_interface(
530  &Asc, &fsc, &x, &itParam, &iluParam, &amgParam, NULL, &order,
531  RESET_CONST);
532 #else
533  status = fasp_solver_dbsr_krylov_FASP4_share_interface(
534  &Asc, &fsc, &x, &itParam, &iluParam, &amgParam, NULL, &order,
535  RESET_CONST);
536 #endif
537  break;
538  case PC_FASP5:
539  Decoupling(&A, &b, &Asc, &fsc, &order, Dmat.data(), decoup_type);
540  status = fasp_solver_dbsr_krylov_FASP5(
541  &Asc, &fsc, &x, &itParam, &iluParam, &amgParam, NULL, &order);
542  break;
543 #endif
544  default:
545  OCP_WARNING("fasp4blkoil was not linked correctly!");
546  OCP_ABORT("Preconditioner type " + to_string(precond_type) +
547  " not supported!");
548  }
549  fill(Dmat.begin(), Dmat.end(), 0.0);
550  }
551 
552 #if WITH_MUMPS // use MUMPS directly
553  else if (solver_type == SOLVER_MUMPS) {
554  dCSRmat Acsr = fasp_format_dbsr_dcsr(&A);
555  status = fasp_solver_mumps(&Acsr, &b, &x, print_level);
556  fasp_dcsr_free(&Acsr);
557  if (status >= 0) status = 1; // Direct solver returns 1
558  }
559 #endif
560 
561 #if WITH_SuperLU // use SuperLU directly
562  else if (solver_type == SOLVER_SUPERLU) {
563  dCSRmat Acsr = fasp_format_dbsr_dcsr(&A);
564  status = fasp_solver_superlu(&Acsr, &b, &x, print_level);
565  fasp_dcsr_free(&Acsr);
566  if (status >= 0) status = 1; // Direct solver returns 1
567  }
568 #endif
569 
570 #if WITH_UMFPACK // use UMFPACK directly
571  else if (solver_type == SOLVER_UMFPACK) {
572  // Need to sort the matrix A for UMFPACK to work
573  dCSRmat Acsr = fasp_format_dbsr_dcsr(&A);
574  dCSRmat A_trans = fasp_dcsr_create(Acsr.row, Acsr.col, Acsr.nnz);
575  fasp_dcsr_transz(&Acsr, NULL, &A_trans);
576  fasp_dcsr_sort(&A_trans);
577  status = fasp_solver_umfpack(&A_trans, &b, &x, print_level);
578  fasp_dcsr_free(&A_trans);
579  fasp_dcsr_free(&Acsr);
580  if (status >= 0) status = 1; // Direct solver returns 1
581  }
582 #endif
583 
584 #ifdef WITH_PARDISO // use PARDISO directly
585  else if (solver_type == SOLVER_PARDISO) {
586  dCSRmat Acsr = fasp_format_dbsr_dcsr(&A);
587  fasp_dcsr_sort(&Acsr);
588  status = fasp_solver_pardiso(&Acsr, &b, &x, print_level);
589  fasp_dcsr_free(&Acsr);
590  if (status >= 0) status = 1; // Direct solver returns 1
591  }
592 #endif
593 
594  else {
595  printf("### ERROR: Wrong solver type %d!!!\n", solver_type);
596  status = ERROR_SOLVER_TYPE;
597  }
598 
599  if (print_level > PRINT_MIN) {
600  if (status < 0) {
601  cout << "\n### WARNING: Solver does not converge!\n" << endl;
602  } else {
603  cout << "\nSolver converges successfully!\n" << endl;
604  }
605  }
606 
607  if (output_type) fclose(stdout);
608 
609  return status;
610 }
611 
612 /*----------------------------------------------------------------------------*/
613 /* Brief Change History of This File */
614 /*----------------------------------------------------------------------------*/
615 /* Author Date Actions */
616 /*----------------------------------------------------------------------------*/
617 /* Shizhe Li Nov/22/2021 Create file */
618 /* Chensong Zhang Jan/19/2022 Set FASP4BLKOIL as optional */
619 /* Li Zhao Apr/04/2022 Set FASP4CUDA as optional */
620 /*----------------------------------------------------------------------------*/
Declaration of classes interfacing to the FASP solvers.
#define PC_BILU
BILU: block ILU preconditioner.
Definition: FaspSolver.hpp:64
#define PC_FASP2
FASP2: MSP, experimental only.
Definition: FaspSolver.hpp:59
#define PC_FASP3
FASP3: MSP, monolithic preconditioner.
Definition: FaspSolver.hpp:60
#define PC_FASP5
FASP5: MSP, experimental only.
Definition: FaspSolver.hpp:62
#define RESET_CONST
Sharing threshold for PC_FASP1_SHARE, PC_FASP4_SHARE.
Definition: FaspSolver.hpp:69
#define PC_NULL
None: no preconditioner.
Definition: FaspSolver.hpp:57
#define PC_FASP4
FASP4: MSP, default for FIM from 2015.
Definition: FaspSolver.hpp:61
#define PC_DIAG
DIAG: diagonal preconditioner.
Definition: FaspSolver.hpp:63
#define PC_FASP1_SHARE
Sharing setup stage for PC_FASP1, use with caution.
Definition: FaspSolver.hpp:67
#define PC_FASP4_SHARE
Sharing setup stage for PC_FASP4, use with caution.
Definition: FaspSolver.hpp:68
#define PC_FASP1
FASP1: MSP, default for FIM from 2020.
Definition: FaspSolver.hpp:58
unsigned int USI
Generic unsigned integer.
Definition: OCPConst.hpp:23
unsigned int OCP_USI
Long unsigned integer.
Definition: OCPConst.hpp:25
int OCP_INT
Long integer.
Definition: OCPConst.hpp:26
#define OCP_WARNING(msg)
Log warning messages.
Definition: UtilError.hpp:39
#define OCP_ABORT(msg)
Abort if critical error happens.
Definition: UtilError.hpp:47
AMG_param amgParam
Parameters for AMG method.
Definition: FaspSolver.hpp:86
input_param inParam
Parameters from input files.
Definition: FaspSolver.hpp:84
ILU_param iluParam
Parameters for ILU method.
Definition: FaspSolver.hpp:87
string solveFile
Relative path of fasp file.
Definition: FaspSolver.hpp:83
SWZ_param swzParam
Parameters for Schwarz method.
Definition: FaspSolver.hpp:88
string solveDir
Current work dir.
Definition: FaspSolver.hpp:82
ITS_param itParam
Parameters for iterative method.
Definition: FaspSolver.hpp:85
void SetupParam(const string &dir, const string &file) override
Set FASP parameters.
Definition: FaspSolver.cpp:16
virtual void InitParam()=0
Initialize the params for linear solvers.