ParsedPostprocessor

Computes a parsed expression with post-processors

The "function" to parse may only include other postprocessors, the time variable and constants from the input file.

note

Derivatives and integrals are not natively supported by the parsing operation, unless the postprocessors in the parsed expression are already the derivatives / integrals of interest.

Example Input File Syntax

In this test input file, we compute various quantities using ParsedPostprocessors. First, we look at the ratio of two postprocessors:

[parsed]
  type = ParsedPostprocessor
  function = 'L2_norm / L1_norm'
  pp_names = 'L2_norm L1_norm'
[]

Then we add a time dependence to the postprocessor result:

[parsed_with_t]
  type = ParsedPostprocessor
  function = 'L2_norm + L1_norm + t'
  pp_names = 'L2_norm L1_norm'
  use_t = true
[]

And finally we introduce constants from the input file in the parsed expression:

[parsed_with_constants]
  type = ParsedPostprocessor
  function = 'L2_norm + 3*L1_norm + mu'
  pp_names = 'L2_norm L1_norm'
  constant_names = 'mu'
  constant_expressions = '4'
[]

Input Parameters

functionfunction expression
C++ Type:FunctionExpression
Controllable:No
Description:function expression
pp_namesPost-processors arguments
C++ Type:std::vector<PostprocessorName>
Controllable:No
Description:Post-processors arguments

Required Parameters

constant_expressionsVector of values for the constants in constant_names (can be an FParser expression)
C++ Type:std::vector<std::string>
Controllable:No
Description:Vector of values for the constants in constant_names (can be an FParser expression)
constant_namesVector of constants used in the parsed function (use this for kB etc.)
C++ Type:std::vector<std::string>
Controllable:No
Description:Vector of constants used in the parsed function (use this for kB etc.)
execute_onTIMESTEP_ENDThe list of flag(s) indicating when this object should be executed, the available options include NONE, INITIAL, LINEAR, NONLINEAR, TIMESTEP_END, TIMESTEP_BEGIN, FINAL, CUSTOM, ALWAYS.
Default:TIMESTEP_END
C++ Type:ExecFlagEnum
Options:NONE, INITIAL, LINEAR, NONLINEAR, TIMESTEP_END, TIMESTEP_BEGIN, FINAL, CUSTOM, TRANSFER, ALWAYS
Controllable:No
Description:The list of flag(s) indicating when this object should be executed, the available options include NONE, INITIAL, LINEAR, NONLINEAR, TIMESTEP_END, TIMESTEP_BEGIN, FINAL, CUSTOM, ALWAYS.
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
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_tFalseMake time (t) variables available in the function expression.
Default:False
C++ Type:bool
Controllable:No
Description:Make time (t) variables available in the function expression.

Optional Parameters

