EquilibrationHelpers.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 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*/
29#ifndef OPM_EQUILIBRATION_HELPERS_HPP
30#define OPM_EQUILIBRATION_HELPERS_HPP
31
32#include <opm/material/common/Tabulated1DFunction.hpp>
33
34#include <opm/input/eclipse/EclipseState/InitConfig/Equil.hpp>
35
36#include <memory>
37#include <vector>
38
39/*
40 ---- synopsis of EquilibrationHelpers.hpp ----
41
42 namespace Opm
43 {
44 namespace EQUIL {
45
46 namespace Miscibility {
47 class RsFunction;
48 class NoMixing;
49 template <class FluidSystem>
50 class RsVD;
51 template <class FluidSystem>
52 class RsSatAtContact;
53 }
54
55 class EquilReg;
56
57
58 template <class FluidSystem, class MaterialLawManager>
59 struct PcEq;
60
61 template <class FluidSystem, class MaterialLawManager>
62 double satFromPc(const MaterialLawManager& materialLawManager,
63 const int phase,
64 const int cell,
65 const double targetPc,
66 const bool increasing = false)
67 template <class FluidSystem, class MaterialLawManager>
68 double satFromSumOfPcs(const MaterialLawManager& materialLawManager,
69 const int phase1,
70 const int phase2,
71 const int cell,
72 const double targetPc)
73 } // namespace Equil
74 } // namespace Opm
75
76 ---- end of synopsis of EquilibrationHelpers.hpp ----
77*/
78
79namespace Opm {
80
88namespace EQUIL {
89
94namespace Miscibility {
95
99template<class Scalar>
101{
102public:
103 virtual ~RsFunction() = default;
104
120 virtual Scalar operator()(const Scalar depth,
121 const Scalar press,
122 const Scalar temp,
123 const Scalar sat = 0.0) const = 0;
124};
125
126
130template<class Scalar>
131class NoMixing : public RsFunction<Scalar>
132{
133public:
150 Scalar
151 operator()(const Scalar /* depth */,
152 const Scalar /* press */,
153 const Scalar /* temp */,
154 const Scalar /* sat */ = 0.0) const override
155 {
156 return 0.0;
157 }
158};
159
160
166template <class FluidSystem>
167class RsVD : public RsFunction<typename FluidSystem::Scalar>
168{
169public:
170 using Scalar = typename FluidSystem::Scalar;
178 RsVD(const int pvtRegionIdx,
179 const std::vector<Scalar>& depth,
180 const std::vector<Scalar>& rs);
181
197 Scalar operator()(const Scalar depth,
198 const Scalar press,
199 const Scalar temp,
200 const Scalar satGas = 0.0) const override;
201
202private:
203 using RsVsDepthFunc = Tabulated1DFunction<Scalar>;
204
205 const int pvtRegionIdx_;
206 RsVsDepthFunc rsVsDepth_;
207
208 Scalar satRs(const Scalar press, const Scalar temp) const;
209};
210
211
217template <class FluidSystem>
218class PBVD : public RsFunction<typename FluidSystem::Scalar>
219{
220public:
221 using Scalar = typename FluidSystem::Scalar;
229 PBVD(const int pvtRegionIdx,
230 const std::vector<Scalar>& depth,
231 const std::vector<Scalar>& pbub);
232
247 Scalar operator()(const Scalar depth,
248 const Scalar cellPress,
249 const Scalar temp,
250 const Scalar satGas = 0.0) const override;
251
252private:
253 using PbubVsDepthFunc = Tabulated1DFunction<Scalar>;
254
255 const int pvtRegionIdx_;
256 PbubVsDepthFunc pbubVsDepth_;
257
258 Scalar satRs(const Scalar press, const Scalar temp) const;
259};
260
261
267template <class FluidSystem>
268class PDVD : public RsFunction<typename FluidSystem::Scalar>
269{
270 using Scalar = typename FluidSystem::Scalar;
271public:
279 PDVD(const int pvtRegionIdx,
280 const std::vector<Scalar>& depth,
281 const std::vector<Scalar>& pdew);
282
297 Scalar operator()(const Scalar depth,
298 const Scalar cellPress,
299 const Scalar temp,
300 const Scalar satOil = 0.0) const override;
301
302private:
303 using PdewVsDepthFunc = Tabulated1DFunction<Scalar>;
304
305 const int pvtRegionIdx_;
306 PdewVsDepthFunc pdewVsDepth_;
307
308 Scalar satRv(const Scalar press, const Scalar temp) const;
309};
310
311
317template <class FluidSystem>
318class RvVD : public RsFunction<typename FluidSystem::Scalar>
319{
320 using Scalar = typename FluidSystem::Scalar;
321public:
329 RvVD(const int pvtRegionIdx,
330 const std::vector<Scalar>& depth,
331 const std::vector<Scalar>& rv);
332
348 Scalar operator()(const Scalar depth,
349 const Scalar press,
350 const Scalar temp,
351 const Scalar satOil = 0.0) const override;
352
353private:
354 using RvVsDepthFunc = Tabulated1DFunction<Scalar>;
355
356 const int pvtRegionIdx_;
357 RvVsDepthFunc rvVsDepth_;
358
359 Scalar satRv(const Scalar press, const Scalar temp) const;
360};
361
362
368template <class FluidSystem>
369class RvwVD : public RsFunction<typename FluidSystem::Scalar>
370{
371 using Scalar = typename FluidSystem::Scalar;
372public:
380 RvwVD(const int pvtRegionIdx,
381 const std::vector<Scalar>& depth,
382 const std::vector<Scalar>& rvw);
383
399 Scalar operator()(const Scalar depth,
400 const Scalar press,
401 const Scalar temp,
402 const Scalar satWat = 0.0) const override;
403
404private:
405 using RvwVsDepthFunc = Tabulated1DFunction<Scalar>;
406
407 const int pvtRegionIdx_;
408 RvwVsDepthFunc rvwVsDepth_;
409
410 Scalar satRvw(const Scalar press, const Scalar temp) const;
411};
412
413
428template <class FluidSystem>
429class RsSatAtContact : public RsFunction<typename FluidSystem::Scalar>
430{
431 using Scalar = typename FluidSystem::Scalar;
432public:
440 RsSatAtContact(const int pvtRegionIdx,
441 const Scalar pContact,
442 const Scalar T_contact);
443
459 Scalar operator()(const Scalar /* depth */,
460 const Scalar press,
461 const Scalar temp,
462 const Scalar satGas = 0.0) const override;
463
464private:
465 const int pvtRegionIdx_;
466 Scalar rsSatContact_;
467
468 Scalar satRs(const Scalar press, const Scalar temp) const;
469};
470
471
486template <class FluidSystem>
487class RvSatAtContact : public RsFunction<typename FluidSystem::Scalar>
488{
489 using Scalar = typename FluidSystem::Scalar;
490public:
498 RvSatAtContact(const int pvtRegionIdx,
499 const Scalar pContact,
500 const Scalar T_contact);
501
517 Scalar operator()(const Scalar /*depth*/,
518 const Scalar press,
519 const Scalar temp,
520 const Scalar satOil = 0.0) const override;
521
522private:
523 const int pvtRegionIdx_;
524 Scalar rvSatContact_;
525
526 Scalar satRv(const Scalar press, const Scalar temp) const;
527};
528
543template <class FluidSystem>
544class RvwSatAtContact : public RsFunction<typename FluidSystem::Scalar>
545{
546 using Scalar = typename FluidSystem::Scalar;
547public:
555 RvwSatAtContact(const int pvtRegionIdx,
556 const Scalar pContact,
557 const Scalar T_contact);
558
574 Scalar operator()(const Scalar /*depth*/,
575 const Scalar press,
576 const Scalar temp,
577 const Scalar satWat = 0.0) const override;
578
579private:
580 const int pvtRegionIdx_;
581 Scalar rvwSatContact_;
582
583 Scalar satRvw(const Scalar press, const Scalar temp) const;
584};
585
586} // namespace Miscibility
587
607template<class Scalar>
609{
610 using TabulatedFunction = Tabulated1DFunction<Scalar>;
611
612public:
622 EquilReg(const EquilRecord& rec,
623 std::shared_ptr<Miscibility::RsFunction<Scalar>> rs,
624 std::shared_ptr<Miscibility::RsFunction<Scalar>> rv,
625 std::shared_ptr<Miscibility::RsFunction<Scalar>> rvw,
626 const TabulatedFunction& tempVdTable,
627 const TabulatedFunction& saltVdTable,
628 const int pvtIdx);
629
634
639
644
645
649 Scalar datum() const;
650
654 Scalar pressure() const;
655
659 Scalar zwoc() const;
660
666 Scalar pcowWoc() const;
667
671 Scalar zgoc() const;
672
678 Scalar pcgoGoc() const;
679
687 int equilibrationAccuracy() const;
688
694
700
706
707 const TabulatedFunction& saltVdTable() const;
708 const TabulatedFunction& tempVdTable() const;
709
713 int pvtIdx() const;
714
715private:
716 EquilRecord rec_;
717 std::shared_ptr<Miscibility::RsFunction<Scalar>> rs_;
718 std::shared_ptr<Miscibility::RsFunction<Scalar>> rv_;
719 std::shared_ptr<Miscibility::RsFunction<Scalar>> rvw_;
720 const TabulatedFunction& tempVdTable_;
721 const TabulatedFunction& saltVdTable_;
722 const int pvtIdx_;
723};
724
725
726
730template <class FluidSystem, class MaterialLawManager>
731struct PcEq
732{
733 using Scalar = typename FluidSystem::Scalar;
734 PcEq(const MaterialLawManager& materialLawManager,
735 const int phase,
736 const int cell,
737 const Scalar targetPc);
738
739 Scalar operator()(Scalar s) const;
740
741private:
742 const MaterialLawManager& materialLawManager_;
743 const int phase_;
744 const int cell_;
745 const Scalar targetPc_;
746};
747
748template <class FluidSystem, class MaterialLawManager>
749typename FluidSystem::Scalar
750minSaturations(const MaterialLawManager& materialLawManager,
751 const int phase, const int cell);
752
753template <class FluidSystem, class MaterialLawManager>
754typename FluidSystem::Scalar
755maxSaturations(const MaterialLawManager& materialLawManager,
756 const int phase, const int cell);
757
760template <class FluidSystem, class MaterialLawManager>
761typename FluidSystem::Scalar
762satFromPc(const MaterialLawManager& materialLawManager,
763 const int phase,
764 const int cell,
765 const typename FluidSystem::Scalar targetPc,
766 const bool increasing = false);
767
771template <class FluidSystem, class MaterialLawManager>
773{
774 using Scalar = typename FluidSystem::Scalar;
775 PcEqSum(const MaterialLawManager& materialLawManager,
776 const int phase1,
777 const int phase2,
778 const int cell,
779 const Scalar targetPc);
780
781 Scalar operator()(Scalar s) const;
782
783private:
784 const MaterialLawManager& materialLawManager_;
785 const int phase1_;
786 const int phase2_;
787 const int cell_;
788 const Scalar targetPc_;
789};
790
794template <class FluidSystem, class MaterialLawManager>
795typename FluidSystem::Scalar
796satFromSumOfPcs(const MaterialLawManager& materialLawManager,
797 const int phase1,
798 const int phase2,
799 const int cell,
800 const typename FluidSystem::Scalar targetPc);
801
803template <class FluidSystem, class MaterialLawManager>
804typename FluidSystem::Scalar
805satFromDepth(const MaterialLawManager& materialLawManager,
806 const typename FluidSystem::Scalar cellDepth,
807 const typename FluidSystem::Scalar contactDepth,
808 const int phase,
809 const int cell,
810 const bool increasing = false);
811
813template <class FluidSystem, class MaterialLawManager>
814bool isConstPc(const MaterialLawManager& materialLawManager,
815 const int phase,
816 const int cell);
817
818} // namespace Equil
819} // namespace Opm
820
821#endif // OPM_EQUILIBRATION_HELPERS_HPP
Definition: EquilibrationHelpers.hpp:609
Scalar datum() const
Definition: EquilibrationHelpers_impl.hpp:373
int pvtIdx() const
Definition: EquilibrationHelpers_impl.hpp:450
EquilReg(const EquilRecord &rec, std::shared_ptr< Miscibility::RsFunction< Scalar > > rs, std::shared_ptr< Miscibility::RsFunction< Scalar > > rv, std::shared_ptr< Miscibility::RsFunction< Scalar > > rvw, const TabulatedFunction &tempVdTable, const TabulatedFunction &saltVdTable, const int pvtIdx)
Definition: EquilibrationHelpers_impl.hpp:355
Scalar pcgoGoc() const
Definition: EquilibrationHelpers_impl.hpp:403
int equilibrationAccuracy() const
Definition: EquilibrationHelpers_impl.hpp:409
const CalcEvaporation & evaporationCalculator() const
Definition: EquilibrationHelpers_impl.hpp:423
Scalar zgoc() const
Definition: EquilibrationHelpers_impl.hpp:397
Scalar pressure() const
Definition: EquilibrationHelpers_impl.hpp:379
const TabulatedFunction & saltVdTable() const
Definition: EquilibrationHelpers_impl.hpp:437
Scalar pcowWoc() const
Definition: EquilibrationHelpers_impl.hpp:391
const TabulatedFunction & tempVdTable() const
Definition: EquilibrationHelpers_impl.hpp:444
Scalar zwoc() const
Definition: EquilibrationHelpers_impl.hpp:385
const CalcDissolution & dissolutionCalculator() const
Definition: EquilibrationHelpers_impl.hpp:416
const CalcWaterEvaporation & waterEvaporationCalculator() const
Definition: EquilibrationHelpers_impl.hpp:430
Definition: EquilibrationHelpers.hpp:132
Scalar operator()(const Scalar, const Scalar, const Scalar, const Scalar=0.0) const override
Definition: EquilibrationHelpers.hpp:151
Definition: EquilibrationHelpers.hpp:219
Scalar operator()(const Scalar depth, const Scalar cellPress, const Scalar temp, const Scalar satGas=0.0) const override
Definition: EquilibrationHelpers_impl.hpp:107
PBVD(const int pvtRegionIdx, const std::vector< Scalar > &depth, const std::vector< Scalar > &pbub)
Definition: EquilibrationHelpers_impl.hpp:96
typename FluidSystem::Scalar Scalar
Definition: EquilibrationHelpers.hpp:221
Definition: EquilibrationHelpers.hpp:269
Scalar operator()(const Scalar depth, const Scalar cellPress, const Scalar temp, const Scalar satOil=0.0) const override
Definition: EquilibrationHelpers_impl.hpp:144
PDVD(const int pvtRegionIdx, const std::vector< Scalar > &depth, const std::vector< Scalar > &pdew)
Definition: EquilibrationHelpers_impl.hpp:133
Definition: EquilibrationHelpers.hpp:101
virtual Scalar operator()(const Scalar depth, const Scalar press, const Scalar temp, const Scalar sat=0.0) const =0
Definition: EquilibrationHelpers.hpp:430
RsSatAtContact(const int pvtRegionIdx, const Scalar pContact, const Scalar T_contact)
Definition: EquilibrationHelpers_impl.hpp:258
Scalar operator()(const Scalar, const Scalar press, const Scalar temp, const Scalar satGas=0.0) const override
Definition: EquilibrationHelpers_impl.hpp:267
Definition: EquilibrationHelpers.hpp:168
Scalar operator()(const Scalar depth, const Scalar press, const Scalar temp, const Scalar satGas=0.0) const override
Definition: EquilibrationHelpers_impl.hpp:69
typename FluidSystem::Scalar Scalar
Definition: EquilibrationHelpers.hpp:170
RsVD(const int pvtRegionIdx, const std::vector< Scalar > &depth, const std::vector< Scalar > &rs)
Definition: EquilibrationHelpers_impl.hpp:58
Definition: EquilibrationHelpers.hpp:488
Scalar operator()(const Scalar, const Scalar press, const Scalar temp, const Scalar satOil=0.0) const override
Definition: EquilibrationHelpers_impl.hpp:299
RvSatAtContact(const int pvtRegionIdx, const Scalar pContact, const Scalar T_contact)
Definition: EquilibrationHelpers_impl.hpp:290
Definition: EquilibrationHelpers.hpp:319
Scalar operator()(const Scalar depth, const Scalar press, const Scalar temp, const Scalar satOil=0.0) const override
Definition: EquilibrationHelpers_impl.hpp:181
RvVD(const int pvtRegionIdx, const std::vector< Scalar > &depth, const std::vector< Scalar > &rv)
Definition: EquilibrationHelpers_impl.hpp:170
Definition: EquilibrationHelpers.hpp:545
Scalar operator()(const Scalar, const Scalar press, const Scalar temp, const Scalar satWat=0.0) const override
Definition: EquilibrationHelpers_impl.hpp:331
RvwSatAtContact(const int pvtRegionIdx, const Scalar pContact, const Scalar T_contact)
Definition: EquilibrationHelpers_impl.hpp:322
Definition: EquilibrationHelpers.hpp:370
Scalar operator()(const Scalar depth, const Scalar press, const Scalar temp, const Scalar satWat=0.0) const override
Definition: EquilibrationHelpers_impl.hpp:224
RvwVD(const int pvtRegionIdx, const std::vector< Scalar > &depth, const std::vector< Scalar > &rvw)
Definition: EquilibrationHelpers_impl.hpp:213
FluidSystem::Scalar satFromDepth(const MaterialLawManager &materialLawManager, const typename FluidSystem::Scalar cellDepth, const typename FluidSystem::Scalar 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:650
FluidSystem::Scalar satFromPc(const MaterialLawManager &materialLawManager, const int phase, const int cell, const typename FluidSystem::Scalar targetPc, const bool increasing=false)
Definition: EquilibrationHelpers_impl.hpp:586
FluidSystem::Scalar minSaturations(const MaterialLawManager &materialLawManager, const int phase, const int cell)
Definition: EquilibrationHelpers_impl.hpp:536
FluidSystem::Scalar maxSaturations(const MaterialLawManager &materialLawManager, const int phase, const int cell)
Definition: EquilibrationHelpers_impl.hpp:561
bool isConstPc(const MaterialLawManager &materialLawManager, const int phase, const int cell)
Return true if capillary pressure function is constant.
Definition: EquilibrationHelpers_impl.hpp:670
FluidSystem::Scalar satFromSumOfPcs(const MaterialLawManager &materialLawManager, const int phase1, const int phase2, const int cell, const typename FluidSystem::Scalar targetPc)
Definition: EquilibrationHelpers_impl.hpp:619
Definition: blackoilboundaryratevector.hh:39
Definition: EquilibrationHelpers.hpp:732
Scalar operator()(Scalar s) const
Definition: EquilibrationHelpers_impl.hpp:471
typename FluidSystem::Scalar Scalar
Definition: EquilibrationHelpers.hpp:733
PcEq(const MaterialLawManager &materialLawManager, const int phase, const int cell, const Scalar targetPc)
Definition: EquilibrationHelpers_impl.hpp:457
Definition: EquilibrationHelpers.hpp:773
PcEqSum(const MaterialLawManager &materialLawManager, const int phase1, const int phase2, const int cell, const Scalar targetPc)
Definition: EquilibrationHelpers_impl.hpp:497
Scalar operator()(Scalar s) const
Definition: EquilibrationHelpers_impl.hpp:513
typename FluidSystem::Scalar Scalar
Definition: EquilibrationHelpers.hpp:774