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
123 virtual Scalar operator()(const Scalar depth,
124 const Scalar press,
125 const Scalar temp,
126 const Scalar sat = 0.0) const = 0;
127};
128
129
133template<class Scalar>
134class NoMixing : public RsFunction<Scalar>
135{
136public:
144 Scalar
145 operator()(const Scalar /* depth */,
146 const Scalar /* press */,
147 const Scalar /* temp */,
148 const Scalar /* sat */ = 0.0) const override
149 {
150 return 0.0;
151 }
152};
153
154
160template <class FluidSystem>
161class RsVD : public RsFunction<typename FluidSystem::Scalar>
162{
163public:
164 using Scalar = typename FluidSystem::Scalar;
172 RsVD(const int pvtRegionIdx,
173 const std::vector<Scalar>& depth,
174 const std::vector<Scalar>& rs);
175
194 Scalar operator()(const Scalar depth,
195 const Scalar press,
196 const Scalar temp,
197 const Scalar satGas = 0.0) const override;
198
199private:
200 using RsVsDepthFunc = Tabulated1DFunction<Scalar>;
201
202 const int pvtRegionIdx_;
203 RsVsDepthFunc rsVsDepth_;
204
205 Scalar satRs(const Scalar press, const Scalar temp) const;
206};
207
208
214template <class FluidSystem>
215class PBVD : public RsFunction<typename FluidSystem::Scalar>
216{
217public:
218 using Scalar = typename FluidSystem::Scalar;
226 PBVD(const int pvtRegionIdx,
227 const std::vector<Scalar>& depth,
228 const std::vector<Scalar>& pbub);
229
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
300 Scalar operator()(const Scalar depth,
301 const Scalar cellPress,
302 const Scalar temp,
303 const Scalar satOil = 0.0) const override;
304
305private:
306 using PdewVsDepthFunc = Tabulated1DFunction<Scalar>;
307
308 const int pvtRegionIdx_;
309 PdewVsDepthFunc pdewVsDepth_;
310
311 Scalar satRv(const Scalar press, const Scalar temp) const;
312};
313
314
320template <class FluidSystem>
321class RvVD : public RsFunction<typename FluidSystem::Scalar>
322{
323 using Scalar = typename FluidSystem::Scalar;
324public:
332 RvVD(const int pvtRegionIdx,
333 const std::vector<Scalar>& depth,
334 const std::vector<Scalar>& rv);
335
354 Scalar operator()(const Scalar depth,
355 const Scalar press,
356 const Scalar temp,
357 const Scalar satOil = 0.0) const override;
358
359private:
360 using RvVsDepthFunc = Tabulated1DFunction<Scalar>;
361
362 const int pvtRegionIdx_;
363 RvVsDepthFunc rvVsDepth_;
364
365 Scalar satRv(const Scalar press, const Scalar temp) const;
366};
367
368
374template <class FluidSystem>
375class RvwVD : public RsFunction<typename FluidSystem::Scalar>
376{
377 using Scalar = typename FluidSystem::Scalar;
378public:
386 RvwVD(const int pvtRegionIdx,
387 const std::vector<Scalar>& depth,
388 const std::vector<Scalar>& rvw);
389
408 Scalar operator()(const Scalar depth,
409 const Scalar press,
410 const Scalar temp,
411 const Scalar satWat = 0.0) const override;
412
413private:
414 using RvwVsDepthFunc = Tabulated1DFunction<Scalar>;
415
416 const int pvtRegionIdx_;
417 RvwVsDepthFunc rvwVsDepth_;
418
419 Scalar satRvw(const Scalar press, const Scalar temp) const;
420};
421
422
437template <class FluidSystem>
438class RsSatAtContact : public RsFunction<typename FluidSystem::Scalar>
439{
440 using Scalar = typename FluidSystem::Scalar;
441public:
449 RsSatAtContact(const int pvtRegionIdx,
450 const Scalar pContact,
451 const Scalar T_contact);
452
471 Scalar operator()(const Scalar /* depth */,
472 const Scalar press,
473 const Scalar temp,
474 const Scalar satGas = 0.0) const override;
475
476private:
477 const int pvtRegionIdx_;
478 Scalar rsSatContact_;
479
480 Scalar satRs(const Scalar press, const Scalar temp) const;
481};
482
483
498template <class FluidSystem>
499class RvSatAtContact : public RsFunction<typename FluidSystem::Scalar>
500{
501 using Scalar = typename FluidSystem::Scalar;
502public:
510 RvSatAtContact(const int pvtRegionIdx,
511 const Scalar pContact,
512 const Scalar T_contact);
513
532 Scalar operator()(const Scalar /*depth*/,
533 const Scalar press,
534 const Scalar temp,
535 const Scalar satOil = 0.0) const override;
536
537private:
538 const int pvtRegionIdx_;
539 Scalar rvSatContact_;
540
541 Scalar satRv(const Scalar press, const Scalar temp) const;
542};
543
558template <class FluidSystem>
559class RvwSatAtContact : public RsFunction<typename FluidSystem::Scalar>
560{
561 using Scalar = typename FluidSystem::Scalar;
562public:
570 RvwSatAtContact(const int pvtRegionIdx,
571 const Scalar pContact,
572 const Scalar T_contact);
573
592 Scalar operator()(const Scalar /*depth*/,
593 const Scalar press,
594 const Scalar temp,
595 const Scalar satWat = 0.0) const override;
596
597private:
598 const int pvtRegionIdx_;
599 Scalar rvwSatContact_;
600
601 Scalar satRvw(const Scalar press, const Scalar temp) const;
602};
603
604} // namespace Miscibility
605
625template<class Scalar>
627{
628 using TabulatedFunction = Tabulated1DFunction<Scalar>;
629
630public:
642 EquilReg(const EquilRecord& rec,
643 std::shared_ptr<Miscibility::RsFunction<Scalar>> rs,
644 std::shared_ptr<Miscibility::RsFunction<Scalar>> rv,
645 std::shared_ptr<Miscibility::RsFunction<Scalar>> rvw,
646 const TabulatedFunction& tempVdTable,
647 const TabulatedFunction& saltVdTable,
648 const int pvtIdx);
649
654
659
664
665
669 Scalar datum() const;
670
674 Scalar pressure() const;
675
679 Scalar zwoc() const;
680
686 Scalar pcowWoc() const;
687
691 Scalar zgoc() const;
692
698 Scalar pcgoGoc() const;
699
707 int equilibrationAccuracy() const;
708
714
720
726
727 const TabulatedFunction& saltVdTable() const;
728 const TabulatedFunction& tempVdTable() const;
729
733 int pvtIdx() const;
734
735private:
736 EquilRecord rec_;
737 std::shared_ptr<Miscibility::RsFunction<Scalar>> rs_;
738 std::shared_ptr<Miscibility::RsFunction<Scalar>> rv_;
739 std::shared_ptr<Miscibility::RsFunction<Scalar>> rvw_;
740 const TabulatedFunction& tempVdTable_;
741 const TabulatedFunction& saltVdTable_;
742 const int pvtIdx_;
743};
744
745
746
750template <class FluidSystem, class MaterialLawManager>
751struct PcEq
752{
753 using Scalar = typename FluidSystem::Scalar;
754 PcEq(const MaterialLawManager& materialLawManager,
755 const int phase,
756 const int cell,
757 const Scalar targetPc);
758
759 Scalar operator()(Scalar s) const;
760
761private:
762 const MaterialLawManager& materialLawManager_;
763 const int phase_;
764 const int cell_;
765 const Scalar targetPc_;
766};
767
768template <class FluidSystem, class MaterialLawManager>
769typename FluidSystem::Scalar
770minSaturations(const MaterialLawManager& materialLawManager,
771 const int phase, const int cell);
772
773template <class FluidSystem, class MaterialLawManager>
774typename FluidSystem::Scalar
775maxSaturations(const MaterialLawManager& materialLawManager,
776 const int phase, const int cell);
777
780template <class FluidSystem, class MaterialLawManager>
781typename FluidSystem::Scalar
782satFromPc(const MaterialLawManager& materialLawManager,
783 const int phase,
784 const int cell,
785 const typename FluidSystem::Scalar targetPc,
786 const bool increasing = false);
787
791template <class FluidSystem, class MaterialLawManager>
793{
794 using Scalar = typename FluidSystem::Scalar;
795 PcEqSum(const MaterialLawManager& materialLawManager,
796 const int phase1,
797 const int phase2,
798 const int cell,
799 const Scalar targetPc);
800
801 Scalar operator()(Scalar s) const;
802
803private:
804 const MaterialLawManager& materialLawManager_;
805 const int phase1_;
806 const int phase2_;
807 const int cell_;
808 const Scalar targetPc_;
809};
810
814template <class FluidSystem, class MaterialLawManager>
815typename FluidSystem::Scalar
816satFromSumOfPcs(const MaterialLawManager& materialLawManager,
817 const int phase1,
818 const int phase2,
819 const int cell,
820 const typename FluidSystem::Scalar targetPc);
821
823template <class FluidSystem, class MaterialLawManager>
824typename FluidSystem::Scalar
825satFromDepth(const MaterialLawManager& materialLawManager,
826 const typename FluidSystem::Scalar cellDepth,
827 const typename FluidSystem::Scalar contactDepth,
828 const int phase,
829 const int cell,
830 const bool increasing = false);
831
833template <class FluidSystem, class MaterialLawManager>
834bool isConstPc(const MaterialLawManager& materialLawManager,
835 const int phase,
836 const int cell);
837
838} // namespace Equil
839} // namespace Opm
840
841#endif // OPM_EQUILIBRATION_HELPERS_HPP
Definition: EquilibrationHelpers.hpp:627
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:135
Scalar operator()(const Scalar, const Scalar, const Scalar, const Scalar=0.0) const override
Definition: EquilibrationHelpers.hpp:145
Definition: EquilibrationHelpers.hpp:216
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:218
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:439
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:162
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:164
RsVD(const int pvtRegionIdx, const std::vector< Scalar > &depth, const std::vector< Scalar > &rs)
Definition: EquilibrationHelpers_impl.hpp:58
Definition: EquilibrationHelpers.hpp:500
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:322
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:560
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:376
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:752
Scalar operator()(Scalar s) const
Definition: EquilibrationHelpers_impl.hpp:471
typename FluidSystem::Scalar Scalar
Definition: EquilibrationHelpers.hpp:753
PcEq(const MaterialLawManager &materialLawManager, const int phase, const int cell, const Scalar targetPc)
Definition: EquilibrationHelpers_impl.hpp:457
Definition: EquilibrationHelpers.hpp:793
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:794