Start | Previous | Next

Porous Flow Tutorial Page 05. Using a realistic equation of state for the fluid

It is trivial to add a realistic equation of state to any PorousFlow simulation. For instance, the high-precision Water97 equation of state may be used:

# Darcy flow with heat advection and conduction, using Water97 properties
[Mesh<<<{"href": "../../syntax/Mesh/index.html"}>>>]
  [annular]
    type = AnnularMeshGenerator<<<{"description": "For rmin>0: creates an annular mesh of QUAD4 elements. For rmin=0: creates a disc mesh of QUAD4 and TRI3 elements. Boundary sidesets are created at rmax and rmin, and given these names. If dmin!=0 and dmax!=360, a sector of an annulus or disc is created. In this case boundary sidesets are also created at dmin and dmax, and given these names", "href": "../../source/meshgenerators/AnnularMeshGenerator.html"}>>>
    nr<<<{"description": "Number of elements in the radial direction"}>>> = 10
    rmin<<<{"description": "Inner radius.  If rmin=0 then a disc mesh (with no central hole) will be created."}>>> = 1.0
    rmax<<<{"description": "Outer radius"}>>> = 10
    growth_r<<<{"description": "The ratio of radial sizes of successive rings of elements"}>>> = 1.4
    nt<<<{"description": "Number of elements in the angular direction"}>>> = 4
    dmin<<<{"description": "Minimum degree, measured in degrees anticlockwise from x axis"}>>> = 0
    dmax<<<{"description": "Maximum angle, measured in degrees anticlockwise from x axis. If dmin=0 and dmax=360 an annular mesh is created. Otherwise, only a sector of an annulus is created"}>>> = 90
  []
  [make3D]
    type = MeshExtruderGenerator<<<{"description": "Takes a 1D or 2D mesh and extrudes the entire structure along the specified axis increasing the dimensionality of the mesh.", "href": "../../source/meshgenerators/MeshExtruderGenerator.html"}>>>
    extrusion_vector<<<{"description": "The direction and length of the extrusion"}>>> = '0 0 12'
    num_layers<<<{"description": "The number of layers in the extruded mesh"}>>> = 3
    bottom_sideset<<<{"description": "The boundary that will be applied to the bottom of the extruded mesh"}>>> = 'bottom'
    top_sideset<<<{"description": "The boundary that will be to the top of the extruded mesh"}>>> = 'top'
    input<<<{"description": "the mesh we want to extrude"}>>> = annular
  []
  [shift_down]
    type = TransformGenerator<<<{"description": "Applies a linear transform to the entire mesh.", "href": "../../source/meshgenerators/TransformGenerator.html"}>>>
    transform<<<{"description": "The type of transformation to perform (TRANSLATE, TRANSLATE_CENTER_ORIGIN, TRANSLATE_MIN_ORIGIN, ROTATE, SCALE)"}>>> = TRANSLATE
    vector_value<<<{"description": "The value to use for the transformation. When using TRANSLATE or SCALE, the xyz coordinates are applied in each direction respectively. When using ROTATE, the values are interpreted as the Euler angles phi, theta and psi given in degrees."}>>> = '0 0 -6'
    input<<<{"description": "The mesh we want to modify"}>>> = make3D
  []
  [aquifer]
    type = SubdomainBoundingBoxGenerator<<<{"description": "Changes the subdomain ID of elements either (XOR) inside or outside the specified box to the specified ID.", "href": "../../source/meshgenerators/SubdomainBoundingBoxGenerator.html"}>>>
    block_id<<<{"description": "Subdomain id to set for inside/outside the bounding box"}>>> = 1
    bottom_left<<<{"description": "The bottom left point (in x,y,z with spaces in-between)."}>>> = '0 0 -2'
    top_right<<<{"description": "The bottom left point (in x,y,z with spaces in-between)."}>>> = '10 10 2'
    input<<<{"description": "The mesh we want to modify"}>>> = shift_down
  []
  [injection_area]
    type = ParsedGenerateSideset<<<{"description": "A MeshGenerator that adds element sides to a sideset if the centroid of the side satisfies the `combinatorial_geometry` expression.", "href": "../../source/meshgenerators/ParsedGenerateSideset.html"}>>>
    combinatorial_geometry<<<{"description": "Function expression encoding a combinatorial geometry"}>>> = 'x*x+y*y<1.01'
    included_subdomains<<<{"description": "A set of subdomain names or ids whose sides will be included in the new sidesets. A side is only added if the subdomain id of the corresponding element is in this set."}>>> = 1
    new_sideset_name<<<{"description": "The name of the new sideset"}>>> = 'injection_area'
    input<<<{"description": "The mesh we want to modify"}>>> = 'aquifer'
  []
  [rename]
    type = RenameBlockGenerator<<<{"description": "Changes the block IDs and/or block names for a given set of blocks defined by either block ID or block name. The changes are independent of ordering. The merging of blocks is supported.", "href": "../../source/meshgenerators/RenameBlockGenerator.html"}>>>
    old_block<<<{"description": "Elements with these block ID(s)/name(s) will be given the new block information specified in 'new_block'"}>>> = '0 1'
    new_block<<<{"description": "The new block ID(s)/name(s) to be given by the elements defined in 'old_block'."}>>> = 'caps aquifer'
    input<<<{"description": "The mesh we want to modify"}>>> = 'injection_area'
  []
