www.mooseframework.org
PorousFlowBrooksCorey.h
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://www.mooseframework.org
3 //*
4 //* All rights reserved, see COPYRIGHT for full restrictions
5 //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 //*
7 //* Licensed under LGPL 2.1, please see LICENSE for details
8 //* https://www.gnu.org/licenses/lgpl-2.1.html
9 
10 #pragma once
11 
12 #include "MooseTypes.h"
13 
28 {
37 Real effectiveSaturation(Real pc, Real pe, Real lambda);
38 
46 Real dEffectiveSaturation(Real pc, Real pe, Real lambda);
47 
55 Real d2EffectiveSaturation(Real pc, Real pe, Real lambda);
56 
66 Real capillaryPressure(Real seff, Real pe, Real lambda, Real pc_max);
67 
77 Real dCapillaryPressure(Real seff, Real pe, Real lambda, Real pc_max);
78 
88 Real d2CapillaryPressure(Real seff, Real pe, Real lambda, Real pc_max);
89 
96 template <typename T>
97 T
98 relativePermeabilityW(const T & seff, Real lambda)
99 {
100  if (MetaPhysicL::raw_value(seff) <= 0.0)
101  return 0.0;
102  else if (MetaPhysicL::raw_value(seff) >= 1.0)
103  return 1.0;
104 
105  return std::pow(seff, (2.0 + 3.0 * lambda) / lambda);
106 }
107 
114 Real dRelativePermeabilityW(Real seff, Real lambda);
115 
122 template <typename T>
123 T
124 relativePermeabilityNW(const T & seff, Real lambda)
125 {
126  if (MetaPhysicL::raw_value(seff) <= 0.0)
127  return 0.0;
128  else if (MetaPhysicL::raw_value(seff) >= 1.0)
129  return 1.0;
130 
131  return seff * seff * (1.0 - std::pow(1.0 - seff, (2.0 + lambda) / lambda));
132 }
133 
140 Real dRelativePermeabilityNW(Real seff, Real lambda);
141 }
Real capillaryPressure(Real seff, Real pe, Real lambda, Real pc_max)
Capillary pressure as a function of effective saturation.
Real dRelativePermeabilityNW(Real seff, Real lambda)
Derivative of relative permeability of the non-wetting phase wrt to effective saturation.
Real effectiveSaturation(Real pc, Real pe, Real lambda)
Effective saturation as a function of capillary pressure Note: seff = 1 for p >= 0.
T relativePermeabilityW(const T &seff, Real lambda)
Relative permeability of the wetting phase as a function of effective saturation. ...
Real dCapillaryPressure(Real seff, Real pe, Real lambda, Real pc_max)
Derivative of capillary pressure wrt effective saturation.
T relativePermeabilityNW(const T &seff, Real lambda)
Relative permeability of the non-wetting phase as a function of effective saturation.
auto raw_value(const Eigen::Map< T > &in)
Real d2CapillaryPressure(Real seff, Real pe, Real lambda, Real pc_max)
Second derivative of capillary pressure wrt effective saturation.
Real dEffectiveSaturation(Real pc, Real pe, Real lambda)
Derivative of effective saturation wrt porepressure.
Real dRelativePermeabilityW(Real seff, Real lambda)
Derivative of relative permeability of the wetting phase wrt to effective saturation.
Real d2EffectiveSaturation(Real pc, Real pe, Real lambda)
Second derivative of effective saturation wrt porepressure.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
MooseUnits pow(const MooseUnits &, int)
Brooks-Corey effective saturation, capillary pressure and relative permeability functions.