Step 2 - Mortar Contact

Mortar based mechanical contact is a fairly recent addition to MOOSE. While it offers benefits such as tighter convergence and better enforcement of contact constraints as well as a better solution quality for contact interface quantities like the contact pressure, it is not quite as mature as the node on face contact used in the previous step.

Mortar based contact for 3D problems is currently in an experimental stage. Edge to edge contact and edge dropping are not fully implemented and can lead to artifacts at sharp transitions form in-contact to out of contact.

#
# Switching to mortar based mechanical contact
# https://mooseframework.inl.gov/modules/contact/tutorials/introduction/step02.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"}>>> = 5
    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"}>>> = 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"}>>> = 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
[]

[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
    correct_edge_dropping<<<{"description": "Whether to enable correct edge dropping treatment for mortar constraints. When disabled any Lagrange Multiplier degree of freedom on a secondary element without full primary contributions will be set (strongly) to 0."}>>> = true
  []
[]

[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'
      # we square time here to get a more progressive loading curve
      # (more pressure later on once contact is established)
      function<<<{"description": "The function that describes the pressure"}>>> = 1e4*t^2
    []
  []
[]

[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"}>>>
  []
[]

[Executioner<<<{"href": "../../../../syntax/Executioner/index.html"}>>>]
  type = Transient
  solve_type = NEWTON
  line_search = none
  petsc_options_iname = '-pc_type -pc_factor_shift_type'
  petsc_options_value = 'lu       NONZERO'
  end_time = 5
  dt = 0.5
  [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
[]

Here we show the steps for migrating a node on face contact problem (the penalty based contact problem from the previous step) to a mortar formulation. Most of the extra work is performed by the Contact action.

MOOSE's mortar based contact uses Lagrange multipliers for the enforcement of the contact constraints. These Lagrange multipliers are additional solution variables (added by the Contact action) that live on lower dimensional subdomains along the contact interfaces.

The Contact action automatically adds these lower dimensional subdomains, but now that they exist in the simulation we have to be careful not to add any physics on them beyond the Lagrange multipliers. This means we have to add block restrictions for all our kernels, materials, and variables to only add them to the volume subdomains.

Input file

`GlobalParams`

We add block = 0 here to make sure all block restrictable objects explicitly listed in the input (and not generated by actions) are by default restricted to subdomain 0 (the two cantilever volumes) and not unrestricted (meaning they would also live on the lower dimensional meshes).

Try removing those lines to see what happens. Running mechanics on the lower dimensional subdomains will introduce an unphysical artificial stiffness and will change the results significantly.

Besides the added block restriction we only change the Contact block.

`Contact`

We change the "formulation" to mortar and remove the "penalty" and "normalize_penalty" parameters. They have no function in mortar contact.

Tasks and questions

Penetration

Let's look at the penetration variable again:

Visualize the penetration variable. You may have to rescale the visualization to start the color scale at 0. Negative penetrations are not of interest here (they effectively are the gap width)

You should see a maximum interpenetration of about 5.1e-9. That's quite an improvement compared to penalty contact!

Contact pressure

Note that the contact_pressure variable is added but not used with the mortar formulation. The physical meaning of the Lagrange multiplier is the normal contact pressure.

Visualize the contact pressure by plotting the pillars_normal_lm variable.

Extracting contact pressure data

To get a better picture of the contact pressure distribution along a contact surface it can be helpful to create a line plot of the pressure.

Look at the NodalValueSamplervectorpostprocessor and see if you can use its "block" parameter to output the pillars_normal_lm variable on the pillars_secondary_subdomain subdomain.

Click here for the answer.

Once you've answered the questions and run this example we will move on to add thermal transport and later on thermal contact. You will need to build the combined-opt executable to follow along.

Click here to continue the tutorial

(modules/contact/tutorials/introduction/step02.i)

#
# Switching to mortar based mechanical contact
# https://mooseframework.inl.gov/modules/contact/tutorials/introduction/step02.html
#

[GlobalParams]
  displacements = 'disp_x disp_y'
  block = 0
[]

[Mesh]
  [generated1]
    type = GeneratedMeshGenerator
    dim = 2
    nx = 5
    ny = 15
    xmin = -0.6
    xmax = -0.1
    ymax = 5
    bias_y = 0.9
    boundary_name_prefix = pillar1
  []
  [generated2]
    type = GeneratedMeshGenerator
    dim = 2
    nx = 5
    ny = 15
    xmin = 0.1
    xmax = 0.6
    ymax = 5
    bias_y = 0.9
    boundary_name_prefix = pillar2
    boundary_id_offset = 4
  []
  [collect_meshes]
    type = MeshCollectionGenerator
    inputs = 'generated1 generated2'
  []

  patch_update_strategy = iteration
[]

[Physics/SolidMechanics/QuasiStatic]
  [all]
    add_variables = true
    strain = FINITE
    generate_output = 'vonmises_stress'
  []
[]

[Contact]
  [pillars]
    primary = pillar1_right
    secondary = pillar2_left
    model = frictionless
    formulation = mortar
    correct_edge_dropping = true
  []
[]

[BCs]
  [bottom_x]
    type = DirichletBC
    variable = disp_x
    boundary = 'pillar1_bottom pillar2_bottom'
    value = 0
  []
  [bottom_y]
    type = DirichletBC
    variable = disp_y
    boundary = 'pillar1_bottom pillar2_bottom'
    value = 0
  []
  [Pressure]
    [sides]
      boundary = 'pillar1_left pillar2_right'
      # we square time here to get a more progressive loading curve
      # (more pressure later on once contact is established)
      function = 1e4*t^2
    []
  []
[]

[Materials]
  [elasticity]
    type = ComputeIsotropicElasticityTensor
    youngs_modulus = 1e9
    poissons_ratio = 0.3
  []
  [stress]
    type = ComputeFiniteStrainElasticStress
  []
[]

[Executioner]
  type = Transient
  solve_type = NEWTON
  line_search = none
  petsc_options_iname = '-pc_type -pc_factor_shift_type'
  petsc_options_value = 'lu       NONZERO'
  end_time = 5
  dt = 0.5
  [Predictor]
    type = SimplePredictor
    scale = 1
  []
[]

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

Input file
Tasks and questions