Step 1 - First Thermal Mortar Contact
Continuing from the final step in the mechanical contact tutorial, for the first time we incorporate thermal contact, staying within the mortar method.
#
# A first attempt at thermo mechanical contact
# https://mooseframework.inl.gov/modules/combined/tutorials/introduction/step01.html
#
[GlobalParams<<<{"href": "../../../../syntax/GlobalParams/index.html"}>>>]
displacements = 'disp_x disp_y'
block = 0
[]
[Mesh<<<{"href": "../../../../syntax/Mesh/index.html"}>>>]
[generated1]
type = GeneratedMeshGenerator<<<{"description": "Create a line, square, or cube mesh with uniformly spaced or biased elements.", "href": "../../../../source/meshgenerators/GeneratedMeshGenerator.html"}>>>
dim<<<{"description": "The dimension of the mesh to be generated"}>>> = 2
nx<<<{"description": "Number of elements in the X direction"}>>> = 5
ny<<<{"description": "Number of elements in the Y direction"}>>> = 15
xmin<<<{"description": "Lower X Coordinate of the generated mesh"}>>> = -0.6
xmax<<<{"description": "Upper X Coordinate of the generated mesh"}>>> = -0.1
ymax<<<{"description": "Upper Y Coordinate of the generated mesh"}>>> = 5
bias_y<<<{"description": "The amount by which to grow (or shrink) the cells in the y-direction."}>>> = 0.9
boundary_name_prefix<<<{"description": "If provided, prefix the built in boundary names with this string"}>>> = pillar1
[]
[generated2]
type = GeneratedMeshGenerator<<<{"description": "Create a line, square, or cube mesh with uniformly spaced or biased elements.", "href": "../../../../source/meshgenerators/GeneratedMeshGenerator.html"}>>>
dim<<<{"description": "The dimension of the mesh to be generated"}>>> = 2
nx<<<{"description": "Number of elements in the X direction"}>>> = 6
ny<<<{"description": "Number of elements in the Y direction"}>>> = 15
xmin<<<{"description": "Lower X Coordinate of the generated mesh"}>>> = 0.1
xmax<<<{"description": "Upper X Coordinate of the generated mesh"}>>> = 0.6
ymax<<<{"description": "Upper Y Coordinate of the generated mesh"}>>> = 4.999
bias_y<<<{"description": "The amount by which to grow (or shrink) the cells in the y-direction."}>>> = 0.9
boundary_name_prefix<<<{"description": "If provided, prefix the built in boundary names with this string"}>>> = pillar2
boundary_id_offset<<<{"description": "This offset is added to the generated boundary IDs"}>>> = 4
[]
[collect_meshes]
type = MeshCollectionGenerator<<<{"description": "Collects multiple meshes into a single (unconnected) mesh.", "href": "../../../../source/meshgenerators/MeshCollectionGenerator.html"}>>>
inputs<<<{"description": "The input MeshGenerators."}>>> = 'generated1 generated2'
[]
patch_update_strategy = iteration
[]
[Variables<<<{"href": "../../../../syntax/Variables/index.html"}>>>]
# temperature field variable
[T]
# initialize to an average temperature
initial_condition<<<{"description": "Specifies a constant initial condition for this variable"}>>> = 50
order<<<{"description": "Specifies the order of the FE shape function to use for this variable (additional orders not listed are allowed)"}>>> = FIRST
family<<<{"description": "Specifies the family of FE shape functions to use for this variable"}>>> = LAGRANGE
[]
# temperature lagrange multiplier
[Tlm]
block = 'pillars_secondary_subdomain'
order<<<{"description": "Specifies the order of the FE shape function to use for this variable (additional orders not listed are allowed)"}>>> = FIRST
family<<<{"description": "Specifies the family of FE shape functions to use for this variable"}>>> = LAGRANGE
[]
[]
[Kernels<<<{"href": "../../../../syntax/Kernels/index.html"}>>>]
[heat_conduction]
type = HeatConduction<<<{"description": "Diffusive heat conduction term $-\\nabla\\cdot(k\\nabla T)$ of the thermal energy conservation equation", "href": "../../../../source/kernels/HeatConduction.html"}>>>
variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = T
[]
[dTdt]
type = HeatConductionTimeDerivative<<<{"description": "Time derivative term $\\rho c_p \\frac{\\partial T}{\\partial t}$ of the thermal energy conservation equation.", "href": "../../../../source/kernels/HeatConductionTimeDerivative.html"}>>>
variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = T
[]
[]
[Physics<<<{"href": "../../../../syntax/Physics/index.html"}>>>/SolidMechanics<<<{"href": "../../../../syntax/Physics/SolidMechanics/index.html"}>>>/QuasiStatic<<<{"href": "../../../../syntax/Physics/SolidMechanics/QuasiStatic/index.html"}>>>]
[all]
add_variables<<<{"description": "Add the displacement variables"}>>> = true
strain<<<{"description": "Strain formulation"}>>> = FINITE
generate_output<<<{"description": "Add scalar quantity output for stress and/or strain"}>>> = 'vonmises_stress'
[]
[]
[Contact<<<{"href": "../../../../syntax/Contact/index.html"}>>>]
[pillars]
primary<<<{"description": "The list of boundary IDs referring to primary sidesets"}>>> = pillar1_right
secondary<<<{"description": "The list of boundary IDs referring to secondary sidesets"}>>> = pillar2_left
model<<<{"description": "The contact model to use"}>>> = frictionless
formulation<<<{"description": "The contact formulation"}>>> = mortar
[]
[]
[Constraints<<<{"href": "../../../../syntax/Constraints/index.html"}>>>]
# thermal contact constraint
[Tlm]
type = GapConductanceConstraint<<<{"description": "Computes the residual and Jacobian contributions for the 'Lagrange Multiplier' implementation of the thermal contact problem. For more information, see the detailed description here: http://tinyurl.com/gmmhbe9", "href": "../../../../source/constraints/GapConductanceConstraint.html"}>>>
variable<<<{"description": "The name of the lagrange multiplier variable that this constraint is applied to. This parameter may not be supplied in the case of using penalty methods for example"}>>> = Tlm
secondary_variable<<<{"description": "Primal variable on secondary surface."}>>> = T
use_displaced_mesh<<<{"description": "Whether or not this object should use the displaced mesh for computation. Note that in the case this is true but no displacements are provided in the Mesh block the undisplaced mesh will still be used."}>>> = true
k<<<{"description": "Gap conductance"}>>> = 1e-1
primary_boundary<<<{"description": "The name of the primary boundary sideset."}>>> = pillar1_right
primary_subdomain<<<{"description": "The name of the primary subdomain."}>>> = pillars_primary_subdomain
secondary_boundary<<<{"description": "The name of the secondary boundary sideset."}>>> = pillar2_left
secondary_subdomain<<<{"description": "The name of the secondary subdomain."}>>> = pillars_secondary_subdomain
[]
[]
[BCs<<<{"href": "../../../../syntax/BCs/index.html"}>>>]
[bottom_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
boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = 'pillar1_bottom pillar2_bottom'
value<<<{"description": "Value of the BC"}>>> = 0
[]
[bottom_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
boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = 'pillar1_bottom pillar2_bottom'
value<<<{"description": "Value of the BC"}>>> = 0
[]
[Pressure<<<{"href": "../../../../syntax/BCs/Pressure/index.html"}>>>]
[sides]
boundary<<<{"description": "The list of boundaries (ids or names) from the mesh where this object applies"}>>> = 'pillar1_left pillar2_right'
function<<<{"description": "The function that describes the pressure"}>>> = 1e4*t^2
[]
[]
# thermal boundary conditions (pillars are heated/cooled from the bottom)
[heat_left]
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"}>>> = T
boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = pillar1_bottom
value<<<{"description": "Value of the BC"}>>> = 100
[]
[cool_right]
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"}>>> = T
boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = pillar2_bottom
value<<<{"description": "Value of the BC"}>>> = 0
[]
[]
[Materials<<<{"href": "../../../../syntax/Materials/index.html"}>>>]
[elasticity]
type = ComputeIsotropicElasticityTensor<<<{"description": "Compute a constant isotropic elasticity tensor.", "href": "../../../../source/materials/ComputeIsotropicElasticityTensor.html"}>>>
youngs_modulus<<<{"description": "Young's modulus of the material."}>>> = 1e9
poissons_ratio<<<{"description": "Poisson's ratio for the material."}>>> = 0.3
[]
[stress]
type = ComputeFiniteStrainElasticStress<<<{"description": "Compute stress using elasticity for finite strains", "href": "../../../../source/materials/ComputeFiniteStrainElasticStress.html"}>>>
[]
# thermal properties
[thermal_conductivity]
type = HeatConductionMaterial<<<{"description": "General-purpose material model for heat conduction", "href": "../../../../source/materials/HeatConductionMaterial.html"}>>>
thermal_conductivity<<<{"description": "The thermal conductivity value"}>>> = 100
specific_heat<<<{"description": "The specific heat value"}>>> = 1
[]
[density]
type = Density<<<{"description": "Creates density material property. This class is deprecated, and its functionalityis replaced by StrainAdjustedDensity for cases when the density should be adjustedto account for material deformation. If it is not desired to adjust the density fordeformation, a variety of general-purpose Materials, such as GenericConstantMaterialor ParsedMaterial can be used to define the density.", "href": "../../../../source/materials/Density.html"}>>>
density<<<{"description": "Density"}>>> = 1
[]
[]
[Executioner<<<{"href": "../../../../syntax/Executioner/index.html"}>>>]
type = Transient
solve_type = NEWTON
line_search = none
# we deal with the saddle point structure of the system by adding a small shift
petsc_options_iname = '-pc_type -pc_factor_shift_type'
petsc_options_value = 'lu nonzero'
end_time = 5
dt = 0.1
[Predictor<<<{"href": "../../../../syntax/Executioner/Predictor/index.html"}>>>]
type = SimplePredictor
scale = 1
[]
[]
[Outputs<<<{"href": "../../../../syntax/Outputs/index.html"}>>>]
exodus<<<{"description": "Output the results using the default settings for Exodus output."}>>> = true
print_linear_residuals<<<{"description": "Enable printing of linear residuals to the screen (Console)"}>>> = false
perf_graph<<<{"description": "Enable printing of the performance graph to the screen (Console)"}>>> = true
[]Note that we do not have a "thermal contact action", so we need to add all objects required for mortar based thermal contact manually.
Input file
Variables
For the thermal part of the problem we add the temperature T variable and supply an initial_condition of 50 temperature units (this value is right between the 0 and 100 values we're setting below for the cooling and heating BCs). Note that like the other physics variables and objects this variable is restricted to block 0 (see GlobalParams).
The second variable we add is the Lagrange multiplier for the thermal contact. Its physical meaning is the heat flux! This variable will be block restricted to the lower dimensional mortar subdomain the [Contact] action created for us (the name of that subdomain is the name of the contact action subblock with _secondary_subdomain appended).
The "family" and "order" parameters define the shape function family and element order. First order Lagrange (no connection to "Lagrange multiplier!") are the MOOSE default settings. In future inputs you might see those parameters omitted. Lagrange shape functions are purely nodal interpolary shape functions.
Constraints
The GapConductanceConstraint is the bare bones gap conductance constraint object available in MOOSE. It could serve as a starting point for the development of more complex gap conductance models (taking into account surface roughness and contact pressure for example).
The constraint object acts on the Lagrange multiplier "variable" Tlm, and it affects the "secondary_variable" T, the temperature, enforcing the energy balance.
We make sure to enable the evaluation of the constraint on the displaced mesh to capture large displacement and deformation by setting the "use_displaced_mesh" parameter to true. The displaced mesh is a copy of the mesh that has all nodes moved by their displacement variable values. Try setting this parameter to false and compare the results (hint the gap distance the constraint sees will be computed from the original mesh and will not change over time).
We chose a gap conductance "k". The heat flux across the gap will be computed as where is the width of the gap.
Then we set the primary and secondary boundaries and subdomains. Again for the subdomains we choose the lower dimensional meshes that the mechanical contact action set up for us. If you have a purely thermal contact problem you will need to create those lower dimensional subdomains manually using the LowerDBlockFromSidesetGenerator mesh generator.
BCs
heat_left and cool_right are two Dirichlet boundary conditions that set the bottom of the left cantilever to an elevated temperature and the bottom of the right cantilever to a lowered temperature, so we can establish a temperature difference and a resulting heat flux between the two cantilevers.
Materials
We need to add two material blocks for thermal properties. HeatConductionMaterial provides thermal conductivity and specific heat. And the Density material computes an updated material density on the displaced mesh, taking volumetric deformation into account. The specific_heat and density properties are used by the HeatConductionTimeDerivative kernel, and the thermal_conductivity property is used by the HeatConduction kernel.
Executioner
We make a small change to the PETSc options and add the -pc_factor_shift_type
nonzero option in "petsc_options_iname" and "petsc_options_value" respectively. This helps the direct solver deal with the saddle point structure of the Jacobian matrix (with zeroes on the diagonal).
Note that we also lower the timestep to better see what is going on. You can set it back to 0.5 to convince yourself that the problem still converges just as well.
Once you've answered the questions and run this example we will move on to Step 2 in which we introduce a third block participating in the contact problem.