DryCaskHeatFlux

Applies a boundary condition that models fuel rod in a dry cask storage system. The rod is assumed to be the center rod in an assembly of identical rods so that the peak cladding temperature is reached. This uses the Manteufel and Trodreas correlations inside the assembly and models the assembly-to-ambient flux using a single parameter.

Description

Lifecycle analyses for Zircaloy-clad fuel rods includes estimation of cladding integrity during dry storage and transport in storage casks. Decay heat from radioactive nuclides increases the rod internal pressure and hoop stress and causes the cladding to reach temperatures up to 400 $var element = document.getElementById("moose-equation-c6468fec-97f7-49a0-9d1f-944e36939d71");katex.render("^\\circ", element, {displayMode:false,throwOnError:false});$ C. These conditions can lead to precipitation of zirconium hydride in radial orientations, reducing the ductility of the cladding. DryCaskHeatFlux calculates the worst-case heat flux and peak clad temperatures.

Predicting the heat loss from a rod located inside an assembly which is packaged with many other assemblies inside of a single dry storage cask is difficult for several reasons. The emissivities and axial power profiles of the fuel rods have large uncertainties. In addition, the spatial distribution of power from the multiple rod assemblies is often unknown ahead of loading. The composition of the fill gas in the cask can also be difficult to predict.

These considerations lead to simplified calculations to predict worst-case heating and peak clad temperatures using homogenized models of the interior of the assembly (Manteufel and Todreas, 1994). The calculation includes terms for heat flux from the (homogenized) interior of the assembly to the edge of the rod bundle, from the edge of the rod bundle to the wall of the assembly, and from a single assembly to the exterior of the cask: $var element = document.getElementById("moose-equation-56b275e0-bb7d-4292-b3e1-228ecac41d83");katex.render("\\begin{aligned} Q &= C_1 \\left(T_m - T_e \\right) + C_2 \\left(T^4_m - T^4_e \\right) \\\\ Q &= C_3 \\left(T_e - T_w \\right) + C_4 \\left(T^4_e - T^4_w \\right) \\\\ Q &= C_5 \\left(T_w - T_a \\right) \\end{aligned}", element, {displayMode:true,throwOnError:false});$ where $var element = document.getElementById("moose-equation-03de2687-0c1f-4ad8-a463-f94f41865922");katex.render("Q", element, {displayMode:false,throwOnError:false});$ is the heat flux, $var element = document.getElementById("moose-equation-56929210-3bef-48d4-89c3-b9eefddd7599");katex.render("C_1, \\ldots, C_5", element, {displayMode:false,throwOnError:false});$ are the lumped heat conductivities, $var element = document.getElementById("moose-equation-1f2d448d-86f8-4b36-8c4b-243404dc5413");katex.render("T_m", element, {displayMode:false,throwOnError:false});$ is the temperature at middle of assembly (hottest rod), $var element = document.getElementById("moose-equation-8fe794c0-818e-486f-aa35-f4ff1ec8b5a4");katex.render("T_e", element, {displayMode:false,throwOnError:false});$ is the temperature at edge of rod bundle, $var element = document.getElementById("moose-equation-7f0726ad-41e9-46ba-8703-c24a6fd89cb0");katex.render("T_w", element, {displayMode:false,throwOnError:false});$ is the temperature at of assembly wall, and $var element = document.getElementById("moose-equation-5d6b493a-c30d-477d-aa14-11a02d82843a");katex.render("T_a", element, {displayMode:false,throwOnError:false});$ is the temperature outside cask (ambient). Manteufel and Todreas (1994) tabulate values of $var element = document.getElementById("moose-equation-e98cb6c8-014a-444d-adbb-61106d3e7e5c");katex.render("C_1, \\ldots, C_4", element, {displayMode:false,throwOnError:false});$ for common geometries for pressurized water reactor (PWR) and BWR fuel assemblies using He and N $var element = document.getElementById("moose-equation-e9fdd8da-2220-442c-8e50-06a7679b6cde");katex.render("_2", element, {displayMode:false,throwOnError:false});$ fill gases. $var element = document.getElementById("moose-equation-9e1c8050-5865-462d-982f-63fcf78b9086");katex.render("C_5", element, {displayMode:false,throwOnError:false});$ depends on the cask type and loading and is specified by the user in BISON.

