THM-Rehbinder test description
Introduction
Rehbinder (Rehbinder, 1995) derives analytical solutions of cylindrically-symmetric and spherically-symmetric THM problems in certain limits. The solution is a steady-state THM solution of a pressurised and heated cavity of radius in a finite domain of radius . In the cylindrically-symmetric case the cavity is a cylinder of radius and the domain is a cylinder of radius . In this document only the cylindrically-symmetric situation is explored. The cylindrical coordinates are denoted by .
Initial and boundary conditions
The initial conditions are zero porepressure, ; zero temperature, ; and zero displacement . The boundary conditions are at the cavity wall are
The latter condition implies the total radial stress is at the cavity wall, corresponding to the fluid in the cavity pushing on the cavity wall. The boundary conditions at the outer radius are
and either or . In this document the former is called the free outer'' boundary condition, while the latter is called the fixed outer'' boundary condition
Simplifying assumptions
In order to derive the solution, Rehbinder makes various simplifications. Translated to the language of the other PorousFlow documents these are as follows.
There is no gravity.
All quantities depend only on the radial coordinate, , and not the angular coordinate or the axial coordinate .
Plane-strain is assumed, so .
There is only one fully-saturated fluid phase that contains one component.
The fluid relative permeability is unity.
The fluid density is constant. In the numerical simulation below, the density is assumed to obey , with , and .
The fluid dynamic viscosity is constant. In the numerical simulation below, the dynamic viscosity is .
The Biot coefficient is unity .
The fluid internal energy is assumed to be linear in the temperature. In the numerical simulation below, where .
The fluid enthalpy is assumed to be equal to the fluid internal energy.
The velocity of the solid skeleton is zero: . To model this using PorousFlow, the
{\tt VolumetricExpansion}'' Kernels are not included.
The rock-matrix density is constant. In the numerical simulation below, the value is chosen.
The rock-matrix specific heat capacity is constant. In the numerical simulation below, the value is chosen.
The rock deforms elastically (so there is no plastic heating). In the numerical simulation below, the Young's modulus is and the Poisson's ratio is .
The porosity is constant. In the simulation below is chosen.
The permeability is constant. In the simulation below is chosen.
The thermal conductivity of the rock-fluid system is constant. In the simulation below is chosen.
Steady-state heat flow, fluid-flow and mechanical deformation has been reached (Rehbinder equation (41)). Using the aove assumptions, this is true if (the left-hand-side is measured in seconds, the right-hand-side in metres). In the simulation below, steady-state is achieved by not including any time-derivative Kernels.
The liquid flow is quasi-stationary, in the sense that the diffusion of heat is much slower than the diffusion of pore pressure (Rehbinder equation (32)). Using the above assumptions this may be written as where is the thermal conductivity of the rock-fluid system, and . Using the above values, this condition reads which is clearly satisfied.
The liquid flow is quasi-stationary, in the sense that a pressure change in the cavity is transmitted instantaneously through the pores to the outer boundary (Rehbinder equation (33)). Using the above assumptions this may be written as . In the simulation below, and are chosen.
The heat is mainly conducted through the matrix, and a negligible part is convected with the flux. This is true if the Peclet number is very small (Rehbinder equation (35))
(Rehbinder writes this in terms of a reference temperature, , which is chosen here to be .) Using the above values this is . In the simulation below, and are chosen, yielding .
Rehbinder's derives the solution of this THM problem as an expansion in the Peclet number. I shall not write the solution here as it is fairly lengthy. The paper contains a few minor typos and they are corrected in the accompanying Python file:
#!/usr/bin/env python3
import os
import sys
import numpy as np
import matplotlib.pyplot as plt
def rehbinder(r):
# Results from Rehbinder with parameters used in the MOOSE simulation.
# Rehbinder's manuscript contains a few typos - I've corrected them here.
# G Rehbinder "Analytic solutions of stationary coupled thermo-hydro-mechanical problems" Int J Rock Mech Min Sci & Geomech Abstr 32 (1995) 453-463
poisson = 0.2
thermal_expansion = 1E-6
young = 1E10
fluid_density = 1000
fluid_specific_heat = 1000
permeability = 1E-12
fluid_viscosity = 1E-3
thermal_conductivity = 1E6
P0 = 1E6
T0 = 1E3
Tref = T0
r0 = 0.1
r1 = 1.0
xi = r / r0
xi1 = r1 / r0
Peclet = fluid_density * fluid_specific_heat * thermal_expansion * Tref * young * permeability / fluid_viscosity / thermal_conductivity / (1 - poisson)
That0 = T0 / Tref
sigmahat0 = -P0 * (1 - poisson) / thermal_expansion / Tref / young
Tzeroth = That0 * (1 - np.log(xi) / np.log(xi1))
Tfirst_pr = 2 * sigmahat0 * That0 * xi * (np.log(xi) - np.log(xi1)) / np.log(xi1)**2
Cone = 2 * That0 * sigmahat0 * (2 + np.log(xi1)) / np.log(xi1)**2
Cone = 2 * That0 * sigmahat0 / np.log(xi1) # Corrected Eqn(87)
Done = 2 * That0 * sigmahat0 * (2 * (xi1 - 1) / np.log(xi1) - 1) / np.log(xi1)**2
Done = 2 * That0 * sigmahat0 * (- 1) / np.log(xi1)**2 # Corrected Eqn(87)
Tfirst_hm = Cone + Done * np.log(xi)
Tfirst = Tfirst_pr + Tfirst_hm
That = Tzeroth + Peclet * Tfirst
T = Tref * That
Pzeroth = -sigmahat0 * (1 - np.log(xi) / np.log(xi1))
Pfirst = 0
Phat = Pzeroth + Peclet * Pfirst
P = thermal_expansion * Tref * young * Phat / (1 - poisson)
g0 = Tzeroth + (1 - 2 * poisson) * Pzeroth / (1 - poisson)
uzeroth_pr = (That0 - sigmahat0 * (1 - 2 * poisson) / (1 - poisson)) * (0.5 * (xi**2 - 1) - 0.25 * (1 - xi**2 + 2 * xi**2 * np.log(xi)) / np.log(xi1)) / xi
uzeroth_pr_xi1 = (That0 - sigmahat0 * (1 - 2 * poisson) / (1 - poisson)) * (0.5 * (xi1**2 - 1) - 0.25 * (1 - xi1**2 + 2 * xi1**2 * np.log(xi1)) / np.log(xi1)) / xi1
# fixed outer boundary
Bzeroth = - ((1 - 2 * poisson) * sigmahat0 + uzeroth_pr_xi1 / xi1) / (1 - 2 * poisson + 1.0 / xi1)
Azeroth = - Bzeroth / xi1**2 - uzeroth_pr_xi1 / xi1
fixed_uzeroth_hm = Azeroth * xi + Bzeroth / xi
fixed_uzeroth = uzeroth_pr + fixed_uzeroth_hm
# free outer boundary
Bzeroth = (xi1**2 * sigmahat0 - xi1 * uzeroth_pr_xi1) / (1 - xi1**2)
Azeroth = (1 - 2 * poisson) * (Bzeroth + sigmahat0)
free_uzeroth_hm = Azeroth * xi + Bzeroth / xi
free_uzeroth = uzeroth_pr + free_uzeroth_hm
ufirst_pr = (1.0 / xi) * (0.5 * (xi**2 - 1) * (2 * Cone - Done) + 0.5 * Done * xi**2 * np.log(xi) + 2 * sigmahat0 * That0 / np.log(xi1)**2 * (xi**3 * np.log(xi) / 3 + (1 - xi**3) / 9 + 0.5 * np.log(xi1) * (1 - xi**2)))
ufirst_pr_xi1 = (1.0 / xi1) * (0.5 * (xi1**2 - 1) * (2 * Cone - Done) + 0.5 * Done * xi1**2 * np.log(xi1) + 2 * sigmahat0 * That0 / np.log(xi1)**2 * (xi1**3 * np.log(xi1) / 3 + (1 - xi1**3) / 9 + 0.5 * np.log(xi1) * (1 - xi1**2)))
# fixed outer boundary
Bfirst = - ufirst_pr_xi1 / xi1 / (1 - 2 * poisson + 1.0 / xi1**2)
Afirst = - Bfirst / xi1**2 - ufirst_pr_xi1 / xi1
fixed_ufirst_hm = Afirst * xi + Bfirst / xi
fixed_ufirst = ufirst_pr + fixed_ufirst_hm
# free outer boundary
Bfirst = xi1 * ufirst_pr_xi1 / (1 - xi1**2)
Afirst = (1 - 2 * poisson) * Bfirst
free_ufirst_hm = Afirst * xi + Bfirst / xi
free_ufirst = ufirst_pr + free_ufirst_hm
fixed_uhat = fixed_uzeroth + Peclet * fixed_ufirst
fixed_u = thermal_expansion * Tref * r0 * fixed_uhat * (1 + poisson) / (1 - poisson) # Corrected Eqn(16)
free_uhat = free_uzeroth + Peclet * free_ufirst
free_u = thermal_expansion * Tref * r0 * free_uhat * (1 + poisson) / (1 - poisson) # Corrected Eqn(16)
return (T, P, fixed_u, free_u)
def moose(fn, rcol, datacol):
try:
f = open(fn)
data = f.readlines()[1:-1]
data = [list(map(float, d.strip().split(","))) for d in data]
data = ([d[rcol] for d in data], [d[datacol] for d in data])
f.close()
except:
sys.stderr.write("Cannot read " + fn + ", or it contains erroneous data\n")
sys.exit(1)
return data
mooser = [0.1 * i for i in range(1, 11)]
fixedT = moose("gold/fixed_outer_T_0001.csv", 0, 4)
fixedP = moose("gold/fixed_outer_P_0001.csv", 0, 4)
fixedu = moose("gold/fixed_outer_U_0001.csv", 0, 4)
freeu = moose("gold/free_outer_U_0001.csv", 0, 4)
rpoints = np.arange(0.1, 1.0, 0.01)
expected = list(zip(*[rehbinder(r) for r in rpoints]))
plt.figure()
plt.plot(rpoints, expected[0], 'k-', linewidth = 3.0, label = 'expected')
plt.plot(fixedT[0], fixedT[1], 'rs', markersize = 10.0, label = 'MOOSE')
plt.legend(loc = 'upper right')
plt.xlabel("r (m)")
plt.ylabel("Temperature (K)")
plt.title("Temperature around cavity")
plt.savefig("temperature_fig.png")
plt.figure()
plt.plot(rpoints, [1E-6 * p for p in expected[1]], 'k-', linewidth = 3.0, label = 'expected')
plt.plot(fixedP[0], [1E-6 * p for p in fixedP[1]], 'rs', markersize = 10.0, label = 'MOOSE')
plt.legend(loc = 'upper right')
plt.xlabel("r (m)")
plt.ylabel("Porepressure (MPa)")
plt.title("Porepressure around cavity")
plt.savefig("porepressure_fig.png")
plt.figure()
plt.plot(rpoints, [1000 * u for u in expected[2]], 'k-', linewidth = 3.0, label = 'expected (fixed)')
plt.plot(fixedu[0], [1000 * u for u in fixedu[1]], 'rs', markersize = 10.0, label = 'MOOSE (fixed)')
plt.plot(rpoints, [1000 * u for u in expected[3]], 'b-', linewidth = 2.0, label = 'expected (free)')
plt.plot(freeu[0], [1000 * u for u in freeu[1]], 'g*', markersize = 13.0, label = 'MOOSE (free)')
plt.legend(loc = 'center right')
plt.xlabel("r (m)")
plt.ylabel("displacement (mm)")
plt.title("Radial displacement around cavity")
plt.savefig("displacement_fig.png")
sys.exit(0)
Comparison with PorousFlow
The PorousFlow module is designed to simulate THM problems. The output compares favourably with Rehbinder's analytical solution, as demonstrated in the figures below.
Figure 1: Comparison between the MOOSE result (squares) and the analytic expression derived by Rehbinder for the porepressure.
The associated input file:
[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"}>>> = 40
nt<<<{"description": "Number of elements in the angular direction"}>>> = 16
rmin<<<{"description": "Inner radius. If rmin=0 then a disc mesh (with no central hole) will be created."}>>> = 0.1
rmax<<<{"description": "Outer radius"}>>> = 1
dmin<<<{"description": "Minimum degree, measured in degrees anticlockwise from x axis"}>>> = 0.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
growth_r<<<{"description": "The ratio of radial sizes of successive rings of elements"}>>> = 1.1
[]
[make3D]
input<<<{"description": "the mesh we want to extrude"}>>> = annular
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"}>>>
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
extrusion_vector<<<{"description": "The direction and length of the extrusion"}>>> = '0 0 1'
num_layers<<<{"description": "The number of layers in the extruded mesh"}>>> = 1
[]
# To get consistent ordering of results with distributed meshes
allow_renumbering = false
[]
[GlobalParams<<<{"href": "../../../../syntax/GlobalParams/index.html"}>>>]
displacements = 'disp_x disp_y disp_z'
PorousFlowDictator = dictator
biot_coefficient = 1.0
[]
[Variables<<<{"href": "../../../../syntax/Variables/index.html"}>>>]
[disp_x]
[]
[disp_y]
[]
[disp_z]
[]
[porepressure]
[]
[temperature]
[]
[]
[BCs<<<{"href": "../../../../syntax/BCs/index.html"}>>>]
[plane_strain]
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"}>>> = disp_z
value<<<{"description": "Value of the BC"}>>> = 0
boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = 'top bottom'
[]
[ymin]
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"}>>> = disp_y
value<<<{"description": "Value of the BC"}>>> = 0
boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = dmin
[]
[xmin]
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"}>>> = disp_x
value<<<{"description": "Value of the BC"}>>> = 0
boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = dmax
[]
[cavity_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"}>>> = 1000
boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = rmin
[]
[cavity_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"}>>> = 1E6
boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = rmin
[]
[cavity_zero_effective_stress_x]
type = Pressure<<<{"description": "Applies a pressure on a given boundary in a given direction", "href": "../../../../source/bcs/Pressure.html"}>>>
variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = disp_x
function<<<{"description": "The function that describes the pressure"}>>> = 1E6
boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = rmin
use_displaced_mesh<<<{"description": "Whether to use the displaced mesh."}>>> = false
[]
[cavity_zero_effective_stress_y]
type = Pressure<<<{"description": "Applies a pressure on a given boundary in a given direction", "href": "../../../../source/bcs/Pressure.html"}>>>
variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = disp_y
function<<<{"description": "The function that describes the pressure"}>>> = 1E6
boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = rmin
use_displaced_mesh<<<{"description": "Whether to use the displaced mesh."}>>> = false
[]
[outer_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"}>>> = 0
boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = rmax
[]
[outer_pressure]
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"}>>> = 0
boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = rmax
[]
[fixed_outer_x]
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"}>>> = disp_x
value<<<{"description": "Value of the BC"}>>> = 0
boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = rmax
[]
[fixed_outer_y]
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"}>>> = disp_y
value<<<{"description": "Value of the BC"}>>> = 0
boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = rmax
[]
[]
[AuxVariables<<<{"href": "../../../../syntax/AuxVariables/index.html"}>>>]
[stress_rr]
family<<<{"description": "Specifies the family of FE shape functions to use for this variable"}>>> = MONOMIAL
order<<<{"description": "Specifies the order of the FE shape function to use for this variable (additional orders not listed are allowed)"}>>> = CONSTANT
[]
[stress_pp]
family<<<{"description": "Specifies the family of FE shape functions to use for this variable"}>>> = MONOMIAL
order<<<{"description": "Specifies the order of the FE shape function to use for this variable (additional orders not listed are allowed)"}>>> = CONSTANT
[]
[]
[AuxKernels<<<{"href": "../../../../syntax/AuxKernels/index.html"}>>>]
[stress_rr]
type = RankTwoScalarAux<<<{"description": "Compute a scalar property of a RankTwoTensor", "href": "../../../../source/auxscalarkernels/RankTwoScalarAux.html"}>>>
rank_two_tensor<<<{"description": "The rank two material tensor name"}>>> = stress
variable<<<{"description": "The name of the variable that this object applies to"}>>> = stress_rr
scalar_type<<<{"description": "Type of scalar output"}>>> = RadialStress
point1<<<{"description": "Start point for axis used to calculate some cylindrical material tensor quantities"}>>> = '0 0 0'
point2<<<{"description": "End point for axis used to calculate some material tensor quantities"}>>> = '0 0 1'
[]
[stress_pp]
type = RankTwoScalarAux<<<{"description": "Compute a scalar property of a RankTwoTensor", "href": "../../../../source/auxscalarkernels/RankTwoScalarAux.html"}>>>
rank_two_tensor<<<{"description": "The rank two material tensor name"}>>> = stress
variable<<<{"description": "The name of the variable that this object applies to"}>>> = stress_pp
scalar_type<<<{"description": "Type of scalar output"}>>> = HoopStress
point1<<<{"description": "Start point for axis used to calculate some cylindrical material tensor quantities"}>>> = '0 0 0'
point2<<<{"description": "End point for axis used to calculate some material tensor quantities"}>>> = '0 0 1'
[]
[]
[FluidProperties<<<{"href": "../../../../syntax/FluidProperties/index.html"}>>>]
[the_simple_fluid]
type = SimpleFluidProperties<<<{"description": "Fluid properties for a simple fluid with a constant bulk density", "href": "../../../../source/fluidproperties/SimpleFluidProperties.html"}>>>
thermal_expansion<<<{"description": "Constant coefficient of thermal expansion (1/K)"}>>> = 0.0
bulk_modulus<<<{"description": "Constant bulk modulus (Pa)"}>>> = 1E12
viscosity<<<{"description": "Constant dynamic viscosity (Pa.s)"}>>> = 1.0E-3
density0<<<{"description": "Density at zero pressure and zero temperature"}>>> = 1000.0
cv<<<{"description": "Constant specific heat capacity at constant volume (J/kg/K)"}>>> = 1000.0
cp<<<{"description": "Constant specific heat capacity at constant pressure (J/kg/K)"}>>> = 1000.0
porepressure_coefficient<<<{"description": "The enthalpy is internal_energy + P / density * porepressure_coefficient. Physically this should be 1.0, but analytic solutions are simplified when it is zero"}>>> = 0.0
[]
[]
[PorousFlowBasicTHM<<<{"href": "../../../../syntax/PorousFlowBasicTHM/index.html"}>>>]
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"}>>> = ThermoHydroMechanical
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
add_stress_aux<<<{"description": "Add AuxVariables that record effective stress"}>>> = true
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
eigenstrain_names<<<{"description": "List of all eigenstrain models used in mechanics calculations. Typically the eigenstrain_name used in ComputeThermalExpansionEigenstrain. Only needed for thermally-coupled simulations with thermal expansion."}>>> = thermal_contribution
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
[]
[Materials<<<{"href": "../../../../syntax/Materials/index.html"}>>>]
[elasticity_tensor]
type = ComputeIsotropicElasticityTensor<<<{"description": "Compute a constant isotropic elasticity tensor.", "href": "../../../../source/materials/ComputeIsotropicElasticityTensor.html"}>>>
youngs_modulus<<<{"description": "Young's modulus of the material."}>>> = 1E10
poissons_ratio<<<{"description": "Poisson's ratio for the material."}>>> = 0.2
[]
[strain]
type = ComputeSmallStrain<<<{"description": "Compute a small strain.", "href": "../../../../source/materials/ComputeSmallStrain.html"}>>>
eigenstrain_names<<<{"description": "List of eigenstrains to be applied in this strain calculation"}>>> = thermal_contribution
[]
[thermal_contribution]
type = ComputeThermalExpansionEigenstrain<<<{"description": "Computes eigenstrain due to thermal expansion with a constant coefficient", "href": "../../../../source/materials/ComputeThermalExpansionEigenstrain.html"}>>>
temperature<<<{"description": "Coupled temperature"}>>> = temperature
thermal_expansion_coeff<<<{"description": "Thermal expansion coefficient"}>>> = 1E-6
eigenstrain_name<<<{"description": "Material property name for the eigenstrain tensor computed by this model. IMPORTANT: The name of this property must also be provided to the strain calculator."}>>> = thermal_contribution
stress_free_temperature<<<{"description": "Reference temperature at which there is no thermal expansion for thermal eigenstrain calculation"}>>> = 0.0
[]
[stress]
type = ComputeLinearElasticStress<<<{"description": "Compute stress using elasticity for small strains", "href": "../../../../source/materials/ComputeLinearElasticStress.html"}>>>
[]
[porosity]
type = PorousFlowPorosityConst<<<{"description": "This Material calculates the porosity assuming it is constant", "href": "../../../../source/materials/PorousFlowPorosityConst.html"}>>> # only the initial value of this is ever used
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
[]
[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"}>>>
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."}>>> = 1E-10
fluid_bulk_modulus<<<{"description": "Fluid bulk modulus"}>>> = 1E12
[]
[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'
[]
[thermal_expansion]
type = PorousFlowConstantThermalExpansionCoefficient<<<{"description": "Computes the effective thermal expansion coefficient, (biot_coeff - porosity) * drained_coefficient + porosity * fluid_coefficient.", "href": "../../../../source/materials/PorousFlowConstantThermalExpansionCoefficient.html"}>>>
fluid_coefficient<<<{"description": "Volumetric coefficient of thermal expansion for the fluid"}>>> = 1E-6
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)"}>>> = 1E-6
[]
[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"}>>> = '1E6 0 0 0 1E6 0 0 0 1E6'
[]
[]
[VectorPostprocessors<<<{"href": "../../../../syntax/VectorPostprocessors/index.html"}>>>]
[P]
type = LineValueSampler<<<{"description": "Samples variable(s) along a specified line", "href": "../../../../source/vectorpostprocessors/LineValueSampler.html"}>>>
start_point<<<{"description": "The beginning of the line"}>>> = '0.1 0 0'
end_point<<<{"description": "The ending of the line"}>>> = '1.0 0 0'
num_points<<<{"description": "The number of points to sample along the line"}>>> = 10
sort_by<<<{"description": "What to sort the samples by"}>>> = x
variable<<<{"description": "The names of the variables that this VectorPostprocessor operates on"}>>> = porepressure
[]
[T]
type = LineValueSampler<<<{"description": "Samples variable(s) along a specified line", "href": "../../../../source/vectorpostprocessors/LineValueSampler.html"}>>>
start_point<<<{"description": "The beginning of the line"}>>> = '0.1 0 0'
end_point<<<{"description": "The ending of the line"}>>> = '1.0 0 0'
num_points<<<{"description": "The number of points to sample along the line"}>>> = 10
sort_by<<<{"description": "What to sort the samples by"}>>> = x
variable<<<{"description": "The names of the variables that this VectorPostprocessor operates on"}>>> = temperature
[]
[U]
type = LineValueSampler<<<{"description": "Samples variable(s) along a specified line", "href": "../../../../source/vectorpostprocessors/LineValueSampler.html"}>>>
start_point<<<{"description": "The beginning of the line"}>>> = '0.1 0 0'
end_point<<<{"description": "The ending of the line"}>>> = '1.0 0 0'
num_points<<<{"description": "The number of points to sample along the line"}>>> = 10
sort_by<<<{"description": "What to sort the samples by"}>>> = x
variable<<<{"description": "The names of the variables that this VectorPostprocessor operates on"}>>> = disp_x
[]
[]
[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_iname<<<{"description": "Names of PETSc name/value pairs"}>>> = '-ksp_type -pc_type -sub_pc_type -snes_rtol'
petsc_options_value<<<{"description": "Values of PETSc name/value pairs (must correspond with \"petsc_options_iname\""}>>> = 'gmres asm lu 1E-10'
[]
[]
[Executioner<<<{"href": "../../../../syntax/Executioner/index.html"}>>>]
type = Steady
solve_type = Newton
[]
[Outputs<<<{"href": "../../../../syntax/Outputs/index.html"}>>>]
file_base<<<{"description": "Common file base name to be utilized with all output objects"}>>> = fixed_outer
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
csv<<<{"description": "Output the scalar variable and postprocessors to a *.csv file using the default CSV output."}>>> = true
[]
Figure 2: Comparison between the MOOSE result (squares) and the analytic expression derived by Rehbinder for the porepressure.
The associated input file:
[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"}>>> = 40
nt<<<{"description": "Number of elements in the angular direction"}>>> = 16
rmin<<<{"description": "Inner radius. If rmin=0 then a disc mesh (with no central hole) will be created."}>>> = 0.1
rmax<<<{"description": "Outer radius"}>>> = 1
dmin<<<{"description": "Minimum degree, measured in degrees anticlockwise from x axis"}>>> = 0.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
growth_r<<<{"description": "The ratio of radial sizes of successive rings of elements"}>>> = 1.1
[]
[make3D]
input<<<{"description": "the mesh we want to extrude"}>>> = annular
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"}>>>
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
extrusion_vector<<<{"description": "The direction and length of the extrusion"}>>> = '0 0 1'
num_layers<<<{"description": "The number of layers in the extruded mesh"}>>> = 1
[]
# To get consistent ordering of results with distributed meshes
allow_renumbering = false
[]
[GlobalParams<<<{"href": "../../../../syntax/GlobalParams/index.html"}>>>]
displacements = 'disp_x disp_y disp_z'
PorousFlowDictator = dictator
biot_coefficient = 1.0
[]
[Variables<<<{"href": "../../../../syntax/Variables/index.html"}>>>]
[disp_x]
[]
[disp_y]
[]
[disp_z]
[]
[porepressure]
[]
[temperature]
[]
[]
[BCs<<<{"href": "../../../../syntax/BCs/index.html"}>>>]
[plane_strain]
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"}>>> = disp_z
value<<<{"description": "Value of the BC"}>>> = 0
boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = 'top bottom'
[]
[ymin]
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"}>>> = disp_y
value<<<{"description": "Value of the BC"}>>> = 0
boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = dmin
[]
[xmin]
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"}>>> = disp_x
value<<<{"description": "Value of the BC"}>>> = 0
boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = dmax
[]
[cavity_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"}>>> = 1000
boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = rmin
[]
[cavity_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"}>>> = 1E6
boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = rmin
[]
[cavity_zero_effective_stress_x]
type = Pressure<<<{"description": "Applies a pressure on a given boundary in a given direction", "href": "../../../../source/bcs/Pressure.html"}>>>
variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = disp_x
function<<<{"description": "The function that describes the pressure"}>>> = 1E6
boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = rmin
use_displaced_mesh<<<{"description": "Whether to use the displaced mesh."}>>> = false
[]
[cavity_zero_effective_stress_y]
type = Pressure<<<{"description": "Applies a pressure on a given boundary in a given direction", "href": "../../../../source/bcs/Pressure.html"}>>>
variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = disp_y
function<<<{"description": "The function that describes the pressure"}>>> = 1E6
boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = rmin
use_displaced_mesh<<<{"description": "Whether to use the displaced mesh."}>>> = false
[]
[outer_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"}>>> = 0
boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = rmax
[]
[outer_pressure]
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"}>>> = 0
boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = rmax
[]
[fixed_outer_x]
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"}>>> = disp_x
value<<<{"description": "Value of the BC"}>>> = 0
boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = rmax
[]
[fixed_outer_y]
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"}>>> = disp_y
value<<<{"description": "Value of the BC"}>>> = 0
boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = rmax
[]
[]
[AuxVariables<<<{"href": "../../../../syntax/AuxVariables/index.html"}>>>]
[stress_rr]
family<<<{"description": "Specifies the family of FE shape functions to use for this variable"}>>> = MONOMIAL
order<<<{"description": "Specifies the order of the FE shape function to use for this variable (additional orders not listed are allowed)"}>>> = CONSTANT
[]
[stress_pp]
family<<<{"description": "Specifies the family of FE shape functions to use for this variable"}>>> = MONOMIAL
order<<<{"description": "Specifies the order of the FE shape function to use for this variable (additional orders not listed are allowed)"}>>> = CONSTANT
[]
[]
[AuxKernels<<<{"href": "../../../../syntax/AuxKernels/index.html"}>>>]
[stress_rr]
type = RankTwoScalarAux<<<{"description": "Compute a scalar property of a RankTwoTensor", "href": "../../../../source/auxscalarkernels/RankTwoScalarAux.html"}>>>
rank_two_tensor<<<{"description": "The rank two material tensor name"}>>> = stress
variable<<<{"description": "The name of the variable that this object applies to"}>>> = stress_rr
scalar_type<<<{"description": "Type of scalar output"}>>> = RadialStress
point1<<<{"description": "Start point for axis used to calculate some cylindrical material tensor quantities"}>>> = '0 0 0'
point2<<<{"description": "End point for axis used to calculate some material tensor quantities"}>>> = '0 0 1'
[]
[stress_pp]
type = RankTwoScalarAux<<<{"description": "Compute a scalar property of a RankTwoTensor", "href": "../../../../source/auxscalarkernels/RankTwoScalarAux.html"}>>>
rank_two_tensor<<<{"description": "The rank two material tensor name"}>>> = stress
variable<<<{"description": "The name of the variable that this object applies to"}>>> = stress_pp
scalar_type<<<{"description": "Type of scalar output"}>>> = HoopStress
point1<<<{"description": "Start point for axis used to calculate some cylindrical material tensor quantities"}>>> = '0 0 0'
point2<<<{"description": "End point for axis used to calculate some material tensor quantities"}>>> = '0 0 1'
[]
[]
[FluidProperties<<<{"href": "../../../../syntax/FluidProperties/index.html"}>>>]
[the_simple_fluid]
type = SimpleFluidProperties<<<{"description": "Fluid properties for a simple fluid with a constant bulk density", "href": "../../../../source/fluidproperties/SimpleFluidProperties.html"}>>>
thermal_expansion<<<{"description": "Constant coefficient of thermal expansion (1/K)"}>>> = 0.0
bulk_modulus<<<{"description": "Constant bulk modulus (Pa)"}>>> = 1E12
viscosity<<<{"description": "Constant dynamic viscosity (Pa.s)"}>>> = 1.0E-3
density0<<<{"description": "Density at zero pressure and zero temperature"}>>> = 1000.0
cv<<<{"description": "Constant specific heat capacity at constant volume (J/kg/K)"}>>> = 1000.0
cp<<<{"description": "Constant specific heat capacity at constant pressure (J/kg/K)"}>>> = 1000.0
porepressure_coefficient<<<{"description": "The enthalpy is internal_energy + P / density * porepressure_coefficient. Physically this should be 1.0, but analytic solutions are simplified when it is zero"}>>> = 0.0
[]
[]
[PorousFlowBasicTHM<<<{"href": "../../../../syntax/PorousFlowBasicTHM/index.html"}>>>]
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"}>>> = ThermoHydroMechanical
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
add_stress_aux<<<{"description": "Add AuxVariables that record effective stress"}>>> = true
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
eigenstrain_names<<<{"description": "List of all eigenstrain models used in mechanics calculations. Typically the eigenstrain_name used in ComputeThermalExpansionEigenstrain. Only needed for thermally-coupled simulations with thermal expansion."}>>> = thermal_contribution
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
[]
[Materials<<<{"href": "../../../../syntax/Materials/index.html"}>>>]
[elasticity_tensor]
type = ComputeIsotropicElasticityTensor<<<{"description": "Compute a constant isotropic elasticity tensor.", "href": "../../../../source/materials/ComputeIsotropicElasticityTensor.html"}>>>
youngs_modulus<<<{"description": "Young's modulus of the material."}>>> = 1E10
poissons_ratio<<<{"description": "Poisson's ratio for the material."}>>> = 0.2
[]
[strain]
type = ComputeSmallStrain<<<{"description": "Compute a small strain.", "href": "../../../../source/materials/ComputeSmallStrain.html"}>>>
eigenstrain_names<<<{"description": "List of eigenstrains to be applied in this strain calculation"}>>> = thermal_contribution
[]
[thermal_contribution]
type = ComputeThermalExpansionEigenstrain<<<{"description": "Computes eigenstrain due to thermal expansion with a constant coefficient", "href": "../../../../source/materials/ComputeThermalExpansionEigenstrain.html"}>>>
temperature<<<{"description": "Coupled temperature"}>>> = temperature
thermal_expansion_coeff<<<{"description": "Thermal expansion coefficient"}>>> = 1E-6
eigenstrain_name<<<{"description": "Material property name for the eigenstrain tensor computed by this model. IMPORTANT: The name of this property must also be provided to the strain calculator."}>>> = thermal_contribution
stress_free_temperature<<<{"description": "Reference temperature at which there is no thermal expansion for thermal eigenstrain calculation"}>>> = 0.0
[]
[stress]
type = ComputeLinearElasticStress<<<{"description": "Compute stress using elasticity for small strains", "href": "../../../../source/materials/ComputeLinearElasticStress.html"}>>>
[]
[porosity]
type = PorousFlowPorosityConst<<<{"description": "This Material calculates the porosity assuming it is constant", "href": "../../../../source/materials/PorousFlowPorosityConst.html"}>>> # only the initial value of this is ever used
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
[]
[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"}>>>
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."}>>> = 1E-10
fluid_bulk_modulus<<<{"description": "Fluid bulk modulus"}>>> = 1E12
[]
[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'
[]
[thermal_expansion]
type = PorousFlowConstantThermalExpansionCoefficient<<<{"description": "Computes the effective thermal expansion coefficient, (biot_coeff - porosity) * drained_coefficient + porosity * fluid_coefficient.", "href": "../../../../source/materials/PorousFlowConstantThermalExpansionCoefficient.html"}>>>
fluid_coefficient<<<{"description": "Volumetric coefficient of thermal expansion for the fluid"}>>> = 1E-6
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)"}>>> = 1E-6
[]
[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"}>>> = '1E6 0 0 0 1E6 0 0 0 1E6'
[]
[]
[VectorPostprocessors<<<{"href": "../../../../syntax/VectorPostprocessors/index.html"}>>>]
[P]
type = LineValueSampler<<<{"description": "Samples variable(s) along a specified line", "href": "../../../../source/vectorpostprocessors/LineValueSampler.html"}>>>
start_point<<<{"description": "The beginning of the line"}>>> = '0.1 0 0'
end_point<<<{"description": "The ending of the line"}>>> = '1.0 0 0'
num_points<<<{"description": "The number of points to sample along the line"}>>> = 10
sort_by<<<{"description": "What to sort the samples by"}>>> = x
variable<<<{"description": "The names of the variables that this VectorPostprocessor operates on"}>>> = porepressure
[]
[T]
type = LineValueSampler<<<{"description": "Samples variable(s) along a specified line", "href": "../../../../source/vectorpostprocessors/LineValueSampler.html"}>>>
start_point<<<{"description": "The beginning of the line"}>>> = '0.1 0 0'
end_point<<<{"description": "The ending of the line"}>>> = '1.0 0 0'
num_points<<<{"description": "The number of points to sample along the line"}>>> = 10
sort_by<<<{"description": "What to sort the samples by"}>>> = x
variable<<<{"description": "The names of the variables that this VectorPostprocessor operates on"}>>> = temperature
[]
[U]
type = LineValueSampler<<<{"description": "Samples variable(s) along a specified line", "href": "../../../../source/vectorpostprocessors/LineValueSampler.html"}>>>
start_point<<<{"description": "The beginning of the line"}>>> = '0.1 0 0'
end_point<<<{"description": "The ending of the line"}>>> = '1.0 0 0'
num_points<<<{"description": "The number of points to sample along the line"}>>> = 10
sort_by<<<{"description": "What to sort the samples by"}>>> = x
variable<<<{"description": "The names of the variables that this VectorPostprocessor operates on"}>>> = disp_x
[]
[]
[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_iname<<<{"description": "Names of PETSc name/value pairs"}>>> = '-ksp_type -pc_type -sub_pc_type -snes_rtol'
petsc_options_value<<<{"description": "Values of PETSc name/value pairs (must correspond with \"petsc_options_iname\""}>>> = 'gmres asm lu 1E-10'
[]
[]
[Executioner<<<{"href": "../../../../syntax/Executioner/index.html"}>>>]
type = Steady
solve_type = Newton
[]
[Outputs<<<{"href": "../../../../syntax/Outputs/index.html"}>>>]
file_base<<<{"description": "Common file base name to be utilized with all output objects"}>>> = fixed_outer
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
csv<<<{"description": "Output the scalar variable and postprocessors to a *.csv file using the default CSV output."}>>> = true
[]
Figure 3: Comparison between the MOOSE result (squares) and the analytic expression derived by Rehbinder for the temperature.
The associated input file for the fixed case:
[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"}>>> = 40
nt<<<{"description": "Number of elements in the angular direction"}>>> = 16
rmin<<<{"description": "Inner radius. If rmin=0 then a disc mesh (with no central hole) will be created."}>>> = 0.1
rmax<<<{"description": "Outer radius"}>>> = 1
dmin<<<{"description": "Minimum degree, measured in degrees anticlockwise from x axis"}>>> = 0.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
growth_r<<<{"description": "The ratio of radial sizes of successive rings of elements"}>>> = 1.1
[]
[make3D]
input<<<{"description": "the mesh we want to extrude"}>>> = annular
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"}>>>
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
extrusion_vector<<<{"description": "The direction and length of the extrusion"}>>> = '0 0 1'
num_layers<<<{"description": "The number of layers in the extruded mesh"}>>> = 1
[]
# To get consistent ordering of results with distributed meshes
allow_renumbering = false
[]
[GlobalParams<<<{"href": "../../../../syntax/GlobalParams/index.html"}>>>]
displacements = 'disp_x disp_y disp_z'
PorousFlowDictator = dictator
biot_coefficient = 1.0
[]
[Variables<<<{"href": "../../../../syntax/Variables/index.html"}>>>]
[disp_x]
[]
[disp_y]
[]
[disp_z]
[]
[porepressure]
[]
[temperature]
[]
[]
[BCs<<<{"href": "../../../../syntax/BCs/index.html"}>>>]
[plane_strain]
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"}>>> = disp_z
value<<<{"description": "Value of the BC"}>>> = 0
boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = 'top bottom'
[]
[ymin]
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"}>>> = disp_y
value<<<{"description": "Value of the BC"}>>> = 0
boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = dmin
[]
[xmin]
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"}>>> = disp_x
value<<<{"description": "Value of the BC"}>>> = 0
boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = dmax
[]
[cavity_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"}>>> = 1000
boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = rmin
[]
[cavity_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"}>>> = 1E6
boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = rmin
[]
[cavity_zero_effective_stress_x]
type = Pressure<<<{"description": "Applies a pressure on a given boundary in a given direction", "href": "../../../../source/bcs/Pressure.html"}>>>
variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = disp_x
function<<<{"description": "The function that describes the pressure"}>>> = 1E6
boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = rmin
use_displaced_mesh<<<{"description": "Whether to use the displaced mesh."}>>> = false
[]
[cavity_zero_effective_stress_y]
type = Pressure<<<{"description": "Applies a pressure on a given boundary in a given direction", "href": "../../../../source/bcs/Pressure.html"}>>>
variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = disp_y
function<<<{"description": "The function that describes the pressure"}>>> = 1E6
boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = rmin
use_displaced_mesh<<<{"description": "Whether to use the displaced mesh."}>>> = false
[]
[outer_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"}>>> = 0
boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = rmax
[]
[outer_pressure]
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"}>>> = 0
boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = rmax
[]
[fixed_outer_x]
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"}>>> = disp_x
value<<<{"description": "Value of the BC"}>>> = 0
boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = rmax
[]
[fixed_outer_y]
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"}>>> = disp_y
value<<<{"description": "Value of the BC"}>>> = 0
boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = rmax
[]
[]
[AuxVariables<<<{"href": "../../../../syntax/AuxVariables/index.html"}>>>]
[stress_rr]
family<<<{"description": "Specifies the family of FE shape functions to use for this variable"}>>> = MONOMIAL
order<<<{"description": "Specifies the order of the FE shape function to use for this variable (additional orders not listed are allowed)"}>>> = CONSTANT
[]
[stress_pp]
family<<<{"description": "Specifies the family of FE shape functions to use for this variable"}>>> = MONOMIAL
order<<<{"description": "Specifies the order of the FE shape function to use for this variable (additional orders not listed are allowed)"}>>> = CONSTANT
[]
[]
[AuxKernels<<<{"href": "../../../../syntax/AuxKernels/index.html"}>>>]
[stress_rr]
type = RankTwoScalarAux<<<{"description": "Compute a scalar property of a RankTwoTensor", "href": "../../../../source/auxscalarkernels/RankTwoScalarAux.html"}>>>
rank_two_tensor<<<{"description": "The rank two material tensor name"}>>> = stress
variable<<<{"description": "The name of the variable that this object applies to"}>>> = stress_rr
scalar_type<<<{"description": "Type of scalar output"}>>> = RadialStress
point1<<<{"description": "Start point for axis used to calculate some cylindrical material tensor quantities"}>>> = '0 0 0'
point2<<<{"description": "End point for axis used to calculate some material tensor quantities"}>>> = '0 0 1'
[]
[stress_pp]
type = RankTwoScalarAux<<<{"description": "Compute a scalar property of a RankTwoTensor", "href": "../../../../source/auxscalarkernels/RankTwoScalarAux.html"}>>>
rank_two_tensor<<<{"description": "The rank two material tensor name"}>>> = stress
variable<<<{"description": "The name of the variable that this object applies to"}>>> = stress_pp
scalar_type<<<{"description": "Type of scalar output"}>>> = HoopStress
point1<<<{"description": "Start point for axis used to calculate some cylindrical material tensor quantities"}>>> = '0 0 0'
point2<<<{"description": "End point for axis used to calculate some material tensor quantities"}>>> = '0 0 1'
[]
[]
[FluidProperties<<<{"href": "../../../../syntax/FluidProperties/index.html"}>>>]
[the_simple_fluid]
type = SimpleFluidProperties<<<{"description": "Fluid properties for a simple fluid with a constant bulk density", "href": "../../../../source/fluidproperties/SimpleFluidProperties.html"}>>>
thermal_expansion<<<{"description": "Constant coefficient of thermal expansion (1/K)"}>>> = 0.0
bulk_modulus<<<{"description": "Constant bulk modulus (Pa)"}>>> = 1E12
viscosity<<<{"description": "Constant dynamic viscosity (Pa.s)"}>>> = 1.0E-3
density0<<<{"description": "Density at zero pressure and zero temperature"}>>> = 1000.0
cv<<<{"description": "Constant specific heat capacity at constant volume (J/kg/K)"}>>> = 1000.0
cp<<<{"description": "Constant specific heat capacity at constant pressure (J/kg/K)"}>>> = 1000.0
porepressure_coefficient<<<{"description": "The enthalpy is internal_energy + P / density * porepressure_coefficient. Physically this should be 1.0, but analytic solutions are simplified when it is zero"}>>> = 0.0
[]
[]
[PorousFlowBasicTHM<<<{"href": "../../../../syntax/PorousFlowBasicTHM/index.html"}>>>]
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"}>>> = ThermoHydroMechanical
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
add_stress_aux<<<{"description": "Add AuxVariables that record effective stress"}>>> = true
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
eigenstrain_names<<<{"description": "List of all eigenstrain models used in mechanics calculations. Typically the eigenstrain_name used in ComputeThermalExpansionEigenstrain. Only needed for thermally-coupled simulations with thermal expansion."}>>> = thermal_contribution
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
[]
[Materials<<<{"href": "../../../../syntax/Materials/index.html"}>>>]
[elasticity_tensor]
type = ComputeIsotropicElasticityTensor<<<{"description": "Compute a constant isotropic elasticity tensor.", "href": "../../../../source/materials/ComputeIsotropicElasticityTensor.html"}>>>
youngs_modulus<<<{"description": "Young's modulus of the material."}>>> = 1E10
poissons_ratio<<<{"description": "Poisson's ratio for the material."}>>> = 0.2
[]
[strain]
type = ComputeSmallStrain<<<{"description": "Compute a small strain.", "href": "../../../../source/materials/ComputeSmallStrain.html"}>>>
eigenstrain_names<<<{"description": "List of eigenstrains to be applied in this strain calculation"}>>> = thermal_contribution
[]
[thermal_contribution]
type = ComputeThermalExpansionEigenstrain<<<{"description": "Computes eigenstrain due to thermal expansion with a constant coefficient", "href": "../../../../source/materials/ComputeThermalExpansionEigenstrain.html"}>>>
temperature<<<{"description": "Coupled temperature"}>>> = temperature
thermal_expansion_coeff<<<{"description": "Thermal expansion coefficient"}>>> = 1E-6
eigenstrain_name<<<{"description": "Material property name for the eigenstrain tensor computed by this model. IMPORTANT: The name of this property must also be provided to the strain calculator."}>>> = thermal_contribution
stress_free_temperature<<<{"description": "Reference temperature at which there is no thermal expansion for thermal eigenstrain calculation"}>>> = 0.0
[]
[stress]
type = ComputeLinearElasticStress<<<{"description": "Compute stress using elasticity for small strains", "href": "../../../../source/materials/ComputeLinearElasticStress.html"}>>>
[]
[porosity]
type = PorousFlowPorosityConst<<<{"description": "This Material calculates the porosity assuming it is constant", "href": "../../../../source/materials/PorousFlowPorosityConst.html"}>>> # only the initial value of this is ever used
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
[]
[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"}>>>
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."}>>> = 1E-10
fluid_bulk_modulus<<<{"description": "Fluid bulk modulus"}>>> = 1E12
[]
[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'
[]
[thermal_expansion]
type = PorousFlowConstantThermalExpansionCoefficient<<<{"description": "Computes the effective thermal expansion coefficient, (biot_coeff - porosity) * drained_coefficient + porosity * fluid_coefficient.", "href": "../../../../source/materials/PorousFlowConstantThermalExpansionCoefficient.html"}>>>
fluid_coefficient<<<{"description": "Volumetric coefficient of thermal expansion for the fluid"}>>> = 1E-6
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)"}>>> = 1E-6
[]
[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"}>>> = '1E6 0 0 0 1E6 0 0 0 1E6'
[]
[]
[VectorPostprocessors<<<{"href": "../../../../syntax/VectorPostprocessors/index.html"}>>>]
[P]
type = LineValueSampler<<<{"description": "Samples variable(s) along a specified line", "href": "../../../../source/vectorpostprocessors/LineValueSampler.html"}>>>
start_point<<<{"description": "The beginning of the line"}>>> = '0.1 0 0'
end_point<<<{"description": "The ending of the line"}>>> = '1.0 0 0'
num_points<<<{"description": "The number of points to sample along the line"}>>> = 10
sort_by<<<{"description": "What to sort the samples by"}>>> = x
variable<<<{"description": "The names of the variables that this VectorPostprocessor operates on"}>>> = porepressure
[]
[T]
type = LineValueSampler<<<{"description": "Samples variable(s) along a specified line", "href": "../../../../source/vectorpostprocessors/LineValueSampler.html"}>>>
start_point<<<{"description": "The beginning of the line"}>>> = '0.1 0 0'
end_point<<<{"description": "The ending of the line"}>>> = '1.0 0 0'
num_points<<<{"description": "The number of points to sample along the line"}>>> = 10
sort_by<<<{"description": "What to sort the samples by"}>>> = x
variable<<<{"description": "The names of the variables that this VectorPostprocessor operates on"}>>> = temperature
[]
[U]
type = LineValueSampler<<<{"description": "Samples variable(s) along a specified line", "href": "../../../../source/vectorpostprocessors/LineValueSampler.html"}>>>
start_point<<<{"description": "The beginning of the line"}>>> = '0.1 0 0'
end_point<<<{"description": "The ending of the line"}>>> = '1.0 0 0'
num_points<<<{"description": "The number of points to sample along the line"}>>> = 10
sort_by<<<{"description": "What to sort the samples by"}>>> = x
variable<<<{"description": "The names of the variables that this VectorPostprocessor operates on"}>>> = disp_x
[]
[]
[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_iname<<<{"description": "Names of PETSc name/value pairs"}>>> = '-ksp_type -pc_type -sub_pc_type -snes_rtol'
petsc_options_value<<<{"description": "Values of PETSc name/value pairs (must correspond with \"petsc_options_iname\""}>>> = 'gmres asm lu 1E-10'
[]
[]
[Executioner<<<{"href": "../../../../syntax/Executioner/index.html"}>>>]
type = Steady
solve_type = Newton
[]
[Outputs<<<{"href": "../../../../syntax/Outputs/index.html"}>>>]
file_base<<<{"description": "Common file base name to be utilized with all output objects"}>>> = fixed_outer
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
csv<<<{"description": "Output the scalar variable and postprocessors to a *.csv file using the default CSV output."}>>> = true
[]The associated input file for the free case:
[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"}>>> = 40
nt<<<{"description": "Number of elements in the angular direction"}>>> = 16
rmin<<<{"description": "Inner radius. If rmin=0 then a disc mesh (with no central hole) will be created."}>>> = 0.1
rmax<<<{"description": "Outer radius"}>>> = 1
dmin<<<{"description": "Minimum degree, measured in degrees anticlockwise from x axis"}>>> = 0.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
growth_r<<<{"description": "The ratio of radial sizes of successive rings of elements"}>>> = 1.1
[]
[make3D]
input<<<{"description": "the mesh we want to extrude"}>>> = annular
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"}>>>
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
extrusion_vector<<<{"description": "The direction and length of the extrusion"}>>> = '0 0 1'
num_layers<<<{"description": "The number of layers in the extruded mesh"}>>> = 1
[]
[]
[GlobalParams<<<{"href": "../../../../syntax/GlobalParams/index.html"}>>>]
displacements = 'disp_x disp_y disp_z'
PorousFlowDictator = dictator
biot_coefficient = 1.0
[]
[Variables<<<{"href": "../../../../syntax/Variables/index.html"}>>>]
[disp_x]
[]
[disp_y]
[]
[disp_z]
[]
[porepressure]
[]
[temperature]
[]
[]
[BCs<<<{"href": "../../../../syntax/BCs/index.html"}>>>]
# sideset 1 = outer
# sideset 2 = cavity
# sideset 3 = ymin
# sideset 4 = xmin
[plane_strain]
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"}>>> = disp_z
value<<<{"description": "Value of the BC"}>>> = 0
boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = 'top bottom'
[]
[ymin]
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"}>>> = disp_y
value<<<{"description": "Value of the BC"}>>> = 0
boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = dmin
[]
[xmin]
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"}>>> = disp_x
value<<<{"description": "Value of the BC"}>>> = 0
boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = dmax
[]
[cavity_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"}>>> = 1000
boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = rmin
[]
[cavity_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"}>>> = 1E6
boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = rmin
[]
[cavity_zero_effective_stress_x]
type = Pressure<<<{"description": "Applies a pressure on a given boundary in a given direction", "href": "../../../../source/bcs/Pressure.html"}>>>
variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = disp_x
function<<<{"description": "The function that describes the pressure"}>>> = 1E6
boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = rmin
use_displaced_mesh<<<{"description": "Whether to use the displaced mesh."}>>> = false
[]
[cavity_zero_effective_stress_y]
type = Pressure<<<{"description": "Applies a pressure on a given boundary in a given direction", "href": "../../../../source/bcs/Pressure.html"}>>>
variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = disp_y
function<<<{"description": "The function that describes the pressure"}>>> = 1E6
boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = rmin
use_displaced_mesh<<<{"description": "Whether to use the displaced mesh."}>>> = false
[]
[outer_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"}>>> = 0
boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = rmax
[]
[outer_pressure]
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"}>>> = 0
boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = rmax
[]
[]
[AuxVariables<<<{"href": "../../../../syntax/AuxVariables/index.html"}>>>]
[stress_rr]
family<<<{"description": "Specifies the family of FE shape functions to use for this variable"}>>> = MONOMIAL
order<<<{"description": "Specifies the order of the FE shape function to use for this variable (additional orders not listed are allowed)"}>>> = CONSTANT
[]
[stress_pp]
family<<<{"description": "Specifies the family of FE shape functions to use for this variable"}>>> = MONOMIAL
order<<<{"description": "Specifies the order of the FE shape function to use for this variable (additional orders not listed are allowed)"}>>> = CONSTANT
[]
[]
[AuxKernels<<<{"href": "../../../../syntax/AuxKernels/index.html"}>>>]
[stress_rr]
type = RankTwoScalarAux<<<{"description": "Compute a scalar property of a RankTwoTensor", "href": "../../../../source/auxscalarkernels/RankTwoScalarAux.html"}>>>
rank_two_tensor<<<{"description": "The rank two material tensor name"}>>> = stress
variable<<<{"description": "The name of the variable that this object applies to"}>>> = stress_rr
scalar_type<<<{"description": "Type of scalar output"}>>> = RadialStress
point1<<<{"description": "Start point for axis used to calculate some cylindrical material tensor quantities"}>>> = '0 0 0'
point2<<<{"description": "End point for axis used to calculate some material tensor quantities"}>>> = '0 0 1'
[]
[stress_pp]
type = RankTwoScalarAux<<<{"description": "Compute a scalar property of a RankTwoTensor", "href": "../../../../source/auxscalarkernels/RankTwoScalarAux.html"}>>>
rank_two_tensor<<<{"description": "The rank two material tensor name"}>>> = stress
variable<<<{"description": "The name of the variable that this object applies to"}>>> = stress_pp
scalar_type<<<{"description": "Type of scalar output"}>>> = HoopStress
point1<<<{"description": "Start point for axis used to calculate some cylindrical material tensor quantities"}>>> = '0 0 0'
point2<<<{"description": "End point for axis used to calculate some material tensor quantities"}>>> = '0 0 1'
[]
[]
[FluidProperties<<<{"href": "../../../../syntax/FluidProperties/index.html"}>>>]
[the_simple_fluid]
type = SimpleFluidProperties<<<{"description": "Fluid properties for a simple fluid with a constant bulk density", "href": "../../../../source/fluidproperties/SimpleFluidProperties.html"}>>>
thermal_expansion<<<{"description": "Constant coefficient of thermal expansion (1/K)"}>>> = 0.0
bulk_modulus<<<{"description": "Constant bulk modulus (Pa)"}>>> = 1E12
viscosity<<<{"description": "Constant dynamic viscosity (Pa.s)"}>>> = 1.0E-3
density0<<<{"description": "Density at zero pressure and zero temperature"}>>> = 1000.0
cv<<<{"description": "Constant specific heat capacity at constant volume (J/kg/K)"}>>> = 1000.0
cp<<<{"description": "Constant specific heat capacity at constant pressure (J/kg/K)"}>>> = 1000.0
porepressure_coefficient<<<{"description": "The enthalpy is internal_energy + P / density * porepressure_coefficient. Physically this should be 1.0, but analytic solutions are simplified when it is zero"}>>> = 0.0
[]
[]
[PorousFlowBasicTHM<<<{"href": "../../../../syntax/PorousFlowBasicTHM/index.html"}>>>]
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"}>>> = ThermoHydroMechanical
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
add_stress_aux<<<{"description": "Add AuxVariables that record effective stress"}>>> = true
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
eigenstrain_names<<<{"description": "List of all eigenstrain models used in mechanics calculations. Typically the eigenstrain_name used in ComputeThermalExpansionEigenstrain. Only needed for thermally-coupled simulations with thermal expansion."}>>> = thermal_contribution
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
[]
[Materials<<<{"href": "../../../../syntax/Materials/index.html"}>>>]
[elasticity_tensor]
type = ComputeIsotropicElasticityTensor<<<{"description": "Compute a constant isotropic elasticity tensor.", "href": "../../../../source/materials/ComputeIsotropicElasticityTensor.html"}>>>
youngs_modulus<<<{"description": "Young's modulus of the material."}>>> = 1E10
poissons_ratio<<<{"description": "Poisson's ratio for the material."}>>> = 0.2
[]
[strain]
type = ComputeSmallStrain<<<{"description": "Compute a small strain.", "href": "../../../../source/materials/ComputeSmallStrain.html"}>>>
eigenstrain_names<<<{"description": "List of eigenstrains to be applied in this strain calculation"}>>> = thermal_contribution
[]
[thermal_contribution]
type = ComputeThermalExpansionEigenstrain<<<{"description": "Computes eigenstrain due to thermal expansion with a constant coefficient", "href": "../../../../source/materials/ComputeThermalExpansionEigenstrain.html"}>>>
temperature<<<{"description": "Coupled temperature"}>>> = temperature
thermal_expansion_coeff<<<{"description": "Thermal expansion coefficient"}>>> = 1E-6
eigenstrain_name<<<{"description": "Material property name for the eigenstrain tensor computed by this model. IMPORTANT: The name of this property must also be provided to the strain calculator."}>>> = thermal_contribution
stress_free_temperature<<<{"description": "Reference temperature at which there is no thermal expansion for thermal eigenstrain calculation"}>>> = 0.0
[]
[stress]
type = ComputeLinearElasticStress<<<{"description": "Compute stress using elasticity for small strains", "href": "../../../../source/materials/ComputeLinearElasticStress.html"}>>>
[]
[porosity]
type = PorousFlowPorosityConst<<<{"description": "This Material calculates the porosity assuming it is constant", "href": "../../../../source/materials/PorousFlowPorosityConst.html"}>>> # only the initial value of this is ever used
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
[]
[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"}>>>
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."}>>> = 1E-10
fluid_bulk_modulus<<<{"description": "Fluid bulk modulus"}>>> = 1E12
[]
[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'
[]
[thermal_expansion]
type = PorousFlowConstantThermalExpansionCoefficient<<<{"description": "Computes the effective thermal expansion coefficient, (biot_coeff - porosity) * drained_coefficient + porosity * fluid_coefficient.", "href": "../../../../source/materials/PorousFlowConstantThermalExpansionCoefficient.html"}>>>
fluid_coefficient<<<{"description": "Volumetric coefficient of thermal expansion for the fluid"}>>> = 1E-6
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)"}>>> = 1E-6
[]
[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"}>>> = '1E6 0 0 0 1E6 0 0 0 1E6'
[]
[]
[VectorPostprocessors<<<{"href": "../../../../syntax/VectorPostprocessors/index.html"}>>>]
[P]
type = LineValueSampler<<<{"description": "Samples variable(s) along a specified line", "href": "../../../../source/vectorpostprocessors/LineValueSampler.html"}>>>
start_point<<<{"description": "The beginning of the line"}>>> = '0.1 0 0'
end_point<<<{"description": "The ending of the line"}>>> = '1.0 0 0'
num_points<<<{"description": "The number of points to sample along the line"}>>> = 10
sort_by<<<{"description": "What to sort the samples by"}>>> = x
variable<<<{"description": "The names of the variables that this VectorPostprocessor operates on"}>>> = porepressure
[]
[T]
type = LineValueSampler<<<{"description": "Samples variable(s) along a specified line", "href": "../../../../source/vectorpostprocessors/LineValueSampler.html"}>>>
start_point<<<{"description": "The beginning of the line"}>>> = '0.1 0 0'
end_point<<<{"description": "The ending of the line"}>>> = '1.0 0 0'
num_points<<<{"description": "The number of points to sample along the line"}>>> = 10
sort_by<<<{"description": "What to sort the samples by"}>>> = x
variable<<<{"description": "The names of the variables that this VectorPostprocessor operates on"}>>> = temperature
[]
[U]
type = LineValueSampler<<<{"description": "Samples variable(s) along a specified line", "href": "../../../../source/vectorpostprocessors/LineValueSampler.html"}>>>
start_point<<<{"description": "The beginning of the line"}>>> = '0.1 0 0'
end_point<<<{"description": "The ending of the line"}>>> = '1.0 0 0'
num_points<<<{"description": "The number of points to sample along the line"}>>> = 10
sort_by<<<{"description": "What to sort the samples by"}>>> = x
variable<<<{"description": "The names of the variables that this VectorPostprocessor operates on"}>>> = disp_x
[]
[]
[Preconditioning<<<{"href": "../../../../syntax/Preconditioning/index.html"}>>>]
[andy]
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"}>>> = '-ksp_type -pc_type -sub_pc_type -snes_rtol'
petsc_options_value<<<{"description": "Values of PETSc name/value pairs (must correspond with \"petsc_options_iname\""}>>> = 'gmres asm lu 1E-8'
[]
[]
[Executioner<<<{"href": "../../../../syntax/Executioner/index.html"}>>>]
type = Steady
solve_type = Newton
[]
[Outputs<<<{"href": "../../../../syntax/Outputs/index.html"}>>>]
file_base<<<{"description": "Common file base name to be utilized with all output objects"}>>> = free_outer
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
csv<<<{"description": "Output the scalar variable and postprocessors to a *.csv file using the default CSV output."}>>> = true
[]References
- G. Rehbinder.
Analytical solutions of stationary coupled thermo-hydro-mechanical problems.
International Journal of Rock Mechanics and Mining Sciences & Geomechanics Abstracts, 32(5):453–463, 1995.
Thermo-Hydro-Mechanical Coupling in Rock Mechanics.
URL: https://www.sciencedirect.com/science/article/pii/014890629500035F, doi:https://doi.org/10.1016/0148-9062(95)00035-F.[BibTeX]