My Project
FlowLinearSolverParameters.hpp
1 /*
2  Copyright 2015, 2020 SINTEF Digital, Mathematics and Cybernetics.
3  Copyright 2015 IRIS AS
4  Copyright 2015 Dr. Blatt - HPC-Simulation-Software & Services
5  Copyright 2015 NTNU
6  Copyright 2015 Statoil AS
7 
8  This file is part of the Open Porous Media project (OPM).
9 
10  OPM is free software: you can redistribute it and/or modify
11  it under the terms of the GNU General Public License as published by
12  the Free Software Foundation, either version 3 of the License, or
13  (at your option) any later version.
14 
15  OPM is distributed in the hope that it will be useful,
16  but WITHOUT ANY WARRANTY; without even the implied warranty of
17  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
18  GNU General Public License for more details.
19 
20  You should have received a copy of the GNU General Public License
21  along with OPM. If not, see <http://www.gnu.org/licenses/>.
22 */
23 
24 #ifndef OPM_FLOWLINEARSOLVERPARAMETERS_HEADER_INCLUDED
25 #define OPM_FLOWLINEARSOLVERPARAMETERS_HEADER_INCLUDED
26 
27 #include <opm/simulators/linalg/MILU.hpp>
28 
29 #include <opm/simulators/linalg/linalgproperties.hh>
30 #include <opm/models/utils/parametersystem.hh>
31 
32 namespace Opm {
33 template <class TypeTag>
34 class ISTLSolverEbos;
35 }
36 
37 
38 namespace Opm::Properties {
39 
40 namespace TTag {
42 }
43 
44 template<class TypeTag, class MyTypeTag>
46  using type = UndefinedProperty;
47 };
48 template<class TypeTag, class MyTypeTag>
49 struct IluRelaxation {
50  using type = UndefinedProperty;
51 };
52 template<class TypeTag, class MyTypeTag>
54  using type = UndefinedProperty;
55 };
56 template<class TypeTag, class MyTypeTag>
58  using type = UndefinedProperty;
59 };
60 //
61 // LinearSolverVerbosity defined in opm-models
62 //
63 template<class TypeTag, class MyTypeTag>
65  using type = UndefinedProperty;
66 };
67 template<class TypeTag, class MyTypeTag>
68 struct MiluVariant {
69  using type = UndefinedProperty;
70 };
71 template<class TypeTag, class MyTypeTag>
72 struct IluRedblack {
73  using type = UndefinedProperty;
74 };
75 template<class TypeTag, class MyTypeTag>
77  using type = UndefinedProperty;
78 };
79 template<class TypeTag, class MyTypeTag>
80 struct UseGmres {
81  using type = UndefinedProperty;
82 };
83 template<class TypeTag, class MyTypeTag>
85  using type = UndefinedProperty;
86 };
87 template<class TypeTag, class MyTypeTag>
89  using type = UndefinedProperty;
90 };
91 template<class TypeTag, class MyTypeTag>
93  using type = UndefinedProperty;
94 };
95 template<class TypeTag, class MyTypeTag>
96 struct CprMaxEllIter {
97  using type = UndefinedProperty;
98 };
99 template<class TypeTag, class MyTypeTag>
101  using type = UndefinedProperty;
102 };
103 template<class TypeTag, class MyTypeTag>
105  using type = UndefinedProperty;
106 };
107 template<class TypeTag, class MyTypeTag>
109  using type = UndefinedProperty;
110 };
111 template<class TypeTag, class MyTypeTag>
112 struct LinearSolver {
113  using type = UndefinedProperty;
114 };
115 template<class TypeTag, class MyTypeTag>
117  using type = UndefinedProperty;
118 };
119 template<class TypeTag, class MyTypeTag>
120 struct BdaDeviceId {
121  using type = UndefinedProperty;
122 };
123 template<class TypeTag, class MyTypeTag>
125  using type = UndefinedProperty;
126 };
127 template<class TypeTag, class MyTypeTag>
129  using type = UndefinedProperty;
130 };
131 template<class TypeTag, class MyTypeTag>
133  using type = UndefinedProperty;
134 };
135 template<class TypeTag>
136 struct LinearSolverReduction<TypeTag, TTag::FlowIstlSolverParams> {
137  using type = GetPropType<TypeTag, Scalar>;
138  static constexpr type value = 1e-2;
139 };
140 template<class TypeTag>
141 struct IluRelaxation<TypeTag, TTag::FlowIstlSolverParams> {
142  using type = GetPropType<TypeTag, Scalar>;
143  static constexpr type value = 0.9;
144 };
145 template<class TypeTag>
146 struct LinearSolverMaxIter<TypeTag, TTag::FlowIstlSolverParams> {
147  static constexpr int value = 200;
148 };
149 template<class TypeTag>
150 struct LinearSolverRestart<TypeTag, TTag::FlowIstlSolverParams> {
151  static constexpr int value = 40;
152 };
153 template<class TypeTag>
154 struct LinearSolverVerbosity<TypeTag, TTag::FlowIstlSolverParams> {
155  static constexpr int value = 0;
156 };
157 template<class TypeTag>
158 struct IluFillinLevel<TypeTag, TTag::FlowIstlSolverParams> {
159  static constexpr int value = 0;
160 };
161 template<class TypeTag>
162 struct MiluVariant<TypeTag, TTag::FlowIstlSolverParams> {
163  static constexpr auto value = "ILU";
164 };
165 template<class TypeTag>
166 struct IluRedblack<TypeTag, TTag::FlowIstlSolverParams> {
167  static constexpr bool value = false;
168 };
169 template<class TypeTag>
170 struct IluReorderSpheres<TypeTag, TTag::FlowIstlSolverParams> {
171  static constexpr bool value = false;
172 };
173 template<class TypeTag>
174 struct UseGmres<TypeTag, TTag::FlowIstlSolverParams> {
175  static constexpr bool value = false;
176 };
177 template<class TypeTag>
178 struct LinearSolverIgnoreConvergenceFailure<TypeTag, TTag::FlowIstlSolverParams> {
179  static constexpr bool value = false;
180 };
181 template<class TypeTag>
182 struct LinearSolverBackend<TypeTag, TTag::FlowIstlSolverParams> {
184 };
185 template<class TypeTag>
186 struct PreconditionerAddWellContributions<TypeTag, TTag::FlowIstlSolverParams> {
187  static constexpr bool value = false;
188 };
189 template<class TypeTag>
190 struct ScaleLinearSystem<TypeTag, TTag::FlowIstlSolverParams> {
191  static constexpr bool value = false;
192 };
193 template<class TypeTag>
194 struct CprMaxEllIter<TypeTag, TTag::FlowIstlSolverParams> {
195  static constexpr int value = 20;
196 };
197 template<class TypeTag>
198 struct CprEllSolvetype<TypeTag, TTag::FlowIstlSolverParams> {
199  static constexpr int value = 0;
200 };
201 template<class TypeTag>
202 struct CprReuseSetup<TypeTag, TTag::FlowIstlSolverParams> {
203  static constexpr int value = 3;
204 };
205 template<class TypeTag>
206 struct CprReuseInterval<TypeTag, TTag::FlowIstlSolverParams> {
207  static constexpr int value = 10;
208 };
209 template<class TypeTag>
210 struct LinearSolver<TypeTag, TTag::FlowIstlSolverParams> {
211  static constexpr auto value = "ilu0";
212 };
213 template<class TypeTag>
214 struct AcceleratorMode<TypeTag, TTag::FlowIstlSolverParams> {
215  static constexpr auto value = "none";
216 };
217 template<class TypeTag>
218 struct BdaDeviceId<TypeTag, TTag::FlowIstlSolverParams> {
219  static constexpr int value = 0;
220 };
221 template<class TypeTag>
222 struct OpenclPlatformId<TypeTag, TTag::FlowIstlSolverParams> {
223  static constexpr int value = 0;
224 };
225 template<class TypeTag>
226 struct OpenclIluReorder<TypeTag, TTag::FlowIstlSolverParams> {
227  static constexpr auto value = ""; // note: default value is chosen depending on the solver used
228 };
229 template<class TypeTag>
230 struct FpgaBitstream<TypeTag, TTag::FlowIstlSolverParams> {
231  static constexpr auto value = "";
232 };
233 
234 } // namespace Opm::Properties
235 
236 namespace Opm
237 {
238 
239 
242  {
243  double linear_solver_reduction_;
244  double ilu_relaxation_;
245  int linear_solver_maxiter_;
246  int linear_solver_restart_;
247  int linear_solver_verbosity_;
248  int ilu_fillin_level_;
249  MILU_VARIANT ilu_milu_;
250  bool ilu_redblack_;
251  bool ilu_reorder_sphere_;
252  bool newton_use_gmres_;
253  bool require_full_sparsity_pattern_;
254  bool ignoreConvergenceFailure_;
255  bool scale_linear_system_;
256  std::string linsolver_;
257  std::string accelerator_mode_;
258  int bda_device_id_;
259  int opencl_platform_id_;
260  int cpr_max_ell_iter_;
261  int cpr_reuse_setup_;
262  int cpr_reuse_interval_;
263  std::string opencl_ilu_reorder_;
264  std::string fpga_bitstream_;
265 
266  template <class TypeTag>
267  void init()
268  {
269  // TODO: these parameters have undocumented non-trivial dependencies
270  linear_solver_reduction_ = EWOMS_GET_PARAM(TypeTag, double, LinearSolverReduction);
271  ilu_relaxation_ = EWOMS_GET_PARAM(TypeTag, double, IluRelaxation);
272  linear_solver_maxiter_ = EWOMS_GET_PARAM(TypeTag, int, LinearSolverMaxIter);
273  linear_solver_restart_ = EWOMS_GET_PARAM(TypeTag, int, LinearSolverRestart);
274  linear_solver_verbosity_ = EWOMS_GET_PARAM(TypeTag, int, LinearSolverVerbosity);
275  ilu_fillin_level_ = EWOMS_GET_PARAM(TypeTag, int, IluFillinLevel);
276  ilu_milu_ = convertString2Milu(EWOMS_GET_PARAM(TypeTag, std::string, MiluVariant));
277  ilu_redblack_ = EWOMS_GET_PARAM(TypeTag, bool, IluRedblack);
278  ilu_reorder_sphere_ = EWOMS_GET_PARAM(TypeTag, bool, IluReorderSpheres);
279  newton_use_gmres_ = EWOMS_GET_PARAM(TypeTag, bool, UseGmres);
280  ignoreConvergenceFailure_ = EWOMS_GET_PARAM(TypeTag, bool, LinearSolverIgnoreConvergenceFailure);
281  scale_linear_system_ = EWOMS_GET_PARAM(TypeTag, bool, ScaleLinearSystem);
282  cpr_max_ell_iter_ = EWOMS_GET_PARAM(TypeTag, int, CprMaxEllIter);
283  cpr_reuse_setup_ = EWOMS_GET_PARAM(TypeTag, int, CprReuseSetup);
284  cpr_reuse_interval_ = EWOMS_GET_PARAM(TypeTag, int, CprReuseInterval);
285  linsolver_ = EWOMS_GET_PARAM(TypeTag, std::string, LinearSolver);
286  accelerator_mode_ = EWOMS_GET_PARAM(TypeTag, std::string, AcceleratorMode);
287  bda_device_id_ = EWOMS_GET_PARAM(TypeTag, int, BdaDeviceId);
288  opencl_platform_id_ = EWOMS_GET_PARAM(TypeTag, int, OpenclPlatformId);
289  opencl_ilu_reorder_ = EWOMS_GET_PARAM(TypeTag, std::string, OpenclIluReorder);
290  fpga_bitstream_ = EWOMS_GET_PARAM(TypeTag, std::string, FpgaBitstream);
291  }
292 
293  template <class TypeTag>
294  static void registerParameters()
295  {
296  EWOMS_REGISTER_PARAM(TypeTag, double, LinearSolverReduction, "The minimum reduction of the residual which the linear solver must achieve");
297  EWOMS_REGISTER_PARAM(TypeTag, double, IluRelaxation, "The relaxation factor of the linear solver's ILU preconditioner");
298  EWOMS_REGISTER_PARAM(TypeTag, int, LinearSolverMaxIter, "The maximum number of iterations of the linear solver");
299  EWOMS_REGISTER_PARAM(TypeTag, int, LinearSolverRestart, "The number of iterations after which GMRES is restarted");
300  EWOMS_REGISTER_PARAM(TypeTag, int, LinearSolverVerbosity, "The verbosity level of the linear solver (0: off, 2: all)");
301  EWOMS_REGISTER_PARAM(TypeTag, int, IluFillinLevel, "The fill-in level of the linear solver's ILU preconditioner");
302  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");
303  EWOMS_REGISTER_PARAM(TypeTag, bool, IluRedblack, "Use red-black partitioning for the ILU preconditioner");
304  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).");
305  EWOMS_REGISTER_PARAM(TypeTag, bool, UseGmres, "Use GMRES as the linear solver");
306  EWOMS_REGISTER_PARAM(TypeTag, bool, LinearSolverIgnoreConvergenceFailure, "Continue with the simulation like nothing happened after the linear solver did not converge");
307  EWOMS_REGISTER_PARAM(TypeTag, bool, ScaleLinearSystem, "Scale linear system according to equation scale and primary variable types");
308  EWOMS_REGISTER_PARAM(TypeTag, int, CprMaxEllIter, "MaxIterations of the elliptic pressure part of the cpr solver");
309  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");
310  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.");
311  EWOMS_REGISTER_PARAM(TypeTag, std::string, LinearSolver, "Configuration of solver. Valid options are: ilu0 (default), cpr (an alias for cpr_trueimpes), 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.'");
312  EWOMS_REGISTER_PARAM(TypeTag, std::string, AcceleratorMode, "Use GPU (cusparseSolver or openclSolver) or FPGA (fpgaSolver) as the linear solver, usage: '--accelerator-mode=[none|cusparse|opencl|fpga|amgcl]'");
313  EWOMS_REGISTER_PARAM(TypeTag, int, BdaDeviceId, "Choose device ID for cusparseSolver or openclSolver, use 'nvidia-smi' or 'clinfo' to determine valid IDs");
314  EWOMS_REGISTER_PARAM(TypeTag, int, OpenclPlatformId, "Choose platform ID for openclSolver, use 'clinfo' to determine valid platform IDs");
315  EWOMS_REGISTER_PARAM(TypeTag, std::string, OpenclIluReorder, "Choose the reordering strategy for ILU for openclSolver and fpgaSolver, usage: '--opencl-ilu-reorder=[level_scheduling|graph_coloring], level_scheduling behaves like Dune and cusparse, graph_coloring is more aggressive and likely to be faster, but is random-based and generally increases the number of linear solves and linear iterations significantly.");
316  EWOMS_REGISTER_PARAM(TypeTag, std::string, FpgaBitstream, "Specify the bitstream file for fpgaSolver (including path), usage: '--fpga-bitstream=<filename>'");
317  }
318 
319  FlowLinearSolverParameters() { reset(); }
320 
321  // set default values
322  void reset()
323  {
324  newton_use_gmres_ = false;
325  linear_solver_reduction_ = 1e-2;
326  linear_solver_maxiter_ = 150;
327  linear_solver_restart_ = 40;
328  linear_solver_verbosity_ = 0;
329  require_full_sparsity_pattern_ = false;
330  ignoreConvergenceFailure_ = false;
331  ilu_fillin_level_ = 0;
332  ilu_relaxation_ = 0.9;
333  ilu_milu_ = MILU_VARIANT::ILU;
334  ilu_redblack_ = false;
335  ilu_reorder_sphere_ = true;
336  accelerator_mode_ = "none";
337  bda_device_id_ = 0;
338  opencl_platform_id_ = 0;
339  opencl_ilu_reorder_ = ""; // note: the default value is chosen depending on the solver used
340  fpga_bitstream_ = "";
341  }
342  };
343 
344 
345 } // namespace Opm
346 
347 
348 
349 
350 #endif // OPM_FLOWLINEARSOLVERPARAMETERS_HEADER_INCLUDED
This class solves the fully implicit black-oil system by solving the reduced system (after eliminatin...
Definition: ISTLSolverEbos.hpp:197
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:242
Definition: FlowLinearSolverParameters.hpp:116
Definition: FlowLinearSolverParameters.hpp:120
Definition: FlowLinearSolverParameters.hpp:100
Definition: FlowLinearSolverParameters.hpp:96
Definition: FlowLinearSolverParameters.hpp:108
Definition: FlowLinearSolverParameters.hpp:104
Definition: FlowLinearSolverParameters.hpp:132
Definition: FlowLinearSolverParameters.hpp:64
Definition: FlowLinearSolverParameters.hpp:72
Definition: FlowLinearSolverParameters.hpp:49
Definition: FlowLinearSolverParameters.hpp:76
Definition: FlowLinearSolverParameters.hpp:84
Definition: FlowLinearSolverParameters.hpp:53
Definition: FlowLinearSolverParameters.hpp:45
Definition: FlowLinearSolverParameters.hpp:57
Definition: FlowLinearSolverParameters.hpp:112
Definition: FlowLinearSolverParameters.hpp:68
Definition: FlowLinearSolverParameters.hpp:128
Definition: FlowLinearSolverParameters.hpp:124
Definition: FlowLinearSolverParameters.hpp:88
Definition: FlowLinearSolverParameters.hpp:92
Definition: FlowLinearSolverParameters.hpp:41
Definition: FlowLinearSolverParameters.hpp:80