Example Input Syntax

[BCs<<<{"href": "../../syntax/NuclearMaterials/BCs/index.html"}>>>]
  [cask_cooling]
    type = DryCaskHeatFlux<<<{"description": "Applies a boundary condition that models fuel rod in a dry cask storage system. The rod is assumed to be the center rod in an assembly of identical rods so that the peak cladding temperature is reached. This uses the Manteufel and Trodreas correlations inside the assembly and models the assembly-to-ambient flux using a single parameter.", "href": "DryCaskHeatFlux.html"}>>>
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = temperature
    boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = right
    bwr_or_pwr<<<{"description": "Cask assembly type: BWR or PWR"}>>> = 'pwr'
    fill_gas<<<{"description": "Fill gas in cask after drying process."}>>> = 'nitrogen'
    ambient_temperature<<<{"description": "Temperature outside of cask [K]."}>>> = 298
    cask_effective_htc<<<{"description": "Effective heat transfer coefficient from assembly to ambient [W/K]."}>>> = 6.857 # W/K per assembly, from: 1200 W/assembly, 298 K ambient, 200 C wall
    start_time<<<{"description": "Time to begin the drying process."}>>> = 0
    drying_duration<<<{"description": "Duration of the drying process before storage."}>>> = 0 # no drying
  []
[]

Input Parameters

boundaryThe list of boundary IDs from the mesh where this object applies
C++ Type:std::vector<BoundaryName>
Controllable:No
Description:The list of boundary IDs from the mesh where this object applies
bwr_or_pwrCask assembly type: BWR or PWR
C++ Type:MooseEnum
Options:bwr, pwr
Controllable:No
Description:Cask assembly type: BWR or PWR
cask_effective_htcEffective heat transfer coefficient from assembly to ambient [W/K].
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:Effective heat transfer coefficient from assembly to ambient [W/K].
fill_gasFill gas in cask after drying process.
C++ Type:MooseEnum
Options:helium, nitrogen, vacuum
Controllable:No
Description:Fill gas in cask after drying process.
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

ambient_temperature298Temperature outside of cask [K].
Default:298
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:Temperature outside of cask [K].
displacementsThe displacements
C++ Type:std::vector<VariableName>
Unit:(no unit assumed)
Controllable:No
Description:The displacements
drying_duration0Duration of the drying process before storage.
Default:0
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:Duration of the drying process before storage.
start_time0Time to begin the drying process.
Default:0
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:Time to begin the drying process.

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 BC'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 BC'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 BC'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 BC'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
skip_execution_outside_variable_domainFalseWhether to skip execution of this boundary condition when the variable it applies to is not defined on the boundary. This can facilitate setups with moving variable domains and fixed boundaries. Note that the FEProblem boundary-restricted integrity checks will also need to be turned off if using this option
Default:False
C++ Type:bool
Controllable:No
Description:Whether to skip execution of this boundary condition when the variable it applies to is not defined on the boundary. This can facilitate setups with moving variable domains and fixed boundaries. Note that the FEProblem boundary-restricted integrity checks will also need to be turned off if using this option
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

(examples/spent_fuel/full_life_cycle_coarse/discrete.i)
(test/tests/dryCask/input.i)

References

Randall D Manteufel and Neil E Todreas. Effective thermal conductivity and edge conductance model for a spent-fuel assembly. Nuclear Technology, 105(3):421–440, 1994.

@article{manteufel1994,
    author = "Manteufel, Randall D and Todreas, Neil E",
    title = "Effective thermal conductivity and edge conductance model for a spent-fuel assembly",
    journal = "Nuclear Technology",
    volume = "105",
    number = "3",
    pages = "421-440",
    year = "1994",
    publisher = "American Nuclear Society"
}

(test/tests/dryCask/input.i)