[]

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

[Variables<<<{"href": "../../syntax/Variables/index.html"}>>>]
  [porepressure]
    initial_condition<<<{"description": "Specifies a constant initial condition for this variable"}>>> = 1E6
  []
  [temperature]
    initial_condition<<<{"description": "Specifies a constant initial condition for this variable"}>>> = 313
    scaling<<<{"description": "Specifies a scaling factor to apply to this variable"}>>> = 1E-8
  []
[]

[PorousFlowBasicTHM<<<{"href": "../../syntax/PorousFlowBasicTHM/index.html"}>>>]
  porepressure<<<{"description": "The name of the porepressure variable"}>>> = porepressure
  temperature<<<{"description": "For isothermal simulations, this is the temperature at which fluid properties (and stress-free strains) are evaluated at.  Otherwise, this is the name of the temperature variable.  Units = Kelvin"}>>> = temperature
  coupling_type<<<{"description": "The type of simulation.  For simulations involving Mechanical deformations, you will need to supply the correct Biot coefficient.  For simulations involving Thermal flows, you will need an associated ConstantThermalExpansionCoefficient Material"}>>> = ThermoHydro
  gravity<<<{"description": "Gravitational acceleration vector downwards (m/s^2)"}>>> = '0 0 0'
  fp<<<{"description": "The name of the user object for fluid properties. Only needed if fluid_properties_type = PorousFlowSingleComponentFluid"}>>> = the_simple_fluid
[]

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

[FluidProperties<<<{"href": "../../syntax/FluidProperties/index.html"}>>>]
  [the_simple_fluid]
    type = Water97FluidProperties<<<{"description": "Fluid properties for water and steam (H2O) using IAPWS-IF97", "href": "../../source/fluidproperties/Water97FluidProperties.html"}>>>
  []
[]
(modules/porous_flow/examples/tutorial/05.i)

(The name "the_simple_fluid" could also be changed to something more appropriate.)

In this case, the initial and boundary conditions must also be changed because the Water97 equation of state is not defined for (an absolute vacuum). For example:

[Variables<<<{"href": "../../syntax/Variables/index.html"}>>>]
  [porepressure]
    initial_condition<<<{"description": "Specifies a constant initial condition for this variable"}>>> = 1E6
  []
  [temperature]
    initial_condition<<<{"description": "Specifies a constant initial condition for this variable"}>>> = 313
    scaling<<<{"description": "Specifies a scaling factor to apply to this variable"}>>> = 1E-8
  []
