PorousFlow1PhaseHysP

This Material computes saturation given porepressure in single-phase, partially-saturated, hysteretic situations, assuming a van Genuchten relationship. Detailed documentation about hysteresis and the van Genuchten relationship can be found on the hysteresis page. It is the hysteretic cousin of PorousFlow1PhaseP

This Material requires a PorousFlowHysteresisOrder Material, which will compute the hysteresis order. PorousFlow1PhaseHysP then uses the hysteresis order along with the single-phase porepressure to compute the saturation.

A simple example usage is:

# Simple example of a 1-phase situation with hysteretic capillary pressure.  Water is removed and added to the system in order to observe the hysteresis
[Mesh<<<{"href": "../../syntax/Mesh/index.html"}>>>]
  [mesh]
    type = GeneratedMeshGenerator<<<{"description": "Create a line, square, or cube mesh with uniformly spaced or biased elements.", "href": "../meshgenerators/GeneratedMeshGenerator.html"}>>>
    dim<<<{"description": "The dimension of the mesh to be generated"}>>> = 1
  []
[]

[GlobalParams<<<{"href": "../../syntax/GlobalParams/index.html"}>>>]
  PorousFlowDictator = dictator
[]

[UserObjects<<<{"href": "../../syntax/UserObjects/index.html"}>>>]
  [dictator]
    type = PorousFlowDictator<<<{"description": "Holds information on the PorousFlow variable names", "href": "../userobjects/PorousFlowDictator.html"}>>>
    number_fluid_phases<<<{"description": "The number of fluid phases in the simulation"}>>> = 1
    number_fluid_components<<<{"description": "The number of fluid components in the simulation"}>>> = 1
    porous_flow_vars<<<{"description": "List of primary variables that are used in the PorousFlow simulation.  Jacobian entries involving derivatives wrt these variables will be computed.  In single-phase models you will just have one (eg 'pressure'), in two-phase models you will have two (eg 'p_water p_gas', or 'p_water s_water'), etc."}>>> = 'pp'
  []
[]

[Variables<<<{"href": "../../syntax/Variables/index.html"}>>>]
  [pp]
    initial_condition<<<{"description": "Specifies a constant initial condition for this variable"}>>> = 0
  []
[]

[Kernels<<<{"href": "../../syntax/Kernels/index.html"}>>>]
  [mass_conservation]
    type = PorousFlowMassTimeDerivative<<<{"description": "Derivative of fluid-component mass with respect to time.  Mass lumping to the nodes is used.", "href": "../kernels/PorousFlowMassTimeDerivative.html"}>>>
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = pp
  []
[]

[DiracKernels<<<{"href": "../../syntax/DiracKernels/index.html"}>>>]
  [pump]
    type = PorousFlowPointSourceFromPostprocessor<<<{"description": "Point source (or sink) that adds (or removes) fluid at a mass flux rate specified by a postprocessor.", "href": "../dirackernels/PorousFlowPointSourceFromPostprocessor.html"}>>>
    mass_flux<<<{"description": "The postprocessor name holding the mass flux at this point in kg/s (positive is flux in, negative is flux out)"}>>> = flux
    point<<<{"description": "The x,y,z coordinates of the point source (or sink)"}>>> = '0.5 0 0'
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = pp
  []
[]

[AuxVariables<<<{"href": "../../syntax/AuxVariables/index.html"}>>>]
  [sat]
    family<<<{"description": "Specifies the family of FE shape functions to use for this variable"}>>> = MONOMIAL
    order<<<{"description": "Specifies the order of the FE shape function to use for this variable (additional orders not listed are allowed)"}>>> = CONSTANT
  []
  [hys_order]
    family<<<{"description": "Specifies the family of FE shape functions to use for this variable"}>>> = MONOMIAL
    order<<<{"description": "Specifies the order of the FE shape function to use for this variable (additional orders not listed are allowed)"}>>> = CONSTANT
  []
[]

[AuxKernels<<<{"href": "../../syntax/AuxKernels/index.html"}>>>]
  [sat]
    type = PorousFlowPropertyAux<<<{"description": "AuxKernel to provide access to properties evaluated at quadpoints. Note that elemental AuxVariables must be used, so that these properties are integrated over each element.", "href": "../auxkernels/PorousFlowPropertyAux.html"}>>>
    variable<<<{"description": "The name of the variable that this object applies to"}>>> = sat
    property<<<{"description": "The fluid property that this auxillary kernel is to calculate"}>>> = saturation
  []
  [hys_order]
    type = PorousFlowPropertyAux<<<{"description": "AuxKernel to provide access to properties evaluated at quadpoints. Note that elemental AuxVariables must be used, so that these properties are integrated over each element.", "href": "../auxkernels/PorousFlowPropertyAux.html"}>>>
    variable<<<{"description": "The name of the variable that this object applies to"}>>> = hys_order
    property<<<{"description": "The fluid property that this auxillary kernel is to calculate"}>>> = hysteresis_order
  []
