MassMatrix

This object is meant to build mass matrices for preconditioning techniques that require them. It is only meant for filling matrices, so all vector tag parameters are suppressed. The "matrix_tags" parameter default is cleared such that the user should provide some non-empty parameter value for "matrix_tags" or "extra_matrix_tags".

Input Parameters

variableThe name of the variable that this residual object operates on
C++ Type:NonlinearVariableName
Unit:(no unit assumed)
Controllable:No
Description:The name of the variable that this residual object operates on

Required Parameters

blockThe list of blocks (ids or names) that this object will be applied
C++ Type:std::vector<SubdomainName>
Controllable:No
Description:The list of blocks (ids or names) that this object will be applied
density1The material property defining the density
Default:1
C++ Type:MaterialPropertyName
Unit:(no unit assumed)
Controllable:No
Description:The material property defining the density
displacementsThe displacements
C++ Type:std::vector<VariableName>
Unit:(no unit assumed)
Controllable:No
Description:The displacements
matrix_onlyTrueWhether this object is only doing assembly to matrices (no vectors)
Default:True
C++ Type:bool
Controllable:No
Description:Whether this object is only doing assembly to matrices (no vectors)

Optional Parameters

control_tagsAdds user-defined labels for accessing object parameters via control logic.
C++ Type:std::vector<std::string>
Controllable:No
Description:Adds user-defined labels for accessing object parameters via control logic.
diag_save_inThe name of auxiliary variables to save this Kernel's diagonal Jacobian contributions to. Everything about that variable must match everything about this variable (the type, what blocks it's on, etc.)
C++ Type:std::vector<AuxVariableName>
Unit:(no unit assumed)
Controllable:No
Description:The name of auxiliary variables to save this Kernel's diagonal Jacobian contributions to. Everything about that variable must match everything about this variable (the type, what blocks it's on, etc.)
enableTrueSet the enabled status of the MooseObject.
Default:True
C++ Type:bool
Controllable:Yes
Description:Set the enabled status of the MooseObject.
implicitTrueDetermines whether this object is calculated using an implicit or explicit form
Default:True
C++ Type:bool
Controllable:No
Description:Determines whether this object is calculated using an implicit or explicit form
save_inThe name of auxiliary variables to save this Kernel's residual contributions to. Everything about that variable must match everything about this variable (the type, what blocks it's on, etc.)
C++ Type:std::vector<AuxVariableName>
Unit:(no unit assumed)
Controllable:No
Description:The name of auxiliary variables to save this Kernel's residual contributions to. Everything about that variable must match everything about this variable (the type, what blocks it's on, etc.)
search_methodnearest_node_connected_sidesChoice of search algorithm. All options begin by finding the nearest node in the primary boundary to a query point in the secondary boundary. In the default nearest_node_connected_sides algorithm, primary boundary elements are searched iff that nearest node is one of their nodes. This is fast to determine via a pregenerated node-to-elem map and is robust on conforming meshes. In the optional all_proximate_sides algorithm, primary boundary elements are searched iff they touch that nearest node, even if they are not topologically connected to it. This is more CPU-intensive but is necessary for robustness on any boundary surfaces which has disconnections (such as Flex IGA meshes) or non-conformity (such as hanging nodes in adaptively h-refined meshes).
Default:nearest_node_connected_sides
C++ Type:MooseEnum
Options:nearest_node_connected_sides, all_proximate_sides
Controllable:No
Description:Choice of search algorithm. All options begin by finding the nearest node in the primary boundary to a query point in the secondary boundary. In the default nearest_node_connected_sides algorithm, primary boundary elements are searched iff that nearest node is one of their nodes. This is fast to determine via a pregenerated node-to-elem map and is robust on conforming meshes. In the optional all_proximate_sides algorithm, primary boundary elements are searched iff they touch that nearest node, even if they are not topologically connected to it. This is more CPU-intensive but is necessary for robustness on any boundary surfaces which has disconnections (such as Flex IGA meshes) or non-conformity (such as hanging nodes in adaptively h-refined meshes).
seed0The seed for the master random number generator
Default:0
C++ Type:unsigned int
Controllable:No
Description:The seed for the master random number generator
use_displaced_meshFalseWhether 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.
Default:False
C++ Type:bool
Controllable:No
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.

Advanced Parameters

extra_matrix_tagsThe extra tags for the matrices this Kernel should fill
C++ Type:std::vector<TagName>
Controllable:No
Description:The extra tags for the matrices this Kernel should fill
matrix_tagsThe tag for the matrices this Kernel should fill
C++ Type:MultiMooseEnum
Options:nontime, system
Controllable:No
Description:The tag for the matrices this Kernel should fill

Contribution To Tagged Field Data Parameters

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

Material Property Retrieval Parameters

Input Files

(modules/navier_stokes/test/tests/finite_element/ins/hdg/ip/lid-driven/lid-driven-fsp.i)
(modules/navier_stokes/test/tests/finite_element/ins/lid_driven/transient_fsp.i)
(modules/solid_mechanics/test/tests/dynamics/time_integration/direct_central_difference_varied_dt.i)
(modules/navier_stokes/test/tests/finite_element/ins/lid_driven/steady_fsp.i)
(modules/navier_stokes/test/tests/finite_element/ins/lid_driven/steady_vector_fsp_stokes.i)
(modules/solid_mechanics/test/tests/dynamics/explicit_mms/mms_direct_combined.i)
(modules/contact/test/tests/explicit_dynamics/deep_impact.i)
(modules/contact/test/tests/explicit_dynamics/settlement.i)
(test/tests/tag/fe-mass-matrix.i)
(modules/contact/test/tests/explicit_dynamics/test_balance_optimized.i)
(modules/navier_stokes/test/tests/finite_element/ins/hdg/ip/lid-driven/matrix-testing.i)
(modules/solid_mechanics/test/tests/dynamics/time_integration/direct_central_difference_multiVarBC.i)
(test/tests/geomsearch/3d_penetration_locator/3d_disconnected_tet.i)
(modules/navier_stokes/test/tests/finite_element/ins/lid_driven/steady_vector_fsp.i)
(modules/solid_mechanics/test/tests/dynamics/time_integration/direct_central_difference.i)
(modules/solid_mechanics/test/tests/dynamics/explicit_mms/mms_direct_second_order_with_ic.i)
(modules/navier_stokes/test/tests/finite_element/ins/lid_driven/steady_vector_fsp_al.i)
(modules/contact/test/tests/explicit_dynamics/highvel.i)
(test/tests/geomsearch/penetration_locator/disconnected_penetration.i)
(modules/contact/test/tests/explicit_dynamics/test_balance.i)

(modules/navier_stokes/test/tests/finite_element/ins/hdg/ip/lid-driven/lid-driven-fsp.i)

mu = 1
rho = 1
l = 1
U = 1e3
n = 200
gamma = 1e5

[Mesh]
  [gen]
    type = GeneratedMeshGenerator
    dim = 2
    xmin = 0
    xmax = ${l}
    ymin = 0
    ymax = ${l}
    nx = ${n}
    ny = ${n}
    elem_type = TRI6
  []
[]

[Problem]
  type = NavierStokesProblem
  extra_tag_matrices = 'mass'
  mass_matrix = 'mass'
  use_pressure_mass_matrix = true
[]

[Variables]
  [vel_x]
    family = MONOMIAL
    order = FIRST
  []
  [vel_y]
    family = MONOMIAL
    order = FIRST
  []
  [pressure]
    family = MONOMIAL
    order = CONSTANT
  []
  [vel_bar_x]
    family = SIDE_HIERARCHIC
    order = FIRST
  []
  [vel_bar_y]
    family = SIDE_HIERARCHIC
    order = FIRST
  []
  [pressure_bar]
    family = SIDE_HIERARCHIC
    order = FIRST
  []
[]

[HDGKernels]
  [momentum_x_convection]
    type = AdvectionIPHDGKernel
    variable = vel_x
    face_variable = vel_bar_x
    velocity = 'velocity'
    coeff = ${rho}
  []
  [momentum_x_diffusion]
    type = NavierStokesStressIPHDGKernel
    variable = vel_x
    face_variable = vel_bar_x
    diffusivity = 'mu'
    alpha = 6
    pressure_variable = pressure
    pressure_face_variable = pressure_bar
    component = 0
  []
  [momentum_y_convection]
    type = AdvectionIPHDGKernel
    variable = vel_y
    face_variable = vel_bar_y
    velocity = 'velocity'
    coeff = ${rho}
  []
  [momentum_y_diffusion]
    type = NavierStokesStressIPHDGKernel
    variable = vel_y
    face_variable = vel_bar_y
    diffusivity = 'mu'
    alpha = 6
    pressure_variable = pressure
    pressure_face_variable = pressure_bar
    component = 1
  []
  [pressure_convection]
    type = AdvectionIPHDGKernel
    variable = pressure
    face_variable = pressure_bar
    velocity = 'velocity'
    coeff = '${fparse -rho}'
    self_advection = false
  []
[]

[Kernels]
  [mass_matrix_pressure]
    type = MassMatrix
    variable = pressure
    matrix_tags = 'mass'
    density = '${fparse -1/gamma}'
  []

  [grad_div_x]
    type = GradDiv
    variable = vel_x
    u = vel_x
    v = vel_y
    gamma = ${gamma}
    component = 0
  []
  [grad_div_y]
    type = GradDiv
    variable = vel_y
    u = vel_x
    v = vel_y
    gamma = ${gamma}
    component = 1
  []
[]

[DGKernels]
  [pb_mass]
    type = MassMatrixDGKernel
    variable = pressure_bar
    matrix_tags = 'mass'
    density = '${fparse -1/gamma}'
  []

  [u_jump]
    type = MassFluxPenalty
    variable = vel_x
    u = vel_x
    v = vel_y
    component = 0
    gamma = ${gamma}
  []
  [v_jump]
    type = MassFluxPenalty
    variable = vel_y
    u = vel_x
    v = vel_y
    component = 1
    gamma = ${gamma}
  []
[]

[BCs]
  [momentum_x_diffusion_walls]
    type = NavierStokesStressIPHDGDirichletBC
    boundary = 'left bottom right'
    variable = vel_x
    face_variable = vel_bar_x
    pressure_variable = pressure
    pressure_face_variable = pressure_bar
    alpha = 6
    functor = '0'
    diffusivity = 'mu'
    component = 0
  []
  [momentum_x_diffusion_top]
    type = NavierStokesStressIPHDGDirichletBC
    boundary = 'top'
    variable = vel_x
    face_variable = vel_bar_x
    pressure_variable = pressure
    pressure_face_variable = pressure_bar
    alpha = 6
    functor = '${U}'
    diffusivity = 'mu'
    component = 0
  []
  [momentum_y_diffusion_all]
    type = NavierStokesStressIPHDGDirichletBC
    boundary = 'left bottom right top'
    variable = vel_y
    face_variable = vel_bar_y
    pressure_variable = pressure
    pressure_face_variable = pressure_bar
    alpha = 6
    functor = '0'
    diffusivity = 'mu'
    component = 1
  []

  [mass_convection]
    type = AdvectionIPHDGPrescribedFluxBC
    face_variable = pressure_bar
    variable = pressure
    velocity = 'velocity'
    coeff = '${fparse -rho}'
    self_advection = false
    boundary = 'left bottom top right'
    prescribed_normal_flux = 0
  []

  [pb_mass]
    type = MassMatrixIntegratedBC
    variable = pressure_bar
    matrix_tags = 'mass'
    boundary = 'left right bottom top'
    density = '${fparse -1/gamma}'
  []

  [u_jump_walls]
    type = MassFluxPenaltyBC
    variable = vel_x
    u = vel_x
    v = vel_y
    component = 0
    boundary = 'left right bottom'
    gamma = ${gamma}
    dirichlet_value = 'walls'
  []
  [v_jump_walls]
    type = MassFluxPenaltyBC
    variable = vel_y
    u = vel_x
    v = vel_y
    component = 1
    boundary = 'left right bottom'
    gamma = ${gamma}
    dirichlet_value = 'walls'
  []
  [u_jump_top]
    type = MassFluxPenaltyBC
    variable = vel_x
    u = vel_x
    v = vel_y
    component = 0
    boundary = 'top'
    gamma = ${gamma}
    dirichlet_value = 'top'
  []
  [v_jump_top]
    type = MassFluxPenaltyBC
    variable = vel_y
    u = vel_x
    v = vel_y
    component = 1
    boundary = 'top'
    gamma = ${gamma}
    dirichlet_value = 'top'
  []
[]

[Functions]
  [top]
    type = ParsedVectorFunction
    value_x = ${U}
    value_y = 0
  []
  [walls]
    type = ParsedVectorFunction
    value_x = 0
    value_y = 0
  []
[]

