www.mooseframework.org
RichardsRelPermBW.C
Go to the documentation of this file.
1 /****************************************************************/
2 /* MOOSE - Multiphysics Object Oriented Simulation Environment */
3 /* */
4 /* All contents are licensed under LGPL V2.1 */
5 /* See LICENSE for full restrictions */
6 /****************************************************************/
7 
8 // "Broadbridge-White" form of relative permeability (P Broadbridge and I White ``Constant rate
9 // rainfall infiltration: A versatile nonlinear model 1. Analytic Solution'', Water Resources
10 // Research 24 (1988) 145-154)
11 //
12 #include "RichardsRelPermBW.h"
13 #include "libmesh/utility.h"
14 
15 template <>
16 InputParameters
18 {
19  InputParameters params = validParams<RichardsRelPerm>();
20  params.addRequiredRangeCheckedParam<Real>(
21  "Sn",
22  "Sn >= 0",
23  "Low saturation. This must be < Ss, and non-negative. This is BW's "
24  "initial effective saturation, below which effective saturation never goes "
25  "in their simulations/models. If Kn=0 then Sn is the immobile saturation.");
26  params.addRangeCheckedParam<Real>(
27  "Ss",
28  1.0,
29  "Ss <= 1",
30  "High saturation. This must be > Sn and <= 1. Effective saturation "
31  "where porepressure = 0. Effective saturation never exceeds this "
32  "value in BW's simulations/models.");
33  params.addRangeCheckedParam<Real>(
34  "Kn", 0.0, "Kn >= 0", "Relative permeability at Seff = Sn. Must be < Ks");
35  params.addRangeCheckedParam<Real>(
36  "Ks", 1.0, "Ks <= 1", "Relative permeability at Seff = Ss. Must be > Kn");
37  params.addRequiredRangeCheckedParam<Real>(
38  "C",
39  "C > 1",
40  "BW's C parameter. Must be > 1. Define s = (seff - Sn)/(Ss - Sn). Then "
41  "relperm = Kn + s^2(c-1)(Kn-Ks)/(c-s) if 0<s<1, otherwise relperm = Kn if "
42  "s<=0, otherwise relperm = Ks if s>=1.");
43  params.addClassDescription("Broadbridge-White form of relative permeability. Define s = (seff - "
44  "Sn)/(Ss - Sn). Then relperm = Kn + s^2(c-1)(Kn-Ks)/(c-s) if 0<s<1, "
45  "otherwise relperm = Kn if s<=0, otherwise relperm = Ks if s>=1.");
46  return params;
47 }
48 
49 RichardsRelPermBW::RichardsRelPermBW(const InputParameters & parameters)
50  : RichardsRelPerm(parameters),
51  _sn(getParam<Real>("Sn")),
52  _ss(getParam<Real>("Ss")),
53  _kn(getParam<Real>("Kn")),
54  _ks(getParam<Real>("Ks")),
55  _c(getParam<Real>("C"))
56 {
57  if (_ss <= _sn)
58  mooseError("In BW relative permeability Sn set to ",
59  _sn,
60  " and Ss set to ",
61  _ss,
62  " but these must obey Ss > Sn");
63  if (_ks <= _kn)
64  mooseError("In BW relative permeability Kn set to ",
65  _kn,
66  " and Ks set to ",
67  _ks,
68  " but these must obey Ks > Kn");
69  _coef = (_ks - _kn) * (_c - 1); // shorthand coefficient
70 }
71 
72 Real
74 {
75  if (seff <= _sn)
76  return _kn;
77 
78  if (seff >= _ss)
79  return _ks;
80 
81  const Real s_internal = (seff - _sn) / (_ss - _sn);
82  const Real krel = _kn + _coef * Utility::pow<2>(s_internal) / (_c - s_internal);
83  return krel;
84 }
85 
86 Real
88 {
89  if (seff <= _sn)
90  return 0.0;
91 
92  if (seff >= _ss)
93  return 0.0;
94 
95  const Real s_internal = (seff - _sn) / (_ss - _sn);
96  const Real krelp = _coef * (2.0 * s_internal / (_c - s_internal) +
97  Utility::pow<2>(s_internal) / Utility::pow<2>(_c - s_internal));
98  return krelp / (_ss - _sn);
99 }
100 
101 Real
103 {
104  if (seff <= _sn)
105  return 0.0;
106 
107  if (seff >= _ss)
108  return 0.0;
109 
110  const Real s_internal = (seff - _sn) / (_ss - _sn);
111  const Real krelpp =
112  _coef * (2.0 / (_c - s_internal) + 4.0 * s_internal / Utility::pow<2>(_c - s_internal) +
113  2.0 * Utility::pow<2>(s_internal) / Utility::pow<3>(_c - s_internal));
114  return krelpp / Utility::pow<2>(_ss - _sn);
115 }
Real relperm(Real seff) const
relative permeability as a function of effective saturation
Base class for Richards relative permeability classes that provide relative permeability as a functio...
RichardsRelPermBW(const InputParameters &parameters)
Real drelperm(Real seff) const
derivative of relative permeability wrt effective saturation
InputParameters validParams< RichardsRelPerm >()
InputParameters validParams< RichardsRelPermBW >()
Real d2relperm(Real seff) const
second derivative of relative permeability wrt effective saturation