MethaneFluidProperties

Fluid properties for methane (CH4)

Fluid properties for methane are mainly calculated using the Setzmann and Wagner equation of state (Setzmann and Wagner, 1991). This formulation uses density and temperature as the primary variables with which to calculate properties such as density, enthalpy and internal energy.

When used with the pressure and temperature interface, which is the case in the Porous Flow module, methane properties are typically calculated by first calculating density iteratively for a given pressure and temperature. This density is then used to calculate the other properties, such as internal energy, directly. The computational expense associated with the iterative calculation can be mitigated using TabulatedFluidProperties.

For low pressures (typically less than 10 MPa), the properties of methane can be approximated using an ideal gas, which are much faster to calculate. However, at higher pressures, this approximation can lead to large differences, see Figure 1.

Figure 1: Methane density at 350K for various pressures.

Transport properties such as viscosity and thermal conductivity are calculated using the formulations provided in Irvine Jr and Liley (1984).

Dissolution of methane into water is calculated using Henry's law (IAPWS, 2004).

Properties of methane

Property	value
Molar mass	0.0160425 kg/mol
Critical temperature	190.564 K
Critical pressure	4.5992 MPa
Critical density	162.66 kg/m $var element = document.getElementById("moose-equation-4caa8020-3b85-48c8-94c0-25da97f8a40f");katex.render("^3", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$
Triple point temperature	90.6941 K
Triple point pressure	0.01169 MPa

Range of validity

The MethaneFluidProperties UserObject is valid for:

90.69 K $var element = document.getElementById("moose-equation-9134df84-3c6f-40e7-96cc-af069ef001b6");katex.render("\\le", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ T $var element = document.getElementById("moose-equation-a1cc70b8-e529-43c8-91ed-0f863cb1b7c2");katex.render("\\le", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ 625 K

and pressures up to 100 MPa.

Input Parameters

T_initial_guess400Temperature initial guess for Newton Method variable set conversion
Default:400
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:Temperature initial guess for Newton Method variable set conversion
max_newton_its100Maximum number of Newton iterations for variable set conversions
Default:100
C++ Type:unsigned int
Controllable:No
Description:Maximum number of Newton iterations for variable set conversions
p_initial_guess200000Pressure initial guess for Newton Method variable set conversion
Default:200000
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:Pressure initial guess for Newton Method variable set conversion
tolerance1e-08Tolerance for 2D Newton variable set conversion
Default:1e-08
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:Tolerance for 2D Newton variable set conversion

Variable Set Conversions Newton Solve Parameters

allow_imperfect_jacobiansFalsetrue to allow unimplemented property derivative terms to be set to zero for the AD API
Default:False
C++ Type:bool
Controllable:No
Description:true to allow unimplemented property derivative terms to be set to zero for the AD API
control_tagsAdds user-defined labels for accessing object parameters via control logic.
C++ Type:std::vector<std::string>
Controllable:No
Description:Adds user-defined labels for accessing object parameters via control logic.
enableTrueSet the enabled status of the MooseObject.
Default:True
C++ Type:bool
Controllable:Yes
Description:Set the enabled status of the MooseObject.
fp_typesingle-phase-fpType of the fluid property object
Default:single-phase-fp
C++ Type:FPType
Controllable:No
Description:Type of the fluid property object

Advanced Parameters

prop_getter_suffixAn optional suffix parameter that can be appended to any attempt to retrieve/get material properties. The suffix will be prepended with a '_' character.
C++ Type:MaterialPropertyName
Unit:(no unit assumed)
Controllable:No
Description:An optional suffix parameter that can be appended to any attempt to retrieve/get material properties. The suffix will be prepended with a '_' character.
use_interpolated_stateFalseFor the old and older state use projected material properties interpolated at the quadrature points. To set up projection use the ProjectedStatefulMaterialStorageAction.
Default:False
C++ Type:bool
Controllable:No
Description:For the old and older state use projected material properties interpolated at the quadrature points. To set up projection use the ProjectedStatefulMaterialStorageAction.

Material Property Retrieval Parameters

Input Files

(modules/porous_flow/examples/restart/gas_injection_new_mesh.i)
(modules/porous_flow/examples/restart/gas_injection.i)
(modules/porous_flow/test/tests/jacobian/basic_advection4.i)
(modules/porous_flow/test/tests/fluids/methane.i)
(modules/porous_flow/test/tests/fluidstate/theis_nonisothermal.i)
(modules/fluid_properties/test/tests/methane/methane.i)

References

IAPWS. Guidelines on the Henry's constant and vapour liquid distribution constant for gases in H$_2$O and D$_2$O at high temperatures. Technical Report, IAPWS, 2004.

@techreport{iapws2004,
    author = "IAPWS",
    title = "{Guidelines on the Henry's constant and vapour liquid distribution constant for gases in H$\_2$O and D$\_2$O at high temperatures}",
    institution = "IAPWS",
    year = "2004"
}

T. F. Irvine Jr and P. E. Liley. Steam and gas tables with computer equations. Academic Press Inc., Orlando, 1984.

@book{irvine1984,
    author = "Irvine Jr, T. F. and Liley, P. E.",
    title = "Steam and gas tables with computer equations",
    publisher = "Academic Press Inc.",
    address = "Orlando",
    year = "1984"
}

Ulrich Setzmann and Wolfgang Wagner. A new equation of state and tables of thermodynamic properties for methane covering the range from the melting line to 625 k at pressures up to 100 mpa. J. Phys. Chem. Ref. Data, 20(6):1061–1155, 1991.

@article{setzmann1991,
    author = "Setzmann, Ulrich and Wagner, Wolfgang",
    title = "A new equation of state and tables of thermodynamic properties for methane covering the range from the melting line to 625 K at pressures up to 100 MPa",
    journal = "J. Phys. Chem. Ref. Data",
    volume = "20",
    number = "6",
    pages = "1061--1155",
    year = "1991"
}

(modules/porous_flow/examples/restart/gas_injection_new_mesh.i)

# Using the results from the equilibrium run to provide the initial condition for
# porepressure, we now inject a gas phase into the brine-saturated reservoir. In this
# example, the mesh is not identical to the mesh used in gravityeq.i. Rather, it is
# generated so that it is more refined near the injection boundary and at the top of
# the model, as that is where the gas plume will be present.
#
# To use the hydrostatic pressure calculated using the gravity equilibrium run as the initial
# condition for the pressure, a SolutionUserObject is used, along with a SolutionFunction to
# interpolate the pressure from the gravity equilibrium run to the initial condition for liqiud
# porepressure in this example.
#
# Even though the gravity equilibrium is established using a 2D mesh, in this example,
# we use a mesh shifted 0.1 m to the right and rotate it about the Y axis to make a 2D radial
# model.
#
# Methane injection takes place over the surface of the hole created by rotating the mesh,
# and hence the injection area is 2 pi r h. We can calculate this using an AreaPostprocessor,
# and then use this in a ParsedFunction to calculate the injection rate so that 10 kg/s of
# methane is injected.
#
# Note: as this example uses the results from a previous simulation, gravityeq.i MUST be
# run before running this input file.

[Mesh]
  type = GeneratedMesh
  dim = 2
  ny = 25
  nx = 50
  ymax = 100
  xmin = 0.1
  xmax = 5000
  bias_x = 1.05
  bias_y = 0.95
  coord_type = RZ
  rz_coord_axis = Y
[]

[GlobalParams]
  PorousFlowDictator = dictator
  gravity = '0 -9.81 0'
  temperature_unit = Celsius
[]

[Variables]
  [pp_liq]
  []
  [sat_gas]
    initial_condition = 0
  []
[]

[ICs]
  [ppliq_ic]
    type = FunctionIC
    variable = pp_liq
    function = ppliq_ic
  []
[]

[AuxVariables]
  [temperature]
    initial_condition = 50
  []
  [xnacl]
    initial_condition = 0.1
  []
  [brine_density]
    family = MONOMIAL
    order = CONSTANT
  []
  [methane_density]
    family = MONOMIAL
    order = CONSTANT
  []
  [massfrac_ph0_sp0]
    initial_condition = 1
  []
  [massfrac_ph1_sp0]
    initial_condition = 0
  []
  [pp_gas]
    family = MONOMIAL
    order = CONSTANT
  []
  [sat_liq]
    family = MONOMIAL
    order = CONSTANT
  []
[]

[Kernels]
  [mass0]
    type = PorousFlowMassTimeDerivative
    variable = pp_liq
  []
  [flux0]
    type = PorousFlowAdvectiveFlux
    variable = pp_liq
  []
  [mass1]
    type = PorousFlowMassTimeDerivative
    variable = sat_gas
    fluid_component = 1
  []
  [flux1]
    type = PorousFlowAdvectiveFlux
    variable = sat_gas
    fluid_component = 1
  []
[]

[AuxKernels]
  [brine_density]
    type = PorousFlowPropertyAux
    property = density
    variable = brine_density
    execute_on = 'initial timestep_end'
  []
  [methane_density]
    type = PorousFlowPropertyAux
    property = density
    variable = methane_density
    phase = 1
    execute_on = 'initial timestep_end'
  []
  [pp_gas]
    type = PorousFlowPropertyAux
    property = pressure
    phase = 1
    variable = pp_gas
    execute_on = 'initial timestep_end'
  []
  [sat_liq]
    type = PorousFlowPropertyAux
    property = saturation
    variable = sat_liq
    execute_on = 'initial timestep_end'
  []
[]

[BCs]
  [gas_injection]
    type = PorousFlowSink
    boundary = left
    variable = sat_gas
    flux_function = injection_rate
    fluid_phase = 1
  []
  [brine_out]
    type = PorousFlowPiecewiseLinearSink
    boundary = right
    variable = pp_liq
    multipliers = '0 1e9'
    pt_vals = '0 1e9'
    fluid_phase = 0
    flux_function = 1e-6
    use_mobility = true
    use_relperm = true
    mass_fraction_component = 0
  []
[]

[Functions]
  [injection_rate]
    type = ParsedFunction
    symbol_values = injection_area
    symbol_names = area
    expression = '-1/area'
  []
  [ppliq_ic]
    type = SolutionFunction
    solution = soln
  []
[]

[UserObjects]
  [dictator]
    type = PorousFlowDictator
    porous_flow_vars = 'pp_liq sat_gas'
    number_fluid_phases = 2
    number_fluid_components = 2
  []
  [pc]
    type = PorousFlowCapillaryPressureVG
    alpha = 1e-5
    m = 0.5
    sat_lr = 0.2
    pc_max = 1e7
  []
  [soln]
    type = SolutionUserObject
    mesh = gravityeq_out.e
    system_variables = porepressure
  []
[]

[FluidProperties]
  [brine]
    type = BrineFluidProperties
  []
  [methane]
    type = MethaneFluidProperties
  []
  [methane_tab]
    type = TabulatedBicubicFluidProperties
    fp = methane
    save_file = false
  []
[]

[Materials]
  [temperature]
    type = PorousFlowTemperature
    temperature = temperature
  []
  [ps]
    type = PorousFlow2PhasePS
    phase0_porepressure = pp_liq
    phase1_saturation = sat_gas
    capillary_pressure = pc
  []
  [massfrac]
    type = PorousFlowMassFraction
    mass_fraction_vars = 'massfrac_ph0_sp0 massfrac_ph1_sp0'
  []
  [brine]
    type = PorousFlowBrine
    compute_enthalpy = false
    compute_internal_energy = false
    xnacl = xnacl
    phase = 0
  []
  [methane]
    type = PorousFlowSingleComponentFluid
    compute_enthalpy = false
    compute_internal_energy = false
    fp = methane_tab
    phase = 1
  []
  [porosity]
    type = PorousFlowPorosityConst
    porosity = 0.1
  []
  [permeability]
    type = PorousFlowPermeabilityConst
    permeability = '1e-13 0 0 0 5e-14 0  0 0 1e-13'
  []
  [relperm_liq]
    type = PorousFlowRelativePermeabilityCorey
    n = 2
    phase = 0
    s_res = 0.2
    sum_s_res = 0.3
  []
  [relperm_gas]
    type = PorousFlowRelativePermeabilityCorey
    n = 2
    phase = 1
    s_res = 0.1
    sum_s_res = 0.3
  []
[]

[Preconditioning]
  [smp]
    type = SMP
    full = true
    petsc_options_iname = '-pc_type -sub_pc_type -sub_pc_factor_shift_type'
    petsc_options_value = ' asm      lu           NONZERO'
  []
[]

[Executioner]
  type = Transient
  solve_type = Newton
  end_time = 1e8
  nl_abs_tol = 1e-12
  nl_rel_tol = 1e-06
  nl_max_its = 20
  dtmax = 1e6
  [TimeStepper]
    type = IterationAdaptiveDT
    dt = 1e1
    growth_factor = 1.5
  []
[]

[Postprocessors]
  [mass_ph0]
    type = PorousFlowFluidMass
    fluid_component = 0
    execute_on = 'initial timestep_end'
  []
  [mass_ph1]
    type = PorousFlowFluidMass
    fluid_component = 1
    execute_on = 'initial timestep_end'
  []
  [injection_area]
    type = AreaPostprocessor
    boundary = left
    execute_on = initial
  []
[]

[Outputs]
  execute_on = 'initial timestep_end'
  exodus = true
  perf_graph = true
[]

(modules/porous_flow/examples/restart/gas_injection.i)

# Using the results from the equilibrium run to provide the initial condition for
# porepressure, we now inject a gas phase into the brine-saturated reservoir. In this
# example, where the mesh used is identical to the mesh used in gravityeq.i, we can use
# the basic restart capability by simply setting the initial condition for porepressure
# using the results from gravityeq.i.
#
# Even though the gravity equilibrium is established using a 2D mesh, in this example,
# we shift the mesh 0.1 m to the right and rotate it about the Y axis to make a 2D radial
# model.
#
# Methane injection takes place over the surface of the hole created by rotating the mesh,
# and hence the injection area is 2 pi r h. We can calculate this using an AreaPostprocessor,
# and then use this in a ParsedFunction to calculate the injection rate so that 10 kg/s of
# methane is injected.
#
# Results can be improved by uniformly refining the initial mesh.
#
# Note: as this example uses the results from a previous simulation, gravityeq.i MUST be
# run before running this input file.

[Mesh]
  uniform_refine = 1
  [file]
    type = FileMeshGenerator
    file = gravityeq_out.e
  []
  [translate]
    type = TransformGenerator
    transform = TRANSLATE
    vector_value = '0.1 0 0'
    input = file
  []
  coord_type = RZ
  rz_coord_axis = Y
[]

[GlobalParams]
  PorousFlowDictator = dictator
  gravity = '0 -9.81 0'
  temperature_unit = Celsius
[]

[Variables]
  [pp_liq]
    initial_from_file_var = porepressure
  []
  [sat_gas]
    initial_condition = 0
  []
[]

[AuxVariables]
  [temperature]
    initial_condition = 50
  []
  [xnacl]
    initial_condition = 0.1
  []
  [brine_density]
    family = MONOMIAL
    order = CONSTANT
  []
  [methane_density]
    family = MONOMIAL
    order = CONSTANT
  []
  [massfrac_ph0_sp0]
    initial_condition = 1
  []
  [massfrac_ph1_sp0]
    initial_condition = 0
  []
  [pp_gas]
    family = MONOMIAL
    order = CONSTANT
  []
  [sat_liq]
    family = MONOMIAL
    order = CONSTANT
  []
[]

[Kernels]
  [mass0]
    type = PorousFlowMassTimeDerivative
    variable = pp_liq
  []
  [flux0]
    type = PorousFlowAdvectiveFlux
    variable = pp_liq
  []
  [mass1]
    type = PorousFlowMassTimeDerivative
    variable = sat_gas
    fluid_component = 1
  []
  [flux1]
    type = PorousFlowAdvectiveFlux
    variable = sat_gas
    fluid_component = 1
  []
[]

[AuxKernels]
  [brine_density]
    type = PorousFlowPropertyAux
    property = density
    variable = brine_density
    execute_on = 'initial timestep_end'
  []
  [methane_density]
    type = PorousFlowPropertyAux
    property = density
    variable = methane_density
    phase = 1
    execute_on = 'initial timestep_end'
  []
  [pp_gas]
    type = PorousFlowPropertyAux
    property = pressure
    phase = 1
    variable = pp_gas
    execute_on = 'initial timestep_end'
  []
  [sat_liq]
    type = PorousFlowPropertyAux
    property = saturation
    variable = sat_liq
    execute_on = 'initial timestep_end'
  []
[]

[BCs]
  [gas_injection]
    type = PorousFlowSink
    boundary = left
    variable = sat_gas
    flux_function = injection_rate
    fluid_phase = 1
  []
  [brine_out]
    type = PorousFlowPiecewiseLinearSink
    boundary = right
    variable = pp_liq
    multipliers = '0 1e9'
    pt_vals = '0 1e9'
    fluid_phase = 0
    flux_function = 1e-6
    use_mobility = true
  []
[]

[Functions]
  [injection_rate]
    type = ParsedFunction
    symbol_values = injection_area
    symbol_names = area
    expression = '-10/area'
  []
[]

[UserObjects]
  [dictator]
    type = PorousFlowDictator
    porous_flow_vars = 'pp_liq sat_gas'
    number_fluid_phases = 2
    number_fluid_components = 2
  []
  [pc]
    type = PorousFlowCapillaryPressureVG
    alpha = 1e-5
    m = 0.5
    sat_lr = 0.2
  []
[]

[FluidProperties]
  [brine]
    type = BrineFluidProperties
  []
  [methane]
    type = MethaneFluidProperties
  []
  [methane_tab]
    type = TabulatedBicubicFluidProperties
    fp = methane
    save_file = false
  []
[]

[Materials]
  [temperature]
    type = PorousFlowTemperature
    temperature = temperature
  []
  [ps]
    type = PorousFlow2PhasePS
    phase0_porepressure = pp_liq
    phase1_saturation = sat_gas
    capillary_pressure = pc
  []
  [massfrac]
    type = PorousFlowMassFraction
    mass_fraction_vars = 'massfrac_ph0_sp0 massfrac_ph1_sp0'
  []
  [brine]
    type = PorousFlowBrine
    compute_enthalpy = false
    compute_internal_energy = false
    xnacl = xnacl
    phase = 0
  []
  [methane]
    type = PorousFlowSingleComponentFluid
    compute_enthalpy = false
    compute_internal_energy = false
    fp = methane_tab
    phase = 1
  []
  [porosity]
    type = PorousFlowPorosityConst
    porosity = 0.1
  []
  [permeability]
    type = PorousFlowPermeabilityConst
    permeability = '1e-13 0 0 0 1e-13 0  0 0 1e-13'
  []
  [relperm_liq]
    type = PorousFlowRelativePermeabilityCorey
    n = 2
    phase = 0
    s_res = 0.2
    sum_s_res = 0.3
  []
  [relperm_gas]
    type = PorousFlowRelativePermeabilityCorey
    n = 2
    phase = 1
    s_res = 0.1
    sum_s_res = 0.3
  []
[]

[Preconditioning]
  [smp]
    type = SMP
    full = true
    petsc_options_iname = '-pc_type -sub_pc_type -sub_pc_factor_shift_type'
    petsc_options_value = ' asm      lu           NONZERO'
  []
[]

[Executioner]
  type = Transient
  solve_type = Newton
  end_time = 1e8
  nl_abs_tol = 1e-12
  nl_rel_tol = 1e-06
  nl_max_its = 20
  dtmax = 1e6
  [TimeStepper]
    type = IterationAdaptiveDT
    dt = 1e1
  []
[]

[Postprocessors]
  [mass_ph0]
    type = PorousFlowFluidMass
    fluid_component = 0
    execute_on = 'initial timestep_end'
  []
  [mass_ph1]
    type = PorousFlowFluidMass
    fluid_component = 1
    execute_on = 'initial timestep_end'
  []
  [injection_area]
    type = AreaPostprocessor
    boundary = left
    execute_on = initial
  []
[]

[Outputs]
  execute_on = 'initial timestep_end'
  exodus = true
  perf_graph = true
  checkpoint = true
[]

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

# Basic advection with 1 porepressure and temperature as PorousFlow variables
# Constant permeability
[Mesh]
  type = GeneratedMesh
  dim = 1
  nx = 1
[]

[GlobalParams]
  PorousFlowDictator = dictator
[]

[Variables]
  [u]
  []
  [T]
  []
  [P]
  []
[]

[ICs]
  [P]
    type = RandomIC
    variable = P
    min = 2E5
    max = 4E5
  []
  [T]
    type = RandomIC
    variable = T
    min = 300
    max = 900
  []
  [u]
    type = RandomIC
    variable = u
  []
[]

[Kernels]
  [dummy_T]
    type = NullKernel
    variable = T
  []
  [dummy_P]
    type = NullKernel
    variable = P
  []
  [u_advection]
    type = PorousFlowBasicAdvection
    variable = u
    phase = 0
  []
[]

[UserObjects]
  [dictator]
    type = PorousFlowDictator
    porous_flow_vars = 'P T'
    number_fluid_phases = 1
    number_fluid_components = 1
  []
  [pc]
    type = PorousFlowCapillaryPressureVG
    alpha = 1E-5
    m = 0.6
    sat_lr = 0.1
  []
[]

[FluidProperties]
  [methane]
    type = MethaneFluidProperties
  []
[]

[Materials]
  [temperature_qp]
    type = PorousFlowTemperature
    temperature = T
  []
  [ppss_qp]
    type = PorousFlow1PhaseP
    porepressure = P
    capillary_pressure = pc
  []
  [fluid_qp]
    type = PorousFlowSingleComponentFluid
    fp = methane
    phase = 0
  []
  [permeability]
    type = PorousFlowPermeabilityConst
    permeability = '5 0 0 0 5 0 0 0 5'
  []
  [relperm_qp]
    type = PorousFlowRelativePermeabilityCorey
    n = 2
    phase = 0
    s_res = 0.1
    sum_s_res = 0.1
  []
  [darcy_velocity_qp]
    type = PorousFlowDarcyVelocityMaterial
    gravity = '0.25 0 0'
  []
[]

[Preconditioning]
  [check]
    type = SMP
    full = true
    #petsc_options = '-snes_test_display'
    petsc_options_iname = '-snes_type'
    petsc_options_value = ' test'
  []
[]

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

(modules/porous_flow/test/tests/fluids/methane.i)

# Test MethaneFluidProperties
# Reference data from Irvine Jr, T. F. and Liley, P. E. (1984) Steam and
# Gas Tables with Computer Equations
#
# For temperature = 350K, the fluid properties should be:
# density = 55.13 kg/m^3
# viscosity = 0.01276 mPa.s
# h = 708.5 kJ/kg

[Mesh]
  type = GeneratedMesh
  dim = 1
  nx = 1
[]

[GlobalParams]
  PorousFlowDictator = dictator
[]

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

[Variables]
  [pp]
    initial_condition = 10e6
  []
[]

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

[AuxVariables]
  [temp]
    initial_condition = 350.0
  []
[]

[Materials]
  [temperature]
    type = PorousFlowTemperature
    temperature = 'temp'
  []
  [ppss]
    type = PorousFlow1PhaseFullySaturated
    porepressure = pp
  []
  [methane]
    type = PorousFlowSingleComponentFluid
    temperature_unit = Kelvin
    fp = methane
    phase = 0
  []
[]

[FluidProperties]
  [methane]
    type = MethaneFluidProperties
  []
[]

[Executioner]
  type = Steady
  solve_type = Newton
[]

[Postprocessors]
  [pressure]
    type = ElementIntegralVariablePostprocessor
    variable = pp
  []
  [temperature]
    type = ElementIntegralVariablePostprocessor
    variable = temp
  []
  [density]
    type = ElementIntegralMaterialProperty
    mat_prop = 'PorousFlow_fluid_phase_density_qp0'
  []
  [viscosity]
    type = ElementIntegralMaterialProperty
    mat_prop = 'PorousFlow_viscosity_qp0'
  []
  [enthalpy]
    type = ElementIntegralMaterialProperty
    mat_prop = 'PorousFlow_fluid_phase_enthalpy_qp0'
  []
[]

[Outputs]
  execute_on = 'timestep_end'
  file_base = methane
  csv = true
[]

(modules/porous_flow/test/tests/fluidstate/theis_nonisothermal.i)

# Two-phase nonisothermal Theis problem: Flow from single source using WaterNCG fluidstate.
# Constant rate injection 2 kg/s of cold gas into warm reservoir
# 1D cylindrical mesh
# Initially, system has only a liquid phase, until enough gas is injected
# to form a gas phase, in which case the system becomes two phase.

[Mesh]
  [mesh]
    type = GeneratedMeshGenerator
    dim = 1
    nx = 40
    xmin = 0.1
    xmax = 200
    bias_x = 1.05
  []
  coord_type = RZ
  rz_coord_axis = Y
[]

[GlobalParams]
  PorousFlowDictator = dictator
  gravity = '0 0 0'
[]

[AuxVariables]
  [saturation_gas]
    order = CONSTANT
    family = MONOMIAL
  []
  [x1]
    order = CONSTANT
    family = MONOMIAL
  []
  [y0]
    order = CONSTANT
    family = MONOMIAL
  []
[]

[AuxKernels]
  [saturation_gas]
    type = PorousFlowPropertyAux
    variable = saturation_gas
    property = saturation
    phase = 1
    execute_on = timestep_end
  []
  [x1]
    type = PorousFlowPropertyAux
    variable = x1
    property = mass_fraction
    phase = 0
    fluid_component = 1
    execute_on = timestep_end
  []
  [y0]
    type = PorousFlowPropertyAux
    variable = y0
    property = mass_fraction
    phase = 1
    fluid_component = 0
    execute_on = timestep_end
  []
[]

[Variables]
  [pgas]
    initial_condition = 20e6
  []
  [zi]
    initial_condition = 0
  []
  [temperature]
    initial_condition = 70
    scaling = 1e-4
  []
[]

[Kernels]
  [mass0]
    type = PorousFlowMassTimeDerivative
    fluid_component = 0
    variable = pgas
  []
  [flux0]
    type = PorousFlowAdvectiveFlux
    fluid_component = 0
    variable = pgas
  []
  [mass1]
    type = PorousFlowMassTimeDerivative
    fluid_component = 1
    variable = zi
  []
  [flux1]
    type = PorousFlowAdvectiveFlux
    fluid_component = 1
    variable = zi
  []
  [energy]
    type = PorousFlowEnergyTimeDerivative
    variable = temperature
  []
  [heatadv]
    type = PorousFlowHeatAdvection
    variable = temperature
  []
  [conduction]
    type = PorousFlowHeatConduction
    variable = temperature
  []
[]

[UserObjects]
  [dictator]
    type = PorousFlowDictator
    porous_flow_vars = 'pgas zi temperature'
    number_fluid_phases = 2
    number_fluid_components = 2
  []
  [pc]
    type = PorousFlowCapillaryPressureConst
    pc = 0
  []
  [fs]
    type = PorousFlowWaterNCG
    water_fp = water
    gas_fp = methane
    capillary_pressure = pc
  []
[]

[FluidProperties]
  [methane]
    type = MethaneFluidProperties
  []
  [water]
    type = Water97FluidProperties
  []
[]

[Materials]
  [temperature]
    type = PorousFlowTemperature
    temperature = temperature
  []
  [waterncg]
    type = PorousFlowFluidState
    gas_porepressure = pgas
    z = zi
    temperature = temperature
    temperature_unit = Celsius
    capillary_pressure = pc
    fluid_state = fs
  []
  [porosity]
    type = PorousFlowPorosityConst
    porosity = 0.2
  []
  [permeability]
    type = PorousFlowPermeabilityConst
    permeability = '1e-12 0 0 0 1e-12 0 0 0 1e-12'
  []
  [relperm_water]
    type = PorousFlowRelativePermeabilityCorey
    n = 2
    phase = 0
    s_res = 0.1
    sum_s_res = 0.1
  []
  [relperm_gas]
    type = PorousFlowRelativePermeabilityCorey
    n = 2
    phase = 1
  []
  [rockheat]
    type = PorousFlowMatrixInternalEnergy
    specific_heat_capacity = 1000
    density = 2500
  []
  [rock_thermal_conductivity]
    type = PorousFlowThermalConductivityIdeal
    dry_thermal_conductivity = '50 0 0  0 50 0  0 0 50'
  []
[]

[BCs]
  [cold_gas]
    type = DirichletBC
    boundary = left
    variable = temperature
    value = 20
  []
  [gas_injecton]
    type = PorousFlowSink
    boundary = left
    variable = zi
    flux_function = -0.159155
  []
  [rightwater]
    type = DirichletBC
    boundary = right
    value = 20e6
    variable = pgas
  []
  [righttemp]
    type = DirichletBC
    boundary = right
    value = 70
    variable = temperature
  []
[]

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

[Executioner]
  type = Transient
  solve_type = NEWTON
  end_time = 1e4
  nl_abs_tol = 1e-7
  nl_rel_tol = 1e-5
  # Avoids failing first time step in parallel
  line_search = 'none'
  [TimeStepper]
    type = IterationAdaptiveDT
    dt = 1
    growth_factor = 1.5
  []
[]

[Postprocessors]
  [pgas]
    type = PointValue
    point = '2 0 0'
    variable = pgas
  []
  [sgas]
    type = PointValue
    point = '2 0 0'
    variable = saturation_gas
  []
  [zi]
    type = PointValue
    point = '2 0 0'
    variable = zi
  []
  [temperature]
    type = PointValue
    point = '2 0 0'
    variable = temperature
  []
  [massgas]
    type = PorousFlowFluidMass
    fluid_component = 1
  []
  [x1]
    type = PointValue
    point = '2 0 0'
    variable = x1
  []
  [y0]
    type = PointValue
    point = '2 0 0'
    variable = y0
  []
[]

[Outputs]
  print_linear_residuals = false
  perf_graph = true
  csv = true
[]

(modules/fluid_properties/test/tests/methane/methane.i)

# Test MethaneFluidProperties
# Reference data from Irvine Jr, T. F. and Liley, P. E. (1984) Steam and
# Gas Tables with Computer Equations
#
# For temperature = 350K, the fluid properties should be:
# density = 55.13 kg/m^3
# viscosity = 0.01276 mPa.s
# cp = 2.375 kJ/kg/K
# h = 708.5 kJ/kg
# s = 11.30 kJ/kg/K
# c = 481.7 m/s
# k = 0.04113 W/m/K

[Mesh]
  type = GeneratedMesh
  dim = 2
  nx = 1
  ny = 1
[]

[Variables]
  [./dummy]
  [../]
[]

[AuxVariables]
  [./pressure]
    family = MONOMIAL
    order = CONSTANT
    initial_condition = 10.0e6
  [../]
  [./temperature]
    family = MONOMIAL
    order = CONSTANT
    initial_condition = 350
  [../]
  [./density]
    family = MONOMIAL
    order = CONSTANT
  [../]
  [./viscosity]
    family = MONOMIAL
    order = CONSTANT
  [../]
  [./cp]
    family = MONOMIAL
    order = CONSTANT
  [../]
  [./cv]
    family = MONOMIAL
    order = CONSTANT
  [../]
  [./internal_energy]
    family = MONOMIAL
    order = CONSTANT
  [../]
  [./enthalpy]
    family = MONOMIAL
    order = CONSTANT
  [../]
  [./entropy]
    family = MONOMIAL
    order = CONSTANT
  [../]
  [./thermal_cond]
    family = MONOMIAL
    order = CONSTANT
  [../]
  [./c]
    family = MONOMIAL
    order = CONSTANT
  [../]
[]

[AuxKernels]
  [./density]
    type = MaterialRealAux
     variable = density
     property = density
  [../]
  [./viscosity]
    type = MaterialRealAux
     variable = viscosity
     property = viscosity
  [../]
  [./cp]
    type = MaterialRealAux
     variable = cp
     property = cp
  [../]
  [./cv]
    type = MaterialRealAux
     variable = cv
     property = cv
  [../]
  [./e]
    type = MaterialRealAux
     variable = internal_energy
     property = e
  [../]
  [./enthalpy]
    type = MaterialRealAux
     variable = enthalpy
     property = h
  [../]
  [./entropy]
    type = MaterialRealAux
     variable = entropy
     property = s
  [../]
  [./thermal_cond]
    type = MaterialRealAux
     variable = thermal_cond
     property = k
  [../]
  [./c]
    type = MaterialRealAux
     variable = c
     property = c
  [../]
[]

[FluidProperties]
  [./methane]
    type = MethaneFluidProperties
  [../]
[]

[Materials]
  [./fp_mat]
    type = FluidPropertiesMaterialPT
    pressure = pressure
    temperature = temperature
    fp = methane
  [../]
[]

[Kernels]
  [./diff]
    type = Diffusion
    variable = dummy
  [../]
[]

[Executioner]
  type = Steady
  solve_type = NEWTON
[]

[Outputs]
  exodus = true
[]

Properties of methane
Range of validity
Input Parameters
Input Files
References