ConstitutiveHeatConduction

Computes the thermal diffusion component/kernel for $var element = document.getElementById("moose-equation-360f4655-ae19-4246-92cd-3028045cc2ee");katex.render("(k \\nabla T, \\nabla \\psi)", element, {displayMode:false,throwOnError:false});$ term with Jacobians contributions introduced by the thermal conductivity calculated via optional material properties.

Description

ConstitutiveHeatConduction calculates the thermal flux term $var element = document.getElementById("moose-equation-2d342c18-aa5b-499f-9546-bece9f6e7c81");katex.render("\\mathbf{q}", element, {displayMode:false,throwOnError:false});$ of the heat equation, $var element = document.getElementById("moose-equation-8b8d02af-5110-4a7a-ac5d-a8e22824737a");katex.render(" \\mathbf{q} = - k(u_i, u_j,...) \\nabla T", element, {displayMode:true,throwOnError:false});$ where $var element = document.getElementById("moose-equation-d14bb2d9-2992-498a-917d-a4f082778c87");katex.render("T", element, {displayMode:false,throwOnError:false});$ is the temperature, and $var element = document.getElementById("moose-equation-3d630ba0-14b1-4e6b-9448-975bd20e0e45");katex.render("k", element, {displayMode:false,throwOnError:false});$ is the thermal conductivity dependent on any number of primal variables $var element = document.getElementById("moose-equation-2e83ec08-55b6-4304-8b80-6d3bfda17629");katex.render("u", element, {displayMode:false,throwOnError:false});$ . The derivatives of $var element = document.getElementById("moose-equation-3fba6ed6-8f39-4350-a7b6-50468285bf38");katex.render("k", element, {displayMode:false,throwOnError:false});$ with respect to the primal variables $var element = document.getElementById("moose-equation-ee1e7341-5498-4fc7-b29e-ceda974e19f6");katex.render("u", element, {displayMode:false,throwOnError:false});$ can be optionally passed to ConstitutiveHeatConduction via the thermal_conductivity_derivs vector. A corresponding vector of coupled primal variables is also required to be passed via thermal_conductivity_args and is used to calculate the Jacobian for both on- and off-diagonal terms.

Example Input Syntax

[Kernels<<<{"href": "../../syntax/Kernels/index.html"}>>>]
  [HeatDiff]
    type = ConstitutiveHeatConduction<<<{"description": "Computes the thermal diffusion component/kernel for $(k \\nabla T, \\nabla \\psi)$ term with Jacobians contributions introduced by the thermal conductivity calculated via optional material properties.", "href": "ConstitutiveHeatConduction.html"}>>>
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = T
    thermal_conductivity<<<{"description": "Fission rate material name"}>>> = thermal_conductivity
    thermal_conductivity_args<<<{"description": "Thermal conductivity formulation arguments"}>>> = u
    thermal_conductivity_derivs<<<{"description": "Material property names of thermal conductivity derivatives with respect to thermal_conductivity_args"}>>> = thermal_conductivity_du
  []
[]

Input Parameters

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

Required Parameters

blockThe list of blocks (ids or names) that this object will be applied
C++ Type:std::vector<SubdomainName>
Controllable:No
Description:The list of blocks (ids or names) that this object will be applied
displacementsThe displacements
C++ Type:std::vector<VariableName>
Unit:(no unit assumed)
Controllable:No
Description:The displacements
thermal_conductivitythermal_conductivityFission rate material name
Default:thermal_conductivity
C++ Type:MaterialPropertyName
Unit:(no unit assumed)
Controllable:No
Description:Fission rate material name
thermal_conductivity_argsThermal conductivity formulation arguments
C++ Type:std::vector<VariableName>
Unit:(no unit assumed)
Controllable:No
Description:Thermal conductivity formulation arguments
thermal_conductivity_derivsMaterial property names of thermal conductivity derivatives with respect to thermal_conductivity_args
C++ Type:std::vector<MaterialPropertyName>
Unit:(no unit assumed)
Controllable:No
Description:Material property names of thermal conductivity derivatives with respect to thermal_conductivity_args

Optional Parameters

absolute_value_vector_tagsThe tags for the vectors this residual object should fill with the absolute value of the residual contribution
C++ Type:std::vector<TagName>
Controllable:No
Description:The tags for the vectors this residual object should fill with the absolute value of the residual contribution
extra_matrix_tagsThe extra tags for the matrices this Kernel should fill
C++ Type:std::vector<TagName>
Controllable:No
Description:The extra tags for the matrices this Kernel should fill
extra_vector_tagsThe extra tags for the vectors this Kernel should fill
C++ Type:std::vector<TagName>
Controllable:No
Description:The extra tags for the vectors this Kernel should fill
matrix_tagssystemThe tag for the matrices this Kernel should fill
Default:system
C++ Type:MultiMooseEnum
Options:nontime, system
Controllable:No
Description:The tag for the matrices this Kernel should fill
vector_tagsnontimeThe tag for the vectors this Kernel should fill
Default:nontime
C++ Type:MultiMooseEnum
Options:nontime, time
Controllable:No
Description:The tag for the vectors this Kernel should fill

