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,
94 const SamplingPoints& samplePoints)
95 { saturatedGasDissolutionFactorTable_[regionIdx].setContainerOfTuples(samplePoints); }
96
106 void setSaturatedOilFormationVolumeFactor(unsigned regionIdx,
107 const SamplingPoints& samplePoints);
108
121 void setInverseOilFormationVolumeFactor(unsigned regionIdx,
122 const TabulatedTwoDFunction& invBo)
123 { inverseOilBTable_[regionIdx] = invBo; }
124
130 void setOilViscosity(unsigned regionIdx, const TabulatedTwoDFunction& muo)
131 { oilMuTable_[regionIdx] = muo; }
132
140 void setSaturatedOilViscosity(unsigned regionIdx,
141 const SamplingPoints& samplePoints);
142
146 void initEnd();
147
151 unsigned numRegions() const
152 { return inverseOilBMuTable_.size(); }
153
157 template <class Evaluation>
158 Evaluation internalEnergy(unsigned,
159 const Evaluation&,
160 const Evaluation&,
161 const Evaluation&) const
162 {
163 throw std::runtime_error("Requested the enthalpy of oil but the thermal "
164 "option is not enabled");
165 }
166
167 Scalar hVap(unsigned) const
168 {
169 throw std::runtime_error("Requested the hvap of oil but the thermal "
170 "option is not enabled");
171 }
172
176 template <class Evaluation>
177 Evaluation viscosity(unsigned regionIdx,
178 const Evaluation& /*temperature*/,
179 const Evaluation& pressure,
180 const Evaluation& Rs) const
181 {
182 // ATTENTION: Rs is the first axis!
183 const Evaluation& invBo = inverseOilBTable_[regionIdx].eval(Rs, pressure, /*extrapolate=*/true);
184 const Evaluation& invMuoBo = inverseOilBMuTable_[regionIdx].eval(Rs, pressure, /*extrapolate=*/true);
185
186 return invBo / invMuoBo;
187 }
188
192 template <class Evaluation>
193 Evaluation saturatedViscosity(unsigned regionIdx,
194 const Evaluation& /*temperature*/,
195 const Evaluation& pressure) const
196 {
197 // ATTENTION: Rs is the first axis!
198 const Evaluation& invBo = inverseSaturatedOilBTable_[regionIdx].eval(pressure, /*extrapolate=*/true);
199 const Evaluation& invMuoBo = inverseSaturatedOilBMuTable_[regionIdx].eval(pressure, /*extrapolate=*/true);
200
201 return invBo / invMuoBo;
202 }
203
207 template <class Evaluation>
208 Evaluation inverseFormationVolumeFactor(unsigned regionIdx,
209 const Evaluation& /*temperature*/,
210 const Evaluation& pressure,
211 const Evaluation& Rs) const
212 {
213 // ATTENTION: Rs is represented by the _first_ axis!
214 return inverseOilBTable_[regionIdx].eval(Rs, pressure, /*extrapolate=*/true);
215 }
216
220 template <class Evaluation>
221 Evaluation saturatedInverseFormationVolumeFactor(unsigned regionIdx,
222 const Evaluation& /*temperature*/,
223 const Evaluation& pressure) const
224 {
225 // ATTENTION: Rs is represented by the _first_ axis!
226 return inverseSaturatedOilBTable_[regionIdx].eval(pressure, /*extrapolate=*/true);
227 }
228
232 template <class Evaluation>
233 Evaluation saturatedGasDissolutionFactor(unsigned regionIdx,
234 const Evaluation& /*temperature*/,
235 const Evaluation& pressure) const
236 { return saturatedGasDissolutionFactorTable_[regionIdx].eval(pressure, /*extrapolate=*/true); }
237
245 template <class Evaluation>
246 Evaluation saturatedGasDissolutionFactor(unsigned regionIdx,
247 const Evaluation& /*temperature*/,
248 const Evaluation& pressure,
249 const Evaluation& oilSaturation,
250 Evaluation maxOilSaturation) const
251 {
252 Evaluation tmp =
253 saturatedGasDissolutionFactorTable_[regionIdx].eval(pressure, /*extrapolate=*/true);
254
255 // apply the vaporization parameters for the gas phase (cf. the Eclipse VAPPARS
256 // keyword)
257 maxOilSaturation = min(maxOilSaturation, Scalar(1.0));
258 if (vapPar2_ > 0.0 && maxOilSaturation > 0.01 && oilSaturation < maxOilSaturation) {
259 constexpr const Scalar eps = 0.001;
260 const Evaluation& So = max(oilSaturation, eps);
261 tmp *= max(1e-3, pow(So / maxOilSaturation, vapPar2_));
262 }
263
264 return tmp;
265 }
266
273 template <class Evaluation>
274 Evaluation saturationPressure(unsigned regionIdx,
275 const Evaluation&,
276 const Evaluation& Rs) const
277 {
278 using Toolbox = MathToolbox<Evaluation>;
279
280 const auto& RsTable = saturatedGasDissolutionFactorTable_[regionIdx];
281 constexpr const Scalar eps = std::numeric_limits<typename Toolbox::Scalar>::epsilon() * 1e6;
282
283 // use the saturation pressure function to get a pretty good initial value
284 Evaluation pSat = saturationPressure_[regionIdx].eval(Rs, /*extrapolate=*/true);
285
286 // Newton method to do the remaining work. If the initial
287 // value is good, this should only take two to three
288 // iterations...
289 bool onProbation = false;
290 for (int i = 0; i < 20; ++i) {
291 const Evaluation& f = RsTable.eval(pSat, /*extrapolate=*/true) - Rs;
292 const Evaluation& fPrime = RsTable.evalDerivative(pSat, /*extrapolate=*/true);
293
294 // If the derivative is "zero" Newton will not converge,
295 // so simply return our initial guess.
296 if (std::abs(scalarValue(fPrime)) < 1.0e-30) {
297 return pSat;
298 }
299
300 const Evaluation& delta = f / fPrime;
301
302 pSat -= delta;
303
304 if (pSat < 0.0) {
305 // if the pressure is lower than 0 Pascals, we set it back to 0. if this
306 // happens twice, we give up and just return 0 Pa...
307 if (onProbation) {
308 return 0.0;
309 }
310
311 onProbation = true;
312 pSat = 0.0;
313 }
314
315 if (std::abs(scalarValue(delta)) < std::abs(scalarValue(pSat))*eps) {
316 return pSat;
317 }
318 }
319
320 const std::string msg =
321 "Finding saturation pressure did not converge: "
322 " pSat = " + std::to_string(getValue(pSat)) +
323 ", Rs = " + std::to_string(getValue(Rs));
324 OpmLog::debug("Live oil saturation pressure", msg);
325 throw NumericalProblem(msg);
326 }
327
328 template <class Evaluation>
329 Evaluation diffusionCoefficient(const Evaluation& /*temperature*/,
330 const Evaluation& /*pressure*/,
331 unsigned /*compIdx*/) const
332 {
333 throw std::runtime_error("Not implemented: The PVT model does not provide "
334 "a diffusionCoefficient()");
335 }
336
337 Scalar gasReferenceDensity(unsigned regionIdx) const
338 { return gasReferenceDensity_[regionIdx]; }
339
340 Scalar oilReferenceDensity(unsigned regionIdx) const
341 { return oilReferenceDensity_[regionIdx]; }
342
343 const std::vector<TabulatedTwoDFunction>& inverseOilBTable() const
344 { return inverseOilBTable_; }
345
346 const std::vector<TabulatedTwoDFunction>& oilMuTable() const
347 { return oilMuTable_; }
348
349 const std::vector<TabulatedTwoDFunction>& inverseOilBMuTable() const
350 { return inverseOilBMuTable_; }
351
352 const std::vector<TabulatedOneDFunction>& saturatedOilMuTable() const
353 { return saturatedOilMuTable_; }
354
355 const std::vector<TabulatedOneDFunction>& inverseSaturatedOilBTable() const
356 { return inverseSaturatedOilBTable_; }
357
358 const std::vector<TabulatedOneDFunction>& inverseSaturatedOilBMuTable() const
359 { return inverseSaturatedOilBMuTable_; }
360
361 const std::vector<TabulatedOneDFunction>& saturatedGasDissolutionFactorTable() const
362 { return saturatedGasDissolutionFactorTable_; }
363
364 const std::vector<TabulatedOneDFunction>& saturationPressure() const
365 { return saturationPressure_; }
366
367 Scalar vapPar2() const
368 { return vapPar2_; }
369
370private:
371 void updateSaturationPressure_(unsigned regionIdx);
372
373 std::vector<Scalar> gasReferenceDensity_{};
374 std::vector<Scalar> oilReferenceDensity_{};
375 std::vector<TabulatedTwoDFunction> inverseOilBTable_{};
376 std::vector<TabulatedTwoDFunction> oilMuTable_{};
377 std::vector<TabulatedTwoDFunction> inverseOilBMuTable_{};
378 std::vector<TabulatedOneDFunction> saturatedOilMuTable_{};
379 std::vector<TabulatedOneDFunction> inverseSaturatedOilBTable_{};
380 std::vector<TabulatedOneDFunction> inverseSaturatedOilBMuTable_{};
381 std::vector<TabulatedOneDFunction> saturatedGasDissolutionFactorTable_{};
382 std::vector<TabulatedOneDFunction> saturationPressure_{};
383
384 Scalar vapPar2_ = 0.0;
385};
386
387} // namespace Opm
388
389#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:62
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:158
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:208
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:246
Evaluation saturatedInverseFormationVolumeFactor(unsigned regionIdx, const Evaluation &, const Evaluation &pressure) const
Returns the formation volume factor [-] of the fluid phase.
Definition LiveOilPvt.hpp:221
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:274
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:177
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:130
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:233
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:193
unsigned numRegions() const
Return the number of PVT regions which are considered by this PVT-object.
Definition LiveOilPvt.hpp:151
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:121
Definition Exceptions.hpp:40
Definition Schedule.hpp:101
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:55
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30
Definition MathToolbox.hpp:51