www.mooseframework.org
Functions
PorousFlowBroadbridgeWhite Namespace Reference

Broadbridge-White version of relative permeability, and effective saturation as a function of capillary pressure. More...

Functions

Real LambertW (Real z)
 Provides the Lambert W function, which satisfies W(z) * exp(W(z)) = z. More...
 
Real effectiveSaturation (Real pc, Real c, Real sn, Real ss, Real las)
 Effective saturation as a function of capillary pressure If pc>=0 this will yield 1, otherwise it will yield <1. More...
 
Real dEffectiveSaturation (Real pc, Real c, Real sn, Real ss, Real las)
 Derivative of effective saturation wrt capillary pressure. More...
 
Real d2EffectiveSaturation (Real pc, Real c, Real sn, Real ss, Real las)
 Second derivative of effective saturation wrt capillary pressure. More...
 
Real relativePermeability (Real s, Real c, Real sn, Real ss, Real kn, Real ks)
 Relative permeability as a function of saturation. More...
 
Real dRelativePermeability (Real s, Real c, Real sn, Real ss, Real kn, Real ks)
 Derivative of relative permeability with respect to saturation. More...
 
Real d2RelativePermeability (Real s, Real c, Real sn, Real ss, Real kn, Real ks)
 Second derivative of relative permeability with respect to saturation. More...
 

Detailed Description

Broadbridge-White version of relative permeability, and effective saturation as a function of capillary pressure.

Valid for Kn small. P Broadbridge, I White ``Constant rate rainfall infiltration: A versatile nonlinear model, 1 Analytical solution''. Water Resources Research 24 (1988) 145–154.

Function Documentation

Real PorousFlowBroadbridgeWhite::d2EffectiveSaturation ( Real  pc,
Real  c,
Real  sn,
Real  ss,
Real  las 
)

Second derivative of effective saturation wrt capillary pressure.

Parameters
pccapillary pressure
cBW's C parameter
snBW's S_n parameter
ssBW's S_s parameter
lasBW's lambda_s parameter
Returns
second derivative of effective saturation wrt capillary pressure

Definition at line 112 of file PorousFlowBroadbridgeWhite.C.

Referenced by PorousFlowCapillaryPressureBW::d2EffectiveSaturation(), and PorousFlow1PhaseP_BW::d2EffectiveSaturation_dP2().

113 {
114  if (pressure >= 0.0)
115  return 0.0;
116  const Real x = (c - 1.0) * std::exp(c - 1.0 - c * pressure / las);
117  const Real lamw = LambertW(x);
118  return -(ss - sn) * Utility::pow<3>(c) / Utility::pow<2>(las) * lamw * (1.0 - 2.0 * lamw) /
119  Utility::pow<5>(1.0 + lamw);
120 }
Real LambertW(Real z)
Provides the Lambert W function, which satisfies W(z) * exp(W(z)) = z.
const std::string pressure
Definition: NS.h:24
Real PorousFlowBroadbridgeWhite::d2RelativePermeability ( Real  s,
Real  c,
Real  sn,
Real  ss,
Real  kn,
Real  ks 
)

Second derivative of relative permeability with respect to saturation.

Parameters
ssaturation
cBW's C parameter
snBW's S_n parameter
ssBW's S_s parameter
knBW's K_n parameter
ksBW's K_s parameter
Returns
second derivative of relative permeability wrt saturation

Definition at line 153 of file PorousFlowBroadbridgeWhite.C.

154 {
155  if (s <= sn)
156  return 0.0;
157 
158  if (s >= ss)
159  return 0.0;
160 
161  const Real coef = (ks - kn) * (c - 1.0);
162  const Real th = (s - sn) / (ss - sn);
163  const Real krelpp = coef * (2.0 / (c - th) + 4.0 * th / Utility::pow<2>(c - th) +
164  2.0 * Utility::pow<2>(th) / Utility::pow<3>(c - th));
165  return krelpp / Utility::pow<2>(ss - sn);
166 }
Real PorousFlowBroadbridgeWhite::dEffectiveSaturation ( Real  pc,
Real  c,
Real  sn,
Real  ss,
Real  las 
)