Contribution To Tagged Field Data Parameters

control_tagsAdds user-defined labels for accessing object parameters via control logic.
C++ Type:std::vector<std::string>
Controllable:No
Description:Adds user-defined labels for accessing object parameters via control logic.
diag_save_inThe name of auxiliary variables to save this Kernel's diagonal Jacobian contributions to. Everything about that variable must match everything about this variable (the type, what blocks it's on, etc.)
C++ Type:std::vector<AuxVariableName>
Unit:(no unit assumed)
Controllable:No
Description:The name of auxiliary variables to save this Kernel's diagonal Jacobian contributions to. Everything about that variable must match everything about this variable (the type, what blocks it's on, etc.)
enableTrueSet the enabled status of the MooseObject.
Default:True
C++ Type:bool
Controllable:Yes
Description:Set the enabled status of the MooseObject.
implicitTrueDetermines whether this object is calculated using an implicit or explicit form
Default:True
C++ Type:bool
Controllable:No
Description:Determines whether this object is calculated using an implicit or explicit form
save_inThe name of auxiliary variables to save this Kernel's residual contributions to. Everything about that variable must match everything about this variable (the type, what blocks it's on, etc.)
C++ Type:std::vector<AuxVariableName>
Unit:(no unit assumed)
Controllable:No
Description:The name of auxiliary variables to save this Kernel's residual contributions to. Everything about that variable must match everything about this variable (the type, what blocks it's on, etc.)
seed0The seed for the master random number generator
Default:0
C++ Type:unsigned int
Controllable:No
Description:The seed for the master random number generator
use_displaced_meshFalseWhether or not this object should use the displaced mesh for computation. Note that in the case this is true but no displacements are provided in the Mesh block the undisplaced mesh will still be used.
Default:False
C++ Type:bool
Controllable:No
Description:Whether or not this object should use the displaced mesh for computation. Note that in the case this is true but no displacements are provided in the Mesh block the undisplaced mesh will still be used.

Advanced Parameters

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

Material Property Retrieval Parameters

Input Files

(test/tests/constitutive_heat_conduction/2d_coupled.i)
(examples/metal_fuel/uzr/pin.i)
(test/tests/constitutive_heat_conduction/2d_steady_state_final_prob.i)
(assessment/metallic_fuel/TREAT/M7/analysis/binary_t427.i)
(test/tests/constitutive_heat_conduction/1D_transient.i)

References

No citations exist within this document.

(test/tests/constitutive_heat_conduction/2d_coupled.i)

# This test solves a 2D steady state heat equation
# The error is found by comparing to the analytical solution

# Note that the thermal conductivity, specific heat, and density in this problem
# Are set to 1, and need to be changed to the constants of the material being
# Analyzed

[Mesh]
  [mesh]
    type = GeneratedMeshGenerator
    dim = 2
    nx = 30
    ny = 30
    xmax = 2
    ymax = 2
  []
[]

[Variables]
  [T]
    initial_condition = 300
  []
  [u]
    initial_condition = 5.5
  []
[]

[Kernels]
  [u_diff]
    type = Diffusion
    variable = u
  []
  [u_dt]
    type = TimeDerivative
    variable = u
  []
  [HeatDiff]
    type = ConstitutiveHeatConduction
    variable = T
    thermal_conductivity = thermal_conductivity
    thermal_conductivity_args = u
    thermal_conductivity_derivs = thermal_conductivity_du
  []
  [HeatTdot]
    type = ConstitutiveHeatConductionTimeDerivative
    variable = T
    density_name = 1
    specific_heat = specific_heat
    specific_heat_args = u
    specific_heat_derivs = specific_heat_du
  []
[]

[BCs]
  [zero]
    type = DirichletBC
    variable = T
    boundary = 'left right bottom'
    value = 300
  []
  [top]
    type = DirichletBC
    variable = T
    boundary = top
    value = 500
  []
  [u_left]
    type = DirichletBC
    variable = u
    value = 10
    boundary = left
  []
  [u_right]
    type = DirichletBC
    variable = u
    value = 1
    boundary = right
  []
[]

[Materials]
  [thcond]
    type = ParsedMaterial
    property_name = thermal_conductivity
    coupled_variables = u
    expression = 'u^2'
  []
  [thcond_du]
    type = ParsedMaterial
    property_name = thermal_conductivity_du
    coupled_variables = u
    expression = '2*u'
  []
  [spheat]
    type = ParsedMaterial
    property_name = specific_heat
    coupled_variables = u
    expression = 'u^3'
  []
  [spheat_du]
    type = ParsedMaterial
    property_name = specific_heat_du
    coupled_variables = u
    expression = '3*u^2'
  []
[]

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

[Executioner]
  type = Transient

  num_steps = 1
[]

(test/tests/constitutive_heat_conduction/2d_coupled.i)

# This test solves a 2D steady state heat equation
# The error is found by comparing to the analytical solution

# Note that the thermal conductivity, specific heat, and density in this problem
# Are set to 1, and need to be changed to the constants of the material being
# Analyzed

