PorousFlowHalfCubicSink

Applies a flux sink to a boundary. The base flux defined by PorousFlowSink is multiplied by a cubic.

The basic sink $var element = document.getElementById("moose-equation-006b40b7-d010-4a18-b873-588e9fda490d");katex.render("f(x,t)", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ is scaled by a cubic flux multiplier of the pressure of a fluid phase or the temperature. Label the pressure or temperature by $var element = document.getElementById("moose-equation-eb1b349f-1f12-4ea0-be07-f4c21faa419f");katex.render("v", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ . Then the sink strength is $var element = document.getElementById("moose-equation-add7bc19-3e74-4bb0-80e0-34c72c9bb60e");katex.render(" s = f(t, x) \\times g(v) \\ .", element, {displayMode:true,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$

The function $var element = document.getElementById("moose-equation-7b42cb26-54db-4c15-ba5e-7ac4adb805ce");katex.render("g", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ interpolates smoothly between $var element = document.getElementById("moose-equation-80beb161-9ac3-4aae-a6eb-94db457f96e2");katex.render("g(v_{\\mathrm{cutoff}} + v_{\\mathrm{center}}) = 0", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ and $var element = document.getElementById("moose-equation-e2193d09-e9eb-47ec-9097-cf3ecdba1837");katex.render("g(v_{\\mathrm{center}}) = g_{\\mathrm{max}}", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ . Since a cubic interpolation is used, the function $var element = document.getElementById("moose-equation-56bc850e-6304-4b44-b8d7-0382aea64015");katex.render("g", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ is continuous and its derivative is continuous everywhere. Explicitly: $var element = document.getElementById("moose-equation-9ef6c035-e389-4e0e-b94b-ef71133455c7");katex.render(" g(v) = \\begin{cases} g_{\\mathrm{max}} & \\text{if}\\ v > v_{\\mathrm{center}} \\\\ \\frac{g_{\\mathrm{max}} \\left(2[v - v_{\\mathrm{center}}] + v_{\\mathrm{cutoff}}\\right) (v - v_{\\mathrm{center}} - v_{\\mathrm{cutoff}})^2}{v_{\\mathrm{cutoff}}^3}, & v_{\\mathrm{cutoff}} + v_{\\mathrm{center}} \\leq v \\leq v_{\\mathrm{center}} \\\\ 0 & \\text{if}\\ v < v_{\\mathrm{cutoff}} + v_{\\mathrm{center}} \\\\ \\end{cases}", element, {displayMode:true,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ Here the units of $var element = document.getElementById("moose-equation-64e7e62a-15b7-4d22-a249-830298a827af");katex.render("f\\times g", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ are kg.m $var element = document.getElementById("moose-equation-a4f6bdbe-2555-4f61-b8a8-f03bacc63bee");katex.render("^{-2}", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ .s $var element = document.getElementById("moose-equation-8177bb24-ce9f-4da2-aad0-61104bcad871");katex.render("^{-1}", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ (for fluids) or J.m $var element = document.getElementById("moose-equation-ed8c3f47-eca5-414e-947d-e41065852a4a");katex.render("^{-1}", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ .s $var element = document.getElementById("moose-equation-87a547e5-b2ff-4543-b376-cfc396a5a423");katex.render("^{-1}", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ (for heat). The parameters $var element = document.getElementById("moose-equation-c8fd7267-5e32-4fe5-b867-8fd026b626c2");katex.render("g_{\\mathrm{max}}", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ , $var element = document.getElementById("moose-equation-834cc2ff-335e-4a57-82a0-c6f854fdf07a");katex.render("v_{\\mathrm{center}}", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ and $var element = document.getElementById("moose-equation-8ba67a35-ed41-42b7-85ed-62abc8f6f08e");katex.render("v_{\\mathrm{cutoff}}", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ are given in the input file using the max, center and cutoff options, respectively.

This sink is often used in groundwater modelling of evapotranspiration, where:

$var element = document.getElementById("moose-equation-e7168230-38db-4318-9860-c1691ac0e730");katex.render("v_{\\mathrm{center}} = 0", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ , since evapotranspiration is maximum when the groundwater table is at or above the topography;
$var element = document.getElementById("moose-equation-f935dc81-dfbb-488c-bbcd-065159cadce3");katex.render("v_{\\mathrm{cutoff}} = -10^{4}d", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ , where $var element = document.getElementById("moose-equation-fe03ab75-662d-45db-9980-fe15f04ebfd8");katex.render("d", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ is the root depth at which there is no evapotranspiration from the groundwater system;
$var element = document.getElementById("moose-equation-9ab2ff48-7048-448c-85ac-ef7507e36b1d");katex.render("g_{\\mathrm{max}}", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ is the value of the pan evaporation, in kg.m $var element = document.getElementById("moose-equation-0c69fc90-7c31-4e31-ae56-455cf31a34a3");katex.render("^{-2}", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ .s $var element = document.getElementById("moose-equation-4ef49c33-b34f-41ec-a488-f17890749001");katex.render("^{-1}", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ .

If $var element = document.getElementById("moose-equation-2866ae32-c6f3-4996-a45e-baa591ceca72");katex.render("f>0", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ then the boundary condition will act as a sink, while if $var element = document.getElementById("moose-equation-b6d5b297-fc08-4467-8592-f22a4d02f8ba");katex.render("f<0", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ the boundary condition acts as a source. If applied to a fluid-component equation, the function $var element = document.getElementById("moose-equation-86f3ac55-473c-46b6-9fa0-ceb68f61ea71");katex.render("f", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ has units kg.m $var element = document.getElementById("moose-equation-6e9db660-daaa-478a-8bf6-863709d2ca86");katex.render("^{-2}", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ .s $var element = document.getElementById("moose-equation-2113c655-0a23-443c-8657-146d0cf9b01c");katex.render("^{-1}", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ . If applied to the heat equation, the function $var element = document.getElementById("moose-equation-568fed9c-1737-4fd0-bd98-c1aff45f91a0");katex.render("f", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ has units J.m $var element = document.getElementById("moose-equation-32d42db6-73d3-4fe3-b6b5-f740291cfdff");katex.render("^{-2}", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ .s $var element = document.getElementById("moose-equation-f73288af-b735-45cb-a8ae-bffc673e827f");katex.render("^{-1}", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ . These units are potentially modified if the extra building blocks enumerated below are used.

In addition, the sink may be multiplied by any or all of the following quantities through the optional parameters list.

Fluid relative permeability
Fluid mobility ( $var element = document.getElementById("moose-equation-863e4f11-d8e5-471e-a237-bf8cfcc7f0da");katex.render("k_{ij}n_{i}n_{j}k_{r} \\rho / \\nu", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ , where $var element = document.getElementById("moose-equation-d574fb2f-d3c2-4cf2-a1cc-7823e26dc983");katex.render("n", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ is the normal vector to the boundary)
Fluid mass fraction
Fluid internal energy
Thermal conductivity

See boundary conditions for many more details and discussion.

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
boundaryThe list of boundary IDs from the mesh where this object applies
C++ Type:std::vector<BoundaryName>
Unit:(no unit assumed)
Controllable:No
Description:The list of boundary IDs from the mesh where this object applies
centerCenter of the cubic flux multiplier (measured in Pa (or K for temperature BCs)).
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:Center of the cubic flux multiplier (measured in Pa (or K for temperature BCs)).
cutoffCutoff of the cubic (measured in Pa (or K for temperature BCs)). This needs to be less than zero.
C++ Type:FunctionName
Unit:(no unit assumed)
Controllable:No
Description:Cutoff of the cubic (measured in Pa (or K for temperature BCs)). This needs to be less than zero.
maxMaximum of the cubic flux multiplier. Denote x = porepressure - center (or in the case of a heat flux with no fluid, x = temperature - center). Then Flux out is multiplied by (max/cutoff^3)*(2x + cutoff)(x - cutoff)^2 for cutoff < x < 0. Flux out is multiplied by max for x >= 0. Flux out is multiplied by 0 for x <= cutoff.
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:Maximum of the cubic flux multiplier. Denote x = porepressure - center (or in the case of a heat flux with no fluid, x = temperature - center). Then Flux out is multiplied by (max/cutoff^3)*(2x + cutoff)(x - cutoff)^2 for cutoff < x < 0. Flux out is multiplied by max for x >= 0. Flux out is multiplied by 0 for x <= cutoff.
variableThe name of the variable that this residual object operates on
C++ Type:NonlinearVariableName
Unit:(no unit assumed)
Controllable:No
Description:The name of the variable that this residual object operates on

Required Parameters

PT_shiftWhenever the sink is an explicit function of porepressure (such as a PiecewiseLinear function) the argument of the function is set to P - PT_shift instead of simply P. Similarly for temperature. PT_shift does not enter into any use_* calculations.
C++ Type:std::vector<VariableName>
Unit:(no unit assumed)
Controllable:No
Description:Whenever the sink is an explicit function of porepressure (such as a PiecewiseLinear function) the argument of the function is set to P - PT_shift instead of simply P. Similarly for temperature. PT_shift does not enter into any use_* calculations.
displacementsThe displacements
C++ Type:std::vector<VariableName>
Unit:(no unit assumed)
Controllable:No
Description:The displacements
fluid_phaseIf supplied, then this BC will potentially be a function of fluid pressure, and you can use mass_fraction_component, use_mobility, use_relperm, use_enthalpy and use_energy. If not supplied, then this BC can only be a function of temperature
C++ Type:unsigned int
Unit:(no unit assumed)
Controllable:No
Description:If supplied, then this BC will potentially be a function of fluid pressure, and you can use mass_fraction_component, use_mobility, use_relperm, use_enthalpy and use_energy. If not supplied, then this BC can only be a function of temperature
flux_function1The flux. The flux is OUT of the medium: hence positive values of this function means this BC will act as a SINK, while negative values indicate this flux will be a SOURCE. The functional form is useful for spatially or temporally varying sinks. Without any use_*, this function is measured in kg.m^-2.s^-1 (or J.m^-2.s^-1 for the case with only heat and no fluids)
Default:1
C++ Type:FunctionName
Unit:(no unit assumed)
Controllable:No
Description:The flux. The flux is OUT of the medium: hence positive values of this function means this BC will act as a SINK, while negative values indicate this flux will be a SOURCE. The functional form is useful for spatially or temporally varying sinks. Without any use_*, this function is measured in kg.m^-2.s^-1 (or J.m^-2.s^-1 for the case with only heat and no fluids)
mass_fraction_componentThe index corresponding to a fluid component. If supplied, the flux will be multiplied by the nodal mass fraction for the component
C++ Type:unsigned int
Unit:(no unit assumed)
Controllable:No
Description:The index corresponding to a fluid component. If supplied, the flux will be multiplied by the nodal mass fraction for the component
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_enthalpyFalseIf true, then fluxes are multiplied by enthalpy. In this case bare_flux is measured in kg.m^-2.s^-1 / (J.kg). This can be used in conjunction with other use_*
Default:False
C++ Type:bool
Unit:(no unit assumed)
Controllable:No
Description:If true, then fluxes are multiplied by enthalpy. In this case bare_flux is measured in kg.m^-2.s^-1 / (J.kg). This can be used in conjunction with other use_*
use_internal_energyFalseIf true, then fluxes are multiplied by fluid internal energy. In this case bare_flux is measured in kg.m^-2.s^-1 / (J.kg). This can be used in conjunction with other use_*
Default:False
C++ Type:bool
Unit:(no unit assumed)
Controllable:No
Description:If true, then fluxes are multiplied by fluid internal energy. In this case bare_flux is measured in kg.m^-2.s^-1 / (J.kg). This can be used in conjunction with other use_*
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.
use_mobilityFalseIf true, then fluxes are multiplied by (density*permeability_nn/viscosity), where the '_nn' indicates the component normal to the boundary. In this case bare_flux is measured in Pa.m^-1. This can be used in conjunction with other use_*
Default:False
C++ Type:bool
Unit:(no unit assumed)
Controllable:No
Description:If true, then fluxes are multiplied by (density*permeability_nn/viscosity), where the '_nn' indicates the component normal to the boundary. In this case bare_flux is measured in Pa.m^-1. This can be used in conjunction with other use_*
use_relpermFalseIf true, then fluxes are multiplied by relative permeability. This can be used in conjunction with other use_*
Default:False
C++ Type:bool
Unit:(no unit assumed)
Controllable:No
Description:If true, then fluxes are multiplied by relative permeability. This can be used in conjunction with other use_*
use_thermal_conductivityFalseIf true, then fluxes are multiplied by thermal conductivity projected onto the normal direction. This can be used in conjunction with other use_*
Default:False
C++ Type:bool
Unit:(no unit assumed)
Controllable:No
Description:If true, then fluxes are multiplied by thermal conductivity projected onto the normal direction. This can be used in conjunction with other use_*

Optional Parameters

absolute_value_vector_tagsThe tags for the vectors this residual object should fill with the absolute value of the residual contribution
C++ Type:std::vector<TagName>
Unit:(no unit assumed)
Controllable:No
Description:The tags for the vectors this residual object should fill with the absolute value of the residual contribution
extra_matrix_tagsThe extra tags for the matrices this Kernel should fill
C++ Type:std::vector<TagName>
Unit:(no unit assumed)
Controllable:No
Description:The extra tags for the matrices this Kernel should fill
extra_vector_tagsThe extra tags for the vectors this Kernel should fill
C++ Type:std::vector<TagName>
Unit:(no unit assumed)
Controllable:No
Description:The extra tags for the vectors this Kernel should fill
matrix_tagssystemThe tag for the matrices this Kernel should fill
Default:system
C++ Type:MultiMooseEnum
Unit:(no unit assumed)
Options:nontime, system
Controllable:No
Description:The tag for the matrices this Kernel should fill
vector_tagsnontimeThe tag for the vectors this Kernel should fill
Default:nontime
C++ Type:MultiMooseEnum
Unit:(no unit assumed)
Options:nontime, time
Controllable:No
Description:The tag for the vectors this Kernel should fill

Tagging 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.
diag_save_inThe name of auxiliary variables to save this BC's diagonal jacobian contributions to. Everything about that variable must match everything about this variable (the type, what blocks it's on, etc.)
C++ Type:std::vector<AuxVariableName>
Unit:(no unit assumed)
Controllable:No
Description:The name of auxiliary variables to save this BC's diagonal jacobian contributions to. Everything about that variable must match everything about this variable (the type, what blocks it's on, etc.)
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
save_inThe name of auxiliary variables to save this BC's residual contributions to. Everything about that variable must match everything about this variable (the type, what blocks it's on, etc.)
C++ Type:std::vector<AuxVariableName>
Unit:(no unit assumed)
Controllable:No
Description:The name of auxiliary variables to save this BC's residual contributions to. Everything about that variable must match everything about this variable (the type, what blocks it's on, etc.)
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

Input Files

(modules/porous_flow/examples/groundwater/ex02_steady_state.i)
(modules/porous_flow/test/tests/jacobian/hcs01.i)
(modules/porous_flow/test/tests/jacobian/hcs02.i)
(modules/porous_flow/test/tests/sinks/s06.i)
(modules/porous_flow/examples/groundwater/ex02_abstraction.i)

(modules/porous_flow/examples/groundwater/ex02_steady_state.i)

# Steady-state groundwater model.  See groundwater_models.md for a detailed description
[Mesh]
  [basic_mesh]
    # mesh create by external program: lies within -500<=x<=500 and -200<=y<=200, with varying z
    type = FileMeshGenerator
    file = ex02_mesh.e
  []
  [name_blocks]
    type = RenameBlockGenerator
    input = basic_mesh
    old_block = '2 3 4'
    new_block = 'bot_aquifer aquitard top_aquifer'
  []
  [zmax]
    type = SideSetsFromNormalsGenerator
    input = name_blocks
    normal_tol = 0.1
    new_boundary = zmax
    normals = '0 0 1'
  []
  [xmin_bot_aquifer]
    type = ParsedGenerateSideset
    input = zmax
    included_subdomains = 2
    normal = '-1 0 0'
    combinatorial_geometry = 'x <= -500.0'
    new_sideset_name = xmin_bot_aquifer
  []
[]

[GlobalParams]
  PorousFlowDictator = dictator
[]

[Variables]
  [pp]
  []
[]

[ICs]
  [pp]
    type = FunctionIC
    variable = pp
    function = initial_pp
  []
[]

[BCs]
  [rainfall_recharge]
    type = PorousFlowSink
    boundary = zmax
    variable = pp
    flux_function = -1E-6  # recharge of 0.1mm/day = 1E-4m3/m2/day = 0.1kg/m2/day ~ 1E-6kg/m2/s
  []
  [evapotranspiration]
    type = PorousFlowHalfCubicSink
    boundary = zmax
    variable = pp
    center = 0.0
    cutoff = -5E4 # roots of depth 5m.  5m of water = 5E4 Pa
    use_mobility = true
    fluid_phase = 0
    # Assume pan evaporation of 4mm/day = 4E-3m3/m2/day = 4kg/m2/day ~ 4E-5kg/m2/s
    # Assume that if permeability was 1E-10m^2 and water table at topography then ET acts as pan strength
    # Because use_mobility = true, then 4E-5 = maximum_flux = max * perm * density / visc = max * 1E-4, so max = 40
    max = 40
  []
[]

[DiracKernels]
  [river]
    type = PorousFlowPolyLineSink
    SumQuantityUO = baseflow
    point_file = ex02_river.bh
    # Assume a perennial river.
    # Assume the river has an incision depth of 1m and a stage height of 1.5m, and these are constant in time and uniform over the whole model.  Hence, if groundwater head is 0.5m (5000Pa) there will be no baseflow and leakage.
    p_or_t_vals = '-999995000 5000 1000005000'
    # Assume the riverbed conductance, k_zz*density*river_segment_length*river_width/riverbed_thickness/viscosity = 1E-6*river_segment_length kg/Pa/s
    fluxes = '-1E3 0 1E3'
    variable = pp
  []
[]

[Functions]
  [initial_pp]
    type = SolutionFunction
    scale_factor = 1E4
    from_variable = cosflow_depth
    solution = initial_mesh
  []
  [baseflow_rate]
    type = ParsedFunction
    symbol_names = 'baseflow_kg dt'
    symbol_values = 'baseflow_kg dt'
    expression = 'baseflow_kg / dt * 24.0 * 3600.0 / 400.0'
  []
[]

[PorousFlowUnsaturated]
  fp = simple_fluid
  porepressure = pp
[]

[FluidProperties]
  [simple_fluid]
    type = SimpleFluidProperties
  []
[]

[Materials]
  [porosity_everywhere]
    type = PorousFlowPorosityConst
    porosity = 0.05
  []
  [permeability_aquifers]
    type = PorousFlowPermeabilityConst
    block = 'top_aquifer bot_aquifer'
    permeability = '1E-12 0 0 0 1E-12 0 0 0 1E-13'
  []
  [permeability_aquitard]
    type = PorousFlowPermeabilityConst
    block = aquitard
    permeability = '1E-16 0 0 0 1E-16 0 0 0 1E-17'
  []
[]

[UserObjects]
  [initial_mesh]
    type = SolutionUserObject
    execute_on = INITIAL
    mesh = ex02_mesh.e
    timestep = LATEST
    system_variables = cosflow_depth
  []
  [baseflow]
    type = PorousFlowSumQuantity
  []
[]

[Postprocessors]
  [baseflow_kg]
    type = PorousFlowPlotQuantity
    uo = baseflow
    outputs = 'none'
  []
  [dt]
    type = TimestepSize
    outputs = 'none'
  []
  [baseflow_l_per_m_per_day]
    type = FunctionValuePostprocessor
    function = baseflow_rate
    indirect_dependencies = 'baseflow_kg dt'
  []
[]

[Preconditioning]
  [smp]
    type = SMP
    full = true
    # following 2 lines are not mandatory, but illustrate a popular preconditioner choice in groundwater models
    petsc_options_iname = '-pc_type -sub_pc_type  -pc_asm_overlap'
    petsc_options_value = ' asm      ilu           2              '
  []
[]

[Executioner]
  type = Transient
  solve_type = Newton
  dt = 1E6
  [TimeStepper]
    type = FunctionDT
    function = 'max(1E6, t)'
  []
  end_time = 1E12
  nl_abs_tol = 1E-13
[]

[Outputs]
  print_linear_residuals = false
  [ex]
    type = Exodus
    execute_on = final
  []
  [csv]
    type = CSV
  []
[]

(modules/porous_flow/test/tests/jacobian/hcs01.i)

# apply a half-cubic sink flux and observe the correct behavior
[Mesh]
  type = GeneratedMesh
  dim = 3
  nx = 1
  ny = 1
  nz = 2
  xmin = -1
  xmax = 1
  ymin = -1
  ymax = 1
[]

[GlobalParams]
  PorousFlowDictator = dictator
[]

[Variables]
  [ppwater]
  []
  [ppgas]
  []
  [massfrac_ph0_sp0]
  []
  [massfrac_ph1_sp0]
  []
[]

[UserObjects]
  [dictator]
    type = PorousFlowDictator
    porous_flow_vars = 'ppwater ppgas massfrac_ph0_sp0 massfrac_ph1_sp0'
    number_fluid_phases = 2
    number_fluid_components = 2
  []
  [pc]
    type = PorousFlowCapillaryPressureVG
    m = 0.5
    alpha = 1
  []
[]

[ICs]
  [ppwater]
    type = RandomIC
    variable = ppwater
    min = -1
    max = 0
  []
  [ppgas]
    type = RandomIC
    variable = ppgas
    min = 0
    max = 1
  []
  [massfrac_ph0_sp0]
    type = RandomIC
    variable = massfrac_ph0_sp0
    min = 0
    max = 1
  []
  [massfrac_ph1_sp0]
    type = RandomIC
    variable = massfrac_ph1_sp0
    min = 0
    max = 1
  []
[]

[Kernels]
  [dummy_ppwater]
    type = TimeDerivative
    variable = ppwater
  []
  [dummy_ppgas]
    type = TimeDerivative
    variable = ppgas
  []
  [dummy_m00]
    type = TimeDerivative
    variable = massfrac_ph0_sp0
  []
  [dummy_m10]
    type = TimeDerivative
    variable = massfrac_ph1_sp0
  []
[]

[FluidProperties]
  [simple_fluid0]
    type = SimpleFluidProperties
    bulk_modulus = 1.5
    density0 = 1
    thermal_expansion = 0
    viscosity = 1
  []
  [simple_fluid1]
    type = SimpleFluidProperties
    bulk_modulus = 0.5
    density0 = 0.5
    thermal_expansion = 0
    viscosity = 1.4
  []
[]

[Materials]
  [temperature]
    type = PorousFlowTemperature
  []
  [ppss]
    type = PorousFlow2PhasePP
    phase0_porepressure = ppwater
    phase1_porepressure = ppgas
    capillary_pressure = pc
  []
  [massfrac]
    type = PorousFlowMassFraction
    mass_fraction_vars = 'massfrac_ph0_sp0 massfrac_ph1_sp0'
  []
  [simple_fluid0]
    type = PorousFlowSingleComponentFluid
    fp = simple_fluid0
    phase = 0
  []
  [simple_fluid1]
    type = PorousFlowSingleComponentFluid
    fp = simple_fluid1
    phase = 1
  []
  [permeability]
    type = PorousFlowPermeabilityConst
    permeability = '1 0 0 0 2 0 0 0 3'
  []
  [relperm0]
    type = PorousFlowRelativePermeabilityCorey
    n = 2
    phase = 0
  []
  [relperm1]
    type = PorousFlowRelativePermeabilityCorey
    n = 3
    phase = 1
  []
[]

[BCs]
  [flux_w]
    type = PorousFlowHalfCubicSink
    boundary = 'left'
    center = 0.1
    cutoff = -1.1
    max = 2.2
    variable = ppwater
    mass_fraction_component = 0
    fluid_phase = 0
    use_relperm = true
    use_mobility = true
    flux_function = 'x*y'
  []
  [flux_g]
    type = PorousFlowHalfCubicSink
    boundary = 'top left front'
    center = 0.5
    cutoff = -1.1
    max = -2.2
    mass_fraction_component = 0
    variable = ppgas
    fluid_phase = 1
    use_relperm = true
    use_mobility = true
    flux_function = '-x*y'
  []
  [flux_1]
    type = PorousFlowHalfCubicSink
    boundary = 'right'
    center = -0.1
    cutoff = -1.1
    max = 1.2
    mass_fraction_component = 1
    variable = massfrac_ph0_sp0
    fluid_phase = 1
    use_relperm = true
    use_mobility = true
    flux_function = '-1.1*x*y'
  []
  [flux_2]
    type = PorousFlowHalfCubicSink
    boundary = 'bottom'
    center = 3.2
    cutoff = -1.1
    max = 1.2
    mass_fraction_component = 1
    variable = massfrac_ph1_sp0
    fluid_phase = 1
    use_relperm = true
    use_mobility = true
    flux_function = '0.5*x*y'
  []
[]

[Preconditioning]
  [check]
    type = SMP
    full = true
    petsc_options_iname = '-ksp_type -pc_type -snes_atol -snes_rtol -snes_max_it -snes_type'
    petsc_options_value = 'bcgs bjacobi 1E-15 1E-10 10000 test'
  []
[]

[Executioner]
  type = Transient
  solve_type = Newton
  dt = 1
  end_time = 2
[]

[Outputs]
  file_base = hcs01
[]

(modules/porous_flow/test/tests/jacobian/hcs02.i)

# apply a half-cubic heat sink flux
[Mesh]
  type = GeneratedMesh
  dim = 3
  nx = 1
  ny = 1
  nz = 2
  xmin = -1
  xmax = 1
  ymin = -1
  ymax = 1
[]

[GlobalParams]
  PorousFlowDictator = dictator
[]

[Variables]
  [temp]
  []
[]

[UserObjects]
  [dictator]
    type = PorousFlowDictator
    porous_flow_vars = temp
    number_fluid_phases = 0
    number_fluid_components = 0
  []
[]

[ICs]
  [temp]
    type = RandomIC
    variable = temp
    min = -1
    max = 0
  []
[]

[Kernels]
  [dummy_temp]
    type = TimeDerivative
    variable = temp
  []
[]

[Materials]
  [temperature]
    type = PorousFlowTemperature
    temperature = temp
  []
[]

[BCs]
  [flux_w]
    type = PorousFlowHalfCubicSink
    boundary = 'left'
    center = 0.1
    cutoff = -1.1
    max = 2.2
    variable = temp
    flux_function = 'x*y'
  []
  [flux_g]
    type = PorousFlowHalfCubicSink
    boundary = 'top left front'
    center = 0.5
    cutoff = -1.1
    max = -2.2
    variable = temp
    flux_function = '-x*y'
  []
  [flux_1]
    type = PorousFlowHalfCubicSink
    boundary = 'right'
    center = -0.1
    cutoff = -1.1
    max = 1.2
    variable = temp
    flux_function = '-1.1*x*y'
  []
  [flux_2]
    type = PorousFlowHalfCubicSink
    boundary = 'bottom'
    center = 3.2
    cutoff = -1.1
    max = 1.2
    variable = temp
    flux_function = '0.5*x*y'
  []
[]


[Preconditioning]
  [check]
    type = SMP
    full = true
    petsc_options_iname = '-ksp_type -pc_type -snes_atol -snes_rtol -snes_max_it -snes_type'
    petsc_options_value = 'bcgs bjacobi 1E-15 1E-10 10000 test'
  []
[]

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

[Outputs]
  file_base = hcs02
[]

(modules/porous_flow/test/tests/sinks/s06.i)

# apply a half-cubic sink flux and observe the correct behavior
[Mesh]
  type = GeneratedMesh
  dim = 3
  nx = 1
  ny = 1
  nz = 1
  xmin = 0
  xmax = 1
  ymin = 0
  ymax = 1
  zmin = 0
  zmax = 2
[]

[GlobalParams]
  PorousFlowDictator = dictator
[]

[UserObjects]
  [dictator]
    type = PorousFlowDictator
    porous_flow_vars = 'pp'
    number_fluid_phases = 1
    number_fluid_components = 1
  []
  [pc]
    type = PorousFlowCapillaryPressureVG
    m = 0.5
    alpha = 1.1
  []
[]

[Variables]
  [pp]
  []
[]

[ICs]
  [pp]
    type = FunctionIC
    variable = pp
    function = x*(y+1)
  []
[]

[Kernels]
  [mass0]
    type = PorousFlowMassTimeDerivative
    fluid_component = 0
    variable = pp
  []
[]

[FluidProperties]
  [simple_fluid]
    type = SimpleFluidProperties
    bulk_modulus = 1.3
    density0 = 1.1
    thermal_expansion = 0
    viscosity = 1.1
  []
[]

[Materials]
  [temperature]
    type = PorousFlowTemperature
  []
  [ppss]
    type = PorousFlow1PhaseP
    porepressure = pp
    capillary_pressure = pc
  []
  [massfrac]
    type = PorousFlowMassFraction
  []
  [simple_fluid]
    type = PorousFlowSingleComponentFluid
    fp = simple_fluid
    phase = 0
  []
  [porosity]
    type = PorousFlowPorosityConst
    porosity = 0.1
  []
  [permeability]
    type = PorousFlowPermeabilityConst
    permeability = '1E-5 0 0 0 1E-5 0 0 0 1E-5'
  []
  [relperm]
    type = PorousFlowRelativePermeabilityCorey
    n = 2
    phase = 0
  []
[]

[AuxVariables]
  [flux_out]
  []
[]

[Functions]
  [mass10]
    type = ParsedFunction
    expression = 'vol*por*dens0*exp(pp/bulk)*if(pp>=0,1,pow(1+pow(-al*pp,1.0/(1-m)),-m))'
    symbol_names = 'vol por dens0 pp bulk al m'
    symbol_values = '0.25 0.1 1.1 p10 1.3 1.1 0.5'
  []
  [rate10]
    type = ParsedFunction
    expression = 'fcn*if(pp>center,m,if(pp<themin,0,m/c/c/c*(2*(pp-center)+c)*((pp-center)-c)*((pp-center)-c)))'
    symbol_names = 'm fcn pp  center sd  themin c'
    symbol_values = '2 3   p10 0.9    0.5 0.1   -0.8'
  []
  [mass10_expect]
    type = ParsedFunction
    expression = 'mass_prev-rate*area*dt'
    symbol_names = 'mass_prev rate     area dt'
    symbol_values = 'm10_prev  m10_rate 0.5 2E-3'
  []
  [mass11]
    type = ParsedFunction
    expression = 'vol*por*dens0*exp(pp/bulk)*if(pp>=0,1,pow(1+pow(-al*pp,1.0/(1-m)),-m))'
    symbol_names = 'vol por dens0 pp bulk al m'
    symbol_values = '0.25 0.1 1.1 p11 1.3 1.1 0.5'
  []
  [rate11]
    type = ParsedFunction
    expression = 'fcn*if(pp>center,m,if(pp<themin,0,m/c/c/c*(2*(pp-center)+c)*((pp-center)-c)*((pp-center)-c)))'
    symbol_names = 'm fcn pp  center sd  themin c'
    symbol_values = '2 3   p11 0.9    0.5 0.1   -0.8'
  []
  [mass11_expect]
    type = ParsedFunction
    expression = 'mass_prev-rate*area*dt'
    symbol_names = 'mass_prev rate     area dt'
    symbol_values = 'm11_prev  m11_rate 0.5 2E-3'
  []
[]

[Postprocessors]
  [flux00]
    type = PointValue
    variable = flux_out
    point = '0 0 0'
  []
  [flux01]
    type = PointValue
    variable = flux_out
    point = '0 1 0'
  []
  [flux10]
    type = PointValue
    variable = flux_out
    point = '1 0 0'
  []
  [flux11]
    type = PointValue
    variable = flux_out
    point = '1 1 0'
  []
  [p00]
    type = PointValue
    point = '0 0 0'
    variable = pp
    execute_on = 'initial timestep_end'
  []
  [p10]
    type = PointValue
    point = '1 0 0'
    variable = pp
    execute_on = 'initial timestep_end'
  []
  [m10]
    type = FunctionValuePostprocessor
    function = mass10
    execute_on = 'initial timestep_end'
  []
  [m10_prev]
    type = FunctionValuePostprocessor
    function = mass10
    execute_on = 'timestep_begin'
    outputs = 'console'
  []
  [m10_rate]
    type = FunctionValuePostprocessor
    function = rate10
    execute_on = 'timestep_end'
  []
  [m10_expect]
    type = FunctionValuePostprocessor
    function = mass10_expect
    execute_on = 'timestep_end'
  []
  [p01]
    type = PointValue
    point = '0 1 0'
    variable = pp
    execute_on = 'initial timestep_end'
  []
  [p11]
    type = PointValue
    point = '1 1 0'
    variable = pp
    execute_on = 'initial timestep_end'
  []
  [m11]
    type = FunctionValuePostprocessor
    function = mass11
    execute_on = 'initial timestep_end'
  []
  [m11_prev]
    type = FunctionValuePostprocessor
    function = mass11
    execute_on = 'timestep_begin'
    outputs = 'console'
  []
  [m11_rate]
    type = FunctionValuePostprocessor
    function = rate11
    execute_on = 'timestep_end'
  []
  [m11_expect]
    type = FunctionValuePostprocessor
    function = mass11_expect
    execute_on = 'timestep_end'
  []
[]

[BCs]
  [flux]
    type = PorousFlowHalfCubicSink
    boundary = 'left right'
    max = 2
    cutoff = -0.8
    center = 0.9
    variable = pp
    use_mobility = false
    use_relperm = false
    fluid_phase = 0
    flux_function = 3
    save_in = flux_out
  []
[]

[Preconditioning]
  [andy]
    type = SMP
    full = true
    petsc_options_iname = '-ksp_type -pc_type -sub_pc_type -snes_max_it -sub_pc_factor_shift_type -pc_asm_overlap'
    petsc_options_value = 'gmres asm lu 10000 NONZERO 2'
  []
[]

[Executioner]
  type = Transient
  solve_type = Newton
  dt = 2E-3
  end_time = 6E-2
  nl_rel_tol = 1E-12
  nl_abs_tol = 1E-12
[]

[Outputs]
  file_base = s06
  [console]
    type = Console
    execute_on = 'nonlinear linear'
    time_step_interval = 5
  []
  [csv]
    type = CSV
    execute_on = 'timestep_end'
    time_step_interval = 3
  []
[]

(modules/porous_flow/examples/groundwater/ex02_abstraction.i)

# Abstraction groundwater model.  See groundwater_models.md for a detailed description
[Mesh]
  [from_steady_state]
    type = FileMeshGenerator
    file = gold/ex02_steady_state_ex.e
  []
[]

[GlobalParams]
  PorousFlowDictator = dictator
[]

[Variables]
  [pp]
  []
[]

[ICs]
  [pp]
    type = FunctionIC
    variable = pp
    function = steady_state_pp
  []
[]

[BCs]
  [rainfall_recharge]
    type = PorousFlowSink
    boundary = zmax
    variable = pp
    flux_function = -1E-6  # recharge of 0.1mm/day = 1E-4m3/m2/day = 0.1kg/m2/day ~ 1E-6kg/m2/s
  []
  [evapotranspiration]
    type = PorousFlowHalfCubicSink
    boundary = zmax
    variable = pp
    center = 0.0
    cutoff = -5E4 # roots of depth 5m.  5m of water = 5E4 Pa
    use_mobility = true
    fluid_phase = 0
    # Assume pan evaporation of 4mm/day = 4E-3m3/m2/day = 4kg/m2/day ~ 4E-5kg/m2/s
    # Assume that if permeability was 1E-10m^2 and water table at topography then ET acts as pan strength
    # Because use_mobility = true, then 4E-5 = maximum_flux = max * perm * density / visc = max * 1E-4, so max = 40
    max = 40
  []
[]

[DiracKernels]
  inactive = polyline_sink_borehole
  [river]
    type = PorousFlowPolyLineSink
    SumQuantityUO = baseflow
    point_file = ex02_river.bh
    # Assume a perennial river.
    # Assume the river has an incision depth of 1m and a stage height of 1.5m, and these are constant in time and uniform over the whole model.  Hence, if groundwater head is 0.5m (5000Pa) there will be no baseflow and leakage.
    p_or_t_vals = '-999995000 5000 1000005000'
    # Assume the riverbed conductance, k_zz*density*river_segment_length*river_width/riverbed_thickness/viscosity = 1E-6*river_segment_length kg/Pa/s
    fluxes = '-1E3 0 1E3'
    variable = pp
  []
  [horizontal_borehole]
    type = PorousFlowPeacemanBorehole
    SumQuantityUO = abstraction
    bottom_p_or_t = -1E5
    unit_weight = '0 0 -1E4'
    character = 1.0
    point_file = ex02.bh
    variable = pp
  []
  [polyline_sink_borehole]
    type = PorousFlowPolyLineSink
    SumQuantityUO = abstraction
    fluxes = '-0.4 0 0.4'
    p_or_t_vals = '-1E8 0 1E8'
    point_file = ex02.bh
    variable = pp
  []
[]

[Functions]
  [steady_state_pp]
    type = SolutionFunction
    from_variable = pp
    solution = steady_state_solution
  []
  [baseflow_rate]
    type = ParsedFunction
    symbol_names = 'baseflow_kg dt'
    symbol_values = 'baseflow_kg dt'
    expression = 'baseflow_kg / dt * 24.0 * 3600.0 / 400.0'
  []
  [abstraction_rate]
    type = ParsedFunction
    symbol_names = 'abstraction_kg dt'
    symbol_values = 'abstraction_kg dt'
    expression = 'abstraction_kg / dt * 24.0 * 3600.0'
  []
[]

[AuxVariables]
  [ini_pp]
  []
  [pp_change]
  []
[]

[AuxKernels]
  [ini_pp]
    type = FunctionAux
    variable = ini_pp
    function = steady_state_pp
    execute_on = INITIAL
  []
  [pp_change]
    type = ParsedAux
    variable = pp_change
    coupled_variables = 'pp ini_pp'
    expression = 'pp - ini_pp'
  []
[]

[PorousFlowUnsaturated]
  fp = simple_fluid
  porepressure = pp
[]

[FluidProperties]
  [simple_fluid]
    type = SimpleFluidProperties
  []
[]

[Materials]
  [porosity_everywhere]
    type = PorousFlowPorosityConst
    porosity = 0.05
  []
  [permeability_aquifers]
    type = PorousFlowPermeabilityConst
    block = 'top_aquifer bot_aquifer'
    permeability = '1E-12 0 0 0 1E-12 0 0 0 1E-13'
  []
  [permeability_aquitard]
    type = PorousFlowPermeabilityConst
    block = aquitard
    permeability = '1E-16 0 0 0 1E-16 0 0 0 1E-17'
  []
[]

[UserObjects]
  [steady_state_solution]
    type = SolutionUserObject
    execute_on = INITIAL
    mesh = gold/ex02_steady_state_ex.e
    timestep = LATEST
    system_variables = pp
  []
  [baseflow]
    type = PorousFlowSumQuantity
  []
  [abstraction]
    type = PorousFlowSumQuantity
  []
[]

[Postprocessors]
  [baseflow_kg]
    type = PorousFlowPlotQuantity
    uo = baseflow
    outputs = 'none'
  []
  [dt]
    type = TimestepSize
    outputs = 'none'
  []
  [baseflow_l_per_m_per_day]
    type = FunctionValuePostprocessor
    function = baseflow_rate
    indirect_dependencies = 'baseflow_kg dt'
  []
  [abstraction_kg]
    type = PorousFlowPlotQuantity
    uo = abstraction
    outputs = 'none'
  []
  [abstraction_kg_per_day]
    type = FunctionValuePostprocessor
    function = abstraction_rate
    indirect_dependencies = 'abstraction_kg dt'
  []
[]

[Preconditioning]
  [smp]
    type = SMP
    full = true
    # following 2 lines are not mandatory, but illustrate a popular preconditioner choice in groundwater models
    petsc_options_iname = '-pc_type -sub_pc_type  -pc_asm_overlap'
    petsc_options_value = ' asm      ilu           2              '
  []
[]

[Executioner]
  type = Transient
  solve_type = Newton
  dt = 100
  [TimeStepper]
    type = FunctionDT
    function = 'max(100, t)'
  []
  end_time = 8.64E5 # 10 days
  nl_abs_tol = 1E-11
[]

[Outputs]
  print_linear_residuals = false
  [ex]
    type = Exodus
    execute_on = final
  []
  [csv]
    type = CSV
  []
[]

Input Parameters
Input Files