www.mooseframework.org
Functions
PorousFlowVanGenuchten Namespace Reference

van Genuchten effective saturation, capillary pressure and relative permeability functions. More...

Functions

Real effectiveSaturation (Real p, Real alpha, Real m)
 Effective saturation as a function of porepressure. More...
 
Real dEffectiveSaturation (Real p, Real alpha, Real m)
 Derivative of effective saturation wrt porepressure. More...
 
Real d2EffectiveSaturation (Real p, Real alpha, Real m)
 Second derivative of effective saturation wrt porepressure. More...
 
Real capillaryPressure (Real seff, Real alpha, Real m, Real pc_max)
 Capillary pressure as a function of effective saturation. More...
 
Real dCapillaryPressure (Real seff, Real alpha, Real m, Real pc_max)
 Derivative of capillary pressure wrt effective saturation. More...
 
Real d2CapillaryPressure (Real seff, Real alpha, Real m, Real pc_max)
 Second derivative of capillary pressure wrt effective saturation. More...
 
Real relativePermeability (Real seff, Real m)
 Relative permeability as a function of effective saturation. More...
 
Real dRelativePermeability (Real seff, Real m)
 Derivative of relative permeability with respect to effective saturation. More...
 
Real d2RelativePermeability (Real seff, Real m)
 Second derivative of relative permeability with respect to effective saturation. More...
 

Detailed Description

van Genuchten effective saturation, capillary pressure and relative permeability functions.

Note: effective saturation is provided as a function of porepressure, not capillary pressure. Note: capillary pressure and relative permeability are functions of effective saturation. The derivatives are therefore given wrt effective saturation. These derivatives must be multiplied by the derivative of effective saturation wrt the true saturation in objects using these relations.

Based on van Genuchten, M. Th., A closed for equation for predicting the hydraulic conductivity of unsaturated soils, Soil Sci. Soc., 44, 892-898 (1980).

Function Documentation

Real PorousFlowVanGenuchten::capillaryPressure ( Real  seff,
Real  alpha,
Real  m,
Real  pc_max 
)

Capillary pressure as a function of effective saturation.

Parameters
seffeffective saturation
alphavan Genuchten parameter
mvan Genuchten exponent
pc_maxmaximum capillary pressure (Pa)
Returns
capillary pressure (Pa)

Definition at line 61 of file PorousFlowVanGenuchten.C.

Referenced by PorousFlow2PhasePS_VG::capillaryPressure(), and PorousFlowCapillaryPressureVG::capillaryPressureCurve().

62 {
63  if (seff >= 1.0)
64  return 0.0;
65  else if (seff <= 0.0)
66  return pc_max;
67  else
68  {
69  Real a = std::pow(seff, -1.0 / m) - 1.0;
70  return std::min(std::pow(a, 1.0 - m) / alpha, pc_max);
71  }
72 }
ExpressionBuilder::EBTerm pow(const ExpressionBuilder::EBTerm &left, T exponent)
Real PorousFlowVanGenuchten::d2CapillaryPressure ( Real  seff,
Real  alpha,
Real  m,
Real  pc_max 
)

Second derivative of capillary pressure wrt effective saturation.

Parameters
seffeffective saturation
alphavan Genuchten parameter
mvan Genuchten exponent
pc_maxmaximum capillary pressure (Pa)
Returns
second derivative of capillary pressure wrt effective saturation

Definition at line 91 of file PorousFlowVanGenuchten.C.

Referenced by PorousFlow2PhasePS_VG::d2CapillaryPressure_dS2(), and PorousFlowCapillaryPressureVG::d2CapillaryPressureCurve().

92 {
93  if (seff <= 0.0 || seff >= 1.0)
94  return 0.0;
95  else
96  {
97  Real a = std::pow(seff, -1.0 / m) - 1.0;
98  // Return 0 if pc > pc_max
99  if (std::pow(a, 1.0 - m) / alpha > pc_max)
100  return 0.0;
101  else
102  {
103  Real d2pc = -std::pow(a, -1.0 - m) * std::pow(seff, -2.0 - 2.0 / m) +
104  ((1.0 + m) / m) * std::pow(a, -m) * std::pow(seff, -1.0 / m - 2.0);
105  d2pc *= (1.0 - m) / m / alpha;
106  return d2pc;
107  }
108  }
109 }
ExpressionBuilder::EBTerm pow(const ExpressionBuilder::EBTerm &left, T exponent)
Real PorousFlowVanGenuchten::d2EffectiveSaturation ( Real  p,
Real  alpha,
Real  m 
)

