All Classes Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
EclHysteresisTwoPhaseLawParams.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 2 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_ECL_HYSTERESIS_TWO_PHASE_LAW_PARAMS_HPP
28 #define OPM_ECL_HYSTERESIS_TWO_PHASE_LAW_PARAMS_HPP
29 
30 #include "EclHysteresisConfig.hpp"
31 #include "EclEpsScalingPoints.hpp"
32 
33 #if HAVE_OPM_PARSER
34 #include <opm/parser/eclipse/Deck/Deck.hpp>
35 #include <opm/parser/eclipse/EclipseState/EclipseState.hpp>
36 #endif
37 
38 #include <string>
39 #include <memory>
40 #include <cassert>
41 #include <algorithm>
42 
44 
45 namespace Opm {
52 template <class EffLawT>
54 {
55  typedef typename EffLawT::Params EffLawParams;
56  typedef typename EffLawParams::Traits::Scalar Scalar;
57 
58 public:
59  typedef typename EffLawParams::Traits Traits;
60 
62  {
63  // These are initialized to two (even though they represent saturations)
64  // to signify that the values are larger than physically possible, and force
65  // using the drainage curve before the first saturation update.
66  pcSwMdc_ = 2.0;
67  krnSwMdc_ = 2.0;
68  // krwSwMdc_ = 2.0;
69 
70  deltaSwImbKrn_ = 0.0;
71  // deltaSwImbKrw_ = 0.0;
72  }
73 
78  void finalize()
79  {
80  if (config().enableHysteresis()) {
81  //C_ = 1.0/(Sncri_ - Sncrd_) + 1.0/(Snmaxd_ - Sncrd_);
82 
83  updateDynamicParams_();
84  }
85 
87  }
88 
92  void setConfig(std::shared_ptr<EclHysteresisConfig> value)
93  { config_ = value; }
94 
98  const EclHysteresisConfig& config() const
99  { return *config_; }
100 
104  void setDrainageParams(std::shared_ptr<EffLawParams> value,
105  const EclEpsScalingPointsInfo<Scalar>& /* info */,
106  EclTwoPhaseSystemType /* twoPhaseSystem */)
107 
108  {
109  drainageParams_ = *value;
110  }
111 
115  const EffLawParams& drainageParams() const
116  { return drainageParams_; }
117 
118  EffLawParams& drainageParams()
119  { return drainageParams_; }
120 
124  void setImbibitionParams(std::shared_ptr<EffLawParams> value,
125  const EclEpsScalingPointsInfo<Scalar>& /* info */,
126  EclTwoPhaseSystemType /* twoPhaseSystem */)
127  {
128  imbibitionParams_ = *value;
129 
130 /*
131  if (twoPhaseSystem == EclGasOilSystem) {
132  Sncri_ = info.Sgcr;
133  }
134  else {
135  assert(twoPhaseSystem == EclOilWaterSystem);
136  Sncri_ = info.Sowcr;
137  }
138 */
139  }
140 
144  const EffLawParams& imbibitionParams() const
145  { return imbibitionParams_; }
146 
147  EffLawParams& imbibitionParams()
148  { return imbibitionParams_; }
149 
154  void setPcSwMdc(Scalar value)
155  { pcSwMdc_ = value; }
156 
161  Scalar pcSwMdc() const
162  { return pcSwMdc_; }
163 
169  void setKrwSwMdc(Scalar /* value */)
170  {}
171  // { krwSwMdc_ = value; };
172 
178  Scalar krwSwMdc() const
179  { return 2.0; }
180  // { return krwSwMdc_; };
181 
187  void setKrnSwMdc(Scalar value)
188  { krnSwMdc_ = value; }
189 
195  Scalar krnSwMdc() const
196  { return krnSwMdc_; }
197 
205  void setDeltaSwImbKrw(Scalar /* value */)
206  {}
207  // { deltaSwImbKrw_ = value; }
208 
216  Scalar deltaSwImbKrw() const
217  { return 0.0; }
218 // { return deltaSwImbKrw_; }
219 
227  void setDeltaSwImbKrn(Scalar value)
228  { deltaSwImbKrn_ = value; }
229 
237  Scalar deltaSwImbKrn() const
238  { return deltaSwImbKrn_; }
239 
246  void update(Scalar pcSw, Scalar /* krwSw */, Scalar krnSw)
247  {
248  bool updateParams = false;
249  if (pcSw < pcSwMdc_) {
250  pcSwMdc_ = pcSw;
251  updateParams = true;
252  }
253 
254 /*
255  // This is quite hacky: Eclipse says that it only uses relperm hysteresis for the
256  // wetting phase (indicated by '0' for the second item of the EHYSTER keyword),
257  // even though this makes about zero sense: one would expect that hysteresis can
258  // be limited to the oil phase, but the oil phase is the wetting phase of the
259  // gas-oil twophase system whilst it is non-wetting for water-oil.
260  if (krwSw < krwSwMdc_)
261  {
262  krwSwMdc_ = krwSw;
263  updateParams = true;
264  }
265 */
266 
267  if (krnSw < krnSwMdc_) {
268  krnSwMdc_ = krnSw;
269  updateParams = true;
270  }
271 
272  if (updateParams)
273  updateDynamicParams_();
274  }
275 
276 private:
277  void updateDynamicParams_()
278  {
279  // HACK: Eclipse seems to disable the wetting-phase relperm even though this is
280  // quite pointless from the physical POV. (see comment above)
281 /*
282  // calculate the saturation deltas for the relative permeabilities
283  Scalar krwMdcDrainage = EffLawT::twoPhaseSatKrw(drainageParams(), krwSwMdc_);
284  Scalar SwKrwMdcImbibition = EffLawT::twoPhaseSatKrwInv(imbibitionParams(), krwMdcDrainage);
285  deltaSwImbKrw_ = SwKrwMdcImbibition - krwSwMdc_;
286 */
287 
288  Scalar krnMdcDrainage = EffLawT::twoPhaseSatKrn(drainageParams(), krnSwMdc_);
289  Scalar SwKrnMdcImbibition = EffLawT::twoPhaseSatKrnInv(imbibitionParams(), krnMdcDrainage);
290  deltaSwImbKrn_ = SwKrnMdcImbibition - krnSwMdc_;
291 
292  // Scalar pcMdcDrainage = EffLawT::twoPhaseSatPcnw(drainageParams(), pcSwMdc_);
293  // Scalar SwPcMdcImbibition = EffLawT::twoPhaseSatPcnwInv(imbibitionParams(), pcMdcDrainage);
294  // deltaSwImbPc_ = SwPcMdcImbibition - pcSwMdc_;
295 
296 // assert(std::abs(EffLawT::twoPhaseSatPcnw(imbibitionParams(), pcSwMdc_ + deltaSwImbPc_)
297 // - EffLawT::twoPhaseSatPcnw(drainageParams(), pcSwMdc_)) < 1e-8);
298  assert(std::abs(EffLawT::twoPhaseSatKrn(imbibitionParams(), krnSwMdc_ + deltaSwImbKrn_)
299  - EffLawT::twoPhaseSatKrn(drainageParams(), krnSwMdc_)) < 1e-8);
300 // assert(std::abs(EffLawT::twoPhaseSatKrw(imbibitionParams(), krwSwMdc_ + deltaSwImbKrw_)
301 // - EffLawT::twoPhaseSatKrw(drainageParams(), krwSwMdc_)) < 1e-8);
302 
303 #if 0
304  Scalar Snhy = 1.0 - SwMdc_;
305 
306  Sncrt_ = Sncrd_ + (Snhy - Sncrd_)/(1 + C_*(Snhy - Sncrd_));
307 #endif
308  }
309 
310  std::shared_ptr<EclHysteresisConfig> config_;
311  EffLawParams imbibitionParams_;
312  EffLawParams drainageParams_;
313 
314  // largest wettinging phase saturation which is on the main-drainage curve. These are
315  // three different values because the sourounding code can choose to use different
316  // definitions for the saturations for different quantities
317  //Scalar krwSwMdc_;
318  Scalar krnSwMdc_;
319  Scalar pcSwMdc_;
320 
321  // offsets added to wetting phase saturation uf using the imbibition curves need to
322  // be used to calculate the wetting phase relperm, the non-wetting phase relperm and
323  // the capillary pressure
324  //Scalar deltaSwImbKrw_;
325  Scalar deltaSwImbKrn_;
326  //Scalar deltaSwImbPc_;
327 
328  // trapped non-wetting phase saturation
329  //Scalar Sncrt_;
330 
331  // the following uses the conventions of the Eclipse technical description:
332  //
333  // Sncrd_: critical non-wetting phase saturation for the drainage curve
334  // Sncri_: critical non-wetting phase saturation for the imbibition curve
335  // Snmaxd_: non-wetting phase saturation where the non-wetting relperm reaches its
336  // maximum on the drainage curve
337  // C_: factor required to calculate the trapped non-wetting phase saturation using
338  // the Killough approach
339  //Scalar Sncrd_;
340  //Scalar Sncri_;
341  //Scalar Snmaxd_;
342  //Scalar C_;
343 };
344 
345 } // namespace Opm
346 
347 #endif
Scalar krwSwMdc() const
Get the saturation of the wetting phase where the last switch from the main drainage curve to imbibit...
Definition: EclHysteresisTwoPhaseLawParams.hpp:178
Scalar pcSwMdc() const
Get the saturation of the wetting phase where the last switch from the main drainage curve to imbibit...
Definition: EclHysteresisTwoPhaseLawParams.hpp:161
void setDrainageParams(std::shared_ptr< EffLawParams > value, const EclEpsScalingPointsInfo< Scalar > &, EclTwoPhaseSystemType)
Sets the parameters used for the drainage curve.
Definition: EclHysteresisTwoPhaseLawParams.hpp:104
Scalar deltaSwImbKrn() const
Returns the saturation value which must be added if krn is calculated using the imbibition curve...
Definition: EclHysteresisTwoPhaseLawParams.hpp:237
A default implementation of the parameters for the material law which implements the ECL relative per...
Definition: EclHysteresisTwoPhaseLawParams.hpp:53
void setDeltaSwImbKrw(Scalar)
Sets the saturation value which must be added if krw is calculated using the imbibition curve...
Definition: EclHysteresisTwoPhaseLawParams.hpp:205
This structure represents all values which can be possibly used as scaling points in the endpoint sca...
Definition: EclEpsScalingPoints.hpp:154
void finalize()
Calculate all dependent quantities once the independent quantities of the parameter object have been ...
Definition: EclHysteresisTwoPhaseLawParams.hpp:78
const EffLawParams & drainageParams() const
Returns the parameters used for the drainage curve.
Definition: EclHysteresisTwoPhaseLawParams.hpp:115
void update(Scalar pcSw, Scalar, Scalar krnSw)
Notify the hysteresis law that a given wetting-phase saturation has been seen.
Definition: EclHysteresisTwoPhaseLawParams.hpp:246
Default implementation for asserting finalization of parameter objects.
void setKrwSwMdc(Scalar)
Set the saturation of the wetting phase where the last switch from the main drainage curve (MDC) to i...
Definition: EclHysteresisTwoPhaseLawParams.hpp:169
Specifies the configuration used by the ECL kr/pC hysteresis code.
Definition: EclHysteresisConfig.hpp:50
void setConfig(std::shared_ptr< EclHysteresisConfig > value)
Set the endpoint scaling configuration object.
Definition: EclHysteresisTwoPhaseLawParams.hpp:92
const EffLawParams & imbibitionParams() const
Returns the parameters used for the imbibition curve.
Definition: EclHysteresisTwoPhaseLawParams.hpp:144
Scalar deltaSwImbKrw() const
Returns the saturation value which must be added if krw is calculated using the imbibition curve...
Definition: EclHysteresisTwoPhaseLawParams.hpp:216
void setImbibitionParams(std::shared_ptr< EffLawParams > value, const EclEpsScalingPointsInfo< Scalar > &, EclTwoPhaseSystemType)
Sets the parameters used for the imbibition curve.
Definition: EclHysteresisTwoPhaseLawParams.hpp:124
const EclHysteresisConfig & config() const
Returns the endpoint scaling configuration object.
Definition: EclHysteresisTwoPhaseLawParams.hpp:98
void setDeltaSwImbKrn(Scalar value)
Sets the saturation value which must be added if krn is calculated using the imbibition curve...
Definition: EclHysteresisTwoPhaseLawParams.hpp:227
void finalize()
Mark the object as finalized.
Definition: EnsureFinalized.hpp:77
void setPcSwMdc(Scalar value)
Set the saturation of the wetting phase where the last switch from the main drainage curve (MDC) to i...
Definition: EclHysteresisTwoPhaseLawParams.hpp:154
Default implementation for asserting finalization of parameter objects.
Definition: EnsureFinalized.hpp:46
void setKrnSwMdc(Scalar value)
Set the saturation of the wetting phase where the last switch from the main drainage curve (MDC) to i...
Definition: EclHysteresisTwoPhaseLawParams.hpp:187
Specifies the configuration used by the ECL kr/pC hysteresis code.
Scalar krnSwMdc() const
Get the saturation of the wetting phase where the last switch from the main drainage curve to imbibit...
Definition: EclHysteresisTwoPhaseLawParams.hpp:195