[Mesh]
  [mesh]
    type = GeneratedMeshGenerator
    dim = 2
    nx = 30
    ny = 30
    xmax = 2
    ymax = 2
  []
[]

[Variables]
  [T]
    initial_condition = 300
  []
  [u]
    initial_condition = 5.5
  []
[]

[Kernels]
  [u_diff]
    type = Diffusion
    variable = u
  []
  [u_dt]
    type = TimeDerivative
    variable = u
  []
  [HeatDiff]
    type = ConstitutiveHeatConduction
    variable = T
    thermal_conductivity = thermal_conductivity
    thermal_conductivity_args = u
    thermal_conductivity_derivs = thermal_conductivity_du
  []
  [HeatTdot]
    type = ConstitutiveHeatConductionTimeDerivative
    variable = T
    density_name = 1
    specific_heat = specific_heat
    specific_heat_args = u
    specific_heat_derivs = specific_heat_du
  []
[]

[BCs]
  [zero]
    type = DirichletBC
    variable = T
    boundary = 'left right bottom'
    value = 300
  []
  [top]
    type = DirichletBC
    variable = T
    boundary = top
    value = 500
  []
  [u_left]
    type = DirichletBC
    variable = u
    value = 10
    boundary = left
  []
  [u_right]
    type = DirichletBC
    variable = u
    value = 1
    boundary = right
  []
[]

[Materials]
  [thcond]
    type = ParsedMaterial
    property_name = thermal_conductivity
    coupled_variables = u
    expression = 'u^2'
  []
  [thcond_du]
    type = ParsedMaterial
    property_name = thermal_conductivity_du
    coupled_variables = u
    expression = '2*u'
  []
  [spheat]
    type = ParsedMaterial
    property_name = specific_heat
    coupled_variables = u
    expression = 'u^3'
  []
  [spheat_du]
    type = ParsedMaterial
    property_name = specific_heat_du
    coupled_variables = u
    expression = '3*u^2'
  []
[]

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

[Executioner]
  type = Transient

  num_steps = 1
[]

(examples/metal_fuel/uzr/pin.i)

initial_fuel_density = 15800

[GlobalParams]
  energy_per_fission = 3.2e-11  # J/fission
  temperature = temp
[]

[Mesh]
  coord_type = RZ
  [mesh]
    type = GeneratedMeshGenerator
    dim = 1
    nx = 300
    xmax = 0.002
  []
[]

[Problem]
  type = ReferenceResidualProblem
  reference_vector = 'ref'
  extra_tag_vectors = 'ref'
[]

[Variables]
  [temp]
    initial_condition = 298
  []
  [X_Zr]
    scaling = 1e10
    initial_condition = 0.225
  []
[]

[Functions]
  [power_history]
    type = PiecewiseLinear
    x = '0 100000 3901600 4333600 6148000 7948800  '
    y = '0 40016 40344 42312 41984 41984  '
  []
  [temp_history]
    type = PiecewiseLinear
    x = '0 100000 3901600 4333600 6148000 7948800  '
    y = '298.0 805.5 803.8 811.1 808.5 807.4'
  []
[]

[Kernels]
  [heat]
    type = ConstitutiveHeatConduction
    variable = temp
    thermal_conductivity = 'thermal_conductivity'
    thermal_conductivity_args = 'temp X_Zr'
    thermal_conductivity_derivs = 'thermal_conductivity_dT thermal_conductivity_dZr'
    extra_vector_tags = 'ref'
  []
  [heat_ie]
    type = ConstitutiveHeatConductionTimeDerivative
    variable = temp
    specific_heat = 'specific_heat'
    specific_heat_args = 'temp'
    specific_heat_derivs = 'specific_heat_dT'
    extra_vector_tags = 'ref'
  []
  [heat_source]
    type = FissionRateHeatSource
    variable = temp
    fission_rate_args = 'X_Zr'
    fission_rate_derivs = 'fission_rate_dZr'
    extra_vector_tags = 'ref'
  []
  [ZrDiffusion]
    type = ZirconiumDiffusion
    variable = X_Zr
    temperature = temp
    extra_vector_tags = 'ref'
  []
  [diffusion_ie]
    type = TimeDerivative
    variable = X_Zr
    extra_vector_tags = 'ref'
  []
[]

[BCs]
  [temp_rightside]
    type = FunctionDirichletBC
    variable = temp
    boundary = right
    function = temp_history
  []
[]