[]
(modules/porous_flow/examples/tutorial/05.i)

and

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

[FluidProperties<<<{"href": "../../syntax/FluidProperties/index.html"}>>>]
  [the_simple_fluid]
    type = Water97FluidProperties<<<{"description": "Fluid properties for water and steam (H2O) using IAPWS-IF97", "href": "../../source/fluidproperties/Water97FluidProperties.html"}>>>
  []
[]

[Materials<<<{"href": "../../syntax/Materials/index.html"}>>>]
  [porosity]
    type = PorousFlowPorosity<<<{"description": "This Material calculates the porosity PorousFlow simulations", "href": "../../source/materials/PorousFlowPorosity.html"}>>>
    porosity_zero<<<{"description": "The porosity at zero volumetric strain and reference temperature and reference effective porepressure and reference chemistry.  This must be a real number or a constant monomial variable (not a linear lagrange or other type of variable)"}>>> = 0.1
  []
  [biot_modulus]
    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"}>>>
    biot_coefficient<<<{"description": "Biot coefficient"}>>> = 0.8
    solid_bulk_compliance<<<{"description": "Reciprocal of the drained bulk modulus of the porous skeleton.  If strain = C * stress, then solid_bulk_compliance = de_ij de_kl C_ijkl.  If the grain bulk modulus is Kg then 1/Kg = (1 - biot_coefficient) * solid_bulk_compliance."}>>> = 2E-7
    fluid_bulk_modulus<<<{"description": "Fluid bulk modulus"}>>> = 1E7
  []
  [permeability_aquifer]
    type = PorousFlowPermeabilityConst<<<{"description": "This Material calculates the permeability tensor assuming it is constant", "href": "../../source/materials/PorousFlowPermeabilityConst.html"}>>>
    block<<<{"description": "The list of blocks (ids or names) that this object will be applied"}>>> = aquifer
    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'
  []
  [permeability_caps]
    type = PorousFlowPermeabilityConst<<<{"description": "This Material calculates the permeability tensor assuming it is constant", "href": "../../source/materials/PorousFlowPermeabilityConst.html"}>>>
    block<<<{"description": "The list of blocks (ids or names) that this object will be applied"}>>> = caps
    permeability<<<{"description": "The permeability tensor (usually in m^2), which is assumed constant for this material"}>>> = '1E-15 0 0   0 1E-15 0   0 0 1E-16'
  []

  [thermal_expansion]
    type = PorousFlowConstantThermalExpansionCoefficient<<<{"description": "Computes the effective thermal expansion coefficient, (biot_coeff - porosity) * drained_coefficient + porosity * fluid_coefficient.", "href": "../../source/materials/PorousFlowConstantThermalExpansionCoefficient.html"}>>>
    biot_coefficient<<<{"description": "Biot coefficient"}>>> = 0.8
    drained_coefficient<<<{"description": "Volumetric coefficient of thermal expansion of the drained porous skeleton (ie the porous rock without fluid, or with a fluid that is free to move in and out of the rock)"}>>> = 0.003
    fluid_coefficient<<<{"description": "Volumetric coefficient of thermal expansion for the fluid"}>>> = 0.0002
  []
  [rock_internal_energy]
    type = PorousFlowMatrixInternalEnergy<<<{"description": "This Material calculates the internal energy of solid rock grains, which is specific_heat_capacity * density * temperature.  Kernels multiply this by (1 - porosity) to find the energy density of the porous rock in a rock-fluid system", "href": "../../source/materials/PorousFlowMatrixInternalEnergy.html"}>>>
    density<<<{"description": "Density of the rock grains"}>>> = 2500.0
    specific_heat_capacity<<<{"description": "Specific heat capacity of the rock grains (J/kg/K)."}>>> = 1200.0
  []
  [thermal_conductivity]
    type = PorousFlowThermalConductivityIdeal<<<{"description": "This Material calculates rock-fluid combined thermal conductivity by using a weighted sum.  Thermal conductivity = dry_thermal_conductivity + S^exponent * (wet_thermal_conductivity - dry_thermal_conductivity), where S is the aqueous saturation", "href": "../../source/materials/PorousFlowThermalConductivityIdeal.html"}>>>
    dry_thermal_conductivity<<<{"description": "The thermal conductivity of the rock matrix when the aqueous saturation is zero"}>>> = '10 0 0  0 10 0  0 0 10'
    block<<<{"description": "The list of blocks (ids or names) that this object will be applied"}>>> = 'caps aquifer'
  []