# This tests the DryCaskHeatFlux heat flux boundary condition. The test uses a uniform slab to
# simulate a piece of cladding. A heat flux is applied to the left side (the fuel side) and the
# test checks the temperature at the right side (the air side) at steady state. The test is
# designed to be used with the chart in Figure 19 of Manteufel and Todreas, which shows that for a
# PWR assembly with assembly wall temperature of 200 C, the temperature rise from the wall to the
# hottest fuel rod should be about 80 C (200 C + 80 C = 280 C = 553 K) for an assembly power of
# 1000 W.
#
# Note that the calculations used in the paper use an axial peaking factor of 1.2. Since we usually
# want to simulate an axially-resolved rod, we drop the peaking factor. Thus, a heat flux of
# 1200 W/assembly in BISON corresponds to 1000 W/assembly in Figure 19.
[Mesh]
  [mesh]
    type = GeneratedMeshGenerator
    dim = 2
    xmin = 0
    xmax = 1
    ymin = 0
    ymax = 1
    nx = 3
    ny = 1
  []
[]

[Variables]
  [temperature]
    initial_condition = 500 # Expecting 553 K (200 C at wall + 80 C delta from assembly wall to hottest rod)
  []
[]

[Kernels]
  [heat]
    type = CoefDiffusion
    variable = temperature
    coef = 1 # The value is irrelevant since we apply a heat flux boundary condition.
  []
[]

[BCs]
  [leftflux]
    type = NeumannBC
    variable = temperature
    boundary = left
    value = 43.35 # W/m^2 from 1200 W/assembly, 15x15 array, rod diam = 0.0107, length = 3.66
  []

  [cask_cooling]
    type = DryCaskHeatFlux
    variable = temperature
    boundary = right
    bwr_or_pwr = 'pwr'
    fill_gas = 'nitrogen'
    ambient_temperature = 298
    cask_effective_htc = 6.857 # W/K per assembly, from: 1200 W/assembly, 298 K ambient, 200 C wall
    start_time = 0
    drying_duration = 0 # no drying
  []
[]

[Postprocessors]
  [left_temperature]
    type = AverageNodalVariableValue
    variable = temperature
    boundary = left
  []
  [right_temperature]
    type = AverageNodalVariableValue
    variable = temperature
    boundary = right
  []
[]

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

[Outputs]
  execute_on = 'timestep_end'
  exodus = true
  [console]
    type = Console
    max_rows = 25
  []
[]

(examples/spent_fuel/full_life_cycle_coarse/discrete.i)

# This model is a linear element, 10 discrete fuel pellet stack (pellet_type_1) with a fine mesh.
# Modifying the base model to simulate the complete fuel life cycle from
# irradiation through dry storage
# Irradiation Time 3 years (6% burnup, ~ 60 MWd/kgU)
# Spent Fuel Pool 3 years
# Vacuum Drying 24 hours
# Dry Cask Storage (DCSS) 5 years
#
irrad_ramp = 8.64e4
irrad_end = 9.46944e7
cool_start = 9.47808e7
cool_end = 4.101408e8
dry_start = 4.101409e8
# dry_end = 4.102272e8   # 24 hour drying
dry_end = 4.101696e8    # 8 hour drying
# store_end = 5.679072e8 # 5 yrs storage
store_end = 4.732416e8 # 2 yrs storage
#


initial_fuel_density = 10431.0

[GlobalParams]
  # Set initial fuel density, other global parameters
  density = ${initial_fuel_density}
  initial_porosity = 0.05
  displacements = 'disp_x disp_y'
  order = FIRST
  family = LAGRANGE
  energy_per_fission = 3.2e-11  # J/fission
  temperature = temp
  volumetric_locking_correction = false
[]

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

[Mesh]
  # Specify coordinate system type
  coord_type = RZ
  # Import mesh file
  patch_update_strategy = auto
  patch_size = 20 # For contact algorithm
  partitioner = centroid
  centroid_partitioner_direction = y
  [mesh]
    type = FileMeshGenerator
    file = coarse10_rz.e
  []
[]

[Variables]
  # Define dependent variables and initial conditions
  [temp]
    initial_condition = 298.0     # set initial temp to coolant inlet
  []
[]