Second derivative of effective saturation wrt porepressure.

Parameters
pporepressure
alphavan Genuchten parameter
mvan Genuchten exponent
Returns
second derivative of effective saturation wrt porepressure

Definition at line 44 of file PorousFlowVanGenuchten.C.

Referenced by PorousFlowCapillaryPressureVG::d2EffectiveSaturation(), PorousFlow1PhaseP_VG::d2EffectiveSaturation_dP2(), and PorousFlow2PhasePP_VG::d2EffectiveSaturation_dP2().

45 {
46  if (p >= 0.0)
47  return 0.0;
48  else
49  {
50  Real n = 1.0 / (1.0 - m);
51  Real inner = 1.0 + std::pow(-alpha * p, n);
52  Real dinner_dp = -n * alpha * std::pow(-alpha * p, n - 1.0);
53  Real d2inner_dp2 = n * (n - 1.0) * alpha * alpha * std::pow(-alpha * p, n - 2.0);
54  Real d2seff_dp2 = m * (m + 1.0) * std::pow(inner, -m - 2.0) * std::pow(dinner_dp, 2.0) -
55  m * std::pow(inner, -m - 1.0) * d2inner_dp2;
56  return d2seff_dp2;
57  }
58 }
ExpressionBuilder::EBTerm pow(const ExpressionBuilder::EBTerm &left, T exponent)
Real PorousFlowVanGenuchten::d2RelativePermeability ( Real  seff,
Real  m 
)

Second derivative of relative permeability with respect to effective saturation.

Parameters
seffeffective saturation
mvan Genuchten exponent
Returns
second derivative of relative permeability wrt effective saturation

Definition at line 141 of file PorousFlowVanGenuchten.C.

142 {
143  // Guard against division by zero
144  if (seff <= 0.0 || seff >= 1.0)
145  return 0.0;
146 
147  const Real a = 1.0 - std::pow(seff, 1.0 / m);
148  const Real da = -1.0 / m * std::pow(seff, 1.0 / m - 1.0);
149  const Real d2a = -(1.0 / m) * (1.0 / m - 1.0) * std::pow(seff, 1.0 / m - 2.0);
150  const Real b = 1.0 - std::pow(a, m);
151  const Real db = -m * std::pow(a, m - 1.0) * da;
152  const Real d2b = -m * (m - 1.0) * std::pow(a, m - 2.0) * da * da - m * std::pow(a, m - 1.0) * d2a;
153 
154  return -0.25 * std::pow(seff, -1.5) * Utility::pow<2>(b) + 2.0 * std::pow(seff, -0.5) * b * db +
155  2.0 * std::sqrt(seff) * db * db + 2.0 * std::sqrt(seff) * b * d2b;
156 }
ExpressionBuilder::EBTerm pow(const ExpressionBuilder::EBTerm &left, T exponent)
Real PorousFlowVanGenuchten::dCapillaryPressure ( Real  seff,
Real  alpha,
Real  m,
Real  pc_max 
)

Derivative of capillary pressure wrt effective saturation.

Parameters
seffeffective saturation
alphavan Genuchten parameter
mvan Genuchten exponent
pc_maxmaximum capillary pressure (Pa)
Returns
derivative of capillary pressure wrt effective saturation

Definition at line 75 of file PorousFlowVanGenuchten.C.

Referenced by PorousFlow2PhasePS_VG::dCapillaryPressure_dS(), and PorousFlowCapillaryPressureVG::dCapillaryPressureCurve().

76 {
77  if (seff <= 0.0 || seff >= 1.0)
78  return 0.0;
79  else
80  {
81  Real a = std::pow(seff, -1.0 / m) - 1.0;
82  // Return 0 if pc > pc_max
83  if (std::pow(a, 1.0 - m) / alpha > pc_max)
84  return 0.0;
85  else
86  return (m - 1.0) * std::pow(a, -m) * std::pow(seff, -1.0 - 1.0 / m) / m / alpha;
87  }
88 }
ExpressionBuilder::EBTerm pow(const ExpressionBuilder::EBTerm &left, T exponent)
Real PorousFlowVanGenuchten::dEffectiveSaturation ( Real  p,
Real  alpha,
Real  m 
)

