All Classes Namespaces Files Functions Variables Typedefs Enumerator Pages
FlowMainPolymer.hpp
1 /*
2  Copyright 2014, 2015 STATOIL ASA.
3  Copyright 2015 SINTEF ICT, Applied Mathematics.
4 
5  This file is part of the Open Porous Media project (OPM).
6 
7  OPM is free software: you can redistribute it and/or modify
8  it under the terms of the GNU General Public License as published by
9  the Free Software Foundation, either version 3 of the License, or
10  (at your option) any later version.
11 
12  OPM is distributed in the hope that it will be useful,
13  but WITHOUT ANY WARRANTY; without even the implied warranty of
14  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15  GNU General Public License for more details.
16 
17  You should have received a copy of the GNU General Public License
18  along with OPM. If not, see <http://www.gnu.org/licenses/>.
19 */
20 
21 #ifndef OPM_FLOWMAINPOLYMER_HEADER_INCLUDED
22 #define OPM_FLOWMAINPOLYMER_HEADER_INCLUDED
23 
24 
25 
26 #include <opm/autodiff/FlowMain.hpp>
27 #include <opm/polymer/PolymerProperties.hpp>
28 #include <opm/polymer/fullyimplicit/PolymerPropsAd.hpp>
29 
30 
31 
32 namespace Opm
33 {
34 
35  // The FlowMainPolymer class is for a black-oil simulator with polymer.
36  template <class Grid, class Simulator>
37  class FlowMainPolymer : public FlowMainBase<FlowMainPolymer<Grid, Simulator>, Grid, Simulator>
38  {
39  protected:
40  using Base = FlowMainBase<FlowMainPolymer<Grid, Simulator>, Grid, Simulator>;
41  using Base::eclipse_state_;
42  using Base::param_;
43  using Base::fis_solver_;
44  using Base::parallel_information_;
45  friend Base;
46 
47  // Set in setupGridAndProps()
48  std::unique_ptr<PolymerProperties> polymer_props_legacy_; // Held by reference in polymer_props_
49  std::unique_ptr<PolymerPropsAd> polymer_props_;
50 
51  // ------------ Methods ------------
52 
53 
54  // Print startup message if on output rank.
55  void printStartupMessage()
56  {
57  if (Base::output_cout_) {
58  const std::string version = moduleVersionName();
59  std::cout << "**********************************************************************\n";
60  std::cout << "* *\n";
61  std::cout << "* This is Flow-Polymer (version " << version << ")"
62  << std::string(18 - version.size(), ' ') << "*\n";
63  std::cout << "* *\n";
64  std::cout << "* Flow-Polymer is a simulator for fully implicit three-phase, *\n";
65  std::cout << "* four-component (black-oil + polymer) flow, and is part of OPM. *\n";
66  std::cout << "* For more information see http://opm-project.org *\n";
67  std::cout << "* *\n";
68  std::cout << "**********************************************************************\n\n";
69 
70  }
71  }
72 
73 
74 
75 
76 
77  // Set up grid and property objects, by calling base class
78  // version and then creating polymer property objects.
79  // Writes to:
80  // polymer_props_legacy_
81  // polymer_props_
82  void setupGridAndProps()
83  {
84  Base::setupGridAndProps();
85 
86  if (Base::deck_->hasKeyword("POLYMER")) {
87  polymer_props_legacy_.reset(new PolymerProperties(*Base::deck_, *Base::eclipse_state_));
88  polymer_props_.reset(new PolymerPropsAd(*polymer_props_legacy_));
89  }
90  }
91 
92  // Setup linear solver.
93  // Writes to:
94  // fis_solver_
95  // Currently, the CPR solver is not ready for polymer solver yet
96  void setupLinearSolver()
97  {
98  const std::string cprSolver = "cpr";
99  const std::string interleavedSolver = "interleaved";
100  const std::string directSolver = "direct";
101  const std::string flowDefaultSolver = interleavedSolver;
102 
103  const Opm::SimulationConfig& simCfg = eclipse_state_->getSimulationConfig();
104  std::string solver_approach = flowDefaultSolver;
105 
106  if (param_.has("solver_approach")) {
107  solver_approach = param_.template get<std::string>("solver_approach");
108  } else {
109  if (simCfg.useCPR()) {
110  solver_approach = cprSolver;
111  }
112  }
113 
114  if (solver_approach == cprSolver) {
115  OPM_THROW( std::runtime_error , "CPR solver is not ready for use with polymer solver yet.");
116  } else if (solver_approach == interleavedSolver) {
117  fis_solver_.reset(new NewtonIterationBlackoilInterleaved(param_, parallel_information_));
118  } else if (solver_approach == directSolver) {
119  fis_solver_.reset(new NewtonIterationBlackoilSimple(param_, parallel_information_));
120  } else {
121  OPM_THROW( std::runtime_error , "Internal error - solver approach " << solver_approach << " not recognized.");
122  }
123  }
124 
125 
126 
127 
128 
129  // Create simulator instance.
130  // Writes to:
131  // simulator_
132  void createSimulator()
133  {
134  // Create the simulator instance.
135  Base::simulator_.reset(new Simulator(Base::param_,
136  Base::grid_init_->grid(),
137  *Base::geoprops_,
138  *Base::fluidprops_,
139  *polymer_props_,
140  Base::rock_comp_->isActive() ? Base::rock_comp_.get() : nullptr,
141  *Base::fis_solver_,
142  Base::gravity_.data(),
143  Base::deck_->hasKeyword("DISGAS"),
144  Base::deck_->hasKeyword("VAPOIL"),
145  Base::deck_->hasKeyword("POLYMER"),
146  Base::deck_->hasKeyword("PLYSHLOG"),
147  Base::deck_->hasKeyword("SHRATE"),
148  Base::eclipse_state_,
149  *Base::output_writer_,
150  Base::deck_,
151  Base::threshold_pressures_));
152  }
153 
154 
155  };
156 
157 
158 } // namespace Opm
159 
160 
161 
162 
163 #endif // OPM_FLOWMAINPOLYMER_HEADER_INCLUDED
This class solves the fully implicit black-oil system by simply concatenating the Jacobian matrices a...
Definition: NewtonIterationBlackoilSimple.hpp:37
Definition: PolymerPropsAd.hpp:32
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
Definition: PolymerProperties.hpp:43
Definition: FlowMainPolymer.hpp:37