All Classes Namespaces Files Functions Variables Typedefs Enumerator Pages
BlackoilSolventModel.hpp
1 /*
2  Copyright 2015 IRIS AS
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_BLACKOILSOLVENTMODEL_HEADER_INCLUDED
21 #define OPM_BLACKOILSOLVENTMODEL_HEADER_INCLUDED
22 
23 #include <opm/autodiff/BlackoilModelBase.hpp>
24 #include <opm/autodiff/BlackoilModelParameters.hpp>
25 #include <opm/autodiff/WellStateFullyImplicitBlackoilSolvent.hpp>
26 #include <opm/autodiff/SolventPropsAdFromDeck.hpp>
27 #include <opm/autodiff/StandardWellsSolvent.hpp>
28 
29 namespace Opm {
30 
37  template<class Grid>
38  class BlackoilSolventModel : public BlackoilModelBase<Grid, StandardWellsSolvent, BlackoilSolventModel<Grid> >
39  {
40  public:
41 
42  // --------- Types and enums ---------
43 
45  typedef typename Base::ReservoirState ReservoirState;
46  typedef typename Base::WellState WellState;
47  // The next line requires C++11 support available in g++ 4.7.
48  // friend Base;
49  friend class BlackoilModelBase<Grid, StandardWellsSolvent, BlackoilSolventModel<Grid> >;
50 
67  BlackoilSolventModel(const typename Base::ModelParameters& param,
68  const Grid& grid,
69  const BlackoilPropsAdFromDeck& fluid,
70  const DerivedGeology& geo,
71  const RockCompressibility* rock_comp_props,
72  const SolventPropsAdFromDeck& solvent_props,
73  const StandardWellsSolvent& well_model,
74  const NewtonIterationBlackoilInterface& linsolver,
75  std::shared_ptr< const EclipseState > eclState,
76  const bool has_disgas,
77  const bool has_vapoil,
78  const bool terminal_output,
79  const bool has_solvent,
80  const bool is_miscible);
81 
86  void updateState(const V& dx,
87  ReservoirState& reservoir_state,
88  WellState& well_state);
89 
90  using Base::wellModel;
91 
92 
93  std::vector<std::vector<double> >
94  computeFluidInPlace(const ReservoirState& x,
95  const std::vector<int>& fipnum);
96 
97  protected:
98 
99  // --------- Types and enums ---------
100 
101  typedef typename Base::SolutionState SolutionState;
102  typedef typename Base::DataBlock DataBlock;
103  enum { Solvent = CanonicalVariablePositions::Next };
104 
105  // --------- Data members ---------
106  const bool has_solvent_;
107  const int solvent_pos_;
108  const SolventPropsAdFromDeck& solvent_props_;
109  const bool is_miscible_;
110  std::vector<ADB> mu_eff_;
111  std::vector<ADB> b_eff_;
112 
113  // Need to declare Base members we want to use here.
114  using Base::grid_;
115  using Base::fluid_;
116  using Base::geo_;
117  using Base::rock_comp_props_;
118  using Base::linsolver_;
119  using Base::active_;
120  using Base::canph_;
121  using Base::cells_;
122  using Base::ops_;
123  using Base::has_disgas_;
124  using Base::has_vapoil_;
125  using Base::param_;
126  using Base::use_threshold_pressure_;
127  using Base::threshold_pressures_by_connection_;
128  using Base::sd_;
129  using Base::phaseCondition_;
130  using Base::residual_;
132  using Base::pvdt_;
133 
134  // --------- Protected methods ---------
135 
136  // Need to declare Base members we want to use here.
137  using Base::wells;
138  using Base::variableState;
139  using Base::computeGasPressure;
140  using Base::applyThresholdPressures;
141  using Base::fluidRsSat;
142  using Base::fluidRvSat;
143  using Base::poroMult;
144  using Base::transMult;
147  using Base::dpMaxRel;
148  using Base::dsMax;
149  using Base::drMaxRel;
150  using Base::maxResidualAllowed;
151  // using Base::updateWellControls;
152  // using Base::computeWellConnectionPressures;
153  // using Base::addWellControlEq;
154  // using Base::computePropertiesForWellConnectionPressures;
155 
156  std::vector<ADB>
157  computeRelPerm(const SolutionState& state) const;
158 
159  ADB
160  fluidViscosity(const int phase,
161  const ADB& p ,
162  const ADB& temp ,
163  const ADB& rs ,
164  const ADB& rv ,
165  const std::vector<PhasePresence>& cond) const;
166 
167  ADB
168  fluidReciprocFVF(const int phase,
169  const ADB& p ,
170  const ADB& temp ,
171  const ADB& rs ,
172  const ADB& rv ,
173  const std::vector<PhasePresence>& cond) const;
174 
175  ADB
176  fluidDensity(const int phase,
177  const ADB& b,
178  const ADB& rs,
179  const ADB& rv) const;
180 
181  void
182  makeConstantState(SolutionState& state) const;
183 
184  std::vector<V>
185  variableStateInitials(const ReservoirState& x,
186  const WellState& xw) const;
187 
188  std::vector<int>
189  variableStateIndices() const;
190 
191  SolutionState
192  variableStateExtractVars(const ReservoirState& x,
193  const std::vector<int>& indices,
194  std::vector<ADB>& vars) const;
195 
196  void
197  computeAccum(const SolutionState& state,
198  const int aix );
199 
200  void
201  assembleMassBalanceEq(const SolutionState& state);
202 
203  void
204  addWellContributionToMassBalanceEq(const std::vector<ADB>& cq_s,
205  const SolutionState& state,
206  WellState& xw);
207 
208  void updateEquationsScaling();
209 
210  void
211  computeMassFlux(const int actph ,
212  const V& transi,
213  const ADB& kr ,
214  const ADB& mu ,
215  const ADB& rho ,
216  const ADB& p ,
217  const SolutionState& state );
218 
219  const std::vector<PhasePresence>
220  phaseCondition() const {return this->phaseCondition_;}
221 
222  // compute effective viscosities (mu_eff_) and effective b factors (b_eff_) using the ToddLongstaff model
223  void computeEffectiveProperties(const SolutionState& state);
224 
225  // compute density and viscosity using the ToddLongstaff mixing model
226  void computeToddLongstaffMixing(std::vector<ADB>& viscosity, std::vector<ADB>& density, const std::vector<ADB>& saturations, const ADB po, const Opm::PhaseUsage pu);
227 
228  // compute phase pressures.
229  std::vector<ADB>
230  computePressures(const ADB& po,
231  const ADB& sw,
232  const ADB& so,
233  const ADB& sg,
234  const ADB& ss) const;
235  };
236 
237 
238 
242  {
243  explicit BlackoilSolventSolutionState(const int np)
245  solvent_saturation( ADB::null())
246  {
247  }
248  ADB solvent_saturation;
249  };
250 
251 
252 
254  template <class Grid>
256  {
257  typedef BlackoilState ReservoirState;
261  };
262 
263 } // namespace Opm
264 
265 #include "BlackoilSolventModel_impl.hpp"
266 
267 #endif // OPM_BLACKOILSOLVENTMODEL_HEADER_INCLUDED
static AutoDiffBlock null()
Construct an empty AutoDiffBlock.
Definition: AutoDiffBlock.hpp:104
Definition: WellStateFullyImplicitBlackoilSolvent.hpp:28
BlackoilSolventModel(const typename Base::ModelParameters &param, const Grid &grid, const BlackoilPropsAdFromDeck &fluid, const DerivedGeology &geo, const RockCompressibility *rock_comp_props, const SolventPropsAdFromDeck &solvent_props, const StandardWellsSolvent &well_model, const NewtonIterationBlackoilInterface &linsolver, std::shared_ptr< const EclipseState > eclState, const bool has_disgas, const bool has_vapoil, const bool terminal_output, const bool has_solvent, const bool is_miscible)
Construct the model.
Definition: BlackoilSolventModel_impl.hpp:73
Need to include concentration in our state variables, otherwise all is as the default blackoil model...
Definition: BlackoilSolventModel.hpp:241
A model implementation for three-phase black oil.
Definition: BlackoilModelBase.hpp:76
const Wells & wells() const
return the Well struct in the WellModel
Definition: BlackoilModelBase.hpp:383
Solver parameters for the BlackoilModel.
Definition: BlackoilModelParameters.hpp:31
void updateState(const V &dx, ReservoirState &reservoir_state, WellState &well_state)
Apply an update to the primary variables, chopped if appropriate.
Definition: BlackoilSolventModel_impl.hpp:389
void updatePrimalVariableFromState(const ReservoirState &state)
update the primal variable for Sg, Rv or Rs.
Traits to encapsulate the types used by classes using or extending this model.
Definition: BlackoilModelBase.hpp:60
StandardWellsSolvent & wellModel()
return the WellModel object
Definition: BlackoilModelBase.hpp:268
bool terminal_output_
Whether we print something to std::cout.
Definition: BlackoilModelBase.hpp:353
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
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: SolventPropsAdFromDeck.hpp:37
void updatePhaseCondFromPrimalVariable(const ReservoirState &state)
Update the phaseCondition_ member based on the primalVariable_ member.
A model implementation for three-phase black oil with one extra component.
Definition: BlackoilSolventModel.hpp:38