[Materials]
  [fission_rate]
    type = UPuZrFissionRate
    rod_linear_power = power_history
    axial_power_profile = 1.0
    pellet_radius = 0.002
    initial_X_Zr = 0.225
    X_Zr = X_Zr
    X_Pu_function = 0.163
  []
  [burnup]
    type = UPuZrBurnup
    initial_X_Zr = 0.225
    initial_X_Pu = 0.163
    density = ${initial_fuel_density}
  []
  [metal_fuel_thermal]
    type = UPuZrThermal
    spheat_model = savage
    thcond_model = lanl
    porosity = porosity
    X_Zr = X_Zr
    X_Pu = 0.163
  []
  [swelling]
    type = UPuZrVolumetricSwellingEigenstrain
    burnup = 0
    hydrostatic_stress = 0
    eigenstrain_name = 'swelling'
    output_properties = 'porosity gaseous_porosity'
  []
  [fuel_density]
    type = GenericConstantMaterial
    prop_names = density
    prop_values = ${initial_fuel_density}
  []
  [phase]
    type = PhaseUPuZr
    X_Zr = X_Zr
    X_Pu = 0.163
    AB_temp = 965.15
    CD_temp = 995.15
  []
  [zr_diff]
    type = ZrDiffusivityUPuZr
    X_Zr = X_Zr
    X_Pu = 0.163
    p_alpha = 0.20
    p_delta = 0.20
    p_beta  = 0.20
    p_gamma = 0.20
    D0_scale_alpha = 15.0
    D0_scale_delta = 1.0
    D0_scale_beta = 2.0
    D0_scale_gamma = 7.0
    outputs = all
  []
[]

[Dampers]
  [bounded_zr]
    type = BoundingValueNodalDamper
    variable = X_Zr
    max_value = 0.83
    min_value = 0
  []
[]

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

[Executioner]
  type = Transient

  solve_type = 'PJFNK'
  petsc_options_iname = '-pc_type'
  petsc_options_value = 'lu'
  line_search = 'none'

  end_time = 3e6
  dtmin = 1e-2
  dtmax = 1e6
  nl_max_its = 30
  l_max_its = 100
  nl_abs_tol = 1e-9
  nl_rel_tol = 1e-9
  l_tol = 1e-05

  [TimeStepper]
    type = IterationAdaptiveDT
    dt = 1e3
    optimal_iterations = 8
    iteration_window = 2
    growth_factor = 2
    cutback_factor = .5
    timestep_limiting_function = power_history
    force_step_every_function_point = true
  []
[]

[Postprocessors]
  [dt]
    type = TimestepSize
  []
  [num_lin]
    type = NumLinearIterations
  []
  [num_nonlin]
    type = NumNonlinearIterations
  []
  [ave_temp]
    type = ElementAverageValue
    variable = temp
    block = 0
  []
  [X_Zr]
    type = NodalExtremeValue
    variable = X_Zr
    block = 0
    value_type = max
  []
[]

[Outputs]
  perf_graph = true
  #exodus = true

  [csv_output]
    type = CSV
    precision = 8
    file_base = pin_out
    hide = 'dt num_lin num_nonlin'
    execute_on = final
  []
[]

(test/tests/constitutive_heat_conduction/2d_steady_state_final_prob.i)

# This test solves a 2D steady state heat equation
# The error is found by comparing to the analytical solution
#
# This test duplicates the test in moose/modules/heat_conduction/test/tests/verify_against_analytical

[Mesh]
  [mesh]
    type = GeneratedMeshGenerator
    dim = 2
    nx = 30
    ny = 30
    xmax = 2
    ymax = 2
  []
[]

[Variables]
  [T]
  []
[]

[Kernels]
  [HeatDiff]
    type = ConstitutiveHeatConduction
    variable = T
  []
[]

[BCs]
  [zero]
    type = DirichletBC
    variable = T
    boundary = 'left right bottom'
    value = 0
  []
  [top]
    type = FunctionDirichletBC
    variable = T
    boundary = top
    function = '10*sin(pi*x*0.5)'
  []
[]

[Materials]
  [properties]
    type = GenericConstantMaterial
    prop_names = 'thermal_conductivity specific_heat density'
    prop_values = '1 1 1'
  []
[]

[Postprocessors]
  [nodal_error]
    type = NodalL2Error
    function = '10/(sinh(pi))*sin(pi*x*0.5)*sinh(pi*y*0.5)'
    variable = T
  []
  [elemental_error]
    type = ElementL2Error
    function = '10/(sinh(pi))*sin(pi*x*0.5)*sinh(pi*y*0.5)'
    variable = T
  []
[]

[Executioner]
  type = Steady
[]

[Outputs]
  exodus = true
[]

(assessment/metallic_fuel/TREAT/M7/analysis/binary_t427.i)

# TREAT M7 U-10Zr pin clad in HT9 transient only simulation.
#   Not using transferred initial condition.
# ANL-IFR-124 Test M7 T-427
# Units: m, W, kg, Pa

[Debug]
  show_var_residual = 'disp_x disp_y temperature'
  show_var_residual_norms = true
[]

[GlobalParams]
  order = SECOND
  family = LAGRANGE
  energy_per_fission = 3.28451e-11
  displacements = 'disp_x disp_y'
  stress_free_temperature = 586.0
  temperature = temperature
  # Some calculated values, which are nice to know.
  # u_weight = 0.90126
  # pu_weight = 0.00115
  # zr_weight = 0.0976
  # X_U = 0.779
  X_Pu = 0.000969
  X_Zr = 0.22
[]

[Problem]
  type = ReferenceResidualProblem
  extra_tag_vectors = 'ref'
  reference_vector = 'ref'
  group_variables = 'disp_x disp_y'
