My Project
Loading...
Searching...
No Matches
LiveOilPvt.hpp
Go to the documentation of this file.
1// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
2// vi: set et ts=4 sw=4 sts=4:
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 2 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 Consult the COPYING file in the top-level source directory of this
20 module for the precise wording of the license and the list of
21 copyright holders.
22*/
27#ifndef OPM_LIVE_OIL_PVT_HPP
28#define OPM_LIVE_OIL_PVT_HPP
29
31#include <opm/common/OpmLog/OpmLog.hpp>
32
36
37namespace Opm {
38
39#if HAVE_ECL_INPUT
40class EclipseState;
41class Schedule;
42class SimpleTable;
43#endif
44
49template <class Scalar>
51{
52 using SamplingPoints = std::vector<std::pair<Scalar, Scalar>>;
53
54public:
57
58#if HAVE_ECL_INPUT
62 void initFromState(const EclipseState& eclState, const Schedule& schedule);
63
64private:
65 void extendPvtoTable_(unsigned regionIdx,
66 unsigned xIdx,
67 const SimpleTable& curTable,
68 const SimpleTable& masterTable);
69
70public:
71#endif // HAVE_ECL_INPUT
72
73 void setNumRegions(size_t numRegions);
74
75 void setVapPars(const Scalar, const Scalar par2)
76 {
77 vapPar2_ = par2;
78 }
79
83 void setReferenceDensities(unsigned regionIdx,
84 Scalar rhoRefOil,
85 Scalar rhoRefGas,
86 Scalar /*rhoRefWater*/);
87
93 void setSaturatedOilGasDissolutionFactor(unsigned regionIdx, const SamplingPoints& samplePoints)
94 { saturatedGasDissolutionFactorTable_[regionIdx].setContainerOfTuples(samplePoints); }
95
105 void setSaturatedOilFormationVolumeFactor(unsigned regionIdx,
106 const SamplingPoints& samplePoints);
107
120 void setInverseOilFormationVolumeFactor(unsigned regionIdx, const TabulatedTwoDFunction& invBo)
121 { inverseOilBTable_[regionIdx] = invBo; }
122
128 void setOilViscosity(unsigned regionIdx, const TabulatedTwoDFunction& muo)
129 { oilMuTable_[regionIdx] = muo; }
130
138 void setSaturatedOilViscosity(unsigned regionIdx,
139 const SamplingPoints& samplePoints);
140
144 void initEnd();
145
149 unsigned numRegions() const
150 { return inverseOilBMuTable_.size(); }
151
155 template <class Evaluation>
156 Evaluation internalEnergy(unsigned,
157 const Evaluation&,
158 const Evaluation&,
159 const Evaluation&) const
160 {
161 throw std::runtime_error("Requested the enthalpy of oil but the thermal option is not enabled");
162 }
163
167 template <class Evaluation>
168 Evaluation viscosity(unsigned regionIdx,
169 const Evaluation& /*temperature*/,
170 const Evaluation& pressure,
171 const Evaluation& Rs) const
172 {
173 // ATTENTION: Rs is the first axis!
174 const Evaluation& invBo = inverseOilBTable_[regionIdx].eval(Rs, pressure, /*extrapolate=*/true);
175 const Evaluation& invMuoBo = inverseOilBMuTable_[regionIdx].eval(Rs, pressure, /*extrapolate=*/true);
176
177 return invBo/invMuoBo;
178 }
179
183 template <class Evaluation>
184 Evaluation saturatedViscosity(unsigned regionIdx,
185 const Evaluation& /*temperature*/,
186 const Evaluation& pressure) const
187 {
188 // ATTENTION: Rs is the first axis!
189 const Evaluation& invBo = inverseSaturatedOilBTable_[regionIdx].eval(pressure, /*extrapolate=*/true);
190 const Evaluation& invMuoBo = inverseSaturatedOilBMuTable_[regionIdx].eval(pressure, /*extrapolate=*/true);
191
192 return invBo/invMuoBo;
193 }
194
198 template <class Evaluation>
199 Evaluation inverseFormationVolumeFactor(unsigned regionIdx,
200 const Evaluation& /*temperature*/,
201 const Evaluation& pressure,
202 const Evaluation& Rs) const
203 {
204 // ATTENTION: Rs is represented by the _first_ axis!
205 return inverseOilBTable_[regionIdx].eval(Rs, pressure, /*extrapolate=*/true);
206 }
207
211 template <class Evaluation>
212 Evaluation saturatedInverseFormationVolumeFactor(unsigned regionIdx,
213 const Evaluation& /*temperature*/,
214 const Evaluation& pressure) const
215 {
216 // ATTENTION: Rs is represented by the _first_ axis!
217 return inverseSaturatedOilBTable_[regionIdx].eval(pressure, /*extrapolate=*/true);
218 }
219
223 template <class Evaluation>
224 Evaluation saturatedGasDissolutionFactor(unsigned regionIdx,
225 const Evaluation& /*temperature*/,
226 const Evaluation& pressure) const
227 { return saturatedGasDissolutionFactorTable_[regionIdx].eval(pressure, /*extrapolate=*/true); }
228
236 template <class Evaluation>
237 Evaluation saturatedGasDissolutionFactor(unsigned regionIdx,
238 const Evaluation& /*temperature*/,
239 const Evaluation& pressure,
240 const Evaluation& oilSaturation,
241 Evaluation maxOilSaturation) const
242 {
243 Evaluation tmp =
244 saturatedGasDissolutionFactorTable_[regionIdx].eval(pressure, /*extrapolate=*/true);
245
246 // apply the vaporization parameters for the gas phase (cf. the Eclipse VAPPARS
247 // keyword)
248 maxOilSaturation = min(maxOilSaturation, Scalar(1.0));
249 if (vapPar2_ > 0.0 && maxOilSaturation > 0.01 && oilSaturation < maxOilSaturation) {
250 constexpr const Scalar eps = 0.001;
251 const Evaluation& So = max(oilSaturation, eps);
252 tmp *= max(1e-3, pow(So/maxOilSaturation, vapPar2_));
253 }
254
255 return tmp;
256 }
257
264 template <class Evaluation>
265 Evaluation saturationPressure(unsigned regionIdx,
266 const Evaluation&,
267 const Evaluation& Rs) const
268 {
269 typedef MathToolbox<Evaluation> Toolbox;
270
271 const auto& RsTable = saturatedGasDissolutionFactorTable_[regionIdx];
272 constexpr const Scalar eps = std::numeric_limits<typename Toolbox::Scalar>::epsilon()*1e6;
273
274 // use the saturation pressure function to get a pretty good initial value
275 Evaluation pSat = saturationPressure_[regionIdx].eval(Rs, /*extrapolate=*/true);
276
277 // Newton method to do the remaining work. If the initial
278 // value is good, this should only take two to three
279 // iterations...
280 bool onProbation = false;
281 for (int i = 0; i < 20; ++i) {
282 const Evaluation& f = RsTable.eval(pSat, /*extrapolate=*/true) - Rs;
283 const Evaluation& fPrime = RsTable.evalDerivative(pSat, /*extrapolate=*/true);
284
285 // If the derivative is "zero" Newton will not converge,
286 // so simply return our initial guess.
287 if (std::abs(scalarValue(fPrime)) < 1.0e-30) {
288 return pSat;
289 }
290
291 const Evaluation& delta = f/fPrime;
292
293 pSat -= delta;
294
295 if (pSat < 0.0) {
296 // if the pressure is lower than 0 Pascals, we set it back to 0. if this
297 // happens twice, we give up and just return 0 Pa...
298 if (onProbation)
299 return 0.0;
300
301 onProbation = true;
302 pSat = 0.0;
303 }
304
305 if (std::abs(scalarValue(delta)) < std::abs(scalarValue(pSat))*eps)
306 return pSat;
307 }
308
309 const std::string msg =
310 "Finding saturation pressure did not converge: "
311 " pSat = " + std::to_string(getValue(pSat)) +
312 ", Rs = " + std::to_string(getValue(Rs));
313 OpmLog::debug("Live oil saturation pressure", msg);
314 throw NumericalProblem(msg);
315 }
316
317 template <class Evaluation>
318 Evaluation diffusionCoefficient(const Evaluation& /*temperature*/,
319 const Evaluation& /*pressure*/,
320 unsigned /*compIdx*/) const
321 {
322 throw std::runtime_error("Not implemented: The PVT model does not provide a diffusionCoefficient()");
323 }
324
325 Scalar gasReferenceDensity(unsigned regionIdx) const
326 { return gasReferenceDensity_[regionIdx]; }
327
328 Scalar oilReferenceDensity(unsigned regionIdx) const
329 { return oilReferenceDensity_[regionIdx]; }
330
331 const std::vector<TabulatedTwoDFunction>& inverseOilBTable() const
332 { return inverseOilBTable_; }
333
334 const std::vector<TabulatedTwoDFunction>& oilMuTable() const
335 { return oilMuTable_; }
336
337 const std::vector<TabulatedTwoDFunction>& inverseOilBMuTable() const
338 { return inverseOilBMuTable_; }
339
340 const std::vector<TabulatedOneDFunction>& saturatedOilMuTable() const
341 { return saturatedOilMuTable_; }
342
343 const std::vector<TabulatedOneDFunction>& inverseSaturatedOilBTable() const
344 { return inverseSaturatedOilBTable_; }
345
346 const std::vector<TabulatedOneDFunction>& inverseSaturatedOilBMuTable() const
347 { return inverseSaturatedOilBMuTable_; }
348
349 const std::vector<TabulatedOneDFunction>& saturatedGasDissolutionFactorTable() const
350 { return saturatedGasDissolutionFactorTable_; }
351
352 const std::vector<TabulatedOneDFunction>& saturationPressure() const
353 { return saturationPressure_; }
354
355 Scalar vapPar2() const
356 { return vapPar2_; }
357
358private:
359 void updateSaturationPressure_(unsigned regionIdx);
360
361 std::vector<Scalar> gasReferenceDensity_;
362 std::vector<Scalar> oilReferenceDensity_;
363 std::vector<TabulatedTwoDFunction> inverseOilBTable_;
364 std::vector<TabulatedTwoDFunction> oilMuTable_;
365 std::vector<TabulatedTwoDFunction> inverseOilBMuTable_;
366 std::vector<TabulatedOneDFunction> saturatedOilMuTable_;
367 std::vector<TabulatedOneDFunction> inverseSaturatedOilBTable_;
368 std::vector<TabulatedOneDFunction> inverseSaturatedOilBMuTable_;
369 std::vector<TabulatedOneDFunction> saturatedGasDissolutionFactorTable_;
370 std::vector<TabulatedOneDFunction> saturationPressure_;
371
372 Scalar vapPar2_ = 0.0;
373};
374
375} // namespace Opm
376
377#endif
Provides the OPM specific exception classes.
A traits class which provides basic mathematical functions for arbitrary scalar floating point values...
Implements a linearly interpolated scalar function that depends on one variable.
Implements a scalar function that depends on two variables and which is sampled uniformly in the X di...
Definition EclipseState.hpp:63
This class represents the Pressure-Volume-Temperature relations of the oil phas with dissolved gas.
Definition LiveOilPvt.hpp:51
Evaluation internalEnergy(unsigned, const Evaluation &, const Evaluation &, const Evaluation &) const
Returns the specific enthalpy [J/kg] of oil given a set of parameters.
Definition LiveOilPvt.hpp:156
void setSaturatedOilFormationVolumeFactor(unsigned regionIdx, const SamplingPoints &samplePoints)
Initialize the function for the oil formation volume factor.
Definition LiveOilPvt.cpp:246
Evaluation inverseFormationVolumeFactor(unsigned regionIdx, const Evaluation &, const Evaluation &pressure, const Evaluation &Rs) const
Returns the formation volume factor [-] of the fluid phase.
Definition LiveOilPvt.hpp:199
Evaluation saturatedGasDissolutionFactor(unsigned regionIdx, const Evaluation &, const Evaluation &pressure, const Evaluation &oilSaturation, Evaluation maxOilSaturation) const
Returns the gas dissolution factor [m^3/m^3] of the oil phase.
Definition LiveOilPvt.hpp:237
Evaluation saturatedInverseFormationVolumeFactor(unsigned regionIdx, const Evaluation &, const Evaluation &pressure) const
Returns the formation volume factor [-] of the fluid phase.
Definition LiveOilPvt.hpp:212
void setSaturatedOilViscosity(unsigned regionIdx, const SamplingPoints &samplePoints)
Initialize the phase viscosity for gas saturated oil.
Definition LiveOilPvt.cpp:273
Evaluation saturationPressure(unsigned regionIdx, const Evaluation &, const Evaluation &Rs) const
Returns the saturation pressure of the oil phase [Pa] depending on its mass fraction of the gas compo...
Definition LiveOilPvt.hpp:265
void setSaturatedOilGasDissolutionFactor(unsigned regionIdx, const SamplingPoints &samplePoints)
Initialize the function for the gas dissolution factor .
Definition LiveOilPvt.hpp:93
Evaluation viscosity(unsigned regionIdx, const Evaluation &, const Evaluation &pressure, const Evaluation &Rs) const
Returns the dynamic viscosity [Pa s] of the fluid phase given a set of parameters.
Definition LiveOilPvt.hpp:168
void setReferenceDensities(unsigned regionIdx, Scalar rhoRefOil, Scalar rhoRefGas, Scalar)
Initialize the reference densities of all fluids for a given PVT region.
Definition LiveOilPvt.cpp:235
void setOilViscosity(unsigned regionIdx, const TabulatedTwoDFunction &muo)
Initialize the viscosity of the oil phase.
Definition LiveOilPvt.hpp:128
Evaluation saturatedGasDissolutionFactor(unsigned regionIdx, const Evaluation &, const Evaluation &pressure) const
Returns the gas dissolution factor [m^3/m^3] of the oil phase.
Definition LiveOilPvt.hpp:224
Evaluation saturatedViscosity(unsigned regionIdx, const Evaluation &, const Evaluation &pressure) const
Returns the dynamic viscosity [Pa s] of the fluid phase given a set of parameters.
Definition LiveOilPvt.hpp:184
unsigned numRegions() const
Return the number of PVT regions which are considered by this PVT-object.
Definition LiveOilPvt.hpp:149
void initEnd()
Finish initializing the oil phase PVT properties.
Definition LiveOilPvt.cpp:300
void setInverseOilFormationVolumeFactor(unsigned regionIdx, const TabulatedTwoDFunction &invBo)
Initialize the function for the oil formation volume factor.
Definition LiveOilPvt.hpp:120
Definition Exceptions.hpp:40
Definition Schedule.hpp:88
Definition SimpleTable.hpp:35
Implements a linearly interpolated scalar function that depends on one variable.
Definition Tabulated1DFunction.hpp:51
Implements a scalar function that depends on two variables and which is sampled uniformly in the X di...
Definition UniformXTabulated2DFunction.hpp:54
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30
Definition MathToolbox.hpp:50