[Materials]
  [const]
    type = ADGenericConstantMaterial
    prop_names = 'rho mu'
    prop_values = '${rho} ${mu}'
  []
  [vel]
    type = ADVectorFromComponentVariablesMaterial
    vector_prop_name = 'velocity'
    u = vel_x
    v = vel_y
  []
[]

[Preconditioning]
  [FSP]
    type = FSP
    topsplit = 'up'
    [up]
      splitting = 'u p'
      splitting_type = schur
      petsc_options_iname = '-pc_fieldsplit_schur_fact_type  -pc_fieldsplit_schur_precondition -ksp_gmres_restart -ksp_type -ksp_pc_side -ksp_rtol'
      petsc_options_value = 'full                            self                              300                fgmres    right        1e-4'
    []
    [u]
      vars = 'vel_x vel_y vel_bar_x vel_bar_y'
      petsc_options = '-ksp_converged_reason'
      petsc_options_iname = '-pc_type -ksp_type -ksp_rtol -ksp_gmres_restart -ksp_pc_side -pc_factor_mat_solver_type'
      petsc_options_value = 'ilu      gmres     1e-2      300                right        strumpack'
    []
    [p]
      vars = 'pressure pressure_bar'
      petsc_options = '-ksp_converged_reason'
      petsc_options_iname = '-ksp_type -ksp_gmres_restart -ksp_rtol -pc_type -ksp_pc_side'
      petsc_options_value = 'gmres     300                1e-2      ilu      right'
    []
  []
[]

[Executioner]
  type = Steady
  solve_type = 'NEWTON'
[]

[Outputs]
  active = ''
  [out]
    type = Exodus
    hide = 'pressure_average vel_bar_x vel_bar_y pressure_bar'
  []
[]

[Postprocessors]
  [Re]
    type = ParsedPostprocessor
    pp_names = ''
    expression = '${rho} * ${U} * ${l} / ${mu}'
  []
  [pressure_average]
    type = ElementAverageValue
    variable = pressure
  []
[]

[Correctors]
  [set_pressure]
    type = NSPressurePin
    pin_type = 'average'
    variable = pressure
    pressure_average = 'pressure_average'
  []
[]

(modules/navier_stokes/test/tests/finite_element/ins/lid_driven/transient_fsp.i)

n=64
mu=2e-3

[GlobalParams]
  gravity = '0 0 0'
  preset = true
  supg = false
[]

[Problem]
  extra_tag_matrices = 'mass'
  previous_nl_solution_required = true
  type = NavierStokesProblem
  mass_matrix = 'mass'
  schur_fs_index = '1'
[]

[Mesh]
  [gen]
    type = GeneratedMeshGenerator
    dim = 2
    xmin = 0
    xmax = 1.0
    ymin = 0
    ymax = 1.0
    nx = ${n}
    ny = ${n}
    elem_type = QUAD9
  []
[]

[Variables]
  [vel_x]
    order = SECOND
    family = LAGRANGE
  []
  [vel_y]
    order = SECOND
    family = LAGRANGE
  []
  [p]
    order = FIRST
    family = LAGRANGE
  []
[]

[Kernels]
  # mass
  [mass]
    type = INSMass
    variable = p
    u = vel_x
    v = vel_y
    pressure = p
  []

  [x_time]
    type = INSMomentumTimeDerivative
    variable = vel_x
  []
  [x_momentum_space]
    type = INSMomentumLaplaceForm
    variable = vel_x
    u = vel_x
    v = vel_y
    pressure = p
    component = 0
  []
  [x_mass]
    type = MassMatrix
    variable = vel_x
    matrix_tags = 'mass'
  []

  [y_time]
    type = INSMomentumTimeDerivative
    variable = vel_y
  []
  [y_momentum_space]
    type = INSMomentumLaplaceForm
    variable = vel_y
    u = vel_x
    v = vel_y
    pressure = p
    component = 1
  []
  [y_mass]
    type = MassMatrix
    variable = vel_y
    matrix_tags = 'mass'
  []
[]

[BCs]
  [x_no_slip]
    type = DirichletBC
    variable = vel_x
    boundary = 'bottom right left'
    value = 0.0
  []
  [lid]
    type = FunctionDirichletBC
    variable = vel_x
    boundary = 'top'
    function = 'lid_function'
  []
  [y_no_slip]
    type = DirichletBC
    variable = vel_y
    boundary = 'bottom right top left'
    value = 0.0
  []
[]

[Materials]
  [const]
    type = GenericConstantMaterial
    block = 0
    prop_names = 'rho mu'
    prop_values = '1  ${mu}'
  []
[]

[Functions]
  [lid_function]
    # We pick a function that is exactly represented in the velocity
    # space so that the Dirichlet conditions are the same regardless
    # of the mesh spacing.
    type = ParsedFunction
    expression = '4*x*(1-x)'
  []
[]

[Preconditioning]
  [FSP]
    type = FSP
    topsplit = 'by_diri_others'
    [by_diri_others]
      splitting = 'diri others'
      splitting_type  = additive
      petsc_options_iname = '-ksp_type'
      petsc_options_value = 'preonly'
    []
      [diri]
        sides = 'left right top bottom'
        vars = 'vel_x vel_y'
        petsc_options_iname = '-pc_type'
        petsc_options_value = 'jacobi'
      []
      [others]
        splitting = 'u p'
        splitting_type  = schur
        petsc_options_iname = '-pc_fieldsplit_schur_fact_type  -pc_fieldsplit_schur_precondition -ksp_gmres_restart -ksp_rtol -ksp_type -ksp_atol'
        petsc_options_value = 'full                            self                              300                1e-5      fgmres    1e-9'
        unside_by_var_boundary_name = 'left top right bottom left top right bottom'
        unside_by_var_var_name = 'vel_x vel_x vel_x vel_x vel_y vel_y vel_y vel_y'
      []
        [u]
          vars = 'vel_x vel_y'
          unside_by_var_boundary_name = 'left top right bottom left top right bottom'
          unside_by_var_var_name = 'vel_x vel_x vel_x vel_x vel_y vel_y vel_y vel_y'
          # petsc_options = '-ksp_converged_reason'
          petsc_options_iname = '-pc_type -ksp_pc_side -ksp_type -ksp_rtol -pc_hypre_type -ksp_gmres_restart'
          petsc_options_value = 'hypre    right        gmres     1e-2      boomeramg      300'
        []
        [p]
          vars = 'p'
          petsc_options = '-pc_lsc_scale_diag -ksp_converged_reason'# -lsc_ksp_converged_reason -lsc_ksp_monitor_true_residual
          petsc_options_iname = '-ksp_type -ksp_gmres_restart -ksp_rtol -pc_type -ksp_pc_side -lsc_pc_type -lsc_pc_hypre_type -lsc_ksp_type -lsc_ksp_rtol -lsc_ksp_pc_side -lsc_ksp_gmres_restart'
          petsc_options_value = 'fgmres    300                1e-2      lsc      right        hypre        boomeramg          gmres         1e-1          right            300'
        []
  []
[]

[Postprocessors]
  [pavg]
    type = ElementAverageValue
    variable = p
  []
[]

[Correctors]
  [set_pressure]
    type = NSPressurePin
    pin_type = 'average'
    variable = p
    pressure_average = 'pavg'
  []
[]

[Executioner]
  solve_type = NEWTON
  type = Transient
  petsc_options_iname = '-snes_max_it'
  petsc_options_value = '100'
  line_search = 'none'
  nl_rel_tol = 1e-8
  nl_abs_tol = 1e-9
  abort_on_solve_fail = true
  normalize_solution_diff_norm_by_dt = false
  [TimeStepper]
    type = IterationAdaptiveDT
    optimal_iterations = 6
    dt = 1e-2
  []
  steady_state_detection = true
[]

[Outputs]
  [exo]
    type = Exodus
    execute_on = 'final'
    hide = 'pavg'
  []
[]

(modules/solid_mechanics/test/tests/dynamics/time_integration/direct_central_difference_varied_dt.i)

###########################################################
# This is a simple test with a time-dependent problem
# demonstrating the use of a central difference with a
# direct calculation of acceleration.
#
# Testing that the first and second time derivatives
# are calculated correctly using the Direct Central Difference
# method

# Testing the accuracy of the timestep averaging method within
# the Direct Central Difference method
###########################################################

[GlobalParams]
  displacements = 'disp_x disp_y'
[]

[Problem]
  extra_tag_matrices = 'mass'
[]

[Mesh]
  type = GeneratedMesh
  dim = 2
  xmin = -1
  xmax = 1
  ymin = -1
  ymax = 1
  nx = 1
  ny = 1
[]

[Variables]
  [disp_x]
  []
  [disp_y]
  []
[]

[Functions]
  #Mid-step velocities
  #0 0.00625 0.015 0.0075 0.25 0
  #Accelerations
  #0.025 0.01944 -0.01 0.48 -2.17
  [forcing_fn]
    type = PiecewiseLinear
    x = '0.0 0.1 0.5    1.0  2.0   2.01 2.23'
    y = '0.0 0.0 0.0025 0.01 0.0175 0.02 0.02'
  []
[]

[Kernels]
  [DynamicSolidMechanics]
    displacements = 'disp_x disp_y'
  []
  [massmatrix]
    type = MassMatrix
    density = density
    matrix_tags = 'mass'
    variable = disp_x
  []
  [massmatrix_y]
    type = MassMatrix
    density = density
    matrix_tags = 'mass'
    variable = disp_y
  []
[]

[Materials]
  [elasticity_tensor_block_one]
    type = ComputeIsotropicElasticityTensor
    youngs_modulus = 1e1
    poissons_ratio = 0.0
  []

  [strain_block]
    type = ComputeFiniteStrain
    displacements = 'disp_x disp_y'
    implicit = false
  []
  [stress_block]
    type = ComputeFiniteStrainElasticStress
    implicit = false
  []
  [density]
    type = GenericConstantMaterial
    prop_names = 'density'
    prop_values = 5
  []
[]

[BCs]
  [left_x]
    type = ExplicitFunctionDirichletBC
    variable = disp_x
    boundary = 'left'
    function = forcing_fn
  []
  [right_x]
    type = ExplicitFunctionDirichletBC
    variable = disp_x
    boundary = 'right'
    function = forcing_fn
  []
[]

[Executioner]
  type = Transient

  [TimeIntegrator]
    type = ExplicitMixedOrder
    mass_matrix_tag = 'mass'
    second_order_vars = 'disp_x disp_y'
  []
  [TimeStepper]
    type = TimeSequenceStepper
    time_sequence = '0.0 0.1 0.5 1.0 2.0 2.01 2.23'
  []
  start_time = 0.0
  num_steps = 6
  dt = 0.1
[]

[Postprocessors]
  [udot]
    type = ElementAverageTimeDerivative
    variable = disp_x
  []
  [udotdot]
    type = ElementAverageSecondTimeDerivative
    variable = disp_x
  []
  [u]
    type = ElementAverageValue
    variable = disp_x
  []
[]

[Outputs]
  exodus = true
[]

(modules/navier_stokes/test/tests/finite_element/ins/lid_driven/steady_fsp.i)

rho=1
mu=2e-3
U=1
l=1
prefactor=${fparse 1/(l/2)^2}
n=64

[GlobalParams]
  gravity = '0 0 0'
[]

[Mesh]
  [gen]
    type = DistributedRectilinearMeshGenerator
    dim = 2
    xmin = 0
    xmax = ${l}
    ymin = 0
    ymax = ${l}
    nx = ${n}
    ny = ${n}
    elem_type = QUAD4
  []
  second_order = true
  parallel_type = distributed
[]

[Variables]
  [vel_x]
    order = SECOND
    family = LAGRANGE
  []
  [vel_y]
    order = SECOND
    family = LAGRANGE
  []
  [p]
    order = FIRST
    family = LAGRANGE
  []
[]

[Kernels]
  [mass]
    type = INSMass
    variable = p
    u = vel_x
    v = vel_y
    pressure = p
  []

  [x_momentum_space]
    type = INSMomentumLaplaceForm
    variable = vel_x
    u = vel_x
    v = vel_y
    pressure = p
    component = 0
  []
  [momentum_x_mass]
    type = MassMatrix
    variable = vel_x
    density = ${rho}
    matrix_tags = 'mass'
  []
  [y_momentum_space]
    type = INSMomentumLaplaceForm
    variable = vel_y
    u = vel_x
    v = vel_y
    pressure = p
    component = 1
  []
  [momentum_y_mass]
    type = MassMatrix
    variable = vel_y
    density = ${rho}
    matrix_tags = 'mass'
  []
[]

[BCs]
  [x_no_slip]
    type = DirichletBC
    variable = vel_x
    boundary = 'bottom right left'
    value = 0.0
  []
  [lid]
    type = FunctionDirichletBC
    variable = vel_x
    boundary = 'top'
    function = 'lid_function'
  []
  [y_no_slip]
    type = DirichletBC
    variable = vel_y
    boundary = 'bottom right top left'
    value = 0.0
  []
