All Classes Namespaces Files Functions Variables Typedefs Enumerator Pages
NewtonIterationBlackoilInterleaved.hpp
1 /*
2  Copyright 2015 SINTEF ICT, Applied Mathematics.
3  Copyright 2015 IRIS AS
4  Copyright 2015 Dr. Blatt - HPC-Simulation-Software & Services
5  Copyright 2015 NTNU
6  Copyright 2015 Statoil AS
7  Copyright 2015 IRIS AS
8 
9  This file is part of the Open Porous Media project (OPM).
10 
11  OPM is free software: you can redistribute it and/or modify
12  it under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OPM is distributed in the hope that it will be useful,
17  but WITHOUT ANY WARRANTY; without even the implied warranty of
18  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
19  GNU General Public License for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OPM. If not, see <http://www.gnu.org/licenses/>.
23 */
24 
25 #ifndef OPM_NEWTONITERATIONBLACKOILINTERLEAVED_HEADER_INCLUDED
26 #define OPM_NEWTONITERATIONBLACKOILINTERLEAVED_HEADER_INCLUDED
27 
28 #include <opm/autodiff/NewtonIterationBlackoilInterface.hpp>
29 #include <opm/core/utility/parameters/ParameterGroup.hpp>
30 
31 #include <array>
32 #include <memory>
33 
34 namespace Opm
35 {
38  {
39  double linear_solver_reduction_;
40  double ilu_relaxation_;
41  int linear_solver_maxiter_;
42  int linear_solver_restart_;
43  int linear_solver_verbosity_;
44  int ilu_fillin_level_;
45  bool newton_use_gmres_;
46  bool require_full_sparsity_pattern_;
47  bool ignoreConvergenceFailure_;
48  bool linear_solver_use_amg_;
49 
51  // read values from parameter class
52  NewtonIterationBlackoilInterleavedParameters( const ParameterGroup& param )
53  {
54  // set default parameters
55  reset();
56 
57  // read parameters (using previsouly set default values)
58  newton_use_gmres_ = param.getDefault("newton_use_gmres", newton_use_gmres_ );
59  linear_solver_reduction_ = param.getDefault("linear_solver_reduction", linear_solver_reduction_ );
60  linear_solver_maxiter_ = param.getDefault("linear_solver_maxiter", linear_solver_maxiter_);
61  linear_solver_restart_ = param.getDefault("linear_solver_restart", linear_solver_restart_);
62  linear_solver_verbosity_ = param.getDefault("linear_solver_verbosity", linear_solver_verbosity_);
63  require_full_sparsity_pattern_ = param.getDefault("require_full_sparsity_pattern", require_full_sparsity_pattern_);
64  ignoreConvergenceFailure_ = param.getDefault("linear_solver_ignoreconvergencefailure", ignoreConvergenceFailure_);
65  linear_solver_use_amg_ = param.getDefault("linear_solver_use_amg", linear_solver_use_amg_ );
66  ilu_relaxation_ = param.getDefault("ilu_relaxation", ilu_relaxation_ );
67  ilu_fillin_level_ = param.getDefault("ilu_fillin_level", ilu_fillin_level_ );
68  }
69 
70  // set default values
71  void reset()
72  {
73  newton_use_gmres_ = false;
74  linear_solver_reduction_ = 1e-2;
75  linear_solver_maxiter_ = 150;
76  linear_solver_restart_ = 40;
77  linear_solver_verbosity_ = 0;
78  require_full_sparsity_pattern_ = false;
79  ignoreConvergenceFailure_ = false;
80  linear_solver_use_amg_ = false;
81  ilu_fillin_level_ = 0;
82  ilu_relaxation_ = 0.9;
83  }
84  };
85 
86 
91  {
92  public:
93 
98  NewtonIterationBlackoilInterleaved(const ParameterGroup& param,
99  const boost::any& parallelInformation=boost::any());
100 
107 
109  virtual int iterations () const { return iterations_; }
110 
112  virtual const boost::any& parallelInformation() const;
113 
114  private:
115  // max number of equations supported, increase if necessary
116  static const int maxNumberEquations_ = 6 ;
117 
118  mutable std::array< std::unique_ptr< NewtonIterationBlackoilInterface >, maxNumberEquations_+1 > newtonIncrementDoublePrecision_;
119  mutable std::array< std::unique_ptr< NewtonIterationBlackoilInterface >, maxNumberEquations_+1 > newtonIncrementSinglePrecision_;
121  boost::any parallelInformation_;
122  mutable int iterations_;
123  };
124 
125 } // namespace Opm
126 
127 
128 #endif // OPM_NEWTONITERATIONBLACKOILINTERLEAVED_HEADER_INCLUDED
LinearisedBlackoilResidual::ADB::V SolutionVector
Return type for linearSolve(). A simple, non-ad vector type.
Definition: NewtonIterationBlackoilInterface.hpp:35
This class solves the fully implicit black-oil system by solving the reduced system (after eliminatin...
Definition: NewtonIterationBlackoilInterleaved.hpp:90
virtual int iterations() const
Definition: NewtonIterationBlackoilInterleaved.hpp:109
NewtonIterationBlackoilInterleaved(const ParameterGroup &param, const boost::any &parallelInformation=boost::any())
Construct a system solver.
Definition: NewtonIterationBlackoilInterleaved.cpp:303
Residual structure of the fully implicit solver.
Definition: LinearisedBlackoilResidual.hpp:47
virtual const boost::any & parallelInformation() const
Get the information about the parallelization of the grid.
Definition: NewtonIterationBlackoilInterleaved.cpp:490
Interface class for (linear) solvers for the fully implicit black-oil system.
Definition: NewtonIterationBlackoilInterface.hpp:31
This class carries all parameters for the NewtonIterationBlackoilInterleaved class.
Definition: NewtonIterationBlackoilInterleaved.hpp:37
virtual SolutionVector computeNewtonIncrement(const LinearisedBlackoilResidual &residual) const
Solve the system of linear equations Ax = b, with A being the combined derivative matrix of the resid...
Definition: NewtonIterationBlackoilInterleaved.cpp:469