allow_duplicate_execution_on_initialFalseIn the case where this UserObject is depended upon by an initial condition, allow it to be executed twice during the initial setup (once before the IC and again after mesh adaptivity (if applicable).
Default:False
C++ Type:bool
Controllable:No
Description:In the case where this UserObject is depended upon by an initial condition, allow it to be executed twice during the initial setup (once before the IC and again after mesh adaptivity (if applicable).
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.
disable_fpoptimizerFalseDisable the function parser algebraic optimizer
Default:False
C++ Type:bool
Controllable:No
Description:Disable the function parser algebraic optimizer
enableTrueSet the enabled status of the MooseObject.
Default:True
C++ Type:bool
Controllable:Yes
Description:Set the enabled status of the MooseObject.
enable_ad_cacheTrueEnable cacheing of function derivatives for faster startup time
Default:True
C++ Type:bool
Controllable:No
Description:Enable cacheing of function derivatives for faster startup time
enable_auto_optimizeTrueEnable automatic immediate optimization of derivatives
Default:True
C++ Type:bool
Controllable:No
Description:Enable automatic immediate optimization of derivatives
enable_jitTrueEnable just-in-time compilation of function expressions for faster evaluation
Default:True
C++ Type:bool
Controllable:No
Description:Enable just-in-time compilation of function expressions for faster evaluation
evalerror_behaviornanWhat to do if evaluation error occurs. Options are to pass a nan, pass a nan with a warning, throw a error, or throw an exception
Default:nan
C++ Type:MooseEnum
Options:nan, nan_warning, error, exception
Controllable:No
Description:What to do if evaluation error occurs. Options are to pass a nan, pass a nan with a warning, throw a error, or throw an exception
force_postauxFalseForces the UserObject to be executed in POSTAUX
Default:False
C++ Type:bool
Controllable:No
Description:Forces the UserObject to be executed in POSTAUX
force_preauxFalseForces the UserObject to be executed in PREAUX
Default:False
C++ Type:bool
Controllable:No
Description:Forces the UserObject to be executed in PREAUX
force_preicFalseForces the UserObject to be executed in PREIC during initial setup
Default:False
C++ Type:bool
Controllable:No
Description:Forces the UserObject to be executed in PREIC during initial setup
outputsVector of output names were you would like to restrict the output of variables(s) associated with this object
C++ Type:std::vector<OutputName>
Controllable:No
Description:Vector of output names were you would like to restrict the output of variables(s) associated with this object
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

Input Files

(modules/heat_conduction/test/tests/gap_heat_transfer_mortar/gap_heat_transfer_sphere3D_mortar.i)
(modules/heat_conduction/test/tests/gap_heat_transfer_mortar/large_gap_heat_transfer_test_cylinder.i)
(modules/heat_conduction/test/tests/gap_heat_transfer_balance/large_gap_heat_transfer_test_cylinder.i)
(modules/heat_conduction/test/tests/gap_heat_transfer_mortar/gap_heat_transfer_sphere3D.i)
(test/tests/postprocessors/parsed_postprocessor/parsed_pp.i)
(modules/heat_conduction/test/tests/gap_heat_transfer_mortar/large_gap_heat_transfer_test_rz_cylinder.i)
(modules/thermal_hydraulics/test/tests/components/simple_turbine_1phase/phy.conservation.i)
(modules/heat_conduction/test/tests/gap_heat_transfer_mortar/gap_heat_transfer_3D_mortar.i)
(modules/heat_conduction/test/tests/gap_heat_transfer_balance/large_gap_heat_transfer_test_rz_cylinder.i)
(modules/thermal_hydraulics/test/tests/problems/brayton_cycle/closed_brayton_cycle.i)
(modules/heat_conduction/test/tests/gap_heat_transfer_mortar/large_gap_heat_transfer_test_sphere_mortar_error.i)
(modules/tensor_mechanics/test/tests/rom_stress_update/ad_verification.i)
(modules/heat_conduction/test/tests/gap_heat_transfer_mortar/large_gap_heat_transfer_test_cylinder_mortar_error.i)
(modules/heat_conduction/test/tests/gap_heat_transfer_mortar/gap_heat_transfer_3D.i)
(modules/tensor_mechanics/test/tests/rom_stress_update/nonad_verification.i)
(modules/heat_conduction/test/tests/gap_heat_transfer_balance/large_gap_heat_transfer_test_sphere.i)
(modules/heat_conduction/test/tests/gap_heat_transfer_mortar/large_gap_heat_transfer_test_sphere.i)
(modules/heat_conduction/test/tests/gap_heat_transfer_mortar/large_gap_heat_transfer_test_cylinder_mortar.i)
(modules/heat_conduction/test/tests/gap_heat_transfer_mortar/large_gap_heat_transfer_test_rz_cylinder_mortar.i)
(modules/heat_conduction/test/tests/gap_heat_transfer_mortar/large_gap_heat_transfer_test_sphere_mortar.i)
(modules/thermal_hydraulics/test/tests/problems/brayton_cycle/open_brayton_cycle.i)

(test/tests/postprocessors/parsed_postprocessor/parsed_pp.i)

[Mesh]
  type = GeneratedMesh
  dim = 2
  nx = 10
  ny = 10
[]

[Variables]
  [./u]
  [../]
[]

[Kernels]
  [./diff]
    type = Diffusion
    variable = u
  [../]
[]

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

[Postprocessors]
  [./L2_norm]
    type = ElementL2Norm
    variable = u
  [../]
  [L1_norm]
    type = ElementL1Error
    function = 0
    variable = u
  []
  [parsed]
    type = ParsedPostprocessor
    function = 'L2_norm / L1_norm'
    pp_names = 'L2_norm L1_norm'
  []
  [parsed_with_t]
    type = ParsedPostprocessor
    function = 'L2_norm + L1_norm + t'
    pp_names = 'L2_norm L1_norm'
    use_t = true
  []
  [parsed_with_constants]
    type = ParsedPostprocessor
    function = 'L2_norm + 3*L1_norm + mu'
    pp_names = 'L2_norm L1_norm'
    constant_names = 'mu'
    constant_expressions = '4'
  []
[]

[Executioner]
  type = Transient
  num_steps = 4
  nl_abs_tol = 1e-8
  solve_type = PJFNK
  petsc_options_iname = '-pc_type -pc_hypre_type'
  petsc_options_value = 'hypre boomeramg'
  line_search = 'none'
[]

[Outputs]
  csv = true
[]

(test/tests/postprocessors/parsed_postprocessor/parsed_pp.i)

[Mesh]
  type = GeneratedMesh
  dim = 2
  nx = 10
  ny = 10
[]

[Variables]
  [./u]
  [../]
[]

[Kernels]
  [./diff]
    type = Diffusion
    variable = u
  [../]
[]

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

[Postprocessors]
  [./L2_norm]
    type = ElementL2Norm
    variable = u
  [../]
  [L1_norm]
    type = ElementL1Error
    function = 0
    variable = u
  []
  [parsed]
    type = ParsedPostprocessor
    function = 'L2_norm / L1_norm'
    pp_names = 'L2_norm L1_norm'
  []
  [parsed_with_t]
    type = ParsedPostprocessor
    function = 'L2_norm + L1_norm + t'
    pp_names = 'L2_norm L1_norm'
    use_t = true
  []
  [parsed_with_constants]
    type = ParsedPostprocessor
    function = 'L2_norm + 3*L1_norm + mu'
    pp_names = 'L2_norm L1_norm'
    constant_names = 'mu'
    constant_expressions = '4'
  []
[]

[Executioner]
  type = Transient
  num_steps = 4
  nl_abs_tol = 1e-8
  solve_type = PJFNK
  petsc_options_iname = '-pc_type -pc_hypre_type'
  petsc_options_value = 'hypre boomeramg'
  line_search = 'none'
[]

[Outputs]
  csv = true
[]

(test/tests/postprocessors/parsed_postprocessor/parsed_pp.i)

[Mesh]
  type = GeneratedMesh
  dim = 2
  nx = 10
  ny = 10
[]

[Variables]
  [./u]
  [../]
[]

[Kernels]
  [./diff]
    type = Diffusion
    variable = u
  [../]
[]

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

[Postprocessors]
  [./L2_norm]
    type = ElementL2Norm
    variable = u
  [../]
  [L1_norm]
    type = ElementL1Error
    function = 0
    variable = u
  []
  [parsed]
    type = ParsedPostprocessor
    function = 'L2_norm / L1_norm'
    pp_names = 'L2_norm L1_norm'
  []
  [parsed_with_t]
    type = ParsedPostprocessor
    function = 'L2_norm + L1_norm + t'
    pp_names = 'L2_norm L1_norm'
    use_t = true
  []
  [parsed_with_constants]
    type = ParsedPostprocessor
    function = 'L2_norm + 3*L1_norm + mu'
    pp_names = 'L2_norm L1_norm'
    constant_names = 'mu'
    constant_expressions = '4'
  []
[]

[Executioner]
  type = Transient
  num_steps = 4
  nl_abs_tol = 1e-8
  solve_type = PJFNK
  petsc_options_iname = '-pc_type -pc_hypre_type'
  petsc_options_value = 'hypre boomeramg'
  line_search = 'none'
[]

[Outputs]
  csv = true
[]

(modules/heat_conduction/test/tests/gap_heat_transfer_mortar/gap_heat_transfer_sphere3D_mortar.i)

sphere_outer_htc = 10 # W/m^2/K
sphere_outer_Tinf = 300 # K

[GlobalParams]
  order = SECOND
  family = LAGRANGE
[]

[Problem]
  kernel_coverage_check = false
  material_coverage_check = false
[]

[Mesh]
  [file]
    type = FileMeshGenerator
    file = sphere3D.e
  []
  [secondary]
    type = LowerDBlockFromSidesetGenerator
    sidesets = '2'
    new_block_id = 10001
    new_block_name = 'secondary_lower'
    input = file
  []
  [primary]
    type = LowerDBlockFromSidesetGenerator
    sidesets = '3'
    new_block_id = 10000
    new_block_name = 'primary_lower'
    input = secondary
  []
[]

[Functions]
  [temp]
    type = PiecewiseLinear
    x = '0   1'
    y = '100 200'
  []
[]

[Variables]
  [temp]
    initial_condition = 500
  []
  [lm]
    order = FIRST
    family = LAGRANGE
    block = 'secondary_lower'
  []
[]

[AuxVariables]
  # [gap_conductance]
  #   order = CONSTANT
  #   family = MONOMIAL
  # []
  [power_density]
    block = 'fuel'
    initial_condition = 50e3
  []
[]

[Kernels]
  [heat_conduction]
    type = HeatConduction
    variable = temp
    block = '1 2'
  []
  [heat_source]
    type = CoupledForce
    variable = temp
    block = 'fuel'
    v = power_density
  []
[]

# [AuxKernels]
#   [gap_cond]
#     type = MaterialRealAux
#     property = gap_conductance
#     variable = gap_conductance
#     boundary = 2
#   []
# []

[Materials]
  [heat1]
    type = HeatConductionMaterial
    block = '1 2'
    specific_heat = 1.0
    thermal_conductivity = 34.6
  []
[]

[UserObjects]
  [radiation]
    type = GapFluxModelRadiation
    T = temp
    boundary = 2
    primary_emissivity = 0.0
    secondary_emissivity = 0.0
  []
  [conduction]
    type = GapFluxModelConduction
    T = temp
    boundary = 2
    gap_conductivity = 5.0
  []
[]

[Constraints]
  [ced]
    type = ModularGapConductanceConstraint
    variable = lm
    secondary_variable = temp
    primary_boundary = 3
    primary_subdomain = 10000
    secondary_boundary = 2
    secondary_subdomain = 10001
    gap_flux_models = 'radiation conduction'
    gap_geometry_type = SPHERE
    sphere_origin = '0 0 0'
  []
[]

[BCs]
  [RPV_out_BC] # k \nabla T = h (T- T_inf) at RPV outer boundary
    type = ConvectiveFluxFunction # (Robin BC)
    variable = temp
    boundary = '4' # outer RPV
    coefficient = ${sphere_outer_htc}
    T_infinity = ${sphere_outer_Tinf}
  []
[]

[Executioner]
  type = Transient
  solve_type = 'PJFNK'

  petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
  petsc_options_value = 'lu       superlu_dist'

  dt = 1
  dtmin = 0.01
  end_time = 1

  nl_rel_tol = 1e-12
  nl_abs_tol = 1e-7
[]

[Outputs]
  exodus = true
  csv = true
  [Console]
    type = Console
  []
[]

[Postprocessors]
  [temp_left]
    type = SideAverageValue
    boundary = 2
    variable = temp
  []

  [temp_right]
    type = SideAverageValue
    boundary = 3
    variable = temp
  []

  [flux_left]
    type = SideDiffusiveFluxIntegral
    variable = temp
    boundary = 2
    diffusivity = thermal_conductivity
  []
  [flux_right]
    type = SideDiffusiveFluxIntegral
    variable = temp
    boundary = 3
    diffusivity = thermal_conductivity
  []
  [ptot]
    type = ElementIntegralVariablePostprocessor
    variable = power_density
    block = 'fuel'
  []
  [sphere_convective_out]
    type = ConvectiveHeatTransferSideIntegral
    T_solid = temp
    boundary = '4' # outer RVP
    T_fluid = ${sphere_outer_Tinf}
    htc = ${sphere_outer_htc}
  []
  [heat_balance] # should be equal to 0 upon convergence
    type = ParsedPostprocessor
    function = '(sphere_convective_out - ptot) / ptot'
    pp_names = 'sphere_convective_out ptot'
  []
[]

[VectorPostprocessors]
  [NodalTemperature]
    type = NodalValueSampler
    sort_by = id
    boundary = '2 3'
    variable = temp
  []
[]

(modules/heat_conduction/test/tests/gap_heat_transfer_mortar/large_gap_heat_transfer_test_cylinder.i)

rpv_core_gap_size = 0.15

core_outer_radius = 2
rpv_inner_radius = '${fparse 2 + rpv_core_gap_size}'
rpv_outer_radius = '${fparse 2.5 + rpv_core_gap_size}'

rpv_outer_htc = 10 # W/m^2/K
rpv_outer_Tinf = 300 # K

core_blocks = '1'
rpv_blocks = '3'

[Mesh]
  [core_gap_rpv]
    type = ConcentricCircleMeshGenerator
    num_sectors = 10
    radii = '${core_outer_radius} ${rpv_inner_radius} ${rpv_outer_radius}'
    rings = '2 1 2'
    has_outer_square = false
    preserve_volumes = true
    portion = full
  []
  [rename_core_bdy]
    type = SideSetsBetweenSubdomainsGenerator
    input = core_gap_rpv
    primary_block = 1
    paired_block = 2
    new_boundary = 'core_outer'
  []
  [rename_inner_rpv_bdy]
    type = SideSetsBetweenSubdomainsGenerator
    input = rename_core_bdy
    primary_block = 3
    paired_block = 2
    new_boundary = 'rpv_inner'
  []
  [2d_mesh]
    type = BlockDeletionGenerator
    input = rename_inner_rpv_bdy
    block = 2
  []
  allow_renumbering = false
[]

[Variables]
  [Tsolid]
    initial_condition = 500
  []
[]

[Kernels]
  [heat_source]
    type = CoupledForce
    variable = Tsolid
    block = '${core_blocks}'
    v = power_density
  []
  [heat_conduction]
    type = HeatConduction
    variable = Tsolid
  []
[]

[BCs]
  [RPV_out_BC] # k \nabla T = h (T- T_inf) at RPV outer boundary
    type = ConvectiveFluxFunction # (Robin BC)
    variable = Tsolid
    boundary = 'outer' # outer RPV
    coefficient = ${rpv_outer_htc}
    T_infinity = ${rpv_outer_Tinf}
  []
[]

[ThermalContact]
  [RPV_gap]
    type = GapHeatTransfer
    gap_geometry_type = 'CYLINDER'
    emissivity_primary = 0.8
    emissivity_secondary = 0.8
    variable = Tsolid
    primary = 'core_outer'
    secondary = 'rpv_inner'
    gap_conductivity = 0.1
    quadrature = true
    cylinder_axis_point_1 = '0 0 0'
    cylinder_axis_point_2 = '0 0 5'
  []
[]

[AuxVariables]
  [power_density]
    block = '${core_blocks}'
    initial_condition = 50e3
  []
[]

[Materials]
  [simple_mat]
    type = HeatConductionMaterial
    thermal_conductivity = 34.6 # W/m/K
  []
[]

[Postprocessors]
  [Tcore_avg]
    type = ElementAverageValue
    variable = Tsolid
    block = '${core_blocks}'
  []
  [Tcore_max]
    type = ElementExtremeValue
    value_type = max
    variable = Tsolid
    block = '${core_blocks}'
  []
  [Tcore_min]
    type = ElementExtremeValue
    value_type = min
    variable = Tsolid
    block = '${core_blocks}'
  []
  [Trpv_avg]
    type = ElementAverageValue
    variable = Tsolid
    block = '${rpv_blocks}'
  []
  [Trpv_max]
    type = ElementExtremeValue
    value_type = max
    variable = Tsolid
    block = '${rpv_blocks}'
  []
  [Trpv_min]
    type = ElementExtremeValue
    value_type = min
    variable = Tsolid
    block = '${rpv_blocks}'
  []
  [ptot]
    type = ElementIntegralVariablePostprocessor
    variable = power_density
    block = '${core_blocks}'
  []
  [rpv_convective_out]
    type = ConvectiveHeatTransferSideIntegral
    T_solid = Tsolid
    boundary = 'outer' # outer RVP
    T_fluid = ${rpv_outer_Tinf}
    htc = ${rpv_outer_htc}
  []
  [heat_balance] # should be equal to 0 upon convergence
    type = ParsedPostprocessor
    function = '(rpv_convective_out - ptot) / ptot'
    pp_names = 'rpv_convective_out ptot'
  []
[]

[VectorPostprocessors]
  [NodalTemperature]
    type = NodalValueSampler
    sort_by = id
    boundary = 'rpv_inner core_outer'
    variable = Tsolid
  []
[]

[Executioner]
  type = Steady
  automatic_scaling = true
  compute_scaling_once = false
  petsc_options_iname = '-pc_type -pc_hypre_type -ksp_gmres_restart '
  petsc_options_value = 'hypre boomeramg 100'

  nl_rel_tol = 1e-10
  nl_abs_tol = 1e-10
  l_max_its = 100

  [Quadrature]
    side_order = seventh
  []
  line_search = none
[]

[Outputs]
  exodus = false
  csv = true
[]

(modules/heat_conduction/test/tests/gap_heat_transfer_balance/large_gap_heat_transfer_test_cylinder.i)

rpv_core_gap_size = 0.15

core_outer_radius = 2
rpv_inner_radius = ${fparse 2 + rpv_core_gap_size}
rpv_outer_radius = ${fparse 2.5 + rpv_core_gap_size}

rpv_outer_htc = 10 # W/m^2/K
rpv_outer_Tinf = 300 # K

core_blocks = '1'
rpv_blocks = '3'

[Mesh]
  [core_gap_rpv]
    type = ConcentricCircleMeshGenerator
    num_sectors = 10
    radii = '${core_outer_radius} ${rpv_inner_radius} ${rpv_outer_radius}'
    rings = '2 1 2'
    has_outer_square = false
    preserve_volumes = true
    portion = full
  []
  [rename_core_bdy]
    type = SideSetsBetweenSubdomainsGenerator
    input = core_gap_rpv
    primary_block = 1
    paired_block = 2
    new_boundary = 'core_outer'
  []
  [rename_inner_rpv_bdy]
    type = SideSetsBetweenSubdomainsGenerator
    input = rename_core_bdy
    primary_block = 3
    paired_block = 2
    new_boundary = 'rpv_inner'
  []
  [2d_mesh]
    type = BlockDeletionGenerator
    input = rename_inner_rpv_bdy
    block = 2
  []
[]

[Variables]
  [Tsolid]
    initial_condition = 500
  []
[]

[Kernels]
  [heat_source]
    type = CoupledForce
    variable = Tsolid
    block = '${core_blocks}'
    v = power_density
  []
  [heat_conduction]
    type = HeatConduction
    variable = Tsolid
  []
[]

[BCs]
  [RPV_out_BC] # k \nabla T = h (T- T_inf) at RPV outer boundary
    type = ConvectiveFluxFunction # (Robin BC)
    variable = Tsolid
    boundary = 'outer' # outer RPV
    coefficient = ${rpv_outer_htc}
    T_infinity = ${rpv_outer_Tinf}
  []
[]

[ThermalContact]
  [RPV_gap]
    type = GapHeatTransfer
    gap_geometry_type = 'CYLINDER'
    emissivity_primary = 0.8
    emissivity_secondary = 0.8
    variable = Tsolid
    primary = 'core_outer'
    secondary = 'rpv_inner'
    gap_conductivity = 0.1
    quadrature = true
    cylinder_axis_point_1 = '0 0 0'
    cylinder_axis_point_2 = '0 0 5'
  []
[]

[AuxVariables]
  [power_density]
    block = '${core_blocks}'
    initial_condition = 50e3
  []
[]

[Materials]
  [simple_mat]
    type = HeatConductionMaterial
    thermal_conductivity = 34.6 # W/m/K
  []
[]

[Postprocessors]
  [Tcore_avg]
    type = ElementAverageValue
    variable = Tsolid
    block = '${core_blocks}'
  []
  [Tcore_max]
    type = ElementExtremeValue
    value_type = max
    variable = Tsolid
    block = '${core_blocks}'
  []
  [Tcore_min]
    type = ElementExtremeValue
    value_type = min
    variable = Tsolid
    block = '${core_blocks}'
  []
  [Trpv_avg]
    type = ElementAverageValue
    variable = Tsolid
    block = '${rpv_blocks}'
  []
  [Trpv_max]
    type = ElementExtremeValue
    value_type = max
    variable = Tsolid
    block = '${rpv_blocks}'
  []
  [Trpv_min]
    type = ElementExtremeValue
    value_type = min
    variable = Tsolid
    block = '${rpv_blocks}'
  []
  [ptot]
    type = ElementIntegralVariablePostprocessor
    variable = power_density
    block = '${core_blocks}'
  []
  [rpv_convective_out]
    type = ConvectiveHeatTransferSideIntegral
    T_solid = Tsolid
    boundary = 'outer' # outer RVP
    T_fluid = ${rpv_outer_Tinf}
    htc = ${rpv_outer_htc}
  []
  [heat_balance] # should be equal to 0 upon convergence
    type = ParsedPostprocessor
    function = '(rpv_convective_out - ptot) / ptot'
    pp_names = 'rpv_convective_out ptot'
  []
[]

[Executioner]
  type = Steady
  automatic_scaling = true
  compute_scaling_once = false
  petsc_options_iname = '-pc_type -pc_hypre_type -ksp_gmres_restart '
  petsc_options_value = 'hypre boomeramg 100'

  nl_rel_tol = 1e-10
  nl_abs_tol = 1e-10
  l_max_its = 100

  [Quadrature]
    side_order = seventh
  []
  line_search = none
[]

[Outputs]
  exodus = false
  csv = true
[]

(modules/heat_conduction/test/tests/gap_heat_transfer_mortar/gap_heat_transfer_sphere3D.i)

sphere_outer_htc = 10 # W/m^2/K
sphere_outer_Tinf = 300 # K

[GlobalParams]
  order = SECOND
  family = LAGRANGE
[]

[Mesh]
  file = sphere3D.e
[]

[Functions]
  [temp]
    type = PiecewiseLinear
    x = '0   1'
    y = '100 200'
  []
[]

[Variables]
  [temp]
    initial_condition = 500
  []
[]

[AuxVariables]
  [gap_conductance]
    order = CONSTANT
    family = MONOMIAL
  []
  [power_density]
    block = 'fuel'
    initial_condition = 50e3
  []
[]

[Kernels]
  [heat_conduction]
    type = HeatConduction
    variable = temp
  []
  [heat_source]
    type = CoupledForce
    variable = temp
    block = 'fuel'
    v = power_density
  []
[]

[AuxKernels]
  [gap_cond]
    type = MaterialRealAux
    property = gap_conductance
    variable = gap_conductance
    boundary = 2
  []
[]

[Materials]
  [heat1]
    type = HeatConductionMaterial
    block = '1 2'
    specific_heat = 1.0
    thermal_conductivity = 34.6
  []
[]

[ThermalContact]
  [thermal_contact]
    type = GapHeatTransfer
    variable = temp
    primary = 3
    secondary = 2
    emissivity_primary = 0
    emissivity_secondary = 0
    gap_conductivity = 5
    gap_geometry_type = SPHERE
    sphere_origin = '0 0 0'
  []
[]

[BCs]
  [RPV_out_BC] # k \nabla T = h (T- T_inf) at RPV outer boundary
    type = ConvectiveFluxFunction # (Robin BC)
    variable = temp
    boundary = '4' # outer RPV
    coefficient = ${sphere_outer_htc}
    T_infinity = ${sphere_outer_Tinf}
  []
[]

[Executioner]
  type = Transient
  solve_type = 'PJFNK'

  petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
  petsc_options_value = 'lu       superlu_dist'

  dt = 1
  dtmin = 0.01
  end_time = 1

  nl_rel_tol = 1e-12
  nl_abs_tol = 1e-7

  [Quadrature]
    order = fifth
    side_order = seventh
  []
[]

[Outputs]
  exodus = true
  csv = true
  [Console]
    type = Console
  []
[]

[Postprocessors]
  [temp_left]
    type = SideAverageValue
    boundary = 2
    variable = temp
  []

  [temp_right]
    type = SideAverageValue
    boundary = 3
    variable = temp
  []

  [flux_left]
    type = SideDiffusiveFluxIntegral
    variable = temp
    boundary = 2
    diffusivity = thermal_conductivity
  []
  [flux_right]
    type = SideDiffusiveFluxIntegral
    variable = temp
    boundary = 3
    diffusivity = thermal_conductivity
  []
  [ptot]
    type = ElementIntegralVariablePostprocessor
    variable = power_density
    block = 'fuel'
  []
  [sphere_convective_out]
    type = ConvectiveHeatTransferSideIntegral
    T_solid = temp
    boundary = '4' # outer RVP
    T_fluid = ${sphere_outer_Tinf}
    htc = ${sphere_outer_htc}
  []
  [heat_balance] # should be equal to 0 upon convergence
    type = ParsedPostprocessor
    function = '(sphere_convective_out - ptot) / ptot'
    pp_names = 'sphere_convective_out ptot'
  []
[]

[VectorPostprocessors]
  [NodalTemperature]
    type = NodalValueSampler
    sort_by = id
    boundary = '2 3'
    variable = temp
  []
[]

(test/tests/postprocessors/parsed_postprocessor/parsed_pp.i)

[Mesh]
  type = GeneratedMesh
  dim = 2
  nx = 10
  ny = 10
[]

[Variables]
  [./u]
  [../]
[]

[Kernels]
  [./diff]
    type = Diffusion
    variable = u
  [../]
[]

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

[Postprocessors]
  [./L2_norm]
    type = ElementL2Norm
    variable = u
  [../]
  [L1_norm]
    type = ElementL1Error
    function = 0
    variable = u
  []
  [parsed]
    type = ParsedPostprocessor
    function = 'L2_norm / L1_norm'
    pp_names = 'L2_norm L1_norm'
  []
  [parsed_with_t]
    type = ParsedPostprocessor
    function = 'L2_norm + L1_norm + t'
    pp_names = 'L2_norm L1_norm'
    use_t = true
  []
  [parsed_with_constants]
    type = ParsedPostprocessor
    function = 'L2_norm + 3*L1_norm + mu'
    pp_names = 'L2_norm L1_norm'
    constant_names = 'mu'
    constant_expressions = '4'
  []
[]

[Executioner]
  type = Transient
  num_steps = 4
  nl_abs_tol = 1e-8
  solve_type = PJFNK
  petsc_options_iname = '-pc_type -pc_hypre_type'
  petsc_options_value = 'hypre boomeramg'
  line_search = 'none'
[]

[Outputs]
  csv = true
[]

(modules/heat_conduction/test/tests/gap_heat_transfer_mortar/large_gap_heat_transfer_test_rz_cylinder.i)

rpv_core_gap_size = 0.2

core_outer_radius = 2
rpv_inner_radius = '${fparse 2 + rpv_core_gap_size}'
rpv_outer_radius = '${fparse 2.5 + rpv_core_gap_size}'
rpv_width = '${fparse rpv_outer_radius - rpv_inner_radius}'

rpv_outer_htc = 10 # W/m^2/K
rpv_outer_Tinf = 300 # K

core_blocks = '1'
rpv_blocks = '3'

[Mesh]
  [gmg]
    type = CartesianMeshGenerator
    dim = 2
    dx = '${core_outer_radius} ${rpv_core_gap_size} ${rpv_width}'
    ix = '400 1 100'
    dy = 1
    iy = '5'
  []
  [set_block_id1]
    type = SubdomainBoundingBoxGenerator
    input = gmg
    bottom_left = '0 0 0'
    top_right = '${core_outer_radius} 1 0'
    block_id = 1
    location = INSIDE
  []
  [rename_core_bdy]
    type = SideSetsBetweenSubdomainsGenerator
    input = set_block_id1
    primary_block = 1
    paired_block = 0
    new_boundary = 'core_outer'
  []
  [set_block_id3]
    type = SubdomainBoundingBoxGenerator
    input = rename_core_bdy
    bottom_left = '${rpv_inner_radius} 0 0'
    top_right = '${rpv_outer_radius} 1 0'
    block_id = 3
    location = INSIDE
  []
  [rename_inner_rpv_bdy]
    type = SideSetsBetweenSubdomainsGenerator
    input = set_block_id3
    primary_block = 3
    paired_block = 0
    new_boundary = 'rpv_inner'
  []
  # comment out for test without gap
  [2d_mesh]
    type = BlockDeletionGenerator
    input = rename_inner_rpv_bdy
    block = 0
  []
  allow_renumbering = false
[]

[Problem]
  coord_type = RZ
[]

[Variables]
  [Tsolid]
    initial_condition = 500
  []
[]

[Kernels]
  [heat_source]
    type = CoupledForce
    variable = Tsolid
    block = '${core_blocks}'
    v = power_density
  []
  [heat_conduction]
    type = HeatConduction
    variable = Tsolid
  []
[]

[BCs]
  [RPV_out_BC] # k \nabla T = h (T- T_inf) at RPV outer boundary
    type = ConvectiveFluxFunction # (Robin BC)
    variable = Tsolid
    boundary = 'right' # outer RPV
    coefficient = ${rpv_outer_htc}
    T_infinity = ${rpv_outer_Tinf}
  []
[]

[ThermalContact]
  [RPV_gap]
    type = GapHeatTransfer
    gap_geometry_type = 'CYLINDER'
    emissivity_primary = 0.8
    emissivity_secondary = 0.8
    variable = Tsolid
    primary = 'core_outer'
    secondary = 'rpv_inner'
    gap_conductivity = 0.1
    quadrature = true
  []
[]

[AuxVariables]
  [power_density]
    block = '${core_blocks}'
    initial_condition = 50e3
  []
[]

[Materials]
  [simple_mat]
    type = HeatConductionMaterial
    thermal_conductivity = 34.6 # W/m/K
  []
[]

[Postprocessors]
  [Tcore_avg]
    type = ElementAverageValue
    variable = Tsolid
    block = '${core_blocks}'
  []
  [Tcore_max]
    type = ElementExtremeValue
    value_type = max
    variable = Tsolid
    block = '${core_blocks}'
  []
  [Tcore_min]
    type = ElementExtremeValue
    value_type = min
    variable = Tsolid
    block = '${core_blocks}'
  []
  [Trpv_avg]
    type = ElementAverageValue
    variable = Tsolid
    block = '${rpv_blocks}'
  []
  [Trpv_max]
    type = ElementExtremeValue
    value_type = max
    variable = Tsolid
    block = '${rpv_blocks}'
  []
  [Trpv_min]
    type = ElementExtremeValue
    value_type = min
    variable = Tsolid
    block = '${rpv_blocks}'
  []
  [ptot]
    type = ElementIntegralVariablePostprocessor
    variable = power_density
    block = '${core_blocks}'
  []
  [rpv_convective_out]
    type = ConvectiveHeatTransferSideIntegral
    T_solid = Tsolid
    boundary = 'right' # outer RVP
    T_fluid = ${rpv_outer_Tinf}
    htc = ${rpv_outer_htc}
  []
  [heat_balance] # should be equal to 0 upon convergence
    type = ParsedPostprocessor
    function = '(rpv_convective_out - ptot) / ptot'
    pp_names = 'rpv_convective_out ptot'
  []
  [flux_from_core] # converges to ptot as the mesh is refined
    type = SideDiffusiveFluxIntegral
    variable = Tsolid
    boundary = core_outer
    diffusivity = thermal_conductivity
  []
  [flux_into_rpv] # converges to rpv_convective_out as the mesh is refined
    type = SideDiffusiveFluxIntegral
    variable = Tsolid
    boundary = rpv_inner
    diffusivity = thermal_conductivity
  []
[]

[VectorPostprocessors]
  [NodalTemperature]
    type = NodalValueSampler
    sort_by = id
    boundary = 'rpv_inner core_outer'
    variable = Tsolid
  []
[]

[Executioner]
  type = Steady
  automatic_scaling = true
  compute_scaling_once = false
  petsc_options_iname = '-pc_type -pc_hypre_type -ksp_gmres_restart '
  petsc_options_value = 'hypre boomeramg 100'
  nl_rel_tol = 1e-10
  nl_abs_tol = 1e-10
  l_max_its = 100
  [Quadrature]
    # order = fifth
    side_order = seventh
  []
  line_search = none
[]

[Outputs]
  exodus = false
  csv = true
[]

(modules/thermal_hydraulics/test/tests/components/simple_turbine_1phase/phy.conservation.i)

[GlobalParams]
  initial_p = 1e6
  initial_T = 517
  initial_vel = 4.3
  initial_vel_x = 4.3
  initial_vel_y = 0
  initial_vel_z = 0

  fp = fp

  closures = simple_closures
  f = 0

  rdg_slope_reconstruction = minmod
  gravity_vector = '0 0 0'
[]

[Modules/FluidProperties]
  [fp]
    type = IdealGasFluidProperties
    gamma = 1.4
    molar_mass = 0.01
  []
[]

[Closures]
  [simple_closures]
    type = Closures1PhaseSimple
  []
[]

[Components]
  [inlet]
    type = InletMassFlowRateTemperature1Phase
    input = 'pipe1:in'
    m_dot = 10
    T = 517
  []

  [pipe1]
    type = FlowChannel1Phase
    position = '0 0 0'
    orientation = '1 0 0'
    length = 1
    n_elems = 10
    A = 1
  []

  [turbine]
    type = SimpleTurbine1Phase
    connections = 'pipe1:out pipe2:in'
    position = '1 0 0'
    volume = 1
    A_ref = 1.0
    K = 0
    on = true
    power = 1000
  []

  [pipe2]
    type = FlowChannel1Phase
    position = '1. 0 0'
    orientation = '1 0 0'
    length = 1
    n_elems = 10
    A = 1
  []

  [outlet]
    type = Outlet1Phase
    input = 'pipe2:out'
    p = 1e6
  []
[]

[Postprocessors]
  [mass_in]
    type = ADFlowBoundaryFlux1Phase
    equation = mass
    boundary = inlet
  []
  [mass_out]
    type = ADFlowBoundaryFlux1Phase
    equation = mass
    boundary = outlet
  []
  [mass_diff]
    type = LinearCombinationPostprocessor
    pp_coefs = '1 -1'
    pp_names = 'mass_in mass_out'
  []
  [p_in]
    type = SideAverageValue
    boundary = pipe1:in
    variable = p
  []
  [vel_in]
    type = SideAverageValue
    boundary = pipe1:in
    variable = vel_x
  []
  [momentum_in]
    type = ADFlowBoundaryFlux1Phase
    equation = momentum
    boundary = inlet
  []
  [momentum_out]
    type = ADFlowBoundaryFlux1Phase
    equation = momentum
    boundary = outlet
  []
  [dP]
    type = ParsedPostprocessor
    pp_names = 'p_in W_dot'
    function = 'p_in * (1 - (1-W_dot/(10*2910.06*517))^(1.4/0.4))'
  []
  [momentum_diff]
    type = LinearCombinationPostprocessor
    pp_coefs = '1 -1 -1'
    pp_names = 'momentum_in momentum_out dP' # momentum source = -dP * A and A=1
  []

  [energy_in]
    type = ADFlowBoundaryFlux1Phase
    equation = energy
    boundary = inlet
  []
  [energy_out]
    type = ADFlowBoundaryFlux1Phase
    equation = energy
    boundary = outlet
  []
  [W_dot]
    type = ScalarVariable
    variable = turbine:W_dot
  []
  [energy_diff]
    type = LinearCombinationPostprocessor
    pp_coefs = '1 -1 -1'
    pp_names = 'energy_in energy_out W_dot'
  []
[]

[Preconditioning]
  [pc]
    type = SMP
    full = true
  []
[]

[Executioner]
  type = Transient
  scheme = bdf2

  start_time = 0
  dt = 1
  num_steps = 30

  abort_on_solve_fail = true

  solve_type = 'newton'
  line_search = 'basic'
  petsc_options_iname = '-pc_type'
  petsc_options_value = ' lu'

  nl_rel_tol = 1e-7
  nl_abs_tol = 2e-6

  nl_max_its = 10
  l_tol = 1e-3
[]

[Outputs]
  [csv]
    type = CSV
    show = 'mass_diff energy_diff momentum_diff'
    execute_on = 'final'
  []
[]

(modules/heat_conduction/test/tests/gap_heat_transfer_mortar/gap_heat_transfer_3D_mortar.i)

outer_htc = 10 # W/m^2/K
outer_Tinf = 300 # K

[GlobalParams]
  order = SECOND
  family = LAGRANGE
[]

[Problem]
  kernel_coverage_check = false
  material_coverage_check = false
[]

[Mesh]
  [left_block]
    type = GeneratedMeshGenerator
    dim = 3
    nx = 3
    ny = 6
    nz = 6
    xmin = -1
    xmax = -0.5
    ymin = -0.5
    ymax = 0.5
    zmin = -0.5
    zmax = 0.5
    elem_type = HEX27
  []
  [left_block_sidesets]
    type = RenameBoundaryGenerator
    input = left_block
    old_boundary = '0 1 2 3 4 5'
    new_boundary = 'left_bottom left_back left_right left_front left_left left_top'
  []
  [left_block_id]
    type = SubdomainIDGenerator
    input = left_block_sidesets
    subdomain_id = 1
  []

  [right_block]
    type = GeneratedMeshGenerator
    dim = 3
    nx = 4
    ny = 8
    nz = 8
    xmin = 0.5
    xmax = 1
    ymin = -0.5
    ymax = 0.5
    zmin = -0.5
    zmax = 0.5
    elem_type = HEX27
  []
  [right_block_sidesets]
    type = RenameBoundaryGenerator
    input = right_block
    old_boundary = '0 1 2 3 4 5'
    # new_boundary = 'right_bottom right_back right_right right_front right_left right_top'
    new_boundary = '100 101 102 103 104 105'
  []
  [right_block_sidesets_rename]
    type = RenameBoundaryGenerator
    input = right_block_sidesets
    old_boundary = '100 101 102 103 104 105'
    new_boundary = 'right_bottom right_back right_right right_front right_left right_top'
  []
  [right_block_id]
    type = SubdomainIDGenerator
    input = right_block_sidesets_rename
    subdomain_id = 2
  []

  [combined_mesh]
    type = MeshCollectionGenerator
    inputs = 'left_block_id right_block_id'
  []

  [left_lower]
    type = LowerDBlockFromSidesetGenerator
    input = combined_mesh
    sidesets = 'left_right'
    new_block_id = '10001'
    new_block_name = 'secondary_lower'
  []
  [right_lower]
    type = LowerDBlockFromSidesetGenerator
    input = left_lower
    sidesets = 'right_left'
    new_block_id = '10000'
    new_block_name = 'primary_lower'
  []
[]

[Functions]
  [temp]
    type = PiecewiseLinear
    x = '0   1'
    y = '100 200'
  []
[]

[Variables]
  [temp]
    initial_condition = 500
  []
  [lm]
    order = SECOND
    family = LAGRANGE
    block = 'secondary_lower'
  []
[]

[AuxVariables]
  [power_density]
    block = 1
    initial_condition = 50e3
  []
[]

[Kernels]
  [heat_conduction]
    type = HeatConduction
    variable = temp
    block = '1 2'
  []
  [heat_source]
    type = CoupledForce
    variable = temp
    block = '1'
    v = power_density
  []
[]

[Materials]
  [heat1]
    type = HeatConductionMaterial
    block = '1 2'
    specific_heat = 1.0
    thermal_conductivity = 34.6
  []
[]

[UserObjects]
  [radiation]
    type = GapFluxModelRadiation
    T = temp
    boundary = 'left_right'
    primary_emissivity = 0.0
    secondary_emissivity = 0.0
  []
  [conduction]
    type = GapFluxModelConduction
    T = temp
    boundary = 'left_right'
    gap_conductivity = 5.0
  []
[]

[Constraints]
  [ced]
    type = ModularGapConductanceConstraint
    variable = lm
    secondary_variable = temp
    primary_boundary = 'right_left'
    primary_subdomain = 'primary_lower'
    secondary_boundary = 'left_right'
    secondary_subdomain = 'secondary_lower'
    gap_flux_models = 'radiation conduction'
    gap_geometry_type = PLATE
  []
[]

[BCs]
  [RPV_out_BC] # k \nabla T = h (T- T_inf) at RPV outer boundary
    type = ConvectiveFluxFunction # (Robin BC)
    variable = temp
    boundary = 'right_right' # outer RPV
    coefficient = ${outer_htc}
    T_infinity = ${outer_Tinf}
  []
[]

[Executioner]
  type = Transient
  solve_type = 'PJFNK'

  petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
  petsc_options_value = 'lu       superlu_dist'

  dt = 1
  dtmin = 0.01
  end_time = 1

  nl_rel_tol = 1e-12
  nl_abs_tol = 1e-8
[]

[Outputs]
  exodus = true
  csv = true
  [Console]
    type = Console
  []
[]

[Postprocessors]
  [temp_left]
    type = SideAverageValue
    boundary = 'left_right'
    variable = temp
  []

  [temp_right]
    type = SideAverageValue
    boundary = 'right_left'
    variable = temp
  []

  [flux_left]
    type = SideDiffusiveFluxIntegral
    variable = temp
    boundary = 'left_right'
    diffusivity = thermal_conductivity
  []
  [flux_right]
    type = SideDiffusiveFluxIntegral
    variable = temp
    boundary = 'right_left'
    diffusivity = thermal_conductivity
  []
  [ptot]
    type = ElementIntegralVariablePostprocessor
    variable = power_density
    block = 1
  []
  [convective_out]
    type = ConvectiveHeatTransferSideIntegral
    T_solid = temp
    boundary = 'right_right' # outer RVP
    T_fluid = ${outer_Tinf}
    htc = ${outer_htc}
  []
  [heat_balance] # should be equal to 0 upon convergence
    type = ParsedPostprocessor
    function = '(convective_out - ptot) / ptot'
    pp_names = 'convective_out ptot'
  []
[]

[VectorPostprocessors]
  [NodalTemperature]
    type = NodalValueSampler
    sort_by = id
    boundary = 'left_right right_left'
    variable = temp
  []
[]

(modules/heat_conduction/test/tests/gap_heat_transfer_balance/large_gap_heat_transfer_test_rz_cylinder.i)

rpv_core_gap_size = 0.2

core_outer_radius = 2
rpv_inner_radius = '${fparse 2 + rpv_core_gap_size}'
rpv_outer_radius = '${fparse 2.5 + rpv_core_gap_size}'
rpv_width = '${fparse rpv_outer_radius - rpv_inner_radius}'

rpv_outer_htc = 10 # W/m^2/K
rpv_outer_Tinf = 300 # K

core_blocks = '1'
rpv_blocks = '3'

[Mesh]
  [gmg]
    type = CartesianMeshGenerator
    dim = 2
    dx = '${core_outer_radius} ${rpv_core_gap_size} ${rpv_width}'
    ix = '400 1 100'
    dy = 1
    iy = '5'
  []
  [set_block_id1]
    type = SubdomainBoundingBoxGenerator
    input = gmg
    bottom_left = '0 0 0'
    top_right = '${core_outer_radius} 1 0'
    block_id = 1
    location = INSIDE
  []
  [rename_core_bdy]
    type = SideSetsBetweenSubdomainsGenerator
    input = set_block_id1
    primary_block = 1
    paired_block = 0
    new_boundary = 'core_outer'
  []
  [set_block_id3]
    type = SubdomainBoundingBoxGenerator
    input = rename_core_bdy
    bottom_left = '${rpv_inner_radius} 0 0'
    top_right = '${rpv_outer_radius} 1 0'
    block_id = 3
    location = INSIDE
  []
  [rename_inner_rpv_bdy]
    type = SideSetsBetweenSubdomainsGenerator
    input = set_block_id3
    primary_block = 3
    paired_block = 0
    new_boundary = 'rpv_inner'
  []
  # comment out for test without gap
  [2d_mesh]
    type = BlockDeletionGenerator
    input = rename_inner_rpv_bdy
    block = 0
  []
[]

[Problem]
  coord_type = RZ
[]

[Variables]
  [Tsolid]
    initial_condition = 500
  []
[]

[Kernels]
  [heat_source]
    type = CoupledForce
    variable = Tsolid
    block = '${core_blocks}'
    v = power_density
  []
  [heat_conduction]
    type = HeatConduction
    variable = Tsolid
  []
[]

[BCs]
  [RPV_out_BC] # k \nabla T = h (T- T_inf) at RPV outer boundary
    type = ConvectiveFluxFunction # (Robin BC)
    variable = Tsolid
    boundary = 'right' # outer RPV
    coefficient = ${rpv_outer_htc}
    T_infinity = ${rpv_outer_Tinf}
  []
[]

[ThermalContact]
  [RPV_gap]
    type = GapHeatTransfer
    gap_geometry_type = 'CYLINDER'
    emissivity_primary = 0.8
    emissivity_secondary = 0.8
    variable = Tsolid
    primary = 'core_outer'
    secondary = 'rpv_inner'
    gap_conductivity = 0.1
    quadrature = true
  []
[]

[AuxVariables]
  [power_density]
    block = '${core_blocks}'
    initial_condition = 50e3
  []
[]

[Materials]
  [simple_mat]
    type = HeatConductionMaterial
    thermal_conductivity = 34.6 # W/m/K
  []
[]

[Postprocessors]
  [Tcore_avg]
    type = ElementAverageValue
    variable = Tsolid
    block = '${core_blocks}'
  []
  [Tcore_max]
    type = ElementExtremeValue
    value_type = max
    variable = Tsolid
    block = '${core_blocks}'
  []
  [Tcore_min]
    type = ElementExtremeValue
    value_type = min
    variable = Tsolid
    block = '${core_blocks}'
  []
  [Trpv_avg]
    type = ElementAverageValue
    variable = Tsolid
    block = '${rpv_blocks}'
  []
  [Trpv_max]
    type = ElementExtremeValue
    value_type = max
    variable = Tsolid
    block = '${rpv_blocks}'
  []
  [Trpv_min]
    type = ElementExtremeValue
    value_type = min
    variable = Tsolid
    block = '${rpv_blocks}'
  []
  [ptot]
    type = ElementIntegralVariablePostprocessor
    variable = power_density
    block = '${core_blocks}'
  []
  [rpv_convective_out]
    type = ConvectiveHeatTransferSideIntegral
    T_solid = Tsolid
    boundary = 'right' # outer RVP
    T_fluid = ${rpv_outer_Tinf}
    htc = ${rpv_outer_htc}
  []
  [heat_balance] # should be equal to 0 upon convergence
    type = ParsedPostprocessor
    function = '(rpv_convective_out - ptot) / ptot'
    pp_names = 'rpv_convective_out ptot'
  []
  [flux_from_core] # converges to ptot as the mesh is refined
    type = SideDiffusiveFluxIntegral
    variable = Tsolid
    boundary = core_outer
    diffusivity = thermal_conductivity
  []
  [flux_into_rpv] # converges to rpv_convective_out as the mesh is refined
    type = SideDiffusiveFluxIntegral
    variable = Tsolid
    boundary = rpv_inner
    diffusivity = thermal_conductivity
  []

[]

[Executioner]
  type = Steady
  automatic_scaling = true
  compute_scaling_once = false
  petsc_options_iname = '-pc_type -pc_hypre_type -ksp_gmres_restart '
  petsc_options_value = 'hypre boomeramg 100'
  nl_rel_tol = 1e-10
  nl_abs_tol = 1e-10
  l_max_its = 100
  [Quadrature]
    # order = fifth
    side_order = seventh
  []
  line_search = none
[]

[Outputs]
  exodus = false
  csv = true
[]

(modules/thermal_hydraulics/test/tests/problems/brayton_cycle/closed_brayton_cycle.i)

# This input file is used to demonstrate a simple closed, air Brayton cycle using
# a compressor, turbine, shaft, motor, and generator.
# The flow length is divided into 6 segments as illustrated below, where
#   - "(C)" denotes the compressor
#   - "(T)" denotes the turbine
#   - "*" denotes a fictitious junction
#
#                Heated section               Cooled section
# *-----(C)-----*--------------*-----(T)-----*--------------*
#    1       2         3          4       5         6
#
# Initially the fluid is at rest at ambient conditions, the shaft speed is zero,
# and no heat transfer occurs with the system.
# The transient is controlled as follows:
#   * 0   - 100 s: motor ramps up torque linearly from zero
#   * 100 - 200 s: motor ramps down torque linearly to zero, HTC ramps up linearly from zero.
#   * 200 - 300 s: (no changes; should approach steady condition)

I_motor = 1.0
motor_torque_max = 400.0

I_generator = 1.0
generator_torque_per_shaft_speed = -0.00025

motor_ramp_up_duration = 100.0
motor_ramp_down_duration = 100.0
post_motor_time = 100.0
t1 = ${motor_ramp_up_duration}
t2 = ${fparse t1 + motor_ramp_down_duration}
t3 = ${fparse t2 + post_motor_time}

D1 = 0.15
D2 = ${D1}
D3 = ${D1}
D4 = ${D1}
D5 = ${D1}
D6 = ${D1}

A1 = ${fparse 0.25 * pi * D1^2}
A2 = ${fparse 0.25 * pi * D2^2}
A3 = ${fparse 0.25 * pi * D3^2}
A4 = ${fparse 0.25 * pi * D4^2}
A5 = ${fparse 0.25 * pi * D5^2}
A6 = ${fparse 0.25 * pi * D6^2}

L1 = 10.0
L2 = ${L1}
L3 = ${L1}
L4 = ${L1}
L5 = ${L1}
L6 = ${L1}

x1 = 0.0
x2 = ${fparse x1 + L1}
x3 = ${fparse x2 + L2}
x4 = ${fparse x3 + L3}
x5 = ${fparse x4 + L4}
x6 = ${fparse x5 + L5}

x2_minus = ${fparse x2 - 0.001}
x2_plus = ${fparse x2 + 0.001}
x5_minus = ${fparse x5 - 0.001}
x5_plus = ${fparse x5 + 0.001}

n_elems1 = 10
n_elems2 = ${n_elems1}
n_elems3 = ${n_elems1}
n_elems4 = ${n_elems1}
n_elems5 = ${n_elems1}
n_elems6 = ${n_elems1}

A_ref_comp = ${fparse 0.5 * (A1 + A2)}
V_comp = ${fparse A_ref_comp * 1.0}
I_comp = 1.0

A_ref_turb = ${fparse 0.5 * (A4 + A5)}
V_turb = ${fparse A_ref_turb * 1.0}
I_turb = 1.0

c0_rated_comp = 351.6925137
rho0_rated_comp = 1.146881112

rated_mfr = 0.25

speed_rated_rpm = 96000
speed_rated = ${fparse speed_rated_rpm * 2 * pi / 60.0}

speed_initial = 0

eff_comp = 0.79
eff_turb = 0.843

T_hot = 1000
T_cold = 300
T_ambient = 300
p_ambient = 1e5

[GlobalParams]
  orientation = '1 0 0'
  gravity_vector = '0 0 0'

  initial_p = ${p_ambient}
  initial_T = ${T_ambient}
  initial_vel = 0
  initial_vel_x = 0
  initial_vel_y = 0
  initial_vel_z = 0

  fp = fp_air
  closures = closures
  f = 0

  scaling_factor_1phase = '1 1 1e-5'
  scaling_factor_rhoV = 1
  scaling_factor_rhouV = 1
  scaling_factor_rhovV = 1
  scaling_factor_rhowV = 1
  scaling_factor_rhoEV = 1e-5

  rdg_slope_reconstruction = none
[]

[Functions]
  [motor_torque_fn]
    type = PiecewiseLinear
    x = '0 ${t1} ${t2}'
    y = '0 ${motor_torque_max} 0'
  []
  [motor_power_fn]
    type = ParsedFunction
    value = 'torque * speed'
    vars = 'torque speed'
    vals = 'motor_torque shaft:omega'
  []
  [generator_torque_fn]
    type = ParsedFunction
    value = 'slope * t'
    vars = 'slope'
    vals = '${generator_torque_per_shaft_speed}'
  []
  [generator_power_fn]
    type = ParsedFunction
    value = 'torque * speed'
    vars = 'torque speed'
    vals = 'generator_torque shaft:omega'
  []
  [htc_wall_fn]
    type = PiecewiseLinear
    x = '0 ${t1} ${t2}'
    y = '0 0 1e3'
  []
[]

[Modules/FluidProperties]
  [fp_air]
    type = IdealGasFluidProperties
    emit_on_nan = none
  []
[]

[Closures]
  [closures]
    type = Closures1PhaseSimple
  []
[]

[Components]
  [shaft]
    type = Shaft
    connected_components = 'motor compressor turbine generator'
    initial_speed = ${speed_initial}
  []
  [motor]
    type = ShaftConnectedMotor
    inertia = ${I_motor}
    torque = 0 # controlled
  []
  [generator]
    type = ShaftConnectedMotor
    inertia = ${I_generator}
    torque = generator_torque_fn
  []

  [pipe1]
    type = FlowChannel1Phase
    position = '${x1} 0 0'
    length = ${L1}
    n_elems = ${n_elems1}
    A = ${A1}
  []
  [compressor]
    type = ShaftConnectedCompressor1Phase
    position = '${x2} 0 0'
    inlet = 'pipe1:out'
    outlet = 'pipe2:in'
    A_ref = ${A_ref_comp}
    volume = ${V_comp}

    omega_rated = ${speed_rated}
    mdot_rated = ${rated_mfr}
    c0_rated = ${c0_rated_comp}
    rho0_rated = ${rho0_rated_comp}

    speeds = '0.5208 0.6250 0.7292 0.8333 0.9375'
    Rp_functions = 'rp_comp1 rp_comp2 rp_comp3 rp_comp4 rp_comp5'
    eff_functions = 'eff_comp1 eff_comp2 eff_comp3 eff_comp4 eff_comp5'

    min_pressure_ratio = 1.0

    speed_cr_I = 0
    inertia_const = ${I_comp}
    inertia_coeff = '${I_comp} 0 0 0'

    # assume no shaft friction
    speed_cr_fr = 0
    tau_fr_const = 0
    tau_fr_coeff = '0 0 0 0'
  []
  [pipe2]
    type = FlowChannel1Phase
    position = '${x2} 0 0'
    length = ${L2}
    n_elems = ${n_elems2}
    A = ${A2}
  []
  [junction2_3]
    type = JunctionOneToOne1Phase
    connections = 'pipe2:out pipe3:in'
  []
  [pipe3]
    type = FlowChannel1Phase
    position = '${x3} 0 0'
    length = ${L3}
    n_elems = ${n_elems3}
    A = ${A3}
  []
  [junction3_4]
    type = JunctionOneToOne1Phase
    connections = 'pipe3:out pipe4:in'
  []
  [pipe4]
    type = FlowChannel1Phase
    position = '${x4} 0 0'
    length = ${L4}
    n_elems = ${n_elems4}
    A = ${A4}
  []
  [turbine]
    type = ShaftConnectedCompressor1Phase
    position = '${x5} 0 0'
    inlet = 'pipe4:out'
    outlet = 'pipe5:in'
    A_ref = ${A_ref_turb}
    volume = ${V_turb}

    treat_as_turbine = true

    omega_rated = ${speed_rated}
    mdot_rated = ${rated_mfr}
    c0_rated = ${c0_rated_comp}
    rho0_rated = ${rho0_rated_comp}

    speeds = '0 0.5208 0.6250 0.7292 0.8333 0.9375'
    Rp_functions = 'rp_turb0 rp_turb1 rp_turb2 rp_turb3 rp_turb4 rp_turb5'
    eff_functions = 'eff_turb1 eff_turb1 eff_turb2 eff_turb3 eff_turb4 eff_turb5'

    min_pressure_ratio = 1.0

    speed_cr_I = 0
    inertia_const = ${I_turb}
    inertia_coeff = '${I_turb} 0 0 0'

    # assume no shaft friction
    speed_cr_fr = 0
    tau_fr_const = 0
    tau_fr_coeff = '0 0 0 0'
  []
  [pipe5]
    type = FlowChannel1Phase
    position = '${x5} 0 0'
    length = ${L5}
    n_elems = ${n_elems5}
    A = ${A5}
  []
  [junction5_6]
    type = JunctionOneToOne1Phase
    connections = 'pipe5:out pipe6:in'
  []
  [pipe6]
    type = FlowChannel1Phase
    position = '${x6} 0 0'
    length = ${L6}
    n_elems = ${n_elems6}
    A = ${A6}
  []
  [junction6_1]
    type = JunctionOneToOne1Phase
    connections = 'pipe6:out pipe1:in'
  []

  [heating]
    type = HeatTransferFromSpecifiedTemperature1Phase
    flow_channel = pipe3
    T_wall = ${T_hot}
    Hw = htc_wall_fn
  []
  [cooling]
    type = HeatTransferFromSpecifiedTemperature1Phase
    flow_channel = pipe6
    T_wall = ${T_cold}
    Hw = htc_wall_fn
  []
[]

[ControlLogic]
  [motor_ctrl]
    type = TimeFunctionComponentControl
    component = motor
    parameter = torque
    function = motor_torque_fn
  []
[]

[Postprocessors]
  [heating_rate]
    type = ADHeatRateConvection1Phase
    block = 'pipe3'
    T = T
    T_wall = T_wall
    Hw = Hw
    P_hf = P_hf
    execute_on = 'INITIAL TIMESTEP_END'
  []
  [cooling_rate]
    type = ADHeatRateConvection1Phase
    block = 'pipe6'
    T = T
    T_wall = T_wall
    Hw = Hw
    P_hf = P_hf
    execute_on = 'INITIAL TIMESTEP_END'
  []

  [motor_torque]
    type = RealComponentParameterValuePostprocessor
    component = motor
    parameter = torque
    execute_on = 'INITIAL TIMESTEP_END'
  []
  [motor_power]
    type = FunctionValuePostprocessor
    function = motor_power_fn
    execute_on = 'INITIAL TIMESTEP_END'
  []

  [generator_torque]
    type = ShaftConnectedComponentPostprocessor
    quantity = torque
    shaft_connected_component_uo = generator:shaftconnected_uo
    execute_on = 'INITIAL TIMESTEP_END'
  []
  [generator_power]
    type = FunctionValuePostprocessor
    function = generator_power_fn
    execute_on = 'INITIAL TIMESTEP_END'
  []

  [shaft_speed]
    type = ScalarVariable
    variable = 'shaft:omega'
    execute_on = 'INITIAL TIMESTEP_END'
  []

  [p_in_comp]
    type = PointValue
    variable = p
    point = '${x2_minus} 0 0'
    execute_on = 'INITIAL TIMESTEP_END'
  []
  [p_out_comp]
    type = PointValue
    variable = p
    point = '${x2_plus} 0 0'
    execute_on = 'INITIAL TIMESTEP_END'
  []
  [p_ratio_comp]
    type = ParsedPostprocessor
    pp_names = 'p_in_comp p_out_comp'
    function = 'p_out_comp / p_in_comp'
    execute_on = 'INITIAL TIMESTEP_END'
  []

  [p_in_turb]
    type = PointValue
    variable = p
    point = '${x5_minus} 0 0'
    execute_on = 'INITIAL TIMESTEP_END'
  []
  [p_out_turb]
    type = PointValue
    variable = p
    point = '${x5_plus} 0 0'
    execute_on = 'INITIAL TIMESTEP_END'
  []
  [p_ratio_turb]
    type = ParsedPostprocessor
    pp_names = 'p_in_turb p_out_turb'
    function = 'p_in_turb / p_out_turb'
    execute_on = 'INITIAL TIMESTEP_END'
  []

  [mfr_comp]
    type = ADFlowJunctionFlux1Phase
    boundary = pipe1:out
    connection_index = 0
    equation = mass
    junction = compressor
  []
  [mfr_turb]
    type = ADFlowJunctionFlux1Phase
    boundary = pipe4:out
    connection_index = 0
    equation = mass
    junction = turbine
  []
[]

[Preconditioning]
  [pc]
    type = SMP
    full = true
  []
[]

[Executioner]
  type = Transient
  scheme = 'bdf2'

  end_time = ${t3}
  [TimeStepper]
    type = IterationAdaptiveDT
    dt = 0.01
    optimal_iterations = 5
    iteration_window = 1
    growth_factor = 1.1
    cutback_factor = 0.9
  []
  dtmin = 1e-5

  steady_state_detection = true
  steady_state_start_time = ${t2}

  solve_type = NEWTON
  nl_rel_tol = 1e-8
  nl_abs_tol = 1e-8
  nl_max_its = 15

  l_tol = 1e-4
  l_max_its = 10
[]

[Outputs]
  exodus = true
  [csv]
    type = CSV
    file_base = 'closed_brayton_cycle'
    execute_vector_postprocessors_on = 'INITIAL'
  []
  [console]
    type = Console
    show = 'shaft_speed p_ratio_comp p_ratio_turb compressor:pressure_ratio turbine:pressure_ratio'
  []
[]

[Functions]
  # compressor pressure ratio
  [rp_comp1]
    type = PiecewiseLinear
    data_file = 'rp_comp1.csv'
    x_index_in_file = 0
    y_index_in_file = 1
    format = columns
    extrap = true
  []
  [rp_comp2]
    type = PiecewiseLinear
    data_file = 'rp_comp2.csv'
    x_index_in_file = 0
    y_index_in_file = 1
    format = columns
    extrap = true
  []
  [rp_comp3]
    type = PiecewiseLinear
    data_file = 'rp_comp3.csv'
    x_index_in_file = 0
    y_index_in_file = 1
    format = columns
    extrap = true
  []
  [rp_comp4]
    type = PiecewiseLinear
    data_file = 'rp_comp4.csv'
    x_index_in_file = 0
    y_index_in_file = 1
    format = columns
    extrap = true
  []
  [rp_comp5]
    type = PiecewiseLinear
    data_file = 'rp_comp5.csv'
    x_index_in_file = 0
    y_index_in_file = 1
    format = columns
    extrap = true
  []

  # compressor efficiency
  [eff_comp1]
    type = ConstantFunction
    value = ${eff_comp}
  []
  [eff_comp2]
    type = ConstantFunction
    value = ${eff_comp}
  []
  [eff_comp3]
    type = ConstantFunction
    value = ${eff_comp}
  []
  [eff_comp4]
    type = ConstantFunction
    value = ${eff_comp}
  []
  [eff_comp5]
    type = ConstantFunction
    value = ${eff_comp}
  []

  # turbine pressure ratio
  [rp_turb0]
    type = ConstantFunction
    value = 1
  []
  [rp_turb1]
    type = PiecewiseLinear
    data_file = 'rp_turb1.csv'
    x_index_in_file = 0
    y_index_in_file = 1
    format = columns
    extrap = true
  []
  [rp_turb2]
    type = PiecewiseLinear
    data_file = 'rp_turb2.csv'
    x_index_in_file = 0
    y_index_in_file = 1
    format = columns
    extrap = true
  []
  [rp_turb3]
    type = PiecewiseLinear
    data_file = 'rp_turb3.csv'
    x_index_in_file = 0
    y_index_in_file = 1
    format = columns
    extrap = true
  []
  [rp_turb4]
    type = PiecewiseLinear
    data_file = 'rp_turb4.csv'
    x_index_in_file = 0
    y_index_in_file = 1
    format = columns
    extrap = true
  []
  [rp_turb5]
    type = PiecewiseLinear
    data_file = 'rp_turb5.csv'
    x_index_in_file = 0
    y_index_in_file = 1
    format = columns
    extrap = true
  []

  # turbine efficiency
  [eff_turb1]
    type = ConstantFunction
    value = ${eff_turb}
  []
  [eff_turb2]
    type = ConstantFunction
    value = ${eff_turb}
  []
  [eff_turb3]
    type = ConstantFunction
    value = ${eff_turb}
  []
  [eff_turb4]
    type = ConstantFunction
    value = ${eff_turb}
  []
  [eff_turb5]
    type = ConstantFunction
    value = ${eff_turb}
  []
[]

(modules/heat_conduction/test/tests/gap_heat_transfer_mortar/large_gap_heat_transfer_test_sphere_mortar_error.i)

sphere_outer_htc = 10 # W/m^2/K
sphere_outer_Tinf = 300 # K

[GlobalParams]
  order = SECOND
  family = LAGRANGE
[]

[Problem]
  coord_type = RZ
  kernel_coverage_check = false
  material_coverage_check = false
[]

[Mesh]
  [file]
    type = FileMeshGenerator
    file = cyl2D.e
  []
  [secondary]
    type = LowerDBlockFromSidesetGenerator
    sidesets = '2'
    new_block_id = 10001
    new_block_name = 'secondary_lower'
    input = file
  []
  [primary]
    type = LowerDBlockFromSidesetGenerator
    sidesets = '3'
    new_block_id = 10000
    new_block_name = 'primary_lower'
    input = secondary
  []
  allow_renumbering = false
[]

[Functions]
  [temp]
    type = PiecewiseLinear
    x = '0   1'
    y = '100 200'
  []
[]

[Variables]
  [temp]
    initial_condition = 500
  []
  [lm]
    order = SECOND
    family = LAGRANGE
    block = 'secondary_lower'
  []
[]

[AuxVariables]
  [power_density]
    block = 'fuel'
    initial_condition = 50e3
  []
[]

[Kernels]
  [heat_conduction]
    type = HeatConduction
    variable = temp
    block = '1 2'
  []
  [heat_source]
    type = CoupledForce
    variable = temp
    block = 'fuel'
    v = power_density
  []
[]

[Materials]
  [heat1]
    type = HeatConductionMaterial
    block = '1 2'
    specific_heat = 1.0
    thermal_conductivity = 34.6
  []
[]

[UserObjects]
  [radiation]
    type = GapFluxModelRadiation
    temperature = temp
    boundary = 2
    primary_emissivity = 0.0
    secondary_emissivity = 0.0
  []
  [conduction]
    type = GapFluxModelConduction
    temperature = temp
    boundary = 2
    gap_conductivity = 5.0
  []
[]

[Constraints]
  [ced]
    type = ModularGapConductanceConstraint
    variable = lm
    secondary_variable = temp
    primary_boundary = 3
    primary_subdomain = 10000
    secondary_boundary = 2
    secondary_subdomain = 10001
    gap_flux_models = 'radiation conduction'
    gap_geometry_type = SPHERE
  []
[]

[BCs]
  [RPV_out_BC] # k \nabla T = h (T- T_inf) at RPV outer boundary
    type = ConvectiveFluxFunction # (Robin BC)
    variable = temp
    boundary = '4' # outer RPV
    coefficient = ${sphere_outer_htc}
    T_infinity = ${sphere_outer_Tinf}
  []
[]

[Executioner]
  type = Transient
  solve_type = 'PJFNK'

  petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
  petsc_options_value = 'lu       superlu_dist'

  dt = 1
  dtmin = 0.01
  end_time = 1

  nl_rel_tol = 1e-12
  nl_abs_tol = 1e-7

[]

[Outputs]
  exodus = true
  csv = true
  [Console]
    type = Console
  []
[]

[VectorPostprocessors]
  [NodalTemperature]
    type = NodalValueSampler
    sort_by = id
    boundary = '2 3'
    variable = temp
  []
[]

[Postprocessors]
  [temp_left]
    type = SideAverageValue
    boundary = 2
    variable = temp
  []

  [temp_right]
    type = SideAverageValue
    boundary = 3
    variable = temp
  []

  [flux_left]
    type = SideDiffusiveFluxIntegral
    variable = temp
    boundary = 2
    diffusivity = thermal_conductivity
  []

  [flux_right]
    type = SideDiffusiveFluxIntegral
    variable = temp
    boundary = 3
    diffusivity = thermal_conductivity
  []
  [ptot]
    type = ElementIntegralVariablePostprocessor
    variable = power_density
    block = 'fuel'
  []
  [sphere_convective_out]
    type = ConvectiveHeatTransferSideIntegral
    T_solid = temp
    boundary = '4' # outer RVP
    T_fluid = ${sphere_outer_Tinf}
    htc = ${sphere_outer_htc}
  []
  [heat_balance] # should be equal to 0 upon convergence
    type = ParsedPostprocessor
    function = '(sphere_convective_out - ptot) / ptot'
    pp_names = 'sphere_convective_out ptot'
  []
[]

(modules/tensor_mechanics/test/tests/rom_stress_update/ad_verification.i)

[Mesh]
  type = GeneratedMesh
  dim = 3
[]

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

[AuxVariables]
  [temperature]
  []
[]

[AuxKernels]
  [temp_aux]
    type = FunctionAux
    variable = temperature
    function = temp_fcn
    execute_on = 'initial timestep_begin'
  []
[]

[Functions]
  [rhom_fcn]
    type = PiecewiseConstant
    data_file = ss316_verification_data.csv
    x_index_in_file = 0
    y_index_in_file = 1
    format = columns
    xy_in_file_only = false
    direction = LEFT_INCLUSIVE
  []
  [rhoi_fcn]
    type = PiecewiseConstant
    data_file = ss316_verification_data.csv
    x_index_in_file = 0
    y_index_in_file = 2
    format = columns
    xy_in_file_only = false
    direction = LEFT_INCLUSIVE
  []
  [vmJ2_fcn]
    type = PiecewiseConstant
    data_file = ss316_verification_data.csv
    x_index_in_file = 0
    y_index_in_file = 3
    format = columns
    xy_in_file_only = false
    direction = LEFT_INCLUSIVE
  []
  [evm_fcn]
    type = PiecewiseConstant
    data_file = ss316_verification_data.csv
    x_index_in_file = 0
    y_index_in_file = 4
    format = columns
    xy_in_file_only = false
    direction = LEFT_INCLUSIVE
  []
  [temp_fcn]
    type = PiecewiseConstant
    data_file = ss316_verification_data.csv
    x_index_in_file = 0
    y_index_in_file = 5
    format = columns
    xy_in_file_only = false
    direction = LEFT_INCLUSIVE
  []

  [rhom_soln_fcn]
    type = PiecewiseConstant
    data_file = ss316_verification_data.csv
    x_index_in_file = 0
    y_index_in_file = 7
    format = columns
    xy_in_file_only = false
    direction = LEFT_INCLUSIVE
  []
  [rhoi_soln_fcn]
    type = PiecewiseConstant
    data_file = ss316_verification_data.csv
    x_index_in_file = 0
    y_index_in_file = 8
    format = columns
    xy_in_file_only = false
    direction = LEFT_INCLUSIVE
  []
  [creep_rate_soln_fcn]
    type = PiecewiseConstant
    data_file = ss316_verification_data.csv
    x_index_in_file = 0
    y_index_in_file = 10
    format = columns
    xy_in_file_only = false
    direction = LEFT_INCLUSIVE
  []
[]

[Modules/TensorMechanics/Master]
  [all]
    strain = FINITE
    add_variables = true
    generate_output = 'vonmises_stress'
    use_automatic_differentiation = true
  []
[]

[BCs]
  [symmx]
    type = ADDirichletBC
    variable = disp_x
    boundary = left
    value = 0
  []
  [symmy]
    type = ADDirichletBC
    variable = disp_y
    boundary = bottom
    value = 0
  []
  [symmz]
    type = ADDirichletBC
    variable = disp_z
    boundary = back
    value = 0
  []
  [pull_x]
    type = ADDirichletBC
    variable = disp_x
    boundary = right
    value = 1e-5 # This is required to make a non-zero effective trial stress so radial return is engaged
  []
[]

[Materials]
  [elasticity_tensor]
    type = ADComputeIsotropicElasticityTensor
    shear_modulus = 1e13
    poissons_ratio = 0.3
  []
  [stress]
    type = ADComputeMultipleInelasticStress
    inelastic_models = rom_stress_prediction
  []
  [rom_stress_prediction]
    type = ADSS316HLAROMANCEStressUpdateTest
    temperature = temperature
    effective_inelastic_strain_name = effective_creep_strain
    internal_solve_full_iteration_history = true
    apply_strain = false
    outputs = all
    wall_dislocation_density_forcing_function = rhoi_fcn
    cell_dislocation_density_forcing_function = rhom_fcn
    old_creep_strain_forcing_function = evm_fcn
    wall_input_window_low_failure = ERROR
    wall_input_window_high_failure = ERROR
    cell_input_window_low_failure = ERROR
    cell_input_window_high_failure = ERROR
    temperature_input_window_low_failure = ERROR
    temperature_input_window_high_failure = ERROR
    stress_input_window_low_failure = ERROR
    stress_input_window_high_failure = ERROR
    old_strain_input_window_low_failure = ERROR
    old_strain_input_window_high_failure = ERROR
    environment_input_window_low_failure = ERROR
    environment_input_window_high_failure = ERROR
    effective_stress_forcing_function = vmJ2_fcn
  []
[]

[Executioner]
  type = Transient

  solve_type = 'NEWTON'

  petsc_options_iname = '-pc_type'
  petsc_options_value = 'lu'
  nl_abs_tol = 1e-1 # Nothing is really being solved here, so loose tolerances are okay

  dt = 1e-3
  end_time = 1e-2
  timestep_tolerance = 1e-3
[]

[Postprocessors]
  [extrapolation]
    type = ElementAverageValue
    variable = ROM_extrapolation
    outputs = console
  []
  [old_strain_in]
    type = FunctionValuePostprocessor
    function = evm_fcn
    execute_on = 'TIMESTEP_END initial'
    outputs = console
  []
  [temperature]
    type = ElementAverageValue
    variable = temperature
    outputs = console
  []
  [rhom]
    type = ElementAverageValue
    variable = cell_dislocations
  []
  [rhoi]
    type = ElementAverageValue
    variable = wall_dislocations
  []
  [creep_rate]
    type = ElementAverageValue
    variable = creep_rate
  []
  [rhom_in]
    type = FunctionValuePostprocessor
    function = rhom_fcn
    execute_on = 'TIMESTEP_END initial'
    outputs = console
  []
  [rhoi_in]
    type = FunctionValuePostprocessor
    function = rhoi_fcn
    execute_on = 'TIMESTEP_END initial'
    outputs = console
  []
  [vmJ2_in]
    type = FunctionValuePostprocessor
    function = vmJ2_fcn
    execute_on = 'TIMESTEP_END initial'
    outputs = console
  []
  [rhom_soln]
    type = FunctionValuePostprocessor
    function = rhom_soln_fcn
    outputs = console
  []
  [rhoi_soln]
    type = FunctionValuePostprocessor
    function = rhoi_soln_fcn
    outputs = console
  []
  [creep_rate_soln]
    type = FunctionValuePostprocessor
    function = creep_rate_soln_fcn
  []

  [rhom_diff]
    type = ParsedPostprocessor
    pp_names = 'rhom_soln rhom'
    function = '(rhom_soln - rhom) / rhom_soln'
    outputs = console
  []
  [rhoi_diff]
    type = ParsedPostprocessor
    pp_names = 'rhoi_soln rhoi'
    function = '(rhoi_soln - rhoi) / rhoi_soln'
    outputs = console
  []
  [creep_rate_diff]
    type = ParsedPostprocessor
    pp_names = 'creep_rate creep_rate_soln'
    function = '(creep_rate_soln - creep_rate) / creep_rate_soln'
    outputs = console
  []

  [z_rhom_max_diff]
    type = TimeExtremeValue
    postprocessor = rhom_diff
    value_type = abs_max
  []
  [z_rhoi_max_diff]
    type = TimeExtremeValue
    postprocessor = rhoi_diff
    value_type = abs_max
  []
  [z_creep_rate_max_diff]
    type = TimeExtremeValue
    postprocessor = creep_rate_diff
    value_type = abs_max
  []
[]

[Outputs]
  csv = true
  execute_on = 'INITIAL TIMESTEP_END FINAL'
[]

(modules/heat_conduction/test/tests/gap_heat_transfer_mortar/large_gap_heat_transfer_test_cylinder_mortar_error.i)

rpv_core_gap_size = 0.15

core_outer_radius = 2
rpv_inner_radius = ${fparse 2 + rpv_core_gap_size}
rpv_outer_radius = ${fparse 2.5 + rpv_core_gap_size}

rpv_outer_htc = 10 # W/m^2/K
rpv_outer_Tinf = 300 # K

core_blocks = '1'
rpv_blocks = '3'

[Mesh]
  [core_gap_rpv]
    type = ConcentricCircleMeshGenerator
    num_sectors = 10
    radii = '${core_outer_radius} ${rpv_inner_radius} ${rpv_outer_radius}'
    rings = '2 1 2'
    has_outer_square = false
    preserve_volumes = true
    portion = full
  []
  [rename_core_bdy]
    type = SideSetsBetweenSubdomainsGenerator
    input = core_gap_rpv
    primary_block = 1
    paired_block = 2
    new_boundary = 'core_outer'
  []
  [rename_inner_rpv_bdy]
    type = SideSetsBetweenSubdomainsGenerator
    input = rename_core_bdy
    primary_block = 3
    paired_block = 2
    new_boundary = 'rpv_inner'
  []
  [2d_mesh]
    type = BlockDeletionGenerator
    input = rename_inner_rpv_bdy
    block = 2
  []
  [secondary]
    type = LowerDBlockFromSidesetGenerator
    sidesets = 'rpv_inner'
    new_block_id = 10001
    new_block_name = 'secondary_lower'
    input = 2d_mesh
  []
  [primary]
    type = LowerDBlockFromSidesetGenerator
    sidesets = 'core_outer'
    new_block_id = 10000
    new_block_name = 'primary_lower'
    input = secondary
  []
  allow_renumbering = false
[]

[Variables]
  [Tsolid]
    initial_condition = 500
  []
  [lm]
    order = FIRST
    family = LAGRANGE
    block = 'secondary_lower'
  []
[]

[Kernels]
  [heat_source]
    type = CoupledForce
    variable = Tsolid
    block = '${core_blocks}'
    v = power_density
  []
  [heat_conduction]
    type = HeatConduction
    variable = Tsolid
  []
[]

[BCs]
  [RPV_out_BC] # k \nabla T = h (T- T_inf) at RPV outer boundary
    type = ConvectiveFluxFunction # (Robin BC)
    variable = Tsolid
    boundary = 'outer' # outer RPV
    coefficient = ${rpv_outer_htc}
    T_infinity = ${rpv_outer_Tinf}
  []
[]

[UserObjects]
  [radiation]
    type = GapFluxModelRadiation
    temperature = Tsolid
    boundary = 'rpv_inner'
    primary_emissivity = 0.8
    secondary_emissivity = 0.8
  []
  [conduction]
    type = GapFluxModelConduction
    temperature = Tsolid
    boundary = 'rpv_inner'
    gap_conductivity = 0.1
  []
[]

[Constraints]
  [ced]
    type = ModularGapConductanceConstraint
    variable = lm
    secondary_variable = Tsolid
    primary_boundary = 'core_outer'
    primary_subdomain = 10000
    secondary_boundary = 'rpv_inner'
    secondary_subdomain = 10001
    gap_flux_models = 'radiation conduction'
    gap_geometry_type = 'CYLINDER'
    cylinder_axis_point_2 = '0 0 5'
  []
[]

[AuxVariables]
  [power_density]
    block = '${core_blocks}'
    initial_condition = 50e3
  []
[]

[Materials]
  [simple_mat]
    type = HeatConductionMaterial
    thermal_conductivity = 34.6 # W/m/K
  []
[]

[Postprocessors]
  [Tcore_avg]
    type = ElementAverageValue
    variable = Tsolid
    block = '${core_blocks}'
  []
  [Tcore_max]
    type = ElementExtremeValue
    value_type = max
    variable = Tsolid
    block = '${core_blocks}'
  []
  [Tcore_min]
    type = ElementExtremeValue
    value_type = min
    variable = Tsolid
    block = '${core_blocks}'
  []
  [Trpv_avg]
    type = ElementAverageValue
    variable = Tsolid
    block = '${rpv_blocks}'
  []
  [Trpv_max]
    type = ElementExtremeValue
    value_type = max
    variable = Tsolid
    block = '${rpv_blocks}'
  []
  [Trpv_min]
    type = ElementExtremeValue
    value_type = min
    variable = Tsolid
    block = '${rpv_blocks}'
  []
  [ptot]
    type = ElementIntegralVariablePostprocessor
    variable = power_density
    block = '${core_blocks}'
  []
  [rpv_convective_out]
    type = ConvectiveHeatTransferSideIntegral
    T_solid = Tsolid
    boundary = 'outer' # outer RVP
    T_fluid = ${rpv_outer_Tinf}
    htc = ${rpv_outer_htc}
  []
  [heat_balance] # should be equal to 0 upon convergence
    type = ParsedPostprocessor
    function = '(rpv_convective_out - ptot) / ptot'
    pp_names = 'rpv_convective_out ptot'
  []
[]

[Preconditioning]
  [smp]
    type = SMP
    full = true
  []
[]

[VectorPostprocessors]
  [NodalTemperature]
    type = NodalValueSampler
    sort_by = id
    boundary = 'rpv_inner core_outer'
    variable = 'Tsolid'
  []
[]

[Executioner]
  type = Steady

  petsc_options = '-snes_converged_reason -pc_svd_monitor'
  petsc_options_iname = '-pc_type -pc_factor_mat_solver_package -mat_mffd_err -pc_factor_shift_type '
                        '-pc_factor_shift_amount'
  petsc_options_value = ' lu       superlu_dist                  1e-5          NONZERO               '
                        '1e-15'
  snesmf_reuse_base = false

  nl_rel_tol = 1e-10
  nl_abs_tol = 1e-10
  l_max_its = 100

  line_search = none
[]

[Outputs]
  exodus = false
  csv = true
[]

(modules/heat_conduction/test/tests/gap_heat_transfer_mortar/gap_heat_transfer_3D.i)

outer_htc = 10 # W/m^2/K
outer_Tinf = 300 # K

[GlobalParams]
  order = SECOND
  family = LAGRANGE
[]

[Mesh]
  [left_block]
    type = GeneratedMeshGenerator
    dim = 3
    nx = 3
    ny = 6
    nz = 6
    xmin = -1
    xmax = -0.5
    ymin = -0.5
    ymax = 0.5
    zmin = -0.5
    zmax = 0.5
    elem_type = HEX27
  []
  [left_block_sidesets]
    type = RenameBoundaryGenerator
    input = left_block
    old_boundary = '0 1 2 3 4 5'
    new_boundary = 'left_bottom left_back left_right left_front left_left left_top'
  []
  [left_block_id]
    type = SubdomainIDGenerator
    input = left_block_sidesets
    subdomain_id = 1
  []

  [right_block]
    type = GeneratedMeshGenerator
    dim = 3
    nx = 4
    ny = 8
    nz = 8
    xmin = 0.5
    xmax = 1
    ymin = -0.5
    ymax = 0.5
    zmin = -0.5
    zmax = 0.5
    elem_type = HEX27
  []
  [right_block_sidesets]
    type = RenameBoundaryGenerator
    input = right_block
    old_boundary = '0 1 2 3 4 5'
    # new_boundary = 'right_bottom right_back right_right right_front right_left right_top'
    new_boundary = '100 101 102 103 104 105'
  []
  [right_block_sidesets_rename]
    type = RenameBoundaryGenerator
    input = right_block_sidesets
    old_boundary = '100 101 102 103 104 105'
    new_boundary = 'right_bottom right_back right_right right_front right_left right_top'
  []
  [right_block_id]
    type = SubdomainIDGenerator
    input = right_block_sidesets_rename
    subdomain_id = 2
  []

  [combined_mesh]
    type = MeshCollectionGenerator
    inputs = 'left_block_id right_block_id'
  []
[]

[Functions]
  [temp]
    type = PiecewiseLinear
    x = '0   1'
    y = '100 200'
  []
[]

[Variables]
  [temp]
    initial_condition = 500
  []
[]

[AuxVariables]
  [gap_conductance]
    order = CONSTANT
    family = MONOMIAL
  []
  [power_density]
    block = 1
    initial_condition = 50e3
  []
[]

[Kernels]
  [heat_conduction]
    type = HeatConduction
    variable = temp
  []
  [heat_source]
    type = CoupledForce
    variable = temp
    block = 1
    v = power_density
  []
[]

[AuxKernels]
  [gap_cond]
    type = MaterialRealAux
    property = gap_conductance
    variable = gap_conductance
    boundary = 'left_right'
  []
[]

[Materials]
  [heat1]
    type = HeatConductionMaterial
    block = '1 2'
    specific_heat = 1.0
    thermal_conductivity = 34.6
  []
[]

[ThermalContact]
  [thermal_contact]
    type = GapHeatTransfer
    variable = temp
    primary = 'right_left'
    secondary = 'left_right'
    emissivity_primary = 0
    emissivity_secondary = 0
    gap_conductivity = 5
    gap_geometry_type = PLATE
  []
[]

[BCs]
  [RPV_out_BC] # k \nabla T = h (T- T_inf) at RPV outer boundary
    type = ConvectiveFluxFunction # (Robin BC)
    variable = temp
    boundary = 'right_right' # outer RPV
    coefficient = ${outer_htc}
    T_infinity = ${outer_Tinf}
  []
[]

[Executioner]
  type = Transient
  solve_type = 'PJFNK'

  petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
  petsc_options_value = 'lu       superlu_dist'

  dt = 1
  dtmin = 0.01
  end_time = 1

  nl_rel_tol = 1e-12
  nl_abs_tol = 1e-8

  [Quadrature]
    order = fifth
    side_order = seventh
  []
[]

[Outputs]
  exodus = true
  csv = true
  [Console]
    type = Console
  []
[]

[Postprocessors]
  [temp_left]
    type = SideAverageValue
    boundary = 'left_right'
    variable = temp
  []

  [temp_right]
    type = SideAverageValue
    boundary = 'right_left'
    variable = temp
  []

  [flux_left]
    type = SideDiffusiveFluxIntegral
    variable = temp
    boundary = 'left_right'
    diffusivity = thermal_conductivity
  []
  [flux_right]
    type = SideDiffusiveFluxIntegral
    variable = temp
    boundary = 'right_left'
    diffusivity = thermal_conductivity
  []
  [ptot]
    type = ElementIntegralVariablePostprocessor
    variable = power_density
    block = 1
  []
  [convective_out]
    type = ConvectiveHeatTransferSideIntegral
    T_solid = temp
    boundary = 'right_right' # outer RVP
    T_fluid = ${outer_Tinf}
    htc = ${outer_htc}
  []
  [heat_balance] # should be equal to 0 upon convergence
    type = ParsedPostprocessor
    function = '(convective_out - ptot) / ptot'
    pp_names = 'convective_out ptot'
  []
[]

[VectorPostprocessors]
  [NodalTemperature]
    type = NodalValueSampler
    sort_by = id
    boundary = 'left_right right_left'
    variable = temp
  []
[]

(modules/tensor_mechanics/test/tests/rom_stress_update/nonad_verification.i)

[Mesh]
  type = GeneratedMesh
  dim = 3
[]

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

[AuxVariables]
  [temperature]
  []
[]

[AuxKernels]
  [temp_aux]
    type = FunctionAux
    variable = temperature
    function = temp_fcn
    execute_on = 'initial timestep_begin'
  []
[]

[Functions]
  [rhom_fcn]
    type = PiecewiseConstant
    data_file = ss316_verification_data.csv
    x_index_in_file = 0
    y_index_in_file = 1
    format = columns
    xy_in_file_only = false
    direction = LEFT_INCLUSIVE
  []
  [rhoi_fcn]
    type = PiecewiseConstant
    data_file = ss316_verification_data.csv
    x_index_in_file = 0
    y_index_in_file = 2
    format = columns
    xy_in_file_only = false
    direction = LEFT_INCLUSIVE
  []
  [vmJ2_fcn]
    type = PiecewiseConstant
    data_file = ss316_verification_data.csv
    x_index_in_file = 0
    y_index_in_file = 3
    format = columns
    xy_in_file_only = false
    direction = LEFT_INCLUSIVE
  []
  [evm_fcn]
    type = PiecewiseConstant
    data_file = ss316_verification_data.csv
    x_index_in_file = 0
    y_index_in_file = 4
    format = columns
    xy_in_file_only = false
    direction = LEFT_INCLUSIVE
  []
  [temp_fcn]
    type = PiecewiseConstant
    data_file = ss316_verification_data.csv
    x_index_in_file = 0
    y_index_in_file = 5
    format = columns
    xy_in_file_only = false
    direction = LEFT_INCLUSIVE
  []

  [rhom_soln_fcn]
    type = PiecewiseConstant
    data_file = ss316_verification_data.csv
    x_index_in_file = 0
    y_index_in_file = 7
    format = columns
    xy_in_file_only = false
    direction = LEFT_INCLUSIVE
  []
  [rhoi_soln_fcn]
    type = PiecewiseConstant
    data_file = ss316_verification_data.csv
    x_index_in_file = 0
    y_index_in_file = 8
    format = columns
    xy_in_file_only = false
    direction = LEFT_INCLUSIVE
  []
  [creep_rate_soln_fcn]
    type = PiecewiseConstant
    data_file = ss316_verification_data.csv
    x_index_in_file = 0
    y_index_in_file = 10
    format = columns
    xy_in_file_only = false
    direction = LEFT_INCLUSIVE
  []
[]

[Modules/TensorMechanics/Master]
  [all]
    strain = FINITE
    add_variables = true
    generate_output = 'vonmises_stress'
  []
[]

[BCs]
  [symmx]
    type = DirichletBC
    variable = disp_x
    boundary = left
    value = 0
  []
  [symmy]
    type = DirichletBC
    variable = disp_y
    boundary = bottom
    value = 0
  []
  [symmz]
    type = DirichletBC
    variable = disp_z
    boundary = back
    value = 0
  []
  [pull_x]
    type = DirichletBC
    variable = disp_x
    boundary = right
    value = 1e-5 # This is required to make a non-zero effective trial stress so radial return is engaged
  []
[]

[Materials]
  [elasticity_tensor]
    type = ComputeIsotropicElasticityTensor
    shear_modulus = 1e13
    poissons_ratio = 0.3
  []
  [stress]
    type = ComputeMultipleInelasticStress
    inelastic_models = rom_stress_prediction
  []
  [rom_stress_prediction]
    type =  SS316HLAROMANCEStressUpdateTest
    temperature = temperature
    effective_inelastic_strain_name = effective_creep_strain
    internal_solve_full_iteration_history = true
    apply_strain = false
    outputs = all
    wall_dislocation_density_forcing_function = rhoi_fcn
    cell_dislocation_density_forcing_function = rhom_fcn
    old_creep_strain_forcing_function = evm_fcn
    wall_input_window_low_failure = ERROR
    wall_input_window_high_failure = ERROR
    cell_input_window_low_failure = ERROR
    cell_input_window_high_failure = ERROR
    temperature_input_window_low_failure = ERROR
    temperature_input_window_high_failure = ERROR
    stress_input_window_low_failure = ERROR
    stress_input_window_high_failure = ERROR
    old_strain_input_window_low_failure = ERROR
    old_strain_input_window_high_failure = ERROR
    environment_input_window_low_failure = ERROR
    environment_input_window_high_failure = ERROR
    effective_stress_forcing_function = vmJ2_fcn
  []
[]

[Executioner]
  type = Transient

  solve_type = 'NEWTON'

  petsc_options_iname = '-pc_type'
  petsc_options_value = 'lu'
  nl_abs_tol = 1e-1 # Nothing is really being solved here, so loose tolerances are okay

  dt = 1e-3
  end_time = 1e-2
  timestep_tolerance = 1e-3
[]

[Postprocessors]
  [extrapolation]
    type = ElementAverageValue
    variable = ROM_extrapolation
    outputs = console
  []
  [old_strain_in]
    type = FunctionValuePostprocessor
    function = evm_fcn
    execute_on = 'TIMESTEP_END initial'
    outputs = console
  []
  [temperature]
    type = ElementAverageValue
    variable = temperature
    outputs = console
  []
  [rhom]
    type = ElementAverageValue
    variable = cell_dislocations
  []
  [rhoi]
    type = ElementAverageValue
    variable = wall_dislocations
  []
  [creep_rate]
    type = ElementAverageValue
    variable = creep_rate
  []
  [rhom_in]
    type = FunctionValuePostprocessor
    function = rhom_fcn
    execute_on = 'TIMESTEP_END initial'
    outputs = console
  []
  [rhoi_in]
    type = FunctionValuePostprocessor
    function = rhoi_fcn
    execute_on = 'TIMESTEP_END initial'
    outputs = console
  []
  [vmJ2_in]
    type = FunctionValuePostprocessor
    function = vmJ2_fcn
    execute_on = 'TIMESTEP_END initial'
    outputs = console
  []
  [rhom_soln]
    type = FunctionValuePostprocessor
    function = rhom_soln_fcn
    outputs = console
  []
  [rhoi_soln]
    type = FunctionValuePostprocessor
    function = rhoi_soln_fcn
    outputs = console
  []
  [creep_rate_soln]
    type = FunctionValuePostprocessor
    function = creep_rate_soln_fcn
  []

  [rhom_diff]
    type = ParsedPostprocessor
    pp_names = 'rhom_soln rhom'
    function = '(rhom_soln - rhom) / rhom_soln'
    outputs = console
  []
  [rhoi_diff]
    type = ParsedPostprocessor
    pp_names = 'rhoi_soln rhoi'
    function = '(rhoi_soln - rhoi) / rhoi_soln'
    outputs = console
  []
  [creep_rate_diff]
    type = ParsedPostprocessor
    pp_names = 'creep_rate creep_rate_soln'
    function = '(creep_rate_soln - creep_rate) / creep_rate_soln'
    outputs = console
  []

  [z_rhom_max_diff]
    type = TimeExtremeValue
    postprocessor = rhom_diff
    value_type = abs_max
  []
  [z_rhoi_max_diff]
    type = TimeExtremeValue
    postprocessor = rhoi_diff
    value_type = abs_max
  []
  [z_creep_rate_max_diff]
    type = TimeExtremeValue
    postprocessor = creep_rate_diff
    value_type = abs_max
  []
[]

[Outputs]
  csv = true
  execute_on = 'INITIAL TIMESTEP_END FINAL'
[]

(modules/heat_conduction/test/tests/gap_heat_transfer_balance/large_gap_heat_transfer_test_sphere.i)

sphere_outer_htc = 10 # W/m^2/K
sphere_outer_Tinf = 300 # K

[GlobalParams]
  order = SECOND
  family = LAGRANGE
[]

[Problem]
  coord_type = RZ
[]

[Mesh]
  [file]
    type = FileMeshGenerator
    file = cyl2D.e
  []
[]

[Functions]
  [temp]
    type = PiecewiseLinear
    x = '0   1'
    y = '100 200'
  []
[]

[Variables]
  [temp]
    initial_condition = 500
  []
[]

[AuxVariables]
  [gap_conductance]
    order = CONSTANT
    family = MONOMIAL
  []
  [power_density]
    block = 'fuel'
    initial_condition = 50e3
  []
[]

[Kernels]
  [heat_conduction]
    type = HeatConduction
    variable = temp
  []
  [heat_source]
    type = CoupledForce
    variable = temp
    block = 'fuel'
    v = power_density
  []
[]

[AuxKernels]
  [gap_cond]
    type = MaterialRealAux
    property = gap_conductance
    variable = gap_conductance
    boundary = 2
  []
[]

[Materials]
  [heat1]
    type = HeatConductionMaterial
    block = '1 2'
    specific_heat = 1.0
    thermal_conductivity = 34.6
  []
[]

[ThermalContact]
  [thermal_contact]
    type = GapHeatTransfer
    variable = temp
    primary = 3
    secondary = 2
    emissivity_primary = 0.8
    emissivity_secondary = 0.8
    gap_conductivity = 0.1
    quadrature = true
    gap_geometry_type = SPHERE
    sphere_origin = '0 0 0'
  []
[]

[BCs]
  [RPV_out_BC] # k \nabla T = h (T- T_inf) at RPV outer boundary
    type = ConvectiveFluxFunction # (Robin BC)
    variable = temp
    boundary = '4' # outer RPV
    coefficient = ${sphere_outer_htc}
    T_infinity = ${sphere_outer_Tinf}
  []
[]

[Executioner]
  type = Transient
  solve_type = 'PJFNK'

  petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
  petsc_options_value = 'lu       superlu_dist'

  dt = 1
  dtmin = 0.01
  end_time = 1

  nl_rel_tol = 1e-12
  nl_abs_tol = 1e-7

  [Quadrature]
    order = fifth
    side_order = seventh
  []
[]

[Outputs]
  exodus = false
  csv = true
  [Console]
    type = Console
  []
[]

[Postprocessors]
  [temp_left]
    type = SideAverageValue
    boundary = 2
    variable = temp
  []

  [temp_right]
    type = SideAverageValue
    boundary = 3
    variable = temp
  []

  [flux_left]
    type = SideDiffusiveFluxIntegral
    variable = temp
    boundary = 2
    diffusivity = thermal_conductivity
  []

  [flux_right]
    type = SideDiffusiveFluxIntegral
    variable = temp
    boundary = 3
    diffusivity = thermal_conductivity
  []
  [ptot]
    type = ElementIntegralVariablePostprocessor
    variable = power_density
    block = 'fuel'
  []
  [sphere_convective_out]
    type = ConvectiveHeatTransferSideIntegral
    T_solid = temp
    boundary = '4' # outer RVP
    T_fluid = ${sphere_outer_Tinf}
    htc = ${sphere_outer_htc}
  []
  [heat_balance] # should be equal to 0 upon convergence
    type = ParsedPostprocessor
    function = '(sphere_convective_out - ptot) / ptot'
    pp_names = 'sphere_convective_out ptot'
  []
[]

(modules/heat_conduction/test/tests/gap_heat_transfer_mortar/large_gap_heat_transfer_test_sphere.i)

sphere_outer_htc = 10 # W/m^2/K
sphere_outer_Tinf = 300 # K

[GlobalParams]
  order = SECOND
  family = LAGRANGE
[]

[Problem]
  coord_type = RZ
[]

[Mesh]
  [file]
    type = FileMeshGenerator
    file = cyl2D.e
  []
  allow_renumbering = false
[]

[Functions]
  [temp]
    type = PiecewiseLinear
    x = '0   1'
    y = '100 200'
  []
[]

[Variables]
  [temp]
    initial_condition = 500
  []
[]

[AuxVariables]
  [gap_conductance]
    order = CONSTANT
    family = MONOMIAL
  []
  [power_density]
    block = 'fuel'
    initial_condition = 50e3
  []
[]

[Kernels]
  [heat_conduction]
    type = HeatConduction
    variable = temp
  []
  [heat_source]
    type = CoupledForce
    variable = temp
    block = 'fuel'
    v = power_density
  []
[]

[AuxKernels]
  [gap_cond]
    type = MaterialRealAux
    property = gap_conductance
    variable = gap_conductance
    boundary = 2
  []
[]

[Materials]
  [heat1]
    type = HeatConductionMaterial
    block = '1 2'
    specific_heat = 1.0
    thermal_conductivity = 34.6
  []
[]

[ThermalContact]
  [thermal_contact]
    type = GapHeatTransfer
    variable = temp
    primary = 3
    secondary = 2
    emissivity_primary = 0.0
    emissivity_secondary = 0.0
    gap_conductivity = 5
  #  quadrature = true
    gap_geometry_type = SPHERE
    sphere_origin = '0 0 0'
  []
[]

[BCs]
  [RPV_out_BC] # k \nabla T = h (T- T_inf) at RPV outer boundary
    type = ConvectiveFluxFunction # (Robin BC)
    variable = temp
    boundary = '4' # outer RPV
    coefficient = ${sphere_outer_htc}
    T_infinity = ${sphere_outer_Tinf}
  []
[]

[VectorPostprocessors]
  [NodalTemperature]
    type = NodalValueSampler
    sort_by = id
    boundary = '2 3'
    variable = temp
  []
[]

[Executioner]
  type = Transient
  solve_type = 'PJFNK'

  petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
  petsc_options_value = 'lu       superlu_dist'

  dt = 1
  dtmin = 0.01
  end_time = 1

  nl_rel_tol = 1e-12
  nl_abs_tol = 1e-7

  [Quadrature]
    order = fifth
    side_order = seventh
  []
[]

[Outputs]
  exodus = true
  csv = true
  [Console]
    type = Console
  []
[]

[Postprocessors]
  [temp_left]
    type = SideAverageValue
    boundary = 2
    variable = temp
  []

  [temp_right]
    type = SideAverageValue
    boundary = 3
    variable = temp
  []

  [flux_left]
    type = SideDiffusiveFluxIntegral
    variable = temp
    boundary = 2
    diffusivity = thermal_conductivity
  []

  [flux_right]
    type = SideDiffusiveFluxIntegral
    variable = temp
    boundary = 3
    diffusivity = thermal_conductivity
  []
  [ptot]
    type = ElementIntegralVariablePostprocessor
    variable = power_density
    block = 'fuel'
  []
  [sphere_convective_out]
    type = ConvectiveHeatTransferSideIntegral
    T_solid = temp
    boundary = '4' # outer RVP
    T_fluid = ${sphere_outer_Tinf}
    htc = ${sphere_outer_htc}
  []
  [heat_balance] # should be equal to 0 upon convergence
    type = ParsedPostprocessor
    function = '(sphere_convective_out - ptot) / ptot'
    pp_names = 'sphere_convective_out ptot'
  []
[]

(modules/heat_conduction/test/tests/gap_heat_transfer_mortar/large_gap_heat_transfer_test_cylinder_mortar.i)

rpv_core_gap_size = 0.15

core_outer_radius = 2
rpv_inner_radius = ${fparse 2 + rpv_core_gap_size}
rpv_outer_radius = ${fparse 2.5 + rpv_core_gap_size}

rpv_outer_htc = 10 # W/m^2/K
rpv_outer_Tinf = 300 # K

core_blocks = '1'
rpv_blocks = '3'

[Mesh]
  [core_gap_rpv]
    type = ConcentricCircleMeshGenerator
    num_sectors = 10
    radii = '${core_outer_radius} ${rpv_inner_radius} ${rpv_outer_radius}'
    rings = '2 1 2'
    has_outer_square = false
    preserve_volumes = true
    portion = full
  []
  [rename_core_bdy]
    type = SideSetsBetweenSubdomainsGenerator
    input = core_gap_rpv
    primary_block = 1
    paired_block = 2
    new_boundary = 'core_outer'
  []
  [rename_inner_rpv_bdy]
    type = SideSetsBetweenSubdomainsGenerator
    input = rename_core_bdy
    primary_block = 3
    paired_block = 2
    new_boundary = 'rpv_inner'
  []
  [2d_mesh]
    type = BlockDeletionGenerator
    input = rename_inner_rpv_bdy
    block = 2
  []
  [secondary]
    type = LowerDBlockFromSidesetGenerator
    sidesets = 'rpv_inner'
    new_block_id = 10001
    new_block_name = 'secondary_lower'
    input = 2d_mesh
  []
  [primary]
    type = LowerDBlockFromSidesetGenerator
    sidesets = 'core_outer'
    new_block_id = 10000
    new_block_name = 'primary_lower'
    input = secondary
  []
  allow_renumbering = false
[]

[Variables]
  [Tsolid]
    initial_condition = 500
  []
  [lm]
    order = FIRST
    family = LAGRANGE
    block = 'secondary_lower'
  []
[]

[Kernels]
  [heat_source]
    type = CoupledForce
    variable = Tsolid
    block = '${core_blocks}'
    v = power_density
  []
  [heat_conduction]
    type = HeatConduction
    variable = Tsolid
  []
[]

[BCs]
  [RPV_out_BC] # k \nabla T = h (T- T_inf) at RPV outer boundary
    type = ConvectiveFluxFunction # (Robin BC)
    variable = Tsolid
    boundary = 'outer' # outer RPV
    coefficient = ${rpv_outer_htc}
    T_infinity = ${rpv_outer_Tinf}
  []
[]

[UserObjects]
  [radiation]
    type = GapFluxModelRadiation
    temperature = Tsolid
    boundary = 'rpv_inner'
    primary_emissivity = 0.8
    secondary_emissivity = 0.8
  []
  [conduction]
    type = GapFluxModelConduction
    temperature = Tsolid
    boundary = 'rpv_inner'
    gap_conductivity = 0.1
  []
[]

[Constraints]
  [ced]
    type = ModularGapConductanceConstraint
    variable = lm
    secondary_variable = Tsolid
    primary_boundary = 'core_outer'
    primary_subdomain = 10000
    secondary_boundary = 'rpv_inner'
    secondary_subdomain = 10001
    gap_flux_models = 'radiation conduction'
    gap_geometry_type = 'CYLINDER'
    cylinder_axis_point_1 = '0 0 0'
    cylinder_axis_point_2 = '0 0 5'
  []
[]

[AuxVariables]
  [power_density]
    block = '${core_blocks}'
    initial_condition = 50e3
  []
[]

[Materials]
  [simple_mat]
    type = HeatConductionMaterial
    thermal_conductivity = 34.6 # W/m/K
  []
[]

[Postprocessors]
  [Tcore_avg]
    type = ElementAverageValue
    variable = Tsolid
    block = '${core_blocks}'
  []
  [Tcore_max]
    type = ElementExtremeValue
    value_type = max
    variable = Tsolid
    block = '${core_blocks}'
  []
  [Tcore_min]
    type = ElementExtremeValue
    value_type = min
    variable = Tsolid
    block = '${core_blocks}'
  []
  [Trpv_avg]
    type = ElementAverageValue
    variable = Tsolid
    block = '${rpv_blocks}'
  []
  [Trpv_max]
    type = ElementExtremeValue
    value_type = max
    variable = Tsolid
    block = '${rpv_blocks}'
  []
  [Trpv_min]
    type = ElementExtremeValue
    value_type = min
    variable = Tsolid
    block = '${rpv_blocks}'
  []
  [ptot]
    type = ElementIntegralVariablePostprocessor
    variable = power_density
    block = '${core_blocks}'
  []
  [rpv_convective_out]
    type = ConvectiveHeatTransferSideIntegral
    T_solid = Tsolid
    boundary = 'outer' # outer RVP
    T_fluid = ${rpv_outer_Tinf}
    htc = ${rpv_outer_htc}
  []
  [heat_balance] # should be equal to 0 upon convergence
    type = ParsedPostprocessor
    function = '(rpv_convective_out - ptot) / ptot'
    pp_names = 'rpv_convective_out ptot'
  []
[]

[Preconditioning]
  [smp]
    type = SMP
    full = true
  []
[]

[VectorPostprocessors]
  [NodalTemperature]
    type = NodalValueSampler
    sort_by = id
    boundary = 'rpv_inner core_outer'
    variable = 'Tsolid'
  []
[]

[Executioner]
  type = Steady

  petsc_options = '-snes_converged_reason -pc_svd_monitor'
  petsc_options_iname = '-pc_type -pc_factor_mat_solver_package -mat_mffd_err -pc_factor_shift_type '
                        '-pc_factor_shift_amount'
  petsc_options_value = ' lu       superlu_dist                  1e-5          NONZERO               '
                        '1e-15'
  snesmf_reuse_base = false

  nl_rel_tol = 1e-10
  nl_abs_tol = 1e-10
  l_max_its = 100

  line_search = none
[]

[Outputs]
  exodus = false
  csv = true
[]

(modules/heat_conduction/test/tests/gap_heat_transfer_mortar/large_gap_heat_transfer_test_rz_cylinder_mortar.i)

rpv_core_gap_size = 0.2

core_outer_radius = 2
rpv_inner_radius = '${fparse 2 + rpv_core_gap_size}'
rpv_outer_radius = '${fparse 2.5 + rpv_core_gap_size}'
rpv_width = '${fparse rpv_outer_radius - rpv_inner_radius}'

rpv_outer_htc = 10 # W/m^2/K
rpv_outer_Tinf = 300 # K

core_blocks = '1'
rpv_blocks = '3'

[Mesh]
  [gmg]
    type = CartesianMeshGenerator
    dim = 2
    dx = '${core_outer_radius} ${rpv_core_gap_size} ${rpv_width}'
    ix = '400 1 100'
    dy = 1
    iy = '5'
  []
  [set_block_id1]
    type = SubdomainBoundingBoxGenerator
    input = gmg
    bottom_left = '0 0 0'
    top_right = '${core_outer_radius} 1 0'
    block_id = 1
    location = INSIDE
  []
  [rename_core_bdy]
    type = SideSetsBetweenSubdomainsGenerator
    input = set_block_id1
    primary_block = 1
    paired_block = 0
    new_boundary = 'core_outer'
  []
  [set_block_id3]
    type = SubdomainBoundingBoxGenerator
    input = rename_core_bdy
    bottom_left = '${rpv_inner_radius} 0 0'
    top_right = '${rpv_outer_radius} 1 0'
    block_id = 3
    location = INSIDE
  []
  [rename_inner_rpv_bdy]
    type = SideSetsBetweenSubdomainsGenerator
    input = set_block_id3
    primary_block = 3
    paired_block = 0
    new_boundary = 'rpv_inner'
  []
  # comment out for test without gap
  [2d_mesh]
    type = BlockDeletionGenerator
    input = rename_inner_rpv_bdy
    block = 0
  []
  [secondary]
    type = LowerDBlockFromSidesetGenerator
    sidesets = 'rpv_inner'
    new_block_id = 10001
    new_block_name = 'secondary_lower'
    input = 2d_mesh
  []
  [primary]
    type = LowerDBlockFromSidesetGenerator
    sidesets = 'core_outer'
    new_block_id = 10000
    new_block_name = 'primary_lower'
    input = secondary
  []
  allow_renumbering = false
[]

[Problem]
  coord_type = RZ
[]

[Variables]
  [Tsolid]
    initial_condition = 500
  []
    [lm]
      order = FIRST
      family = LAGRANGE
      block = 'secondary_lower'
    []
[]

[Kernels]
  [heat_source]
    type = CoupledForce
    variable = Tsolid
    block = '${core_blocks}'
    v = power_density
  []
  [heat_conduction]
    type = HeatConduction
    variable = Tsolid
  []
[]

[BCs]
  [RPV_out_BC] # k \nabla T = h (T- T_inf) at RPV outer boundary
    type = ConvectiveFluxFunction # (Robin BC)
    variable = Tsolid
    boundary = 'right' # outer RPV
    coefficient = ${rpv_outer_htc}
    T_infinity = ${rpv_outer_Tinf}
  []
[]

[UserObjects]
  [radiation]
    type = GapFluxModelRadiation
    temperature = Tsolid
    boundary = 'rpv_inner'
    primary_emissivity = 0.8
    secondary_emissivity = 0.8
  []
  [conduction]
    type = GapFluxModelConduction
    temperature = Tsolid
    boundary = 'rpv_inner'
    gap_conductivity = 0.1
  []
[]

[Constraints]
  [ced]
    type = ModularGapConductanceConstraint
    variable = lm
    secondary_variable = Tsolid
    primary_boundary = 'core_outer'
    primary_subdomain = 10000
    secondary_boundary = 'rpv_inner'
    secondary_subdomain = 10001
    gap_flux_models = 'radiation conduction'
    gap_geometry_type = 'CYLINDER'
  []
[]

[AuxVariables]
  [power_density]
    block = '${core_blocks}'
    initial_condition = 50e3
  []
[]

[Materials]
  [simple_mat]
    type = HeatConductionMaterial
    thermal_conductivity = 34.6 # W/m/K
  []
[]

[Postprocessors]
  [Tcore_avg]
    type = ElementAverageValue
    variable = Tsolid
    block = '${core_blocks}'
  []
  [Tcore_max]
    type = ElementExtremeValue
    value_type = max
    variable = Tsolid
    block = '${core_blocks}'
  []
  [Tcore_min]
    type = ElementExtremeValue
    value_type = min
    variable = Tsolid
    block = '${core_blocks}'
  []
  [Trpv_avg]
    type = ElementAverageValue
    variable = Tsolid
    block = '${rpv_blocks}'
  []
  [Trpv_max]
    type = ElementExtremeValue
    value_type = max
    variable = Tsolid
    block = '${rpv_blocks}'
  []
  [Trpv_min]
    type = ElementExtremeValue
    value_type = min
    variable = Tsolid
    block = '${rpv_blocks}'
  []
  [ptot]
    type = ElementIntegralVariablePostprocessor
    variable = power_density
    block = '${core_blocks}'
  []
  [rpv_convective_out]
    type = ConvectiveHeatTransferSideIntegral
    T_solid = Tsolid
    boundary = 'right' # outer RVP
    T_fluid = ${rpv_outer_Tinf}
    htc = ${rpv_outer_htc}
  []
  [heat_balance] # should be equal to 0 upon convergence
    type = ParsedPostprocessor
    function = '(rpv_convective_out - ptot) / ptot'
    pp_names = 'rpv_convective_out ptot'
  []
  [flux_from_core] # converges to ptot as the mesh is refined
    type = SideDiffusiveFluxIntegral
    variable = Tsolid
    boundary = core_outer
    diffusivity = thermal_conductivity
  []
  [flux_into_rpv] # converges to rpv_convective_out as the mesh is refined
    type = SideDiffusiveFluxIntegral
    variable = Tsolid
    boundary = rpv_inner
    diffusivity = thermal_conductivity
  []

[]

[Preconditioning]
  [smp]
    type = SMP
    full = true
  []
[]

[VectorPostprocessors]
  [NodalTemperature]
    type = NodalValueSampler
    sort_by = id
    boundary = 'rpv_inner core_outer'
    variable = Tsolid
  []
[]

[Executioner]
  type = Steady

  petsc_options = '-snes_converged_reason -pc_svd_monitor'
  petsc_options_iname = '-pc_type -pc_factor_mat_solver_package -mat_mffd_err -pc_factor_shift_type '
                        '-pc_factor_shift_amount'
  petsc_options_value = ' lu       superlu_dist                  1e-5          NONZERO               '
                        '1e-15'
  snesmf_reuse_base = false
  nl_rel_tol = 1e-10
  nl_abs_tol = 1e-10
  l_max_its = 100

  line_search = none
[]

[Outputs]
  exodus = false
  csv = true
[]

(modules/heat_conduction/test/tests/gap_heat_transfer_mortar/large_gap_heat_transfer_test_sphere_mortar.i)

sphere_outer_htc = 10 # W/m^2/K
sphere_outer_Tinf = 300 # K

[GlobalParams]
  order = SECOND
  family = LAGRANGE
[]

[Problem]
  coord_type = RZ
  kernel_coverage_check = false
  material_coverage_check = false
[]

[Mesh]
  [file]
    type = FileMeshGenerator
    file = cyl2D.e
  []
  [secondary]
    type = LowerDBlockFromSidesetGenerator
    sidesets = '2'
    new_block_id = 10001
    new_block_name = 'secondary_lower'
    input = file
  []
  [primary]
    type = LowerDBlockFromSidesetGenerator
    sidesets = '3'
    new_block_id = 10000
    new_block_name = 'primary_lower'
    input = secondary
  []
  allow_renumbering = false
[]

[Functions]
  [temp]
    type = PiecewiseLinear
    x = '0   1'
    y = '100 200'
  []
[]

[Variables]
  [temp]
    initial_condition = 500
  []
  [lm]
    order = SECOND
    family = LAGRANGE
    block = 'secondary_lower'
  []
[]

[AuxVariables]
  [power_density]
    block = 'fuel'
    initial_condition = 50e3
  []
[]

[Kernels]
  [heat_conduction]
    type = HeatConduction
    variable = temp
    block = '1 2'
  []
  [heat_source]
    type = CoupledForce
    variable = temp
    block = 'fuel'
    v = power_density
  []
[]

[Materials]
  [heat1]
    type = HeatConductionMaterial
    block = '1 2'
    specific_heat = 1.0
    thermal_conductivity = 34.6
  []
[]

[UserObjects]
  [radiation]
    type = GapFluxModelRadiation
    temperature = temp
    boundary = 2
    primary_emissivity = 0.0
    secondary_emissivity = 0.0
  []
  [conduction]
    type = GapFluxModelConduction
    temperature = temp
    boundary = 2
    gap_conductivity = 5.0
  []
[]

[Constraints]
  [ced]
    type = ModularGapConductanceConstraint
    variable = lm
    secondary_variable = temp
    primary_boundary = 3
    primary_subdomain = 10000
    secondary_boundary = 2
    secondary_subdomain = 10001
    gap_flux_models = 'radiation conduction'
    gap_geometry_type = SPHERE
    sphere_origin = '0 0 0'
  []
[]

[BCs]
  [RPV_out_BC] # k \nabla T = h (T- T_inf) at RPV outer boundary
    type = ConvectiveFluxFunction # (Robin BC)
    variable = temp
    boundary = '4' # outer RPV
    coefficient = ${sphere_outer_htc}
    T_infinity = ${sphere_outer_Tinf}
  []
[]

[Executioner]
  type = Transient
  solve_type = 'PJFNK'

  petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
  petsc_options_value = 'lu       superlu_dist'

  dt = 1
  dtmin = 0.01
  end_time = 1

  nl_rel_tol = 1e-12
  nl_abs_tol = 1e-7

[]

[Outputs]
  exodus = true
  csv = true
  [Console]
    type = Console
  []
[]

[VectorPostprocessors]
  [NodalTemperature]
    type = NodalValueSampler
    sort_by = id
    boundary = '2 3'
    variable = temp
  []
[]

[Postprocessors]
  [temp_left]
    type = SideAverageValue
    boundary = 2
    variable = temp
  []

  [temp_right]
    type = SideAverageValue
    boundary = 3
    variable = temp
  []

  [flux_left]
    type = SideDiffusiveFluxIntegral
    variable = temp
    boundary = 2
    diffusivity = thermal_conductivity
  []

  [flux_right]
    type = SideDiffusiveFluxIntegral
    variable = temp
    boundary = 3
    diffusivity = thermal_conductivity
  []
  [ptot]
    type = ElementIntegralVariablePostprocessor
    variable = power_density
    block = 'fuel'
  []
  [sphere_convective_out]
    type = ConvectiveHeatTransferSideIntegral
    T_solid = temp
    boundary = '4' # outer RVP
    T_fluid = ${sphere_outer_Tinf}
    htc = ${sphere_outer_htc}
  []
  [heat_balance] # should be equal to 0 upon convergence
    type = ParsedPostprocessor
    function = '(sphere_convective_out - ptot) / ptot'
    pp_names = 'sphere_convective_out ptot'
  []
[]

(modules/thermal_hydraulics/test/tests/problems/brayton_cycle/open_brayton_cycle.i)

# This input file is used to demonstrate a simple open-air Brayton cycle using
# a compressor, turbine, shaft, motor, and generator.
# The flow length is divided into 5 segments as illustrated below, where
#   - "(I)" denotes the inlet
#   - "(C)" denotes the compressor
#   - "(T)" denotes the turbine
#   - "(O)" denotes the outlet
#   - "*" denotes a fictitious junction
#
#                  Heated section
# (I)-----(C)-----*--------------*-----(T)-----(O)
#      1       2         3          4       5
#
# Initially the fluid is at rest at ambient conditions, the shaft speed is zero,
# and no heat transfer occurs with the system.
# The transient is controlled as follows:
#   * 0   - 100 s: motor ramps up torque linearly from zero
#   * 100 - 200 s: motor ramps down torque linearly to zero, HTC ramps up linearly from zero.
#   * 200 - 300 s: (no changes; should approach steady condition)

I_motor = 1.0
motor_torque_max = 400.0

I_generator = 1.0
generator_torque_per_shaft_speed = -0.00025

motor_ramp_up_duration = 100.0
motor_ramp_down_duration = 100.0
post_motor_time = 100.0
t1 = ${motor_ramp_up_duration}
t2 = ${fparse t1 + motor_ramp_down_duration}
t3 = ${fparse t2 + post_motor_time}

D1 = 0.15
D2 = ${D1}
D3 = ${D1}
D4 = ${D1}
D5 = ${D1}

A1 = ${fparse 0.25 * pi * D1^2}
A2 = ${fparse 0.25 * pi * D2^2}
A3 = ${fparse 0.25 * pi * D3^2}
A4 = ${fparse 0.25 * pi * D4^2}
A5 = ${fparse 0.25 * pi * D5^2}

L1 = 10.0
L2 = ${L1}
L3 = ${L1}
L4 = ${L1}
L5 = ${L1}

x1 = 0.0
x2 = ${fparse x1 + L1}
x3 = ${fparse x2 + L2}
x4 = ${fparse x3 + L3}
x5 = ${fparse x4 + L4}

x2_minus = ${fparse x2 - 0.001}
x2_plus = ${fparse x2 + 0.001}
x5_minus = ${fparse x5 - 0.001}
x5_plus = ${fparse x5 + 0.001}

n_elems1 = 10
n_elems2 = ${n_elems1}
n_elems3 = ${n_elems1}
n_elems4 = ${n_elems1}
n_elems5 = ${n_elems1}

A_ref_comp = ${fparse 0.5 * (A1 + A2)}
V_comp = ${fparse A_ref_comp * 1.0}
I_comp = 1.0

A_ref_turb = ${fparse 0.5 * (A4 + A5)}
V_turb = ${fparse A_ref_turb * 1.0}
I_turb = 1.0

c0_rated_comp = 351.6925137
rho0_rated_comp = 1.146881112

rated_mfr = 0.25

speed_rated_rpm = 96000
speed_rated = ${fparse speed_rated_rpm * 2 * pi / 60.0}

speed_initial = 0

eff_comp = 0.79
eff_turb = 0.843

T_hot = 1000
T_ambient = 300
p_ambient = 1e5

[GlobalParams]
  orientation = '1 0 0'
  gravity_vector = '0 0 0'

  initial_p = ${p_ambient}
  initial_T = ${T_ambient}
  initial_vel = 0
  initial_vel_x = 0
  initial_vel_y = 0
  initial_vel_z = 0

  fp = fp_air
  closures = closures
  f = 0

  scaling_factor_1phase = '1 1 1e-5'
  scaling_factor_rhoV = 1
  scaling_factor_rhouV = 1
  scaling_factor_rhovV = 1
  scaling_factor_rhowV = 1
  scaling_factor_rhoEV = 1e-5

  rdg_slope_reconstruction = none
[]

[Functions]
  [motor_torque_fn]
    type = PiecewiseLinear
    x = '0 ${t1} ${t2}'
    y = '0 ${motor_torque_max} 0'
  []
  [motor_power_fn]
    type = ParsedFunction
    value = 'torque * speed'
    vars = 'torque speed'
    vals = 'motor_torque shaft:omega'
  []
  [generator_torque_fn]
    type = ParsedFunction
    value = 'slope * t'
    vars = 'slope'
    vals = '${generator_torque_per_shaft_speed}'
  []
  [generator_power_fn]
    type = ParsedFunction
    value = 'torque * speed'
    vars = 'torque speed'
    vals = 'generator_torque shaft:omega'
  []
  [htc_wall_fn]
    type = PiecewiseLinear
    x = '0 ${t1} ${t2}'
    y = '0 0 1e3'
  []
[]

[Modules/FluidProperties]
  [fp_air]
    type = IdealGasFluidProperties
    emit_on_nan = none
  []
[]

[Closures]
  [closures]
    type = Closures1PhaseSimple
  []
[]

[Components]
  [shaft]
    type = Shaft
    connected_components = 'motor compressor turbine generator'
    initial_speed = ${speed_initial}
  []
  [motor]
    type = ShaftConnectedMotor
    inertia = ${I_motor}
    torque = 0 # controlled
  []
  [generator]
    type = ShaftConnectedMotor
    inertia = ${I_generator}
    torque = generator_torque_fn
  []

  [inlet]
    type = InletStagnationPressureTemperature1Phase
    input = 'pipe1:in'
    p0 = ${p_ambient}
    T0 = ${T_ambient}
  []
  [pipe1]
    type = FlowChannel1Phase
    position = '${x1} 0 0'
    length = ${L1}
    n_elems = ${n_elems1}
    A = ${A1}
  []
  [compressor]
    type = ShaftConnectedCompressor1Phase
    position = '${x2} 0 0'
    inlet = 'pipe1:out'
    outlet = 'pipe2:in'
    A_ref = ${A_ref_comp}
    volume = ${V_comp}

    omega_rated = ${speed_rated}
    mdot_rated = ${rated_mfr}
    c0_rated = ${c0_rated_comp}
    rho0_rated = ${rho0_rated_comp}

    speeds = '0.5208 0.6250 0.7292 0.8333 0.9375'
    Rp_functions = 'rp_comp1 rp_comp2 rp_comp3 rp_comp4 rp_comp5'
    eff_functions = 'eff_comp1 eff_comp2 eff_comp3 eff_comp4 eff_comp5'

    min_pressure_ratio = 1.0

    speed_cr_I = 0
    inertia_const = ${I_comp}
    inertia_coeff = '${I_comp} 0 0 0'

    # assume no shaft friction
    speed_cr_fr = 0
    tau_fr_const = 0
    tau_fr_coeff = '0 0 0 0'
  []
  [pipe2]
    type = FlowChannel1Phase
    position = '${x2} 0 0'
    length = ${L2}
    n_elems = ${n_elems2}
    A = ${A2}
  []
  [junction2_3]
    type = JunctionOneToOne1Phase
    connections = 'pipe2:out pipe3:in'
  []
  [pipe3]
    type = FlowChannel1Phase
    position = '${x3} 0 0'
    length = ${L3}
    n_elems = ${n_elems3}
    A = ${A3}
  []
  [junction3_4]
    type = JunctionOneToOne1Phase
    connections = 'pipe3:out pipe4:in'
  []
  [pipe4]
    type = FlowChannel1Phase
    position = '${x4} 0 0'
    length = ${L4}
    n_elems = ${n_elems4}
    A = ${A4}
  []
  [turbine]
    type = ShaftConnectedCompressor1Phase
    position = '${x5} 0 0'
    inlet = 'pipe4:out'
    outlet = 'pipe5:in'
    A_ref = ${A_ref_turb}
    volume = ${V_turb}

    treat_as_turbine = true

    omega_rated = ${speed_rated}
    mdot_rated = ${rated_mfr}
    c0_rated = ${c0_rated_comp}
    rho0_rated = ${rho0_rated_comp}

    speeds = '0 0.5208 0.6250 0.7292 0.8333 0.9375'
    Rp_functions = 'rp_turb0 rp_turb1 rp_turb2 rp_turb3 rp_turb4 rp_turb5'
    eff_functions = 'eff_turb1 eff_turb1 eff_turb2 eff_turb3 eff_turb4 eff_turb5'

    min_pressure_ratio = 1.0

    speed_cr_I = 0
    inertia_const = ${I_turb}
    inertia_coeff = '${I_turb} 0 0 0'

    # assume no shaft friction
    speed_cr_fr = 0
    tau_fr_const = 0
    tau_fr_coeff = '0 0 0 0'
  []
  [pipe5]
    type = FlowChannel1Phase
    position = '${x5} 0 0'
    length = ${L5}
    n_elems = ${n_elems5}
    A = ${A5}
  []
  [outlet]
    type = Outlet1Phase
    input = 'pipe5:out'
    p = ${p_ambient}
  []

  [heating]
    type = HeatTransferFromSpecifiedTemperature1Phase
    flow_channel = pipe3
    T_wall = ${T_hot}
    Hw = htc_wall_fn
  []
[]

[ControlLogic]
  [motor_ctrl]
    type = TimeFunctionComponentControl
    component = motor
    parameter = torque
    function = motor_torque_fn
  []
[]

[Postprocessors]
  [heating_rate]
    type = ADHeatRateConvection1Phase
    block = 'pipe3'
    T = T
    T_wall = T_wall
    Hw = Hw
    P_hf = P_hf
    execute_on = 'INITIAL TIMESTEP_END'
  []

  [motor_torque]
    type = RealComponentParameterValuePostprocessor
    component = motor
    parameter = torque
    execute_on = 'INITIAL TIMESTEP_END'
  []
  [motor_power]
    type = FunctionValuePostprocessor
    function = motor_power_fn
    execute_on = 'INITIAL TIMESTEP_END'
  []

  [generator_torque]
    type = ShaftConnectedComponentPostprocessor
    quantity = torque
    shaft_connected_component_uo = generator:shaftconnected_uo
    execute_on = 'INITIAL TIMESTEP_END'
  []
  [generator_power]
    type = FunctionValuePostprocessor
    function = generator_power_fn
    execute_on = 'INITIAL TIMESTEP_END'
  []

  [shaft_speed]
    type = ScalarVariable
    variable = 'shaft:omega'
    execute_on = 'INITIAL TIMESTEP_END'
  []

  [p_in_comp]
    type = PointValue
    variable = p
    point = '${x2_minus} 0 0'
    execute_on = 'INITIAL TIMESTEP_END'
  []
  [p_out_comp]
    type = PointValue
    variable = p
    point = '${x2_plus} 0 0'
    execute_on = 'INITIAL TIMESTEP_END'
  []
  [p_ratio_comp]
    type = ParsedPostprocessor
    pp_names = 'p_in_comp p_out_comp'
    function = 'p_out_comp / p_in_comp'
    execute_on = 'INITIAL TIMESTEP_END'
  []

  [p_in_turb]
    type = PointValue
    variable = p
    point = '${x5_minus} 0 0'
    execute_on = 'INITIAL TIMESTEP_END'
  []
  [p_out_turb]
    type = PointValue
    variable = p
    point = '${x5_plus} 0 0'
    execute_on = 'INITIAL TIMESTEP_END'
  []
  [p_ratio_turb]
    type = ParsedPostprocessor
    pp_names = 'p_in_turb p_out_turb'
    function = 'p_in_turb / p_out_turb'
    execute_on = 'INITIAL TIMESTEP_END'
  []

  [mfr_comp]
    type = ADFlowJunctionFlux1Phase
    boundary = pipe1:out
    connection_index = 0
    equation = mass
    junction = compressor
  []
  [mfr_turb]
    type = ADFlowJunctionFlux1Phase
    boundary = pipe4:out
    connection_index = 0
    equation = mass
    junction = turbine
  []
[]

[Preconditioning]
  [pc]
    type = SMP
    full = true
  []
[]

[Executioner]
  type = Transient
  scheme = 'bdf2'

  end_time = ${t3}
  [TimeStepper]
    type = IterationAdaptiveDT
    dt = 0.01
    optimal_iterations = 5
    iteration_window = 1
    growth_factor = 1.1
    cutback_factor = 0.9
  []
  dtmin = 1e-5

  steady_state_detection = true
  steady_state_start_time = ${t2}

  solve_type = NEWTON
  nl_rel_tol = 1e-8
  nl_abs_tol = 1e-8
  nl_max_its = 15

  l_tol = 1e-4
  l_max_its = 10
[]

[Outputs]
  exodus = true
  [csv]
    type = CSV
    file_base = 'open_brayton_cycle'
    execute_vector_postprocessors_on = 'INITIAL'
  []
  [console]
    type = Console
    show = 'shaft_speed p_ratio_comp p_ratio_turb compressor:pressure_ratio turbine:pressure_ratio'
  []
[]

[Functions]
  # compressor pressure ratio
  [rp_comp1]
    type = PiecewiseLinear
    data_file = 'rp_comp1.csv'
    x_index_in_file = 0
    y_index_in_file = 1
    format = columns
    extrap = true
  []
  [rp_comp2]
    type = PiecewiseLinear
    data_file = 'rp_comp2.csv'
    x_index_in_file = 0
    y_index_in_file = 1
    format = columns
    extrap = true
  []
  [rp_comp3]
    type = PiecewiseLinear
    data_file = 'rp_comp3.csv'
    x_index_in_file = 0
    y_index_in_file = 1
    format = columns
    extrap = true
  []
  [rp_comp4]
    type = PiecewiseLinear
    data_file = 'rp_comp4.csv'
    x_index_in_file = 0
    y_index_in_file = 1
    format = columns
    extrap = true
  []
  [rp_comp5]
    type = PiecewiseLinear
    data_file = 'rp_comp5.csv'
    x_index_in_file = 0
    y_index_in_file = 1
    format = columns
    extrap = true
  []

  # compressor efficiency
  [eff_comp1]
    type = ConstantFunction
    value = ${eff_comp}
  []
  [eff_comp2]
    type = ConstantFunction
    value = ${eff_comp}
  []
  [eff_comp3]
    type = ConstantFunction
    value = ${eff_comp}
  []
  [eff_comp4]
    type = ConstantFunction
    value = ${eff_comp}
  []
  [eff_comp5]
    type = ConstantFunction
    value = ${eff_comp}
  []

  # turbine pressure ratio
  [rp_turb0]
    type = ConstantFunction
    value = 1
  []
  [rp_turb1]
    type = PiecewiseLinear
    data_file = 'rp_turb1.csv'
    x_index_in_file = 0
    y_index_in_file = 1
    format = columns
    extrap = true
  []
  [rp_turb2]
    type = PiecewiseLinear
    data_file = 'rp_turb2.csv'
    x_index_in_file = 0
    y_index_in_file = 1
    format = columns
    extrap = true
  []
  [rp_turb3]
    type = PiecewiseLinear
    data_file = 'rp_turb3.csv'
    x_index_in_file = 0
    y_index_in_file = 1
    format = columns
    extrap = true
  []
  [rp_turb4]
    type = PiecewiseLinear
    data_file = 'rp_turb4.csv'
    x_index_in_file = 0
    y_index_in_file = 1
    format = columns
    extrap = true
  []
  [rp_turb5]
    type = PiecewiseLinear
    data_file = 'rp_turb5.csv'
    x_index_in_file = 0
    y_index_in_file = 1
    format = columns
    extrap = true
  []

  # turbine efficiency
  [eff_turb1]
    type = ConstantFunction
    value = ${eff_turb}
  []
  [eff_turb2]
    type = ConstantFunction
    value = ${eff_turb}
  []
  [eff_turb3]
    type = ConstantFunction
    value = ${eff_turb}
  []
  [eff_turb4]
    type = ConstantFunction
    value = ${eff_turb}
  []
  [eff_turb5]
    type = ConstantFunction
    value = ${eff_turb}
  []
[]

Example Input File Syntax
Input Parameters
Input Files

Install MOOSE

New Users

Examples and Tutorials

Application Usage

Physics and Syntax

Application Development

Framework Development

MOOSEDocs

Infrastructure

Questions

Information and Tools

INL Applications and Remote Access

ParsedPostprocessor

Example Input File Syntax

Input Parameters

Required Parameters

Optional Parameters

Advanced Parameters

Input Files