[]

[Mesh]
  coord_type = RZ
  [smeared_pellet_mesh]
    type = FuelPinMeshGenerator
    clad_bot_gap_height = 0.38e-3
    clad_gap_width = 1e-6
    clad_thickness = 0.38e-3
    clad_top_gap_height = 252.35e-3
    pellet_height = 37.3e-2
    pellet_outer_radius = 2.54e-3
    pellet_quantity = 1
    bottom_clad_height = 1.9e-3
    top_clad_height = 1.9e-3

    elem_type = QUAD8
    clad_mesh_density = customize
    pellet_mesh_density = customize
    nx_c = 3
    nx_p = 5
    ny_c = 150
    ny_cl = 3
    ny_cu = 3
    ny_p = 175
  []
  patch_size = 25
  partitioner = centroid
  patch_update_strategy = auto
  centroid_partitioner_direction = y
[]

[Functions]
  [power_history] # Peak power density
    type = PiecewiseLinear
    data_file = 'm7fa_digitized.csv'
    format = columns
    scale_factor = 78.8 # g; grams of fuel in pellet for W/g to W
  []
  [axial_peaking_factors] # Peak factor (1 max)
    type = PiecewiseBilinear
    data_file = 'm7_axial_tranAlone.csv'
    axis = 1
  []
  [q_heat]
    type = CompositeFunction
    functions = 'power_history axial_peaking_factors'
    scale_factor = 132273.8097 # m^-3; volume of fuel slug, 1/V, to make power density
  []
  [dt_fun]
    type = PiecewiseLinear
    x = '0   6   7    9    15   20    25'
    y = '0.5 0.5 0.01 0.01 0.01 0.001 0.1'
  []
[]

[Variables]
  [temperature]
    initial_condition = 586.0
  []
[]

[Physics/SolidMechanics/QuasiStatic]
  add_variables = true
  strain = FINITE
  generate_output = 'stress_xx stress_yy stress_zz vonmises_stress hydrostatic_stress elastic_strain_xx elastic_strain_yy elastic_strain_zz strain_xx strain_yy strain_zz'
  [slug_mech]
    block = pellet
    additional_generate_output = 'creep_strain_xx creep_strain_yy creep_strain_zz'
    eigenstrain_names = 'slug_thermal_strain'
    extra_vector_tags = 'ref'
  []
  [clad_mech]
    block = clad
    additional_generate_output = 'creep_strain_xx creep_strain_yy creep_strain_zz'
    eigenstrain_names = 'clad_thermal_eigenstrain'
    extra_vector_tags = 'ref'
  []
[]

[Kernels]
  [gravity]
    type = Gravity
    variable = disp_y
    value = -9.81
    extra_vector_tags = 'ref'
  []
  [heat]
    type = ConstitutiveHeatConduction
    extra_vector_tags = 'ref'
    variable = temperature
    thermal_conductivity = 'thermal_conductivity'
    thermal_conductivity_args = 'temperature'
    thermal_conductivity_derivs = 'thermal_conductivity_dT'
  []
  [heat_ie]
    type = ConstitutiveHeatConductionTimeDerivative
    extra_vector_tags = 'ref'
    variable = temperature
    specific_heat = 'specific_heat'
    specific_heat_args = 'temperature'
    specific_heat_derivs = 'specific_heat_dT'
  []
  [vol_heat_source]
    type = HeatSource
    extra_vector_tags = 'ref'
    block = pellet
    function = q_heat
    variable = temperature
  []
[]

[AuxVariables]
  [power_density]
    block = pellet
  []
  [fission_rate]
    block = pellet
  []
  [fast_neutron_flux]
    block = clad
  []
  [fast_neutron_fluence]
    block = clad
  []
  [creep_rate_aux]
    order = CONSTANT
    family = MONOMIAL
  []

  # Aux variables for output
  [energy_density]
    block = pellet
    initial_condition = 0.0
  []
  [gap_cond]
    order = CONSTANT
    family = MONOMIAL
  []
  [coolant_htc]
    order = CONSTANT
    family = MONOMIAL
    block = clad
  []
  [coolant_temperature]
    order = CONSTANT
    family = MONOMIAL
    block = clad
  []
  [cumulative_damage_index]
    order = CONSTANT
    family = MONOMIAL
    block = clad
  []
  [element_failed]
    order = CONSTANT
    family = MONOMIAL
  []
  [eutectic_thickness_in]
    order = CONSTANT
    family = MONOMIAL
    block = clad
  []
  [eutectic_thickness_out]
    order = CONSTANT
    family = MONOMIAL
    block = clad
  []
[]

