Point and line sources and sinks

This page may be read in conjunction with the theory behind Dirac Kernels in PorousFlow.

Geometric tests

The test suite contains many core tests that demonstrate:

  • when a point sink is placed at a node, it withdraws fluid (or heat) from only that node, for instance

# fully-saturated situation with a poly-line sink at one
# of the nodes.  Because there is no fluid flow, the
# other nodes should not experience any change in
# porepressure.
# The poly-line sink has length=2 and weight=0.1, and
# extracts fluid at a constant rate of 1 kg.m^-1.s^-1.
# Therefore, in 1 second it will have extracted a total
# of 0.2 kg.
# The porosity is 0.1, and the elemental volume is 2,
# so the fluid mass at the node in question = 0.2 * density / 4,
# where the 4 is the number of nodes in the element.
# In this simulation density = dens0 * exp(P / bulk), with
# dens0 = 100, and bulk = 20 MPa.
# The initial porepressure P0 = 10 MPa, so the final (after
# 1 second of simulation) is
# P(t=1) = 0.950879 MPa
#
[Mesh<<<{"href": "../../../../syntax/Mesh/index.html"}>>>]
  type = GeneratedMesh
  dim = 2
  nx = 1
  ny = 1
  xmin = 0
  xmax = 2
  ymin = 0
  ymax = 1
[]

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

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

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

[UserObjects<<<{"href": "../../../../syntax/UserObjects/index.html"}>>>]
  [dictator]
    type = PorousFlowDictator<<<{"description": "Holds information on the PorousFlow variable names", "href": "../../../../source/userobjects/PorousFlowDictator.html"}>>>
    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'
    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
  []
  [pls_total_outflow_mass]
    type = PorousFlowSumQuantity<<<{"description": "Records total mass flowing into a borehole", "href": "../../../../source/userobjects/PorousFlowSumQuantity.html"}>>>
  []
  [pc]
    type = PorousFlowCapillaryPressureVG<<<{"description": "van Genuchten capillary pressure", "href": "../../../../source/userobjects/PorousFlowCapillaryPressureVG.html"}>>>
    m<<<{"description": "van Genuchten exponent m. Must be between 0 and 1, and optimally should be set to >0.5"}>>> = 0.5
    alpha<<<{"description": "van Genuchten parameter alpha. Must be positive"}>>> = 1e-7
  []
[]

[FluidProperties<<<{"href": "../../../../syntax/FluidProperties/index.html"}>>>]
  [simple_fluid]
    type = SimpleFluidProperties<<<{"description": "Fluid properties for a simple fluid with a constant bulk density", "href": "../../../../source/fluidproperties/SimpleFluidProperties.html"}>>>
    bulk_modulus<<<{"description": "Constant bulk modulus (Pa)"}>>> = 2e7
    density0<<<{"description": "Density at zero pressure and zero temperature"}>>> = 100
    thermal_expansion<<<{"description": "Constant coefficient of thermal expansion (1/K)"}>>> = 0
  []
[]