[]

[Materials]
  [const]
    type = GenericConstantMaterial
    block = 0
    prop_names = 'rho mu'
    prop_values = '${rho} ${mu}'
  []
[]

[Functions]
  [lid_function]
    # We pick a function that is exactly represented in the velocity
    # space so that the Dirichlet conditions are the same regardless
    # of the mesh spacing.
    type = ParsedFunction
    expression = '${prefactor}*${U}*x*(${l}-x)'
  []
[]

[Problem]
  type = NavierStokesProblem
  mass_matrix = 'mass'
  extra_tag_matrices = 'mass'
[]

[Preconditioning]
  [FSP]
    type = FSP
    topsplit = 'up'
    [up]
      splitting = 'u p'
      splitting_type  = schur
      petsc_options_iname = '-pc_fieldsplit_schur_fact_type  -pc_fieldsplit_schur_precondition -ksp_gmres_restart -ksp_type -ksp_pc_side -ksp_rtol'
      petsc_options_value = 'full                            self                              300                fgmres    right        1e-4'
    []
      [u]
        vars = 'vel_x vel_y'
        # petsc_options = '-ksp_converged_reason'
        petsc_options_iname = '-pc_type -pc_hypre_type -ksp_type -ksp_rtol -ksp_gmres_restart -ksp_pc_side'
        petsc_options_value = 'hypre    boomeramg      gmres     1e-2      300                right'
      []
      [p]
        vars = 'p'
        petsc_options = '-pc_lsc_scale_diag -ksp_converged_reason'# -lsc_ksp_converged_reason -lsc_ksp_monitor_true_residual
        petsc_options_iname = '-ksp_type -ksp_gmres_restart -ksp_rtol -pc_type -ksp_pc_side -lsc_pc_type -lsc_pc_hypre_type -lsc_ksp_type -lsc_ksp_rtol -lsc_ksp_pc_side -lsc_ksp_gmres_restart'
        petsc_options_value = 'fgmres    300                1e-2      lsc      right        hypre        boomeramg          gmres         1e-1          right            300'
      []
  []
[]

[Postprocessors]
  [pavg]
    type = ElementAverageValue
    variable = p
  []
[]

[Correctors]
  [set_pressure]
    type = NSPressurePin
    pin_type = 'average'
    variable = p
    pressure_average = 'pavg'
  []
[]

[Executioner]
  type = Steady
  solve_type = NEWTON
[]

[Outputs]
  [exo]
    type = Exodus
    execute_on = 'final'
    hide = 'pavg'
  []
[]

(modules/navier_stokes/test/tests/finite_element/ins/lid_driven/steady_vector_fsp_stokes.i)

rho=1
mu=1
U=1
l=1
prefactor=${fparse 1/(l/2)^2}
n=8

[Mesh]
  [gen]
    type = DistributedRectilinearMeshGenerator
    dim = 2
    xmin = 0
    xmax = ${l}
    ymin = 0
    ymax = ${l}
    nx = ${n}
    ny = ${n}
    elem_type = QUAD4
  []
  second_order = true
  parallel_type = distributed
[]

[Variables]
  [vel]
    order = SECOND
    family = LAGRANGE_VEC
  []
  [p]
    order = FIRST
    family = LAGRANGE
  []
[]

[Kernels]
  [mass]
    type = INSADMass
    variable = p
  []
  [mass_kernel]
    type = MassMatrix
    variable = p
    matrix_tags = 'mass'
  []
  [momentum_viscous]
    type = INSADMomentumViscous
    variable = vel
  []
  [momentum_pressure]
    type = INSADMomentumPressure
    variable = vel
    pressure = p
    integrate_p_by_parts = true
  []
[]

[BCs]
  [no_slip]
    type = ADVectorFunctionDirichletBC
    variable = vel
    boundary = 'bottom right left'
  []
  [lid]
    type = ADVectorFunctionDirichletBC
    variable = vel
    boundary = 'top'
    function_x = 'lid_function'
  []
[]

[Materials]
  [const]
    type = ADGenericConstantMaterial
    prop_names = 'rho mu'
    prop_values = '${rho} ${mu}'
  []
  [insad]
    type = INSADMaterial
    velocity = vel
    pressure = p
  []
[]

[Functions]
  [lid_function]
    # We pick a function that is exactly represented in the velocity
    # space so that the Dirichlet conditions are the same regardless
    # of the mesh spacing.
    type = ParsedFunction
    expression = '${prefactor}*${U}*x*(${l}-x)'
  []
[]

[Problem]
  type = NavierStokesProblem
  mass_matrix = 'mass'
  extra_tag_matrices = 'mass'
  use_pressure_mass_matrix = true
[]

[Preconditioning]
  [FSP]
    type = FSP
    topsplit = 'up'
    [up]
      splitting = 'u p'
      splitting_type  = schur
      petsc_options_iname = '-pc_fieldsplit_schur_fact_type  -pc_fieldsplit_schur_precondition -ksp_gmres_restart -ksp_type -ksp_pc_side -ksp_rtol'
      petsc_options_value = 'full                            self                              300                fgmres    right        1e-4'
    []
      [u]
        vars = 'vel'
        # petsc_options = '-ksp_converged_reason'
        petsc_options_iname = '-pc_type -pc_hypre_type -ksp_type -ksp_rtol -ksp_gmres_restart -ksp_pc_side'
        petsc_options_value = 'hypre    boomeramg      gmres     1e-2      300                right'
      []
      [p]
        vars = 'p'
        petsc_options = '-ksp_converged_reason'
        petsc_options_iname = '-ksp_type -ksp_gmres_restart -ksp_rtol -pc_type -ksp_pc_side -pc_hypre_type'
        petsc_options_value = 'fgmres    300                1e-2      hypre    right        boomeramg'
      []
  []
[]

[Postprocessors]
  [pavg]
    type = ElementAverageValue
    variable = p
  []
[]

[Correctors]
  [set_pressure]
    type = NSPressurePin
    pin_type = 'average'
    variable = p
    pressure_average = 'pavg'
  []
[]

[Executioner]
  type = Steady
  solve_type = NEWTON
  nl_rel_tol = 1e-12
[]

[Outputs]
  print_linear_residuals = false
  [exo]
    type = Exodus
    execute_on = 'final'
    hide = 'pavg'
  []
[]

(modules/solid_mechanics/test/tests/dynamics/explicit_mms/mms_direct_combined.i)

[Mesh]
  type = GeneratedMesh
  dim = 2
  xmin = -1
  xmax = 1
  ymin = -1
  ymax = 1
  nx = 10
  ny = 10
  elem_type = QUAD4
[]

[Problem]
  extra_tag_matrices = 'mass'
[]

[Variables]
  [u_first]
  []
  [u_second]
  []
[]

[AuxVariables]
  [dot_u_second]

  []
[]

[AuxKernels]
  [dot_u_second]
    type = TestNewmarkTI
    variable = dot_u_second
    displacement = u_second
    first = true
    execute_on = 'TIMESTEP_END'
  []
[]

[Functions]
  [ic]
    type = ParsedFunction
    expression = 0
  []

  [forcing_fn_first]
    type = ParsedFunction
    expression = (x+y)
  []

  [exact_fn_first]
    type = ParsedFunction
    expression = t*(x+y)
  []

  [exact_fn_second]
    type = ParsedFunction
    expression = '0.5*t^2'
  []
  [exact_dot_fn_second]
    type = ParsedFunction
    expression = 't'
  []
[]

[Kernels]
  [Mass_x_second]
    type = MassMatrix
    variable = u_second
    density = 1
    matrix_tags = 'mass'
  []
  [ffn_second]
    type = BodyForce
    variable = u_second
    function = 1
  []
  [Mass_x_first]
    type = MassMatrix
    variable = u_first
    density = 1
    matrix_tags = 'mass'
  []
  [diff]
    type = Diffusion
    variable = u_first
  []
  [ffn]
    type = BodyForce
    variable = u_first
    function = forcing_fn_first
  []
[]

[BCs]
  [all]
    type = ExplicitFunctionDirichletBC
    variable = u_first
    boundary = '0 1 2 3'
    function = exact_fn_first
  []
[]

[Postprocessors]
  [l2_err_second]
    type = ElementL2Error
    variable = u_second
    function = exact_fn_second
  []
  [l2_dot_err_second]
    type = ElementL2Error
    variable = dot_u_second
    function = exact_dot_fn_second
  []
  [l2_err_first]
    type = ElementL2Error
    variable = u_first
    function = exact_fn_first
  []
[]

[Executioner]
  type = Transient
  start_time = 0.0
  num_steps = 20
  dt = 0.00005
  l_tol = 1e-12
  [TimeIntegrator]
    type = ExplicitMixedOrder
    mass_matrix_tag = 'mass'
    use_constant_mass = true
    second_order_vars = 'u_second'
    first_order_vars = 'u_first'
  []
[]

[Outputs]
  exodus = true
  [console]
    type = Console
    max_rows = 10
  []
[]

# Tests stateful material properties below

[Materials]
  [u_sqrd]
    type = ParsedMaterial
    expression = 'u_first^2'
    property_name = u_sqrd
    coupled_variables = 'u_first'
  []
[]

[AuxVariables]
  [diff_t_begin]
    family = MONOMIAL
    order = CONSTANT
  []
  [diff_t_end]
    family = MONOMIAL
    order = CONSTANT
  []
[]

[AuxKernels]
  [diff_t_begin]
    type = MaterialRateRealAux
    property = u_sqrd
    variable = diff_t_begin
    execute_on = LINEAR
  []
  [diff_t_end]
    type = MaterialRateRealAux
    property = u_sqrd
    variable = diff_t_end
    execute_on = TIMESTEP_END
  []
[]

[Postprocessors]
  [l2_norm_begin]
    type = ElementL2Norm
    variable = diff_t_begin
  []
  [l2_norm_end]
    type = ElementL2Norm
    variable = diff_t_end
  []
[]

(modules/contact/test/tests/explicit_dynamics/deep_impact.i)

# This test demonstrates explicit contact with MOOSE and includes optimizations
# to enhance performance.
[GlobalParams]
  displacements = 'disp_x disp_y disp_z'
[]

[Problem]
  extra_tag_matrices = 'mass'
[]

[Mesh]
  [block_one]
    type = GeneratedMeshGenerator
    dim = 3
    nx = 4
    ny = 4
    nz = 4
    xmin = 4.5
    xmax = 5.5
    ymin = 4.5
    ymax = 5.5
    zmin = 0.0001
    zmax = 1.0001
    boundary_name_prefix = 'ball'
  []
  [block_two]
    type = GeneratedMeshGenerator
    dim = 3
    nx = 9
    ny = 9
    nz = 4
    xmin = 3
    xmax = 7
    ymin = 3
    ymax = 7
    zmin = -2
    zmax = 0
    boundary_name_prefix = 'base'
    boundary_id_offset = 10
  []
  [block_one_id]
    type = SubdomainIDGenerator
    input = block_one
    subdomain_id = 1
  []
  [block_two_id]
    type = SubdomainIDGenerator
    input = block_two
    subdomain_id = 2
  []
  [combine]
    type = MeshCollectionGenerator
    inputs = ' block_one_id block_two_id'
  []
  allow_renumbering = false
  # patch_update_strategy = always
[]

[Variables]
  [disp_x]
  []
  [disp_y]
  []
  [disp_z]
  []
[]

[AuxVariables]
  [gap_rate]
  []
  [vel_x]
  []
  [accel_x]
  []
  [vel_y]
  []
  [accel_y]
  []
  [vel_z]
  []
  [accel_z]
  []
[]

[AuxKernels]
  [accel_x]
    type = TestNewmarkTI
    variable = accel_x
    displacement = disp_x
    first = false
    execute_on = 'TIMESTEP_END'
  []
  [vel_x]
    type = TestNewmarkTI
    variable = vel_x
    displacement = disp_x
    execute_on = 'TIMESTEP_END'
  []
  [accel_y]
    type = TestNewmarkTI
    variable = accel_y
    displacement = disp_y
    first = false
    execute_on = 'TIMESTEP_END'
  []
  [vel_y]
    type = TestNewmarkTI
    variable = vel_y
    displacement = disp_x
    execute_on = 'TIMESTEP_END'
  []
  [accel_z]
    type = TestNewmarkTI
    variable = accel_z
    displacement = disp_z
    first = false
    execute_on = 'TIMESTEP_END'
  []
  [vel_z]
    type = TestNewmarkTI
    variable = vel_z
    displacement = disp_z
    execute_on = 'TIMESTEP_END'
  []
[]