[AuxKernels]
  [calc_fission_rate]
    type = FissionRateGeneral
    fission_rate_formulation = POWER_DENSITY
    block = pellet
    power_density_function = q_heat
    variable = fission_rate
    execute_on = 'initial timestep_end'
  []
  [fast_neutron_flux]
    type = FastNeutronFluxAux
    variable = fast_neutron_flux
    block = clad
    factor = 0.6e19
    execute_on = 'initial timestep_end'
  []
  [fast_neutron_fluence]
    type = FastNeutronFluenceAux
    variable = fast_neutron_fluence
    block = clad
    fast_neutron_flux = fast_neutron_flux
    execute_on = 'initial timestep_end'
  []
  [calc_power_density]
    type = FunctionAux
    block = pellet
    function = q_heat
    variable = power_density
    execute_on = 'initial timestep_end'
  []
  [calc_energy_density]
    type = VariableTimeIntegrationAux
    block = pellet
    order = 2
    variable_to_integrate = power_density
    variable = energy_density
    execute_on = timestep_end
  []
  # Hoop stress_zz.
  [creep_rate_aux]
    type = MaterialRealAux
    property = creep_rate
    variable = creep_rate_aux
    execute_on = timestep_end
  []
  [conductance]
    type = MaterialRealAux
    property = gap_conductance
    variable = gap_cond
    boundary = 10
  []
  [coolant_temperature]
    type = MaterialRealAux
    property = coolant_temperature
    variable = coolant_temperature
    boundary = 2
  []
  [coolant_htc]
    type = MaterialRealAux
    property = coolant_channel_htc
    variable = coolant_htc
    boundary = 2
  []
  [cdf_amount]
    block = clad
    type = MaterialRealAux
    property = cdf_failure
    variable = cumulative_damage_index
  []
  [failed_element]
    type = MaterialRealAux
    property = failed
    variable = element_failed
  []
  [fcci_eutectic_in]
    execute_on = timestep_end
    type = EutecticThicknessFCCI
    boundary = 5
    variable = eutectic_thickness_in
  []
  [fcci_eutectic_out]
    execute_on = timestep_end
    type = EutecticThicknessFCCI
    boundary = 2
    variable = eutectic_thickness_out
  []
[]

[Postprocessors]
  [_dt]
    type = TimestepSize
    outputs = all
  []
  [fis_gas_released]
    type = ElementIntegralMaterialProperty
    mat_prop = fis_gas_rel
    block = pellet
    execute_on = 'initial timestep_end'
  []
  [ave_temp_interior]
    type = SideAverageValue
    boundary = 9
    variable = temperature
    execute_on = 'initial linear'
  []
  [avg_clad_temp]
    type = ElementAverageValue
    variable = temperature
    outputs = all
    block = clad
  []
  [avg_fuel_temp]
    type = ElementAverageValue
    variable = temperature
    outputs = all
    block = pellet
  []
  [max_outer_clad_temp]
    type = NodalExtremeValue
    variable = temperature
    value_type = max
    boundary = 2
  []
  [peak_clad_temp]
    type = ElementExtremeValue
    variable = temperature
    value_type = max
    block = clad
    outputs = all
  []
  [peak_fuel_temp]
    type = ElementExtremeValue
    variable = temperature
    value_type = max
    block = pellet
    outputs = all
  []
  [peak_outer_fuel_temp]
    type = NodalExtremeValue
    variable = temperature
    value_type = max
    boundary = 10
  []
  [peak_coolant_temperature]
    type = ElementExtremeValue
    variable = coolant_temperature
    value_type = max
    block = clad
    outputs = all
  []
  [max_power_density]
    type = ElementExtremeValue
    variable = power_density
    value_type = max
    block = pellet
  []
  [avg_power_density]
    type = ElementAverageValue
    variable = power_density
    block = pellet
    outputs = all
  []
  [total_deposit_pin_energy]
    type = ElementIntegralVariablePostprocessor
    variable = energy_density
    block = pellet
  []
  [clad_inner_vol]
    type = InternalVolume
    boundary = 7
  []
  [pellet_volume]
    type = InternalVolume
    boundary = 8
    execute_on = 'initial timestep_end'
  []
  [gas_volume]
    type = InternalVolume
    boundary = 9
    execute_on = 'initial linear'
    outputs = all
    addition = -2.0e-6 # m3; assuming all bond sodium displaced by closed gap
  []
  [clad_fuel_gap]
    type = NodalExtremeValue
    variable = penetration
    boundary = 10
    outputs = all
  []
  [max_eutectic_pen_in]
    type = ElementExtremeValue
    variable = eutectic_thickness_in
    block = clad
    value_type = max
    execute_on = timestep_end
    outputs = all
  []
  [max_eutectic_pen_out]
    type = ElementExtremeValue
    variable = eutectic_thickness_out
    block = clad
    value_type = max
    execute_on = timestep_end
    outputs = all
  []
[]

[Contact]
  [pellet_clad_mechanical]
    primary = 5
    secondary = 10
    penalty = 1e12
    model = frictionless
    normalize_penalty = true
    tangential_tolerance = 1e-3
    normal_smoothing_distance = 0.1
  []
[]

[ThermalContact]
  [thermal_contact]
    type = GapHeatTransfer
    variable = temperature
    primary = 5
    secondary = 10
    quadrature = true
    gap_conductivity = 73.70
    min_gap = 1e-3
  []
[]

