My Project
Loading...
Searching...
No Matches
EquilibrationHelpers_impl.hpp
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 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 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*/
23
24#include <opm/common/TimingMacros.hpp>
25
26#include <opm/common/utility/numeric/RootFinders.hpp>
27
28#include <opm/material/fluidmatrixinteractions/EclMaterialLawManager.hpp>
29#include <opm/material/fluidstates/SimpleModularFluidState.hpp>
30#include <opm/material/fluidsystems/BlackOilFluidSystem.hpp>
31
33
34#include <fmt/format.h>
35
36namespace Opm {
37namespace EQUIL {
38
39using FluidSystemSimple = BlackOilFluidSystem<double>;
40
41// Adjust oil pressure according to gas saturation and cap pressure
42using SatOnlyFluidState = SimpleModularFluidState<double,
43 /*numPhases=*/3,
44 /*numComponents=*/3,
45 FluidSystemSimple,
46 /*storePressure=*/false,
47 /*storeTemperature=*/false,
48 /*storeComposition=*/false,
49 /*storeFugacity=*/false,
50 /*storeSaturation=*/true,
51 /*storeDensity=*/false,
52 /*storeViscosity=*/false,
53 /*storeEnthalpy=*/false>;
54
55namespace Miscibility {
56
57template<class FluidSystem>
58RsVD<FluidSystem>::RsVD(const int pvtRegionIdx,
59 const std::vector<double>& depth,
60 const std::vector<double>& rs)
61 : pvtRegionIdx_(pvtRegionIdx)
62 , rsVsDepth_(depth, rs)
63{
64}
65
66template<class FluidSystem>
68operator()(const double depth,
69 const double press,
70 const double temp,
71 const double satGas) const
72{
73 const auto sat_rs = satRs(press, temp);
74 if (satGas > std::sqrt(std::numeric_limits<double>::epsilon())) {
75 return sat_rs;
76 }
77 else {
78 if (rsVsDepth_.xMin() > depth)
79 return std::min(sat_rs, rsVsDepth_.valueAt(0));
80 else if (rsVsDepth_.xMax() < depth)
81 return std::min(sat_rs, rsVsDepth_.valueAt(rsVsDepth_.numSamples() - 1));
82 return std::min(sat_rs, rsVsDepth_.eval(depth, /*extrapolate=*/false));
83 }
84}
85
86template<class FluidSystem>
87double RsVD<FluidSystem>::satRs(const double press, const double temp) const
88{
89 return FluidSystem::oilPvt().saturatedGasDissolutionFactor(pvtRegionIdx_, temp, press);
90}
91
92template<class FluidSystem>
93PBVD<FluidSystem>::PBVD(const int pvtRegionIdx,
94 const std::vector<double>& depth,
95 const std::vector<double>& pbub)
96 : pvtRegionIdx_(pvtRegionIdx)
97 , pbubVsDepth_(depth, pbub)
98{
99}
100
101template<class FluidSystem>
103operator()(const double depth,
104 const double cellPress,
105 const double temp,
106 const double satGas) const
107{
108 double press = cellPress;
109 if (satGas <= 0.0) {
110 if (pbubVsDepth_.xMin() > depth)
111 press = pbubVsDepth_.valueAt(0);
112 else if (pbubVsDepth_.xMax() < depth)
113 press = pbubVsDepth_.valueAt(pbubVsDepth_.numSamples() - 1);
114 else
115 press = pbubVsDepth_.eval(depth, /*extrapolate=*/false);
116 }
117 return satRs(std::min(press, cellPress), temp);
118}
119
120template<class FluidSystem>
122satRs(const double press, const double temp) const
123{
124 return FluidSystem::oilPvt().saturatedGasDissolutionFactor(pvtRegionIdx_, temp, press);
125}
126
127template<class FluidSystem>
128PDVD<FluidSystem>::PDVD(const int pvtRegionIdx,
129 const std::vector<double>& depth,
130 const std::vector<double>& pdew)
131 : pvtRegionIdx_(pvtRegionIdx)
132 , pdewVsDepth_(depth, pdew)
133{
134}
135
136template<class FluidSystem>
138operator()(const double depth,
139 const double cellPress,
140 const double temp,
141 const double satOil) const
142{
143 double press = cellPress;
144 if (satOil <= 0.0) {
145 if (pdewVsDepth_.xMin() > depth)
146 press = pdewVsDepth_.valueAt(0);
147 else if (pdewVsDepth_.xMax() < depth)
148 press = pdewVsDepth_.valueAt(pdewVsDepth_.numSamples() - 1);
149 else
150 press = pdewVsDepth_.eval(depth, /*extrapolate=*/false);
151 }
152 return satRv(std::min(press, cellPress), temp);
153}
154
155template<class FluidSystem>
157satRv(const double press, const double temp) const
158{
159 return FluidSystem::gasPvt().saturatedOilVaporizationFactor(pvtRegionIdx_, temp, press);
160}
161
162
163template<class FluidSystem>
164RvVD<FluidSystem>::RvVD(const int pvtRegionIdx,
165 const std::vector<double>& depth,
166 const std::vector<double>& rv)
167 : pvtRegionIdx_(pvtRegionIdx)
168 , rvVsDepth_(depth, rv)
169{
170}
171
172template<class FluidSystem>
174operator()(const double depth,
175 const double press,
176 const double temp,
177 const double satOil) const
178{
179 if (satOil < - std::sqrt(std::numeric_limits<double>::epsilon())) {
180 throw std::logic_error {
181 "Must not pass negative oil saturation"
182 };
183 }
184 const auto sat_rv = satRv(press, temp);
185 if (satOil > std::sqrt(std::numeric_limits<double>::epsilon())) {
186 return sat_rv;
187 }
188 else {
189 if (rvVsDepth_.xMin() > depth)
190 return std::min(sat_rv, rvVsDepth_.valueAt(0));
191 else if (rvVsDepth_.xMax() < depth)
192 return std::min(sat_rv, rvVsDepth_.valueAt(rvVsDepth_.numSamples() - 1));
193 return std::min(sat_rv, rvVsDepth_.eval(depth, /*extrapolate=*/false));
194 }
195}
196
197template<class FluidSystem>
199satRv(const double press, const double temp) const
200{
201 return FluidSystem::gasPvt().saturatedOilVaporizationFactor(pvtRegionIdx_, temp, press);
202}
203
204
205template<class FluidSystem>
206RvwVD<FluidSystem>::RvwVD(const int pvtRegionIdx,
207 const std::vector<double>& depth,
208 const std::vector<double>& rvw)
209 : pvtRegionIdx_(pvtRegionIdx)
210 , rvwVsDepth_(depth, rvw)
211{
212}
213
214template<class FluidSystem>
216operator()(const double depth,
217 const double press,
218 const double temp,
219 const double satWat) const
220{
221 if (satWat < - std::sqrt(std::numeric_limits<double>::epsilon())) {
222 throw std::logic_error {
223 "Must not pass negative water saturation"
224 };
225 }
226
227 const auto sat_rvw = satRvw(press, temp);
228 if (satWat > std::sqrt(std::numeric_limits<double>::epsilon())) {
229 return sat_rvw; //saturated Rvw
230 }
231 else {
232 if (rvwVsDepth_.xMin() > depth)
233 return std::min(sat_rvw,rvwVsDepth_.valueAt(0));
234 else if (rvwVsDepth_.xMax() < depth)
235 return std::min(sat_rvw, rvwVsDepth_.valueAt(rvwVsDepth_.numSamples() - 1));
236 return std::min(sat_rvw, rvwVsDepth_.eval(depth, /*extrapolate=*/false));
237 }
238}
239
240template<class FluidSystem>
242satRvw(const double press, const double temp) const
243{
244 return FluidSystem::gasPvt().saturatedWaterVaporizationFactor(pvtRegionIdx_, temp, press);
245}
246
247
248template<class FluidSystem>
250RsSatAtContact(const int pvtRegionIdx, const double pContact, const double T_contact)
251 : pvtRegionIdx_(pvtRegionIdx)
252{
253 rsSatContact_ = satRs(pContact, T_contact);
254}
255
256template<class FluidSystem>
258operator()(const double /* depth */,
259 const double press,
260 const double temp,
261 const double satGas) const
262{
263 if (satGas > std::sqrt(std::numeric_limits<double>::epsilon())) {
264 return satRs(press, temp);
265 }
266 else {
267 return std::min(satRs(press, temp), rsSatContact_);
268 }
269}
270
271template<class FluidSystem>
273satRs(const double press, const double temp) const
274{
275 return FluidSystem::oilPvt().saturatedGasDissolutionFactor(pvtRegionIdx_, temp, press);
276}
277
278template<class FluidSystem>
280RvSatAtContact(const int pvtRegionIdx, const double pContact, const double T_contact)
281 : pvtRegionIdx_(pvtRegionIdx)
282{
283 rvSatContact_ = satRv(pContact, T_contact);
284}
285
286template<class FluidSystem>
288operator()(const double /*depth*/,
289 const double press,
290 const double temp,
291 const double satOil) const
292{
293 if (satOil > std::sqrt(std::numeric_limits<double>::epsilon())) {
294 return satRv(press, temp);
295 }
296 else {
297 return std::min(satRv(press, temp), rvSatContact_);
298 }
299}
300
301template<class FluidSystem>
303satRv(const double press, const double temp) const
304{
305 return FluidSystem::gasPvt().saturatedOilVaporizationFactor(pvtRegionIdx_, temp, press);;
306}
307
308template<class FluidSystem>
310RvwSatAtContact(const int pvtRegionIdx, const double pContact, const double T_contact)
311 : pvtRegionIdx_(pvtRegionIdx)
312{
313 rvwSatContact_ = satRvw(pContact, T_contact);
314}
315
316template<class FluidSystem>
318operator()(const double /*depth*/,
319 const double press,
320 const double temp,
321 const double satWat) const
322{
323 if (satWat > std::sqrt(std::numeric_limits<double>::epsilon())) {
324 return satRvw(press, temp);
325 }
326 else {
327 return std::min(satRvw(press, temp), rvwSatContact_);
328 }
329}
330
331template<class FluidSystem>
333satRvw(const double press, const double temp) const
334{
335 return FluidSystem::gasPvt().saturatedWaterVaporizationFactor(pvtRegionIdx_, temp, press);;
336}
337
338} // namespace Miscibility
339
340EquilReg::EquilReg(const EquilRecord& rec,
341 std::shared_ptr<Miscibility::RsFunction> rs,
342 std::shared_ptr<Miscibility::RsFunction> rv,
343 std::shared_ptr<Miscibility::RsFunction> rvw,
344 const TabulatedFunction& tempVdTable,
345 const TabulatedFunction& saltVdTable,
346 const int pvtIdx)
347 : rec_ (rec)
348 , rs_ (rs)
349 , rv_ (rv)
350 , rvw_ (rvw)
351 , tempVdTable_ (tempVdTable)
352 , saltVdTable_ (saltVdTable)
353 , pvtIdx_ (pvtIdx)
354{
355}
356
357double EquilReg::datum() const
358{
359 return this->rec_.datumDepth();
360}
361
362double EquilReg::pressure() const
363{
364 return this->rec_.datumDepthPressure();
365}
366
367double EquilReg::zwoc() const
368{
369 return this->rec_.waterOilContactDepth();
370}
371
372double EquilReg::pcowWoc() const
373{
374 return this->rec_.waterOilContactCapillaryPressure();
375}
376
377double EquilReg::zgoc() const
378{
379 return this->rec_.gasOilContactDepth();
380}
381
382double EquilReg::pcgoGoc() const
383{
384 return this->rec_.gasOilContactCapillaryPressure();
385}
386
388{
389 return this->rec_.initializationTargetAccuracy();
390}
391
394{
395 return *this->rs_;
396}
397
400{
401 return *this->rv_;
402}
403
406{
407 return *this->rvw_;
408}
409
410const EquilReg::TabulatedFunction&
411EquilReg::saltVdTable() const
412{
413 return saltVdTable_;
414}
415
416const EquilReg::TabulatedFunction&
417EquilReg::tempVdTable() const
418{
419 return tempVdTable_;
420}
421
423{
424 return this->pvtIdx_;
425}
426
427template<class FluidSystem, class MaterialLawManager>
429PcEq(const MaterialLawManager& materialLawManager,
430 const int phase,
431 const int cell,
432 const double targetPc)
433 : materialLawManager_(materialLawManager),
434 phase_(phase),
435 cell_(cell),
436 targetPc_(targetPc)
437{
438}
439
440template<class FluidSystem, class MaterialLawManager>
441double PcEq<FluidSystem,MaterialLawManager>::
442operator()(double s) const
443{
444 const auto& matParams = materialLawManager_.materialLawParams(cell_);
445 SatOnlyFluidState fluidState;
446 fluidState.setSaturation(FluidSystem::waterPhaseIdx, 0.0);
447 fluidState.setSaturation(FluidSystem::oilPhaseIdx, 0.0);
448 fluidState.setSaturation(FluidSystem::gasPhaseIdx, 0.0);
449 fluidState.setSaturation(phase_, s);
450
451 std::array<double, FluidSystem::numPhases> pc{0.0};
452 using MaterialLaw = typename MaterialLawManager::MaterialLaw;
453 MaterialLaw::capillaryPressures(pc, matParams, fluidState);
454 double sign = (phase_ == FluidSystem::waterPhaseIdx)? -1.0 : 1.0;
455 double pcPhase = pc[FluidSystem::oilPhaseIdx] + sign * pc[phase_];
456 return pcPhase - targetPc_;
457}
458
459template<class FluidSystem, class MaterialLawManager>
460PcEqSum<FluidSystem,MaterialLawManager>::
461PcEqSum(const MaterialLawManager& materialLawManager,
462 const int phase1,
463 const int phase2,
464 const int cell,
465 const double targetPc)
466 : materialLawManager_(materialLawManager),
467 phase1_(phase1),
468 phase2_(phase2),
469 cell_(cell),
470 targetPc_(targetPc)
471{
472}
473
474template<class FluidSystem, class MaterialLawManager>
475double PcEqSum<FluidSystem,MaterialLawManager>::
476operator()(double s) const
477{
478 const auto& matParams = materialLawManager_.materialLawParams(cell_);
479 SatOnlyFluidState fluidState;
480 fluidState.setSaturation(FluidSystem::waterPhaseIdx, 0.0);
481 fluidState.setSaturation(FluidSystem::oilPhaseIdx, 0.0);
482 fluidState.setSaturation(FluidSystem::gasPhaseIdx, 0.0);
483 fluidState.setSaturation(phase1_, s);
484 fluidState.setSaturation(phase2_, 1.0 - s);
485
486 std::array<double, FluidSystem::numPhases> pc {0.0};
487
488 using MaterialLaw = typename MaterialLawManager::MaterialLaw;
489 MaterialLaw::capillaryPressures(pc, matParams, fluidState);
490 double sign1 = (phase1_ == FluidSystem::waterPhaseIdx)? -1.0 : 1.0;
491 double pc1 = pc[FluidSystem::oilPhaseIdx] + sign1 * pc[phase1_];
492 double sign2 = (phase2_ == FluidSystem::waterPhaseIdx)? -1.0 : 1.0;
493 double pc2 = pc[FluidSystem::oilPhaseIdx] + sign2 * pc[phase2_];
494 return pc1 + pc2 - targetPc_;
495}
496
497template <class FluidSystem, class MaterialLawManager>
498double minSaturations(const MaterialLawManager& materialLawManager,
499 const int phase, const int cell)
500{
501 const auto& scaledDrainageInfo =
502 materialLawManager.oilWaterScaledEpsInfoDrainage(cell);
503
504 // Find minimum and maximum saturations.
505 switch(phase) {
506 case FluidSystem::waterPhaseIdx:
507 return scaledDrainageInfo.Swl;
508
509 case FluidSystem::gasPhaseIdx:
510 return scaledDrainageInfo.Sgl;
511
512 case FluidSystem::oilPhaseIdx:
513 throw std::runtime_error("Min saturation not implemented for oil phase.");
514
515 default:
516 throw std::runtime_error("Unknown phaseIdx .");
517 }
518 return -1.0;
519}
520
521template <class FluidSystem, class MaterialLawManager>
522double maxSaturations(const MaterialLawManager& materialLawManager,
523 const int phase, const int cell)
524{
525 const auto& scaledDrainageInfo =
526 materialLawManager.oilWaterScaledEpsInfoDrainage(cell);
527
528 // Find minimum and maximum saturations.
529 switch(phase) {
530 case FluidSystem::waterPhaseIdx:
531 return scaledDrainageInfo.Swu;
532
533 case FluidSystem::gasPhaseIdx:
534 return scaledDrainageInfo.Sgu;
535
536 case FluidSystem::oilPhaseIdx:
537 throw std::runtime_error("Max saturation not implemented for oil phase.");
538
539 default:
540 throw std::runtime_error("Unknown phaseIdx .");
541 }
542 return -1.0;
543}
544
545template <class FluidSystem, class MaterialLawManager>
546double satFromPc(const MaterialLawManager& materialLawManager,
547 const int phase,
548 const int cell,
549 const double targetPc,
550 const bool increasing)
551{
552 // Find minimum and maximum saturations.
553 double s0 = increasing ? maxSaturations<FluidSystem>(materialLawManager, phase, cell) : minSaturations<FluidSystem>(materialLawManager, phase, cell);
554 double s1 = increasing ? minSaturations<FluidSystem>(materialLawManager, phase, cell) : maxSaturations<FluidSystem>(materialLawManager, phase, cell);
555
556 // Create the equation f(s) = pc(s) - targetPc
557 const PcEq<FluidSystem, MaterialLawManager> f(materialLawManager, phase, cell, targetPc);
558 double f0 = f(s0);
559 double f1 = f(s1);
560 if (!std::isfinite(f0 + f1))
561 throw std::logic_error(fmt::format("The capillary pressure values {} and {} are not finite", f0, f1));
562
563 if (f0 <= 0.0)
564 return s0;
565 else if (f1 >= 0.0)
566 return s1;
567
568 const double tol = 1e-10;
569 // should at least converge in 2 times bisection but some safety here:
570 const int maxIter = -2*static_cast<int>(std::log2(tol)) + 10;
571 int usedIterations = -1;
572 const double root = RegulaFalsiBisection<ThrowOnError>::solve(f, s0, s1, maxIter, tol, usedIterations);
573 return root;
574}
575
576template<class FluidSystem, class MaterialLawManager>
577double satFromSumOfPcs(const MaterialLawManager& materialLawManager,
578 const int phase1,
579 const int phase2,
580 const int cell,
581 const double targetPc)
582{
583 // Find minimum and maximum saturations.
584 double s0 = minSaturations<FluidSystem>(materialLawManager, phase1, cell);
585 double s1 = maxSaturations<FluidSystem>(materialLawManager, phase1, cell);
586
587 // Create the equation f(s) = pc1(s) + pc2(1-s) - targetPc
588 const PcEqSum<FluidSystem, MaterialLawManager> f(materialLawManager, phase1, phase2, cell, targetPc);
589 double f0 = f(s0);
590 double f1 = f(s1);
591 if (f0 <= 0.0)
592 return s0;
593 else if (f1 >= 0.0)
594 return s1;
595
596 assert(f0 > 0.0 && f1 < 0.0);
597 const double tol = 1e-10;
598 // should at least converge in 2 times bisection but some safety here:
599 const int maxIter = -2*static_cast<int>(std::log2(tol)) + 10;
600 int usedIterations = -1;
601 const double root = RegulaFalsiBisection<ThrowOnError>::solve(f, s0, s1, maxIter, tol, usedIterations);
602 return root;
603}
604
605template<class FluidSystem, class MaterialLawManager>
606double satFromDepth(const MaterialLawManager& materialLawManager,
607 const double cellDepth,
608 const double contactDepth,
609 const int phase,
610 const int cell,
611 const bool increasing)
612{
613 const double s0 = increasing ? maxSaturations<FluidSystem>(materialLawManager, phase, cell) : minSaturations<FluidSystem>(materialLawManager, phase, cell);
614 const double s1 = increasing ? minSaturations<FluidSystem>(materialLawManager, phase, cell) : maxSaturations<FluidSystem>(materialLawManager, phase, cell);
615
616 if (cellDepth < contactDepth) {
617 return s0;
618 }
619 else {
620 return s1;
621 }
622}
623
624template<class FluidSystem, class MaterialLawManager>
625bool isConstPc(const MaterialLawManager& materialLawManager,
626 const int phase,
627 const int cell)
628{
629 // Create the equation f(s) = pc(s);
630 const PcEq<FluidSystem, MaterialLawManager> f(materialLawManager, phase, cell, 0);
631 const double f0 = f(minSaturations<FluidSystem>(materialLawManager, phase, cell));
632 const double f1 = f(maxSaturations<FluidSystem>(materialLawManager, phase, cell));
633 return std::abs(f0 - f1) < std::numeric_limits<double>::epsilon();
634}
635
636}
637}
Auxiliary routines that to solve the ODEs that emerge from the hydrostatic equilibrium problem.
const CalcDissolution & dissolutionCalculator() const
Retrieve dissolved gas-oil ratio calculator of current region.
Definition EquilibrationHelpers_impl.hpp:393
double datum() const
Datum depth in current region.
Definition EquilibrationHelpers_impl.hpp:357
double pcgoGoc() const
Gas-oil capillary pressure at gas-oil contact.
Definition EquilibrationHelpers_impl.hpp:382
int equilibrationAccuracy() const
Accuracy/strategy for initial fluid-in-place calculation.
Definition EquilibrationHelpers_impl.hpp:387
double zgoc() const
Depth of gas-oil contact.
Definition EquilibrationHelpers_impl.hpp:377
const CalcWaterEvaporation & waterEvaporationCalculator() const
Retrieve vapourised water-gas ratio calculator of current region.
Definition EquilibrationHelpers_impl.hpp:405
const CalcEvaporation & evaporationCalculator() const
Retrieve vapourised oil-gas ratio calculator of current region.
Definition EquilibrationHelpers_impl.hpp:399
double pcowWoc() const
water-oil capillary pressure at water-oil contact.
Definition EquilibrationHelpers_impl.hpp:372
double zwoc() const
Depth of water-oil contact.
Definition EquilibrationHelpers_impl.hpp:367
double pressure() const
Pressure at datum depth in current region.
Definition EquilibrationHelpers_impl.hpp:362
int pvtIdx() const
Retrieve pvtIdx of the region.
Definition EquilibrationHelpers_impl.hpp:422
EquilReg(const EquilRecord &rec, std::shared_ptr< Miscibility::RsFunction > rs, std::shared_ptr< Miscibility::RsFunction > rv, std::shared_ptr< Miscibility::RsFunction > rvw, const TabulatedFunction &tempVdTable, const TabulatedFunction &saltVdTable, const int pvtIdx)
Constructor.
Definition EquilibrationHelpers_impl.hpp:340
Type that implements "dissolved gas-oil ratio" tabulated as a function of depth policy.
Definition EquilibrationHelpers.hpp:220
double operator()(const double depth, const double cellPress, const double temp, const double satGas=0.0) const
Function call.
Definition EquilibrationHelpers_impl.hpp:103
PBVD(const int pvtRegionIdx, const std::vector< double > &depth, const std::vector< double > &pbub)
Constructor.
Definition EquilibrationHelpers_impl.hpp:93
Type that implements "vaporized oil-gas ratio" tabulated as a function of depth policy.
Definition EquilibrationHelpers.hpp:271
PDVD(const int pvtRegionIdx, const std::vector< double > &depth, const std::vector< double > &pdew)
Constructor.
Definition EquilibrationHelpers_impl.hpp:128
double operator()(const double depth, const double cellPress, const double temp, const double satOil=0.0) const
Function call.
Definition EquilibrationHelpers_impl.hpp:138
Base class for phase mixing functions.
Definition EquilibrationHelpers.hpp:100
Class that implements "dissolved gas-oil ratio" (Rs) as function of depth and pressure as follows:
Definition EquilibrationHelpers.hpp:431
RsSatAtContact(const int pvtRegionIdx, const double pContact, const double T_contact)
Constructor.
Definition EquilibrationHelpers_impl.hpp:250
double operator()(const double, const double press, const double temp, const double satGas=0.0) const
Function call.
Definition EquilibrationHelpers_impl.hpp:258
Type that implements "dissolved gas-oil ratio" tabulated as a function of depth policy.
Definition EquilibrationHelpers.hpp:168
RsVD(const int pvtRegionIdx, const std::vector< double > &depth, const std::vector< double > &rs)
Constructor.
Definition EquilibrationHelpers_impl.hpp:58
double operator()(const double depth, const double press, const double temp, const double satGas=0.0) const
Function call.
Definition EquilibrationHelpers_impl.hpp:68
Class that implements "vaporized oil-gas ratio" (Rv) as function of depth and pressure as follows:
Definition EquilibrationHelpers.hpp:486
RvSatAtContact(const int pvtRegionIdx, const double pContact, const double T_contact)
Constructor.
Definition EquilibrationHelpers_impl.hpp:280
double operator()(const double, const double press, const double temp, const double satOil=0.0) const
Function call.
Definition EquilibrationHelpers_impl.hpp:288
Type that implements "vaporized oil-gas ratio" tabulated as a function of depth policy.
Definition EquilibrationHelpers.hpp:322
double operator()(const double depth, const double press, const double temp, const double satOil=0.0) const
Function call.
Definition EquilibrationHelpers_impl.hpp:174
RvVD(const int pvtRegionIdx, const std::vector< double > &depth, const std::vector< double > &rv)
Constructor.
Definition EquilibrationHelpers_impl.hpp:164
Class that implements "vaporized water-gas ratio" (Rvw) as function of depth and pressure as follows:
Definition EquilibrationHelpers.hpp:540
RvwSatAtContact(const int pvtRegionIdx, const double pContact, const double T_contact)
Constructor.
Definition EquilibrationHelpers_impl.hpp:310
double operator()(const double, const double press, const double temp, const double satWat=0.0) const
Function call.
Definition EquilibrationHelpers_impl.hpp:318
Type that implements "vaporized water-gas ratio" tabulated as a function of depth policy.
Definition EquilibrationHelpers.hpp:372
double operator()(const double depth, const double press, const double temp, const double satWat=0.0) const
Function call.
Definition EquilibrationHelpers_impl.hpp:216
RvwVD(const int pvtRegionIdx, const std::vector< double > &depth, const std::vector< double > &rvw)
Constructor.
Definition EquilibrationHelpers_impl.hpp:206
double satFromDepth(const MaterialLawManager &materialLawManager, const double cellDepth, const double contactDepth, const int phase, const int cell, const bool increasing=false)
Compute saturation from depth. Used for constant capillary pressure function.
Definition EquilibrationHelpers_impl.hpp:606
double satFromSumOfPcs(const MaterialLawManager &materialLawManager, const int phase1, const int phase2, const int cell, const double targetPc)
Compute saturation of some phase corresponding to a given capillary pressure, where the capillary pre...
Definition EquilibrationHelpers_impl.hpp:577
double satFromPc(const MaterialLawManager &materialLawManager, const int phase, const int cell, const double targetPc, const bool increasing=false)
Compute saturation of some phase corresponding to a given capillary pressure.
Definition EquilibrationHelpers_impl.hpp:546
bool isConstPc(const MaterialLawManager &materialLawManager, const int phase, const int cell)
Return true if capillary pressure function is constant.
Definition EquilibrationHelpers_impl.hpp:625
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition BlackoilPhases.hpp:27
Functor for inverting a sum of capillary pressure functions.
Definition EquilibrationHelpers.hpp:760
Functor for inverting capillary pressure function.
Definition EquilibrationHelpers.hpp:723