SolventPropsAdFromDeck.hpp
1 /*
2  Copyright 2015 IRIS
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 SOLVENTPROPSADFROMDECK_HPP
21 #define SOLVENTPROPSADFROMDECK_HPP
22 
23 #include <opm/autodiff/BlackoilPropsAdFromDeck.hpp>
24 #include <opm/autodiff/AutoDiffBlock.hpp>
25 
26 #include <opm/core/utility/NonuniformTableLinear.hpp>
27 
28 #include <opm/parser/eclipse/Deck/Deck.hpp>
29 #include <opm/parser/eclipse/EclipseState/EclipseState.hpp>
30 
31 #include <cmath>
32 #include <vector>
33 #include <opm/common/ErrorMacros.hpp>
34 
35 namespace Opm
36 {
38 {
39 public:
40  SolventPropsAdFromDeck(const Deck& deck,
41  const EclipseState& eclipseState,
42  const int number_of_cells,
43  const int* global_cell);
44 
46  // Fluid interface //
48 
49  typedef AutoDiffBlock<double> ADB;
50  typedef ADB::V V;
51  typedef std::vector<int> Cells;
52 
57  ADB bSolvent(const ADB& pg,
58  const Cells& cells) const;
59 
64  ADB muSolvent(const ADB& pg,
65  const Cells& cells) const;
66 
71  ADB gasRelPermMultiplier(const ADB& solventFraction,
72  const Cells& cells) const;
73 
78  ADB solventRelPermMultiplier(const ADB& solventFraction,
79  const Cells& cells) const;
80 
86  const Cells& cells) const;
87 
94  const Cells& cells) const;
95 
101  const Cells& cells) const;
102 
107  ADB miscibilityFunction(const ADB& solventFraction,
108  const Cells& cells) const;
109 
115  const Cells& cells) const;
116 
122  const Cells& cells) const;
123 
129  const Cells& cells) const;
130 
134  V solventSurfaceDensity(const Cells& cells) const;
135 
139  V mixingParameterViscosity(const Cells& cells) const;
140 
144  V mixingParameterDensity(const Cells& cells) const;
145 
150  ADB pressureMixingParameter(const ADB& po,
151  const Cells& cells) const;
152 
153 
154 private:
160  ADB makeADBfromTables(const ADB& X,
161  const Cells& cells,
162  const std::vector<int>& regionIdx,
163  const std::vector<NonuniformTableLinear<double>>& tables) const;
164 
172  void extractTableIndex(const std::string& keyword,
173  const Opm::EclipseState& eclState,
174  size_t numCompressed,
175  const int* compressedToCartesianCellIdx,
176  std::vector<int>& tableIdx) const;
177 
178  // The PVT region which is to be used for each cell
179  std::vector<int> cellPvtRegionIdx_;
180  std::vector<int> cellMiscRegionIdx_;
181  std::vector<int> cellSatNumRegionIdx_;
182  std::vector<NonuniformTableLinear<double> > b_;
183  std::vector<NonuniformTableLinear<double> > viscosity_;
184  std::vector<NonuniformTableLinear<double> > inverseBmu_;
185  std::vector<double> solvent_surface_densities_;
186  std::vector<NonuniformTableLinear<double> > krg_;
187  std::vector<NonuniformTableLinear<double> > krs_;
188  std::vector<NonuniformTableLinear<double> > krn_;
189  std::vector<NonuniformTableLinear<double> > mkro_;
190  std::vector<NonuniformTableLinear<double> > mkrsg_;
191  std::vector<NonuniformTableLinear<double> > misc_;
192  std::vector<NonuniformTableLinear<double> > pmisc_;
193  std::vector<NonuniformTableLinear<double> > sorwmis_;
194  std::vector<NonuniformTableLinear<double> > sgcwmis_;
195  std::vector<NonuniformTableLinear<double> > tlpmix_param_;
196  std::vector<double> mix_param_viscosity_;
197  std::vector<double> mix_param_density_;
198 };
199 
200 } // namespace OPM
201 
202 #endif // SOLVENTPROPSADFROMDECK_HPP
ADB miscibilityFunction(const ADB &solventFraction, const Cells &cells) const
Miscible function.
Definition: SolventPropsAdFromDeck.cpp:386
ADB miscibleSolventGasRelPermMultiplier(const ADB &Ssg, const Cells &cells) const
Miscible Solvent + Gas relPerm multiplier.
Definition: SolventPropsAdFromDeck.cpp:366
ADB bSolvent(const ADB &pg, const Cells &cells) const
Solvent formation volume factor.
Definition: SolventPropsAdFromDeck.cpp:339
ADB misicibleHydrocarbonWaterRelPerm(const ADB &Sn, const Cells &cells) const
Miscible hydrocrabon relPerm wrt water.
Definition: SolventPropsAdFromDeck.cpp:360
ADB pressureMiscibilityFunction(const ADB &po, const Cells &cells) const
Pressure dependent miscibility function.
Definition: SolventPropsAdFromDeck.cpp:393
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: AdditionalObjectDeleter.hpp:22
ADB muSolvent(const ADB &pg, const Cells &cells) const
Solvent viscosity.
Definition: SolventPropsAdFromDeck.cpp:313
ADB miscibleCriticalGasSaturationFunction(const ADB &Sw, const Cells &cells) const
Miscible critical gas saturation function.
Definition: SolventPropsAdFromDeck.cpp:404
V mixingParameterViscosity(const Cells &cells) const
Todd-Longstaff mixing parameter for viscosity calculation.
Definition: SolventPropsAdFromDeck.cpp:457
V mixingParameterDensity(const Cells &cells) const
Todd-Longstaff mixing parameter for density calculation.
Definition: SolventPropsAdFromDeck.cpp:471
V solventSurfaceDensity(const Cells &cells) const
Solvent surface density.
Definition: SolventPropsAdFromDeck.cpp:447
ADB miscibleOilRelPermMultiplier(const ADB &So, const Cells &cells) const
Miscible Oil relPerm multiplier.
Definition: SolventPropsAdFromDeck.cpp:376
ADB miscibleResidualOilSaturationFunction(const ADB &Sw, const Cells &cells) const
Miscible residual oil saturation function.
Definition: SolventPropsAdFromDeck.cpp:414
ADB pressureMixingParameter(const ADB &po, const Cells &cells) const
Todd-Longstaff pressure dependent mixing parameter.
Definition: SolventPropsAdFromDeck.cpp:485
ADB gasRelPermMultiplier(const ADB &solventFraction, const Cells &cells) const
Gas relPerm multipliers.
Definition: SolventPropsAdFromDeck.cpp:346
Definition: SolventPropsAdFromDeck.hpp:37
ADB solventRelPermMultiplier(const ADB &solventFraction, const Cells &cells) const
Solvent relPerm multipliers.
Definition: SolventPropsAdFromDeck.cpp:353
Eigen::Array< double, Eigen::Dynamic, 1 > V
Underlying type for values.
Definition: AutoDiffBlock.hpp:99