[Kernels]
  [DynamicTensorMechanics]
    displacements = 'disp_x disp_y disp_z'
  []

  [Mass_x]
    type = MassMatrix
    variable = disp_x
    density = density
    matrix_tags = 'mass'
  []
  [Mass_y]
    type = MassMatrix
    variable = disp_y
    density = density
    matrix_tags = 'mass'
  []
  [Mass_z]
    type = MassMatrix
    variable = disp_z
    density = density
    matrix_tags = 'mass'
  []
[]

[Kernels]
  [gravity]
    type = Gravity
    variable = disp_z
    value = -981.0
    block = 1
  []
[]

[BCs]
  [x_front]
    type = ExplicitDirichletBC
    variable = disp_x
    boundary = 'ball_front'
    value = 0.0
  []
  [y_front]
    type = ExplicitDirichletBC
    variable = disp_y
    boundary = 'ball_front'
    value = 0.0
  []
  [x_fixed]
    type = ExplicitDirichletBC
    variable = disp_x
    boundary = 'base_back'
    value = 0.0
  []
  [y_fixed]
    type = ExplicitDirichletBC
    variable = disp_y
    boundary = 'base_back'
    value = 0.0
  []
  [z_fixed]
    type = ExplicitDirichletBC
    variable = disp_z
    boundary = 'base_back'
    value = 0.0
  []

[]

[ExplicitDynamicsContact]
  [my_contact]
    model = frictionless_balance
    primary = 'base_front ball_back'
    secondary = 'ball_back base_front'
    vel_x = 'vel_x'
    vel_y = 'vel_y'
    vel_z = 'vel_z'
  []
[]

[Materials]
  [elasticity_tensor_block_one]
    type = ComputeIsotropicElasticityTensor
    youngs_modulus = 1e1
    poissons_ratio = 0.3
    block = 1
    constant_on = SUBDOMAIN
  []
  [elasticity_tensor_block_two]
    type = ComputeIsotropicElasticityTensor
    youngs_modulus = 1e7
    poissons_ratio = 0.3
    block = 2
    constant_on = SUBDOMAIN
  []
  [strain_block]
    type = ComputeFiniteStrain
    displacements = 'disp_x disp_y disp_z'
    implicit = false
  []
  [stress_block]
    type = ComputeFiniteStrainElasticStress
  []
  [density_one]
    type = GenericConstantMaterial
    prop_names = density
    prop_values = 1e7
    output_properties = 'density'
    block = '1'
  []
  [density_two]
    type = GenericConstantMaterial
    prop_names = density
    prop_values = 1e3
    output_properties = 'density'
    block = '2'
  []
  [wave_speed]
    type = WaveSpeed
    outputs = 'exodus'
    output_properties = 'wave_speed'
  []
[]

[Executioner]
  type = Transient
  start_time = 0
  end_time = 0.04
  dt = 0.0001
  timestep_tolerance = 1e-6
  [TimeIntegrator]
    type = ExplicitMixedOrder
    mass_matrix_tag = 'mass'
    use_constant_mass = true
    second_order_vars = 'disp_x disp_y disp_z'
  []
  skip_exception_check = true
[]
[Outputs]
  time_step_interval = 100
  exodus = true
[]

(modules/contact/test/tests/explicit_dynamics/settlement.i)

# One element test to test the central difference time integrator in 3D.
[GlobalParams]
  displacements = 'disp_x disp_y disp_z'
  volumetric_locking_correction = true
[]

[Problem]
  extra_tag_matrices = 'mass'
[]

[Mesh]
  [block_one]
    type = GeneratedMeshGenerator
    dim = 3
    nx = 3
    ny = 3
    nz = 3
    xmin = 4.5
    xmax = 5.5
    ymin = 4.5
    ymax = 5.5
    zmin = 0.0001
    zmax = 1.0001
    boundary_name_prefix = 'ball'
  []
  [block_two]
    type = GeneratedMeshGenerator
    dim = 3
    nx = 2
    ny = 2
    nz = 2
    xmin = 0.0
    xmax = 10
    ymin = 0.0
    ymax = 10
    zmin = -2
    zmax = 0
    boundary_name_prefix = 'base'
    boundary_id_offset = 10
  []
  [block_one_id]
    type = SubdomainIDGenerator
    input = block_one
    subdomain_id = 1
  []
  [block_two_id]
    type = SubdomainIDGenerator
    input = block_two
    subdomain_id = 2
  []
  [combine]
    type = MeshCollectionGenerator
    inputs = ' block_one_id block_two_id'
  []
[]

[AuxVariables]
  [penetration]
  []
[]

[AuxKernels]
  [penetration]
    type = PenetrationAux
    variable = penetration
    boundary = ball_back
    paired_boundary = base_front
    quantity = distance
  []
[]

[Variables]
  [disp_x]
  []
  [disp_y]
  []
  [disp_z]
  []
[]

[AuxVariables]
  [gap_rate]
  []
  [vel_x]
  []
  [accel_x]
  []
  [vel_y]
  []
  [accel_y]
  []
  [vel_z]
  []
  [accel_z]
  []
  [stress_zz]
    family = MONOMIAL
    order = CONSTANT
  []
  [strain_zz]
    family = MONOMIAL
    order = CONSTANT
  []
  [kinetic_energy_one]
    order = CONSTANT
    family = MONOMIAL
  []
  [elastic_energy_one]
    order = CONSTANT
    family = MONOMIAL
  []
  [kinetic_energy_two]
    order = CONSTANT
    family = MONOMIAL
  []
  [elastic_energy_two]
    order = CONSTANT
    family = MONOMIAL
  []
[]

[AuxKernels]
  [stress_zz]
    type = RankTwoAux
    rank_two_tensor = stress
    index_i = 2
    index_j = 2
    variable = stress_zz
    execute_on = 'TIMESTEP_END'
  []
  [strain_zz]
    type = RankTwoAux
    rank_two_tensor = mechanical_strain
    index_i = 2
    index_j = 2
    variable = strain_zz
  []
  [accel_x]
    type = TestNewmarkTI
    variable = accel_x
    displacement = disp_x
    first = false
    execute_on = 'LINEAR TIMESTEP_BEGIN TIMESTEP_END'
  []
  [vel_x]
    type = TestNewmarkTI
    variable = vel_x
    displacement = disp_x
    execute_on = 'LINEAR TIMESTEP_BEGIN TIMESTEP_END'
  []
  [accel_y]
    type = TestNewmarkTI
    variable = accel_y
    displacement = disp_y
    first = false
    execute_on = 'LINEAR TIMESTEP_BEGIN TIMESTEP_END'
  []
  [vel_y]
    type = TestNewmarkTI
    variable = vel_y
    displacement = disp_x
    execute_on = 'LINEAR TIMESTEP_BEGIN TIMESTEP_END'
  []
  [accel_z]
    type = TestNewmarkTI
    variable = accel_z
    displacement = disp_z
    first = false
    execute_on = 'LINEAR TIMESTEP_BEGIN TIMESTEP_END'
  []
  [vel_z]
    type = TestNewmarkTI
    variable = vel_z
    displacement = disp_z
    execute_on = 'LINEAR TIMESTEP_BEGIN TIMESTEP_END'
  []
  [kinetic_energy_one]
    type = KineticEnergyAux
    block = '1'
    variable = kinetic_energy_one
    newmark_velocity_x = vel_x
    newmark_velocity_y = vel_y
    newmark_velocity_z = vel_z
    density = density
  []
  [elastic_energy_one]
    type = ElasticEnergyAux
    variable = elastic_energy_one
    block = '1'
  []
  [kinetic_energy_two]
    type = KineticEnergyAux
    block = '2'
    variable = kinetic_energy_two
    newmark_velocity_x = vel_x
    newmark_velocity_y = vel_y
    newmark_velocity_z = vel_z
    density = density
  []
  [elastic_energy_two]
    type = ElasticEnergyAux
    variable = elastic_energy_two
    block = '2'
  []
[]

[Kernels]
  [DynamicTensorMechanics]
    displacements = 'disp_x disp_y disp_z'
    stiffness_damping_coefficient = 9.5e-4
    generate_output = 'stress_zz strain_zz'
  []
  [Mass_x]
    type = MassMatrix
    variable = disp_x
    density = density
    matrix_tags = 'mass'
  []
  [Mass_y]
    type = MassMatrix
    variable = disp_y
    density = density
    matrix_tags = 'mass'
  []
  [Mass_z]
    type = MassMatrix
    variable = disp_z
    density = density
    matrix_tags = 'mass'
  []
  [gravity]
    type = Gravity
    variable = disp_z
    value = -98.10
    block = 1
  []
[]

[BCs]
  [x_front]
    type = ExplicitDirichletBC
    variable = disp_x
    boundary = 'ball_front'
    value = 0.0
  []
  [y_front]
    type = ExplicitDirichletBC
    variable = disp_y
    boundary = 'ball_front'
    value = 0.0
  []
  [x_fixed]
    type = ExplicitDirichletBC
    variable = disp_x
    boundary = 'base_back'
    value = 0.0
  []
  [y_fixed]
    type = ExplicitDirichletBC
    variable = disp_y
    boundary = 'base_back'
    value = 0.0
  []
  [z_fixed]
    type = ExplicitDirichletBC
    variable = disp_z
    boundary = 'base_back'
    value = 0.0
  []
  [z_fixed_front]
    type = ExplicitDirichletBC
    variable = disp_z
    boundary = 'base_front'
    value = 0.0
  []
[]

[ExplicitDynamicsContact]
  [my_contact]
    model = frictionless_balance
    primary = base_front
    secondary = ball_back
    vel_x = 'vel_x'
    vel_y = 'vel_y'
    vel_z = 'vel_z'
    verbose = true
  []
[]

[Materials]
  [elasticity_tensor_block_one]
    type = ComputeIsotropicElasticityTensor
    youngs_modulus = 1e6
    poissons_ratio = 0.0
    block = 1
    outputs = 'exodus'
    output_properties = __all__
  []
  [elasticity_tensor_block_two]
    type = ComputeIsotropicElasticityTensor
    youngs_modulus = 1e10
    poissons_ratio = 0.0
    block = 2
    outputs = 'exodus'
    output_properties = __all__
  []
  [strain_block]
    type = ComputeFiniteStrain # ComputeIncrementalStrain
    displacements = 'disp_x disp_y disp_z'
    implicit = false
  []
  [stress_block]
    type = ComputeFiniteStrainElasticStress
  []
  [density_one]
    type = GenericConstantMaterial
    prop_names = density
    prop_values = 1e1
    outputs = 'exodus'
    output_properties = 'density'
    block = '1'
  []
  [density_two]
    type = GenericConstantMaterial
    prop_names = density
    prop_values = 1e6
    outputs = 'exodus'
    output_properties = 'density'
    block = '2'
  []
  [wave_speed]
    type = WaveSpeed
    outputs = 'exodus'
    output_properties = 'wave_speed'
  []
[]

[Executioner]
  type = Transient
  start_time = 0
  end_time = 0.05
  dt = 1.0e-4
  timestep_tolerance = 1e-6

  [TimeIntegrator]
    type = ExplicitMixedOrder
    mass_matrix_tag = 'mass'
    second_order_vars = 'disp_x disp_y disp_z'
  []
[]

[Outputs]
  time_step_interval = 1
  exodus = true
  csv = true
  execute_on = 'TIMESTEP_END'
[]

[Postprocessors]
  [accel_58z]
    type = NodalVariableValue
    nodeid = 1
    variable = accel_z
  []
  [vel_58z]
    type = NodalVariableValue
    nodeid = 1
    variable = vel_z
  []
  [critical_time_step]
    type = CriticalTimeStep
  []
  [contact_pressure_max]
    type = NodalExtremeValue
    variable = contact_pressure
    block = '1 2'
    value_type = max
  []
  [penetration_max]
    type = NodalExtremeValue
    variable = penetration
    block = '1 2'
    value_type = max
  []
  [total_kinetic_energy_one]
    type = ElementIntegralVariablePostprocessor
    variable = kinetic_energy_one
    block = '1'
  []
  [total_elastic_energy_one]
    type = ElementIntegralVariablePostprocessor
    variable = elastic_energy_one
    block = '1'
  []
  [total_kinetic_energy_two]
    type = ElementIntegralVariablePostprocessor
    variable = kinetic_energy_two
    block = '2'
  []
  [total_elastic_energy_two]
    type = ElementIntegralVariablePostprocessor
    variable = elastic_energy_two
    block = '2'
  []
[]

(test/tests/tag/fe-mass-matrix.i)

[Mesh]
  [square]
    type = GeneratedMeshGenerator
    nx = 2
    ny = 2
    dim = 2
  []
[]

[Variables]
  [u]
  []
[]

[AuxVariables]
  [mass]
  []
[]

[Kernels]
  [diff]
    type = Diffusion
    variable = u
  []
  [mass]
    type = MassMatrix
    variable = u
    matrix_tags = 'mass'
  []
[]