[BCs]
  [no_x_all]
    type = DirichletBC
    variable = disp_x
    boundary = 12
    value = 0.0
  []
  [no_y_fuel]
    type = PenaltyDirichletBC
    penalty = 1e10
    variable = disp_y
    boundary = 20
    value = 0.0
  []
  [no_y_clad]
    type = DirichletBC
    variable = disp_y
    boundary = 1
    value = 0.0
  []
  [Pressure]
    [coolantPressure]
      boundary = '1 2 3'
      factor = 455054.0
    []
  []
  [PlenumPressure]
    [plenumPressure]
      boundary = 9
      initial_pressure = 1.41343e6
      initial_temperature = 309.0
      startup_time = 0.0
      R = 8.3143
      temperature = ave_temp_interior
      volume = gas_volume
      output = plenum_pressure
      material_input = fis_gas_released
      output_initial_moles = plenum_moles
    []
  []
[]

[CoolantChannel]
  [convective_clad_surface]
    boundary = '1 2 3'
    variable = temperature
    inlet_temperature = 586.0
    inlet_pressure = 455054.0
    inlet_massflux = 4520.72
    coolant_material = sodium
    flow_area = 2.22e-5
    hydraulic_diameter = 2.057e-3
    heated_perimeter = 1.835e-2
    heated_diameter = 4.84e-3
    number_axial_zone = 50
    htc_correlation_type = 3
    compute_enthalpy = true
  []
[]

[Materials]
  # Fuel Slug Properties
  [set_porosity]
    type = GenericConstantMaterial
    block = pellet
    prop_names = porosity
    prop_values = 0.31
  []
  [set_mat_fission_rate]
    type = ParsedMaterial
    block = pellet
    coupled_variables = 'fission_rate'
    expression = 'fission_rate * 1.0'
    property_name = 'fission_rate'
  []
  [melted]
    type = GenericMaterialFailure
    block = pellet
    compared = greater_than
    variable_check = temperature
    constant_criteria = 1600.0
  []
  [for_youngs]
    # Pulled equation from Bison manual and added temperature limit.
    type = ParsedMaterial
    block = pellet
    outputs = all
    output_properties = 'youngs_modulus'
    property_name = 'youngs_modulus'
    coupled_variables = 'temperature'
    material_property_names = 'porosity'
    constant_names = 'T_limit E_U T_meltU W_Zr W_Pu B_E Ta_start Ta_end'
    constant_expressions = '1550.0 1.6e11 1405.0 0.0976 0.00115 1.2 923.0 973.0'
    expression = 'E_p := 1.0 - B_E * porosity; E_W := (1.0 + 0.17 * W_Zr) / (1.0 + 1.34 * W_Zr) - W_Pu; T_act := if(temperature > T_limit, T_limit, temperature); x_smooth := (T_act - Ta_start) / (Ta_end - Ta_start); f_smooth := if(T_act < Ta_start, 0.0, if(T_act > Ta_end, 1.0, 6.0 * pow(x_smooth, 5) - 15.0 * pow(x_smooth, 4) + 10.0 * pow(x_smooth, 3))); E_T := 1.0 - 1.06 * (T_act - 588.0) / T_meltU - f_smooth * 0.3 * (1.0 - 1.06 * (Ta_end - 588.0) / T_meltU); E_U * E_T * E_p * E_W'
  []
  [for_poissons]
    # Pulled equation from Bison manual and added temperature limit.
    type = ParsedMaterial
    block = pellet
    outputs = all
    output_properties = 'poissons_ratio'
    property_name = 'poissons_ratio'
    coupled_variables = 'temperature'
    material_property_names = 'porosity'
    constant_names = 'T_limit nu_U T_meltU W_Zr B_nu'
    constant_expressions = '1550.0 0.24 1405.0 0.0976 0.8'
    expression = 'nu_p := 1.0 - B_nu * porosity; nu_W := (1.0 + 3.4 * W_Zr) / (1.0 + 1.9 * W_Zr); T_act := if(temperature > T_limit, T_limit, temperature); nu_T := 1.0 + 1.2 * (T_act - 588.0) / T_meltU; nu_U *nu_p * nu_W * nu_T'
  []
  [slug_elasticity_tensor]
    type = ComputeVariableIsotropicElasticityTensor
    block = pellet
    args = temperature
    youngs_modulus = youngs_modulus
    poissons_ratio = poissons_ratio
  []
  [slug_stress]
    type = ComputeMultipleInelasticStress
    tangent_operator = nonlinear
    inelastic_models = 'slug_creep'
    block = pellet
  []
  [slug_creep]
    type = UPuZrCreepUpdate
    block = pellet
    porosity = porosity
    max_inelastic_increment = 1e-3
  []
  [slug_thermal_expansion]
    type = ComputeThermalExpansionEigenstrain
    block = pellet
    thermal_expansion_coeff = 1.18e-5
    eigenstrain_name = slug_thermal_strain
  []
  [slug_thermal]
    type = UPuZrThermal
    block = pellet
    Na_depth = 0.0673
    fuel_outer_radius = 2.54e-3
    spheat_model = savage
    thcond_model = billone
    porosity_model = partially_logged
    initiating_porosity = 0.24
    outputs = all
    output_properties = thermal_conductivity
    porosity = porosity
  []
  [slug_density]
    type = StrainAdjustedDensity
    block = pellet
    strain_free_density = 15700.0
  []
  [slug_fission_gas_release]
    block = pellet
    type = UPuZrFissionGasRelease
    fission_rate = fission_rate
  []

  # Cladding Properties
  [clad_elasticity_tensor]
    type = ComputeIsotropicElasticityTensor
    youngs_modulus = 1.88e11
    poissons_ratio = 0.236
    block = clad
  []
  [clad_stress]
    type = ComputeMultipleInelasticStress
    block = clad
    tangent_operator = nonlinear
    inelastic_models = 'clad_creep'
  []
  [clad_fast_flux]
    type = FastNeutronFlux
    block = clad
    factor = 0.6e19
  []
  [clad_creep]
    type = HT9CreepUpdate
    block = clad
  []
  [clad_thermal_expansion]
    type = ComputeThermalExpansionEigenstrain
    block = clad
    thermal_expansion_coeff = 1.2e-5
    eigenstrain_name = clad_thermal_eigenstrain
  []
  [clad_thermal]
    type = HT9Thermal
    block = clad
  []
  [clad_density]
    type = StrainAdjustedDensity
    block = clad
    strain_free_density = 7874.0
  []
  [failclad]
    block = clad
    type = HT9FailureClad
    hoop_stress = stress_zz
    method = cdf_short
  []