[AuxVariables]
  # Define auxilary variables
  [oxide_thickness]
    order = CONSTANT
    family = MONOMIAL
  []
  [fast_neutron_flux]
    block = clad
  []
  [fast_neutron_fluence]
    block = clad
  []
  [max_fission_rate]
    order = CONSTANT
    family = MONOMIAL
  []
  [grain_radius]
    block = pellet_type_1
    initial_condition = 10e-6
  []
  [creep_strain_mag]
    order = CONSTANT
    family = MONOMIAL
  []
  [gap_cond]
    order = CONSTANT
    family = MONOMIAL
  []
[]

[Functions]
  # Define functions to control power and boundary conditions
  [power_history]
    type = PiecewiseLinear
    x = '0 ${irrad_ramp}  ${irrad_end}  ${cool_start}'
    y = '0  25e3    25e3      0'
  []
  [axial_peaking_factors]
    type = PiecewiseLinear
    axis = y
    x = '0.00324 0.0151 0.10998 0.12184'
    y = '1.0     1.0     1.0    1.0'
  []
  [pressure_ramp]              # reads and interpolates input data defining amplitude curve for fill gas pressure
    type = PiecewiseLinear
    x = '-200 0'
    y = '0 1'
  []
  [coolant_pressure]
    type = PiecewiseLinear
    #             ---- irrad ----             --- pool ---           - storage -
    x = '0    ${irrad_ramp} ${irrad_end} ${cool_start} ${cool_end} ${dry_start} ${dry_end} ${store_end}'
    y = '1e5     15.5e6        15.5e6         2e5          2e5          1e5        1e5         1e5'
  []
  [coolant_temperature]
    type = PiecewiseLinear
    #             ---- irrad ----             --- pool ---           - storage -
    x = '0    ${irrad_ramp} ${irrad_end} ${cool_start} ${cool_end} ${dry_start} ${dry_end} ${store_end}'
    y = '300      587            587          308          308          308         308        308'
  []
  [coolant_htc]
    type = PiecewiseLinear
    # From CoolantChannel model, HTC falls from 37000 to 22000 as the oxide grows.
    # Coolant flow is maintained until after CZP, then 1 more day. Flow is then reduced until the
    # correct htc for natural convection is achieved (~400 W/m2-K).
    # Drying is handled by DryCaskHeatFlux.
    x = '0    ${irrad_ramp}  7e7   ${irrad_end} ${cool_start} ${cool_end} ${dry_start} ${dry_end} ${store_end}'
    y = '37e3      37e3     25e3       22e3           400         400           0           0           0'
  []
[]

[Physics/SolidMechanics/QuasiStatic]
  [fuel]
    block = pellet_type_1
    add_variables = true
    strain = FINITE
    eigenstrain_names = 'fuel_relocation_strain fuel_thermal_strain fuel_volumetric_strain'
    generate_output = 'vonmises_stress stress_xx stress_yy stress_zz'
    decomposition_method = EigenSolution
    extra_vector_tags = 'ref'
  []
  [clad]
    block = clad
    add_variables = true
    strain = FINITE
    eigenstrain_names = 'clad_thermal_strain clad_irradiation_growth_strain'
    generate_output = 'vonmises_stress stress_xx stress_yy stress_zz'
    decomposition_method = EigenSolution
    extra_vector_tags = 'ref'
  []
[]

[Kernels]
  # Define kernels for the various terms in the PDE system
  [gravity]       # body force term in stress equilibrium equation
    type = Gravity
    variable = disp_y
    value = -9.81
  []
  [heat]         # gradient term in heat conduction equation
    type = HeatConduction
    variable = temp
    extra_vector_tags = 'ref'
  []
  [heat_ie]       # time term in heat conduction equation
    type = HeatConductionTimeDerivative
    variable = temp
    extra_vector_tags = 'ref'
  []
  [heat_source]  # source term in heat conduction equation
    type = NeutronHeatSource
    variable = temp
    extra_vector_tags = 'ref'
    block = pellet_type_1 # fission rate applied to the fuel (block 2) only
    fission_rate = fission_rate
    decay_heat_function = decay_heat_function
    max_fission_rate = max_fission_rate
  []
[]