[]

[Preconditioning<<<{"href": "../../syntax/Preconditioning/index.html"}>>>]
  active<<<{"description": "If specified only the blocks named will be visited and made active"}>>> = basic
  [basic]
    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"}>>> = '-ksp_diagonal_scale -ksp_diagonal_scale_fix'
    petsc_options_iname<<<{"description": "Names of PETSc name/value pairs"}>>> = '-pc_type -sub_pc_type -sub_pc_factor_shift_type -pc_asm_overlap'
    petsc_options_value<<<{"description": "Values of PETSc name/value pairs (must correspond with \"petsc_options_iname\""}>>> = ' asm      lu           NONZERO                   2'
  []
  [preferred_but_might_not_be_installed]
    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_iname<<<{"description": "Names of PETSc name/value pairs"}>>> = '-pc_type -pc_factor_mat_solver_package'
    petsc_options_value<<<{"description": "Values of PETSc name/value pairs (must correspond with \"petsc_options_iname\""}>>> = ' lu       mumps'
  []
[]

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

[Outputs<<<{"href": "../../syntax/Outputs/index.html"}>>>]
  exodus<<<{"description": "Output the results using the default settings for Exodus output."}>>> = true
[]
(modules/porous_flow/examples/tutorial/05.i)

Using realistic high-precision equations of state can cause PorousFlow to run quite slowly, because the equations of state are so complicated to evaluate. It is always recommended to use TabulatedFluidProperties in the following way:

