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/common/utility/parameters/ParameterGroup.hpp>
28 #include <opm/simulators/linalg/ParallelOverlappingILU0.hpp>
29 
30 #include <opm/simulators/linalg/linalgproperties.hh>
31 #include <opm/models/utils/parametersystem.hh>
32 
33 #include <array>
34 #include <memory>
35 
36 namespace Opm {
37 template <class TypeTag>
38 class ISTLSolverEbos;
39 }
40 
41 
42 namespace Opm::Properties {
43 
44 namespace TTag {
46 }
47 
48 template<class TypeTag, class MyTypeTag>
50  using type = UndefinedProperty;
51 };
52 template<class TypeTag, class MyTypeTag>
53 struct IluRelaxation {
54  using type = UndefinedProperty;
55 };
56 template<class TypeTag, class MyTypeTag>
58  using type = UndefinedProperty;
59 };
60 template<class TypeTag, class MyTypeTag>
62  using type = UndefinedProperty;
63 };
64 template<class TypeTag, class MyTypeTag>
66  using type = UndefinedProperty;
67 };
68 template<class TypeTag, class MyTypeTag>
70  using type = UndefinedProperty;
71 };
72 template<class TypeTag, class MyTypeTag>
73 struct MiluVariant {
74  using type = UndefinedProperty;
75 };
76 template<class TypeTag, class MyTypeTag>
77 struct IluRedblack {
78  using type = UndefinedProperty;
79 };
80 template<class TypeTag, class MyTypeTag>
82  using type = UndefinedProperty;
83 };
84 template<class TypeTag, class MyTypeTag>
85 struct UseGmres {
86  using type = UndefinedProperty;
87 };
88 template<class TypeTag, class MyTypeTag>
90  using type = UndefinedProperty;
91 };
92 template<class TypeTag, class MyTypeTag>
94  using type = UndefinedProperty;
95 };
96 template<class TypeTag, class MyTypeTag>
98  using type = UndefinedProperty;
99 };
100 template<class TypeTag, class MyTypeTag>
102  using type = UndefinedProperty;
103 };
104 template<class TypeTag, class MyTypeTag>
106  using type = UndefinedProperty;
107 };
108 template<class TypeTag, class MyTypeTag>
110  using type = UndefinedProperty;
111 };
112 template<class TypeTag, class MyTypeTag>
114  using type = UndefinedProperty;
115 };
116 template<class TypeTag, class MyTypeTag>
117 struct Linsolver {
118  using type = UndefinedProperty;
119 };
120 template<class TypeTag, class MyTypeTag>
122  using type = UndefinedProperty;
123 };
124 template<class TypeTag, class MyTypeTag>
125 struct BdaDeviceId {
126  using type = UndefinedProperty;
127 };
128 template<class TypeTag, class MyTypeTag>
130  using type = UndefinedProperty;
131 };
132 template<class TypeTag, class MyTypeTag>
134  using type = UndefinedProperty;
135 };
136 template<class TypeTag, class MyTypeTag>
138  using type = UndefinedProperty;
139 };
140 
141 template<class TypeTag>
142 struct LinearSolverReduction<TypeTag, TTag::FlowIstlSolverParams> {
143  using type = GetPropType<TypeTag, Scalar>;
144  static constexpr type value = 1e-2;
145 };
146 template<class TypeTag>
147 struct IluRelaxation<TypeTag, TTag::FlowIstlSolverParams> {
148  using type = GetPropType<TypeTag, Scalar>;
149  static constexpr type value = 0.9;
150 };
151 template<class TypeTag>
152 struct LinearSolverMaxIter<TypeTag, TTag::FlowIstlSolverParams> {
153  static constexpr int value = 200;
154 };
155 template<class TypeTag>
156 struct LinearSolverRestart<TypeTag, TTag::FlowIstlSolverParams> {
157  static constexpr int value = 40;
158 };
159 template<class TypeTag>
160 struct FlowLinearSolverVerbosity<TypeTag, TTag::FlowIstlSolverParams> {
161  static constexpr int value = 0;
162 };
163 template<class TypeTag>
164 struct IluFillinLevel<TypeTag, TTag::FlowIstlSolverParams> {
165  static constexpr int value = 0;
166 };
167 template<class TypeTag>
168 struct MiluVariant<TypeTag, TTag::FlowIstlSolverParams> {
169  static constexpr auto value = "ILU";
170 };
171 template<class TypeTag>
172 struct IluRedblack<TypeTag, TTag::FlowIstlSolverParams> {
173  static constexpr bool value = false;
174 };
175 template<class TypeTag>
176 struct IluReorderSpheres<TypeTag, TTag::FlowIstlSolverParams> {
177  static constexpr bool value = false;
178 };
179 template<class TypeTag>
180 struct UseGmres<TypeTag, TTag::FlowIstlSolverParams> {
181  static constexpr bool value = false;
182 };
183 template<class TypeTag>
184 struct LinearSolverRequireFullSparsityPattern<TypeTag, TTag::FlowIstlSolverParams> {
185  static constexpr bool value = false;
186 };
187 template<class TypeTag>
188 struct LinearSolverIgnoreConvergenceFailure<TypeTag, TTag::FlowIstlSolverParams> {
189  static constexpr bool value = false;
190 };
191 template<class TypeTag>
192 struct LinearSolverBackend<TypeTag, TTag::FlowIstlSolverParams> {
194 };
195 template<class TypeTag>
196 struct PreconditionerAddWellContributions<TypeTag, TTag::FlowIstlSolverParams> {
197  static constexpr bool value = false;
198 };
199 template<class TypeTag>
200 struct ScaleLinearSystem<TypeTag, TTag::FlowIstlSolverParams> {
201  static constexpr bool value = false;
202 };
203 template<class TypeTag>
204 struct CprMaxEllIter<TypeTag, TTag::FlowIstlSolverParams> {
205  static constexpr int value = 20;
206 };
207 template<class TypeTag>
208 struct CprEllSolvetype<TypeTag, TTag::FlowIstlSolverParams> {
209  static constexpr int value = 0;
210 };
211 template<class TypeTag>
212 struct CprReuseSetup<TypeTag, TTag::FlowIstlSolverParams> {
213  static constexpr int value = 3;
214 };
215 template<class TypeTag>
216 struct Linsolver<TypeTag, TTag::FlowIstlSolverParams> {
217  static constexpr auto value = "ilu0";
218 };
219 template<class TypeTag>
220 struct AcceleratorMode<TypeTag, TTag::FlowIstlSolverParams> {
221  static constexpr auto value = "none";
222 };
223 template<class TypeTag>
224 struct BdaDeviceId<TypeTag, TTag::FlowIstlSolverParams> {
225  static constexpr int value = 0;
226 };
227 template<class TypeTag>
228 struct OpenclPlatformId<TypeTag, TTag::FlowIstlSolverParams> {
229  static constexpr int value = 0;
230 };
231 template<class TypeTag>
232 struct OpenclIluReorder<TypeTag, TTag::FlowIstlSolverParams> {
233  static constexpr auto value = ""; // note: default value is chosen depending on the solver used
234 };
235 template<class TypeTag>
236 struct FpgaBitstream<TypeTag, TTag::FlowIstlSolverParams> {
237  static constexpr auto value = "";
238 };
239 
240 } // namespace Opm::Properties
241 
242 namespace Opm
243 {
244 
245 
248  {
249  double linear_solver_reduction_;
250  double ilu_relaxation_;
251  int linear_solver_maxiter_;
252  int linear_solver_restart_;
253  int linear_solver_verbosity_;
254  int ilu_fillin_level_;
255  MILU_VARIANT ilu_milu_;
256  bool ilu_redblack_;
257  bool ilu_reorder_sphere_;
258  bool newton_use_gmres_;
259  bool require_full_sparsity_pattern_;
260  bool ignoreConvergenceFailure_;
261  bool scale_linear_system_;
262  std::string linsolver_;
263  std::string accelerator_mode_;
264  int bda_device_id_;
265  int opencl_platform_id_;
266  int cpr_max_ell_iter_ = 20;
267  int cpr_reuse_setup_ = 0;
268  std::string opencl_ilu_reorder_;
269  std::string fpga_bitstream_;
270 
271  template <class TypeTag>
272  void init()
273  {
274  // TODO: these parameters have undocumented non-trivial dependencies
275  linear_solver_reduction_ = EWOMS_GET_PARAM(TypeTag, double, LinearSolverReduction);
276  ilu_relaxation_ = EWOMS_GET_PARAM(TypeTag, double, IluRelaxation);
277  linear_solver_maxiter_ = EWOMS_GET_PARAM(TypeTag, int, LinearSolverMaxIter);
278  linear_solver_restart_ = EWOMS_GET_PARAM(TypeTag, int, LinearSolverRestart);
279  linear_solver_verbosity_ = EWOMS_GET_PARAM(TypeTag, int, FlowLinearSolverVerbosity);
280  ilu_fillin_level_ = EWOMS_GET_PARAM(TypeTag, int, IluFillinLevel);
281  ilu_milu_ = convertString2Milu(EWOMS_GET_PARAM(TypeTag, std::string, MiluVariant));
282  ilu_redblack_ = EWOMS_GET_PARAM(TypeTag, bool, IluRedblack);
283  ilu_reorder_sphere_ = EWOMS_GET_PARAM(TypeTag, bool, IluReorderSpheres);
284  newton_use_gmres_ = EWOMS_GET_PARAM(TypeTag, bool, UseGmres);
285  require_full_sparsity_pattern_ = EWOMS_GET_PARAM(TypeTag, bool, LinearSolverRequireFullSparsityPattern);
286  ignoreConvergenceFailure_ = EWOMS_GET_PARAM(TypeTag, bool, LinearSolverIgnoreConvergenceFailure);
287  scale_linear_system_ = EWOMS_GET_PARAM(TypeTag, bool, ScaleLinearSystem);
288  cpr_max_ell_iter_ = EWOMS_GET_PARAM(TypeTag, int, CprMaxEllIter);
289  cpr_reuse_setup_ = EWOMS_GET_PARAM(TypeTag, int, CprReuseSetup);
290  linsolver_ = EWOMS_GET_PARAM(TypeTag, std::string, Linsolver);
291  accelerator_mode_ = EWOMS_GET_PARAM(TypeTag, std::string, AcceleratorMode);
292  bda_device_id_ = EWOMS_GET_PARAM(TypeTag, int, BdaDeviceId);
293  opencl_platform_id_ = EWOMS_GET_PARAM(TypeTag, int, OpenclPlatformId);
294  opencl_ilu_reorder_ = EWOMS_GET_PARAM(TypeTag, std::string, OpenclIluReorder);
295  fpga_bitstream_ = EWOMS_GET_PARAM(TypeTag, std::string, FpgaBitstream);
296  }
297 
298  template <class TypeTag>
299  static void registerParameters()
300  {
301  EWOMS_REGISTER_PARAM(TypeTag, double, LinearSolverReduction, "The minimum reduction of the residual which the linear solver must achieve");
302  EWOMS_REGISTER_PARAM(TypeTag, double, IluRelaxation, "The relaxation factor of the linear solver's ILU preconditioner");
303  EWOMS_REGISTER_PARAM(TypeTag, int, LinearSolverMaxIter, "The maximum number of iterations of the linear solver");
304  EWOMS_REGISTER_PARAM(TypeTag, int, LinearSolverRestart, "The number of iterations after which GMRES is restarted");
305  EWOMS_REGISTER_PARAM(TypeTag, int, FlowLinearSolverVerbosity, "The verbosity level of the linear solver (0: off, 2: all)");
306  EWOMS_REGISTER_PARAM(TypeTag, int, IluFillinLevel, "The fill-in level of the linear solver's ILU preconditioner");
307  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");
308  EWOMS_REGISTER_PARAM(TypeTag, bool, IluRedblack, "Use red-black partitioning for the ILU preconditioner");
309  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).");
310  EWOMS_REGISTER_PARAM(TypeTag, bool, UseGmres, "Use GMRES as the linear solver");
311  EWOMS_REGISTER_PARAM(TypeTag, bool, LinearSolverRequireFullSparsityPattern, "Produce the full sparsity pattern for the linear solver");
312  EWOMS_REGISTER_PARAM(TypeTag, bool, LinearSolverIgnoreConvergenceFailure, "Continue with the simulation like nothing happened after the linear solver did not converge");
313  EWOMS_REGISTER_PARAM(TypeTag, bool, ScaleLinearSystem, "Scale linear system according to equation scale and primary variable types");
314  EWOMS_REGISTER_PARAM(TypeTag, int, CprMaxEllIter, "MaxIterations of the elliptic pressure part of the cpr solver");
315  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");
316  EWOMS_REGISTER_PARAM(TypeTag, std::string, Linsolver, "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.'");
317  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]'");
318  EWOMS_REGISTER_PARAM(TypeTag, int, BdaDeviceId, "Choose device ID for cusparseSolver or openclSolver, use 'nvidia-smi' or 'clinfo' to determine valid IDs");
319  EWOMS_REGISTER_PARAM(TypeTag, int, OpenclPlatformId, "Choose platform ID for openclSolver, use 'clinfo' to determine valid platform IDs");
320  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.");
321  EWOMS_REGISTER_PARAM(TypeTag, std::string, FpgaBitstream, "Specify the bitstream file for fpgaSolver (including path), usage: '--fpga-bitstream=<filename>'");
322  }
323 
324  FlowLinearSolverParameters() { reset(); }
325 
326  // set default values
327  void reset()
328  {
329  newton_use_gmres_ = false;
330  linear_solver_reduction_ = 1e-2;
331  linear_solver_maxiter_ = 150;
332  linear_solver_restart_ = 40;
333  linear_solver_verbosity_ = 0;
334  require_full_sparsity_pattern_ = false;
335  ignoreConvergenceFailure_ = false;
336  ilu_fillin_level_ = 0;
337  ilu_relaxation_ = 0.9;
338  ilu_milu_ = MILU_VARIANT::ILU;
339  ilu_redblack_ = false;
340  ilu_reorder_sphere_ = true;
341  accelerator_mode_ = "none";
342  bda_device_id_ = 0;
343  opencl_platform_id_ = 0;
344  opencl_ilu_reorder_ = ""; // note: the default value is chosen depending on the solver used
345  fpga_bitstream_ = "";
346  }
347  };
348 
349 
350 } // namespace Opm
351 
352 
353 
354 
355 #endif // OPM_FLOWLINEARSOLVERPARAMETERS_HEADER_INCLUDED
This class solves the fully implicit black-oil system by solving the reduced system (after eliminatin...
Definition: ISTLSolverEbos.hpp:78
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: BlackoilPhases.hpp:27
MILU_VARIANT
Definition: ParallelOverlappingILU0.hpp:47
@ ILU
Do not perform modified ILU.
This class carries all parameters for the NewtonIterationBlackoilInterleaved class.
Definition: FlowLinearSolverParameters.hpp:248
Definition: FlowLinearSolverParameters.hpp:121
Definition: FlowLinearSolverParameters.hpp:125
Definition: FlowLinearSolverParameters.hpp:109
Definition: FlowLinearSolverParameters.hpp:105
Definition: FlowLinearSolverParameters.hpp:113
Definition: FlowLinearSolverParameters.hpp:65
Definition: FlowLinearSolverParameters.hpp:137
Definition: FlowLinearSolverParameters.hpp:69
Definition: FlowLinearSolverParameters.hpp:77
Definition: FlowLinearSolverParameters.hpp:53
Definition: FlowLinearSolverParameters.hpp:81
Definition: FlowLinearSolverParameters.hpp:93
Definition: FlowLinearSolverParameters.hpp:57
Definition: FlowLinearSolverParameters.hpp:49
Definition: FlowLinearSolverParameters.hpp:89
Definition: FlowLinearSolverParameters.hpp:61
Definition: FlowLinearSolverParameters.hpp:117
Definition: FlowLinearSolverParameters.hpp:73
Definition: FlowLinearSolverParameters.hpp:133
Definition: FlowLinearSolverParameters.hpp:129
Definition: FlowLinearSolverParameters.hpp:97
Definition: FlowLinearSolverParameters.hpp:101
Definition: FlowLinearSolverParameters.hpp:45
Definition: FlowLinearSolverParameters.hpp:85