[AuxKernels]
  [TagMatrixAux1]
    type = TagMatrixAux
    variable = mass
    v = u
    matrix_tag = mass
  []
[]

[BCs]
  [left]
    type = DirichletBC
    variable = u
    boundary = left
    value = 0
  []
  [right]
    type = DirichletBC
    variable = u
    boundary = right
    value = 1
  []
[]

[Problem]
  type = FEProblem
  extra_tag_matrices = 'mass'
[]

[Executioner]
  type = Steady
  solve_type = 'NEWTON'
[]

[Outputs]
  exodus = true
[]

(modules/contact/test/tests/explicit_dynamics/test_balance_optimized.i)

[GlobalParams]
  displacements = 'disp_x disp_y disp_z'
  volumetric_locking_correction = true
[]

[Problem]
  extra_tag_matrices = 'mass'
[]

[Mesh]
  [block_one]
    type = GeneratedMeshGenerator
    dim = 3
    nx = 3
    ny = 3
    nz = 3
    xmin = 4.5
    xmax = 5.5
    ymin = 4.5
    ymax = 5.5
    zmin = 0.0001
    zmax = 1.0001
    boundary_name_prefix = 'ball'
  []
  [block_two]
    type = GeneratedMeshGenerator
    dim = 3
    nx = 2
    ny = 2
    nz = 2
    xmin = 0.0
    xmax = 10
    ymin = 0.0
    ymax = 10
    zmin = -2
    zmax = 0
    boundary_name_prefix = 'base'
    boundary_id_offset = 10
  []
  [block_one_id]
    type = SubdomainIDGenerator
    input = block_one
    subdomain_id = 1
  []
  [block_two_id]
    type = SubdomainIDGenerator
    input = block_two
    subdomain_id = 2
  []
  [combine]
    type = MeshCollectionGenerator
    inputs = ' block_one_id block_two_id'
  []
  allow_renumbering = false
  patch_update_strategy = auto
[]

[Variables]
  [disp_x]
  []
  [disp_y]
  []
  [disp_z]
  []
[]

[AuxVariables]
  [gap_rate]
  []
  [vel_x]
  []
  [accel_x]
  []
  [vel_y]
  []
  [accel_y]
  []
  [vel_z]
  []
  [accel_z]
  []
[]

[AuxKernels]
  [accel_x]
    type = TestNewmarkTI
    variable = accel_x
    displacement = disp_x
    first = false
    execute_on = 'TIMESTEP_END'
  []
  [vel_x]
    type = TestNewmarkTI
    variable = vel_x
    displacement = disp_x
    execute_on = 'TIMESTEP_END'
  []
  [accel_y]
    type = TestNewmarkTI
    variable = accel_y
    displacement = disp_y
    first = false
    execute_on = 'TIMESTEP_END'
  []
  [vel_y]
    type = TestNewmarkTI
    variable = vel_y
    displacement = disp_x
    execute_on = 'TIMESTEP_END'
  []
  [accel_z]
    type = TestNewmarkTI
    variable = accel_z
    displacement = disp_z
    first = false
    execute_on = 'TIMESTEP_END'
  []
  [vel_z]
    type = TestNewmarkTI
    variable = vel_z
    displacement = disp_z
    execute_on = 'TIMESTEP_END'
  []
[]

[Kernels]
  [DynamicTensorMechanics]
    displacements = 'disp_x disp_y disp_z'
    volumetric_locking_correction = true
    stiffness_damping_coefficient = 0.001
    decomposition_method = EigenSolution
  []

  [Mass_x]
    type = MassMatrix
    variable = disp_x
    density = density
    matrix_tags = 'mass'
  []
  [Mass_y]
    type = MassMatrix
    variable = disp_y
    density = density
    matrix_tags = 'mass'
  []
  [Mass_z]
    type = MassMatrix
    variable = disp_z
    density = density
    matrix_tags = 'mass'
  []
[]

[Kernels]
  [gravity]
    type = Gravity
    variable = disp_z
    value = -981.0
  []
[]

[BCs]
  [x_front]
    type = ExplicitDirichletBC
    variable = disp_x
    boundary = 'ball_front'
    value = 0.0
  []
  [y_front]
    type = ExplicitDirichletBC
    variable = disp_y
    boundary = 'ball_front'
    value = 0.0
  []
  [x_fixed]
    type = ExplicitDirichletBC
    variable = disp_x
    boundary = 'base_back'
    value = 0.0
  []
  [y_fixed]
    type = ExplicitDirichletBC
    variable = disp_y
    boundary = 'base_back'
    value = 0.0
  []
  [z_fixed]
    type = ExplicitDirichletBC
    variable = disp_z
    boundary = 'base_back'
    value = 0.0
  []
  [z_fixed_front]
    type = ExplicitDirichletBC
    variable = disp_z
    boundary = 'base_front'
    value = 0.0
  []
[]

[ExplicitDynamicsContact]
  [my_contact]
    model = frictionless_balance
    primary = base_front
    secondary = ball_back
    vel_x = 'vel_x'
    vel_y = 'vel_y'
    vel_z = 'vel_z'
  []
[]

[Materials]
  [elasticity_tensor_block_one]
    type = ComputeIsotropicElasticityTensor
    youngs_modulus = 1e6
    poissons_ratio = 0.0
    block = 1
    constant_on = SUBDOMAIN
  []
  [elasticity_tensor_block_two]
    type = ComputeIsotropicElasticityTensor
    youngs_modulus = 1e10
    poissons_ratio = 0.0
    block = 2
    constant_on = SUBDOMAIN
  []
  [strain_block]
    type = ComputeFiniteStrain
    displacements = 'disp_x disp_y disp_z'
    implicit = false
  []
  [stress_block]
    type = ComputeFiniteStrainElasticStress
  []
  [density_one]
    type = GenericConstantMaterial
    prop_names = density
    prop_values = 1e1
    output_properties = 'density'
    block = '1'
  []
  [density_two]
    type = GenericConstantMaterial
    prop_names = density
    prop_values = 1e6
    output_properties = 'density'
    block = '2'
  []
  [wave_speed]
    type = WaveSpeed
    outputs = 'exodus'
    output_properties = 'wave_speed'
  []
[]

[Executioner]
  type = Transient
  start_time = 0
  end_time = 0.0025
  dt = 0.00001
  timestep_tolerance = 1e-6
  [TimeIntegrator]
    type = ExplicitMixedOrder
    mass_matrix_tag = 'mass'
    use_constant_mass = true
    second_order_vars = 'disp_x disp_y disp_z'
  []
[]
[Outputs]
  time_step_interval = 10
  exodus = true
  csv = true
  file_base = test_balance_out
[]

(modules/navier_stokes/test/tests/finite_element/ins/hdg/ip/lid-driven/matrix-testing.i)

mu = 1
rho = 1
U = 1
n = 5
l = 1

[Mesh]
  [gen]
    type = GeneratedMeshGenerator
    dim = 2
    xmin = 0
    xmax = ${l}
    ymin = 0
    ymax = ${l}
    nx = ${n}
    ny = ${n}
    elem_type = TRI6
  []
[]

[Problem]
  type = PrintMatricesNSProblem
  extra_tag_matrices = 'mass jump combined grad_div'
  pressure_mass_matrix = 'mass'
  velocity_mass_matrix = 'mass'
  augmented_lagrange_matrices = 'jump grad_div combined'
  u = vel_x
  v = vel_y
  pressure = pressure
  pressure_bar = pressure_bar
  print = false
[]

[Variables]
  [vel_x]
    family = MONOMIAL
    order = FIRST
  []
  [vel_y]
    family = MONOMIAL
    order = FIRST
  []
  [pressure]
    family = MONOMIAL
    order = CONSTANT
  []
  [vel_bar_x]
    family = SIDE_HIERARCHIC
    order = FIRST
  []
  [vel_bar_y]
    family = SIDE_HIERARCHIC
    order = FIRST
  []
  [pressure_bar]
    family = SIDE_HIERARCHIC
    order = FIRST
  []
  [lambda]
    family = SCALAR
    order = FIRST
  []
[]

[HDGKernels]
  [momentum_x_diffusion]
    type = NavierStokesStressIPHDGKernel
    variable = vel_x
    face_variable = vel_bar_x
    diffusivity = 'mu'
    alpha = 6
    pressure_variable = pressure
    pressure_face_variable = pressure_bar
    component = 0
  []

  [momentum_y_diffusion]
    type = NavierStokesStressIPHDGKernel
    variable = vel_y
    face_variable = vel_bar_y
    diffusivity = 'mu'
    alpha = 6
    pressure_variable = pressure
    pressure_face_variable = pressure_bar
    component = 1
  []

  [pressure_convection]
    type = AdvectionIPHDGKernel
    variable = pressure
    face_variable = pressure_bar
    velocity = 'velocity'
    coeff = '${fparse -rho}'
    self_advection = false
  []
[]

[Kernels]
  [mean_zero_pressure]
    type = ScalarLagrangeMultiplier
    variable = pressure
    lambda = lambda
  []

  [mass_matrix_vel_x]
    type = MassMatrix
    variable = vel_x
    matrix_tags = 'mass'
  []
  [mass_matrix_vel_y]
    type = MassMatrix
    variable = vel_y
    matrix_tags = 'mass'
  []
  [mass_matrix_pressure]
    type = MassMatrix
    variable = pressure
    matrix_tags = 'mass'
  []

  [u_jump]
    type = GradDiv
    matrix_only = true
    variable = vel_x
    u = vel_x
    v = vel_y
    component = 0
    vector_tags = ''
    matrix_tags = 'grad_div combined'
  []
  [v_jump]
    type = GradDiv
    matrix_only = true
    variable = vel_y
    u = vel_x
    v = vel_y
    component = 1
    vector_tags = ''
    matrix_tags = 'grad_div combined'
  []
[]

[ScalarKernels]
  [mean_zero_pressure_lm]
    type = AverageValueConstraint
    variable = lambda
    pp_name = pressure_integral
    value = 0
  []
[]

[DGKernels]
  [pb_mass]
    type = MassMatrixDGKernel
    variable = pressure_bar
    matrix_tags = 'mass'
  []

  [u_jump]
    type = MassFluxPenalty
    matrix_only = true
    variable = vel_x
    u = vel_x
    v = vel_y
    component = 0
    vector_tags = ''
    matrix_tags = 'jump combined'
  []
  [v_jump]
    type = MassFluxPenalty
    matrix_only = true
    variable = vel_y
    u = vel_x
    v = vel_y
    component = 1
    vector_tags = ''
    matrix_tags = 'jump combined'
  []
[]

[BCs]
  [momentum_x_diffusion_walls]
    type = NavierStokesStressIPHDGDirichletBC
    boundary = 'left bottom right'
    variable = vel_x
    face_variable = vel_bar_x
    pressure_variable = pressure
    pressure_face_variable = pressure_bar
    alpha = 6
    functor = '0'
    diffusivity = 'mu'
    component = 0
  []
  [momentum_x_diffusion_top]
    type = NavierStokesStressIPHDGDirichletBC
    boundary = 'top'
    variable = vel_x
    face_variable = vel_bar_x
    pressure_variable = pressure
    pressure_face_variable = pressure_bar
    alpha = 6
    functor = '${U}'
    diffusivity = 'mu'
    component = 0
  []
  [momentum_y_diffusion_all]
    type = NavierStokesStressIPHDGDirichletBC
    boundary = 'left bottom right top'
    variable = vel_y
    face_variable = vel_bar_y
    pressure_variable = pressure
    pressure_face_variable = pressure_bar
    alpha = 6
    functor = '0'
    diffusivity = 'mu'
    component = 1
  []

  [mass_convection]
    type = AdvectionIPHDGPrescribedFluxBC
    face_variable = pressure_bar
    variable = pressure
    velocity = 'velocity'
    coeff = '${fparse -rho}'
    self_advection = false
    boundary = 'left bottom top right'
    prescribed_normal_flux = 0
  []

  [pb_mass]
    type = MassMatrixIntegratedBC
    variable = pressure_bar
    matrix_tags = 'mass'
    boundary = 'left right bottom top'
  []

  [u_jump_walls]
    type = MassFluxPenaltyBC
    matrix_only = true
    variable = vel_x
    u = vel_x
    v = vel_y
    component = 0
    vector_tags = ''
    matrix_tags = 'jump combined'
    boundary = 'left right bottom'
    dirichlet_value = 'walls'
  []
  [v_jump_walls]
    type = MassFluxPenaltyBC
    matrix_only = true
    variable = vel_y
    u = vel_x
    v = vel_y
    component = 1
    vector_tags = ''
    matrix_tags = 'jump combined'
    boundary = 'left right bottom'
    dirichlet_value = 'walls'
  []
  [u_jump_top]
    type = MassFluxPenaltyBC
    matrix_only = true
    variable = vel_x
    u = vel_x
    v = vel_y
    component = 0
    vector_tags = ''
    matrix_tags = 'jump combined'
    boundary = 'top'
    dirichlet_value = 'top'
  []
  [v_jump_top]
    type = MassFluxPenaltyBC
    matrix_only = true
    variable = vel_y
    u = vel_x
    v = vel_y
    component = 1
    vector_tags = ''
    matrix_tags = 'jump combined'
    boundary = 'top'
    dirichlet_value = 'top'
  []