[Materials<<<{"href": "../../../../syntax/Materials/index.html"}>>>]
  [temperature]
    type = PorousFlowTemperature<<<{"description": "Material to provide temperature at the quadpoints or nodes and derivatives of it with respect to the PorousFlow variables", "href": "../../../../source/materials/PorousFlowTemperature.html"}>>>
  []
  [ppss]
    type = PorousFlow1PhaseP<<<{"description": "This Material is used for the partially saturated single-phase situation where porepressure is the primary variable", "href": "../../../../source/materials/PorousFlow1PhaseP.html"}>>>
    porepressure<<<{"description": "Variable that represents the porepressure of the single phase"}>>> = pp
    capillary_pressure<<<{"description": "Name of the UserObject defining the capillary pressure"}>>> = pc
  []
  [massfrac]
    type = PorousFlowMassFraction<<<{"description": "This Material forms a std::vector<std::vector ...> of mass-fractions out of the individual mass fractions", "href": "../../../../source/materials/PorousFlowMassFraction.html"}>>>
  []
  [simple_fluid]
    type = PorousFlowSingleComponentFluid<<<{"description": "This Material calculates fluid properties at the quadpoints or nodes for a single component fluid", "href": "../../../../source/materials/PorousFlowSingleComponentFluid.html"}>>>
    fp<<<{"description": "The name of the user object for fluid properties"}>>> = simple_fluid
    phase<<<{"description": "The phase number"}>>> = 0
  []
  [porosity]
    type = PorousFlowPorosityConst<<<{"description": "This Material calculates the porosity assuming it is constant", "href": "../../../../source/materials/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
  []
[]

[DiracKernels<<<{"href": "../../../../syntax/DiracKernels/index.html"}>>>]
  [pls]
    type = PorousFlowPolyLineSink<<<{"description": "Approximates a polyline sink by using a number of point sinks with given weighting whose positions are read from a file.  NOTE: if you are using PorousFlowPorosity that depends on volumetric strain, you should set strain_at_nearest_qp=true in your GlobalParams, to ensure the nodal Porosity Material uses the volumetric strain at the Dirac quadpoints, and can therefore be computed", "href": "../../../../source/dirackernels/PorousFlowPolyLineSink.html"}>>>
    fluid_phase<<<{"description": "The fluid phase whose pressure (and potentially mobility, enthalpy, etc) controls the flux to the line sink.  For p_or_t=temperature, and without any use_*, this parameter is irrelevant"}>>> = 0
    point_file<<<{"description": "The file containing the coordinates of the points and their weightings that approximate the line sink.  The physical meaning of the weightings depend on the scenario, eg, they may be borehole radii.  Each line in the file must contain a space-separated weight and coordinate, viz r x y z.  For boreholes, the last point in the file is defined as the borehole bottom, where the borehole pressure is bottom_pressure.  If your file contains just one point, you must also specify the line_length and line_direction parameters.  Note that you will get segementation faults if your points do not lie within your mesh!"}>>> = pls01_21.bh
    line_length<<<{"description": "Line length.  Note this is only used if there is only one point in the point_file."}>>> = 2
    SumQuantityUO<<<{"description": "User Object of type=PorousFlowSumQuantity in which to place the total outflow from the line sink for each time step."}>>> = pls_total_outflow_mass
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = pp
    p_or_t_vals<<<{"description": "Tuple of pressure (or temperature) values.  Must be monotonically increasing."}>>> = '0 1E7'
    fluxes<<<{"description": "Tuple of flux values (measured in kg.m^-1.s^-1 if no 'use_*' are employed).  These flux values are multiplied by the line-segment length to achieve a flux in kg.s^-1.  A piecewise-linear fit is performed to the (p_or_t_vals,flux) pairs to obtain the flux at any arbitrary pressure (or temperature).  If a quad-point pressure is less than the first pressure value, the first flux value is used.  If quad-point pressure exceeds the final pressure value, the final flux value is used.  This flux is OUT of the medium: hence positive values of flux means this will be a SINK, while negative values indicate this flux will be a SOURCE."}>>> = '1 1'
  []
[]

[Postprocessors<<<{"href": "../../../../syntax/Postprocessors/index.html"}>>>]
  [pls_report]
    type = PorousFlowPlotQuantity<<<{"description": "Extracts the value from the PorousFlowSumQuantity UserObject", "href": "../../../../source/postprocessors/PorousFlowPlotQuantity.html"}>>>
    uo<<<{"description": "PorousFlowSumQuantity user object name that holds the required information"}>>> = pls_total_outflow_mass
  []

  [fluid_mass0]
    type = PorousFlowFluidMass<<<{"description": "Calculates the mass of a fluid component in a region", "href": "../../../../source/postprocessors/PorousFlowFluidMass.html"}>>>
    execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = timestep_begin
  []

  [fluid_mass1]
    type = PorousFlowFluidMass<<<{"description": "Calculates the mass of a fluid component in a region", "href": "../../../../source/postprocessors/PorousFlowFluidMass.html"}>>>
    execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = timestep_end
  []

  [zmass_error]
    type = FunctionValuePostprocessor<<<{"description": "Computes the value of a supplied function at a single point (scalable)", "href": "../../../../source/postprocessors/FunctionValuePostprocessor.html"}>>>
    function<<<{"description": "The function which supplies the postprocessor value."}>>> = mass_bal_fcn
    execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = timestep_end
    indirect_dependencies<<<{"description": "If the evaluated function depends on other postprocessors they must be listed here to ensure proper dependency resolution"}>>> = 'fluid_mass1 fluid_mass0 pls_report'
  []

  [p00]
    type = PointValue<<<{"description": "Compute the value of a variable at a specified location", "href": "../../../../source/postprocessors/PointValue.html"}>>>
    variable<<<{"description": "The name of the variable that this postprocessor operates on."}>>> = pp
    point<<<{"description": "The physical point where the solution will be evaluated."}>>> = '0 0 0'
    execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = timestep_end
  []
  [p01]
    type = PointValue<<<{"description": "Compute the value of a variable at a specified location", "href": "../../../../source/postprocessors/PointValue.html"}>>>
    variable<<<{"description": "The name of the variable that this postprocessor operates on."}>>> = pp
    point<<<{"description": "The physical point where the solution will be evaluated."}>>> = '0 1 0'
    execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = timestep_end
  []
  [p20]
    type = PointValue<<<{"description": "Compute the value of a variable at a specified location", "href": "../../../../source/postprocessors/PointValue.html"}>>>
    variable<<<{"description": "The name of the variable that this postprocessor operates on."}>>> = pp
    point<<<{"description": "The physical point where the solution will be evaluated."}>>> = '2 0 0'
    execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = timestep_end
  []
  [p21]
    type = PointValue<<<{"description": "Compute the value of a variable at a specified location", "href": "../../../../source/postprocessors/PointValue.html"}>>>
    variable<<<{"description": "The name of the variable that this postprocessor operates on."}>>> = pp
    point<<<{"description": "The physical point where the solution will be evaluated."}>>> = '2 1 0'
    execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = timestep_end
  []
[]

[Functions<<<{"href": "../../../../syntax/Functions/index.html"}>>>]
  [mass_bal_fcn]
    type = ParsedFunction<<<{"description": "Function created by parsing a string", "href": "../../../../source/functions/MooseParsedFunction.html"}>>>
    expression<<<{"description": "The user defined function."}>>> = abs((a-c+d)/2/(a+c))
    symbol_names<<<{"description": "Symbols (excluding t,x,y,z) that are bound to the values provided by the corresponding items in the vals vector."}>>> = 'a c d'
    symbol_values<<<{"description": "Constant numeric values, postprocessor names, function names, and scalar variables corresponding to the symbols in symbol_names."}>>> = 'fluid_mass1 fluid_mass0 pls_report'
  []
[]

[Preconditioning<<<{"href": "../../../../syntax/Preconditioning/index.html"}>>>]
  [usual]
    type = SMP<<<{"description": "Single matrix preconditioner (SMP) builds a preconditioner using user defined off-diagonal parts of the Jacobian.", "href": "../../../../source/preconditioners/SingleMatrixPreconditioner.html"}>>>
    full<<<{"description": "Set to true if you want the full set of couplings between variables simply for convenience so you don't have to set every off_diag_row and off_diag_column combination."}>>> = true
    petsc_options<<<{"description": "Singleton PETSc options"}>>> = '-snes_converged_reason'
    petsc_options_iname<<<{"description": "Names of PETSc name/value pairs"}>>> = '-ksp_type -pc_type -snes_atol -snes_rtol -snes_max_it -ksp_max_it'
    petsc_options_value<<<{"description": "Values of PETSc name/value pairs (must correspond with \"petsc_options_iname\""}>>> = 'bcgs bjacobi 1E-10 1E-10 10000 30'
  []
[]

[Executioner<<<{"href": "../../../../syntax/Executioner/index.html"}>>>]
  type = Transient
  end_time = 1
  dt = 1
  solve_type = NEWTON
[]

[Outputs<<<{"href": "../../../../syntax/Outputs/index.html"}>>>]
  file_base<<<{"description": "Common file base name to be utilized with all output objects"}>>> = pls01
  exodus<<<{"description": "Output the results using the default settings for Exodus output."}>>> = false
  csv<<<{"description": "Output the scalar variable and postprocessors to a *.csv file using the default CSV output."}>>> = true
  execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = timestep_end
[]
(modules/porous_flow/test/tests/dirackernels/pls01.i)
  • when a point sink that is proportional to mobility (or relative permeability, etc) is placed in an element where some nodes have zero mobility (or relative permeability, etc), then fluid (or heat) is not extracted from those nodes, for instance

# Test that the upwinding works correctly.
#
# A poly-line sink sits at the centre of the element.
# It has length=4 and weight=0.5, and extracts fluid
# at a constant rate of
# (1 * relative_permeability) kg.m^-1.s^-1
# Since it sits at the centre of the element, it extracts
# equally from each node, so the rate of extraction from
# each node is
# (0.5 * relative_permeability) kg.s^-1
# including the length and weight effects.
#
# There is no fluid flow.
#
# The initial conditions are such that all nodes have
# relative_permeability=0, except for one which has
# relative_permeaility = 1.  Therefore, all nodes should
# remain at their initial porepressure, except the one.
#
# The porosity is 0.1, and the elemental volume is 2,
# so the fluid mass at the node in question = 0.2 * density / 4,
# where the 4 is the number of nodes in the element.
# In this simulation density = dens0 * exp(P / bulk), with
# dens0 = 100, and bulk = 20 MPa.
# The initial porepressure P0 = 10 MPa, so the final (after
# 1 second of simulation) is
# P(t=1) = 8.748592 MPa
[Mesh<<<{"href": "../../../../syntax/Mesh/index.html"}>>>]
  type = GeneratedMesh
  dim = 2
  nx = 1
  ny = 1
  xmin = 0
  xmax = 2
  ymin = 0
  ymax = 1
[]

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

[Variables<<<{"href": "../../../../syntax/Variables/index.html"}>>>]
  [pp]
  []
[]

[ICs<<<{"href": "../../../../syntax/ICs/index.html"}>>>]
  [pp]
    type = FunctionIC<<<{"description": "An initial condition that uses a normal function of x, y, z to produce values (and optionally gradients) for a field variable.", "href": "../../../../source/ics/FunctionIC.html"}>>>
    variable<<<{"description": "The variable this initial condition is supposed to provide values for."}>>> = pp
    #function = if((x<1)&(y<0.5),1E7,-1E7)
    function<<<{"description": "The initial condition function."}>>> = if((x<1)&(y>0.5),1E7,-1E7)
    #function = if((x>1)&(y<0.5),1E7,-1E7)
    #function = if((x>1)&(y>0.5),1E7,-1E7)
  []
[]

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

[UserObjects<<<{"href": "../../../../syntax/UserObjects/index.html"}>>>]
  [dictator]
    type = PorousFlowDictator<<<{"description": "Holds information on the PorousFlow variable names", "href": "../../../../source/userobjects/PorousFlowDictator.html"}>>>
    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'
    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
  []
  [pls_total_outflow_mass]
    type = PorousFlowSumQuantity<<<{"description": "Records total mass flowing into a borehole", "href": "../../../../source/userobjects/PorousFlowSumQuantity.html"}>>>
  []
  [pc]
    type = PorousFlowCapillaryPressureVG<<<{"description": "van Genuchten capillary pressure", "href": "../../../../source/userobjects/PorousFlowCapillaryPressureVG.html"}>>>
    m<<<{"description": "van Genuchten exponent m. Must be between 0 and 1, and optimally should be set to >0.5"}>>> = 0.5
    alpha<<<{"description": "van Genuchten parameter alpha. Must be positive"}>>> = 1e-7
  []
[]

[FluidProperties<<<{"href": "../../../../syntax/FluidProperties/index.html"}>>>]
  [simple_fluid]
    type = SimpleFluidProperties<<<{"description": "Fluid properties for a simple fluid with a constant bulk density", "href": "../../../../source/fluidproperties/SimpleFluidProperties.html"}>>>
    bulk_modulus<<<{"description": "Constant bulk modulus (Pa)"}>>> = 2e7
    density0<<<{"description": "Density at zero pressure and zero temperature"}>>> = 100
    thermal_expansion<<<{"description": "Constant coefficient of thermal expansion (1/K)"}>>> = 0
  []
[]

[Materials<<<{"href": "../../../../syntax/Materials/index.html"}>>>]
  [temperature]
    type = PorousFlowTemperature<<<{"description": "Material to provide temperature at the quadpoints or nodes and derivatives of it with respect to the PorousFlow variables", "href": "../../../../source/materials/PorousFlowTemperature.html"}>>>
  []
  [ppss]
    type = PorousFlow1PhaseP<<<{"description": "This Material is used for the partially saturated single-phase situation where porepressure is the primary variable", "href": "../../../../source/materials/PorousFlow1PhaseP.html"}>>>
    porepressure<<<{"description": "Variable that represents the porepressure of the single phase"}>>> = pp
    capillary_pressure<<<{"description": "Name of the UserObject defining the capillary pressure"}>>> = pc
  []
  [massfrac]
    type = PorousFlowMassFraction<<<{"description": "This Material forms a std::vector<std::vector ...> of mass-fractions out of the individual mass fractions", "href": "../../../../source/materials/PorousFlowMassFraction.html"}>>>
  []
  [simple_fluid]
    type = PorousFlowSingleComponentFluid<<<{"description": "This Material calculates fluid properties at the quadpoints or nodes for a single component fluid", "href": "../../../../source/materials/PorousFlowSingleComponentFluid.html"}>>>
    fp<<<{"description": "The name of the user object for fluid properties"}>>> = simple_fluid
    phase<<<{"description": "The phase number"}>>> = 0
  []
  [porosity]
    type = PorousFlowPorosityConst<<<{"description": "This Material calculates the porosity assuming it is constant", "href": "../../../../source/materials/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
  []
  [relperm]
    type = PorousFlowRelativePermeabilityFLAC<<<{"description": "This Material calculates relative permeability of a phase using a model inspired by FLAC", "href": "../../../../source/materials/PorousFlowRelativePermeabilityFLAC.html"}>>>
    phase<<<{"description": "The phase number"}>>> = 0
    m<<<{"description": "relperm = (1 + m)seff^m - m seff^(m+1)"}>>> = 2
    s_res<<<{"description": "The residual saturation of the phase j. Must be between 0 and 1"}>>> = 0.99
    sum_s_res<<<{"description": "Sum of residual saturations over all phases.  Must be between 0 and 1"}>>> = 0.99
  []
[]

[DiracKernels<<<{"href": "../../../../syntax/DiracKernels/index.html"}>>>]
  [pls]
    type = PorousFlowPolyLineSink<<<{"description": "Approximates a polyline sink by using a number of point sinks with given weighting whose positions are read from a file.  NOTE: if you are using PorousFlowPorosity that depends on volumetric strain, you should set strain_at_nearest_qp=true in your GlobalParams, to ensure the nodal Porosity Material uses the volumetric strain at the Dirac quadpoints, and can therefore be computed", "href": "../../../../source/dirackernels/PorousFlowPolyLineSink.html"}>>>
    fluid_phase<<<{"description": "The fluid phase whose pressure (and potentially mobility, enthalpy, etc) controls the flux to the line sink.  For p_or_t=temperature, and without any use_*, this parameter is irrelevant"}>>> = 0
    point_file<<<{"description": "The file containing the coordinates of the points and their weightings that approximate the line sink.  The physical meaning of the weightings depend on the scenario, eg, they may be borehole radii.  Each line in the file must contain a space-separated weight and coordinate, viz r x y z.  For boreholes, the last point in the file is defined as the borehole bottom, where the borehole pressure is bottom_pressure.  If your file contains just one point, you must also specify the line_length and line_direction parameters.  Note that you will get segementation faults if your points do not lie within your mesh!"}>>> = pls03.bh
    use_relative_permeability<<<{"description": "Multiply the flux by the fluid relative permeability"}>>> = true
    line_length<<<{"description": "Line length.  Note this is only used if there is only one point in the point_file."}>>> = 4
    SumQuantityUO<<<{"description": "User Object of type=PorousFlowSumQuantity in which to place the total outflow from the line sink for each time step."}>>> = pls_total_outflow_mass
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = pp
    p_or_t_vals<<<{"description": "Tuple of pressure (or temperature) values.  Must be monotonically increasing."}>>> = '0 1E7'
    fluxes<<<{"description": "Tuple of flux values (measured in kg.m^-1.s^-1 if no 'use_*' are employed).  These flux values are multiplied by the line-segment length to achieve a flux in kg.s^-1.  A piecewise-linear fit is performed to the (p_or_t_vals,flux) pairs to obtain the flux at any arbitrary pressure (or temperature).  If a quad-point pressure is less than the first pressure value, the first flux value is used.  If quad-point pressure exceeds the final pressure value, the final flux value is used.  This flux is OUT of the medium: hence positive values of flux means this will be a SINK, while negative values indicate this flux will be a SOURCE."}>>> = '1 1'
  []
[]

[Postprocessors<<<{"href": "../../../../syntax/Postprocessors/index.html"}>>>]
  [pls_report]
    type = PorousFlowPlotQuantity<<<{"description": "Extracts the value from the PorousFlowSumQuantity UserObject", "href": "../../../../source/postprocessors/PorousFlowPlotQuantity.html"}>>>
    uo<<<{"description": "PorousFlowSumQuantity user object name that holds the required information"}>>> = pls_total_outflow_mass
  []

  [fluid_mass0]
    type = PorousFlowFluidMass<<<{"description": "Calculates the mass of a fluid component in a region", "href": "../../../../source/postprocessors/PorousFlowFluidMass.html"}>>>
    execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = timestep_begin
  []

  [fluid_mass1]
    type = PorousFlowFluidMass<<<{"description": "Calculates the mass of a fluid component in a region", "href": "../../../../source/postprocessors/PorousFlowFluidMass.html"}>>>
    execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = timestep_end
  []

  [zmass_error]
    type = FunctionValuePostprocessor<<<{"description": "Computes the value of a supplied function at a single point (scalable)", "href": "../../../../source/postprocessors/FunctionValuePostprocessor.html"}>>>
    function<<<{"description": "The function which supplies the postprocessor value."}>>> = mass_bal_fcn
    execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = timestep_end
    indirect_dependencies<<<{"description": "If the evaluated function depends on other postprocessors they must be listed here to ensure proper dependency resolution"}>>> = 'fluid_mass1 fluid_mass0 pls_report'
  []

  [p00]
    type = PointValue<<<{"description": "Compute the value of a variable at a specified location", "href": "../../../../source/postprocessors/PointValue.html"}>>>
    variable<<<{"description": "The name of the variable that this postprocessor operates on."}>>> = pp
    point<<<{"description": "The physical point where the solution will be evaluated."}>>> = '0 0 0'
    execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = timestep_end
  []
  [p01]
    type = PointValue<<<{"description": "Compute the value of a variable at a specified location", "href": "../../../../source/postprocessors/PointValue.html"}>>>
    variable<<<{"description": "The name of the variable that this postprocessor operates on."}>>> = pp
    point<<<{"description": "The physical point where the solution will be evaluated."}>>> = '0 1 0'
    execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = timestep_end
  []
  [p20]
    type = PointValue<<<{"description": "Compute the value of a variable at a specified location", "href": "../../../../source/postprocessors/PointValue.html"}>>>
    variable<<<{"description": "The name of the variable that this postprocessor operates on."}>>> = pp
    point<<<{"description": "The physical point where the solution will be evaluated."}>>> = '2 0 0'
    execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = timestep_end
  []
  [p21]
    type = PointValue<<<{"description": "Compute the value of a variable at a specified location", "href": "../../../../source/postprocessors/PointValue.html"}>>>
    variable<<<{"description": "The name of the variable that this postprocessor operates on."}>>> = pp
    point<<<{"description": "The physical point where the solution will be evaluated."}>>> = '2 1 0'
    execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = timestep_end
  []
[]

[Functions<<<{"href": "../../../../syntax/Functions/index.html"}>>>]
  [mass_bal_fcn]
    type = ParsedFunction<<<{"description": "Function created by parsing a string", "href": "../../../../source/functions/MooseParsedFunction.html"}>>>
    expression<<<{"description": "The user defined function."}>>> = abs((a-c+d)/2/(a+c))
    symbol_names<<<{"description": "Symbols (excluding t,x,y,z) that are bound to the values provided by the corresponding items in the vals vector."}>>> = 'a c d'
    symbol_values<<<{"description": "Constant numeric values, postprocessor names, function names, and scalar variables corresponding to the symbols in symbol_names."}>>> = 'fluid_mass1 fluid_mass0 pls_report'
  []
[]

[Preconditioning<<<{"href": "../../../../syntax/Preconditioning/index.html"}>>>]
  [usual]
    type = SMP<<<{"description": "Single matrix preconditioner (SMP) builds a preconditioner using user defined off-diagonal parts of the Jacobian.", "href": "../../../../source/preconditioners/SingleMatrixPreconditioner.html"}>>>
    full<<<{"description": "Set to true if you want the full set of couplings between variables simply for convenience so you don't have to set every off_diag_row and off_diag_column combination."}>>> = true
    petsc_options<<<{"description": "Singleton PETSc options"}>>> = '-snes_converged_reason'
    petsc_options_iname<<<{"description": "Names of PETSc name/value pairs"}>>> = '-ksp_type -pc_type -snes_atol -snes_rtol -snes_max_it -ksp_max_it'
    petsc_options_value<<<{"description": "Values of PETSc name/value pairs (must correspond with \"petsc_options_iname\""}>>> = 'bcgs bjacobi 1E-10 1E-10 10000 30'
  []
[]

[Executioner<<<{"href": "../../../../syntax/Executioner/index.html"}>>>]
  type = Transient
  end_time = 1
  dt = 1
  solve_type = NEWTON
[]

[Outputs<<<{"href": "../../../../syntax/Outputs/index.html"}>>>]
  file_base<<<{"description": "Common file base name to be utilized with all output objects"}>>> = pls03
  exodus<<<{"description": "Output the results using the default settings for Exodus output."}>>> = false
  csv<<<{"description": "Output the scalar variable and postprocessors to a *.csv file using the default CSV output."}>>> = true
  execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = timestep_end
[]
(modules/porous_flow/test/tests/dirackernels/pls03.i)

Basic point sources/sinks

The following input file extracts fluid at the rate (kg.s)

# Test PorousFlowSquarePulsePointSource DiracKernel

[Mesh<<<{"href": "../../../../syntax/Mesh/index.html"}>>>]
  type = GeneratedMesh
  dim = 2
  bias_x = 1.1
  bias_y = 1.1
  ymax = 1
  xmax = 1
[]

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

[Variables<<<{"href": "../../../../syntax/Variables/index.html"}>>>]
  [pp]
  []
[]

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

[UserObjects<<<{"href": "../../../../syntax/UserObjects/index.html"}>>>]
  [dictator]
    type = PorousFlowDictator<<<{"description": "Holds information on the PorousFlow variable names", "href": "../../../../source/userobjects/PorousFlowDictator.html"}>>>
    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
    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
  []
  [pc]
    type = PorousFlowCapillaryPressureVG<<<{"description": "van Genuchten capillary pressure", "href": "../../../../source/userobjects/PorousFlowCapillaryPressureVG.html"}>>>
    m<<<{"description": "van Genuchten exponent m. Must be between 0 and 1, and optimally should be set to >0.5"}>>> = 0.5
    alpha<<<{"description": "van Genuchten parameter alpha. Must be positive"}>>> = 1e-7
  []
[]

[FluidProperties<<<{"href": "../../../../syntax/FluidProperties/index.html"}>>>]
  [simple_fluid]
    type = SimpleFluidProperties<<<{"description": "Fluid properties for a simple fluid with a constant bulk density", "href": "../../../../source/fluidproperties/SimpleFluidProperties.html"}>>>
    bulk_modulus<<<{"description": "Constant bulk modulus (Pa)"}>>> = 2e9
    density0<<<{"description": "Density at zero pressure and zero temperature"}>>> = 1000
    thermal_expansion<<<{"description": "Constant coefficient of thermal expansion (1/K)"}>>> = 0
  []
[]

[Materials<<<{"href": "../../../../syntax/Materials/index.html"}>>>]
  [temperature]
    type = PorousFlowTemperature<<<{"description": "Material to provide temperature at the quadpoints or nodes and derivatives of it with respect to the PorousFlow variables", "href": "../../../../source/materials/PorousFlowTemperature.html"}>>>
  []
  [ppss]
    type = PorousFlow1PhaseP<<<{"description": "This Material is used for the partially saturated single-phase situation where porepressure is the primary variable", "href": "../../../../source/materials/PorousFlow1PhaseP.html"}>>>
    porepressure<<<{"description": "Variable that represents the porepressure of the single phase"}>>> = pp
    capillary_pressure<<<{"description": "Name of the UserObject defining the capillary pressure"}>>> = pc
  []
  [massfrac]
    type = PorousFlowMassFraction<<<{"description": "This Material forms a std::vector<std::vector ...> of mass-fractions out of the individual mass fractions", "href": "../../../../source/materials/PorousFlowMassFraction.html"}>>>
  []
  [simple_fluid]
    type = PorousFlowSingleComponentFluid<<<{"description": "This Material calculates fluid properties at the quadpoints or nodes for a single component fluid", "href": "../../../../source/materials/PorousFlowSingleComponentFluid.html"}>>>
    fp<<<{"description": "The name of the user object for fluid properties"}>>> = simple_fluid
    phase<<<{"description": "The phase number"}>>> = 0
  []
  [porosity]
    type = PorousFlowPorosityConst<<<{"description": "This Material calculates the porosity assuming it is constant", "href": "../../../../source/materials/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.2
  []
[]

[Postprocessors<<<{"href": "../../../../syntax/Postprocessors/index.html"}>>>]
  [total_mass]
    type = PorousFlowFluidMass<<<{"description": "Calculates the mass of a fluid component in a region", "href": "../../../../source/postprocessors/PorousFlowFluidMass.html"}>>>
    execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = 'initial timestep_end'
  []
[]

[Preconditioning<<<{"href": "../../../../syntax/Preconditioning/index.html"}>>>]
  [smp]
    type = SMP<<<{"description": "Single matrix preconditioner (SMP) builds a preconditioner using user defined off-diagonal parts of the Jacobian.", "href": "../../../../source/preconditioners/SingleMatrixPreconditioner.html"}>>>
    full<<<{"description": "Set to true if you want the full set of couplings between variables simply for convenience so you don't have to set every off_diag_row and off_diag_column combination."}>>> = true
  []
[]

[Executioner<<<{"href": "../../../../syntax/Executioner/index.html"}>>>]
  type = Transient
  solve_type = Newton
  nl_abs_tol = 1e-14
  dt = 200
  end_time = 2000
[]

[Outputs<<<{"href": "../../../../syntax/Outputs/index.html"}>>>]
  perf_graph<<<{"description": "Enable printing of the performance graph to the screen (Console)"}>>> = true
  file_base<<<{"description": "Common file base name to be utilized with all output objects"}>>> = squarepulse1
  csv<<<{"description": "Output the scalar variable and postprocessors to a *.csv file using the default CSV output."}>>> = true
  execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = 'initial timestep_end'
  [con]
    output_linear<<<{"description": "Specifies whether output occurs on each PETSc linear residual evaluation"}>>> = true
    type = Console<<<{"description": "Object for screen output.", "href": "../../../../source/outputs/Console.html"}>>>
  []
[]

[ICs<<<{"href": "../../../../syntax/ICs/index.html"}>>>]
  [PressureIC]
    variable<<<{"description": "The variable this initial condition is supposed to provide values for."}>>> = pp
    type = ConstantIC<<<{"description": "Sets a constant field value.", "href": "../../../../source/ics/ConstantIC.html"}>>>
    value<<<{"description": "The value to be set in IC"}>>> = 20e6
  []
[]

[DiracKernels<<<{"href": "../../../../syntax/DiracKernels/index.html"}>>>]
  [sink1]
    type = PorousFlowSquarePulsePointSource<<<{"description": "Point source (or sink) that adds (removes) fluid at a constant mass flux rate for times between the specified start and end times.", "href": "../../../../source/dirackernels/PorousFlowSquarePulsePointSource.html"}>>>
    start_time<<<{"description": "The time at which the source will start (Default is 0)"}>>> = 100
    end_time<<<{"description": "The time at which the source will end (Default is 1e30)"}>>> = 300
    point<<<{"description": "The x,y,z coordinates of the point source (sink)"}>>> = '0.5 0.5 0'
    mass_flux<<<{"description": "The mass flux at this point in kg/s (positive is flux in, negative is flux out)"}>>> = -0.1
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = pp
  []
  [sink]
    type = PorousFlowSquarePulsePointSource<<<{"description": "Point source (or sink) that adds (removes) fluid at a constant mass flux rate for times between the specified start and end times.", "href": "../../../../source/dirackernels/PorousFlowSquarePulsePointSource.html"}>>>
    start_time<<<{"description": "The time at which the source will start (Default is 0)"}>>> = 600
    end_time<<<{"description": "The time at which the source will end (Default is 1e30)"}>>> = 1400
    point<<<{"description": "The x,y,z coordinates of the point source (sink)"}>>> = '0.5 0.5 0'
    mass_flux<<<{"description": "The mass flux at this point in kg/s (positive is flux in, negative is flux out)"}>>> = -0.1
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = pp
  []
  [source]
    point<<<{"description": "The x,y,z coordinates of the point source (sink)"}>>> = '0.5 0.5 0'
    start_time<<<{"description": "The time at which the source will start (Default is 0)"}>>> = 1500
    mass_flux<<<{"description": "The mass flux at this point in kg/s (positive is flux in, negative is flux out)"}>>> = 0.2
    end_time<<<{"description": "The time at which the source will end (Default is 1e30)"}>>> = 2000
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = pp
    type = PorousFlowSquarePulsePointSource<<<{"description": "Point source (or sink) that adds (removes) fluid at a constant mass flux rate for times between the specified start and end times.", "href": "../../../../source/dirackernels/PorousFlowSquarePulsePointSource.html"}>>>
  []
[]
(modules/porous_flow/test/tests/dirackernels/squarepulse1.i)

MOOSE produces the expected result:

Figure 1: Results of a "square-pulsed" extraction and injection of fluid from a fluid-flow simulation.

Theis tests

In the fully-saturated, isothermal case with no mechanical coupling and constant, large fluid bulk modulus, the fluid flow equations reduce to the form conventionally used in groundwater flow: where

  • is the "head", . Here is the fluid density.

  • is the "specific storage". If the rock compressibility is ignored then , where is the fluid bulk modulus

  • is the hydraulic conductivity: , where is the permeability and is the fluid viscosity

This equation is identical to the physics described by the PorousFlowFullySaturatedMassTimeDerivative and the PorousFlowFullySaturatedDarcyBase Kernels, which may be used explicitly, or through the PorousFlowBasicTHM Action. When using these, the specific storage in the above equation is replaced by (reciprocal of the Biot modulus), and by , and by . Finally, remember that this formulation is volume based, not mass based like the remainder of PorousFlow.

Place a constant volumetric sink of strength (m.m.s) acting along an infinite line in an isotropic 3D medium. The situation is therefore two dimensional. Theis provided the solution for the head (1) which is frequently used in the groundwater literature. "Ei" is the exponential integral function, with values for small (where is the Euler number ), and for large .

This is checked using a number of PorousFlow tests (the extraction specified using mass rate, with volume rate, and from single and 2-phase systems). For instance:

# Theis problem: Flow to single sink using BasicTHM
# SinglePhase
# RZ mesh

[Mesh<<<{"href": "../../../../syntax/Mesh/index.html"}>>>]
  type = GeneratedMesh
  dim = 1
  nx = 20
  xmax = 100
  bias_x = 1.05
  coord_type = RZ
[]

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

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

[PorousFlowBasicTHM<<<{"href": "../../../../syntax/PorousFlowBasicTHM/index.html"}>>>]
  dictator_name<<<{"description": "The name of the dictator user object that is created by this Action"}>>> = dictator
  add_darcy_aux<<<{"description": "Add AuxVariables that record Darcy velocity"}>>> = false
  fp<<<{"description": "The name of the user object for fluid properties. Only needed if fluid_properties_type = PorousFlowSingleComponentFluid"}>>> = simple_fluid
  gravity<<<{"description": "Gravitational acceleration vector downwards (m/s^2)"}>>> = '0 0 0'
  multiply_by_density<<<{"description": "If true, then the Kernels for fluid flow are multiplied by the fluid density.  If false, this multiplication is not performed, which means the problem linearises, but that care must be taken when using other PorousFlow objects."}>>> = false
  porepressure<<<{"description": "The name of the porepressure variable"}>>> = pp
[]

[FluidProperties<<<{"href": "../../../../syntax/FluidProperties/index.html"}>>>]
  [simple_fluid]
    type = SimpleFluidProperties<<<{"description": "Fluid properties for a simple fluid with a constant bulk density", "href": "../../../../source/fluidproperties/SimpleFluidProperties.html"}>>>
    viscosity<<<{"description": "Constant dynamic viscosity (Pa.s)"}>>> = 0.001
  []
[]

[Materials<<<{"href": "../../../../syntax/Materials/index.html"}>>>]
  [porosity]
    type = PorousFlowPorosityConst<<<{"description": "This Material calculates the porosity assuming it is constant", "href": "../../../../source/materials/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.05
  []
  [biot_mod]
    type = PorousFlowConstantBiotModulus<<<{"description": "Computes the Biot Modulus, which is assumed to be constant for all time.  Sometimes 1 / BiotModulus is called storativity", "href": "../../../../source/materials/PorousFlowConstantBiotModulus.html"}>>>
    fluid_bulk_modulus<<<{"description": "Fluid bulk modulus"}>>> = 2E9
    biot_coefficient<<<{"description": "Biot coefficient"}>>> = 1.0
  []
  [permeability]
    type = PorousFlowPermeabilityConst<<<{"description": "This Material calculates the permeability tensor assuming it is constant", "href": "../../../../source/materials/PorousFlowPermeabilityConst.html"}>>>
    permeability<<<{"description": "The permeability tensor (usually in m^2), which is assumed constant for this material"}>>> = '1E-14 0 0 0 1E-14 0 0 0 1E-14'
  []
[]

[DiracKernels<<<{"href": "../../../../syntax/DiracKernels/index.html"}>>>]
  [sink]
    type = PorousFlowSquarePulsePointSource<<<{"description": "Point source (or sink) that adds (removes) fluid at a constant mass flux rate for times between the specified start and end times.", "href": "../../../../source/dirackernels/PorousFlowSquarePulsePointSource.html"}>>>
    point<<<{"description": "The x,y,z coordinates of the point source (sink)"}>>> = '0 0 0'
    mass_flux<<<{"description": "The mass flux at this point in kg/s (positive is flux in, negative is flux out)"}>>> = -0.16E-3 # recall this is a volumetric flux because multiply_by_density = false in the Action, so this corresponds to a mass_flux of 0.16 kg/s/m because density=1000
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = pp
  []
[]

[VectorPostprocessors<<<{"href": "../../../../syntax/VectorPostprocessors/index.html"}>>>]
  [pp]
    type = LineValueSampler<<<{"description": "Samples variable(s) along a specified line", "href": "../../../../source/vectorpostprocessors/LineValueSampler.html"}>>>
    num_points<<<{"description": "The number of points to sample along the line"}>>> = 25
    start_point<<<{"description": "The beginning of the line"}>>> = '0 0 0'
    end_point<<<{"description": "The ending of the line"}>>> = '100 0 0'
    sort_by<<<{"description": "What to sort the samples by"}>>> = x
    variable<<<{"description": "The names of the variables that this VectorPostprocessor operates on"}>>> = pp
  []
[]

[Preconditioning<<<{"href": "../../../../syntax/Preconditioning/index.html"}>>>]
  [smp]
    type = SMP<<<{"description": "Single matrix preconditioner (SMP) builds a preconditioner using user defined off-diagonal parts of the Jacobian.", "href": "../../../../source/preconditioners/SingleMatrixPreconditioner.html"}>>>
    full<<<{"description": "Set to true if you want the full set of couplings between variables simply for convenience so you don't have to set every off_diag_row and off_diag_column combination."}>>> = true
  []
[]

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

[Outputs<<<{"href": "../../../../syntax/Outputs/index.html"}>>>]
  perf_graph<<<{"description": "Enable printing of the performance graph to the screen (Console)"}>>> = true
  [csv]
    type = CSV<<<{"description": "Output for postprocessors, vector postprocessors, and scalar variables using comma seperated values (CSV).", "href": "../../../../source/outputs/CSV.html"}>>>
    execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = final
  []
[]
(modules/porous_flow/test/tests/dirackernels/theis_rz.i)

MOOSE agrees with the Theis solution:

Figure 2: Results of the fully-saturated, single-phase Theis simulation.

Flow from a source in a radial 1D problem admits a similarity solution, where the fluid pressure (and saturation, mass fraction, etc) is a function of the similarity variable . This may be observed in Eq. (1). Using this fact, two-phase immiscible flow may be tested. The following input file describes injection of a gas phase into a fully liquid-saturated model at a constant rate:

# Two phase Theis problem: Flow from single source
# Constant rate injection 0.5 kg/s
# 1D cylindrical mesh

[Mesh<<<{"href": "../../../../syntax/Mesh/index.html"}>>>]
  type = GeneratedMesh
  dim = 1
  nx = 100
  xmax = 2000
  bias_x = 1.05
  coord_type = RZ
  rz_coord_axis = Y
[]

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

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

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

[Kernels<<<{"href": "../../../../syntax/Kernels/index.html"}>>>]
  [mass0]
    type = PorousFlowMassTimeDerivative<<<{"description": "Derivative of fluid-component mass with respect to time.  Mass lumping to the nodes is used.", "href": "../../../../source/kernels/PorousFlowMassTimeDerivative.html"}>>>
    fluid_component<<<{"description": "The index corresponding to the component for this kernel"}>>> = 0
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = ppwater
  []
  [flux0]
    type = PorousFlowAdvectiveFlux<<<{"description": "Fully-upwinded advective flux of the component given by fluid_component", "href": "../../../../source/kernels/PorousFlowAdvectiveFlux.html"}>>>
    fluid_component<<<{"description": "The index corresponding to the fluid component for this kernel"}>>> = 0
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = ppwater
  []
  [mass1]
    type = PorousFlowMassTimeDerivative<<<{"description": "Derivative of fluid-component mass with respect to time.  Mass lumping to the nodes is used.", "href": "../../../../source/kernels/PorousFlowMassTimeDerivative.html"}>>>
    fluid_component<<<{"description": "The index corresponding to the component for this kernel"}>>> = 1
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = sgas
  []
  [flux1]
    type = PorousFlowAdvectiveFlux<<<{"description": "Fully-upwinded advective flux of the component given by fluid_component", "href": "../../../../source/kernels/PorousFlowAdvectiveFlux.html"}>>>
    fluid_component<<<{"description": "The index corresponding to the fluid component for this kernel"}>>> = 1
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = sgas
  []
[]

[UserObjects<<<{"href": "../../../../syntax/UserObjects/index.html"}>>>]
  [dictator]
    type = PorousFlowDictator<<<{"description": "Holds information on the PorousFlow variable names", "href": "../../../../source/userobjects/PorousFlowDictator.html"}>>>
    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."}>>> = 'ppwater sgas'
    number_fluid_phases<<<{"description": "The number of fluid phases in the simulation"}>>> = 2
    number_fluid_components<<<{"description": "The number of fluid components in the simulation"}>>> = 2
  []
  [pc]
    type = PorousFlowCapillaryPressureConst<<<{"description": "Constant capillary pressure", "href": "../../../../source/userobjects/PorousFlowCapillaryPressureConst.html"}>>>
    pc<<<{"description": "Constant capillary pressure (Pa). Default is 0"}>>> = 1e5
  []
[]

[FluidProperties<<<{"href": "../../../../syntax/FluidProperties/index.html"}>>>]
  [simple_fluid0]
    type = SimpleFluidProperties<<<{"description": "Fluid properties for a simple fluid with a constant bulk density", "href": "../../../../source/fluidproperties/SimpleFluidProperties.html"}>>>
    bulk_modulus<<<{"description": "Constant bulk modulus (Pa)"}>>> = 2e9
    density0<<<{"description": "Density at zero pressure and zero temperature"}>>> = 1000
    viscosity<<<{"description": "Constant dynamic viscosity (Pa.s)"}>>> = 1e-3
    thermal_expansion<<<{"description": "Constant coefficient of thermal expansion (1/K)"}>>> = 0
  []
  [simple_fluid1]
    type = SimpleFluidProperties<<<{"description": "Fluid properties for a simple fluid with a constant bulk density", "href": "../../../../source/fluidproperties/SimpleFluidProperties.html"}>>>
    bulk_modulus<<<{"description": "Constant bulk modulus (Pa)"}>>> = 2e9
    density0<<<{"description": "Density at zero pressure and zero temperature"}>>> = 10
    viscosity<<<{"description": "Constant dynamic viscosity (Pa.s)"}>>> = 1e-4
    thermal_expansion<<<{"description": "Constant coefficient of thermal expansion (1/K)"}>>> = 0
  []
[]

[Materials<<<{"href": "../../../../syntax/Materials/index.html"}>>>]
  [temperature]
    type = PorousFlowTemperature<<<{"description": "Material to provide temperature at the quadpoints or nodes and derivatives of it with respect to the PorousFlow variables", "href": "../../../../source/materials/PorousFlowTemperature.html"}>>>
  []
  [ppss]
    type = PorousFlow2PhasePS<<<{"description": "This Material calculates the 2 porepressures and the 2 saturations in a 2-phase situation, and derivatives of these with respect to the PorousFlowVariables.", "href": "../../../../source/materials/PorousFlow2PhasePS.html"}>>>
    phase0_porepressure<<<{"description": "Variable that is the porepressure of phase 0 (the liquid phase)"}>>> = ppwater
    phase1_saturation<<<{"description": "Variable that is the saturation of phase 1 (the gas phase)"}>>> = sgas
    capillary_pressure<<<{"description": "Name of the UserObject defining the capillary pressure"}>>> = pc
  []
  [massfrac]
    type = PorousFlowMassFraction<<<{"description": "This Material forms a std::vector<std::vector ...> of mass-fractions out of the individual mass fractions", "href": "../../../../source/materials/PorousFlowMassFraction.html"}>>>
    mass_fraction_vars<<<{"description": "List of variables that represent the mass fractions.  Format is 'f_ph0^c0 f_ph0^c1 f_ph0^c2 ... f_ph0^c(N-2) f_ph1^c0 f_ph1^c1 fph1^c2 ... fph1^c(N-2) ... fphP^c0 f_phP^c1 fphP^c2 ... fphP^c(N-2)' where N=num_components and P=num_phases, and it is assumed that f_ph^c(N-1)=1-sum(f_ph^c,{c,0,N-2}) so that f_ph^c(N-1) need not be given.  If no variables are provided then num_phases=1=num_components."}>>> = 'massfrac_ph0_sp0 massfrac_ph1_sp0'
  []
  [simple_fluid0]
    type = PorousFlowSingleComponentFluid<<<{"description": "This Material calculates fluid properties at the quadpoints or nodes for a single component fluid", "href": "../../../../source/materials/PorousFlowSingleComponentFluid.html"}>>>
    fp<<<{"description": "The name of the user object for fluid properties"}>>> = simple_fluid0
    phase<<<{"description": "The phase number"}>>> = 0
    compute_enthalpy<<<{"description": "Compute the fluid enthalpy"}>>> = false
    compute_internal_energy<<<{"description": "Compute the fluid internal energy"}>>> = false
  []
  [simple_fluid1]
    type = PorousFlowSingleComponentFluid<<<{"description": "This Material calculates fluid properties at the quadpoints or nodes for a single component fluid", "href": "../../../../source/materials/PorousFlowSingleComponentFluid.html"}>>>
    fp<<<{"description": "The name of the user object for fluid properties"}>>> = simple_fluid1
    phase<<<{"description": "The phase number"}>>> = 1
    compute_enthalpy<<<{"description": "Compute the fluid enthalpy"}>>> = false
    compute_internal_energy<<<{"description": "Compute the fluid internal energy"}>>> = false
  []
  [porosity]
    type = PorousFlowPorosityConst<<<{"description": "This Material calculates the porosity assuming it is constant", "href": "../../../../source/materials/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.2
  []
  [permeability]
    type = PorousFlowPermeabilityConst<<<{"description": "This Material calculates the permeability tensor assuming it is constant", "href": "../../../../source/materials/PorousFlowPermeabilityConst.html"}>>>
    permeability<<<{"description": "The permeability tensor (usually in m^2), which is assumed constant for this material"}>>> = '1e-12 0 0 0 1e-12 0 0 0 1e-12'
  []
  [relperm_water]
    type = PorousFlowRelativePermeabilityCorey<<<{"description": "This Material calculates relative permeability of the fluid phase, using the simple Corey model ((S-S_res)/(1-sum(S_res)))^n", "href": "../../../../source/materials/PorousFlowRelativePermeabilityCorey.html"}>>>
    n<<<{"description": "The Corey exponent of the phase."}>>> = 1
    phase<<<{"description": "The phase number"}>>> = 0
  []
  [relperm_gas]
    type = PorousFlowRelativePermeabilityCorey<<<{"description": "This Material calculates relative permeability of the fluid phase, using the simple Corey model ((S-S_res)/(1-sum(S_res)))^n", "href": "../../../../source/materials/PorousFlowRelativePermeabilityCorey.html"}>>>
    n<<<{"description": "The Corey exponent of the phase."}>>> = 1
    phase<<<{"description": "The phase number"}>>> = 1
  []
[]

[BCs<<<{"href": "../../../../syntax/BCs/index.html"}>>>]
  [rightwater]
    type = DirichletBC<<<{"description": "Imposes the essential boundary condition $u=g$, where $g$ is a constant, controllable value.", "href": "../../../../source/bcs/DirichletBC.html"}>>>
    boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = right
    value<<<{"description": "Value of the BC"}>>> = 20e6
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = ppwater
  []
[]

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

[Preconditioning<<<{"href": "../../../../syntax/Preconditioning/index.html"}>>>]
  [smp]
    type = SMP<<<{"description": "Single matrix preconditioner (SMP) builds a preconditioner using user defined off-diagonal parts of the Jacobian.", "href": "../../../../source/preconditioners/SingleMatrixPreconditioner.html"}>>>
    full<<<{"description": "Set to true if you want the full set of couplings between variables simply for convenience so you don't have to set every off_diag_row and off_diag_column combination."}>>> = true
    petsc_options<<<{"description": "Singleton PETSc options"}>>> = '-snes_converged_reason -ksp_diagonal_scale -ksp_diagonal_scale_fix -ksp_gmres_modifiedgramschmidt -snes_linesearch_monitor'
    petsc_options_iname<<<{"description": "Names of PETSc name/value pairs"}>>> = '-ksp_type -pc_type -sub_pc_type -sub_pc_factor_shift_type -pc_asm_overlap -snes_atol -snes_rtol -snes_max_it'
    petsc_options_value<<<{"description": "Values of PETSc name/value pairs (must correspond with \"petsc_options_iname\""}>>> = 'gmres      asm      lu           NONZERO                   2               1E-8       1E-10 20'
  []
[]

[Executioner<<<{"href": "../../../../syntax/Executioner/index.html"}>>>]
  type = Transient
  solve_type = Newton
  end_time = 1e4
  [TimeStepper<<<{"href": "../../../../syntax/Executioner/TimeStepper/index.html"}>>>]
    type = IterationAdaptiveDT
    dt = 10
    growth_factor = 2
  []
[]

[VectorPostprocessors<<<{"href": "../../../../syntax/VectorPostprocessors/index.html"}>>>]
  [line]
    type = NodalValueSampler<<<{"description": "Samples values of nodal variable(s).", "href": "../../../../source/vectorpostprocessors/NodalValueSampler.html"}>>>
    sort_by<<<{"description": "What to sort the samples by"}>>> = x
    variable<<<{"description": "The names of the variables that this VectorPostprocessor operates on"}>>> = 'ppwater sgas'
    execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = 'timestep_end'
  []
[]

[Postprocessors<<<{"href": "../../../../syntax/Postprocessors/index.html"}>>>]
  [ppwater]
    type = PointValue<<<{"description": "Compute the value of a variable at a specified location", "href": "../../../../source/postprocessors/PointValue.html"}>>>
    point<<<{"description": "The physical point where the solution will be evaluated."}>>> = '4 0 0'
    variable<<<{"description": "The name of the variable that this postprocessor operates on."}>>> = ppwater
  []
  [sgas]
    type = PointValue<<<{"description": "Compute the value of a variable at a specified location", "href": "../../../../source/postprocessors/PointValue.html"}>>>
    point<<<{"description": "The physical point where the solution will be evaluated."}>>> = '4 0 0'
    variable<<<{"description": "The name of the variable that this postprocessor operates on."}>>> = sgas
  []
  [massgas]
    type = PorousFlowFluidMass<<<{"description": "Calculates the mass of a fluid component in a region", "href": "../../../../source/postprocessors/PorousFlowFluidMass.html"}>>>
    fluid_component<<<{"description": "The index corresponding to the fluid component that this Postprocessor acts on"}>>> = 1
  []
[]

[Outputs<<<{"href": "../../../../syntax/Outputs/index.html"}>>>]
  file_base<<<{"description": "Common file base name to be utilized with all output objects"}>>> = theis3
  print_linear_residuals<<<{"description": "Enable printing of linear residuals to the screen (Console)"}>>> = false
  perf_graph<<<{"description": "Enable printing of the performance graph to the screen (Console)"}>>> = true
  [csv]
    type = CSV<<<{"description": "Output for postprocessors, vector postprocessors, and scalar variables using comma seperated values (CSV).", "href": "../../../../source/outputs/CSV.html"}>>>
    execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = timestep_end
    execute_vector_postprocessors_on<<<{"description": "Enable/disable the output of VectorPostprocessors"}>>> = final
  []
[]
(modules/porous_flow/test/tests/dirackernels/theis3.i)

Figure 3 shows the comparison of similarity solutions calculated with either fixed radial distance or fixed time. In this case, good agreement is observed between the two results for both liquid pressure and gas saturation.

Figure 3: Results of the 2-phase radial injection simulation.

Further tests of this kind, involving multi-component, multi-phase flow, including mutual dissolution of gas and liquid phases may be found here.

Peaceman borehole fluxes

The test suite that checks that the Peaceman flux (2) is correctly implemented. A vertical borehole is placed through the centre of a single element, and fluid flow to the borehole as a function of porepressure is measured. The tests are

  • A production borehole with , with a fully-saturated medium.

# fully-saturated
# production
[Mesh<<<{"href": "../../../../syntax/Mesh/index.html"}>>>]
  type = GeneratedMesh
  dim = 3
  nx = 1
  ny = 1
  nz = 1
  xmin = -1
  xmax = 1
  ymin = -1
  ymax = 1
  zmin = -1
  zmax = 1
[]

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

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

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

[UserObjects<<<{"href": "../../../../syntax/UserObjects/index.html"}>>>]
  [borehole_total_outflow_mass]
    type = PorousFlowSumQuantity<<<{"description": "Records total mass flowing into a borehole", "href": "../../../../source/userobjects/PorousFlowSumQuantity.html"}>>>
  []
  [dictator]
    type = PorousFlowDictator<<<{"description": "Holds information on the PorousFlow variable names", "href": "../../../../source/userobjects/PorousFlowDictator.html"}>>>
    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'
    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
  []
  [pc]
    type = PorousFlowCapillaryPressureVG<<<{"description": "van Genuchten capillary pressure", "href": "../../../../source/userobjects/PorousFlowCapillaryPressureVG.html"}>>>
    m<<<{"description": "van Genuchten exponent m. Must be between 0 and 1, and optimally should be set to >0.5"}>>> = 0.5
    alpha<<<{"description": "van Genuchten parameter alpha. Must be positive"}>>> = 1e-7
  []
[]

[FluidProperties<<<{"href": "../../../../syntax/FluidProperties/index.html"}>>>]
  [simple_fluid]
    type = SimpleFluidProperties<<<{"description": "Fluid properties for a simple fluid with a constant bulk density", "href": "../../../../source/fluidproperties/SimpleFluidProperties.html"}>>>
    bulk_modulus<<<{"description": "Constant bulk modulus (Pa)"}>>> = 2e9
    viscosity<<<{"description": "Constant dynamic viscosity (Pa.s)"}>>> = 1e-3
    density0<<<{"description": "Density at zero pressure and zero temperature"}>>> = 1000
    thermal_expansion<<<{"description": "Constant coefficient of thermal expansion (1/K)"}>>> = 0
  []
[]

[Materials<<<{"href": "../../../../syntax/Materials/index.html"}>>>]
  [temperature]
    type = PorousFlowTemperature<<<{"description": "Material to provide temperature at the quadpoints or nodes and derivatives of it with respect to the PorousFlow variables", "href": "../../../../source/materials/PorousFlowTemperature.html"}>>>
  []
  [ppss]
    type = PorousFlow1PhaseP<<<{"description": "This Material is used for the partially saturated single-phase situation where porepressure is the primary variable", "href": "../../../../source/materials/PorousFlow1PhaseP.html"}>>>
    porepressure<<<{"description": "Variable that represents the porepressure of the single phase"}>>> = pp
    capillary_pressure<<<{"description": "Name of the UserObject defining the capillary pressure"}>>> = pc
  []
  [massfrac]
    type = PorousFlowMassFraction<<<{"description": "This Material forms a std::vector<std::vector ...> of mass-fractions out of the individual mass fractions", "href": "../../../../source/materials/PorousFlowMassFraction.html"}>>>
  []
  [simple_fluid]
    type = PorousFlowSingleComponentFluid<<<{"description": "This Material calculates fluid properties at the quadpoints or nodes for a single component fluid", "href": "../../../../source/materials/PorousFlowSingleComponentFluid.html"}>>>
    fp<<<{"description": "The name of the user object for fluid properties"}>>> = simple_fluid
    phase<<<{"description": "The phase number"}>>> = 0
  []
  [porosity]
    type = PorousFlowPorosityConst<<<{"description": "This Material calculates the porosity assuming it is constant", "href": "../../../../source/materials/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
  []
  [permeability]
    type = PorousFlowPermeabilityConst<<<{"description": "This Material calculates the permeability tensor assuming it is constant", "href": "../../../../source/materials/PorousFlowPermeabilityConst.html"}>>>
    permeability<<<{"description": "The permeability tensor (usually in m^2), which is assumed constant for this material"}>>> = '1E-12 0 0 0 1E-12 0 0 0 1E-12'
  []
  [relperm]
    type = PorousFlowRelativePermeabilityCorey<<<{"description": "This Material calculates relative permeability of the fluid phase, using the simple Corey model ((S-S_res)/(1-sum(S_res)))^n", "href": "../../../../source/materials/PorousFlowRelativePermeabilityCorey.html"}>>>
    n<<<{"description": "The Corey exponent of the phase."}>>> = 2
    phase<<<{"description": "The phase number"}>>> = 0
  []
[]

[DiracKernels<<<{"href": "../../../../syntax/DiracKernels/index.html"}>>>]
  [bh]
    type = PorousFlowPeacemanBorehole<<<{"description": "Approximates a borehole in the mesh using the Peaceman approach, ie using a number of point sinks with given radii whose positions are read from a file.  NOTE: if you are using PorousFlowPorosity that depends on volumetric strain, you should set strain_at_nearest_qp=true in your GlobalParams, to ensure the nodal Porosity Material uses the volumetric strain at the Dirac quadpoints, and can therefore be computed", "href": "../../../../source/dirackernels/PorousFlowPeacemanBorehole.html"}>>>

    # Because the Variable for this Sink is pp, and pp is associated
    # with the fluid-mass conservation equation, this sink is extracting
    # fluid mass (and not heat energy or something else)
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = pp

    # The following specfies that the total fluid mass coming out of
    # the porespace via this sink in this timestep should be recorded
    # in the pls_total_outflow_mass UserObject
    SumQuantityUO<<<{"description": "User Object of type=PorousFlowSumQuantity in which to place the total outflow from the line sink for each time step."}>>> = borehole_total_outflow_mass

    # The following file defines the polyline geometry
    # which is just two points in this particular example
    point_file<<<{"description": "The file containing the coordinates of the points and their weightings that approximate the line sink.  The physical meaning of the weightings depend on the scenario, eg, they may be borehole radii.  Each line in the file must contain a space-separated weight and coordinate, viz r x y z.  For boreholes, the last point in the file is defined as the borehole bottom, where the borehole pressure is bottom_pressure.  If your file contains just one point, you must also specify the line_length and line_direction parameters.  Note that you will get segementation faults if your points do not lie within your mesh!"}>>> = bh02.bh

    # First, we want Peacemans f to be a function of porepressure (and not
    # temperature or something else).  So bottom_p_or_t is actually porepressure
    function_of<<<{"description": "Modifying functions will be a function of either pressure and permeability (eg, for boreholes that pump fluids) or temperature and thermal conductivity (eg, for boreholes that pump pure heat with no fluid flow)"}>>> = pressure
    fluid_phase<<<{"description": "The fluid phase whose pressure (and potentially mobility, enthalpy, etc) controls the flux to the line sink.  For p_or_t=temperature, and without any use_*, this parameter is irrelevant"}>>> = 0

    # The bottomhole pressure
    bottom_p_or_t<<<{"description": "For function_of=pressure, this function is the pressure at the bottom of the borehole, otherwise it is the temperature at the bottom of the borehole."}>>> = 0

    # In this example there is no increase of the wellbore pressure
    # due to gravity:
    unit_weight<<<{"description": "(fluid_density*gravitational_acceleration) as a vector pointing downwards.  Note that the borehole pressure at a given z position is bottom_p_or_t + unit_weight*(q - q_bottom), where q=(x,y,z) and q_bottom=(x,y,z) of the bottom point of the borehole.  The analogous formula holds for function_of=temperature.  If you don't want bottomhole pressure (or temperature) to vary in the borehole just set unit_weight=0.  Typical value is = (0,0,-1E4), for water"}>>> = '0 0 0'

    # PeacemanBoreholes should almost always have use_mobility = true
    use_mobility<<<{"description": "Multiply the flux by the fluid mobility"}>>> = true

    # This is a production wellbore (a sink of fluid that removes fluid from porespace)
    character<<<{"description": "If zero then borehole does nothing.  If positive the borehole acts as a sink (production well) for porepressure > borehole pressure, and does nothing otherwise.  If negative the borehole acts as a source (injection well) for porepressure < borehole pressure, and does nothing otherwise.  The flow rate to/from the borehole is multiplied by |character|, so usually character = +/- 1, but you can specify other quantities to provide an overall scaling to the flow if you like."}>>> = 1
  []
[]

[Postprocessors<<<{"href": "../../../../syntax/Postprocessors/index.html"}>>>]
  [bh_report]
    type = PorousFlowPlotQuantity<<<{"description": "Extracts the value from the PorousFlowSumQuantity UserObject", "href": "../../../../source/postprocessors/PorousFlowPlotQuantity.html"}>>>
    uo<<<{"description": "PorousFlowSumQuantity user object name that holds the required information"}>>> = borehole_total_outflow_mass
  []
  [fluid_mass0]
    type = PorousFlowFluidMass<<<{"description": "Calculates the mass of a fluid component in a region", "href": "../../../../source/postprocessors/PorousFlowFluidMass.html"}>>>
    execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = timestep_begin
  []

  [fluid_mass1]
    type = PorousFlowFluidMass<<<{"description": "Calculates the mass of a fluid component in a region", "href": "../../../../source/postprocessors/PorousFlowFluidMass.html"}>>>
    execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = timestep_end
  []

  [zmass_error]
    type = FunctionValuePostprocessor<<<{"description": "Computes the value of a supplied function at a single point (scalable)", "href": "../../../../source/postprocessors/FunctionValuePostprocessor.html"}>>>
    function<<<{"description": "The function which supplies the postprocessor value."}>>> = mass_bal_fcn
    execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = timestep_end
    indirect_dependencies<<<{"description": "If the evaluated function depends on other postprocessors they must be listed here to ensure proper dependency resolution"}>>> = 'fluid_mass1 fluid_mass0 bh_report'
  []

  [p0]
    type = PointValue<<<{"description": "Compute the value of a variable at a specified location", "href": "../../../../source/postprocessors/PointValue.html"}>>>
    variable<<<{"description": "The name of the variable that this postprocessor operates on."}>>> = pp
    point<<<{"description": "The physical point where the solution will be evaluated."}>>> = '0 0 0'
    execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = timestep_end
  []
[]

[Functions<<<{"href": "../../../../syntax/Functions/index.html"}>>>]
  [mass_bal_fcn]
    type = ParsedFunction<<<{"description": "Function created by parsing a string", "href": "../../../../source/functions/MooseParsedFunction.html"}>>>
    expression<<<{"description": "The user defined function."}>>> = abs((a-c+d)/2/(a+c))
    symbol_names<<<{"description": "Symbols (excluding t,x,y,z) that are bound to the values provided by the corresponding items in the vals vector."}>>> = 'a c d'
    symbol_values<<<{"description": "Constant numeric values, postprocessor names, function names, and scalar variables corresponding to the symbols in symbol_names."}>>> = 'fluid_mass1 fluid_mass0 bh_report'
  []
[]

[Preconditioning<<<{"href": "../../../../syntax/Preconditioning/index.html"}>>>]
  [usual]
    type = SMP<<<{"description": "Single matrix preconditioner (SMP) builds a preconditioner using user defined off-diagonal parts of the Jacobian.", "href": "../../../../source/preconditioners/SingleMatrixPreconditioner.html"}>>>
    full<<<{"description": "Set to true if you want the full set of couplings between variables simply for convenience so you don't have to set every off_diag_row and off_diag_column combination."}>>> = true
    petsc_options<<<{"description": "Singleton PETSc options"}>>> = '-snes_converged_reason'
    petsc_options_iname<<<{"description": "Names of PETSc name/value pairs"}>>> = '-ksp_type -pc_type -snes_atol -snes_rtol -snes_max_it -ksp_max_it'
    petsc_options_value<<<{"description": "Values of PETSc name/value pairs (must correspond with \"petsc_options_iname\""}>>> = 'bcgs bjacobi 1E-10 1E-10 10000 30'
  []
[]

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

[Outputs<<<{"href": "../../../../syntax/Outputs/index.html"}>>>]
  file_base<<<{"description": "Common file base name to be utilized with all output objects"}>>> = bh02
  exodus<<<{"description": "Output the results using the default settings for Exodus output."}>>> = false
  csv<<<{"description": "Output the scalar variable and postprocessors to a *.csv file using the default CSV output."}>>> = true
  execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = timestep_end
[]
(modules/porous_flow/test/tests/dirackernels/bh02.i)
  • An injection borehole with MPa, with a fully-saturated medium.

# fully-saturated
# injection
[Mesh<<<{"href": "../../../../syntax/Mesh/index.html"}>>>]
  type = GeneratedMesh
  dim = 3
  nx = 1
  ny = 1
  nz = 1
  xmin = -1
  xmax = 1
  ymin = -1
  ymax = 1
  zmin = -1
  zmax = 1
[]

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

[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"}>>>]
  [mass0]
    type = PorousFlowMassTimeDerivative<<<{"description": "Derivative of fluid-component mass with respect to time.  Mass lumping to the nodes is used.", "href": "../../../../source/kernels/PorousFlowMassTimeDerivative.html"}>>>
    fluid_component<<<{"description": "The index corresponding to the component for this kernel"}>>> = 0
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = pp
  []
[]

[UserObjects<<<{"href": "../../../../syntax/UserObjects/index.html"}>>>]
  [dictator]
    type = PorousFlowDictator<<<{"description": "Holds information on the PorousFlow variable names", "href": "../../../../source/userobjects/PorousFlowDictator.html"}>>>
    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'
    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
  []
  [borehole_total_outflow_mass]
    type = PorousFlowSumQuantity<<<{"description": "Records total mass flowing into a borehole", "href": "../../../../source/userobjects/PorousFlowSumQuantity.html"}>>>
  []
  [pc]
    type = PorousFlowCapillaryPressureVG<<<{"description": "van Genuchten capillary pressure", "href": "../../../../source/userobjects/PorousFlowCapillaryPressureVG.html"}>>>
    m<<<{"description": "van Genuchten exponent m. Must be between 0 and 1, and optimally should be set to >0.5"}>>> = 0.5
    alpha<<<{"description": "van Genuchten parameter alpha. Must be positive"}>>> = 1e-7
  []
[]

[FluidProperties<<<{"href": "../../../../syntax/FluidProperties/index.html"}>>>]
  [simple_fluid]
    type = SimpleFluidProperties<<<{"description": "Fluid properties for a simple fluid with a constant bulk density", "href": "../../../../source/fluidproperties/SimpleFluidProperties.html"}>>>
    bulk_modulus<<<{"description": "Constant bulk modulus (Pa)"}>>> = 2e9
    viscosity<<<{"description": "Constant dynamic viscosity (Pa.s)"}>>> = 1e-3
    density0<<<{"description": "Density at zero pressure and zero temperature"}>>> = 1000
    thermal_expansion<<<{"description": "Constant coefficient of thermal expansion (1/K)"}>>> = 0
  []
[]

[Materials<<<{"href": "../../../../syntax/Materials/index.html"}>>>]
  [temperature]
    type = PorousFlowTemperature<<<{"description": "Material to provide temperature at the quadpoints or nodes and derivatives of it with respect to the PorousFlow variables", "href": "../../../../source/materials/PorousFlowTemperature.html"}>>>
  []
  [ppss]
    type = PorousFlow1PhaseP<<<{"description": "This Material is used for the partially saturated single-phase situation where porepressure is the primary variable", "href": "../../../../source/materials/PorousFlow1PhaseP.html"}>>>
    porepressure<<<{"description": "Variable that represents the porepressure of the single phase"}>>> = pp
    capillary_pressure<<<{"description": "Name of the UserObject defining the capillary pressure"}>>> = pc
  []
  [massfrac]
    type = PorousFlowMassFraction<<<{"description": "This Material forms a std::vector<std::vector ...> of mass-fractions out of the individual mass fractions", "href": "../../../../source/materials/PorousFlowMassFraction.html"}>>>
  []
  [simple_fluid]
    type = PorousFlowSingleComponentFluid<<<{"description": "This Material calculates fluid properties at the quadpoints or nodes for a single component fluid", "href": "../../../../source/materials/PorousFlowSingleComponentFluid.html"}>>>
    fp<<<{"description": "The name of the user object for fluid properties"}>>> = simple_fluid
    phase<<<{"description": "The phase number"}>>> = 0
  []
  [porosity]
    type = PorousFlowPorosityConst<<<{"description": "This Material calculates the porosity assuming it is constant", "href": "../../../../source/materials/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
  []
  [permeability]
    type = PorousFlowPermeabilityConst<<<{"description": "This Material calculates the permeability tensor assuming it is constant", "href": "../../../../source/materials/PorousFlowPermeabilityConst.html"}>>>
    permeability<<<{"description": "The permeability tensor (usually in m^2), which is assumed constant for this material"}>>> = '1E-12 0 0 0 1E-12 0 0 0 1E-12'
  []
  [relperm]
    type = PorousFlowRelativePermeabilityCorey<<<{"description": "This Material calculates relative permeability of the fluid phase, using the simple Corey model ((S-S_res)/(1-sum(S_res)))^n", "href": "../../../../source/materials/PorousFlowRelativePermeabilityCorey.html"}>>>
    n<<<{"description": "The Corey exponent of the phase."}>>> = 2
    phase<<<{"description": "The phase number"}>>> = 0
  []
[]

[DiracKernels<<<{"href": "../../../../syntax/DiracKernels/index.html"}>>>]
  [bh]
    type = PorousFlowPeacemanBorehole<<<{"description": "Approximates a borehole in the mesh using the Peaceman approach, ie using a number of point sinks with given radii whose positions are read from a file.  NOTE: if you are using PorousFlowPorosity that depends on volumetric strain, you should set strain_at_nearest_qp=true in your GlobalParams, to ensure the nodal Porosity Material uses the volumetric strain at the Dirac quadpoints, and can therefore be computed", "href": "../../../../source/dirackernels/PorousFlowPeacemanBorehole.html"}>>>
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = pp
    SumQuantityUO<<<{"description": "User Object of type=PorousFlowSumQuantity in which to place the total outflow from the line sink for each time step."}>>> = borehole_total_outflow_mass
    point_file<<<{"description": "The file containing the coordinates of the points and their weightings that approximate the line sink.  The physical meaning of the weightings depend on the scenario, eg, they may be borehole radii.  Each line in the file must contain a space-separated weight and coordinate, viz r x y z.  For boreholes, the last point in the file is defined as the borehole bottom, where the borehole pressure is bottom_pressure.  If your file contains just one point, you must also specify the line_length and line_direction parameters.  Note that you will get segementation faults if your points do not lie within your mesh!"}>>> = bh03.bh
    function_of<<<{"description": "Modifying functions will be a function of either pressure and permeability (eg, for boreholes that pump fluids) or temperature and thermal conductivity (eg, for boreholes that pump pure heat with no fluid flow)"}>>> = pressure
    fluid_phase<<<{"description": "The fluid phase whose pressure (and potentially mobility, enthalpy, etc) controls the flux to the line sink.  For p_or_t=temperature, and without any use_*, this parameter is irrelevant"}>>> = 0
    bottom_p_or_t<<<{"description": "For function_of=pressure, this function is the pressure at the bottom of the borehole, otherwise it is the temperature at the bottom of the borehole."}>>> = 'insitu_pp'
    unit_weight<<<{"description": "(fluid_density*gravitational_acceleration) as a vector pointing downwards.  Note that the borehole pressure at a given z position is bottom_p_or_t + unit_weight*(q - q_bottom), where q=(x,y,z) and q_bottom=(x,y,z) of the bottom point of the borehole.  The analogous formula holds for function_of=temperature.  If you don't want bottomhole pressure (or temperature) to vary in the borehole just set unit_weight=0.  Typical value is = (0,0,-1E4), for water"}>>> = '0 0 0'
    use_mobility<<<{"description": "Multiply the flux by the fluid mobility"}>>> = true
    character<<<{"description": "If zero then borehole does nothing.  If positive the borehole acts as a sink (production well) for porepressure > borehole pressure, and does nothing otherwise.  If negative the borehole acts as a source (injection well) for porepressure < borehole pressure, and does nothing otherwise.  The flow rate to/from the borehole is multiplied by |character|, so usually character = +/- 1, but you can specify other quantities to provide an overall scaling to the flow if you like."}>>> = -1
  []
[]

[Postprocessors<<<{"href": "../../../../syntax/Postprocessors/index.html"}>>>]
  [bh_report]
    type = PorousFlowPlotQuantity<<<{"description": "Extracts the value from the PorousFlowSumQuantity UserObject", "href": "../../../../source/postprocessors/PorousFlowPlotQuantity.html"}>>>
    uo<<<{"description": "PorousFlowSumQuantity user object name that holds the required information"}>>> = borehole_total_outflow_mass
  []

  [fluid_mass0]
    type = PorousFlowFluidMass<<<{"description": "Calculates the mass of a fluid component in a region", "href": "../../../../source/postprocessors/PorousFlowFluidMass.html"}>>>
    execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = timestep_begin
  []

  [fluid_mass1]
    type = PorousFlowFluidMass<<<{"description": "Calculates the mass of a fluid component in a region", "href": "../../../../source/postprocessors/PorousFlowFluidMass.html"}>>>
    execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = timestep_end
  []

  [zmass_error]
    type = FunctionValuePostprocessor<<<{"description": "Computes the value of a supplied function at a single point (scalable)", "href": "../../../../source/postprocessors/FunctionValuePostprocessor.html"}>>>
    function<<<{"description": "The function which supplies the postprocessor value."}>>> = mass_bal_fcn
    execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = timestep_end
    indirect_dependencies<<<{"description": "If the evaluated function depends on other postprocessors they must be listed here to ensure proper dependency resolution"}>>> = 'fluid_mass1 fluid_mass0 bh_report'
  []

  [p0]
    type = PointValue<<<{"description": "Compute the value of a variable at a specified location", "href": "../../../../source/postprocessors/PointValue.html"}>>>
    variable<<<{"description": "The name of the variable that this postprocessor operates on."}>>> = pp
    point<<<{"description": "The physical point where the solution will be evaluated."}>>> = '0 0 0'
    execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = timestep_end
  []
[]

[Functions<<<{"href": "../../../../syntax/Functions/index.html"}>>>]
  [mass_bal_fcn]
    type = ParsedFunction<<<{"description": "Function created by parsing a string", "href": "../../../../source/functions/MooseParsedFunction.html"}>>>
    expression<<<{"description": "The user defined function."}>>> = abs((a-c+d)/2/(a+c))
    symbol_names<<<{"description": "Symbols (excluding t,x,y,z) that are bound to the values provided by the corresponding items in the vals vector."}>>> = 'a c d'
    symbol_values<<<{"description": "Constant numeric values, postprocessor names, function names, and scalar variables corresponding to the symbols in symbol_names."}>>> = 'fluid_mass1 fluid_mass0 bh_report'
  []
  [insitu_pp]
    type = ParsedFunction<<<{"description": "Function created by parsing a string", "href": "../../../../source/functions/MooseParsedFunction.html"}>>>
    expression<<<{"description": "The user defined function."}>>> = '-2e7*z' #bh_bottom is located at z=-0.5
  []
[]

[Preconditioning<<<{"href": "../../../../syntax/Preconditioning/index.html"}>>>]
  [usual]
    type = SMP<<<{"description": "Single matrix preconditioner (SMP) builds a preconditioner using user defined off-diagonal parts of the Jacobian.", "href": "../../../../source/preconditioners/SingleMatrixPreconditioner.html"}>>>
    full<<<{"description": "Set to true if you want the full set of couplings between variables simply for convenience so you don't have to set every off_diag_row and off_diag_column combination."}>>> = true
    petsc_options<<<{"description": "Singleton PETSc options"}>>> = '-snes_converged_reason'
    petsc_options_iname<<<{"description": "Names of PETSc name/value pairs"}>>> = '-ksp_type -pc_type -snes_atol -snes_rtol -snes_max_it -ksp_max_it'
    petsc_options_value<<<{"description": "Values of PETSc name/value pairs (must correspond with \"petsc_options_iname\""}>>> = 'bcgs bjacobi 1E-10 1E-10 10000 30'
  []
[]

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

[Outputs<<<{"href": "../../../../syntax/Outputs/index.html"}>>>]
  file_base<<<{"description": "Common file base name to be utilized with all output objects"}>>> = bh03
  exodus<<<{"description": "Output the results using the default settings for Exodus output."}>>> = false
  csv<<<{"description": "Output the scalar variable and postprocessors to a *.csv file using the default CSV output."}>>> = true
  execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = timestep_end
[]
(modules/porous_flow/test/tests/dirackernels/bh03.i)
  • A production borehole with MPa with a fully-saturated medium, with use_mobility=true

# fully-saturated
# production
[Mesh<<<{"href": "../../../../syntax/Mesh/index.html"}>>>]
  type = GeneratedMesh
  dim = 3
  nx = 1
  ny = 1
  nz = 1
  xmin = -1
  xmax = 1
  ymin = -1
  ymax = 1
  zmin = -1
  zmax = 1
[]

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

[Functions<<<{"href": "../../../../syntax/Functions/index.html"}>>>]
  [dts]
    type = PiecewiseLinear<<<{"description": "Linearly interpolates between pairs of x-y data", "href": "../../../../source/functions/PiecewiseLinear.html"}>>>
    y<<<{"description": "The ordinate values"}>>> = '1E-2 1E-1 1 1E1 1E2 1E3'
    x<<<{"description": "The abscissa values"}>>> = '0 1E-1 1 1E1 1E2 1E3'
  []
[]

[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"}>>>]
  [mass0]
    type = PorousFlowMassTimeDerivative<<<{"description": "Derivative of fluid-component mass with respect to time.  Mass lumping to the nodes is used.", "href": "../../../../source/kernels/PorousFlowMassTimeDerivative.html"}>>>
    fluid_component<<<{"description": "The index corresponding to the component for this kernel"}>>> = 0
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = pp
  []
[]

[UserObjects<<<{"href": "../../../../syntax/UserObjects/index.html"}>>>]
  [dictator]
    type = PorousFlowDictator<<<{"description": "Holds information on the PorousFlow variable names", "href": "../../../../source/userobjects/PorousFlowDictator.html"}>>>
    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'
    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
  []
  [borehole_total_outflow_mass]
    type = PorousFlowSumQuantity<<<{"description": "Records total mass flowing into a borehole", "href": "../../../../source/userobjects/PorousFlowSumQuantity.html"}>>>
  []
  [pc]
    type = PorousFlowCapillaryPressureVG<<<{"description": "van Genuchten capillary pressure", "href": "../../../../source/userobjects/PorousFlowCapillaryPressureVG.html"}>>>
    m<<<{"description": "van Genuchten exponent m. Must be between 0 and 1, and optimally should be set to >0.5"}>>> = 0.8
    alpha<<<{"description": "van Genuchten parameter alpha. Must be positive"}>>> = 1e-5
  []
[]

[FluidProperties<<<{"href": "../../../../syntax/FluidProperties/index.html"}>>>]
  [simple_fluid]
    type = SimpleFluidProperties<<<{"description": "Fluid properties for a simple fluid with a constant bulk density", "href": "../../../../source/fluidproperties/SimpleFluidProperties.html"}>>>
    bulk_modulus<<<{"description": "Constant bulk modulus (Pa)"}>>> = 2e9
    viscosity<<<{"description": "Constant dynamic viscosity (Pa.s)"}>>> = 1e-3
    density0<<<{"description": "Density at zero pressure and zero temperature"}>>> = 1000
    thermal_expansion<<<{"description": "Constant coefficient of thermal expansion (1/K)"}>>> = 0
  []
[]

[Materials<<<{"href": "../../../../syntax/Materials/index.html"}>>>]
  [temperature]
    type = PorousFlowTemperature<<<{"description": "Material to provide temperature at the quadpoints or nodes and derivatives of it with respect to the PorousFlow variables", "href": "../../../../source/materials/PorousFlowTemperature.html"}>>>
  []
  [ppss]
    type = PorousFlow1PhaseP<<<{"description": "This Material is used for the partially saturated single-phase situation where porepressure is the primary variable", "href": "../../../../source/materials/PorousFlow1PhaseP.html"}>>>
    porepressure<<<{"description": "Variable that represents the porepressure of the single phase"}>>> = pp
    capillary_pressure<<<{"description": "Name of the UserObject defining the capillary pressure"}>>> = pc
  []
  [massfrac]
    type = PorousFlowMassFraction<<<{"description": "This Material forms a std::vector<std::vector ...> of mass-fractions out of the individual mass fractions", "href": "../../../../source/materials/PorousFlowMassFraction.html"}>>>
  []
  [simple_fluid]
    type = PorousFlowSingleComponentFluid<<<{"description": "This Material calculates fluid properties at the quadpoints or nodes for a single component fluid", "href": "../../../../source/materials/PorousFlowSingleComponentFluid.html"}>>>
    fp<<<{"description": "The name of the user object for fluid properties"}>>> = simple_fluid
    phase<<<{"description": "The phase number"}>>> = 0
  []
  [porosity]
    type = PorousFlowPorosityConst<<<{"description": "This Material calculates the porosity assuming it is constant", "href": "../../../../source/materials/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
  []
  [permeability]
    type = PorousFlowPermeabilityConst<<<{"description": "This Material calculates the permeability tensor assuming it is constant", "href": "../../../../source/materials/PorousFlowPermeabilityConst.html"}>>>
    permeability<<<{"description": "The permeability tensor (usually in m^2), which is assumed constant for this material"}>>> = '1E-12 0 0 0 1E-12 0 0 0 1E-12'
  []
  [relperm]
    type = PorousFlowRelativePermeabilityFLAC<<<{"description": "This Material calculates relative permeability of a phase using a model inspired by FLAC", "href": "../../../../source/materials/PorousFlowRelativePermeabilityFLAC.html"}>>>
    m<<<{"description": "relperm = (1 + m)seff^m - m seff^(m+1)"}>>> = 2
    phase<<<{"description": "The phase number"}>>> = 0
  []
[]

[DiracKernels<<<{"href": "../../../../syntax/DiracKernels/index.html"}>>>]
  [bh]
    type = PorousFlowPeacemanBorehole<<<{"description": "Approximates a borehole in the mesh using the Peaceman approach, ie using a number of point sinks with given radii whose positions are read from a file.  NOTE: if you are using PorousFlowPorosity that depends on volumetric strain, you should set strain_at_nearest_qp=true in your GlobalParams, to ensure the nodal Porosity Material uses the volumetric strain at the Dirac quadpoints, and can therefore be computed", "href": "../../../../source/dirackernels/PorousFlowPeacemanBorehole.html"}>>>
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = pp
    SumQuantityUO<<<{"description": "User Object of type=PorousFlowSumQuantity in which to place the total outflow from the line sink for each time step."}>>> = borehole_total_outflow_mass
    point_file<<<{"description": "The file containing the coordinates of the points and their weightings that approximate the line sink.  The physical meaning of the weightings depend on the scenario, eg, they may be borehole radii.  Each line in the file must contain a space-separated weight and coordinate, viz r x y z.  For boreholes, the last point in the file is defined as the borehole bottom, where the borehole pressure is bottom_pressure.  If your file contains just one point, you must also specify the line_length and line_direction parameters.  Note that you will get segementation faults if your points do not lie within your mesh!"}>>> = bh02.bh
    fluid_phase<<<{"description": "The fluid phase whose pressure (and potentially mobility, enthalpy, etc) controls the flux to the line sink.  For p_or_t=temperature, and without any use_*, this parameter is irrelevant"}>>> = 0
    bottom_p_or_t<<<{"description": "For function_of=pressure, this function is the pressure at the bottom of the borehole, otherwise it is the temperature at the bottom of the borehole."}>>> = -1E6
    unit_weight<<<{"description": "(fluid_density*gravitational_acceleration) as a vector pointing downwards.  Note that the borehole pressure at a given z position is bottom_p_or_t + unit_weight*(q - q_bottom), where q=(x,y,z) and q_bottom=(x,y,z) of the bottom point of the borehole.  The analogous formula holds for function_of=temperature.  If you don't want bottomhole pressure (or temperature) to vary in the borehole just set unit_weight=0.  Typical value is = (0,0,-1E4), for water"}>>> = '0 0 0'
    use_mobility<<<{"description": "Multiply the flux by the fluid mobility"}>>> = true
    character<<<{"description": "If zero then borehole does nothing.  If positive the borehole acts as a sink (production well) for porepressure > borehole pressure, and does nothing otherwise.  If negative the borehole acts as a source (injection well) for porepressure < borehole pressure, and does nothing otherwise.  The flow rate to/from the borehole is multiplied by |character|, so usually character = +/- 1, but you can specify other quantities to provide an overall scaling to the flow if you like."}>>> = 1
  []
[]

[Postprocessors<<<{"href": "../../../../syntax/Postprocessors/index.html"}>>>]
  [bh_report]
    type = PorousFlowPlotQuantity<<<{"description": "Extracts the value from the PorousFlowSumQuantity UserObject", "href": "../../../../source/postprocessors/PorousFlowPlotQuantity.html"}>>>
    uo<<<{"description": "PorousFlowSumQuantity user object name that holds the required information"}>>> = borehole_total_outflow_mass
  []

  [fluid_mass0]
    type = PorousFlowFluidMass<<<{"description": "Calculates the mass of a fluid component in a region", "href": "../../../../source/postprocessors/PorousFlowFluidMass.html"}>>>
    execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = timestep_begin
  []

  [fluid_mass1]
    type = PorousFlowFluidMass<<<{"description": "Calculates the mass of a fluid component in a region", "href": "../../../../source/postprocessors/PorousFlowFluidMass.html"}>>>
    execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = timestep_end
  []

  [zmass_error]
    type = FunctionValuePostprocessor<<<{"description": "Computes the value of a supplied function at a single point (scalable)", "href": "../../../../source/postprocessors/FunctionValuePostprocessor.html"}>>>
    function<<<{"description": "The function which supplies the postprocessor value."}>>> = mass_bal_fcn
    execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = timestep_end
    indirect_dependencies<<<{"description": "If the evaluated function depends on other postprocessors they must be listed here to ensure proper dependency resolution"}>>> = 'fluid_mass1 fluid_mass0 bh_report'
  []

  [p0]
    type = PointValue<<<{"description": "Compute the value of a variable at a specified location", "href": "../../../../source/postprocessors/PointValue.html"}>>>
    variable<<<{"description": "The name of the variable that this postprocessor operates on."}>>> = pp
    point<<<{"description": "The physical point where the solution will be evaluated."}>>> = '0 0 0'
    execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = timestep_end
  []
[]

[Functions]
  [mass_bal_fcn]
    type = ParsedFunction<<<{"description": "Function created by parsing a string", "href": "../../../../source/functions/MooseParsedFunction.html"}>>>
    expression<<<{"description": "The user defined function."}>>> = abs((a-c+d)/2/(a+c))
    symbol_names<<<{"description": "Symbols (excluding t,x,y,z) that are bound to the values provided by the corresponding items in the vals vector."}>>> = 'a c d'
    symbol_values<<<{"description": "Constant numeric values, postprocessor names, function names, and scalar variables corresponding to the symbols in symbol_names."}>>> = 'fluid_mass1 fluid_mass0 bh_report'
  []
[]

[Preconditioning<<<{"href": "../../../../syntax/Preconditioning/index.html"}>>>]
  [usual]
    type = SMP<<<{"description": "Single matrix preconditioner (SMP) builds a preconditioner using user defined off-diagonal parts of the Jacobian.", "href": "../../../../source/preconditioners/SingleMatrixPreconditioner.html"}>>>
    full<<<{"description": "Set to true if you want the full set of couplings between variables simply for convenience so you don't have to set every off_diag_row and off_diag_column combination."}>>> = true
    petsc_options<<<{"description": "Singleton PETSc options"}>>> = '-snes_converged_reason'
    petsc_options_iname<<<{"description": "Names of PETSc name/value pairs"}>>> = '-ksp_type -pc_type -snes_atol -snes_rtol -snes_max_it -ksp_max_it'
    petsc_options_value<<<{"description": "Values of PETSc name/value pairs (must correspond with \"petsc_options_iname\""}>>> = 'bcgs bjacobi 1E-10 1E-10 10000 30'
  []
[]

[Executioner<<<{"href": "../../../../syntax/Executioner/index.html"}>>>]
  type = Transient
  end_time = 1E3
  solve_type = NEWTON
  [TimeStepper<<<{"href": "../../../../syntax/Executioner/TimeStepper/index.html"}>>>]
    type = FunctionDT
    function = dts
  []
[]

[Outputs<<<{"href": "../../../../syntax/Outputs/index.html"}>>>]
  file_base<<<{"description": "Common file base name to be utilized with all output objects"}>>> = bh04
  exodus<<<{"description": "Output the results using the default settings for Exodus output."}>>> = false
  csv<<<{"description": "Output the scalar variable and postprocessors to a *.csv file using the default CSV output."}>>> = true
  execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = timestep_end
[]
(modules/porous_flow/test/tests/dirackernels/bh04.i)
  • An injection borehole with with an unsaturated medium, with use_mobility=true

# unsaturated
# injection
[Mesh<<<{"href": "../../../../syntax/Mesh/index.html"}>>>]
  type = GeneratedMesh
  dim = 3
  nx = 1
  ny = 1
  nz = 1
  xmin = -1
  xmax = 1
  ymin = -1
  ymax = 1
  zmin = -1
  zmax = 1
[]

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

[Functions<<<{"href": "../../../../syntax/Functions/index.html"}>>>]
  [dts]
    type = PiecewiseLinear<<<{"description": "Linearly interpolates between pairs of x-y data", "href": "../../../../source/functions/PiecewiseLinear.html"}>>>
    y<<<{"description": "The ordinate values"}>>> = '500 500 1E1'
    x<<<{"description": "The abscissa values"}>>> = '4000 5000 6500'
  []
[]

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

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

[UserObjects<<<{"href": "../../../../syntax/UserObjects/index.html"}>>>]
  [dictator]
    type = PorousFlowDictator<<<{"description": "Holds information on the PorousFlow variable names", "href": "../../../../source/userobjects/PorousFlowDictator.html"}>>>
    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'
    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
  []
  [borehole_total_outflow_mass]
    type = PorousFlowSumQuantity<<<{"description": "Records total mass flowing into a borehole", "href": "../../../../source/userobjects/PorousFlowSumQuantity.html"}>>>
  []
  [pc]
    type = PorousFlowCapillaryPressureVG<<<{"description": "van Genuchten capillary pressure", "href": "../../../../source/userobjects/PorousFlowCapillaryPressureVG.html"}>>>
    m<<<{"description": "van Genuchten exponent m. Must be between 0 and 1, and optimally should be set to >0.5"}>>> = 0.8
    alpha<<<{"description": "van Genuchten parameter alpha. Must be positive"}>>> = 1e-5
  []
[]

[FluidProperties<<<{"href": "../../../../syntax/FluidProperties/index.html"}>>>]
  [simple_fluid]
    type = SimpleFluidProperties<<<{"description": "Fluid properties for a simple fluid with a constant bulk density", "href": "../../../../source/fluidproperties/SimpleFluidProperties.html"}>>>
    bulk_modulus<<<{"description": "Constant bulk modulus (Pa)"}>>> = 2e9
    viscosity<<<{"description": "Constant dynamic viscosity (Pa.s)"}>>> = 1e-3
    density0<<<{"description": "Density at zero pressure and zero temperature"}>>> = 1000
    thermal_expansion<<<{"description": "Constant coefficient of thermal expansion (1/K)"}>>> = 0
  []
[]

[Materials<<<{"href": "../../../../syntax/Materials/index.html"}>>>]
  [temperature]
    type = PorousFlowTemperature<<<{"description": "Material to provide temperature at the quadpoints or nodes and derivatives of it with respect to the PorousFlow variables", "href": "../../../../source/materials/PorousFlowTemperature.html"}>>>
  []
  [ppss]
    type = PorousFlow1PhaseP<<<{"description": "This Material is used for the partially saturated single-phase situation where porepressure is the primary variable", "href": "../../../../source/materials/PorousFlow1PhaseP.html"}>>>
    porepressure<<<{"description": "Variable that represents the porepressure of the single phase"}>>> = pp
    capillary_pressure<<<{"description": "Name of the UserObject defining the capillary pressure"}>>> = pc
  []
  [massfrac]
    type = PorousFlowMassFraction<<<{"description": "This Material forms a std::vector<std::vector ...> of mass-fractions out of the individual mass fractions", "href": "../../../../source/materials/PorousFlowMassFraction.html"}>>>
  []
  [simple_fluid]
    type = PorousFlowSingleComponentFluid<<<{"description": "This Material calculates fluid properties at the quadpoints or nodes for a single component fluid", "href": "../../../../source/materials/PorousFlowSingleComponentFluid.html"}>>>
    fp<<<{"description": "The name of the user object for fluid properties"}>>> = simple_fluid
    phase<<<{"description": "The phase number"}>>> = 0
  []
  [porosity]
    type = PorousFlowPorosityConst<<<{"description": "This Material calculates the porosity assuming it is constant", "href": "../../../../source/materials/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
  []
  [permeability]
    type = PorousFlowPermeabilityConst<<<{"description": "This Material calculates the permeability tensor assuming it is constant", "href": "../../../../source/materials/PorousFlowPermeabilityConst.html"}>>>
    permeability<<<{"description": "The permeability tensor (usually in m^2), which is assumed constant for this material"}>>> = '1E-12 0 0 0 1E-12 0 0 0 1E-12'
  []
  [relperm]
    type = PorousFlowRelativePermeabilityFLAC<<<{"description": "This Material calculates relative permeability of a phase using a model inspired by FLAC", "href": "../../../../source/materials/PorousFlowRelativePermeabilityFLAC.html"}>>>
    m<<<{"description": "relperm = (1 + m)seff^m - m seff^(m+1)"}>>> = 2
    phase<<<{"description": "The phase number"}>>> = 0
  []
[]

[DiracKernels<<<{"href": "../../../../syntax/DiracKernels/index.html"}>>>]
  [bh]
    type = PorousFlowPeacemanBorehole<<<{"description": "Approximates a borehole in the mesh using the Peaceman approach, ie using a number of point sinks with given radii whose positions are read from a file.  NOTE: if you are using PorousFlowPorosity that depends on volumetric strain, you should set strain_at_nearest_qp=true in your GlobalParams, to ensure the nodal Porosity Material uses the volumetric strain at the Dirac quadpoints, and can therefore be computed", "href": "../../../../source/dirackernels/PorousFlowPeacemanBorehole.html"}>>>
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = pp
    SumQuantityUO<<<{"description": "User Object of type=PorousFlowSumQuantity in which to place the total outflow from the line sink for each time step."}>>> = borehole_total_outflow_mass
    point_file<<<{"description": "The file containing the coordinates of the points and their weightings that approximate the line sink.  The physical meaning of the weightings depend on the scenario, eg, they may be borehole radii.  Each line in the file must contain a space-separated weight and coordinate, viz r x y z.  For boreholes, the last point in the file is defined as the borehole bottom, where the borehole pressure is bottom_pressure.  If your file contains just one point, you must also specify the line_length and line_direction parameters.  Note that you will get segementation faults if your points do not lie within your mesh!"}>>> = bh03.bh
    fluid_phase<<<{"description": "The fluid phase whose pressure (and potentially mobility, enthalpy, etc) controls the flux to the line sink.  For p_or_t=temperature, and without any use_*, this parameter is irrelevant"}>>> = 0
    bottom_p_or_t<<<{"description": "For function_of=pressure, this function is the pressure at the bottom of the borehole, otherwise it is the temperature at the bottom of the borehole."}>>> = 0
    unit_weight<<<{"description": "(fluid_density*gravitational_acceleration) as a vector pointing downwards.  Note that the borehole pressure at a given z position is bottom_p_or_t + unit_weight*(q - q_bottom), where q=(x,y,z) and q_bottom=(x,y,z) of the bottom point of the borehole.  The analogous formula holds for function_of=temperature.  If you don't want bottomhole pressure (or temperature) to vary in the borehole just set unit_weight=0.  Typical value is = (0,0,-1E4), for water"}>>> = '0 0 0'
    use_mobility<<<{"description": "Multiply the flux by the fluid mobility"}>>> = true
    character<<<{"description": "If zero then borehole does nothing.  If positive the borehole acts as a sink (production well) for porepressure > borehole pressure, and does nothing otherwise.  If negative the borehole acts as a source (injection well) for porepressure < borehole pressure, and does nothing otherwise.  The flow rate to/from the borehole is multiplied by |character|, so usually character = +/- 1, but you can specify other quantities to provide an overall scaling to the flow if you like."}>>> = -1
  []
[]

[Postprocessors<<<{"href": "../../../../syntax/Postprocessors/index.html"}>>>]
  [bh_report]
    type = PorousFlowPlotQuantity<<<{"description": "Extracts the value from the PorousFlowSumQuantity UserObject", "href": "../../../../source/postprocessors/PorousFlowPlotQuantity.html"}>>>
    uo<<<{"description": "PorousFlowSumQuantity user object name that holds the required information"}>>> = borehole_total_outflow_mass
  []

  [fluid_mass0]
    type = PorousFlowFluidMass<<<{"description": "Calculates the mass of a fluid component in a region", "href": "../../../../source/postprocessors/PorousFlowFluidMass.html"}>>>
    execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = timestep_begin
  []

  [fluid_mass1]
    type = PorousFlowFluidMass<<<{"description": "Calculates the mass of a fluid component in a region", "href": "../../../../source/postprocessors/PorousFlowFluidMass.html"}>>>
    execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = timestep_end
  []

  [zmass_error]
    type = FunctionValuePostprocessor<<<{"description": "Computes the value of a supplied function at a single point (scalable)", "href": "../../../../source/postprocessors/FunctionValuePostprocessor.html"}>>>
    function<<<{"description": "The function which supplies the postprocessor value."}>>> = mass_bal_fcn
    execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = timestep_end
    indirect_dependencies<<<{"description": "If the evaluated function depends on other postprocessors they must be listed here to ensure proper dependency resolution"}>>> = 'fluid_mass1 fluid_mass0 bh_report'
  []

  [p0]
    type = PointValue<<<{"description": "Compute the value of a variable at a specified location", "href": "../../../../source/postprocessors/PointValue.html"}>>>
    variable<<<{"description": "The name of the variable that this postprocessor operates on."}>>> = pp
    point<<<{"description": "The physical point where the solution will be evaluated."}>>> = '0 0 0'
    execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = timestep_end
  []
[]

[Functions]
  [mass_bal_fcn]
    type = ParsedFunction<<<{"description": "Function created by parsing a string", "href": "../../../../source/functions/MooseParsedFunction.html"}>>>
    expression<<<{"description": "The user defined function."}>>> = abs((a-c+d)/2/(a+c))
    symbol_names<<<{"description": "Symbols (excluding t,x,y,z) that are bound to the values provided by the corresponding items in the vals vector."}>>> = 'a c d'
    symbol_values<<<{"description": "Constant numeric values, postprocessor names, function names, and scalar variables corresponding to the symbols in symbol_names."}>>> = 'fluid_mass1 fluid_mass0 bh_report'
  []
[]

[Preconditioning<<<{"href": "../../../../syntax/Preconditioning/index.html"}>>>]
  [usual]
    type = SMP<<<{"description": "Single matrix preconditioner (SMP) builds a preconditioner using user defined off-diagonal parts of the Jacobian.", "href": "../../../../source/preconditioners/SingleMatrixPreconditioner.html"}>>>
    full<<<{"description": "Set to true if you want the full set of couplings between variables simply for convenience so you don't have to set every off_diag_row and off_diag_column combination."}>>> = true
    petsc_options<<<{"description": "Singleton PETSc options"}>>> = '-snes_converged_reason'
    petsc_options_iname<<<{"description": "Names of PETSc name/value pairs"}>>> = '-ksp_type -pc_type -snes_atol -snes_rtol -snes_max_it -ksp_max_it'
    petsc_options_value<<<{"description": "Values of PETSc name/value pairs (must correspond with \"petsc_options_iname\""}>>> = 'bcgs bjacobi 1E-10 1E-10 10000 30'
  []
[]

[Executioner<<<{"href": "../../../../syntax/Executioner/index.html"}>>>]
  type = Transient
  end_time = 6500
  solve_type = NEWTON
  [TimeStepper<<<{"href": "../../../../syntax/Executioner/TimeStepper/index.html"}>>>]
    type = FunctionDT
    function = dts
  []
[]

[Outputs<<<{"href": "../../../../syntax/Outputs/index.html"}>>>]
  file_base<<<{"description": "Common file base name to be utilized with all output objects"}>>> = bh05
  exodus<<<{"description": "Output the results using the default settings for Exodus output."}>>> = false
  csv<<<{"description": "Output the scalar variable and postprocessors to a *.csv file using the default CSV output."}>>> = true
  execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = timestep_end
[]
(modules/porous_flow/test/tests/dirackernels/bh05.i)

Further commentary, and the results demonstrating that MOOSE is correct may be found in the page discussing the theory behind Dirac Kernels in PorousFlow.

Comparison with a steady-state 2D analytic solution

The test

# Comparison with analytical solution for cylindrically-symmetric situation
[Mesh<<<{"href": "../../../../syntax/Mesh/index.html"}>>>]
  type = FileMesh
  file = bh07_input.e
[]

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

[Functions<<<{"href": "../../../../syntax/Functions/index.html"}>>>]
  [dts]
    type = PiecewiseLinear<<<{"description": "Linearly interpolates between pairs of x-y data", "href": "../../../../source/functions/PiecewiseLinear.html"}>>>
    y<<<{"description": "The ordinate values"}>>> = '1000 10000'
    x<<<{"description": "The abscissa values"}>>> = '100 1000'
  []
[]

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

[Kernels<<<{"href": "../../../../syntax/Kernels/index.html"}>>>]
  [mass0]
    type = PorousFlowMassTimeDerivative<<<{"description": "Derivative of fluid-component mass with respect to time.  Mass lumping to the nodes is used.", "href": "../../../../source/kernels/PorousFlowMassTimeDerivative.html"}>>>
    fluid_component<<<{"description": "The index corresponding to the component for this kernel"}>>> = 0
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = pp
  []
  [fflux]
    type = PorousFlowAdvectiveFlux<<<{"description": "Fully-upwinded advective flux of the component given by fluid_component", "href": "../../../../source/kernels/PorousFlowAdvectiveFlux.html"}>>>
    fluid_component<<<{"description": "The index corresponding to the fluid component for this kernel"}>>> = 0
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = pp
    gravity<<<{"description": "Gravitational acceleration vector downwards (m/s^2)"}>>> = '0 0 0'
  []
[]

[BCs<<<{"href": "../../../../syntax/BCs/index.html"}>>>]
  [fix_outer]
    type = DirichletBC<<<{"description": "Imposes the essential boundary condition $u=g$, where $g$ is a constant, controllable value.", "href": "../../../../source/bcs/DirichletBC.html"}>>>
    boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = perimeter
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = pp
    value<<<{"description": "Value of the BC"}>>> = 1E7
  []
[]

[UserObjects<<<{"href": "../../../../syntax/UserObjects/index.html"}>>>]
  [dictator]
    type = PorousFlowDictator<<<{"description": "Holds information on the PorousFlow variable names", "href": "../../../../source/userobjects/PorousFlowDictator.html"}>>>
    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'
    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
  []
  [borehole_total_outflow_mass]
    type = PorousFlowSumQuantity<<<{"description": "Records total mass flowing into a borehole", "href": "../../../../source/userobjects/PorousFlowSumQuantity.html"}>>>
  []
  [pc]
    type = PorousFlowCapillaryPressureVG<<<{"description": "van Genuchten capillary pressure", "href": "../../../../source/userobjects/PorousFlowCapillaryPressureVG.html"}>>>
    m<<<{"description": "van Genuchten exponent m. Must be between 0 and 1, and optimally should be set to >0.5"}>>> = 0.8
    alpha<<<{"description": "van Genuchten parameter alpha. Must be positive"}>>> = 1e-5
  []
[]

[FluidProperties<<<{"href": "../../../../syntax/FluidProperties/index.html"}>>>]
  [simple_fluid]
    type = SimpleFluidProperties<<<{"description": "Fluid properties for a simple fluid with a constant bulk density", "href": "../../../../source/fluidproperties/SimpleFluidProperties.html"}>>>
    bulk_modulus<<<{"description": "Constant bulk modulus (Pa)"}>>> = 2e9
    viscosity<<<{"description": "Constant dynamic viscosity (Pa.s)"}>>> = 1e-3
    density0<<<{"description": "Density at zero pressure and zero temperature"}>>> = 1000
    thermal_expansion<<<{"description": "Constant coefficient of thermal expansion (1/K)"}>>> = 0
  []
[]

[Materials<<<{"href": "../../../../syntax/Materials/index.html"}>>>]
  [temperature]
    type = PorousFlowTemperature<<<{"description": "Material to provide temperature at the quadpoints or nodes and derivatives of it with respect to the PorousFlow variables", "href": "../../../../source/materials/PorousFlowTemperature.html"}>>>
  []
  [ppss]
    type = PorousFlow1PhaseP<<<{"description": "This Material is used for the partially saturated single-phase situation where porepressure is the primary variable", "href": "../../../../source/materials/PorousFlow1PhaseP.html"}>>>
    porepressure<<<{"description": "Variable that represents the porepressure of the single phase"}>>> = pp
    capillary_pressure<<<{"description": "Name of the UserObject defining the capillary pressure"}>>> = pc
  []
  [massfrac]
    type = PorousFlowMassFraction<<<{"description": "This Material forms a std::vector<std::vector ...> of mass-fractions out of the individual mass fractions", "href": "../../../../source/materials/PorousFlowMassFraction.html"}>>>
  []
  [simple_fluid]
    type = PorousFlowSingleComponentFluid<<<{"description": "This Material calculates fluid properties at the quadpoints or nodes for a single component fluid", "href": "../../../../source/materials/PorousFlowSingleComponentFluid.html"}>>>
    fp<<<{"description": "The name of the user object for fluid properties"}>>> = simple_fluid
    phase<<<{"description": "The phase number"}>>> = 0
  []
  [porosity]
    type = PorousFlowPorosityConst<<<{"description": "This Material calculates the porosity assuming it is constant", "href": "../../../../source/materials/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
  []
  [permeability]
    type = PorousFlowPermeabilityConst<<<{"description": "This Material calculates the permeability tensor assuming it is constant", "href": "../../../../source/materials/PorousFlowPermeabilityConst.html"}>>>
    permeability<<<{"description": "The permeability tensor (usually in m^2), which is assumed constant for this material"}>>> = '1E-11 0 0 0 1E-11 0 0 0 1E-11'
  []
  [relperm]
    type = PorousFlowRelativePermeabilityFLAC<<<{"description": "This Material calculates relative permeability of a phase using a model inspired by FLAC", "href": "../../../../source/materials/PorousFlowRelativePermeabilityFLAC.html"}>>>
    m<<<{"description": "relperm = (1 + m)seff^m - m seff^(m+1)"}>>> = 2
    phase<<<{"description": "The phase number"}>>> = 0
  []
[]

[DiracKernels<<<{"href": "../../../../syntax/DiracKernels/index.html"}>>>]
  [bh]
    type = PorousFlowPeacemanBorehole<<<{"description": "Approximates a borehole in the mesh using the Peaceman approach, ie using a number of point sinks with given radii whose positions are read from a file.  NOTE: if you are using PorousFlowPorosity that depends on volumetric strain, you should set strain_at_nearest_qp=true in your GlobalParams, to ensure the nodal Porosity Material uses the volumetric strain at the Dirac quadpoints, and can therefore be computed", "href": "../../../../source/dirackernels/PorousFlowPeacemanBorehole.html"}>>>
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = pp
    SumQuantityUO<<<{"description": "User Object of type=PorousFlowSumQuantity in which to place the total outflow from the line sink for each time step."}>>> = borehole_total_outflow_mass
    point_file<<<{"description": "The file containing the coordinates of the points and their weightings that approximate the line sink.  The physical meaning of the weightings depend on the scenario, eg, they may be borehole radii.  Each line in the file must contain a space-separated weight and coordinate, viz r x y z.  For boreholes, the last point in the file is defined as the borehole bottom, where the borehole pressure is bottom_pressure.  If your file contains just one point, you must also specify the line_length and line_direction parameters.  Note that you will get segementation faults if your points do not lie within your mesh!"}>>> = bh07.bh
    fluid_phase<<<{"description": "The fluid phase whose pressure (and potentially mobility, enthalpy, etc) controls the flux to the line sink.  For p_or_t=temperature, and without any use_*, this parameter is irrelevant"}>>> = 0
    bottom_p_or_t<<<{"description": "For function_of=pressure, this function is the pressure at the bottom of the borehole, otherwise it is the temperature at the bottom of the borehole."}>>> = 0
    unit_weight<<<{"description": "(fluid_density*gravitational_acceleration) as a vector pointing downwards.  Note that the borehole pressure at a given z position is bottom_p_or_t + unit_weight*(q - q_bottom), where q=(x,y,z) and q_bottom=(x,y,z) of the bottom point of the borehole.  The analogous formula holds for function_of=temperature.  If you don't want bottomhole pressure (or temperature) to vary in the borehole just set unit_weight=0.  Typical value is = (0,0,-1E4), for water"}>>> = '0 0 0'
    use_mobility<<<{"description": "Multiply the flux by the fluid mobility"}>>> = true
    re_constant<<<{"description": "The dimensionless constant used in evaluating the borehole effective radius.  This depends on the meshing scheme.  Peacemann finite-difference calculations give 0.28, while for rectangular finite elements the result is closer to 0.1594.  (See  Eqn(4.13) of Z Chen, Y Zhang, Well flow models for various numerical methods, Int J Num Analysis and Modeling, 3 (2008) 375-388.)"}>>> = 0.1594  # use Chen and Zhang version
    character<<<{"description": "If zero then borehole does nothing.  If positive the borehole acts as a sink (production well) for porepressure > borehole pressure, and does nothing otherwise.  If negative the borehole acts as a source (injection well) for porepressure < borehole pressure, and does nothing otherwise.  The flow rate to/from the borehole is multiplied by |character|, so usually character = +/- 1, but you can specify other quantities to provide an overall scaling to the flow if you like."}>>> = 2 # double the strength because bh07.bh only fills half the mesh
  []
[]

[Postprocessors<<<{"href": "../../../../syntax/Postprocessors/index.html"}>>>]
  [bh_report]
    type = PorousFlowPlotQuantity<<<{"description": "Extracts the value from the PorousFlowSumQuantity UserObject", "href": "../../../../source/postprocessors/PorousFlowPlotQuantity.html"}>>>
    uo<<<{"description": "PorousFlowSumQuantity user object name that holds the required information"}>>> = borehole_total_outflow_mass
    execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = 'initial timestep_end'
  []
  [fluid_mass]
    type = PorousFlowFluidMass<<<{"description": "Calculates the mass of a fluid component in a region", "href": "../../../../source/postprocessors/PorousFlowFluidMass.html"}>>>
    execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = 'initial timestep_end'
  []
[]

[VectorPostprocessors<<<{"href": "../../../../syntax/VectorPostprocessors/index.html"}>>>]
  [pp]
    type = LineValueSampler<<<{"description": "Samples variable(s) along a specified line", "href": "../../../../source/vectorpostprocessors/LineValueSampler.html"}>>>
    variable<<<{"description": "The names of the variables that this VectorPostprocessor operates on"}>>> = pp
    start_point<<<{"description": "The beginning of the line"}>>> = '0 0 0'
    end_point<<<{"description": "The ending of the line"}>>> = '300 0 0'
    sort_by<<<{"description": "What to sort the samples by"}>>> = x
    num_points<<<{"description": "The number of points to sample along the line"}>>> = 300
    execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = timestep_end
  []
[]

[Preconditioning<<<{"href": "../../../../syntax/Preconditioning/index.html"}>>>]
  [usual]
    type = SMP<<<{"description": "Single matrix preconditioner (SMP) builds a preconditioner using user defined off-diagonal parts of the Jacobian.", "href": "../../../../source/preconditioners/SingleMatrixPreconditioner.html"}>>>
    full<<<{"description": "Set to true if you want the full set of couplings between variables simply for convenience so you don't have to set every off_diag_row and off_diag_column combination."}>>> = true
    petsc_options<<<{"description": "Singleton PETSc options"}>>> = '-snes_converged_reason'
    petsc_options_iname<<<{"description": "Names of PETSc name/value pairs"}>>> = '-ksp_type -pc_type -snes_atol -snes_rtol -snes_max_it -ksp_max_it'
    petsc_options_value<<<{"description": "Values of PETSc name/value pairs (must correspond with \"petsc_options_iname\""}>>> = 'bcgs bjacobi 1E-10 1E-10 10000 30'
  []
[]

[Executioner<<<{"href": "../../../../syntax/Executioner/index.html"}>>>]
  type = Transient
  end_time = 1E3
  solve_type = NEWTON
  [TimeStepper<<<{"href": "../../../../syntax/Executioner/TimeStepper/index.html"}>>>]
    # get only marginally better results for smaller time steps
    type = FunctionDT
    function = dts
  []
[]

[Outputs<<<{"href": "../../../../syntax/Outputs/index.html"}>>>]
  file_base<<<{"description": "Common file base name to be utilized with all output objects"}>>> = bh07
  [along_line]
    type = CSV<<<{"description": "Output for postprocessors, vector postprocessors, and scalar variables using comma seperated values (CSV).", "href": "../../../../source/outputs/CSV.html"}>>>
    execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = final
  []
  [exodus]
    type = Exodus<<<{"description": "Object for output data in the Exodus format", "href": "../../../../source/outputs/Exodus.html"}>>>
    execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = 'initial final'
  []
[]
(modules/porous_flow/test/tests/dirackernels/bh07.i)

checks whether the Peaceman borehole extracts fluid at the correct rate and the corresponding reduction in porepressure agrees with an analytic solution of Laplace's equation. Details and results demonstrating the correctness of PorousFlow can be found in the page discussing the theory behind Dirac Kernels in PorousFlow.