# Darcy flow with heat advection and conduction, using Water97 properties
[Mesh<<<{"href": "../../syntax/Mesh/index.html"}>>>]
  [annular]
    type = AnnularMeshGenerator<<<{"description": "For rmin>0: creates an annular mesh of QUAD4 elements. For rmin=0: creates a disc mesh of QUAD4 and TRI3 elements. Boundary sidesets are created at rmax and rmin, and given these names. If dmin!=0 and dmax!=360, a sector of an annulus or disc is created. In this case boundary sidesets are also created at dmin and dmax, and given these names", "href": "../../source/meshgenerators/AnnularMeshGenerator.html"}>>>
    nr<<<{"description": "Number of elements in the radial direction"}>>> = 10
    rmin<<<{"description": "Inner radius.  If rmin=0 then a disc mesh (with no central hole) will be created."}>>> = 1.0
    rmax<<<{"description": "Outer radius"}>>> = 10
    growth_r<<<{"description": "The ratio of radial sizes of successive rings of elements"}>>> = 1.4
    nt<<<{"description": "Number of elements in the angular direction"}>>> = 4
    dmin<<<{"description": "Minimum degree, measured in degrees anticlockwise from x axis"}>>> = 0
    dmax<<<{"description": "Maximum angle, measured in degrees anticlockwise from x axis. If dmin=0 and dmax=360 an annular mesh is created. Otherwise, only a sector of an annulus is created"}>>> = 90
  []
  [make3D]
    type = MeshExtruderGenerator<<<{"description": "Takes a 1D or 2D mesh and extrudes the entire structure along the specified axis increasing the dimensionality of the mesh.", "href": "../../source/meshgenerators/MeshExtruderGenerator.html"}>>>
    extrusion_vector<<<{"description": "The direction and length of the extrusion"}>>> = '0 0 12'
    num_layers<<<{"description": "The number of layers in the extruded mesh"}>>> = 3
    bottom_sideset<<<{"description": "The boundary that will be applied to the bottom of the extruded mesh"}>>> = 'bottom'
    top_sideset<<<{"description": "The boundary that will be to the top of the extruded mesh"}>>> = 'top'
    input<<<{"description": "the mesh we want to extrude"}>>> = annular
  []
  [shift_down]
    type = TransformGenerator<<<{"description": "Applies a linear transform to the entire mesh.", "href": "../../source/meshgenerators/TransformGenerator.html"}>>>
    transform<<<{"description": "The type of transformation to perform (TRANSLATE, TRANSLATE_CENTER_ORIGIN, TRANSLATE_MIN_ORIGIN, ROTATE, SCALE)"}>>> = TRANSLATE
    vector_value<<<{"description": "The value to use for the transformation. When using TRANSLATE or SCALE, the xyz coordinates are applied in each direction respectively. When using ROTATE, the values are interpreted as the Euler angles phi, theta and psi given in degrees."}>>> = '0 0 -6'
    input<<<{"description": "The mesh we want to modify"}>>> = make3D
  []
  [aquifer]
    type = SubdomainBoundingBoxGenerator<<<{"description": "Changes the subdomain ID of elements either (XOR) inside or outside the specified box to the specified ID.", "href": "../../source/meshgenerators/SubdomainBoundingBoxGenerator.html"}>>>
    block_id<<<{"description": "Subdomain id to set for inside/outside the bounding box"}>>> = 1
    bottom_left<<<{"description": "The bottom left point (in x,y,z with spaces in-between)."}>>> = '0 0 -2'
    top_right<<<{"description": "The bottom left point (in x,y,z with spaces in-between)."}>>> = '10 10 2'
    input<<<{"description": "The mesh we want to modify"}>>> = shift_down
  []
  [injection_area]
    type = ParsedGenerateSideset<<<{"description": "A MeshGenerator that adds element sides to a sideset if the centroid of the side satisfies the `combinatorial_geometry` expression.", "href": "../../source/meshgenerators/ParsedGenerateSideset.html"}>>>
    combinatorial_geometry<<<{"description": "Function expression encoding a combinatorial geometry"}>>> = 'x*x+y*y<1.01'
    included_subdomains<<<{"description": "A set of subdomain names or ids whose sides will be included in the new sidesets. A side is only added if the subdomain id of the corresponding element is in this set."}>>> = 1
    new_sideset_name<<<{"description": "The name of the new sideset"}>>> = 'injection_area'
    input<<<{"description": "The mesh we want to modify"}>>> = 'aquifer'
  []
  [rename]
    type = RenameBlockGenerator<<<{"description": "Changes the block IDs and/or block names for a given set of blocks defined by either block ID or block name. The changes are independent of ordering. The merging of blocks is supported.", "href": "../../source/meshgenerators/RenameBlockGenerator.html"}>>>
    old_block<<<{"description": "Elements with these block ID(s)/name(s) will be given the new block information specified in 'new_block'"}>>> = '0 1'
    new_block<<<{"description": "The new block ID(s)/name(s) to be given by the elements defined in 'old_block'."}>>> = 'caps aquifer'
    input<<<{"description": "The mesh we want to modify"}>>> = 'injection_area'
  []
[]

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

[Variables<<<{"href": "../../syntax/Variables/index.html"}>>>]
  [porepressure]
    initial_condition<<<{"description": "Specifies a constant initial condition for this variable"}>>> = 1E6
  []
  [temperature]
    initial_condition<<<{"description": "Specifies a constant initial condition for this variable"}>>> = 313
    scaling<<<{"description": "Specifies a scaling factor to apply to this variable"}>>> = 1E-8
  []
