20 #ifndef OPM_EQUILIBRATIONHELPERS_HEADER_INCLUDED 21 #define OPM_EQUILIBRATIONHELPERS_HEADER_INCLUDED 23 #include <opm/core/props/BlackoilPropertiesInterface.hpp> 24 #include <opm/core/props/BlackoilPhases.hpp> 25 #include <opm/core/utility/linearInterpolation.hpp> 26 #include <opm/core/utility/RegionMapping.hpp> 27 #include <opm/core/utility/RootFinders.hpp> 29 #include <opm/parser/eclipse/EclipseState/InitConfig/Equil.hpp> 90 template <
class Props>
134 const std::vector<double>& z)
const 136 const int np = props_.numPhases();
137 std::vector<double> A(np * np, 0);
139 assert (z.size() == std::vector<double>::size_type(np));
142 props_.matrix(1, &p, &T, &z[0], &c_[0], &A[0], dAdp);
144 std::vector<double> rho(np, 0.0);
145 props_.density(1, &A[0], &c_[0], &rho[0]);
152 const std::vector<int> c_;
160 namespace Miscibility {
186 const double sat = 0.0)
const = 0;
215 const double = 0.0)
const 239 const std::vector<double>& depth,
240 const std::vector<double>& rs)
247 std::fill(z_, z_ + BlackoilPhases::MaxNumPhases, 0.0);
248 z_[pu.phase_pos[BlackoilPhases::Vapour]] = 1e100;
249 z_[pu.phase_pos[BlackoilPhases::Liquid]] = 1.0;
271 const double sat_gas = 0.0)
const 274 return satRs(press, temp);
276 return std::min(satRs(press, temp), linearInterpolationNoExtrapolation(depth_, rs_, depth));
283 std::vector<double> depth_;
284 std::vector<double> rs_;
285 double z_[BlackoilPhases::MaxNumPhases];
286 mutable double A_[BlackoilPhases::MaxNumPhases * BlackoilPhases::MaxNumPhases];
288 double satRs(
const double press,
const double temp)
const 290 props_.
matrix(1, &press, &temp, z_, &cell_, A_, 0);
294 const int opos = props_.
phaseUsage().phase_pos[BlackoilPhases::Liquid];
295 const int gpos = props_.
phaseUsage().phase_pos[BlackoilPhases::Vapour];
297 return A_[np*opos + gpos] / A_[np*opos + opos];
319 const std::vector<double>& depth,
320 const std::vector<double>& rv)
327 std::fill(z_, z_ + BlackoilPhases::MaxNumPhases, 0.0);
328 z_[pu.phase_pos[BlackoilPhases::Vapour]] = 1.0;
329 z_[pu.phase_pos[BlackoilPhases::Liquid]] = 1e100;
351 const double sat_oil = 0.0 )
const 353 if (std::abs(sat_oil) > 1e-16) {
354 return satRv(press, temp);
356 return std::min(satRv(press, temp), linearInterpolationNoExtrapolation(depth_, rv_, depth));
363 std::vector<double> depth_;
364 std::vector<double> rv_;
365 double z_[BlackoilPhases::MaxNumPhases];
366 mutable double A_[BlackoilPhases::MaxNumPhases * BlackoilPhases::MaxNumPhases];
368 double satRv(
const double press,
const double temp)
const 370 props_.
matrix(1, &press, &temp, z_, &cell_, A_, 0);
374 const int opos = props_.
phaseUsage().phase_pos[BlackoilPhases::Liquid];
375 const int gpos = props_.
phaseUsage().phase_pos[BlackoilPhases::Vapour];
377 return A_[np*gpos + opos] / A_[np*gpos + gpos];
407 : props_(props), cell_(cell)
410 std::fill(z_, z_ + BlackoilPhases::MaxNumPhases, 0.0);
411 z_[pu.phase_pos[BlackoilPhases::Vapour]] = 1e100;
412 z_[pu.phase_pos[BlackoilPhases::Liquid]] = 1.0;
413 rs_sat_contact_ = satRs(p_contact, T_contact);
435 const double sat_gas = 0.0)
const 438 return satRs(press, temp);
440 return std::min(satRs(press, temp), rs_sat_contact_);
447 double z_[BlackoilPhases::MaxNumPhases];
448 double rs_sat_contact_;
449 mutable double A_[BlackoilPhases::MaxNumPhases * BlackoilPhases::MaxNumPhases];
451 double satRs(
const double press,
const double temp)
const 453 props_.
matrix(1, &press, &temp, z_, &cell_, A_, 0);
457 const int opos = props_.
phaseUsage().phase_pos[BlackoilPhases::Liquid];
458 const int gpos = props_.
phaseUsage().phase_pos[BlackoilPhases::Vapour];
460 return A_[np*opos + gpos] / A_[np*opos + opos];
490 : props_(props), cell_(cell)
493 std::fill(z_, z_ + BlackoilPhases::MaxNumPhases, 0.0);
494 z_[pu.phase_pos[BlackoilPhases::Vapour]] = 1.0;
495 z_[pu.phase_pos[BlackoilPhases::Liquid]] = 1e100;
496 rv_sat_contact_ = satRv(p_contact, T_contact);
518 const double sat_oil = 0.0)
const 521 return satRv(press, temp);
523 return std::min(satRv(press, temp), rv_sat_contact_);
530 double z_[BlackoilPhases::MaxNumPhases];
531 double rv_sat_contact_;
532 mutable double A_[BlackoilPhases::MaxNumPhases * BlackoilPhases::MaxNumPhases];
534 double satRv(
const double press,
const double temp)
const 536 props_.
matrix(1, &press, &temp, z_, &cell_, A_, 0);
540 const int opos = props_.
phaseUsage().phase_pos[BlackoilPhases::Liquid];
541 const int gpos = props_.
phaseUsage().phase_pos[BlackoilPhases::Vapour];
543 return A_[np*gpos + opos] / A_[np*gpos + gpos];
568 template <
class DensCalc>
581 const DensCalc& density,
582 std::shared_ptr<Miscibility::RsFunction> rs,
583 std::shared_ptr<Miscibility::RsFunction> rv,
611 double datum()
const {
return this->rec_.datumDepth(); }
616 double pressure()
const {
return this->rec_.datumDepthPressure(); }
621 double zwoc()
const {
return this->rec_.waterOilContactDepth(); }
628 double pcow_woc()
const {
return this->rec_.waterOilContactCapillaryPressure(); }
633 double zgoc()
const {
return this->rec_.gasOilContactDepth(); }
640 double pcgo_goc()
const {
return this->rec_.gasOilContactCapillaryPressure(); }
671 std::shared_ptr<Miscibility::RsFunction> rs_;
672 std::shared_ptr<Miscibility::RsFunction> rv_;
686 const double target_pc)
690 target_pc_(target_pc)
692 std::fill(s_, s_ + BlackoilPhases::MaxNumPhases, 0.0);
693 std::fill(pc_, pc_ + BlackoilPhases::MaxNumPhases, 0.0);
695 double operator()(
double s)
const 698 props_.
capPress(1, s_, &cell_, pc_, 0);
699 return pc_[phase_] - target_pc_;
705 const int target_pc_;
706 mutable double s_[BlackoilPhases::MaxNumPhases];
707 mutable double pc_[BlackoilPhases::MaxNumPhases];
717 const double target_pc,
718 const bool increasing =
false)
721 double sminarr[BlackoilPhases::MaxNumPhases];
722 double smaxarr[BlackoilPhases::MaxNumPhases];
723 props.
satRange(1, &cell, sminarr, smaxarr);
724 const double s0 = increasing ? smaxarr[phase] : sminarr[phase];
725 const double s1 = increasing ? sminarr[phase] : smaxarr[phase];
728 const PcEq f(props, phase, cell, target_pc);
729 const double f0 = f(s0);
730 const double f1 = f(s1);
734 }
else if (f1 > 0.0) {
737 const int max_iter = 60;
738 const double tol = 1e-6;
741 const double sol = ScalarSolver::solve(f, std::min(s0, s1), std::max(s0, s1), max_iter, tol, iter_used);
756 const double target_pc)
761 target_pc_(target_pc)
763 std::fill(s_, s_ + BlackoilPhases::MaxNumPhases, 0.0);
764 std::fill(pc_, pc_ + BlackoilPhases::MaxNumPhases, 0.0);
766 double operator()(
double s)
const 769 s_[phase2_] = 1.0 - s;
770 props_.
capPress(1, s_, &cell_, pc_, 0);
771 return pc_[phase1_] + pc_[phase2_] - target_pc_;
778 const int target_pc_;
779 mutable double s_[BlackoilPhases::MaxNumPhases];
780 mutable double pc_[BlackoilPhases::MaxNumPhases];
793 const double target_pc)
796 double sminarr[BlackoilPhases::MaxNumPhases];
797 double smaxarr[BlackoilPhases::MaxNumPhases];
798 props.
satRange(1, &cell, sminarr, smaxarr);
799 const double smin = sminarr[phase1];
800 const double smax = smaxarr[phase1];
803 const PcEqSum f(props, phase1, phase2, cell, target_pc);
804 const double f0 = f(smin);
805 const double f1 = f(smax);
808 }
else if (f1 > 0.0) {
811 const int max_iter = 30;
812 const double tol = 1e-6;
815 const double sol = ScalarSolver::solve(f, smin, smax, max_iter, tol, iter_used);
822 const double cellDepth,
823 const double contactDepth,
826 const bool increasing =
false)
829 double sminarr[BlackoilPhases::MaxNumPhases];
830 double smaxarr[BlackoilPhases::MaxNumPhases];
831 props.
satRange(1, &cell, sminarr, smaxarr);
832 const double s0 = increasing ? smaxarr[phase] : sminarr[phase];
833 const double s1 = increasing ? sminarr[phase] : smaxarr[phase];
835 if (cellDepth < contactDepth){
849 double sminarr[BlackoilPhases::MaxNumPhases];
850 double smaxarr[BlackoilPhases::MaxNumPhases];
851 props.
satRange(1, &cell, sminarr, smaxarr);
854 const PcEq f(props, phase, cell, 0);
855 const double f0 = f(sminarr[phase]);
856 const double f1 = f(smaxarr[phase]);
857 return std::abs(f0 - f1) < std::numeric_limits<double>::epsilon();
864 #endif // OPM_EQUILIBRATIONHELPERS_HEADER_INCLUDED
std::vector< double > operator()(const double p, const double T, const std::vector< double > &z) const
Compute phase densities of all phases at phase point given by (pressure, surface volume) tuple...
Definition: EquilibrationHelpers.hpp:132
Definition: RootFinders.hpp:103
const CalcEvaporation & evaporationCalculator() const
Retrieve vapourised oil-gas ratio calculator of current region.
Definition: EquilibrationHelpers.hpp:660
virtual void capPress(const int n, const double *s, const int *cells, double *pc, double *dpcds) const =0
double satFromPc(const BlackoilPropertiesInterface &props, const int phase, const int cell, const double target_pc, const bool increasing=false)
Compute saturation of some phase corresponding to a given capillary pressure.
Definition: EquilibrationHelpers.hpp:714
Aggregate information base of an equilibration region.
Definition: EquilibrationHelpers.hpp:569
virtual void satRange(const int n, const int *cells, double *smin, double *smax) const =0
Obtain the range of allowable saturation values.
Definition: EquilibrationHelpers.hpp:91
const CalcDissolution & dissolutionCalculator() const
Retrieve dissolved gas-oil ratio calculator of current region.
Definition: EquilibrationHelpers.hpp:653
virtual double operator()(const double depth, const double press, const double temp, const double sat=0.0) const =0
Function call operator.
double operator()(const double depth, const double press, const double temp, const double sat_gas=0.0) const
Function call.
Definition: EquilibrationHelpers.hpp:268
Definition: AnisotropicEikonal.cpp:446
DensityCalculator(const BlackoilPropertiesInterface &props, const int c)
Constructor.
Definition: EquilibrationHelpers.hpp:112
double pressure() const
Pressure at datum depth in current region.
Definition: EquilibrationHelpers.hpp:616
double pcow_woc() const
water-oil capillary pressure at water-oil contact.
Definition: EquilibrationHelpers.hpp:628
double satFromDepth(const BlackoilPropertiesInterface &props, 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.hpp:821
Functor for inverting a sum of capillary pressure functions.
Definition: EquilibrationHelpers.hpp:750
virtual void matrix(const int n, const double *p, const double *T, const double *z, const int *cells, double *A, double *dAdp) const =0
const CalcDensity & densityCalculator() const
Retrieve phase density calculator of current region.
Definition: EquilibrationHelpers.hpp:646
const PhaseUsage & phaseUsage() const
Retrieve active fluid phase summary.
Definition: EquilibrationHelpers.hpp:666
Functor for inverting capillary pressure function.
Definition: EquilibrationHelpers.hpp:681
Miscibility::RsFunction CalcEvaporation
Type of vapourised oil-gas ratio calculator.
Definition: EquilibrationHelpers.hpp:606
Abstract base class for blackoil fluid and reservoir properties.
Definition: BlackoilPropertiesInterface.hpp:37
RsVD(const BlackoilPropertiesInterface &props, const int cell, const std::vector< double > &depth, const std::vector< double > &rs)
Constructor.
Definition: EquilibrationHelpers.hpp:237
Definition: BlackoilPhases.hpp:36
double datum() const
Datum depth in current region.
Definition: EquilibrationHelpers.hpp:611
double zwoc() const
Depth of water-oil contact.
Definition: EquilibrationHelpers.hpp:621
bool isConstPc(const BlackoilPropertiesInterface &props, const int phase, const int cell)
Return true if capillary pressure function is constant.
Definition: EquilibrationHelpers.hpp:844
double operator()(const double, const double, const double, const double=0.0) const
Function call.
Definition: EquilibrationHelpers.hpp:212
Miscibility::RsFunction CalcDissolution
Type of dissolved gas-oil ratio calculator.
Definition: EquilibrationHelpers.hpp:601
virtual int numPhases() const =0
double operator()(const double depth, const double press, const double temp, const double sat_oil=0.0) const
Function call.
Definition: EquilibrationHelpers.hpp:348
RvVD(const BlackoilPropertiesInterface &props, const int cell, const std::vector< double > &depth, const std::vector< double > &rv)
Constructor.
Definition: EquilibrationHelpers.hpp:317
Base class for phase mixing functions.
Definition: EquilibrationHelpers.hpp:165
double pcgo_goc() const
Gas-oil capillary pressure at gas-oil contact.
Definition: EquilibrationHelpers.hpp:640
Type that implements "vaporized oil-gas ratio" tabulated as a function of depth policy.
Definition: EquilibrationHelpers.hpp:307
Type that implements "dissolved gas-oil ratio" tabulated as a function of depth policy.
Definition: EquilibrationHelpers.hpp:227
virtual PhaseUsage phaseUsage() const =0
double zgoc() const
Depth of gas-oil contact.
Definition: EquilibrationHelpers.hpp:633
Type that implements "no phase mixing" policy.
Definition: EquilibrationHelpers.hpp:193
double satFromSumOfPcs(const BlackoilPropertiesInterface &props, const int phase1, const int phase2, const int cell, const double target_pc)
Compute saturation of some phase corresponding to a given capillary pressure, where the capillary pre...
Definition: EquilibrationHelpers.hpp:789
DensCalc CalcDensity
Type of density calculator.
Definition: EquilibrationHelpers.hpp:596
EquilReg(const EquilRecord &rec, const DensCalc &density, std::shared_ptr< Miscibility::RsFunction > rs, std::shared_ptr< Miscibility::RsFunction > rv, const PhaseUsage &pu)
Constructor.
Definition: EquilibrationHelpers.hpp:580