SimulatorFullyImplicitBlackoilPolymer.hpp
1 /*
2  Copyright 2013 SINTEF ICT, Applied Mathematics.
3  Copyright 2014 STATOIL ASA.
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_SIMULATORFULLYIMPLICITBLACKOILPOLYMER_HEADER_INCLUDED
22 #define OPM_SIMULATORFULLYIMPLICITBLACKOILPOLYMER_HEADER_INCLUDED
23 
24 #include <opm/autodiff/SimulatorBase.hpp>
25 #include <opm/autodiff/SimulatorFullyImplicitBlackoilOutput.hpp>
26 #include <opm/polymer/fullyimplicit/BlackoilPolymerModel.hpp>
27 #include <opm/polymer/fullyimplicit/WellStateFullyImplicitBlackoilPolymer.hpp>
28 #include <opm/polymer/PolymerBlackoilState.hpp>
29 #include <opm/polymer/PolymerInflow.hpp>
30 
31 #include <opm/core/utility/parameters/ParameterGroup.hpp>
32 #include <opm/common/ErrorMacros.hpp>
33 
34 #include <opm/autodiff/GeoProps.hpp>
35 #include <opm/autodiff/BlackoilPropsAdFromDeck.hpp>
37 #include <opm/autodiff/NonlinearSolver.hpp>
38 
39 #include <opm/core/grid.h>
40 #include <opm/core/wells.h>
41 #include <opm/core/well_controls.h>
42 #include <opm/core/pressure/flow_bc.h>
43 
44 #include <opm/core/simulator/SimulatorReport.hpp>
45 #include <opm/simulators/timestepping/SimulatorTimer.hpp>
46 //#include <opm/core/simulator/AdaptiveSimulatorTimer.hpp>
47 #include <opm/core/utility/StopWatch.hpp>
48 #include <opm/core/utility/miscUtilities.hpp>
49 #include <opm/core/utility/miscUtilitiesBlackoil.hpp>
50 
51 #include <opm/core/props/rock/RockCompressibility.hpp>
52 
53 //#include <opm/core/simulator/AdaptiveTimeStepping.hpp>
54 #include <opm/core/transport/reorder/TransportSolverCompressibleTwophaseReorder.hpp>
55 
56 #include <opm/parser/eclipse/EclipseState/Schedule/Schedule.hpp>
57 #include <opm/parser/eclipse/EclipseState/Schedule/ScheduleEnums.hpp>
58 #include <opm/parser/eclipse/EclipseState/Schedule/Well.hpp>
59 #include <opm/parser/eclipse/EclipseState/Schedule/WellProductionProperties.hpp>
60 #include <opm/parser/eclipse/Deck/Deck.hpp>
61 
62 #include <boost/filesystem.hpp>
63 #include <boost/lexical_cast.hpp>
64 
65 #include <algorithm>
66 #include <cstddef>
67 #include <cassert>
68 #include <functional>
69 #include <memory>
70 #include <numeric>
71 #include <fstream>
72 #include <iostream>
73 #include <string>
74 #include <unordered_map>
75 #include <utility>
76 #include <vector>
77 
78 namespace Opm
79 {
80  template <class GridT>
82 
83  class StandardWells;
84 
85  template<class GridT>
87  {
91  typedef GridT Grid;
94  typedef StandardWells WellModel;
95  };
96 
99  template <class GridT>
101  : public SimulatorBase<SimulatorFullyImplicitBlackoilPolymer<GridT> >
102  {
104  typedef SimulatorBase<ThisType> BaseType;
105 
106  typedef SimulatorTraits<ThisType> Traits;
107  typedef typename Traits::Solver Solver;
108  typedef typename Traits::WellModel WellModel;
109 
110  public:
111  SimulatorFullyImplicitBlackoilPolymer(const ParameterGroup& param,
112  const GridT& grid,
113  DerivedGeology& geo,
115  const PolymerPropsAd& polymer_props,
116  const RockCompressibility* rock_comp_props,
118  const double* gravity,
119  const bool disgas,
120  const bool vapoil,
121  const bool polymer,
122  const bool plyshlog,
123  const bool shrate,
124  std::shared_ptr<EclipseState> eclipse_state,
125  BlackoilOutputWriter& output_writer,
126  std::shared_ptr< Deck > deck,
127  const std::vector<double>& threshold_pressures_by_face);
128 
129  std::unique_ptr<Solver> createSolver(const WellModel& well_model);
130 
131 
132  void handleAdditionalWellInflow(SimulatorTimer& timer,
133  WellsManager& wells_manager,
134  typename BaseType::WellState& well_state,
135  const Wells* wells);
136 
137 
138  private:
139  const PolymerPropsAd& polymer_props_;
140  bool has_polymer_;
141  // flag for PLYSHLOG keyword
142  bool has_plyshlog_;
143  // flag for SHRATE keyword
144  bool has_shrate_;
145  std::shared_ptr< Deck > deck_;
146 
147  std::vector<double> wells_rep_radius_;
148  std::vector<double> wells_perf_length_;
149  std::vector<double> wells_bore_diameter_;
150 
151  // generate the mapping from Cartesian grid cells to global compressed cells,
152  // copied from opm-core, to be used in function computeRepRadiusPerfLength()
153  static void
154  setupCompressedToCartesian(const int* global_cell, int number_of_cells, std::map<int,int>& cartesian_to_compressed);
155 
156  // calculate the representative radius and length for for well peforations
157  // and store the wellbore diameters
158  // it will be used in the shear-thinning calcluation only.
159  void
160  computeRepRadiusPerfLength(const EclipseState& eclipseState,
161  const size_t timeStep,
162  const GridT& grid,
163  std::vector<double>& wells_rep_radius,
164  std::vector<double>& wells_perf_length,
165  std::vector<double>& wells_bore_diameter);
166 
167  };
168 
169 } // namespace Opm
170 
171 #include "SimulatorFullyImplicitBlackoilPolymer_impl.hpp"
172 
173 #endif // OPM_SIMULATORFULLYIMPLICITBLACKOILPOLYMER_HEADER_INCLUDED
A nonlinear solver class suitable for general fully-implicit models, as well as pressure, transport and sequential models.
Definition: NonlinearSolver.hpp:37
Definition: PolymerPropsAd.hpp:32
Definition: WellStateFullyImplicitBlackoilPolymer.hpp:28
Class collecting all necessary components for a two-phase simulation.
Definition: SimulatorBase.hpp:85
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: AdditionalObjectDeleter.hpp:22
Facility for converting component rates at surface conditions to phase (voidage) rates at reservoir c...
Class for handling the standard well model.
Definition: StandardWells.hpp:51
Simulator state for a compressible two-phase simulator with polymer.
Definition: PolymerBlackoilState.hpp:33
A model implementation for three-phase black oil with polymer.
Definition: BlackoilPolymerModel.hpp:45
Class collecting all necessary components for a blackoil simulation with polymer injection.
Definition: SimulatorFullyImplicitBlackoilPolymer.hpp:81
Class containing static geological properties that are derived from grid and petrophysical properties...
Definition: GeoProps.hpp:59
Wrapper class for VTK, Matlab, and ECL output.
Definition: SimulatorFullyImplicitBlackoilOutput.hpp:206
Interface class for (linear) solvers for the fully implicit black-oil system.
Definition: NewtonIterationBlackoilInterface.hpp:31
This class implements the AD-adapted fluid interface for three-phase black-oil.
Definition: BlackoilPropsAdFromDeck.hpp:61
Definition: SimulatorBase.hpp:81
Definition: SimulatorTimer.hpp:34