[Burnup]
  [burnup]
    block = pellet_type_1
    rod_ave_lin_pow = power_history          # using the power function defined above
    axial_power_profile = axial_peaking_factors     # using the axial power profile function defined above
    num_radial = 80
    num_axial = 11
    a_lower = 0.00324    # mesh dependent!
    a_upper = 0.12184    # mesh dependent!
    fuel_inner_radius = 0
    fuel_outer_radius = .0041
    fuel_volume_ratio = 0.987775 # for use with dished pellets (ratio of actual volume to cylinder volume)
    order = CONSTANT
    family = MONOMIAL
    RPF = RPF
    #N235 = N235 # Activate to write N235 concentration to output file
    #N238 = N238 # Activate to write N238 concentration to output file
    #N239 = N239 # Activate to write N239 concentration to output file
    #N240 = N240 # Activate to write N240 concentration to output file
    #N241 = N241 # Activate to write N241 concentration to output file
    #N242 = N242 # Activate to write N242 concentration to output file
  []
[]

[AuxKernels]
  # Define auxilliary kernels for each of the aux variables
  [oxide_thickness]
    type = MaterialRealAux
    boundary = 2
    variable = oxide_thickness
    property = oxide_scale_thickness
  []
  [fast_neutron_flux]
    type = FastNeutronFluxAux
    variable = fast_neutron_flux
    block = clad
    rod_ave_lin_pow = power_history
    axial_power_profile = axial_peaking_factors
    factor = 3e13
    execute_on = timestep_begin
  []
  [fast_neutron_fluence]
    type = FastNeutronFluenceAux
    variable = fast_neutron_fluence
    block = clad
    fast_neutron_flux = fast_neutron_flux
    execute_on = timestep_begin
  []
  [max_fission_rate]
    type = MaxFissionRateAux
    variable = max_fission_rate
    block = pellet_type_1
    fission_rate = fission_rate
    execute_on = timestep_begin
  []
  [grain_radius]
    type = GrainRadiusAux
    block = pellet_type_1
    variable = grain_radius
    temperature = temp
    execute_on = nonlinear
  []
  [creep_strain_mag]
    type = MaterialRealAux
    block = clad
    property = effective_creep_strain
    variable = creep_strain_mag
    execute_on = timestep_end
  []
  [conductance]
    type = MaterialRealAux
    property = gap_conductance
    variable = gap_cond
    boundary = 10
  []
[]

[Contact]
  # Define mechanical contact between the fuel (sideset=10) and the clad (sideset=5)
  [pellet_clad_mechanical]
    primary = 5
    secondary = 10
    formulation = kinematic
    model = frictionless
  []
[]

[ThermalContact]
  # Define thermal contact between the fuel (sideset=10) and the clad (sideset=5)
  [thermal_contact]
    type = GasGapHeatTransfer
    variable = temp
    primary = 5
    secondary = 10
    initial_moles = initial_moles       # coupling to a postprocessor which supplies the initial plenum/gap gas mass
    gas_released = fission_gas_released     # coupling to a postprocessor which supplies the fission gas addition
    plenum_pressure = plenum_pressure
    contact_pressure = contact_pressure
    quadrature = true
  []
[]

[BCs]
# Define boundary conditions
  [no_x_all] # pin pellets and clad along axis of symmetry (y)
    type = DirichletBC
    variable = disp_x
    boundary = 12
    value = 0.0
  []
  [no_y_clad_bottom] # pin clad bottom in the axial direction (y)
    type = DirichletBC
    variable = disp_y
    boundary = '1'
    value = 0.0
  []
  [no_y_fuel_bottom] # pin fuel bottom in the axial direction (y)
    type = DirichletBC
    variable = disp_y
    boundary = '1020'
    value = 0.0
  []
  [Pressure] # apply coolant pressure on clad outer walls
    [coolantPressure]
      boundary = '1 2 3'
      factor = 1.0
      function = coolant_pressure
    []
  []
  [PlenumPressure] # apply plenum pressure on clad inner walls and pellet surfaces
    [plenumPressure]
      boundary = 9
      initial_pressure = 2.0e6
      startup_time = 0
      R = 8.3143
      output_initial_moles = initial_moles       # coupling to post processor to get initial fill gas mass
      temperature = plenum_temperature            # coupling to post processor to get gas temperature approximation
      volume = plenum_volume                        # coupling to post processor to get gas volume
      material_input = fission_gas_released          # coupling to post processor to get fission gas added
      output = plenum_pressure                   # coupling to post processor to output plenum/gap pressure
      displacements = 'disp_x disp_y'
    []
  []
  [convective_clad_surface] # apply convective boundary to clad outer surface
    type = ConvectiveFluxFunction
    boundary = '1 2 3'
    variable = temp
    coefficient = 'coolant_htc'
    T_infinity = 'coolant_temperature'
  []
  [cask_cooling]
    type = DryCaskHeatFlux
    variable = temp
    boundary = '1 2 3'
    bwr_or_pwr = 'pwr'
    fill_gas = 'helium'
    ambient_temperature = 298
    cask_effective_htc = 3.1 # W/K from each assembly to ambient
    start_time = ${cool_end}
    drying_duration = 86400
  []