Derivative of effective saturation wrt porepressure.

Parameters
pporepressure
alphavan Genuchten parameter
mvan Genuchten exponent
Returns
derivative of effective saturation wrt porepressure

Definition at line 29 of file PorousFlowVanGenuchten.C.

Referenced by PorousFlowCapillaryPressureVG::dEffectiveSaturation(), PorousFlow1PhaseP_VG::dEffectiveSaturation_dP(), and PorousFlow2PhasePP_VG::dEffectiveSaturation_dP().

30 {
31  if (p >= 0.0)
32  return 0.0;
33  else
34  {
35  Real n = 1.0 / (1.0 - m);
36  Real inner = 1.0 + std::pow(-alpha * p, n);
37  Real dinner_dp = -n * alpha * std::pow(-alpha * p, n - 1.0);
38  Real dseff_dp = -m * std::pow(inner, -m - 1) * dinner_dp;
39  return dseff_dp;
40  }
41 }
ExpressionBuilder::EBTerm pow(const ExpressionBuilder::EBTerm &left, T exponent)
Real PorousFlowVanGenuchten::dRelativePermeability ( Real  seff,
Real  m 
)

Derivative of relative permeability with respect to effective saturation.

Parameters
seffeffective saturation
mvan Genuchten exponent
Returns
derivative of relative permeability wrt effective saturation

Definition at line 126 of file PorousFlowVanGenuchten.C.

Referenced by PorousFlowRelativePermeabilityVG::dRelativePermeability().

127 {
128  // Guard against division by zero
129  if (seff <= 0.0 || seff >= 1.0)
130  return 0.0;
131 
132  const Real a = 1.0 - std::pow(seff, 1.0 / m);
133  const Real da = -1.0 / m * std::pow(seff, 1.0 / m - 1.0);
134  const Real b = 1.0 - std::pow(a, m);
135  const Real db = -m * std::pow(a, m - 1.0) * da;
136 
137  return 0.5 * std::pow(seff, -0.5) * Utility::pow<2>(b) + 2.0 * std::sqrt(seff) * b * db;
138 }
ExpressionBuilder::EBTerm pow(const ExpressionBuilder::EBTerm &left, T exponent)
Real PorousFlowVanGenuchten::effectiveSaturation ( Real  p,
Real  alpha,
Real  m 
)

Effective saturation as a function of porepressure.

Note: seff = 1 for p >= 0

Parameters
pporepressure
alphavan Genuchten parameter
mvan Genuchten exponent
Returns
effective saturation

Definition at line 14 of file PorousFlowVanGenuchten.C.

Referenced by PorousFlow1PhaseP_VG::effectiveSaturation(), PorousFlow2PhasePP_VG::effectiveSaturation(), and PorousFlowCapillaryPressureVG::effectiveSaturation().

15 {
16  Real n, seff;
17 
18  if (p >= 0.0)
19  return 1.0;
20  else
21  {
22  n = 1.0 / (1.0 - m);
23  seff = 1.0 + std::pow(-alpha * p, n);
24  return std::pow(seff, -m);
25  }
26 }
ExpressionBuilder::EBTerm pow(const ExpressionBuilder::EBTerm &left, T exponent)
Real PorousFlowVanGenuchten::relativePermeability ( Real  seff,
Real  m 
)

Relative permeability as a function of effective saturation.

Parameters
seffeffective saturation
mvan Genuchten exponent
Returns
relative permeability

Definition at line 112 of file PorousFlowVanGenuchten.C.

Referenced by PorousFlowRelativePermeabilityVG::relativePermeability().

113 {
114  if (seff <= 0.0)
115  return 0.0;
116  else if (seff >= 1.0)
117  return 1.0;
118 
119  const Real a = 1.0 - std::pow(seff, 1.0 / m);
120  const Real b = 1.0 - std::pow(a, m);
121 
122  return std::sqrt(seff) * Utility::pow<2>(b);
123 }
ExpressionBuilder::EBTerm pow(const ExpressionBuilder::EBTerm &left, T exponent)