24#ifndef OPM_FLOWLINEARSOLVERPARAMETERS_HEADER_INCLUDED
25#define OPM_FLOWLINEARSOLVERPARAMETERS_HEADER_INCLUDED
27#include <opm/simulators/linalg/MILU.hpp>
29#include <opm/simulators/linalg/linalgproperties.hh>
30#include <opm/models/utils/parametersystem.hh>
33template <
class TypeTag>
38namespace Opm::Properties {
44template<
class TypeTag,
class MyTypeTag>
46 using type = UndefinedProperty;
49template<
class TypeTag,
class MyTypeTag>
51 using type = UndefinedProperty;
54template<
class TypeTag,
class MyTypeTag>
56 using type = UndefinedProperty;
58template<
class TypeTag,
class MyTypeTag>
60 using type = UndefinedProperty;
62template<
class TypeTag,
class MyTypeTag>
64 using type = UndefinedProperty;
69template<
class TypeTag,
class MyTypeTag>
71 using type = UndefinedProperty;
73template<
class TypeTag,
class MyTypeTag>
75 using type = UndefinedProperty;
77template<
class TypeTag,
class MyTypeTag>
79 using type = UndefinedProperty;
81template<
class TypeTag,
class MyTypeTag>
83 using type = UndefinedProperty;
85template<
class TypeTag,
class MyTypeTag>
87 using type = UndefinedProperty;
89template<
class TypeTag,
class MyTypeTag>
91 using type = UndefinedProperty;
93template<
class TypeTag,
class MyTypeTag>
95 using type = UndefinedProperty;
97template<
class TypeTag,
class MyTypeTag>
99 using type = UndefinedProperty;
101template<
class TypeTag,
class MyTypeTag>
103 using type = UndefinedProperty;
105template<
class TypeTag,
class MyTypeTag>
107 using type = UndefinedProperty;
109template<
class TypeTag,
class MyTypeTag>
111 using type = UndefinedProperty;
113template<
class TypeTag,
class MyTypeTag>
115 using type = UndefinedProperty;
117template<
class TypeTag,
class MyTypeTag>
119 using type = UndefinedProperty;
121template<
class TypeTag,
class MyTypeTag>
123 using type = UndefinedProperty;
125template<
class TypeTag,
class MyTypeTag>
127 using type = UndefinedProperty;
129template<
class TypeTag>
131 using type = GetPropType<TypeTag, Scalar>;
132 static constexpr type value = 1e-2;
134template<
class TypeTag>
136 using type = GetPropType<TypeTag, Scalar>;
137 static constexpr type value = 1e-2;
139template<
class TypeTag>
141 using type = GetPropType<TypeTag, Scalar>;
142 static constexpr type value = 0.9;
144template<
class TypeTag>
146 static constexpr int value = 200;
148template<
class TypeTag>
150 static constexpr int value = 40;
152template<
class TypeTag>
153struct LinearSolverVerbosity<TypeTag, TTag::FlowIstlSolverParams> {
154 static constexpr int value = 0;
156template<
class TypeTag>
158 static constexpr int value = 0;
160template<
class TypeTag>
162 static constexpr auto value =
"ILU";
164template<
class TypeTag>
166 static constexpr bool value =
false;
168template<
class TypeTag>
170 static constexpr bool value =
false;
172template<
class TypeTag>
173struct UseGmres<TypeTag, TTag::FlowIstlSolverParams> {
174 static constexpr bool value =
false;
176template<
class TypeTag>
178 static constexpr bool value =
false;
180template<
class TypeTag>
181struct LinearSolverBackend<TypeTag, TTag::FlowIstlSolverParams> {
184template<
class TypeTag>
186 static constexpr bool value =
false;
188template<
class TypeTag>
190 static constexpr bool value =
false;
192template<
class TypeTag>
194 static constexpr int value = 4;
196template<
class TypeTag>
198 static constexpr int value = 30;
200template<
class TypeTag>
202 static constexpr auto value =
"ilu0";
204template<
class TypeTag>
206 static constexpr auto value =
"none";
208template<
class TypeTag>
210 static constexpr int value = 0;
212template<
class TypeTag>
214 static constexpr int value = 0;
216template<
class TypeTag>
218 static constexpr bool value =
true;
230 double linear_solver_reduction_;
231 double relaxed_linear_solver_reduction_;
232 double ilu_relaxation_;
233 int linear_solver_maxiter_;
234 int linear_solver_restart_;
235 int linear_solver_verbosity_;
236 int ilu_fillin_level_;
239 bool ilu_reorder_sphere_;
240 bool newton_use_gmres_;
241 bool require_full_sparsity_pattern_;
242 bool ignoreConvergenceFailure_;
243 bool scale_linear_system_;
244 std::string linsolver_;
245 std::string accelerator_mode_;
247 int opencl_platform_id_;
248 int cpr_reuse_setup_;
249 int cpr_reuse_interval_;
250 bool opencl_ilu_parallel_;
252 template <
class TypeTag>
253 void init(
bool cprRequestedInDataFile)
256 linear_solver_reduction_ = EWOMS_GET_PARAM(TypeTag,
double, LinearSolverReduction);
257 relaxed_linear_solver_reduction_ = EWOMS_GET_PARAM(TypeTag,
double, RelaxedLinearSolverReduction);
258 ilu_relaxation_ = EWOMS_GET_PARAM(TypeTag,
double, IluRelaxation);
259 linear_solver_maxiter_ = EWOMS_GET_PARAM(TypeTag,
int, LinearSolverMaxIter);
260 linear_solver_restart_ = EWOMS_GET_PARAM(TypeTag,
int, LinearSolverRestart);
261 linear_solver_verbosity_ = EWOMS_GET_PARAM(TypeTag,
int, LinearSolverVerbosity);
262 ilu_fillin_level_ = EWOMS_GET_PARAM(TypeTag,
int, IluFillinLevel);
263 ilu_milu_ = convertString2Milu(EWOMS_GET_PARAM(TypeTag, std::string, MiluVariant));
264 ilu_redblack_ = EWOMS_GET_PARAM(TypeTag,
bool, IluRedblack);
265 ilu_reorder_sphere_ = EWOMS_GET_PARAM(TypeTag,
bool, IluReorderSpheres);
266 newton_use_gmres_ = EWOMS_GET_PARAM(TypeTag,
bool, UseGmres);
267 ignoreConvergenceFailure_ = EWOMS_GET_PARAM(TypeTag,
bool, LinearSolverIgnoreConvergenceFailure);
268 scale_linear_system_ = EWOMS_GET_PARAM(TypeTag,
bool, ScaleLinearSystem);
269 cpr_reuse_setup_ = EWOMS_GET_PARAM(TypeTag,
int, CprReuseSetup);
270 cpr_reuse_interval_ = EWOMS_GET_PARAM(TypeTag,
int, CprReuseInterval);
272 if (!EWOMS_PARAM_IS_SET(TypeTag,
int, LinearSolverMaxIter) && cprRequestedInDataFile) {
275 linsolver_ = EWOMS_GET_PARAM(TypeTag, std::string, LinearSolver);
278 accelerator_mode_ = EWOMS_GET_PARAM(TypeTag, std::string, AcceleratorMode);
279 bda_device_id_ = EWOMS_GET_PARAM(TypeTag,
int, BdaDeviceId);
280 opencl_platform_id_ = EWOMS_GET_PARAM(TypeTag,
int, OpenclPlatformId);
281 opencl_ilu_parallel_ = EWOMS_GET_PARAM(TypeTag,
bool, OpenclIluParallel);
284 template <
class TypeTag>
285 static void registerParameters()
287 EWOMS_REGISTER_PARAM(TypeTag,
double, LinearSolverReduction,
"The minimum reduction of the residual which the linear solver must achieve for accepting solution as converged.");
288 EWOMS_REGISTER_PARAM(TypeTag,
double, RelaxedLinearSolverReduction,
"The minimum reduction of the residual which the linear solver need to achieve for the solution to be accepted");
289 EWOMS_REGISTER_PARAM(TypeTag,
double, IluRelaxation,
"The relaxation factor of the linear solver's ILU preconditioner");
290 EWOMS_REGISTER_PARAM(TypeTag,
int, LinearSolverMaxIter,
"The maximum number of iterations of the linear solver");
291 EWOMS_REGISTER_PARAM(TypeTag,
int, LinearSolverRestart,
"The number of iterations after which GMRES is restarted");
292 EWOMS_REGISTER_PARAM(TypeTag,
int, LinearSolverVerbosity,
"The verbosity level of the linear solver (0: off, 2: all)");
293 EWOMS_REGISTER_PARAM(TypeTag,
int, IluFillinLevel,
"The fill-in level of the linear solver's ILU preconditioner");
294 EWOMS_REGISTER_PARAM(TypeTag, std::string, MiluVariant,
"Specify which variant of the modified-ILU preconditioner ought to be used. Possible variants are: ILU (default, plain ILU), MILU_1 (lump diagonal with dropped row entries), MILU_2 (lump diagonal with the sum of the absolute values of the dropped row entries), MILU_3 (if diagonal is positive add sum of dropped row entrires. Otherwise subtract them), MILU_4 (if diagonal is positive add sum of dropped row entrires. Otherwise do nothing");
295 EWOMS_REGISTER_PARAM(TypeTag,
bool, IluRedblack,
"Use red-black partitioning for the ILU preconditioner");
296 EWOMS_REGISTER_PARAM(TypeTag,
bool, IluReorderSpheres,
"Whether to reorder the entries of the matrix in the red-black ILU preconditioner in spheres starting at an edge. If false the original ordering is preserved in each color. Otherwise why try to ensure D4 ordering (in a 2D structured grid, the diagonal elements are consecutive).");
297 EWOMS_REGISTER_PARAM(TypeTag,
bool, UseGmres,
"Use GMRES as the linear solver");
298 EWOMS_REGISTER_PARAM(TypeTag,
bool, LinearSolverIgnoreConvergenceFailure,
"Continue with the simulation like nothing happened after the linear solver did not converge");
299 EWOMS_REGISTER_PARAM(TypeTag,
bool, ScaleLinearSystem,
"Scale linear system according to equation scale and primary variable types");
300 EWOMS_REGISTER_PARAM(TypeTag,
int, CprReuseSetup,
"Reuse preconditioner setup. Valid options are 0: recreate the preconditioner for every linear solve, 1: recreate once every timestep, 2: recreate if last linear solve took more than 10 iterations, 3: never recreate, 4: recreated every CprReuseInterval");
301 EWOMS_REGISTER_PARAM(TypeTag,
int, CprReuseInterval,
"Reuse preconditioner interval. Used when CprReuseSetup is set to 4, then the preconditioner will be fully recreated instead of reused every N linear solve, where N is this parameter.");
302 EWOMS_REGISTER_PARAM(TypeTag, std::string, LinearSolver,
"Configuration of solver. Valid options are: ilu0 (default), cprw, cpr (an alias for cprw), cpr_quasiimpes, cpr_trueimpes or amg. Alternatively, you can request a configuration to be read from a JSON file by giving the filename here, ending with '.json.'");
303 EWOMS_REGISTER_PARAM(TypeTag, std::string, AcceleratorMode,
"Choose a linear solver, usage: '--accelerator-mode=[none|cusparse|opencl|amgcl|rocalution]'");
304 EWOMS_REGISTER_PARAM(TypeTag,
int, BdaDeviceId,
"Choose device ID for cusparseSolver or openclSolver, use 'nvidia-smi' or 'clinfo' to determine valid IDs");
305 EWOMS_REGISTER_PARAM(TypeTag,
int, OpenclPlatformId,
"Choose platform ID for openclSolver, use 'clinfo' to determine valid platform IDs");
306 EWOMS_REGISTER_PARAM(TypeTag,
bool, OpenclIluParallel,
"Parallelize ILU decomposition and application on GPU");
314 newton_use_gmres_ =
false;
315 linear_solver_reduction_ = 1e-2;
316 relaxed_linear_solver_reduction_ = 1e-2;
317 linear_solver_maxiter_ = 150;
318 linear_solver_restart_ = 40;
319 linear_solver_verbosity_ = 0;
320 require_full_sparsity_pattern_ =
false;
321 ignoreConvergenceFailure_ =
false;
322 ilu_fillin_level_ = 0;
323 ilu_relaxation_ = 0.9;
325 ilu_redblack_ =
false;
326 ilu_reorder_sphere_ =
true;
327 accelerator_mode_ =
"none";
329 opencl_platform_id_ = 0;
330 opencl_ilu_parallel_ =
true;
This class solves the fully implicit black-oil system by solving the reduced system (after eliminatin...
Definition: ISTLSolverEbos.hpp:200
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: BlackoilPhases.hpp:27
MILU_VARIANT
Definition: MILU.hpp:34
@ ILU
Do not perform modified ILU.
This class carries all parameters for the NewtonIterationBlackoilInterleaved class.
Definition: FlowLinearSolverParameters.hpp:229
Definition: FlowLinearSolverParameters.hpp:114
Definition: FlowLinearSolverParameters.hpp:118
Definition: FlowLinearSolverParameters.hpp:106
Definition: FlowLinearSolverParameters.hpp:102
Definition: FlowLinearSolverParameters.hpp:70
Definition: FlowLinearSolverParameters.hpp:78
Definition: FlowLinearSolverParameters.hpp:55
Definition: FlowLinearSolverParameters.hpp:82
Definition: FlowLinearSolverParameters.hpp:90
Definition: FlowLinearSolverParameters.hpp:59
Definition: FlowLinearSolverParameters.hpp:45
Definition: FlowLinearSolverParameters.hpp:63
Definition: FlowLinearSolverParameters.hpp:110
Definition: FlowLinearSolverParameters.hpp:74
Definition: FlowLinearSolverParameters.hpp:126
Definition: FlowLinearSolverParameters.hpp:94
Definition: FlowLinearSolverParameters.hpp:50
Definition: FlowLinearSolverParameters.hpp:98
Definition: FlowLinearSolverParameters.hpp:41
Definition: FlowLinearSolverParameters.hpp:86