[]

[Controls]
  [DCSS]
    type = TimePeriod
    disable_objects = 'BCs/convective_clad_surface'
    start_time = ${cool_end}
    end_time = 1e9
  []
[]

[Materials]
  # Define material behavior models and input material property data
  [fuel_thermal]                       # temperature and burnup dependent thermal properties of UO2 (BISON kernel)
    type = UO2Thermal
    block = 'pellet_type_1'
    thermal_conductivity_model = NFIR
    temperature = temp
    burnup_function = burnup
  []
  [ZryOxidation]
    type = ZryOxidation
    boundary = '2'
    clad_inner_radius = 0.00418
    clad_outer_radius = 0.00474
    normal_operating_temperature_model = epri_kwu_ce
    high_temperature_model = leistikow
    use_coolant_channel = true
    fast_neutron_flux = fast_neutron_flux
    outputs = all
  []
  [fuel_volumetric_swelling]
    type = UO2VolumetricSwellingEigenstrain
    block = 'pellet_type_1'
    burnup_function = burnup
    initial_fuel_density = 10431.0
    gas_swelling_model_type = SIFGRS
    eigenstrain_name = fuel_volumetric_strain
  []
  [fuel_elasticity_tensor]
    type = UO2ElasticityTensor
    block = 'pellet_type_1'
  []
  [fuel_elastic_stress]
    type = ComputeFiniteStrainElasticStress
    block = 'pellet_type_1'
  []
  [fuel_thermal_expansion]
    type = ComputeThermalExpansionEigenstrain
    block = 'pellet_type_1'
    thermal_expansion_coeff = 10.0e-6
    stress_free_temperature = 298.0
    eigenstrain_name = fuel_thermal_strain
  []
  [fuel_relocation]
    type = UO2RelocationEigenstrain
    block = 'pellet_type_1'
    burnup_function = burnup
    diameter = 0.0082
    rod_ave_lin_pow = power_history
    axial_power_profile = axial_peaking_factors
    diametral_gap =160e-6
    burnup_relocation_stop = 0.3
    relocation_activation1 = 5000
    eigenstrain_name = fuel_relocation_strain
  []
  [clad_thermal]
    type = HeatConductionMaterial
    block = 'clad'
    thermal_conductivity = 16.0
    specific_heat = 330.0
  []
  [clad_elasticity_tensor]
    type = ZryElasticityTensor
    block = 'clad'
  []
  [clad_creep_model]
    type = ZryCreepLimbackHoppeUpdate
    block = 'clad'
    fast_neutron_flux = fast_neutron_flux
    fast_neutron_fluence = fast_neutron_fluence
    zircaloy_material_type = stress_relief_annealed
    model_irradiation_creep = true
    model_primary_creep = true
    model_thermal_creep = true
  []
  [clad_inelastic_stress]
    type = ComputeMultipleInelasticStress
    block = 'clad'
    tangent_operator = elastic
    inelastic_models = 'clad_creep_model'
  []
  [clad_thermal_expansion]
    type = ZryThermalExpansionMATPROEigenstrain
    block = clad
    stress_free_temperature = 298.0
    eigenstrain_name = clad_thermal_strain
  []
  [clad_irradiation_growth]
     type = ZryIrradiationGrowthEigenstrain
     block = clad
     fast_neutron_fluence = fast_neutron_fluence
     zircaloy_material_type = stress_relief_annealed
     eigenstrain_name = clad_irradiation_growth_strain
  []
  [fission_gas_release]
    type = UO2Sifgrs
    block = pellet_type_1
    temperature = temp
    burnup_function = burnup
    grain_radius = grain_radius
    gbs_model = true
  []
  [clad_density]
    type = StrainAdjustedDensity
    block = clad
    strain_free_density = 6551.0
  []
  [fuel_density]
    type = StrainAdjustedDensity
    block = pellet_type_1
    strain_free_density = ${initial_fuel_density}
  []