[]

[FluidProperties<<<{"href": "../../syntax/FluidProperties/index.html"}>>>]
  [simple_fluid]
    type = SimpleFluidProperties<<<{"description": "Fluid properties for a simple fluid with a constant bulk density", "href": "../fluidproperties/SimpleFluidProperties.html"}>>>
  []
[]

[Materials<<<{"href": "../../syntax/Materials/index.html"}>>>]
  [porosity]
    type = PorousFlowPorosityConst<<<{"description": "This Material calculates the porosity assuming it is constant", "href": "PorousFlowPorosityConst.html"}>>>
    porosity<<<{"description": "The porosity (assumed indepenent of porepressure, temperature, strain, etc, for this material).  This should be a real number, or a constant monomial variable (not a linear lagrange or other kind of variable)."}>>> = 0.1
  []
  [temperature]
    type = PorousFlowTemperature<<<{"description": "Material to provide temperature at the quadpoints or nodes and derivatives of it with respect to the PorousFlow variables", "href": "PorousFlowTemperature.html"}>>>
    temperature<<<{"description": "Fluid temperature variable.  Note, the default is suitable if your simulation is using Kelvin units, but probably not for Celsius"}>>> = 20
  []
  [massfrac]
    type = PorousFlowMassFraction<<<{"description": "This Material forms a std::vector<std::vector ...> of mass-fractions out of the individual mass fractions", "href": "PorousFlowMassFraction.html"}>>>
  []
  [simple_fluid]
    type = PorousFlowSingleComponentFluid<<<{"description": "This Material calculates fluid properties at the quadpoints or nodes for a single component fluid", "href": "PorousFlowSingleComponentFluid.html"}>>>
    fp<<<{"description": "The name of the user object for fluid properties"}>>> = simple_fluid
    phase<<<{"description": "The phase number"}>>> = 0
  []
  [hys_order_material]
    type = PorousFlowHysteresisOrder<<<{"description": "Computes hysteresis order for use in hysteretic capillary pressures and relative permeabilities", "href": "PorousFlowHysteresisOrder.html"}>>>
  []
  [pc_calculator]
    type = PorousFlow1PhaseHysP<<<{"description": "This Material is used for unsaturated single-phase situations where porepressure is the primary variable and the capillary pressure is hysteretic.  The hysteretic formulation assumes that the single phase is a liquid", "href": "PorousFlow1PhaseHysP.html"}>>>
    alpha_d<<<{"description": "van Genuchten alpha parameter for the primary drying curve.  If using standard units, this is measured in Pa^-1.  Suggested value is around 1E-5"}>>> = 10.0
    alpha_w<<<{"description": "van Genuchten alpha parameter for the primary wetting curve.  If using standard units, this is measured in Pa^-1.  Suggested value is around 1E-5"}>>> = 7.0
    n_d<<<{"description": "van Genuchten n parameter for the primary drying curve.  Dimensionless.  Suggested value is around 2"}>>> = 1.5
    n_w<<<{"description": "van Genuchten n parameter for the primary wetting curve.  Dimensionless.  Suggested value is around 2"}>>> = 1.9
    S_l_min<<<{"description": "Minimum liquid saturation for which the van Genuchten expression is valid.  If no lower extension is used then Pc = Pc_max for liquid saturation <= S_l_min"}>>> = 0.1
    S_lr<<<{"description": "Liquid residual saturation where the liquid relative permeability is zero.  This is used in the Land expression to find S_gr_del.  Almost definitely you need to set S_lr equal to the quantity used for your relative-permeability curves.  Almost definitely you should set S_lr > S_l_min"}>>> = 0.2
    S_gr_max<<<{"description": "Residual gas saturation.  1 - S_gr_max is the maximum saturation for which the van Genuchten expression is valid for the wetting curve.  You must ensure S_gr_max < 1 - S_l_min.  Often S_gr_max = -0.3136 * ln(porosity) - 0.1334 is used"}>>> = 0.3
    Pc_max<<<{"description": "Value of capillary pressure at which the lower extension commences.  The default value means capillary pressure uses the van Genuchten expression for S > S_l_min and is 'infinity' for S <= S_l_min.  This will result in poor convergence around S = S_l_min"}>>> = 12.0
    high_ratio<<<{"description": "The extension to the wetting curves commences at high_ratio * (1 - S_gr_del), where 1 - S_gr_del is the value of the liquid saturation when Pc = 0 (on the wetting curve)"}>>> = 0.9
    low_extension_type<<<{"description": "Type of extension to use for small liquid saturation values.  The extensions modify the capillary pressure for all saturation values less than S(Pc_max).  That is, if the van Genuchten expression would produce Pc > Pc_max, then the extension is used instead.  NONE: Simply cut-off the capillary-pressure at Pc_max, so that Pc <= Pc_max for all S.  QUADRATIC: Pc is a quadratic in S that is continuous and differentiable at S(Pc_max) and has zero derivative at S = 0 (hence, its value at S = 0 will be greater than Pc_max).  EXPONENTIAL: Pc is an exponential in S that is continuous and differentiable at S(Pc_max) (hence, its value at S = 0 will be much greater than Pc_max"}>>> = quadratic
    high_extension_type<<<{"description": "Type of extension to use for the wetting curves when the liquid saturation is around 1.  The extensions modify the wetting capillary pressure for all saturation values greater than high_ratio * (1 - S_gr_del), where 1 - S_gr_del is the value of liquid saturation when the van Genuchten expression gives Pc = 0.  NONE: use the van Genuchten expression and when S > 1 - S_gr_del, set Pc = 0.  POWER: Pc is proportional to (1 - S)^power, where the coefficient of proportionality and the power are chosen so the resulting curve is continuous and differentiable"}>>> = power
    porepressure<<<{"description": "Variable that represents the porepressure of the single liquid phase"}>>> = pp
  []
