My Project
TwoPhaseLETCurvesParams.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_TWO_PHASE_LET_CURVES_PARAMS_HPP
28 #define OPM_TWO_PHASE_LET_CURVES_PARAMS_HPP
29 
32 
33 namespace Opm {
34 
43 template <class TraitsT>
45 {
46  using Scalar = typename TraitsT::Scalar;
47 
48 public:
49  using Traits = TraitsT;
50 
51  static constexpr int wIdx = 0; //wetting phase index for two phase let
52  static constexpr int nwIdx = 1; //non-wetting phase index for two phase let
53 
55  {
56  Valgrind::SetUndefined(*this);
57  }
58 
59  virtual ~TwoPhaseLETCurvesParams() {}
60 
65  void finalize()
66  {
68 
69  // printLETCoeffs();
70  }
71 
75  Scalar Smin(const unsigned phaseIdx) const
76  { EnsureFinalized::check(); if (phaseIdx<Traits::numPhases) return Smin_[phaseIdx]; return 0.0;}
77 
81  Scalar dS(const unsigned phaseIdx) const
82  { EnsureFinalized::check(); if (phaseIdx<Traits::numPhases) return dS_[phaseIdx]; return 0.0;}
83 
87  Scalar Sminpc() const
88  { EnsureFinalized::check(); return Sminpc_; }
89 
93  Scalar dSpc() const
94  { EnsureFinalized::check(); return dSpc_; }
95 
99  Scalar L(const unsigned phaseIdx) const
100  { EnsureFinalized::check(); if (phaseIdx<Traits::numPhases) return L_[phaseIdx]; return 0.0;}
101 
105  Scalar E(const unsigned phaseIdx) const
106  { EnsureFinalized::check(); if (phaseIdx<Traits::numPhases) return E_[phaseIdx]; return 0.0;}
107 
111  Scalar T(const unsigned phaseIdx) const
112  { EnsureFinalized::check(); if (phaseIdx<Traits::numPhases) return T_[phaseIdx]; return 0.0;}
113 
117  Scalar Krt(const unsigned phaseIdx) const
118  { EnsureFinalized::check(); if (phaseIdx<Traits::numPhases) return Krt_[phaseIdx]; return 0.0;}
119 
123  Scalar Lpc() const
124  { EnsureFinalized::check(); return Lpc_; }
125 
129  Scalar Epc() const
130  { EnsureFinalized::check(); return Epc_; }
131 
135  Scalar Tpc() const
136  { EnsureFinalized::check(); return Tpc_; }
137 
141  Scalar Pcir() const
142  { EnsureFinalized::check(); return Pcir_; }
143 
147  Scalar Pct() const
148  { EnsureFinalized::check(); return Pct_; }
149 
156  template <class Container>
157  void setKrwSamples(const Container& letProp, const Container& )
158  {
159  setLETCoeffs(wIdx, letProp[2], letProp[3], letProp[4], letProp[5]);
160  //Smin_[wIdx] = letProp[1];
161  //dS_[nwIdx] = 1.0 - letProp[0];
162  Smin_[wIdx] = letProp[0];
163  dS_[wIdx] = letProp[1] - letProp[0];
164  }
165 
172  template <class Container>
173  void setKrnSamples(const Container& letProp, const Container& )
174  {
175  setLETCoeffs(nwIdx, letProp[2], letProp[3], letProp[4], letProp[5]);
176  //Smin_[nwIdx] = letProp[1];
177  //dS_[wIdx] = 1.0 - letProp[0];
178  Smin_[nwIdx] = letProp[0];
179  dS_[nwIdx] = letProp[1] - letProp[0];
180  }
181 
182 
189  template <class Container>
190  void setPcnwSamples(const Container& letProp, const Container& )
191  {
192  setLETPcCoeffs(letProp[2], letProp[3], letProp[4], letProp[5], letProp[6]);
193  Sminpc_ = letProp[0];
194  dSpc_ = 1.0 - letProp[0] - letProp[1];
195  }
196 
197 private:
202  void setLETCoeffs(unsigned phaseIdx, Scalar L, Scalar E, Scalar T, Scalar Krt)
203  {
204  if (phaseIdx < Traits::numPhases) {
205  L_[phaseIdx]=L;
206  E_[phaseIdx]=E;
207  T_[phaseIdx]=T;
208  Krt_[phaseIdx]=Krt;
209  }
210  }
211 
216  void setLETPcCoeffs(Scalar L, Scalar E, Scalar T, Scalar Pcir, Scalar Pct)
217  {
218  Lpc_=L;
219  Epc_=E;
220  Tpc_=T;
221  Pcir_=Pcir;
222  Pct_=Pct;
223  }
224 
229  void printLETCoeffs()
230  {
231 /*
232  std::cout << "# LET parameters: "<< std::endl;
233  for (int i=0; i<Traits::numPhases; ++i) {
234  std::cout << "Kr[" << i;
235  std::cout << "]: Smin:" << Smin_[i];
236  std::cout << " dS:" << dS_[i];
237  std::cout << " L:" << L_[i];
238  std::cout << " E:" << E_[i];
239  std::cout << " T:" << T_[i];
240  std::cout << " Krt:" << Krt_[i];
241  std::cout << std::endl;
242  }
243 
244  std::cout << "Pc: Smin:" << Sminpc_;
245  std::cout << " dS:" << dSpc_;
246  std::cout << " L:" << Lpc_;
247  std::cout << " E:" << Epc_;
248  std::cout << " T:" << Tpc_;
249  std::cout << " Pcir:" << Pcir_;
250  std::cout << " Pct:" << Pct_;
251  std::cout << std::endl;
252 
253  std::cout << "================================="<< std::endl;
254 */
255  }
256 
257  Scalar Smin_[Traits::numPhases];
258  Scalar dS_[Traits::numPhases];
259 
260  Scalar L_[Traits::numPhases];
261  Scalar E_[Traits::numPhases];
262  Scalar T_[Traits::numPhases];
263  Scalar Krt_[Traits::numPhases];
264 
265  Scalar Sminpc_;
266  Scalar dSpc_;
267  Scalar Lpc_;
268  Scalar Epc_;
269  Scalar Tpc_;
270  Scalar Pcir_;
271  Scalar Pct_;
272 };
273 
274 } // namespace Opm
275 
276 #endif // OPM_TWO_PHASE_LET_CURVES_PARAMS_HPP
Default implementation for asserting finalization of parameter objects.
Some templates to wrap the valgrind client request macros.
Default implementation for asserting finalization of parameter objects.
Definition: EnsureFinalized.hpp:47
void finalize()
Mark the object as finalized.
Definition: EnsureFinalized.hpp:75
Specification of the material parameters for the LET constitutive relations.
Definition: TwoPhaseLETCurvesParams.hpp:45
Scalar Sminpc() const
Returns the Epc_ parameter.
Definition: TwoPhaseLETCurvesParams.hpp:87
void setPcnwSamples(const Container &letProp, const Container &)
Set the LET-related parameters for the capillary pressure curve of the non-wetting phase.
Definition: TwoPhaseLETCurvesParams.hpp:190
Scalar Smin(const unsigned phaseIdx) const
Returns the Smin_ parameter.
Definition: TwoPhaseLETCurvesParams.hpp:75
void finalize()
Calculate all dependent quantities once the independent quantities of the parameter object have been ...
Definition: TwoPhaseLETCurvesParams.hpp:65
Scalar T(const unsigned phaseIdx) const
Returns the T_ parameter.
Definition: TwoPhaseLETCurvesParams.hpp:111
Scalar dSpc() const
Returns the Epc_ parameter.
Definition: TwoPhaseLETCurvesParams.hpp:93
Scalar Epc() const
Returns the Epc_ parameter.
Definition: TwoPhaseLETCurvesParams.hpp:129
Scalar L(const unsigned phaseIdx) const
Returns the L_ parameter.
Definition: TwoPhaseLETCurvesParams.hpp:99
Scalar dS(const unsigned phaseIdx) const
Returns the dS_ parameter.
Definition: TwoPhaseLETCurvesParams.hpp:81
Scalar Lpc() const
Returns the Lpc_ parameter.
Definition: TwoPhaseLETCurvesParams.hpp:123
Scalar E(const unsigned phaseIdx) const
Returns the E_ parameter.
Definition: TwoPhaseLETCurvesParams.hpp:105
void setKrwSamples(const Container &letProp, const Container &)
Set the LET-related parameters for the relative permeability curve of the wetting phase.
Definition: TwoPhaseLETCurvesParams.hpp:157
Scalar Krt(const unsigned phaseIdx) const
Returns the Krt_ parameter.
Definition: TwoPhaseLETCurvesParams.hpp:117
void setKrnSamples(const Container &letProp, const Container &)
Set the LET-related parameters for the relative permeability curve of the non-wetting phase.
Definition: TwoPhaseLETCurvesParams.hpp:173
Scalar Pct() const
Returns the Pct_ parameter.
Definition: TwoPhaseLETCurvesParams.hpp:147
Scalar Pcir() const
Returns the Pcir_ parameter.
Definition: TwoPhaseLETCurvesParams.hpp:141
Scalar Tpc() const
Returns the Tpc_ parameter.
Definition: TwoPhaseLETCurvesParams.hpp:135