[]

[Dampers]
  [BoundingValueNodalDamper]
    type = BoundingValueNodalDamper
    variable = temp
    max_value = 3200
    min_value = 0
  []
  [limitX]
    type = MaxIncrement
    max_increment = 1e-5
    variable = disp_x
  []
[]

[Executioner]

  # PETSC options:
  #   petsc_options
  #   petsc_options_iname
  #   petsc_options_value
  #
  # controls for linear iterations
  #   l_max_its
  #   l_tol
  #
  # controls for nonlinear iterations
  #   nl_max_its
  #   nl_rel_tol
  #   nl_abs_tol
  #
  # time control
  #   start_time
  #   dt
  #   optimal_iterations
  #   iteration_window
  #   linear_iteration_ratio

  type = Transient
  solve_type = 'PJFNK'

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

  line_search = 'none'

  l_max_its = 50
  l_tol = 1e-3
  nl_max_its = 35
  nl_rel_tol = 1e-4
  nl_abs_tol = 1e-8
  start_time = -200
  n_startup_steps = 1
  end_time = ${store_end}
  dtmax = 2e6
  dtmin = 1

  [TimeStepper]
    type = IterationAdaptiveDT
    dt = 2e2
    optimal_iterations = 20
    iteration_window = 6
    time_t = '0    ${irrad_ramp} ${irrad_end} ${cool_start} ${cool_end} ${dry_start} ${dry_end} ${store_end}'
    time_dt ='1e3       1e4          1e3           100          100          100          100       100'
    growth_factor = 1.5
    cutback_factor = .6
  []

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

[Postprocessors]
  # Define postprocessors (some are required as specified above; others are optional; many others are available)
  [clad_inner_vol]              # volume inside of cladding
    type = InternalVolume
    boundary = 7
    outputs = exodus
    execute_on = 'initial timestep_end'
  []
  [fis_gas_grain]
    type = ElementIntegralFisGasGrainSifgrs
    block = pellet_type_1
    outputs = exodus
    execute_on = linear
  []
  [fis_gas_boundary]
    type = ElementIntegralFisGasBoundarySifgrs
    block = pellet_type_1
    outputs = exodus
    execute_on = linear
  []
  [flux_from_clad]           # area integrated heat flux from the cladding
    type = SideDiffusiveFluxIntegral
    variable = temp
    boundary = 5
    diffusivity = thermal_conductivity
    execute_on = timestep_end
  []
  [flux_from_fuel]          # area integrated heat flux from the fuel
    type = SideDiffusiveFluxIntegral
    variable = temp
    boundary = 10
    diffusivity = thermal_conductivity
    execute_on = timestep_end
  []
  [_dt]                     # time step
    type = TimestepSize
  []
  [num_lin_it]
    type = NumLinearIterations
  []
  [num_nonlin_it]
    type = NumNonlinearIterations
  []
  [tot_lin_it]
    type = CumulativeValuePostprocessor
    postprocessor = num_lin_it
  []
  [tot_nonlin_it]
    type = CumulativeValuePostprocessor
    postprocessor = num_nonlin_it
  []
  [alive_time]
    type = PerfGraphData
    section_name = Root
    data_type = TOTAL
  []
  [decay_heat_function]
    type = DecayHeatFunction
    time_at_shutdown = ${cool_start}
    table_or_sum = sum
  []
  [peak_clad_temp]
    type = NodalExtremeValue
    variable = temp
    block = 'clad'
    execute_on = 'timestep_end'
  []
  [max_clad_hoop_stress]
    type = ElementExtremeValue
    variable = stress_zz
    block = 'clad'
    value_type = 'max'
    execute_on = 'timestep_end'
  []
  [rod_total_power]
    type = ElementIntegralPower
    variable = temp
    burnup_function = burnup
    block = pellet_type_1
  []
  [rod_input_power]
    type = FunctionValuePostprocessor
    function = power_history
    scale_factor = 0.1186 # rod height
  []
  [peak_oxide_thickness]
    type = ElementExtremeValue
    variable = oxide_thickness
    block = 'clad'
    value_type = 'max'
    execute_on = 'timestep_end'
  []