Derivative of effective saturation wrt capillary pressure.

Parameters
pccapillary pressure
cBW's C parameter
snBW's S_n parameter
ssBW's S_s parameter
lasBW's lambda_s parameter
Returns
derivative of effective saturation wrt capillary pressure

Definition at line 102 of file PorousFlowBroadbridgeWhite.C.

Referenced by PorousFlowCapillaryPressureBW::dEffectiveSaturation(), and PorousFlow1PhaseP_BW::dEffectiveSaturation_dP().

103 {
104  if (pressure >= 0.0)
105  return 0.0;
106  const Real x = (c - 1.0) * std::exp(c - 1.0 - c * pressure / las);
107  const Real lamw = LambertW(x);
108  return (ss - sn) * Utility::pow<2>(c) / las * lamw / Utility::pow<3>(1.0 + lamw);
109 }
Real LambertW(Real z)
Provides the Lambert W function, which satisfies W(z) * exp(W(z)) = z.
const std::string pressure
Definition: NS.h:24
Real PorousFlowBroadbridgeWhite::dRelativePermeability ( Real  s,
Real  c,
Real  sn,
Real  ss,
Real  kn,
Real  ks 
)

Derivative of relative permeability with respect to saturation.

Parameters
ssaturation
cBW's C parameter
snBW's S_n parameter
ssBW's S_s parameter
knBW's K_n parameter
ksBW's K_s parameter
Returns
derivative of relative permeability wrt saturation

Definition at line 138 of file PorousFlowBroadbridgeWhite.C.

Referenced by PorousFlowRelativePermeabilityBW::dRelativePermeability().

139 {
140  if (s <= sn)
141  return 0.0;
142 
143  if (s >= ss)
144  return 0.0;
145 
146  const Real coef = (ks - kn) * (c - 1.0);
147  const Real th = (s - sn) / (ss - sn);
148  const Real krelp = coef * (2.0 * th / (c - th) + Utility::pow<2>(th) / Utility::pow<2>(c - th));
149  return krelp / (ss - sn);
150 }
Real PorousFlowBroadbridgeWhite::effectiveSaturation ( Real  pc,
Real  c,
Real  sn,
Real  ss,
Real  las 
)

Effective saturation as a function of capillary pressure If pc>=0 this will yield 1, otherwise it will yield <1.

Parameters
pccapillary pressure
cBW's C parameter
snBW's S_n parameter
ssBW's S_s parameter
lasBW's lambda_s parameter
Returns
effective saturation

Definition at line 92 of file PorousFlowBroadbridgeWhite.C.

Referenced by PorousFlowCapillaryPressureBW::effectiveSaturation(), and PorousFlow1PhaseP_BW::effectiveSaturation().

93 {
94  if (pressure >= 0.0)
95  return 1.0;
96  const Real x = (c - 1.0) * std::exp(c - 1.0 - c * pressure / las);
97  const Real th = c / (1.0 + LambertW(x)); // use branch 0 for positive x
98  return sn + (ss - sn) * th;
99 }
Real LambertW(Real z)
Provides the Lambert W function, which satisfies W(z) * exp(W(z)) = z.
const std::string pressure
Definition: NS.h:24
Real PorousFlowBroadbridgeWhite::LambertW ( Real  z)

Provides the Lambert W function, which satisfies W(z) * exp(W(z)) = z.

Parameters
zz value, which must be positive
Returns
W

Definition at line 14 of file PorousFlowBroadbridgeWhite.C.

Referenced by d2EffectiveSaturation(), dEffectiveSaturation(), and effectiveSaturation().