[]

[Postprocessors<<<{"href": "../../syntax/Postprocessors/index.html"}>>>]
  [flux]
    type = FunctionValuePostprocessor<<<{"description": "Computes the value of a supplied function at a single point (scalable)", "href": "../postprocessors/FunctionValuePostprocessor.html"}>>>
    function<<<{"description": "The function which supplies the postprocessor value."}>>> = 'if(t <= 9, -10, 10)'
  []
  [hys_order]
    type = PointValue<<<{"description": "Compute the value of a variable at a specified location", "href": "../postprocessors/PointValue.html"}>>>
    point<<<{"description": "The physical point where the solution will be evaluated."}>>> = '0 0 0'
    variable<<<{"description": "The name of the variable that this postprocessor operates on."}>>> = hys_order
  []
  [sat]
    type = PointValue<<<{"description": "Compute the value of a variable at a specified location", "href": "../postprocessors/PointValue.html"}>>>
    point<<<{"description": "The physical point where the solution will be evaluated."}>>> = '0 0 0'
    variable<<<{"description": "The name of the variable that this postprocessor operates on."}>>> = sat
  []
  [pp]
    type = PointValue<<<{"description": "Compute the value of a variable at a specified location", "href": "../postprocessors/PointValue.html"}>>>
    point<<<{"description": "The physical point where the solution will be evaluated."}>>> = '0 0 0'
    variable<<<{"description": "The name of the variable that this postprocessor operates on."}>>> = pp
  []
[]

[Executioner<<<{"href": "../../syntax/Executioner/index.html"}>>>]
  type = Transient
  solve_type = Newton
  dt = 0.5
  end_time = 19
  nl_abs_tol = 1E-10
[]

[Outputs<<<{"href": "../../syntax/Outputs/index.html"}>>>]
  csv<<<{"description": "Output the scalar variable and postprocessors to a *.csv file using the default CSV output."}>>> = true
[]

Note the following features:

• The PorousFlowHysteresisOrder Material appears alongside the PorousFlow1PhaseHysP Material:

  [hys_order_material]
    type = PorousFlowHysteresisOrder
  []
  [pc_calculator]
    type = PorousFlow1PhaseHysP
    alpha_d = 10.0
    alpha_w = 7.0
    n_d = 1.5
    n_w = 1.9
    S_l_min = 0.1
    S_lr = 0.2
    S_gr_max = 0.3
    Pc_max = 12.0
    high_ratio = 0.9
    low_extension_type = quadratic
    high_extension_type = power
    porepressure = pp
  []
[]

• Saturation and hysteresis order may be recorded into AuxVariables using a PorousFlowPropertyAux

[AuxKernels<<<{"href": "../../syntax/AuxKernels/index.html"}>>>]
  [sat]
    type = PorousFlowPropertyAux<<<{"description": "AuxKernel to provide access to properties evaluated at quadpoints. Note that elemental AuxVariables must be used, so that these properties are integrated over each element.", "href": "../auxkernels/PorousFlowPropertyAux.html"}>>>
    variable<<<{"description": "The name of the variable that this object applies to"}>>> = sat
    property<<<{"description": "The fluid property that this auxillary kernel is to calculate"}>>> = saturation
  []
  [hys_order]
    type = PorousFlowPropertyAux<<<{"description": "AuxKernel to provide access to properties evaluated at quadpoints. Note that elemental AuxVariables must be used, so that these properties are integrated over each element.", "href": "../auxkernels/PorousFlowPropertyAux.html"}>>>
    variable<<<{"description": "The name of the variable that this object applies to"}>>> = hys_order
    property<<<{"description": "The fluid property that this auxillary kernel is to calculate"}>>> = hysteresis_order
  []
[]