[]

[PorousFlowBasicTHM<<<{"href": "../../syntax/PorousFlowBasicTHM/index.html"}>>>]
  porepressure<<<{"description": "The name of the porepressure variable"}>>> = porepressure
  temperature<<<{"description": "For isothermal simulations, this is the temperature at which fluid properties (and stress-free strains) are evaluated at.  Otherwise, this is the name of the temperature variable.  Units = Kelvin"}>>> = temperature
  coupling_type<<<{"description": "The type of simulation.  For simulations involving Mechanical deformations, you will need to supply the correct Biot coefficient.  For simulations involving Thermal flows, you will need an associated ConstantThermalExpansionCoefficient Material"}>>> = ThermoHydro
  gravity<<<{"description": "Gravitational acceleration vector downwards (m/s^2)"}>>> = '0 0 0'
  fp<<<{"description": "The name of the user object for fluid properties. Only needed if fluid_properties_type = PorousFlowSingleComponentFluid"}>>> = tabulated_water
[]

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

[FluidProperties<<<{"href": "../../syntax/FluidProperties/index.html"}>>>]
  [true_water]
    type = Water97FluidProperties<<<{"description": "Fluid properties for water and steam (H2O) using IAPWS-IF97", "href": "../../source/fluidproperties/Water97FluidProperties.html"}>>>
  []
  [tabulated_water]
    type = TabulatedFluidProperties<<<{"description": "Fluid properties using bicubic interpolation on tabulated values provided", "href": "../../source/fluidproperties/TabulatedBicubicFluidProperties.html"}>>>
    fp<<<{"description": "The name of the FluidProperties UserObject"}>>> = true_water
    temperature_min<<<{"description": "Minimum temperature for tabulated data."}>>> = 275
    interpolated_properties<<<{"description": "Properties to interpolate if no data file is provided"}>>> = 'density viscosity enthalpy internal_energy'
    fluid_property_output_file<<<{"description": "Name of the CSV file which can be output with the tabulation. This file can then be read as a 'fluid_property_file'"}>>> = water97_tabulated.csv
    # Comment out the fp parameter and uncomment below to use the newly generated tabulation
    # fluid_property_file = water97_tabulated.csv
  []
[]
(modules/porous_flow/examples/tutorial/05_tabulated.i)

and

[PorousFlowBasicTHM<<<{"href": "../../syntax/PorousFlowBasicTHM/index.html"}>>>]
  porepressure<<<{"description": "The name of the porepressure variable"}>>> = porepressure
  temperature<<<{"description": "For isothermal simulations, this is the temperature at which fluid properties (and stress-free strains) are evaluated at.  Otherwise, this is the name of the temperature variable.  Units = Kelvin"}>>> = temperature
  coupling_type<<<{"description": "The type of simulation.  For simulations involving Mechanical deformations, you will need to supply the correct Biot coefficient.  For simulations involving Thermal flows, you will need an associated ConstantThermalExpansionCoefficient Material"}>>> = ThermoHydro
  gravity<<<{"description": "Gravitational acceleration vector downwards (m/s^2)"}>>> = '0 0 0'
  fp<<<{"description": "The name of the user object for fluid properties. Only needed if fluid_properties_type = PorousFlowSingleComponentFluid"}>>> = tabulated_water
[]
(modules/porous_flow/examples/tutorial/05_tabulated.i)

Another advantage of using TabulatedfluidProperties is that the bounds on pressure and temperature imposed by the original "true" equation of state can be extrapolated past, which can help to remove problems of MOOSE trying to sample outside the original's region of validity (and thus causing the timestep to be cut).

Finally, in the case at hand, it is worth reminding the reader that PorousFlowBasicTHM is based on the assumptions of a constant, large fluid bulk modulus, and a constant fluid thermal expansion coefficient, which are incorrect for some fluids. User beware!

Start | Previous | Next