[]

[VectorPostprocessors]
  [clad_surf_props]
    type = LineValueSampler
    variable = 'oxide_thickness temp stress_zz'
    start_point = '0.00467 0.0001 0'
    end_point = '0.00467 0.1279 0'
    num_points = 100
    sort_by = y
    outputs = 'outfile_1'
  []
[]

[PerformanceMetricOutputs]
[]

[StandardLWRFuelRodOutputs]
  fuel_pellet_blocks = pellet_type_1
[]

[Outputs]
  perf_graph = true
  exodus = true
  color = false
  csv = true
  [console]
    type = Console
    max_rows = 25
  []
  [outfile_1]
    type = CSV
    execute_on = 'FINAL'
  []
  [chkfile]
    type = CSV
    show = 'peak_clad_temp peak_oxide_thickness max_clad_hoop_stress'
    execute_on = final
  []
[]

(test/tests/dryCask/input.i)

# This tests the DryCaskHeatFlux heat flux boundary condition. The test uses a uniform slab to
# simulate a piece of cladding. A heat flux is applied to the left side (the fuel side) and the
# test checks the temperature at the right side (the air side) at steady state. The test is
# designed to be used with the chart in Figure 19 of Manteufel and Todreas, which shows that for a
# PWR assembly with assembly wall temperature of 200 C, the temperature rise from the wall to the
# hottest fuel rod should be about 80 C (200 C + 80 C = 280 C = 553 K) for an assembly power of
# 1000 W.
#
# Note that the calculations used in the paper use an axial peaking factor of 1.2. Since we usually
# want to simulate an axially-resolved rod, we drop the peaking factor. Thus, a heat flux of
# 1200 W/assembly in BISON corresponds to 1000 W/assembly in Figure 19.
[Mesh]
  [mesh]
    type = GeneratedMeshGenerator
    dim = 2
    xmin = 0
    xmax = 1
    ymin = 0
    ymax = 1
    nx = 3
    ny = 1
  []
[]

[Variables]
  [temperature]
    initial_condition = 500 # Expecting 553 K (200 C at wall + 80 C delta from assembly wall to hottest rod)
  []
[]

[Kernels]
  [heat]
    type = CoefDiffusion
    variable = temperature
    coef = 1 # The value is irrelevant since we apply a heat flux boundary condition.
  []
[]

[BCs]
  [leftflux]
    type = NeumannBC
    variable = temperature
    boundary = left
    value = 43.35 # W/m^2 from 1200 W/assembly, 15x15 array, rod diam = 0.0107, length = 3.66
  []

  [cask_cooling]
    type = DryCaskHeatFlux
    variable = temperature
    boundary = right
    bwr_or_pwr = 'pwr'
    fill_gas = 'nitrogen'
    ambient_temperature = 298
    cask_effective_htc = 6.857 # W/K per assembly, from: 1200 W/assembly, 298 K ambient, 200 C wall
    start_time = 0
    drying_duration = 0 # no drying
  []
[]

[Postprocessors]
  [left_temperature]
    type = AverageNodalVariableValue
    variable = temperature
    boundary = left
  []
  [right_temperature]
    type = AverageNodalVariableValue
    variable = temperature
    boundary = right
  []
[]

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

[Outputs]
  execute_on = 'timestep_end'
  exodus = true
  [console]
    type = Console
    max_rows = 25
  []
[]

Description
Example Input Syntax
Input Parameters
Input Files
References