[]

[Functions]
  [top]
    type = ParsedVectorFunction
    value_x = ${U}
    value_y = 0
  []
  [walls]
    type = ParsedVectorFunction
    value_x = 0
    value_y = 0
  []
[]

[Materials]
  [const]
    type = ADGenericConstantMaterial
    prop_names = 'rho mu'
    prop_values = '${rho} ${mu}'
  []
  [vel]
    type = ADVectorFromComponentVariablesMaterial
    vector_prop_name = 'velocity'
    u = vel_x
    v = vel_y
  []
[]

[Executioner]
  type = Steady
  solve_type = 'NEWTON'
  petsc_options_iname = '-pc_type -pc_factor_shift_type'
  petsc_options_value = 'lu       NONZERO'
  nl_rel_tol = 1e-12
[]

[Outputs]
  [out]
    type = CSV
    hide = 'pressure_integral lambda'
    execute_on = 'timestep_end'
  []
[]

[Postprocessors]
  [symmetric]
    type = MatrixSymmetryCheck
    execute_on = 'timestep_end'
  []
  [pressure_integral]
    type = ElementIntegralVariablePostprocessor
    variable = pressure
    execute_on = linear
  []
  [fe_jump_and_upb_equiv]
    type = MatrixEqualityCheck
    mat1 = 'vel_pb_grad_div.mat'
    mat2 = 'jump.mat'
  []
  [fe_grad_div_and_up_equiv]
    type = MatrixEqualityCheck
    mat1 = 'vel_p_grad_div.mat'
    mat2 = 'grad_div.mat'
  []
  [fe_combined_and_upall_equiv]
    type = MatrixEqualityCheck
    mat1 = 'vel_all_p_grad_div.mat'
    mat2 = 'combined.mat'
  []
  [upb_grad_div_num_zero_eig]
    type = MatrixEigenvalueCheck
    mat = vel_pb_grad_div.mat
  []
  [upb_div_grad_num_zero_eig]
    type = MatrixEigenvalueCheck
    mat = vel_pb_div_grad.mat
  []
  [up_grad_div_num_zero_eig]
    type = MatrixEigenvalueCheck
    mat = vel_p_grad_div.mat
  []
  [up_div_grad_num_zero_eig]
    type = MatrixEigenvalueCheck
    mat = vel_p_div_grad.mat
  []
  [upall_grad_div_num_zero_eig]
    type = MatrixEigenvalueCheck
    mat = vel_all_p_grad_div.mat
  []
  [upall_div_grad_num_zero_eig]
    type = MatrixEigenvalueCheck
    mat = vel_all_p_div_grad.mat
  []
[]

(modules/solid_mechanics/test/tests/dynamics/time_integration/direct_central_difference_multiVarBC.i)

###########################################################
# This is a simple test with a time-dependent problem
# demonstrating the use of a central difference with a
# direct calculation of acceleration.
#
# Testing that the first and second time derivatives
# are calculated correctly using the Central Difference Direct
# method
###########################################################

[GlobalParams]
  displacements = 'disp_x disp_y'
[]

[Problem]
  extra_tag_matrices = 'mass'
[]

[Mesh]
  type = GeneratedMesh
  dim = 2
  xmin = -1
  xmax = 1
  ymin = -1
  ymax = 1
  nx = 4
  ny = 4
[]

[Variables]
  [disp_x]
  []
  [disp_y]
  []
[]

[Functions]
  [forcing_fn]
    type = ParsedFunction
    expression = 't'
  []
[]

[Kernels]
  [DynamicSolidMechanics]
    displacements = 'disp_x disp_y'
  []
  [massmatrix]
    type = MassMatrix
    density = density
    matrix_tags = 'mass'
    variable = disp_x
  []
  [massmatrix_y]
    type = MassMatrix
    density = density
    matrix_tags = 'mass'
    variable = disp_y
  []
[]

[Materials]
  [elasticity_tensor_block_one]
    type = ComputeIsotropicElasticityTensor
    youngs_modulus = 1e1
    poissons_ratio = 0.33
  []

  [strain_block]
    type = ComputeFiniteStrain
    displacements = 'disp_x disp_y'
    implicit = false
  []
  [stress_block]
    type = ComputeFiniteStrainElasticStress
    implicit = false
  []
  [density]
    type = GenericConstantMaterial
    prop_names = 'density'
    prop_values = 1
  []
[]

[BCs]
  [left_x]
    type = ExplicitDirichletBC
    variable = disp_x
    boundary = 'left'
    value = 0
  []
  [left_y]
    type = ExplicitDirichletBC
    variable = disp_y
    boundary = 'left'
    value = 0
  []
  [right_x]
    type = ExplicitFunctionDirichletBC
    variable = disp_x
    boundary = 'right'
    function = forcing_fn
  []
[]

[Executioner]
  type = Transient

  [TimeIntegrator]
    type = ExplicitMixedOrder
    mass_matrix_tag = 'mass'
    second_order_vars = 'disp_x disp_y'
  []
  start_time = 0.0
  num_steps = 6
  dt = 0.1
[]

[Postprocessors]
  [left_x]
    type = AverageNodalVariableValue
    variable = disp_x
    boundary = left
  []
  [right_y]
    type = AverageNodalVariableValue
    variable = disp_x
    boundary = left
  []
[]

[Outputs]
  exodus = true
[]

(test/tests/geomsearch/3d_penetration_locator/3d_disconnected_tet.i)

[Mesh]
  [connected_mesh]
    type = FileMeshGenerator
    file = 3d_thermal_contact_tet.e
  []
  [exploded_mesh]
    type = BreakMeshByElementGenerator
    input = connected_mesh
    interface_name = 'interelement'
  []
  construct_side_list_from_node_list = true
  allow_renumbering = false # fix VPP ordering
[]

[Variables]
  active = 'u'

  [./u]
    order = FIRST
    family = LAGRANGE
  [../]
[]

[AuxVariables]
  [./gap_distance]
    order = FIRST
    family = LAGRANGE
  [../]
[]

[Kernels]
  [./l2]
    type = MassMatrix
    variable = u
    matrix_tags = 'system'
  [../]
[]

[AuxKernels]
  [./distance]
    type = PenetrationAux
    variable = gap_distance
    boundary = 'leftright'
    paired_boundary = 'rightleft'
    search_method = all_proximate_sides
  [../]
[]

[VectorPostprocessors]
  [gap_sampler]
    type = SideValueSampler
    boundary = 'leftright'
    variable = 'gap_distance'
    sort_by = 'id'
  []
[]

[Executioner]
  type = Steady

  solve_type = 'PJFNK'
[]

[Outputs]
  csv = true
[]

(modules/navier_stokes/test/tests/finite_element/ins/lid_driven/steady_vector_fsp.i)

rho=1
mu=1
U=1
l=1
prefactor=${fparse 1/(l/2)^2}
n=8

[Mesh]
  [gen]
    type = DistributedRectilinearMeshGenerator
    dim = 2
    xmin = 0
    xmax = ${l}
    ymin = 0
    ymax = ${l}
    nx = ${n}
    ny = ${n}
    elem_type = QUAD4
  []
  second_order = true
  parallel_type = distributed
[]

[Variables]
  [vel]
    order = SECOND
    family = LAGRANGE_VEC
  []
  [p]
    order = FIRST
    family = LAGRANGE
  []
[]

[Kernels]
  [mass]
    type = INSADMass
    variable = p
  []
  [mass_kernel]
    type = MassMatrix
    variable = p
    matrix_tags = 'mass'
  []
  [momentum_convection]
    type = INSADMomentumAdvection
    variable = vel
  []
  [momentum_viscous]
    type = INSADMomentumViscous
    variable = vel
    extra_matrix_tags = 'L'
  []
  [momentum_pressure]
    type = INSADMomentumPressure
    variable = vel
    pressure = p
    integrate_p_by_parts = true
  []
[]

[BCs]
  [no_slip]
    type = ADVectorFunctionDirichletBC
    variable = vel
    boundary = 'bottom right left'
    extra_matrix_tags = 'L'
  []
  [lid]
    type = ADVectorFunctionDirichletBC
    variable = vel
    boundary = 'top'
    function_x = 'lid_function'
    extra_matrix_tags = 'L'
  []
[]

[Materials]
  [const]
    type = ADGenericConstantMaterial
    prop_names = 'rho mu'
    prop_values = '${rho} ${mu}'
  []
  [insad]
    type = INSADMaterial
    velocity = vel
    pressure = p
  []
[]

[Functions]
  [lid_function]
    # We pick a function that is exactly represented in the velocity
    # space so that the Dirichlet conditions are the same regardless
    # of the mesh spacing.
    type = ParsedFunction
    expression = '${prefactor}*${U}*x*(${l}-x)'
  []
[]

[Problem]
  type = NavierStokesProblem
  mass_matrix = 'mass'
  extra_tag_matrices = 'mass L'
  L_matrix = 'L'
  commute_lsc = true
[]

[Preconditioning]
  [FSP]
    type = FSP
    topsplit = 'up'
    [up]
      splitting = 'u p'
      splitting_type  = schur
      petsc_options_iname = '-pc_fieldsplit_schur_fact_type  -pc_fieldsplit_schur_precondition -ksp_gmres_restart -ksp_type -ksp_pc_side -ksp_rtol'
      petsc_options_value = 'full                            self                              300                fgmres    right        1e-4'
    []
      [u]
        vars = 'vel'
        # petsc_options = '-ksp_converged_reason'
        petsc_options_iname = '-pc_type -pc_hypre_type -ksp_type -ksp_rtol -ksp_gmres_restart -ksp_pc_side'
        petsc_options_value = 'hypre    boomeramg      gmres     1e-2      300                right'
      []
      [p]
        vars = 'p'
        petsc_options = '-ksp_converged_reason -pc_lsc_commute'
        petsc_options_iname = '-ksp_type -ksp_gmres_restart -ksp_rtol -pc_type -ksp_pc_side -lsc_pc_type -lsc_pc_hypre_type -lsc_ksp_type -lsc_ksp_rtol -lsc_ksp_pc_side -lsc_ksp_gmres_restart -lsc_mass_pc_type -lsc_mass_pc_hypre_type -lsc_mass_ksp_rtol -lsc_mass_ksp_type'
        petsc_options_value = 'fgmres    300                1e-2      lsc      right        hypre        boomeramg          fgmres         1e-1         right            300                    hypre             boomeramg               1e-1               gmres'
      []
  []
[]

[Postprocessors]
  [pavg]
    type = ElementAverageValue
    variable = p
  []
[]

[Correctors]
  [set_pressure]
    type = NSPressurePin
    pin_type = 'average'
    variable = p
    pressure_average = 'pavg'
  []
[]

[Executioner]
  type = Steady
  solve_type = NEWTON
[]

[Outputs]
  print_linear_residuals = false
  [exo]
    type = Exodus
    execute_on = 'final'
    hide = 'pavg'
    file_base = 'fsp_steady_low_Re_olshanskii'
  []
[]

(modules/solid_mechanics/test/tests/dynamics/time_integration/direct_central_difference.i)

###########################################################
# This is a simple test with a time-dependent problem
# demonstrating the use of a central difference with a
# direct calculation of acceleration.
#
# Testing that the first and second time derivatives
# are calculated correctly using the Central Difference Direct
# method
###########################################################

[GlobalParams]
  displacements = 'disp_x disp_y'
[]

[Problem]
  extra_tag_matrices = 'mass'
[]

[Mesh]
  type = GeneratedMesh
  dim = 2
  xmin = -1
  xmax = 1
  ymin = -1
  ymax = 1
  nx = 1
  ny = 1
[]

[Variables]
  [disp_x]
  []
  [disp_y]
  []
[]

[Functions]
  [forcing_fn]
    type = PiecewiseLinear
    x = '0.0 0.1 0.2    0.3  0.4    0.5  0.6'
    y = '0.0 0.0 0.0025 0.01 0.0175 0.02 0.02'
  []
[]

[Kernels]
  [DynamicSolidMechanics]
    displacements = 'disp_x disp_y'
  []
  [massmatrix]
    type = MassMatrix
    density = density
    matrix_tags = 'mass'
    variable = disp_x
  []
  [massmatrix_y]
    type = MassMatrix
    density = density
    matrix_tags = 'mass'
    variable = disp_y
  []
[]

