All Classes Namespaces Files Functions Variables Typedefs Enumerator Pages
BlackoilPropsAdFromDeck.hpp
1 /*
2  Copyright 2013 SINTEF ICT, Applied Mathematics.
3  Copyright 2015 Dr. Blatt - HPC-Simulation-Software & Services.
4  Copyright 2015 NTNU.
5 
6  This file is part of the Open Porous Media project (OPM).
7 
8  OPM is free software: you can redistribute it and/or modify
9  it under the terms of the GNU General Public License as published by
10  the Free Software Foundation, either version 3 of the License, or
11  (at your option) any later version.
12 
13  OPM is distributed in the hope that it will be useful,
14  but WITHOUT ANY WARRANTY; without even the implied warranty of
15  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16  GNU General Public License for more details.
17 
18  You should have received a copy of the GNU General Public License
19  along with OPM. If not, see <http://www.gnu.org/licenses/>.
20 */
21 
22 #ifndef OPM_BLACKOILPROPSADFROMDECK_HEADER_INCLUDED
23 #define OPM_BLACKOILPROPSADFROMDECK_HEADER_INCLUDED
24 
25 #include <opm/autodiff/AutoDiffBlock.hpp>
26 #include <opm/autodiff/BlackoilModelEnums.hpp>
27 
28 #include <opm/core/props/satfunc/SaturationPropsFromDeck.hpp>
29 #include <opm/core/props/rock/RockFromDeck.hpp>
30 
31 #include <opm/material/fluidsystems/BlackOilFluidSystem.hpp>
32 #include <opm/material/densead/Math.hpp>
33 #include <opm/material/densead/Evaluation.hpp>
34 
35 #include <opm/parser/eclipse/Deck/Deck.hpp>
36 #include <opm/parser/eclipse/EclipseState/EclipseState.hpp>
37 
38 #include <memory>
39 #include <array>
40 #include <vector>
41 
42 #ifdef HAVE_OPM_GRID
43 #include <opm/common/utility/platform_dependent/disable_warnings.h>
44 #include <dune/grid/CpGrid.hpp>
45 #include <opm/common/utility/platform_dependent/reenable_warnings.h>
46 #endif
47 
48 namespace Opm
49 {
50  class PvtInterface;
51 
62  {
63  friend class BlackoilPropsDataHandle;
64  public:
65  typedef FluidSystems::BlackOil<double> FluidSystem;
66  typedef Opm::GasPvtMultiplexer<double> GasPvt;
67  typedef Opm::OilPvtMultiplexer<double> OilPvt;
68  typedef Opm::WaterPvtMultiplexer<double> WaterPvt;
69 
70  typedef typename SaturationPropsFromDeck::MaterialLawManager MaterialLawManager;
71 
86  BlackoilPropsAdFromDeck(const Opm::Deck& deck,
87  const Opm::EclipseState& eclState,
88  std::shared_ptr<MaterialLawManager> materialLawManager,
89  const UnstructuredGrid& grid,
90  const bool init_rock = true );
91 
92 #ifdef HAVE_OPM_GRID
93  BlackoilPropsAdFromDeck(const Opm::Deck& deck,
108  const Opm::EclipseState& eclState,
109  std::shared_ptr<MaterialLawManager> materialLawManager,
110  const Dune::CpGrid& grid,
111  const bool init_rock = true );
112 #endif
113 
121  BlackoilPropsAdFromDeck(const Opm::Deck& deck,
122  const Opm::EclipseState& eclState,
123  const UnstructuredGrid& grid,
124  const bool init_rock = true );
125 
126 #ifdef HAVE_OPM_GRID
127  BlackoilPropsAdFromDeck(const Opm::Deck& deck,
135  const Opm::EclipseState& eclState,
136  const Dune::CpGrid& grid,
137  const bool init_rock = true );
138 #endif
139 
153  std::shared_ptr<MaterialLawManager> materialLawManager,
154  const int number_of_cells);
155 
156 
158  // Rock interface //
160 
162  int numDimensions() const;
163 
165  int numCells() const;
166 
169  virtual const int* cellPvtRegionIndex() const
170  { return &cellPvtRegionIdx_[0]; }
171 
173  const double* porosity() const;
174 
178  const double* permeability() const;
179 
180 
182  // Fluid interface //
184 
185  typedef AutoDiffBlock<double> ADB;
186  typedef ADB::V V;
187  typedef std::vector<int> Cells;
188 
190  int numPhases() const;
191 
193  PhaseUsage phaseUsage() const;
194 
195  // ------ Density ------
196 
201  V surfaceDensity(const int phaseIdx , const Cells& cells) const;
202 
203 
204  // ------ Viscosity ------
205 
211  ADB muWat(const ADB& pw,
212  const ADB& T,
213  const Cells& cells) const;
214 
222  ADB muOil(const ADB& po,
223  const ADB& T,
224  const ADB& rs,
225  const std::vector<PhasePresence>& cond,
226  const Cells& cells) const;
227 
235  ADB muGas(const ADB& pg,
236  const ADB& T,
237  const ADB& rv,
238  const std::vector<PhasePresence>& cond,
239  const Cells& cells) const;
240 
241  // ------ Formation volume factor (b) ------
242 
248  ADB bWat(const ADB& pw,
249  const ADB& T,
250  const Cells& cells) const;
251 
259  ADB bOil(const ADB& po,
260  const ADB& T,
261  const ADB& rs,
262  const std::vector<PhasePresence>& cond,
263  const Cells& cells) const;
264 
272  ADB bGas(const ADB& pg,
273  const ADB& T,
274  const ADB& rv,
275  const std::vector<PhasePresence>& cond,
276  const Cells& cells) const;
277 
278  // ------ Rs bubble point curve ------
279 
284  V rsSat(const V& po,
285  const Cells& cells) const;
286 
292  V rsSat(const V& po,
293  const V& so,
294  const Cells& cells) const;
295 
300  ADB rsSat(const ADB& po,
301  const Cells& cells) const;
302 
308  ADB rsSat(const ADB& po,
309  const ADB& so,
310  const Cells& cells) const;
311 
312  // ------ Rv condensation curve ------
313 
318  ADB rvSat(const ADB& po,
319  const Cells& cells) const;
320 
326  ADB rvSat(const ADB& po,
327  const ADB& so,
328  const Cells& cells) const;
329 
330  // ------ Relative permeability ------
331 
339  std::vector<ADB> relperm(const ADB& sw,
340  const ADB& so,
341  const ADB& sg,
342  const Cells& cells) const;
343 
352  std::vector<ADB> capPress(const ADB& sw,
353  const ADB& so,
354  const ADB& sg,
355  const Cells& cells) const;
356 
359  void updateSatHyst(const std::vector<double>& saturation,
360  const std::vector<int>& cells);
361 
365  void setGasOilHystParams(const std::vector<double>& pcswmdc,
366  const std::vector<double>& krnswdc,
367  const std::vector<int>& cells);
368 
372  void getGasOilHystParams(std::vector<double>& pcswmdc,
373  std::vector<double>& krnswdc,
374  const std::vector<int>& cells) const;
375 
379  void setOilWaterHystParams(const std::vector<double>& pcswmdc,
380  const std::vector<double>& krnswdc,
381  const std::vector<int>& cells);
382 
386  void getOilWaterHystParams(std::vector<double>& pcswmdc,
387  std::vector<double>& krnswdc,
388  const std::vector<int>& cells) const;
389 
392  void updateSatOilMax(const std::vector<double>& saturation);
393 
395  const std::vector<double>& satOilMax() const;
396 
402  void setSatOilMax(const std::vector<double>& max_sat);
403 
405  std::vector<double> bubblePointPressure(const Cells& cells,
406  const V& T,
407  const V& rs) const;
408 
410  std::vector<double> dewPointPressure(const Cells& cells,
411  const V& T,
412  const V& rv) const;
413 
417  void setSwatInitScaling(const std::vector<double>& saturation,
418  const std::vector<double>& pc);
419 
423  V scaledCriticalOilinGasSaturations(const Cells& cells) const;
424 
428  V scaledCriticalGasSaturations(const Cells& cells) const;
429 
431  const WaterPvt& waterProps() const
432  {
433  return FluidSystem::waterPvt();
434  }
435 
437  const OilPvt& oilProps() const
438  {
439  return FluidSystem::oilPvt();
440  }
441 
443  const GasPvt& gasProps() const
444  {
445  return FluidSystem::gasPvt();
446  }
447 
449  const MaterialLawManager& materialLaws() const
450  {
451  return *materialLawManager_;
452  }
453 
454  // Direct access to pvt region indices.
455  const std::vector<int>& pvtRegions() const
456  {
457  return cellPvtRegionIdx_;
458  }
459 
460 
461  private:
463  void init(const Opm::Deck& deck,
464  const Opm::EclipseState& eclState,
465  std::shared_ptr<MaterialLawManager> materialLawManager,
466  int number_of_cells,
467  const int* global_cell,
468  const int* cart_dims,
469  const bool init_rock);
470 
472  void applyVap(V& r,
473  const V& so,
474  const std::vector<int>& cells,
475  const double vap) const;
476 
477  void applyVap(ADB& r,
478  const ADB& so,
479  const std::vector<int>& cells,
480  const double vap) const;
481 
482  RockFromDeck rock_;
483 
484  // This has to be a shared pointer as we must
485  // be able to make a copy of *this in the parallel case.
486  std::shared_ptr<MaterialLawManager> materialLawManager_;
487  std::shared_ptr<SaturationPropsFromDeck> satprops_;
488 
489  PhaseUsage phase_usage_;
490  // bool has_vapoil_;
491  // bool has_disgas_;
492 
493  // The PVT region which is to be used for each cell
494  std::vector<int> cellPvtRegionIdx_;
495 
496  // VAPPARS
497  double vap1_;
498  double vap2_;
499  std::vector<double> satOilMax_;
500  double vap_satmax_guard_; //Threshold value to promote stability
501  };
502 } // namespace Opm
503 
504 #endif // OPM_BLACKOILPROPSADFROMDECK_HEADER_INCLUDED
std::vector< double > bubblePointPressure(const Cells &cells, const V &T, const V &rs) const
Returns the bubble point pressures.
Definition: BlackoilPropsAdFromDeck.cpp:958
ADB muGas(const ADB &pg, const ADB &T, const ADB &rv, const std::vector< PhasePresence > &cond, const Cells &cells) const
Gas viscosity.
Definition: BlackoilPropsAdFromDeck.cpp:381
std::vector< ADB > relperm(const ADB &sw, const ADB &so, const ADB &sg, const Cells &cells) const
Relative permeabilities for all phases.
Definition: BlackoilPropsAdFromDeck.cpp:748
ADB bWat(const ADB &pw, const ADB &T, const Cells &cells) const
Water formation volume factor.
Definition: BlackoilPropsAdFromDeck.cpp:446
int numDimensions() const
Definition: BlackoilPropsAdFromDeck.cpp:194
ADB muWat(const ADB &pw, const ADB &T, const Cells &cells) const
Water viscosity.
Definition: BlackoilPropsAdFromDeck.cpp:263
V rsSat(const V &po, const Cells &cells) const
Bubble point curve for Rs as function of oil pressure.
void getOilWaterHystParams(std::vector< double > &pcswmdc, std::vector< double > &krnswdc, const std::vector< int > &cells) const
Set oil-water hysteresis parameters.
Definition: BlackoilPropsAdFromDeck.cpp:921
void updateSatOilMax(const std::vector< double > &saturation)
Update for max oil saturation.
Definition: BlackoilPropsAdFromDeck.cpp:932
ADB rvSat(const ADB &po, const Cells &cells) const
Condensation curve for Rv as function of oil pressure.
Definition: BlackoilPropsAdFromDeck.cpp:688
V scaledCriticalGasSaturations(const Cells &cells) const
Obtain the scaled critical gas saturation values.
Definition: BlackoilPropsAdFromDeck.cpp:1102
const OilPvt & oilProps() const
Direct access to lower-level oil pvt props.
Definition: BlackoilPropsAdFromDeck.hpp:437
int numPhases() const
Definition: BlackoilPropsAdFromDeck.cpp:225
void setSatOilMax(const std::vector< double > &max_sat)
Force set max oil saturation (used for restarting)
Definition: BlackoilPropsAdFromDeck.cpp:952
void setGasOilHystParams(const std::vector< double > &pcswmdc, const std::vector< double > &krnswdc, const std::vector< int > &cells)
Set gas-oil hysteresis parameters.
Definition: BlackoilPropsAdFromDeck.cpp:882
ADB bOil(const ADB &po, const ADB &T, const ADB &rs, const std::vector< PhasePresence > &cond, const Cells &cells) const
Oil formation volume factor.
Definition: BlackoilPropsAdFromDeck.cpp:493
int numCells() const
Definition: BlackoilPropsAdFromDeck.cpp:200
void getGasOilHystParams(std::vector< double > &pcswmdc, std::vector< double > &krnswdc, const std::vector< int > &cells) const
Get gas-oil hysteresis parameters.
Definition: BlackoilPropsAdFromDeck.cpp:895
const double * porosity() const
Definition: BlackoilPropsAdFromDeck.cpp:206
std::vector< double > dewPointPressure(const Cells &cells, const V &T, const V &rv) const
Returns the dew point pressures.
Definition: BlackoilPropsAdFromDeck.cpp:982
void updateSatHyst(const std::vector< double > &saturation, const std::vector< int > &cells)
Saturation update for hysteresis behavior.
Definition: BlackoilPropsAdFromDeck.cpp:872
const std::vector< double > & satOilMax() const
Returns the max oil saturation vector.
Definition: BlackoilPropsAdFromDeck.cpp:946
V scaledCriticalOilinGasSaturations(const Cells &cells) const
Obtain the scaled critical oil in gas saturation values.
Definition: BlackoilPropsAdFromDeck.cpp:1082
void setOilWaterHystParams(const std::vector< double > &pcswmdc, const std::vector< double > &krnswdc, const std::vector< int > &cells)
Set oil-water hysteresis parameters.
Definition: BlackoilPropsAdFromDeck.cpp:908
ADB bGas(const ADB &pg, const ADB &T, const ADB &rv, const std::vector< PhasePresence > &cond, const Cells &cells) const
Gas formation volume factor.
Definition: BlackoilPropsAdFromDeck.cpp:566
const WaterPvt & waterProps() const
Direct access to lower-level water pvt props.
Definition: BlackoilPropsAdFromDeck.hpp:431
const MaterialLawManager & materialLaws() const
Direct access to lower-level saturation functions.
Definition: BlackoilPropsAdFromDeck.hpp:449
V surfaceDensity(const int phaseIdx, const Cells &cells) const
Densities of stock components at surface conditions.
Definition: BlackoilPropsAdFromDeck.cpp:242
const double * permeability() const
Definition: BlackoilPropsAdFromDeck.cpp:214
std::vector< ADB > capPress(const ADB &sw, const ADB &so, const ADB &sg, const Cells &cells) const
Capillary pressure for all phases.
Definition: BlackoilPropsAdFromDeck.cpp:809
void setSwatInitScaling(const std::vector< double > &saturation, const std::vector< double > &pc)
Set capillary pressure scaling according to pressure diff.
Definition: BlackoilPropsAdFromDeck.cpp:1008
virtual const int * cellPvtRegionIndex() const
Return an array containing the PVT table index for each grid cell.
Definition: BlackoilPropsAdFromDeck.hpp:169
const GasPvt & gasProps() const
Direct access to lower-level gas pvt props.
Definition: BlackoilPropsAdFromDeck.hpp:443
PhaseUsage phaseUsage() const
Definition: BlackoilPropsAdFromDeck.cpp:231
This class implements the AD-adapted fluid interface for three-phase black-oil.
Definition: BlackoilPropsAdFromDeck.hpp:61
Eigen::Array< Scalar, Eigen::Dynamic, 1 > V
Underlying type for values.
Definition: AutoDiffBlock.hpp:99
BlackoilPropsAdFromDeck(const Opm::Deck &deck, const Opm::EclipseState &eclState, std::shared_ptr< MaterialLawManager > materialLawManager, const UnstructuredGrid &grid, const bool init_rock=true)
Constructor to create a blackoil properties from an ECL deck.
Definition: BlackoilPropsAdFromDeck.cpp:48
ADB muOil(const ADB &po, const ADB &T, const ADB &rs, const std::vector< PhasePresence > &cond, const Cells &cells) const
Oil viscosity.
Definition: BlackoilPropsAdFromDeck.cpp:314