www.mooseframework.org
PorousFlowRelativePermeabilityBase.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 
9 
10 template <>
11 InputParameters
13 {
14  InputParameters params = validParams<PorousFlowMaterialBase>();
15  params.addRangeCheckedParam<Real>(
16  "scaling", 1.0, "scaling>=0", "Relative permeability is multiplied by this factor");
17  params.addRangeCheckedParam<Real>(
18  "s_res",
19  0,
20  "s_res >= 0 & s_res < 1",
21  "The residual saturation of the phase j. Must be between 0 and 1");
22  params.addRangeCheckedParam<Real>(
23  "sum_s_res",
24  0,
25  "sum_s_res >= 0 & sum_s_res < 1",
26  "Sum of residual saturations over all phases. Must be between 0 and 1");
27  params.addClassDescription("Base class for PorousFlow relative permeability materials");
28  return params;
29 }
30 
32  const InputParameters & parameters)
33  : PorousFlowMaterialBase(parameters),
34  _scaling(getParam<Real>("scaling")),
35  _saturation_variable_name(_dictator.saturationVariableNameDummy()),
36  _saturation(_nodal_material
37  ? getMaterialProperty<std::vector<Real>>("PorousFlow_saturation_nodal")
38  : getMaterialProperty<std::vector<Real>>("PorousFlow_saturation_qp")),
39  _relative_permeability(
40  _nodal_material ? declareProperty<Real>("PorousFlow_relative_permeability_nodal" + _phase)
41  : declareProperty<Real>("PorousFlow_relative_permeability_qp" + _phase)),
42  _drelative_permeability_ds(
43  _nodal_material
44  ? declarePropertyDerivative<Real>("PorousFlow_relative_permeability_nodal" + _phase,
45  _saturation_variable_name)
46  : declarePropertyDerivative<Real>("PorousFlow_relative_permeability_qp" + _phase,
47  _saturation_variable_name)),
48  _s_res(getParam<Real>("s_res")),
49  _sum_s_res(getParam<Real>("sum_s_res")),
50  _dseff_ds(1.0 / (1.0 - _sum_s_res))
51 {
52  if (_sum_s_res < _s_res)
53  mooseError("Sum of residual saturations sum_s_res cannot be smaller than s_res in ", name());
54 }
55 
56 void
58 {
59  // Effective saturation
60  Real seff = effectiveSaturation(_saturation[_qp][_phase_num]);
61  Real relperm, drelperm;
62 
63  if (seff < 0.0)
64  {
65  // Relative permeability is 0 for saturation less than residual
66  relperm = 0.0;
67  drelperm = 0.0;
68  }
69  else if (seff >= 0.0 && seff <= 1)
70  {
71  relperm = relativePermeability(seff);
72  drelperm = dRelativePermeability(seff);
73  }
74  else // seff > 1
75  {
76  // Relative permeability is 1 when fully saturated
77  relperm = 1.0;
78  drelperm = 0.0;
79  }
80 
81  _relative_permeability[_qp] = relperm * _scaling;
82  _drelative_permeability_ds[_qp] = drelperm * _dseff_ds * _scaling;
83 }
84 
85 Real
87 {
88  return (saturation - _s_res) / (1.0 - _sum_s_res);
89 }
const Real _scaling
Relative permeability is multiplied by this quantity.
MaterialProperty< Real > & _relative_permeability
Relative permeability material property.
virtual Real dRelativePermeability(Real seff) const =0
Derivative of relative permeability with respect to effective saturation.
const Real _dseff_ds
Derivative of effective saturation with respect to saturation.
PorousFlowRelativePermeabilityBase(const InputParameters &parameters)
const MaterialProperty< std::vector< Real > > & _saturation
Saturation material property.
InputParameters validParams< PorousFlowRelativePermeabilityBase >()
const unsigned int _phase_num
Phase number of fluid.
Base class for all PorousFlow materials that provide phase-dependent properties.
InputParameters validParams< PorousFlowMaterialBase >()
virtual Real effectiveSaturation(Real saturation) const
Effective saturation of fluid phase.
MaterialProperty< Real > & _drelative_permeability_ds
Derivative of relative permeability wrt phase saturation.
virtual Real relativePermeability(Real seff) const =0
Relative permeability equation (must be overriden in derived class)
const Real _sum_s_res
Sum of residual saturations over all phases.
const Real _s_res
Residual saturation of specified phase.
void FORTRAN_CALL() saturation(double &P, double &T, int &N, int &nerr)