[Materials]
  [elasticity_tensor_block_one]
    type = ComputeIsotropicElasticityTensor
    youngs_modulus = 1e1
    poissons_ratio = 0.0
  []

  [strain_block]
    type = ComputeFiniteStrain
    displacements = 'disp_x disp_y'
    implicit = false
  []
  [stress_block]
    type = ComputeFiniteStrainElasticStress
    implicit = false
  []
  [density]
    type = GenericConstantMaterial
    prop_names = 'density'
    prop_values = 1
  []
[]

[BCs]
  [left_x]
    type = ExplicitFunctionDirichletBC
    variable = disp_x
    boundary = 'left'
    function = forcing_fn
  []
  [right_x]
    type = ExplicitFunctionDirichletBC
    variable = disp_x
    boundary = 'right'
    function = forcing_fn
  []
[]

[Executioner]
  type = Transient

  [TimeIntegrator]
    type = ExplicitMixedOrder
    mass_matrix_tag = 'mass'
    second_order_vars = 'disp_x disp_y'
  []

  start_time = 0.0
  num_steps = 6
  dt = 0.1
[]

[Postprocessors]
  [udot]
    type = ElementAverageTimeDerivative
    variable = disp_x
  []
  [udotdot]
    type = ElementAverageSecondTimeDerivative
    variable = disp_x
  []
  [u]
    type = ElementAverageValue
    variable = disp_x
  []
[]

[Outputs]
  exodus = true
[]

(modules/solid_mechanics/test/tests/dynamics/explicit_mms/mms_direct_second_order_with_ic.i)

[Mesh]
  type = GeneratedMesh
  dim = 2
  xmin = -1
  xmax = 1
  ymin = -1
  ymax = 1
  nx = 10
  ny = 10
  elem_type = QUAD4
[]

dt_and_v0 = 0.00005

[Problem]
  extra_tag_matrices = 'mass'
[]

[Variables]
  [u]
  []
[]

[AuxVariables]
  [dot_u]
  []
[]

[ICs]
  [u_old]
    type = ConstantIC
    variable = u
    state = OLD
    # set's v_0 to 1
    value = -${dt_and_v0}
  []
[]

[AuxKernels]
  [dot_u]
    type = TestNewmarkTI
    variable = dot_u
    displacement = u
    first = true
    execute_on = 'TIMESTEP_END'
  []
[]

[Functions]
  [exact_fn]
    type = ParsedFunction
    expression = 't + 0.5*t^2'
  []
  [exact_dot_fn]
    type = ParsedFunction
    expression = '1 + t'
  []
[]

[Kernels]
  [Mass_x]
    type = MassMatrix
    variable = u
    density = 1
    matrix_tags = 'mass'
  []
  [ffn]
    type = BodyForce
    variable = u
    function = 1
  []
[]

[Postprocessors]
  [l2_err]
    type = ElementL2Error
    variable = u
    function = exact_fn
  []
  [l2_dot_err]
    type = ElementL2Error
    variable = dot_u
    function = exact_dot_fn
  []
[]

[Executioner]
  type = Transient
  start_time = 0.0
  num_steps = 20
  dt = ${dt_and_v0}
  l_tol = 1e-12
  [TimeIntegrator]
    type = ExplicitMixedOrder
    mass_matrix_tag = 'mass'
    use_constant_mass = true
    second_order_vars = 'u'
  []
[]

[Outputs]
  exodus = true
  [console]
    type = Console
    max_rows = 10
  []
[]

(modules/navier_stokes/test/tests/finite_element/ins/lid_driven/steady_vector_fsp_al.i)

rho=1
mu=1e-3
U=1
l=1
prefactor=${fparse 1/(l/2)^2}
n=8
gamma=${U}

[Mesh]
  [gen]
    type = GeneratedMeshGenerator
    dim = 2
    xmin = 0
    xmax = ${l}
    ymin = 0
    ymax = ${l}
    nx = ${n}
    ny = ${n}
    elem_type = QUAD4
  []
  second_order = true
[]

[Variables]
  [vel]
    order = SECOND
    family = LAGRANGE_VEC
  []
  [p]
    order = FIRST
    family = LAGRANGE
  []
[]

[Kernels]
  [mass]
    type = INSADMass
    variable = p
  []
  [mass_kernel]
    type = MassMatrix
    variable = p
    matrix_tags = 'mass'
    density = ${fparse -1/(gamma + mu)}
  []
  [momentum_advection]
    type = INSADMomentumAdvection
    variable = vel
  []
  [momentum_viscous]
    type = INSADMomentumViscous
    variable = vel
  []
  [momentum_pressure]
    type = INSADMomentumPressure
    variable = vel
    pressure = p
    integrate_p_by_parts = true
  []
  [momentum_graddiv]
    type = INSADMomentumGradDiv
    variable = vel
    gamma = ${gamma}
  []
[]

[BCs]
  [no_slip]
    type = ADVectorFunctionDirichletBC
    variable = vel
    boundary = 'bottom right left'
    preset = true
  []
  [lid]
    type = ADVectorFunctionDirichletBC
    variable = vel
    boundary = 'top'
    function_x = 'lid_function'
    preset = true
  []
[]

[Materials]
  [const]
    type = ADGenericConstantMaterial
    prop_names = 'rho mu'
    prop_values = '${rho} ${mu}'
  []
  [insad]
    type = INSADTauMaterial
    velocity = vel
    pressure = p
  []
[]

[Functions]
  [lid_function]
    # We pick a function that is exactly represented in the velocity
    # space so that the Dirichlet conditions are the same regardless
    # of the mesh spacing.
    type = ParsedFunction
    expression = '${prefactor}*${U}*x*(${l}-x)'
  []
[]

[Problem]
  type = NavierStokesProblem
  mass_matrix = 'mass'
  extra_tag_matrices = 'mass'
  use_pressure_mass_matrix = true
[]

[Preconditioning]
  [FSP]
    type = FSP
    topsplit = 'up'
    [up]
      splitting = 'u p'
      splitting_type  = schur
      petsc_options_iname = '-pc_fieldsplit_schur_fact_type  -pc_fieldsplit_schur_precondition -ksp_gmres_restart -ksp_type -ksp_pc_side -ksp_rtol'
      petsc_options_value = 'full                            self                              300                fgmres    right        1e-4'
    []
      [u]
        vars = 'vel'
        petsc_options = '-ksp_converged_reason'
        petsc_options_iname = '-pc_type -ksp_type -ksp_rtol -ksp_gmres_restart -ksp_pc_side -pc_factor_mat_solver_type'
        petsc_options_value = 'ilu      gmres     1e-2      300                right        strumpack'
      []
      [p]
        vars = 'p'
        petsc_options = '-ksp_converged_reason'
        petsc_options_iname = '-ksp_type -ksp_gmres_restart -ksp_rtol -pc_type -ksp_pc_side'
        petsc_options_value = 'gmres     300                1e-2      ilu      right'
      []
  []
[]

[Postprocessors]
  [pavg]
    type = ElementAverageValue
    variable = p
  []
[]

[Correctors]
  [set_pressure]
    type = NSPressurePin
    pin_type = 'average'
    variable = p
    pressure_average = 'pavg'
  []
[]

[Executioner]
  type = Steady
  solve_type = NEWTON
  nl_rel_tol = 1e-12
[]

[Outputs]
  print_linear_residuals = false
  [exo]
    type = Exodus
    execute_on = 'final'
    hide = 'pavg'
  []
[]

(modules/contact/test/tests/explicit_dynamics/highvel.i)

# One element test to test the central difference time integrator in 3D.
[GlobalParams]
  displacements = 'disp_x disp_y disp_z'
  volumetric_locking_correction = true
[]

[Problem]
  extra_tag_matrices = 'mass'
[]

[Mesh]
  [block_one]
    type = GeneratedMeshGenerator
    dim = 3
    nx = 3
    ny = 3
    nz = 3
    xmin = 4.5
    xmax = 5.5
    ymin = 4.5
    ymax = 5.5
    zmin = 0.06
    zmax = 1.06
    boundary_name_prefix = 'ball'
  []
  [block_two]
    type = GeneratedMeshGenerator
    dim = 3
    nx = 2
    ny = 2
    nz = 2
    xmin = 0.0
    xmax = 10
    ymin = 0.0
    ymax = 10
    zmin = -2
    zmax = 0
    boundary_name_prefix = 'base'
    boundary_id_offset = 10
  []
  [block_one_id]
    type = SubdomainIDGenerator
    input = block_one
    subdomain_id = 1
  []
  [block_two_id]
    type = SubdomainIDGenerator
    input = block_two
    subdomain_id = 2
  []
  [combine]
    type = MeshCollectionGenerator
    inputs = ' block_one_id block_two_id'
  []
[]

[AuxVariables]
  [penetration]
  []
[]

[AuxKernels]
  [penetration]
    type = PenetrationAux
    variable = penetration
    boundary = ball_back
    paired_boundary = base_front
    quantity = distance
  []
[]

[Variables]
  [disp_x]
  []
  [disp_y]
  []
  [disp_z]
  []
[]

[AuxVariables]
  [gap_rate]
  []
  [vel_x]
  []
  [accel_x]
  []
  [vel_y]
  []
  [accel_y]
  []
  [vel_z]
  []
  [accel_z]
  []
  [stress_zz]
    family = MONOMIAL
    order = CONSTANT
  []
  [strain_zz]
    family = MONOMIAL
    order = CONSTANT
  []
  [kinetic_energy_one]
    order = CONSTANT
    family = MONOMIAL
  []
  [elastic_energy_one]
    order = CONSTANT
    family = MONOMIAL
  []
  [kinetic_energy_two]
    order = CONSTANT
    family = MONOMIAL
  []
  [elastic_energy_two]
    order = CONSTANT
    family = MONOMIAL
  []
[]

[AuxKernels]
  [stress_zz]
    type = RankTwoAux
    rank_two_tensor = stress
    index_i = 2
    index_j = 2
    variable = stress_zz
    execute_on = 'TIMESTEP_END'
  []
  [strain_zz]
    type = RankTwoAux
    rank_two_tensor = mechanical_strain
    index_i = 2
    index_j = 2
    variable = strain_zz
  []
  [accel_x]
    type = TestNewmarkTI
    variable = accel_x
    displacement = disp_x
    first = false
    execute_on = 'LINEAR TIMESTEP_BEGIN TIMESTEP_END'
  []
  [vel_x]
    type = TestNewmarkTI
    variable = vel_x
    displacement = disp_x
    execute_on = 'LINEAR TIMESTEP_BEGIN TIMESTEP_END'
  []
  [accel_y]
    type = TestNewmarkTI
    variable = accel_y
    displacement = disp_y
    first = false
    execute_on = 'LINEAR TIMESTEP_BEGIN TIMESTEP_END'
  []
  [vel_y]
    type = TestNewmarkTI
    variable = vel_y
    displacement = disp_x
    execute_on = 'LINEAR TIMESTEP_BEGIN TIMESTEP_END'
  []
  [accel_z]
    type = TestNewmarkTI
    variable = accel_z
    displacement = disp_z
    first = false
    execute_on = 'LINEAR TIMESTEP_BEGIN TIMESTEP_END'
  []
  [vel_z]
    type = TestNewmarkTI
    variable = vel_z
    displacement = disp_z
    execute_on = 'LINEAR TIMESTEP_BEGIN TIMESTEP_END'
  []
  [kinetic_energy_one]
    type = KineticEnergyAux
    block = '1'
    variable = kinetic_energy_one
    newmark_velocity_x = vel_x
    newmark_velocity_y = vel_y
    newmark_velocity_z = vel_z
    density = density
  []
  [elastic_energy_one]
    type = ElasticEnergyAux
    variable = elastic_energy_one
    block = '1'
  []
  [kinetic_energy_two]
    type = KineticEnergyAux
    block = '2'
    variable = kinetic_energy_two
    newmark_velocity_x = vel_x
    newmark_velocity_y = vel_y
    newmark_velocity_z = vel_z
    density = density
  []
  [elastic_energy_two]
    type = ElasticEnergyAux
    variable = elastic_energy_two
    block = '2'
  []
[]

[Kernels]
  [DynamicTensorMechanics]
    displacements = 'disp_x disp_y disp_z'
    stiffness_damping_coefficient = 1.0e-3
    generate_output = 'stress_zz strain_zz'
  []
  [Mass_x]
    type = MassMatrix
    variable = disp_x
    density = density
    matrix_tags = 'mass'
  []
  [Mass_y]
    type = MassMatrix
    variable = disp_y
    density = density
    matrix_tags = 'mass'
  []
  [Mass_z]
    type = MassMatrix
    variable = disp_z
    density = density
    matrix_tags = 'mass'
  []
  [gravity]
    type = Gravity
    variable = disp_z
    value = -981
    block = 1
  []
[]

