My Project
Loading...
Searching...
No Matches
SatCurveMultiplexerParams.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*/
27#ifndef OPM_SAT_CURVE_MULTIPLEXER_PARAMS_HPP
28#define OPM_SAT_CURVE_MULTIPLEXER_PARAMS_HPP
29
30
31#include "TwoPhaseLETCurves.hpp"
35
37
38#include <type_traits>
39#include <cassert>
40#include <memory>
41
42namespace Opm {
43
44enum class SatCurveMultiplexerApproach {
45 PiecewiseLinear,
46 LET
47};
48
57template <class TraitsT>
59{
60public:
61 using Traits = TraitsT;
62 using Scalar = typename TraitsT::Scalar;
63 enum { numPhases = 2 };
64
65private:
68
69 using LETParams = typename LETTwoPhaseLaw::Params;
70 using PLParams = typename PLTwoPhaseLaw::Params;
71
72 template <class ParamT>
73 struct Deleter
74 {
75 inline void operator () ( void* ptr )
76 {
77 delete static_cast< ParamT* > (ptr);
78 }
79 };
80
81 using ParamPointerType = std::shared_ptr<void>;
82
83public:
84
88 SatCurveMultiplexerParams() : realParams_()
89 {
90 }
91
93 : realParams_()
94 {
95 setApproach( other.approach() );
96 }
97
99 {
100 realParams_.reset();
101 setApproach( other.approach() );
102 return *this;
103 }
104
105 void setApproach(SatCurveMultiplexerApproach newApproach)
106 {
107 assert(realParams_ == 0);
108 approach_ = newApproach;
109
110 switch (approach()) {
111 case SatCurveMultiplexerApproach::LET:
112 realParams_ = ParamPointerType(new LETParams, Deleter< LETParams > () );
113 break;
114
115 case SatCurveMultiplexerApproach::PiecewiseLinear:
116 realParams_ = ParamPointerType(new PLParams, Deleter< PLParams > () );
117 break;
118 }
119 }
120
121 SatCurveMultiplexerApproach approach() const
122 { return approach_; }
123
124 // get the parameter object for the LET curve
125 template <SatCurveMultiplexerApproach approachV>
126 typename std::enable_if<approachV == SatCurveMultiplexerApproach::LET, LETParams>::type&
127 getRealParams()
128 {
129 assert(approach() == approachV);
130 return this->template castTo<LETParams>();
131 }
132
133 template <SatCurveMultiplexerApproach approachV>
134 typename std::enable_if<approachV == SatCurveMultiplexerApproach::LET, const LETParams>::type&
135 getRealParams() const
136 {
137 assert(approach() == approachV);
138 return this->template castTo<LETParams>();
139 }
140
141 // get the parameter object for the PL curve
142 template <SatCurveMultiplexerApproach approachV>
143 typename std::enable_if<approachV == SatCurveMultiplexerApproach::PiecewiseLinear, PLParams>::type&
144 getRealParams()
145 {
146 assert(approach() == approachV);
147 return this->template castTo<PLParams>();
148 }
149
150 template <SatCurveMultiplexerApproach approachV>
151 typename std::enable_if<approachV == SatCurveMultiplexerApproach::PiecewiseLinear, const PLParams>::type&
152 getRealParams() const
153 {
154 assert(approach() == approachV);
155 return this->template castTo<PLParams>();
156 }
157
158 template<class Serializer>
159 void serializeOp(Serializer& serializer)
160 {
161 switch (approach()) {
162 case SatCurveMultiplexerApproach::LET:
163 serializer(castTo<LETParams>());
164 break;
165
166 case SatCurveMultiplexerApproach::PiecewiseLinear:
167 serializer(castTo<PLParams>());
168 break;
169 }
170 }
171
172private:
173 template <class ParamT>
174 ParamT& castTo()
175 {
176 return *(static_cast<ParamT *> (realParams_.operator->()));
177 }
178
179 template <class ParamT>
180 const ParamT& castTo() const
181 {
182 return *(static_cast<const ParamT *> (realParams_.operator->()));
183 }
184
185 SatCurveMultiplexerApproach approach_;
186 ParamPointerType realParams_;
187};
188
189} // namespace Opm
190
191#endif // OPM_SAT_CURVE_MULTIPLEXER_PARAMS_HPP
Default implementation for asserting finalization of parameter objects.
Specification of the material parameters for a two-phase material law which uses a table and piecewis...
Implementation of a tabulated, piecewise linear capillary pressure law.
Default implementation for asserting finalization of parameter objects.
Definition EnsureFinalized.hpp:47
Implementation of a tabulated, piecewise linear capillary pressure law.
Definition PiecewiseLinearTwoPhaseMaterial.hpp:50
ParamsT Params
The type of the parameter objects for this law.
Definition PiecewiseLinearTwoPhaseMaterial.hpp:58
Specification of the material parameters for the saturation function multiplexer.
Definition SatCurveMultiplexerParams.hpp:59
SatCurveMultiplexerParams()
The multiplexer constructor.
Definition SatCurveMultiplexerParams.hpp:88
Implementation of the LET curve saturation functions.
Definition TwoPhaseLETCurves.hpp:50
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30