BlackoilPolymerModel.hpp
1 /*
2  Copyright 2013, 2015 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_BLACKOILPOLYMERMODEL_HEADER_INCLUDED
22 #define OPM_BLACKOILPOLYMERMODEL_HEADER_INCLUDED
23 
24 #include <opm/autodiff/BlackoilModelBase.hpp>
25 #include <opm/autodiff/BlackoilModelParameters.hpp>
26 #include <opm/polymer/PolymerProperties.hpp>
27 #include <opm/polymer/fullyimplicit/PolymerPropsAd.hpp>
28 #include <opm/polymer/PolymerBlackoilState.hpp>
29 #include <opm/polymer/fullyimplicit/WellStateFullyImplicitBlackoilPolymer.hpp>
30 #include <opm/autodiff/StandardWells.hpp>
31 #include <opm/simulators/timestepping/SimulatorTimerInterface.hpp>
32 
33 namespace Opm {
34 
44  template<class Grid>
45  class BlackoilPolymerModel : public BlackoilModelBase<Grid, StandardWells, BlackoilPolymerModel<Grid> >
46  {
47  public:
48 
49  // --------- Types and enums ---------
50 
52  typedef typename Base::ReservoirState ReservoirState;
53  typedef typename Base::WellState WellState;
54  // The next line requires C++11 support available in g++ 4.7.
55  // friend Base;
56  friend class BlackoilModelBase<Grid, StandardWells, BlackoilPolymerModel<Grid> >;
57 
77  BlackoilPolymerModel(const typename Base::ModelParameters& param,
78  const Grid& grid,
79  const BlackoilPropsAdFromDeck& fluid,
80  const DerivedGeology& geo,
81  const RockCompressibility* rock_comp_props,
82  const PolymerPropsAd& polymer_props_ad,
83  const StandardWells& well_model,
84  const NewtonIterationBlackoilInterface& linsolver,
85  std::shared_ptr< const EclipseState > eclipseState,
86  const bool has_disgas,
87  const bool has_vapoil,
88  const bool has_polymer,
89  const bool has_plyshlog,
90  const bool has_shrate,
91  const std::vector<double>& wells_rep_radius,
92  const std::vector<double>& wells_perf_length,
93  const std::vector<double>& wells_bore_diameter,
94  const bool terminal_output);
95 
100  void prepareStep(const SimulatorTimerInterface& timer,
101  const ReservoirState& reservoir_state,
102  const WellState& well_state);
103 
108  void afterStep(const SimulatorTimerInterface& timer,
109  ReservoirState& reservoir_state,
110  WellState& well_state);
111 
116  void updateState(const V& dx,
117  ReservoirState& reservoir_state,
118  WellState& well_state);
119 
124  SimulatorReport
125  assemble(const ReservoirState& reservoir_state,
126  WellState& well_state,
127  const bool initial_assembly);
128 
129  using Base::wellModel;
130 
131  protected:
132 
133  // --------- Types and enums ---------
134 
135  typedef typename Base::SolutionState SolutionState;
136  typedef typename Base::DataBlock DataBlock;
137  enum { Concentration = CanonicalVariablePositions::Next };
138 
139  // --------- Data members ---------
140 
141  const PolymerPropsAd& polymer_props_ad_;
142  const bool has_polymer_;
143  const bool has_plyshlog_;
144  const bool has_shrate_;
145  const int poly_pos_;
146  V cmax_;
147 
148  // representative radius and perforation length of well perforations
149  // to be used in shear-thinning computation.
150  std::vector<double> wells_rep_radius_;
151  std::vector<double> wells_perf_length_;
152  // wellbore diameters
153  std::vector<double> wells_bore_diameter_;
154 
155  // shear-thinning factor for cell faces
156  std::vector<double> shear_mult_faces_;
157  // shear-thinning factor for well perforations
158  std::vector<double> shear_mult_wells_;
159 
160  // Need to declare Base members we want to use here.
161  using Base::grid_;
162  using Base::fluid_;
163  using Base::geo_;
164  using Base::rock_comp_props_;
165  using Base::linsolver_;
166  using Base::active_;
167  using Base::canph_;
168  using Base::cells_;
169  using Base::ops_;
170  using Base::has_disgas_;
171  using Base::has_vapoil_;
172  using Base::param_;
173  using Base::use_threshold_pressure_;
174  using Base::threshold_pressures_by_connection_;
175  using Base::sd_;
176  using Base::phaseCondition_;
177  using Base::residual_;
179  using Base::pvdt_;
180  using Base::vfp_properties_;
181 
182  // --------- Protected methods ---------
183 
184  // Need to declare Base members we want to use here.
185  using Base::wells;
186  using Base::wellsActive;
187  using Base::variableState;
188  using Base::computePressures;
189  using Base::computeGasPressure;
190  using Base::applyThresholdPressures;
191  using Base::fluidViscosity;
192  using Base::fluidReciprocFVF;
193  using Base::fluidDensity;
194  using Base::fluidRsSat;
195  using Base::fluidRvSat;
196  using Base::poroMult;
197  using Base::transMult;
200  using Base::dpMaxRel;
201  using Base::dsMax;
202  using Base::drMaxRel;
203  using Base::maxResidualAllowed;
204 
205  // using Base::updateWellControls;
206  // using Base::computeWellConnectionPressures;
207  // using Base::addWellControlEq;
208  using Base::computeRelPerm;
209 
210 
211  void
212  makeConstantState(SolutionState& state) const;
213 
214  std::vector<V>
215  variableStateInitials(const ReservoirState& x,
216  const WellState& xw) const;
217 
218  std::vector<int>
219  variableStateIndices() const;
220 
221  SolutionState
222  variableStateExtractVars(const ReservoirState& x,
223  const std::vector<int>& indices,
224  std::vector<ADB>& vars) const;
225 
226  void
227  computeAccum(const SolutionState& state,
228  const int aix );
229 
230  void
231  computeInjectionMobility(const SolutionState& state,
232  std::vector<ADB>& mob_perfcells);
233 
234  void
235  assembleMassBalanceEq(const SolutionState& state);
236 
237  void
238  addWellContributionToMassBalanceEq(const std::vector<ADB>& cq_s,
239  const SolutionState& state,
240  WellState& xw);
241 
242  void updateEquationsScaling();
243 
244  void
245  computeMassFlux(const int actph ,
246  const V& transi,
247  const ADB& kr ,
248  const ADB& mu ,
249  const ADB& rho ,
250  const ADB& p ,
251  const SolutionState& state );
252 
253  void
254  computeCmax(ReservoirState& state);
255 
256  ADB
257  computeMc(const SolutionState& state) const;
258 
259  const std::vector<PhasePresence>
260  phaseCondition() const {return this->phaseCondition_;}
261 
264  void computeWaterShearVelocityFaces(const V& transi,
265  const std::vector<ADB>& phasePressure, const SolutionState& state,
266  std::vector<double>& water_vel, std::vector<double>& visc_mult);
267 
270  void computeWaterShearVelocityWells(const SolutionState& state, WellState& xw, const ADB& cq_sw,
271  std::vector<double>& water_vel_wells, std::vector<double>& visc_mult_wells);
272 
273  };
274 
275 
276 
280  {
281  explicit BlackoilPolymerSolutionState(const int np)
283  concentration( ADB::null())
284  {
285  }
286  ADB concentration;
287  };
288 
289 
290 
292  template <class Grid>
294  {
299  };
300 
301 } // namespace Opm
302 
303 #include "BlackoilPolymerModel_impl.hpp"
304 
305 
306 #endif // OPM_BLACKOILPOLYMERMODEL_HEADER_INCLUDED
static AutoDiffBlock null()
Construct an empty AutoDiffBlock.
Definition: AutoDiffBlock.hpp:104
Need to include concentration in our state variables, otherwise all is as the default blackoil model...
Definition: BlackoilPolymerModel.hpp:279
void computeWaterShearVelocityFaces(const V &transi, const std::vector< ADB > &phasePressure, const SolutionState &state, std::vector< double > &water_vel, std::vector< double > &visc_mult)
Computing the water velocity without shear-thinning for the cell faces.
Definition: BlackoilPolymerModel_impl.hpp:599
Definition: PolymerPropsAd.hpp:32
SimulatorReport assemble(const ReservoirState &reservoir_state, WellState &well_state, const bool initial_assembly)
Assemble the residual and Jacobian of the nonlinear system.
Definition: BlackoilPolymerModel_impl.hpp:496
Definition: WellStateFullyImplicitBlackoilPolymer.hpp:28
A model implementation for three-phase black oil.
Definition: BlackoilModelBase.hpp:76
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: AdditionalObjectDeleter.hpp:22
Solver parameters for the BlackoilModel.
Definition: BlackoilModelParameters.hpp:31
Class for handling the standard well model.
Definition: StandardWells.hpp:51
void updateState(const V &dx, ReservoirState &reservoir_state, WellState &well_state)
Apply an update to the primary variables, chopped if appropriate.
Definition: BlackoilPolymerModel_impl.hpp:409
void prepareStep(const SimulatorTimerInterface &timer, const ReservoirState &reservoir_state, const WellState &well_state)
Called once before each time step.
Definition: BlackoilPolymerModel_impl.hpp:124
void afterStep(const SimulatorTimerInterface &timer, ReservoirState &reservoir_state, WellState &well_state)
Called once after each time step.
Definition: BlackoilPolymerModel_impl.hpp:141
void updatePrimalVariableFromState(const ReservoirState &state)
update the primal variable for Sg, Rv or Rs.
Definition: BlackoilModelBase_impl.hpp:2193
Simulator state for a compressible two-phase simulator with polymer.
Definition: PolymerBlackoilState.hpp:33
Traits to encapsulate the types used by classes using or extending this model.
Definition: BlackoilModelBase.hpp:60
StandardWells & wellModel()
return the WellModel object
Definition: BlackoilModelBase.hpp:268
BlackoilPolymerModel(const typename Base::ModelParameters &param, const Grid &grid, const BlackoilPropsAdFromDeck &fluid, const DerivedGeology &geo, const RockCompressibility *rock_comp_props, const PolymerPropsAd &polymer_props_ad, const StandardWells &well_model, const NewtonIterationBlackoilInterface &linsolver, std::shared_ptr< const EclipseState > eclipseState, const bool has_disgas, const bool has_vapoil, const bool has_polymer, const bool has_plyshlog, const bool has_shrate, const std::vector< double > &wells_rep_radius, const std::vector< double > &wells_perf_length, const std::vector< double > &wells_bore_diameter, const bool terminal_output)
Construct the model.
Definition: BlackoilPolymerModel_impl.hpp:77
A model implementation for three-phase black oil with polymer.
Definition: BlackoilPolymerModel.hpp:45
bool terminal_output_
Whether we print something to std::cout.
Definition: BlackoilModelBase.hpp:353
Interface class for SimulatorTimer objects, to be improved.
Definition: SimulatorTimerInterface.hpp:35
const Wells & wells() const
return the Well struct in the WellModel
Definition: BlackoilModelBase.hpp:383
Struct for containing iteration variables.
Definition: DefaultBlackoilSolutionState.hpp:29
Class containing static geological properties that are derived from grid and petrophysical properties...
Definition: GeoProps.hpp:59
bool wellsActive() const
return true if wells are available in the reservoir
Definition: BlackoilModelBase.hpp:386
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
void updatePhaseCondFromPrimalVariable(const ReservoirState &state)
Update the phaseCondition_ member based on the primalVariable_ member.
Definition: BlackoilModelBase_impl.hpp:2206
void computeWaterShearVelocityWells(const SolutionState &state, WellState &xw, const ADB &cq_sw, std::vector< double > &water_vel_wells, std::vector< double > &visc_mult_wells)
Computing the water velocity without shear-thinning for the well perforations based on the water flux...
Definition: BlackoilPolymerModel_impl.hpp:715