[BCs]
  [x_front]
    type = ExplicitDirichletBC
    variable = disp_x
    boundary = 'ball_front'
    value = 0.0
  []
  [y_front]
    type = ExplicitDirichletBC
    variable = disp_y
    boundary = 'ball_front'
    value = 0.0
  []
  [x_fixed]
    type = ExplicitDirichletBC
    variable = disp_x
    boundary = 'base_back'
    value = 0.0
  []
  [y_fixed]
    type = ExplicitDirichletBC
    variable = disp_y
    boundary = 'base_back'
    value = 0.0
  []
  [z_fixed]
    type = ExplicitDirichletBC
    variable = disp_z
    boundary = 'base_back'
    value = 0.0
  []
  [z_fixed_front]
    type = ExplicitDirichletBC
    variable = disp_z
    boundary = 'base_front'
    value = 0.0
  []
[]

[ExplicitDynamicsContact]
  [my_contact]
    model = frictionless_balance
    primary = base_front
    secondary = ball_back
    vel_x = 'vel_x'
    vel_y = 'vel_y'
    vel_z = 'vel_z'
    verbose = true
  []
[]

[Materials]
  [elasticity_tensor_block_one]
    type = ComputeIsotropicElasticityTensor
    youngs_modulus = 1e6
    poissons_ratio = 0.0
    block = 1
    outputs = 'exodus'
    output_properties = __all__
  []
  [elasticity_tensor_block_two]
    type = ComputeIsotropicElasticityTensor
    youngs_modulus = 1e10
    poissons_ratio = 0.0
    block = 2
    outputs = 'exodus'
    output_properties = __all__
  []
  [strain_block]
    type = ComputeFiniteStrain # ComputeIncrementalSmallStrain
    displacements = 'disp_x disp_y disp_z'
    implicit = false
  []
  [stress_block]
    type = ComputeFiniteStrainElasticStress
  []
  [density_one]
    type = GenericConstantMaterial
    prop_names = density
    prop_values = 1e1
    outputs = 'exodus'
    output_properties = 'density'
    block = '1'
  []
  [density_two]
    type = GenericConstantMaterial
    prop_names = density
    prop_values = 1e6
    outputs = 'exodus'
    output_properties = 'density'
    block = '2'
  []
  [wave_speed]
    type = WaveSpeed
    outputs = 'exodus'
    output_properties = 'wave_speed'
  []
[]

[Executioner]
  type = Transient

  end_time = 0.03
  dt = 2e-4
  timestep_tolerance = 1e-6

  [TimeIntegrator]
    type = ExplicitMixedOrder
    mass_matrix_tag = mass
    second_order_vars = 'disp_x disp_y disp_z'
  []
[]

[Outputs]
  time_step_interval = 2
  exodus = true
  csv = true
  execute_on = 'TIMESTEP_END'
  file_base = highvel_out
[]

[Postprocessors]
  [accel_58z]
    type = NodalVariableValue
    nodeid = 1
    variable = accel_z
  []
  [vel_58z]
    type = NodalVariableValue
    nodeid = 1
    variable = vel_z
  []
  [critical_time_step]
    type = CriticalTimeStep
  []
  [contact_pressure_max]
    type = NodalExtremeValue
    variable = contact_pressure
    block = '1 2'
    value_type = max
  []
  [penetration_max]
    type = NodalExtremeValue
    variable = penetration
    block = '1 2'
    value_type = max
  []
  [total_kinetic_energy_one]
    type = ElementIntegralVariablePostprocessor
    variable = kinetic_energy_one
    block = '1'
  []
  [total_elastic_energy_one]
    type = ElementIntegralVariablePostprocessor
    variable = elastic_energy_one
    block = '1'
  []
  [total_kinetic_energy_two]
    type = ElementIntegralVariablePostprocessor
    variable = kinetic_energy_two
    block = '2'
  []
  [total_elastic_energy_two]
    type = ElementIntegralVariablePostprocessor
    variable = elastic_energy_two
    block = '2'
  []
[]

(test/tests/geomsearch/penetration_locator/disconnected_penetration.i)

###########################################################
# This is a test of the Geometric Search System. This test
# uses the penetration location object through the
# PenetrationAux Auxilary Kernel to detect overlaps of
# specified interfaces (boundaries) in the domain.
#
# @Requirement F6.50
###########################################################


[Mesh]
  [connected_mesh]
    type = FileMeshGenerator
    file = 2dcontact_collide.e
  []
  [exploded_mesh]
    type = BreakMeshByElementGenerator
    input = connected_mesh
    subdomains = '1 2'
    interface_name = 'interelement'
  []
[]

[Variables]
  active = 'u'

  [./u]
    order = FIRST
    family = LAGRANGE
  [../]
[]

[AuxVariables]
  [./penetration]
    order = FIRST
    family = LAGRANGE
  [../]
[]

[Kernels]
  [./l2]
    type = MassMatrix
    variable = u
    matrix_tags = 'system'
  [../]
[]

[AuxKernels]
  [./penetrate]
    type = PenetrationAux
    variable = penetration
    boundary = 2
    paired_boundary = 3
    search_method = all_proximate_sides
  [../]
[]

[Executioner]
  type = Steady

  solve_type = 'PJFNK'
[]

[Outputs]
  exodus = true
[]

(modules/contact/test/tests/explicit_dynamics/test_balance.i)

# One element test to test the central difference time integrator in 3D.
[GlobalParams]
  displacements = 'disp_x disp_y disp_z'
  volumetric_locking_correction = true
[]

[Problem]
  extra_tag_matrices = 'mass'
[]

[Mesh]
  [block_one]
    type = GeneratedMeshGenerator
    dim = 3
    nx = 3
    ny = 3
    nz = 3
    xmin = 4.5
    xmax = 5.5
    ymin = 4.5
    ymax = 5.5
    zmin = 0.0001
    zmax = 1.0001
    boundary_name_prefix = 'ball'
  []
  [block_two]
    type = GeneratedMeshGenerator
    dim = 3
    nx = 2
    ny = 2
    nz = 2
    xmin = 0.0
    xmax = 10
    ymin = 0.0
    ymax = 10
    zmin = -2
    zmax = 0
    boundary_name_prefix = 'base'
    boundary_id_offset = 10
  []
  [block_one_id]
    type = SubdomainIDGenerator
    input = block_one
    subdomain_id = 1
  []
  [block_two_id]
    type = SubdomainIDGenerator
    input = block_two
    subdomain_id = 2
  []
  [combine]
    type = MeshCollectionGenerator
    inputs = ' block_one_id block_two_id'
  []
  allow_renumbering = false
[]

[Variables]
  [disp_x]
  []
  [disp_y]
  []
  [disp_z]
  []
[]

[AuxVariables]
  [gap_rate]
  []
  [vel_x]
  []
  [accel_x]
  []
  [vel_y]
  []
  [accel_y]
  []
  [vel_z]
  []
  [accel_z]
  []
  [stress_zz]
    family = MONOMIAL
    order = CONSTANT
  []
  [strain_zz]
    family = MONOMIAL
    order = CONSTANT
  []
[]

[AuxKernels]
  [stress_zz]
    type = RankTwoAux
    rank_two_tensor = stress
    index_i = 2
    index_j = 2
    variable = stress_zz
    execute_on = 'TIMESTEP_END'
  []
  [strain_zz]
    type = RankTwoAux
    rank_two_tensor = mechanical_strain
    index_i = 2
    index_j = 2
    variable = strain_zz
  []
  [accel_x]
    type = TestNewmarkTI
    variable = accel_x
    displacement = disp_x
    first = false
    execute_on = 'LINEAR TIMESTEP_BEGIN TIMESTEP_END'
  []
  [vel_x]
    type = TestNewmarkTI
    variable = vel_x
    displacement = disp_x
    execute_on = 'LINEAR TIMESTEP_BEGIN TIMESTEP_END'
  []
  [accel_y]
    type = TestNewmarkTI
    variable = accel_y
    displacement = disp_y
    first = false
    execute_on = 'LINEAR TIMESTEP_BEGIN TIMESTEP_END'
  []
  [vel_y]
    type = TestNewmarkTI
    variable = vel_y
    displacement = disp_x
    execute_on = 'LINEAR TIMESTEP_BEGIN TIMESTEP_END'
  []
  [accel_z]
    type = TestNewmarkTI
    variable = accel_z
    displacement = disp_z
    first = false
    execute_on = 'LINEAR TIMESTEP_BEGIN TIMESTEP_END'
  []
  [vel_z]
    type = TestNewmarkTI
    variable = vel_z
    displacement = disp_z
    execute_on = 'LINEAR TIMESTEP_BEGIN TIMESTEP_END'
  []
[]

[AuxVariables]
  [penetration]
  []
[]

[AuxKernels]
  [penetration]
    type = PenetrationAux
    variable = penetration
    boundary = ball_back
    paired_boundary = base_front
    quantity = distance
  []
[]

[Kernels]
  [DynamicTensorMechanics]
    displacements = 'disp_x disp_y disp_z'
    volumetric_locking_correction = true
    stiffness_damping_coefficient = 0.001
    generate_output = 'stress_zz strain_zz'
  []
  [Mass_x]
    type = MassMatrix
    variable = disp_x
    density = density
    matrix_tags = 'mass'
  []
  [Mass_y]
    type = MassMatrix
    variable = disp_y
    density = density
    matrix_tags = 'mass'
  []
  [Mass_z]
    type = MassMatrix
    variable = disp_z
    density = density
    matrix_tags = 'mass'
  []
[]

[Kernels]
  [gravity]
    type = Gravity
    variable = disp_z
    value = -981.0
  []
[]

[BCs]
  [x_front]
    type = ExplicitDirichletBC
    variable = disp_x
    boundary = 'ball_front'
    value = 0.0
  []
  [y_front]
    type = ExplicitDirichletBC
    variable = disp_y
    boundary = 'ball_front'
    value = 0.0
  []
  [x_fixed]
    type = ExplicitDirichletBC
    variable = disp_x
    boundary = 'base_back'
    value = 0.0
  []
  [y_fixed]
    type = ExplicitDirichletBC
    variable = disp_y
    boundary = 'base_back'
    value = 0.0
  []
  [z_fixed]
    type = ExplicitDirichletBC
    variable = disp_z
    boundary = 'base_back'
    value = 0.0
  []
  [z_fixed_front]
    type = ExplicitDirichletBC
    variable = disp_z
    boundary = 'base_front'
    value = 0.0
  []
[]

[ExplicitDynamicsContact]
  [my_contact]
    model = frictionless_balance
    primary = base_front
    secondary = ball_back
    vel_x = 'vel_x'
    vel_y = 'vel_y'
    vel_z = 'vel_z'
    verbose = true
  []
[]

[Materials]
  [elasticity_tensor_block_one]
    type = ComputeIsotropicElasticityTensor
    youngs_modulus = 1e6
    poissons_ratio = 0.0
    block = 1
    outputs = 'exodus'
    output_properties = __all__
  []
  [elasticity_tensor_block_two]
    type = ComputeIsotropicElasticityTensor
    youngs_modulus = 1e10
    poissons_ratio = 0.0
    block = 2
    outputs = 'exodus'
    output_properties = __all__
  []
  [strain_block]
    type = ComputeFiniteStrain
    displacements = 'disp_x disp_y disp_z'
    implicit = false
  []
  [stress_block]
    type = ComputeFiniteStrainElasticStress
  []
  [density_one]
    type = GenericConstantMaterial
    prop_names = density
    prop_values = 1e1
    outputs = 'exodus'
    output_properties = 'density'
    block = '1'
  []
  [density_two]
    type = GenericConstantMaterial
    prop_names = density
    prop_values = 1e6
    outputs = 'exodus'
    output_properties = 'density'
    block = '2'
  []
  [wave_speed]
    type = WaveSpeed
    outputs = 'exodus'
    output_properties = 'wave_speed'
  []
[]

[Executioner]
  type = Transient
  start_time = 0
  end_time = 0.0025
  dt = 0.00001
  timestep_tolerance = 1e-6

  [TimeIntegrator]
    type = ExplicitMixedOrder
    mass_matrix_tag = 'mass'
    second_order_vars = 'disp_x disp_y disp_z'
  []
[]

[Outputs]
  time_step_interval = 10
  exodus = true
  csv = true
[]

[Postprocessors]
  [accel_58z]
    type = NodalVariableValue
    nodeid = 1
    variable = accel_z
  []
  [vel_58z]
    type = NodalVariableValue
    nodeid = 1
    variable = vel_z
  []
  [disp_58z]
    type = NodalVariableValue
    nodeid = 1
    variable = disp_z
  []
  [critical_time_step]
    type = CriticalTimeStep
  []
  [contact_pressure_max]
    type = NodalExtremeValue
    variable = contact_pressure
    block = '1 2'
    value_type = max
  []
[]

Input Parameters
Input Files