All Classes Namespaces Files Functions Variables Typedefs Enumerator Pages
FlowMainSequential.hpp
1 /*
2  Copyright 2016 SINTEF ICT, Applied Mathematics.
3 
4  This file is part of the Open Porous Media project (OPM).
5 
6  OPM is free software: you can redistribute it and/or modify
7  it under the terms of the GNU General Public License as published by
8  the Free Software Foundation, either version 3 of the License, or
9  (at your option) any later version.
10 
11  OPM is distributed in the hope that it will be useful,
12  but WITHOUT ANY WARRANTY; without even the implied warranty of
13  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14  GNU General Public License for more details.
15 
16  You should have received a copy of the GNU General Public License
17  along with OPM. If not, see <http://www.gnu.org/licenses/>.
18 */
19 
20 #ifndef OPM_FLOWMAINSEQUENTIAL_HEADER_INCLUDED
21 #define OPM_FLOWMAINSEQUENTIAL_HEADER_INCLUDED
22 
23 
24 
25 #include <opm/autodiff/FlowMain.hpp>
26 
27 
28 
29 namespace Opm
30 {
31 
32  // The FlowMainSequential class is for a black-oil simulator using the sequential models.
33  template <class Grid, class Simulator>
34  class FlowMainSequential : public FlowMainBase<FlowMainSequential<Grid, Simulator>, Grid, Simulator>
35  {
36  protected:
38  using Base::eclipse_state_;
39  using Base::param_;
40  using Base::fis_solver_;
41  using Base::parallel_information_;
42  friend Base;
43 
44  // ------------ Methods ------------
45 
46 
47  // Print startup message if on output rank.
48  void printStartupMessage()
49  {
50  if (Base::output_cout_) {
51  const std::string version = moduleVersionName();
52  std::cout << "**********************************************************************\n";
53  std::cout << "* *\n";
54  std::cout << "* This is Flow-Sequential (version " << version << ")"
55  << std::string(17 - version.size(), ' ') << "*\n";
56  std::cout << "* *\n";
57  std::cout << "* Flow-Sequential is a simulator for fully implicit three-phase, *\n";
58  std::cout << "* black-oil flow, and is part of OPM. *\n";
59  std::cout << "* For more information see http://opm-project.org *\n";
60  std::cout << "* *\n";
61  std::cout << "**********************************************************************\n\n";
62 
63  }
64  }
65 
66 
67 
68 
69 
70  // Setup linear solver.
71  // Writes to:
72  // fis_solver_
73  // param_ (conditionally)
74  // The CPR solver cannot be used with the sequential model.
75  // Also, the interleaved solver requires the full sparsity pattern option.
76  void setupLinearSolver()
77  {
78  const std::string cprSolver = "cpr";
79  const std::string interleavedSolver = "interleaved";
80  const std::string directSolver = "direct";
81  std::string flowDefaultSolver = interleavedSolver;
82 
83  if (!param_.has("solver_approach")) {
84  if (eclipse_state_->getSimulationConfig().useCPR()) {
85  flowDefaultSolver = cprSolver;
86  }
87  }
88 
89  const std::string solver_approach = param_.getDefault("solver_approach", flowDefaultSolver);
90 
91  if (solver_approach == cprSolver) {
92  OPM_THROW( std::runtime_error , "CPR solver is not ready for use with sequential simulator.");
93  } else if (solver_approach == interleavedSolver) {
94  if (!param_.has("require_full_sparsity_pattern")) {
95  param_.insertParameter("require_full_sparsity_pattern", "true");
96  }
97  fis_solver_.reset(new NewtonIterationBlackoilInterleaved(param_, parallel_information_));
98  } else if (solver_approach == directSolver) {
99  fis_solver_.reset(new NewtonIterationBlackoilSimple(param_, parallel_information_));
100  } else {
101  OPM_THROW( std::runtime_error , "Internal error - solver approach " << solver_approach << " not recognized.");
102  }
103  }
104 
105 
106 
107 
108 
109  // Create simulator instance.
110  // Writes to:
111  // simulator_
112  void createSimulator()
113  {
114  // We must override the min_iter argument unless it was already supplied, to avoid requiring iteration.
115  if (!param_.has("min_iter")) {
116  param_.insertParameter("min_iter", "0");
117  }
118 
119  // Create the simulator instance.
120  Base::simulator_.reset(new Simulator(Base::param_,
121  Base::grid_init_->grid(),
122  *Base::geoprops_,
123  *Base::fluidprops_,
124  Base::rock_comp_->isActive() ? Base::rock_comp_.get() : nullptr,
125  *Base::fis_solver_,
126  Base::gravity_.data(),
127  Base::deck_->hasKeyword("DISGAS"),
128  Base::deck_->hasKeyword("VAPOIL"),
129  Base::eclipse_state_,
130  *Base::output_writer_,
131  Base::threshold_pressures_));
132  }
133 
134 
135  };
136 
137 
138 } // namespace Opm
139 
140 
141 #endif // OPM_FLOWMAINSEQUENTIAL_HEADER_INCLUDED
Definition: FlowMainSequential.hpp:34
This class solves the fully implicit black-oil system by simply concatenating the Jacobian matrices a...
Definition: NewtonIterationBlackoilSimple.hpp:37
std::string moduleVersionName()
Return the version name of the module, for example &quot;2015.10&quot; (for a release branch) or &quot;2016...
Definition: moduleVersion.cpp:28
This class solves the fully implicit black-oil system by solving the reduced system (after eliminatin...
Definition: NewtonIterationBlackoilInterleaved.hpp:90
This class encapsulates the setup and running of a simulator based on an input deck.
Definition: FlowMain.hpp:122