15 {
16  /* Lambert W function.
17  Was ~/C/LambertW.c written K M Briggs Keith dot Briggs at bt dot com 97 May 21.
18  Revised KMB 97 Nov 20; 98 Feb 11, Nov 24, Dec 28; 99 Jan 13; 00 Feb 23; 01 Apr 09
19 
20  Computes Lambert W function, principal branch.
21  See LambertW1.c for -1 branch.
22 
23  Returned value W(z) satisfies W(z)*exp(W(z))=z
24  test data...
25  W(1)= 0.5671432904097838730
26  W(2)= 0.8526055020137254914
27  W(20)=2.2050032780240599705
28  To solve (a+b*R)*exp(-c*R)-d=0 for R, use
29  R=-(b*W(-exp(-a*c/b)/b*d*c)+a*c)/b/c
30 
31  Test:
32  gcc -DTESTW LambertW.c -o LambertW -lm && LambertW
33  Library:
34  gcc -O3 -c LambertW.c
35 
36  Modified trivially by Andy to use MOOSE things
37  */
38  mooseAssert(z > 0, "LambertW function in RichardsSeff1BWsmall called with negative argument");
39 
40  const Real eps = 4.0e-16; //, em1=0.3678794411714423215955237701614608;
41  Real p, e, t, w;
42 
43  /* Uncomment this stuff is you ever need to call with a negative argument
44  if (z < -em1)
45  mooseError("LambertW: bad argument ", z, "\n");
46 
47  if (0.0 == z)
48  return 0.0;
49  if (z < -em1+1e-4)
50  {
51  // series near -em1 in sqrt(q)
52  Real q=z+em1,r=std::sqrt(q),q2=q*q,q3=q2*q;
53  return
54  -1.0
55  +2.331643981597124203363536062168*r
56  -1.812187885639363490240191647568*q
57  +1.936631114492359755363277457668*r*q
58  -2.353551201881614516821543561516*q2
59  +3.066858901050631912893148922704*r*q2
60  -4.175335600258177138854984177460*q3
61  +5.858023729874774148815053846119*r*q3
62  -8.401032217523977370984161688514*q3*q; // error approx 1e-16
63  }
64  */
65  /* initial approx for iteration... */
66  if (z < 1.0)
67  {
68  /* series near 0 */
69  p = std::sqrt(2.0 * (2.7182818284590452353602874713526625 * z + 1.0));
70  w = -1.0 + p * (1.0 + p * (-0.333333333333333333333 + p * 0.152777777777777777777777));
71  }
72  else
73  w = std::log(z); /* asymptotic */
74  if (z > 3.0)
75  w -= std::log(w); /* useful? */
76  for (unsigned int i = 0; i < 10; ++i)
77  {
78  /* Halley iteration */
79  e = std::exp(w);
80  t = w * e - z;
81  p = w + 1.0;
82  t /= e * p - 0.5 * (p + 1.0) * t / p;
83  w -= t;
84  if (std::abs(t) < eps * (1.0 + std::abs(w)))
85  return w; /* rel-abs error */
86  }
87  /* should never get here */
88  mooseError("LambertW: No convergence at z= ", z, "\n");
89 }
Real PorousFlowBroadbridgeWhite::relativePermeability ( Real  s,
Real  c,
Real  sn,
Real  ss,
Real  kn,
Real  ks 
)

Relative permeability as a function of saturation.

Parameters
ssaturation
cBW's C parameter
snBW's S_n parameter
ssBW's S_s parameter
knBW's K_n parameter
ksBW's K_s parameter
Returns
relative permeability

Definition at line 123 of file PorousFlowBroadbridgeWhite.C.

Referenced by PorousFlowRelativePermeabilityBW::relativePermeability().

124 {
125  if (s <= sn)
126  return kn;
127 
128  if (s >= ss)
129  return ks;
130 
131  const Real coef = (ks - kn) * (c - 1.0);
132  const Real th = (s - sn) / (ss - sn);
133  const Real krel = kn + coef * Utility::pow<2>(th) / (c - th);
134  return krel;
135 }