PorousFlowHystereticInfo

This Material computes quantities relevant to hysteretic capillary pressure. The quantity computed is written into the PorousFlow_hysteretic_info_qp or PorousFlow_hysteretic_nodal Property (depending on the at_nodes parameter).

warning

This Material does not compute porepressures or saturations, so cannot be used in usual PorousFlow simulations. Instead, it enables preliminary exploration of hysteresis in models before hysteretic PorousFlow simulations are run: informative plots may be generated using this Material.

The quantity computed depends on the info_required parameter, which may be one of the following.

pc capillary pressure, given the saturation
sat saturation. The saturation is provided as an input, the capillary pressure is computed, and then the saturation is computed using the inverse relationship. This is useful when exploring the non-invertibility of the hysteretic relationships
sat_given_pc saturation, given the capillary pressure
various derivative information, such as the numerical error in $var element = document.getElementById("moose-equation-798cf36f-71a4-4d79-a60a-272ee4af7ffa");katex.render("dP_{c}/dS", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ , etc.

An example input file containing PorousFlowHystereticInfo is

# The saturation is varied with time and the capillary pressure is computed
[Mesh]
  [mesh]
    type = GeneratedMeshGenerator
    dim = 1
  []
[]

[GlobalParams]
  PorousFlowDictator = dictator
[]

[UserObjects]
  [dictator]
    type = PorousFlowDictator
    number_fluid_phases = 1
    number_fluid_components = 1
    porous_flow_vars = ''
  []
[]

[Variables]
  [dummy]
  []
[]

[Kernels]
  [dummy]
    type = TimeDerivative
    variable = dummy
  []
[]

[AuxVariables]
  [sat]
    initial_condition = 1
  []
  [hys_order]
    family = MONOMIAL
    order = CONSTANT
  []
  [pc]
    family = MONOMIAL
    order = CONSTANT
  []
[]

[AuxKernels]
  [sat_aux]
    type = FunctionAux
    variable = sat
    function = '1 - t'
  []
  [hys_order]
    type = PorousFlowPropertyAux
    variable = hys_order
    property = hysteresis_order
  []
  [pc]
    type = PorousFlowPropertyAux
    variable = pc
    property = hysteretic_info
  []
[]

[Materials]
  [hys_order]
    type = PorousFlowHysteresisOrder
  []
  [pc_calculator]
    type = PorousFlowHystereticInfo
    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
    sat_var = sat
  []
[]

[Postprocessors]
  [hys_order]
    type = PointValue
    point = '0 0 0'
    variable = hys_order
  []
  [sat]
    type = PointValue
    point = '0 0 0'
    variable = sat
  []
  [pc]
    type = PointValue
    point = '0 0 0'
    variable = pc
  []
[]

[Executioner]
  type = Transient
  solve_type = Linear
  dt = 0.1
  end_time = 1
[]

[Outputs]
  csv = true
[]

In this input file, saturation is an AuxVariable which is varied with time through a FunctionAux. The PorousFlowHystereticInfo Material computes capillary pressure, which is written into an AuxVariable using the PorousFlowPropertyAux AuxKernel. By choosing different FunctionAux, various capillary-pressure curves may be generated such as the one shown below.

The results of two hysteretic simulations employing different FunctionAux. The lines show the expected result while the crosses and asterisks show the MOOSE result.

Input Parameters

PorousFlowDictatorThe UserObject that holds the list of PorousFlow variable names
C++ Type:UserObjectName
Unit:(no unit assumed)
Controllable:No
Description:The UserObject that holds the list of PorousFlow variable names
S_gr_maxResidual 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
C++ Type:double
Unit:(no unit assumed)
Controllable:No
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
S_l_minMinimum 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
C++ Type:double
Unit:(no unit assumed)
Controllable:No
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
alpha_dvan Genuchten alpha parameter for the primary drying curve. If using standard units, this is measured in Pa^-1. Suggested value is around 1E-5
C++ Type:double
Unit:(no unit assumed)
Controllable:No
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
alpha_wvan Genuchten alpha parameter for the primary wetting curve. If using standard units, this is measured in Pa^-1. Suggested value is around 1E-5
C++ Type:double
Unit:(no unit assumed)
Controllable:No
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
n_dvan Genuchten n parameter for the primary drying curve. Dimensionless. Suggested value is around 2
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:van Genuchten n parameter for the primary drying curve. Dimensionless. Suggested value is around 2
n_wvan Genuchten n parameter for the primary wetting curve. Dimensionless. Suggested value is around 2
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:van Genuchten n parameter for the primary wetting curve. Dimensionless. Suggested value is around 2
sat_varVariable that represent liquid saturation. This is always needed to ensure the hysteretic order is computed correctly
C++ Type:std::vector<VariableName>
Unit:(no unit assumed)
Controllable:No
Description:Variable that represent liquid saturation. This is always needed to ensure the hysteretic order is computed correctly

Required Parameters

Pc_max1.79769e+308Value 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
Default:1.79769e+308
C++ Type:double
Unit:(no unit assumed)
Controllable:No
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
S_lr0Liquid 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
Default:0
C++ Type:double
Unit:(no unit assumed)
Controllable:No
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
at_nodesFalseEvaluate Material properties at nodes instead of quadpoints
Default:False
C++ Type:bool
Unit:(no unit assumed)
Controllable:No
Description:Evaluate Material properties at nodes instead of quadpoints
blockThe list of blocks (ids or names) that this object will be applied
C++ Type:std::vector<SubdomainName>
Unit:(no unit assumed)
Controllable:No
Description:The list of blocks (ids or names) that this object will be applied
boundaryThe list of boundaries (ids or names) from the mesh where this object applies
C++ Type:std::vector<BoundaryName>
Unit:(no unit assumed)
Controllable:No
Description:The list of boundaries (ids or names) from the mesh where this object applies
computeTrueWhen false, MOOSE will not call compute methods on this material. The user must call computeProperties() after retrieving the MaterialBase via MaterialBasePropertyInterface::getMaterialBase(). Non-computed MaterialBases are not sorted for dependencies.
Default:True
C++ Type:bool
Unit:(no unit assumed)
Controllable:No
Description:When false, MOOSE will not call compute methods on this material. The user must call computeProperties() after retrieving the MaterialBase via MaterialBasePropertyInterface::getMaterialBase(). Non-computed MaterialBases are not sorted for dependencies.
constant_onNONEWhen ELEMENT, MOOSE will only call computeQpProperties() for the 0th quadrature point, and then copy that value to the other qps.When SUBDOMAIN, MOOSE will only call computeQpProperties() for the 0th quadrature point, and then copy that value to the other qps. Evaluations on element qps will be skipped
Default:NONE
C++ Type:MooseEnum
Unit:(no unit assumed)
Options:NONE, ELEMENT, SUBDOMAIN
Controllable:No
Description:When ELEMENT, MOOSE will only call computeQpProperties() for the 0th quadrature point, and then copy that value to the other qps.When SUBDOMAIN, MOOSE will only call computeQpProperties() for the 0th quadrature point, and then copy that value to the other qps. Evaluations on element qps will be skipped
declare_suffixAn optional suffix parameter that can be appended to any declared properties. The suffix will be prepended with a '_' character.
C++ Type:MaterialPropertyName
Unit:(no unit assumed)
Controllable:No
Description:An optional suffix parameter that can be appended to any declared properties. The suffix will be prepended with a '_' character.
fd_eps1e-08Small quantity used in computing the finite-difference approximations to derivatives
Default:1e-08
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:Small quantity used in computing the finite-difference approximations to derivatives
high_extension_typepowerType 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
Default:power
C++ Type:MooseEnum
Unit:(no unit assumed)
Options:none, power
Controllable:No
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
high_ratio0.9The 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)
Default:0.9
C++ Type:double
Unit:(no unit assumed)
Controllable:No
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)
info_requiredpcThe type of information required. pc: capillary pressure given the saturation (pc_var is not used in this case). sat: given the liquid saturation, compute the capillary pressure, then invert the relationship to yield liquid saturation again (pc_var is not used in this case). This is useful to understand the non-invertibility of the hysteretic relationships. sat_given_pc: given the capillary pressure, compute the saturation. dS_dPc_err: relative error in d(sat)/d(pc) calculation, ie (S(Pc + fd_eps) - S(Pc - fd_eps))/(2 * eps * S'(Pc)) - 1, where S' is the coded derivative (Pc is computed from sat in this case: pc_var is not used). This is useful for checking derivatives. dPc_dS_err: relative error in d(pc)/d(sat) calculation, ie (Pc(S + fd_eps) - Pc(S - fd_eps)) / (2 * eps * Pc'(S)) - 1, where Pc' is the coded derative (pc_var is not used in this case). d2S_dPc2_err: relative error in d^2(sat)/d(pc)^2 (Pc is computed from sat in this case: pc_var is not used). d2Pc_dS2_err: relative error in d^2(pc)/d(sat)^2.
Default:pc
C++ Type:MooseEnum
Unit:(no unit assumed)
Options:pc, sat, sat_given_pc, dS_dPc_err, dPc_dS_err, d2S_dPc2_err, d2Pc_dS2_err
Controllable:No
Description:The type of information required. pc: capillary pressure given the saturation (pc_var is not used in this case). sat: given the liquid saturation, compute the capillary pressure, then invert the relationship to yield liquid saturation again (pc_var is not used in this case). This is useful to understand the non-invertibility of the hysteretic relationships. sat_given_pc: given the capillary pressure, compute the saturation. dS_dPc_err: relative error in d(sat)/d(pc) calculation, ie (S(Pc + fd_eps) - S(Pc - fd_eps))/(2 * eps * S'(Pc)) - 1, where S' is the coded derivative (Pc is computed from sat in this case: pc_var is not used). This is useful for checking derivatives. dPc_dS_err: relative error in d(pc)/d(sat) calculation, ie (Pc(S + fd_eps) - Pc(S - fd_eps)) / (2 * eps * Pc'(S)) - 1, where Pc' is the coded derative (pc_var is not used in this case). d2S_dPc2_err: relative error in d^2(sat)/d(pc)^2 (Pc is computed from sat in this case: pc_var is not used). d2Pc_dS2_err: relative error in d^2(pc)/d(sat)^2.
low_extension_typeexponentialType 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
Default:exponential
C++ Type:MooseEnum
Unit:(no unit assumed)
Options:none, quadratic, exponential
Controllable:No
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
pc_varVariable that represents capillary pressure. Depending on info_required, this may not be used to compute the info
C++ Type:std::vector<VariableName>
Unit:(no unit assumed)
Controllable:No
Description:Variable that represents capillary pressure. Depending on info_required, this may not be used to compute the info
prop_getter_suffixAn optional suffix parameter that can be appended to any attempt to retrieve/get material properties. The suffix will be prepended with a '_' character.
C++ Type:MaterialPropertyName
Unit:(no unit assumed)
Controllable:No
Description:An optional suffix parameter that can be appended to any attempt to retrieve/get material properties. The suffix will be prepended with a '_' character.
use_interpolated_stateFalseFor the old and older state use projected material properties interpolated at the quadrature points. To set up projection use the ProjectedStatefulMaterialStorageAction.
Default:False
C++ Type:bool
Unit:(no unit assumed)
Controllable:No
Description:For the old and older state use projected material properties interpolated at the quadrature points. To set up projection use the ProjectedStatefulMaterialStorageAction.

Optional Parameters

control_tagsAdds user-defined labels for accessing object parameters via control logic.
C++ Type:std::vector<std::string>
Unit:(no unit assumed)
Controllable:No
Description:Adds user-defined labels for accessing object parameters via control logic.
enableTrueSet the enabled status of the MooseObject.
Default:True
C++ Type:bool
Unit:(no unit assumed)
Controllable:Yes
Description:Set the enabled status of the MooseObject.
implicitTrueDetermines whether this object is calculated using an implicit or explicit form
Default:True
C++ Type:bool
Unit:(no unit assumed)
Controllable:No
Description:Determines whether this object is calculated using an implicit or explicit form
seed0The seed for the master random number generator
Default:0
C++ Type:unsigned int
Unit:(no unit assumed)
Controllable:No
Description:The seed for the master random number generator
use_displaced_meshFalseWhether or not this object should use the displaced mesh for computation. Note that in the case this is true but no displacements are provided in the Mesh block the undisplaced mesh will still be used.
Default:False
C++ Type:bool
Unit:(no unit assumed)
Controllable:No
Description:Whether or not this object should use the displaced mesh for computation. Note that in the case this is true but no displacements are provided in the Mesh block the undisplaced mesh will still be used.

Advanced Parameters

output_propertiesList of material properties, from this material, to output (outputs must also be defined to an output type)
C++ Type:std::vector<std::string>
Unit:(no unit assumed)
Controllable:No
Description:List of material properties, from this material, to output (outputs must also be defined to an output type)
outputsnone Vector of output names where you would like to restrict the output of variables(s) associated with this object
Default:none
C++ Type:std::vector<OutputName>
Unit:(no unit assumed)
Controllable:No
Description:Vector of output names where you would like to restrict the output of variables(s) associated with this object

Outputs Parameters

Input Files

(modules/porous_flow/test/tests/hysteresis/hys_pc_03.i)
(modules/porous_flow/test/tests/hysteresis/hys_pc_01.i)
(modules/porous_flow/test/tests/hysteresis/hys_pc_3.i)
(modules/porous_flow/test/tests/hysteresis/hys_pc_1.i)
(modules/porous_flow/test/tests/hysteresis/vary_sat_1.i)
(modules/porous_flow/test/tests/hysteresis/hys_pc_2.i)
(modules/porous_flow/test/tests/hysteresis/hys_pc_02.i)

(modules/porous_flow/test/tests/hysteresis/vary_sat_1.i)

# The saturation is varied with time and the capillary pressure is computed
[Mesh]
  [mesh]
    type = GeneratedMeshGenerator
    dim = 1
  []
[]

[GlobalParams]
  PorousFlowDictator = dictator
[]

[UserObjects]
  [dictator]
    type = PorousFlowDictator
    number_fluid_phases = 1
    number_fluid_components = 1
    porous_flow_vars = ''
  []
[]

[Variables]
  [dummy]
  []
[]

[Kernels]
  [dummy]
    type = TimeDerivative
    variable = dummy
  []
[]

[AuxVariables]
  [sat]
    initial_condition = 1
  []
  [hys_order]
    family = MONOMIAL
    order = CONSTANT
  []
  [pc]
    family = MONOMIAL
    order = CONSTANT
  []
[]

[AuxKernels]
  [sat_aux]
    type = FunctionAux
    variable = sat
    function = '1 - t'
  []
  [hys_order]
    type = PorousFlowPropertyAux
    variable = hys_order
    property = hysteresis_order
  []
  [pc]
    type = PorousFlowPropertyAux
    variable = pc
    property = hysteretic_info
  []
[]

[Materials]
  [hys_order]
    type = PorousFlowHysteresisOrder
  []
  [pc_calculator]
    type = PorousFlowHystereticInfo
    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
    sat_var = sat
  []
[]

[Postprocessors]
  [hys_order]
    type = PointValue
    point = '0 0 0'
    variable = hys_order
  []
  [sat]
    type = PointValue
    point = '0 0 0'
    variable = sat
  []
  [pc]
    type = PointValue
    point = '0 0 0'
    variable = pc
  []
[]

[Executioner]
  type = Transient
  solve_type = Linear
  dt = 0.1
  end_time = 1
[]

[Outputs]
  csv = true
[]

(modules/porous_flow/test/tests/hysteresis/hys_pc_03.i)

# Capillary-pressure calculation.  Primary drying curve with low_extension_type = exponential
# When comparing the results with a by-hand computation, remember the MOOSE results are averaged over an element
[Mesh]
  [mesh]
    type = GeneratedMeshGenerator
    dim = 1
    xmin = 0
    xmax = 1
    nx = 100
  []
[]

[GlobalParams]
  PorousFlowDictator = dictator
[]

[UserObjects]
  [dictator]
    type = PorousFlowDictator
    number_fluid_phases = 1
    number_fluid_components = 1
    porous_flow_vars = ''
  []
[]

[Variables]
  [sat]
  []
[]

[ICs]
  [sat]
    type = FunctionIC
    variable = sat
    function = 'x'
  []
[]

[BCs]
  [sat]
    type = FunctionDirichletBC
    variable = sat
    function = 'x'
    boundary = 'left right'
  []
[]


[Kernels]
  [dummy]
    type = Diffusion
    variable = sat
  []
[]

[Materials]
  [hys_order]
    type = PorousFlowHysteresisOrder
  []
  [pc_calculator]
    type = PorousFlowHystereticInfo
    alpha_d = 10.0
    alpha_w = 10.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
    low_extension_type = exponential
    sat_var = sat
  []
[]

[AuxVariables]
  [hys_order]
    family = MONOMIAL
    order = CONSTANT
  []
  [pc]
    family = MONOMIAL
    order = CONSTANT
  []
[]

[AuxKernels]
  [hys_order]
    type = PorousFlowPropertyAux
    variable = hys_order
    property = hysteresis_order
  []
  [pc]
    type = PorousFlowPropertyAux
    variable = pc
    property = hysteretic_info
  []
[]

[VectorPostprocessors]
  [pc]
    type = LineValueSampler
    warn_discontinuous_face_values = false
    start_point = '0 0 0'
    end_point = '1 0 0'
    num_points = 10
    sort_by = x
    variable = 'sat pc'
  []
[]

[Executioner]
  type = Transient
  solve_type = Linear
  dt = 1
  end_time = 1
[]

[Outputs]
  csv = true
[]

(modules/porous_flow/test/tests/hysteresis/hys_pc_01.i)

# Capillary-pressure calculation.  Primary drying curve with low_extension_type = none
# When comparing the results with a by-hand computation, remember the MOOSE results are averaged over an element
[Mesh]
  [mesh]
    type = GeneratedMeshGenerator
    dim = 1
    xmin = 0
    xmax = 1
    nx = 100
  []
[]

[GlobalParams]
  PorousFlowDictator = dictator
[]

[UserObjects]
  [dictator]
    type = PorousFlowDictator
    number_fluid_phases = 1
    number_fluid_components = 1
    porous_flow_vars = ''
  []
[]

[Variables]
  [sat]
  []
[]

[ICs]
  [sat]
    type = FunctionIC
    variable = sat
    function = 'x'
  []
[]

[BCs]
  [sat]
    type = FunctionDirichletBC
    variable = sat
    function = 'x'
    boundary = 'left right'
  []
[]


[Kernels]
  [dummy]
    type = Diffusion
    variable = sat
  []
[]

[Materials]
  [hys_order]
    type = PorousFlowHysteresisOrder
  []
  [pc_calculator]
    type = PorousFlowHystereticInfo
    alpha_d = 10.0
    alpha_w = 10.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
    low_extension_type = none
    sat_var = sat
  []
[]

[AuxVariables]
  [hys_order]
    family = MONOMIAL
    order = CONSTANT
  []
  [pc]
    family = MONOMIAL
    order = CONSTANT
  []
[]

[AuxKernels]
  [hys_order]
    type = PorousFlowPropertyAux
    variable = hys_order
    property = hysteresis_order
  []
  [pc]
    type = PorousFlowPropertyAux
    variable = pc
    property = hysteretic_info
  []
[]

[VectorPostprocessors]
  [pc]
    type = LineValueSampler
    warn_discontinuous_face_values = false
    start_point = '0 0 0'
    end_point = '1 0 0'
    num_points = 10
    sort_by = x
    variable = 'sat pc'
  []
[]

[Executioner]
  type = Transient
  solve_type = Linear
  dt = 1
  end_time = 1
[]

[Outputs]
  csv = true
[]

(modules/porous_flow/test/tests/hysteresis/hys_pc_3.i)

# Capillary-pressure calculation.  Third-order curve
# When comparing the results with a by-hand computation, remember the MOOSE results are averaged over an element
[Mesh]
  [mesh]
    type = GeneratedMeshGenerator
    dim = 1
    xmin = 0.4
    xmax = 0.9
    nx = 50
  []
[]

[GlobalParams]
  PorousFlowDictator = dictator
[]

[UserObjects]
  [dictator]
    type = PorousFlowDictator
    number_fluid_phases = 1
    number_fluid_components = 1
    porous_flow_vars = ''
  []
[]

[Variables]
  [sat]
  []
[]

[ICs]
  [sat]
    type = FunctionIC
    variable = sat
    function = 'x'
  []
[]

[BCs]
  [sat]
    type = FunctionDirichletBC
    variable = sat
    function = 'x'
    boundary = 'left right'
  []
[]


[Kernels]
  [dummy]
    type = Diffusion
    variable = sat
  []
[]

[Materials]
  [hys_order]
    type = PorousFlowHysteresisOrder
    initial_order = 3
    previous_turning_points = '0.1 0.9 0.4'
  []
  [pc_calculator]
    type = PorousFlowHystereticInfo
    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 = none
    high_extension_type = none
    sat_var = sat
  []
[]

[AuxVariables]
  [hys_order]
    family = MONOMIAL
    order = CONSTANT
  []
  [pc]
    family = MONOMIAL
    order = CONSTANT
  []
[]

[AuxKernels]
  [hys_order]
    type = PorousFlowPropertyAux
    variable = hys_order
    property = hysteresis_order
  []
  [pc]
    type = PorousFlowPropertyAux
    variable = pc
    property = hysteretic_info
  []
[]

[VectorPostprocessors]
  [pc]
    type = LineValueSampler
    warn_discontinuous_face_values = false
    start_point = '0.4 0 0'
    end_point = '0.9 0 0'
    num_points = 8
    sort_by = x
    variable = 'sat pc'
  []
[]

[Executioner]
  type = Transient
  solve_type = Linear
  dt = 1
  end_time = 1
[]

[Outputs]
  csv = true
[]

(modules/porous_flow/test/tests/hysteresis/hys_pc_1.i)

# Capillary-pressure calculation.  First-order wetting curve
# When comparing the results with a by-hand computation, remember the MOOSE results are averaged over an element
# Also, when using info_required=sat, remember that: (1) the hysteretic capillary pressure is not invertible if no high extension is used; (2) if saturation exceeds the turning point (eg sat <= 0.1) then the drying curve will be used
[Mesh]
  [mesh]
    type = GeneratedMeshGenerator
    dim = 1
    xmin = 0
    xmax = 1
    nx = 100
  []
[]

[GlobalParams]
  PorousFlowDictator = dictator
[]

[UserObjects]
  [dictator]
    type = PorousFlowDictator
    number_fluid_phases = 1
    number_fluid_components = 1
    porous_flow_vars = ''
  []
[]

[Variables]
  [sat]
  []
[]

[ICs]
  [sat]
    type = FunctionIC
    variable = sat
    function = 'x'
  []
[]

[BCs]
  [sat]
    type = FunctionDirichletBC
    variable = sat
    function = 'x'
    boundary = 'left right'
  []
[]


[Kernels]
  [dummy]
    type = Diffusion
    variable = sat
  []
[]

[Materials]
  [hys_order]
    type = PorousFlowHysteresisOrder
    initial_order = 1
    previous_turning_points = 0.1
  []
  [pc_calculator]
    type = PorousFlowHystereticInfo
    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 = none
    high_extension_type = none
    sat_var = sat
  []
[]

[AuxVariables]
  [hys_order]
    family = MONOMIAL
    order = CONSTANT
  []
  [pc]
    family = MONOMIAL
    order = CONSTANT
  []
[]

[AuxKernels]
  [hys_order]
    type = PorousFlowPropertyAux
    variable = hys_order
    property = hysteresis_order
  []
  [pc]
    type = PorousFlowPropertyAux
    variable = pc
    property = hysteretic_info
  []
[]

[VectorPostprocessors]
  [pc]
    type = LineValueSampler
    warn_discontinuous_face_values = false
    start_point = '0 0 0'
    end_point = '1 0 0'
    num_points = 10
    sort_by = x
    variable = 'sat pc'
  []
[]

[Executioner]
  type = Transient
  solve_type = Linear
  dt = 1
  end_time = 1
[]

[Outputs]
  csv = true
[]

(modules/porous_flow/test/tests/hysteresis/vary_sat_1.i)

# The saturation is varied with time and the capillary pressure is computed
[Mesh]
  [mesh]
    type = GeneratedMeshGenerator
    dim = 1
  []
[]

[GlobalParams]
  PorousFlowDictator = dictator
[]

[UserObjects]
  [dictator]
    type = PorousFlowDictator
    number_fluid_phases = 1
    number_fluid_components = 1
    porous_flow_vars = ''
  []
[]

[Variables]
  [dummy]
  []
[]

[Kernels]
  [dummy]
    type = TimeDerivative
    variable = dummy
  []
[]

[AuxVariables]
  [sat]
    initial_condition = 1
  []
  [hys_order]
    family = MONOMIAL
    order = CONSTANT
  []
  [pc]
    family = MONOMIAL
    order = CONSTANT
  []
[]

[AuxKernels]
  [sat_aux]
    type = FunctionAux
    variable = sat
    function = '1 - t'
  []
  [hys_order]
    type = PorousFlowPropertyAux
    variable = hys_order
    property = hysteresis_order
  []
  [pc]
    type = PorousFlowPropertyAux
    variable = pc
    property = hysteretic_info
  []
[]

[Materials]
  [hys_order]
    type = PorousFlowHysteresisOrder
  []
  [pc_calculator]
    type = PorousFlowHystereticInfo
    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
    sat_var = sat
  []
[]

[Postprocessors]
  [hys_order]
    type = PointValue
    point = '0 0 0'
    variable = hys_order
  []
  [sat]
    type = PointValue
    point = '0 0 0'
    variable = sat
  []
  [pc]
    type = PointValue
    point = '0 0 0'
    variable = pc
  []
[]

[Executioner]
  type = Transient
  solve_type = Linear
  dt = 0.1
  end_time = 1
[]

[Outputs]
  csv = true
[]

(modules/porous_flow/test/tests/hysteresis/hys_pc_2.i)

# Capillary-pressure calculation.  Second-order drying curve
# When comparing the results with a by-hand computation, remember the MOOSE results are averaged over an element
[Mesh]
  [mesh]
    type = GeneratedMeshGenerator
    dim = 1
    xmin = 0.1
    xmax = 0.9
    nx = 80
  []
[]

[GlobalParams]
  PorousFlowDictator = dictator
[]

[UserObjects]
  [dictator]
    type = PorousFlowDictator
    number_fluid_phases = 1
    number_fluid_components = 1
    porous_flow_vars = ''
  []
[]

[Variables]
  [sat]
  []
[]

[ICs]
  [sat]
    type = FunctionIC
    variable = sat
    function = 'x'
  []
[]

[BCs]
  [sat]
    type = FunctionDirichletBC
    variable = sat
    function = 'x'
    boundary = 'left right'
  []
[]


[Kernels]
  [dummy]
    type = Diffusion
    variable = sat
  []
[]

[Materials]
  [hys_order]
    type = PorousFlowHysteresisOrder
    initial_order = 2
    previous_turning_points = '0.1 0.9'
  []
  [pc_calculator]
    type = PorousFlowHystereticInfo
    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 = none
    high_extension_type = none
    sat_var = sat
  []
[]

[AuxVariables]
  [hys_order]
    family = MONOMIAL
    order = CONSTANT
  []
  [pc]
    family = MONOMIAL
    order = CONSTANT
  []
[]

[AuxKernels]
  [hys_order]
    type = PorousFlowPropertyAux
    variable = hys_order
    property = hysteresis_order
  []
  [pc]
    type = PorousFlowPropertyAux
    variable = pc
    property = hysteretic_info
  []
[]

[VectorPostprocessors]
  [pc]
    type = LineValueSampler
    warn_discontinuous_face_values = false
    start_point = '0.1 0 0'
    end_point = '0.9 0 0'
    num_points = 8
    sort_by = x
    variable = 'sat pc'
  []
[]

[Executioner]
  type = Transient
  solve_type = Linear
  dt = 1
  end_time = 1
[]

[Outputs]
  csv = true
[]

(modules/porous_flow/test/tests/hysteresis/hys_pc_02.i)

# Capillary-pressure calculation.  Primary drying curve with low_extension_type = quadratic
# When comparing the results with a by-hand computation, remember the MOOSE results are averaged over an element
[Mesh]
  [mesh]
    type = GeneratedMeshGenerator
    dim = 1
    xmin = 0
    xmax = 1
    nx = 100
  []
[]

[GlobalParams]
  PorousFlowDictator = dictator
[]

[UserObjects]
  [dictator]
    type = PorousFlowDictator
    number_fluid_phases = 1
    number_fluid_components = 1
    porous_flow_vars = ''
  []
[]

[Variables]
  [sat]
  []
[]

[ICs]
  [sat]
    type = FunctionIC
    variable = sat
    function = 'x'
  []
[]

[BCs]
  [sat]
    type = FunctionDirichletBC
    variable = sat
    function = 'x'
    boundary = 'left right'
  []
[]


[Kernels]
  [dummy]
    type = Diffusion
    variable = sat
  []
[]

[Materials]
  [hys_order]
    type = PorousFlowHysteresisOrder
  []
  [pc_calculator]
    type = PorousFlowHystereticInfo
    alpha_d = 10.0
    alpha_w = 10.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
    low_extension_type = quadratic
    sat_var = sat
  []
[]

[AuxVariables]
  [hys_order]
    family = MONOMIAL
    order = CONSTANT
  []
  [pc]
    family = MONOMIAL
    order = CONSTANT
  []
[]

[AuxKernels]
  [hys_order]
    type = PorousFlowPropertyAux
    variable = hys_order
    property = hysteresis_order
  []
  [pc]
    type = PorousFlowPropertyAux
    variable = pc
    property = hysteretic_info
  []
[]

[VectorPostprocessors]
  [pc]
    type = LineValueSampler
    warn_discontinuous_face_values = false
    start_point = '0 0 0'
    end_point = '1 0 0'
    num_points = 10
    sort_by = x
    variable = 'sat pc'
  []
[]

[Executioner]
  type = Transient
  solve_type = Linear
  dt = 1
  end_time = 1
[]

[Outputs]
  csv = true
[]

Input Parameters
Input Files