21 #ifndef OPM_FLOWMAINPOLYMER_HEADER_INCLUDED 22 #define OPM_FLOWMAINPOLYMER_HEADER_INCLUDED 26 #include <opm/autodiff/FlowMain.hpp> 27 #include <opm/polymer/PolymerProperties.hpp> 28 #include <opm/polymer/fullyimplicit/PolymerPropsAd.hpp> 36 template <
class Gr
id,
class Simulator>
41 using Base::eclipse_state_;
43 using Base::fis_solver_;
44 using Base::parallel_information_;
48 std::unique_ptr<PolymerProperties> polymer_props_legacy_;
49 std::unique_ptr<PolymerPropsAd> polymer_props_;
55 void printStartupMessage()
57 if (Base::output_cout_) {
59 std::cout <<
"**********************************************************************\n";
61 std::cout <<
"* This is Flow-Polymer (version " << version <<
")" 62 << std::string(18 - version.size(),
' ') <<
"*\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";
68 std::cout <<
"**********************************************************************\n\n";
82 void setupGridAndProps()
84 Base::setupGridAndProps();
86 if (Base::deck_->hasKeyword(
"POLYMER")) {
87 polymer_props_legacy_.reset(
new PolymerProperties(*Base::deck_, *Base::eclipse_state_));
96 void setupLinearSolver()
98 const std::string cprSolver =
"cpr";
99 const std::string interleavedSolver =
"interleaved";
100 const std::string directSolver =
"direct";
101 const std::string flowDefaultSolver = interleavedSolver;
103 const Opm::SimulationConfig& simCfg = eclipse_state_->getSimulationConfig();
104 std::string solver_approach = flowDefaultSolver;
106 if (param_.has(
"solver_approach")) {
107 solver_approach = param_.template get<std::string>(
"solver_approach");
109 if (simCfg.useCPR()) {
110 solver_approach = cprSolver;
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) {
118 }
else if (solver_approach == directSolver) {
121 OPM_THROW( std::runtime_error ,
"Internal error - solver approach " << solver_approach <<
" not recognized.");
132 void createSimulator()
135 Base::simulator_.reset(
new Simulator(Base::param_,
136 Base::grid_init_->grid(),
140 Base::rock_comp_->isActive() ? Base::rock_comp_.
get() :
nullptr,
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_,
151 Base::threshold_pressures_));
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 "2015.10" (for a release branch) or "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
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: AdditionalObjectDeleter.hpp:22
Definition: FlowMainPolymer.hpp:37