[]

[Preconditioning]
  [SMP]
    type = SMP
    full = true
    solve_type = 'PJFNK'
    petsc_options = '-snes_ksp_ew -snes_converged_reason'
    petsc_options_iname = '-pc_type -pc_factor_mat_solver_package -ksp_gmres_restart'
    petsc_options_value = 'lu       superlu_dist                  51'
  []
[]

[Executioner]
  type = Transient
  automatic_scaling = true
  compute_scaling_once = true
  line_search = 'none'
  l_max_its = 100
  l_tol = 1e-3
  nl_max_its = 200
  nl_rel_tol = 1e-6
  nl_abs_tol = 1e-9
  start_time = 0.0
  end_time = 25.0
  dtmin = 0.0001
  dtmax = 1.0
  [TimeStepper]
    type = FunctionDT
    function = dt_fun
  []
  [Quadrature]
    order = FIFTH
    side_order = SEVENTH
  []
[]

[PerformanceMetricOutputs]
[]

[Outputs]
  perf_graph = true
  time_step_interval =  1
  exodus = true
  csv = true # Figures comparing to experiment results are based off this csv output.
  [console]
    type = Console
    output_linear = true
    output_nonlinear = true
  []
  [check]
    # Numerical assessment test check at points near important features.
    #  Checking near the feature but not during it hopefully will not affect the feature computation.
    type = CSV
    execute_on = TIMESTEP_END
    sync_only = true
    sync_times = '10 16.5 20'
  []
[]

(test/tests/constitutive_heat_conduction/1D_transient.i)

# This test solves a 1D transient heat equation
# The error is caclulated by comparing to the analytical solution
# The problem setup and analytical solution are taken from "Advanced Engineering
# Mathematics, 10th edition" by Erwin Kreyszig.
# http://www.amazon.com/Advanced-Engineering-Mathematics-Erwin-Kreyszig/dp/0470458364
# It is Example 1 in section 12.6 on page 561
#
# This test duplicates the test in moose/modules/heat_conduction/test/tests/verify_against_analytical

[Mesh]
  [mesh]
    type = GeneratedMeshGenerator
    dim = 1
    nx = 160
    xmax = 80
  []
[]

[Variables]
  [T]
  []
[]

[ICs]
  [T_IC]
    type = FunctionIC
    variable = T
    function = '100*sin(pi*x/80)'
  []
[]

[Kernels]
  [HeatDiff]
    type = ConstitutiveHeatConduction
    variable = T
  []
  [HeatTdot]
    type = ConstitutiveHeatConductionTimeDerivative
    variable = T
  []
[]

[BCs]
  [sides]
    type = DirichletBC
    variable = T
    boundary = 'left right'
    value = 0
  []
[]

[Materials]
  [k]
    type = GenericConstantMaterial
    prop_names = 'thermal_conductivity'
    prop_values = '0.95' #copper in cal/(cm sec C)
  []
  [cp]
    type = GenericConstantMaterial
    prop_names = 'specific_heat'
    prop_values = '0.092' #copper in cal/(g C)
  []
  [rho]
    type = GenericConstantMaterial
    prop_names = 'density'
    prop_values = '8.92' #copper in g/(cm^3)
  []
[]

[Postprocessors]
  [error]
    type = NodalL2Error
    function = '100*sin(pi*x/80)*exp(-0.95/(0.092*8.92)*pi^2/80^2*t)'
    variable = T
  []
[]

[Executioner]
  type = Transient
  scheme = bdf2
  l_tol = 1e-6
  dt = 2
  end_time = 100
[]

[Outputs]
  exodus = true
[]

Description
Example Input Syntax
Input Parameters
Input Files
References