CumulativeDamageIndex

Computes the cumulative damage to cladding for low temperature steady-state conditions.

Description

Intergranular stress corrosion cracking, typically caused by pellet-clad interaction (PCI), can lead to clad failure during normal operation. A cumulative damage model can be used to predict failure in the form of a cumulative damage index. Two models are available, the FALCON model (Rashid et al., 1988) and the ASTM model (Roberts et al., 1979).

FALCON Model

A cumulative damage model is used in BISON to estimate cladding damage for low temperature steady-state conditions. The model is based on the notion of a cumulative damage index, which has the following form. $(1)var element = document.getElementById("moose-equation-4c2f4745-d9c3-499e-ac67-160809925a21");katex.render(" D = \\int_{0}^{t_n} \\frac{dt}{t_f (\\sigma_{\\text{hoop}}, B, T)}", element, {displayMode:true,throwOnError:false});$ In Eq. (1), $var element = document.getElementById("moose-equation-94960206-f831-4c60-b6c8-4ebafe7f3b04");katex.render("D", element, {displayMode:false,throwOnError:false});$ is the amount of damage at time $var element = document.getElementById("moose-equation-551247cb-bb7a-4b12-ba1c-288b7b7f8789");katex.render("t_n", element, {displayMode:false,throwOnError:false});$ (all time in seconds), $var element = document.getElementById("moose-equation-48c5e946-9043-4d56-a5b9-9ace9fe57c56");katex.render("t_f", element, {displayMode:false,throwOnError:false});$ is the failure time at stress $var element = document.getElementById("moose-equation-0cad54f7-ee1e-490b-bd33-b343c7fc4b72");katex.render("\\sigma_{\\text{hoop}}", element, {displayMode:false,throwOnError:false});$ (all stress in units of MPa), $var element = document.getElementById("moose-equation-a9a1a084-e47b-4a5b-916a-a6f695b71a38");katex.render("T", element, {displayMode:false,throwOnError:false});$ is temperature (K), and $var element = document.getElementById("moose-equation-4b0051be-b1d4-406d-82e7-b746e66ea3c5");katex.render("B", element, {displayMode:false,throwOnError:false});$ is burnup (MWd/MTU).

The variable $var element = document.getElementById("moose-equation-f957c484-2033-435b-bf31-3f662b12775b");katex.render("t_f", element, {displayMode:false,throwOnError:false});$ has the form: $(2)var element = document.getElementById("moose-equation-56cc82f3-55e7-4de5-98e9-f22f09ea3327");katex.render(" t_f = \\bar{t} \\exp{\\left[(1.015 \\sigma_y + 1.74 \\sigma_{\\text{ref}} - 2.755 \\sigma_{\\text{hoop}} ) 0.01\\right]}", element, {displayMode:true,throwOnError:false});$ where $var element = document.getElementById("moose-equation-f017372b-04e2-434e-be49-46fd5a9e59e8");katex.render("\\sigma_y", element, {displayMode:false,throwOnError:false});$ is the clad yield stress, and $(3)var element = document.getElementById("moose-equation-31d723f2-1d21-4621-aa7b-95acfe88c2df");katex.render(" \\bar{t} = 5\\times10^5\\left(1.13\\times10^{-4}B-0.13\\right)^{-0.75}\\exp{\\left[-30\\left(1-\\frac{611}{T}\\right)\\right]}", element, {displayMode:true,throwOnError:false});$ where the stress $var element = document.getElementById("moose-equation-a44a46c7-7077-4d8f-ada1-d28ec7a4ff13");katex.render("\\sigma_{\\text{ref}}", element, {displayMode:false,throwOnError:false});$ is a threshold stress, which has the form: $(4)var element = document.getElementById("moose-equation-64f42c20-65ea-4c82-b5fc-b380cc118a92");katex.render(" \\sigma_{\\text{ref}} = \\begin{cases} 336.476 \\left(B-5000\\right)^{-0.07262} &\\text{ for Zr2}\\\\ 310.275 \\left(B-5000\\right)^{-0.0440} &\\text{ for Zr4}. \\end{cases}", element, {displayMode:true,throwOnError:false});$ The model for cumulative damage index (Eq. (1)) activates only when $var element = document.getElementById("moose-equation-970f2db8-e035-4fa4-8fd3-6b1548df5a0a");katex.render("\\sigma_{\\text{hoop}} > \\sigma_{\\text{ref}}", element, {displayMode:false,throwOnError:false});$ and $var element = document.getElementById("moose-equation-dc1a8987-2763-4f1f-a05c-86a410b89481");katex.render("B", element, {displayMode:false,throwOnError:false});$ $var element = document.getElementById("moose-equation-025d3217-6ace-447b-a443-b6b9f3ac557c");katex.render(">", element, {displayMode:false,throwOnError:false});$ 5000 MWd/MTU.

ASTM Model

An alternative model is developed for the PCI screening in BISON. The equation for failure time is written as $var element = document.getElementById("moose-equation-7340de1e-3068-43d6-ad76-430c4feabbf1");katex.render("t_f = -\\frac{14.7}{\\max{\\left(C_{I_2}, C_1\\right)}} \\left[ \\ln{\\left( 1.0 - \\frac{1.0 - \\left[ \\dfrac{\\sigma_{\\text{hoop}}}{\\sigma_{y}} \\right]}{1.633\\cos{\\phi}\\sin{\\phi}} \\right)} \\right]^3", element, {displayMode:true,throwOnError:false});$ where $var element = document.getElementById("moose-equation-4b866728-ed04-4c05-8d60-d9119c7a8fc1");katex.render("\\sigma_{\\text{hoop}}", element, {displayMode:false,throwOnError:false});$ is the clad hoop stress, $var element = document.getElementById("moose-equation-b754c2b1-6553-4e0e-90a2-cddacc6c40e7");katex.render("\\phi", element, {displayMode:false,throwOnError:false});$ is the texture angle, $var element = document.getElementById("moose-equation-10f45409-d4be-4f54-81f5-69d28bb8624a");katex.render("C_{I_2}", element, {displayMode:false,throwOnError:false});$ is iodine concentration in $var element = document.getElementById("moose-equation-71f77ae3-78ac-4553-abc4-264db6489a48");katex.render("g/cm^2", element, {displayMode:false,throwOnError:false});$ . $var element = document.getElementById("moose-equation-58f87dd0-1f7c-4812-b8a0-a70bfed5a3b3");katex.render("C_1", element, {displayMode:false,throwOnError:false});$ = $var element = document.getElementById("moose-equation-24e691e1-b2b6-44f4-9c9d-362b41351b53");katex.render("1.0\\times10^{-6}g/cm^2", element, {displayMode:false,throwOnError:false});$ . Assuming $var element = document.getElementById("moose-equation-705339c0-d468-4f41-9ee8-6d847956f2f5");katex.render("1\\%", element, {displayMode:false,throwOnError:false});$ fission gas release, the iodine concentration is computed as $var element = document.getElementById("moose-equation-e7a704ae-e8a0-4836-bae2-2276a4130869");katex.render("7\\times10^{-4}\\text{Bu}", element, {displayMode:false,throwOnError:false});$ using an empirical equation $var element = document.getElementById("moose-equation-4d509327-f733-44f1-af83-b4bb127d247e");katex.render("C_{I_2}=\\left( 0.02\\times{\\text{Fraction of U atoms}}\\times{\\text{Bu}}\\times{\\text{FGR}} \\right)", element, {displayMode:false,throwOnError:false});$ in Roberts et al. (1979). This equation also indicates a threshold stress for the accumulation of damage. The threshold stress is $var element = document.getElementById("moose-equation-adaf2cb7-b299-4269-b887-ba14974c1311");katex.render(" \\sigma_{th} = \\sigma_y \\left(1.0-1.633\\cos{\\phi}\\sin{\\phi}\\right)", element, {displayMode:true,throwOnError:false});$ The accumulation of damage is limited by a maximum time $var element = document.getElementById("moose-equation-75e778da-a639-427c-bdb2-8d8e498dc791");katex.render("\\tau", element, {displayMode:false,throwOnError:false});$ after reaching the peak stress, and the parameter $var element = document.getElementById("moose-equation-b1a3aafd-3c5a-4772-8e30-00f02db66d70");katex.render("\\tau", element, {displayMode:false,throwOnError:false});$ is calibrated from ramp test data and is set as 1000 sec.

Example Input Syntax

[AuxKernels<<<{"href": "../../syntax/AuxKernels/index.html"}>>>]
  [cumulative_damage_index]
    type = CumulativeDamageIndex<<<{"description": "Computes the cumulative damage to cladding for low temperature steady-state conditions.", "href": "CumulativeDamageIndex.html"}>>>
    variable<<<{"description": "The name of the variable that this object applies to"}>>> = cumulative_damage_index
    burnup<<<{"description": "Name of Postprocessor holding the average burnup."}>>> = burnup
    temperature<<<{"description": "Coupled temperature (K)"}>>> = temp
    hoop_stress<<<{"description": "Coupled hoop stress."}>>> = stress_zz
    yield_stress<<<{"description": "Constant value of the yield stress of the clad (in Pa)."}>>> = 240e6
    clad_type<<<{"description": "The cladding type."}>>> = 'Zr2'
    block<<<{"description": "The list of blocks (ids or names) that this object will be applied"}>>> = 1
    execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = timestep_end
  []
[]

Input Parameters

hoop_stressCoupled hoop stress.
C++ Type:std::vector<VariableName>
Unit:(no unit assumed)
Controllable:No
Description:Coupled hoop stress.
temperatureCoupled temperature (K)
C++ Type:std::vector<VariableName>
Unit:(no unit assumed)
Controllable:No
Description:Coupled temperature (K)
variableThe name of the variable that this object applies to
C++ Type:AuxVariableName
Unit:(no unit assumed)
Controllable:No
Description:The name of the variable that this object applies to

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
boundaryThe list of boundaries (ids or names) from the mesh where this object applies
C++ Type:std::vector<BoundaryName>
Controllable:No
Description:The list of boundaries (ids or names) from the mesh where this object applies
burnupburnupName of Postprocessor holding the average burnup.
Default:burnup
C++ Type:PostprocessorName
Unit:(no unit assumed)
Controllable:No
Description:Name of Postprocessor holding the average burnup.
check_boundary_restrictedTrueWhether to check for multiple element sides on the boundary in the case of a boundary restricted, element aux variable. Setting this to false will allow contribution to a single element's elemental value(s) from multiple boundary sides on the same element (example: when the restricted boundary exists on two or more sides of an element, such as at a corner of a mesh
Default:True
C++ Type:bool
Controllable:No
Description:Whether to check for multiple element sides on the boundary in the case of a boundary restricted, element aux variable. Setting this to false will allow contribution to a single element's elemental value(s) from multiple boundary sides on the same element (example: when the restricted boundary exists on two or more sides of an element, such as at a corner of a mesh
clad_typeZr4The cladding type.
Default:Zr4
C++ Type:MooseEnum
Options:Zr2, Zr4
Controllable:No
Description:The cladding type.
execute_onTIMESTEP_ENDThe list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html.
Default:TIMESTEP_END
C++ Type:ExecFlagEnum
Options:XFEM_MARK, NONE, INITIAL, LINEAR, NONLINEAR_CONVERGENCE, NONLINEAR, POSTCHECK, TIMESTEP_END, TIMESTEP_BEGIN, MULTIAPP_FIXED_POINT_END, MULTIAPP_FIXED_POINT_BEGIN, FINAL, CUSTOM, PRE_DISPLACE
Controllable:No
Description:The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html.
initial_burnup0Burnup in FIMA.
Default:0
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:Burnup in FIMA.
iodine_concentration0.03Iodine concentration in g/cm^2.
Default:0.03
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:Iodine concentration in g/cm^2.
max_relaxation_time1000Maximum time for limiting CDI calculation.
Default:1000
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:Maximum time for limiting CDI calculation.
model_optionFALCONThe model to use for computing the CDI
Default:FALCON
C++ Type:MooseEnum
Options:FALCON, ASTM
Controllable:No
Description:The model to use for computing the CDI
stress_concentration_factor1Stress concentration factor.
Default:1
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:Stress concentration factor.
texture_angle0.523599Texture angle.
Default:0.523599
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:Texture angle.
yield_stressConstant value of the yield stress of the clad (in Pa).
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:Constant value of the yield stress of the clad (in Pa).
yield_stress_propertyThe name of the material property containing the cladding yield stress (in Pa).
C++ Type:MaterialPropertyName
Unit:(no unit assumed)
Controllable:No
Description:The name of the material property containing the cladding yield stress (in Pa).

Optional Parameters

control_tagsAdds user-defined labels for accessing object parameters via control logic.
C++ Type:std::vector<std::string>
Controllable:No
Description:Adds user-defined labels for accessing object parameters via control logic.
enableTrueSet the enabled status of the MooseObject.
Default:True
C++ Type:bool
Controllable:Yes
Description:Set the enabled status of the MooseObject.
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/cumulative_damage_index/cumulative_damage_index_astm_model.i)
(test/tests/cumulative_damage_index/cumulative_damage_index_yield_prop_test.i)
(test/tests/cumulative_damage_index/cumulative_damage_index_test.i)

References

Y. R. Rashid, A. J. Zangari, and C. L. Lin. Modeling of PCI Under Steady State and Transient Operating Conditions. In Proceedings of a technical committee meeting organized by the IAEA: Water Reactor Fuel Element Computer Modelling in Steady State, Transient, and Accident Conditions. Preston, United Kingdom, September 18-22, 1988.

@inproceedings{joerashidIAEA1988,
    author = "Rashid, Y. R. and Zangari, A. J. and Lin, C. L.",
    address = "Preston, United Kingdom",
    booktitle = "Proceedings of a technical committee meeting organized by the IAEA: Water Reactor Fuel Element Computer Modelling in Steady State, Transient, and Accident Conditions",
    month = "September 18-22,",
    title = "{Modeling of PCI Under Steady State and Transient Operating Conditions}",
    year = "1988"
}

J. T. A. Roberts, R. L. Jones, D. Cubicciotti, A. K. Miller, H. F. Wachob, E. Smith, and F. L. Yaggee. A stress corrosion cracking model for pellet-cladding interaction failures in light-water reactor fuel rods. In Zirconium in the Nuclear Industry: Fourth Conference, 285–305. ASTM STP 681, American Society for Testing and Materials, 1979, 1979.

@inproceedings{roberts1979,
    author = "Roberts, J. T. A. and Jones, R. L. and Cubicciotti, D. and Miller, A. K. and Wachob, H. F. and Smith, E. and Yaggee, F. L.",
    title = "A Stress Corrosion Cracking Model for Pellet-Cladding Interaction Failures in Light-Water Reactor Fuel Rods",
    booktitle = "Zirconium in the Nuclear Industry: Fourth Conference",
    year = "1979",
    pages = "285-305",
    publisher = "ASTM STP 681, American Society for Testing and Materials, 1979"
}

(test/tests/cumulative_damage_index/cumulative_damage_index_test.i)

# This tests the AuxKernel CumulativeDamageIndex.  This is a work in progress.
# until we decide what the value of CDI should be.
# The test consists of a ring with an internal and external pressure.
# The delta pressure (inside - outside) is ramped until the hoop stress
# exceeds a reference stress.  When hoop stress > reference stress and
# bunup > 5000 MWd/MTU, the equation for CDI is evaluated and accumulated.
# See the following reference for the CDI equations:
# Y. R. Rashid, A. J. Zangari, C. L. Lin, Water reactor fuel element computer modelling in steady state, transient and accident conditions, in: Proceedings of a technical committee meeting organized by the International Atomic Energy Agency, Preston, 18-22 September, 1988.
# CHANGES
# tkl, 10/17/12 fixed the input mesh so that a pressure boundary condition can be applied to the water side of the clad
# added Sideset 3 - the water side of the clad
# added Sideset 4 - the 'top' of the element
# added a pressure boundary condition to Sideset 3
#
# Updates as of 7/15/2015. Reactivated the test. Re-golded with execute_on = timestep_end for all auxkernels, which gave results consistent with excel calculations (see file).
#
# Note that it seems wrong to use burnup (a fuel property) to calculate damage in the cladding. That's why a post processor was chosen to take an average of the burnup in all the fuel, and use that as input. The regression test was disabled because of an inconsistency of post processing execution. Specifying execute_on = timestep_end should fix the regression test issue and allow further development of this model. Also added more postprocessors to output.
[GlobalParams]
  displacements = 'disp_x disp_y'
  density = 10431.
  temperature = temp
  volumetric_locking_correction = true
[]

[Mesh]
  coord_type = RZ
  [mesh]
    type = FileMeshGenerator
    file = single_element.e     # this one has a sideset (3) so we can apply a pressure BC
  []
[]

[Variables]
  [temp]
    order = FIRST
    family = LAGRANGE
    initial_condition = 600
  []
[]

[AuxVariables]
  [cumulative_damage_index]
    order = CONSTANT
    family = MONOMIAL
    block = 1
  []
  [burnup]
    order = FIRST
    family = LAGRANGE
  []
  [fission_rate]
    order = FIRST
    family = LAGRANGE
  []
[]

[Functions]
  [pressure_ramp_coolant]
    type = PiecewiseLinear
    x = '0 7.5e6 1.5e7'
    y = '0 1     0.1'
  []
  [pressure_ramp_internal]
    type = PiecewiseLinear
    x = '0 7.5e6 1.5e7'
    y = '0 1     2'
    scale_factor = 1
  []
[]

[Physics/SolidMechanics/QuasiStatic]
  [solid]
    block = 1
    strain = FINITE
    add_variables = true
    decomposition_method = EigenSolution
    generate_output = 'stress_zz strain_zz'
  []
[]

[Kernels]
  [heat]
    type = HeatConduction
    variable = temp
  []
  [heat_ie]
    type = HeatConductionTimeDerivative
    variable = temp
  []
[]

[AuxKernels]
  [cumulative_damage_index]
    type = CumulativeDamageIndex
    variable = cumulative_damage_index
    burnup = burnup
    temperature = temp
    hoop_stress = stress_zz
    yield_stress = 240e6
    clad_type = 'Zr2'
    block = 1
    execute_on = timestep_end
  []
  [burnup]
    type = BurnupAux
    variable = burnup
    fission_rate = fission_rate
    execute_on = timestep_end
  []
  [fission_rate]
    type = ConstantAux
    variable = fission_rate
    value = 1e19
    execute_on = timestep_end
  []
[]


[BCs]
  [no_y_clad_bottom]
    type = DirichletBC
    variable = disp_y
    boundary = '2'
    value = 0.0
  []
  [Pressure]
    [internalPressure]
      boundary = 1
      factor = 11.0e6
      function = pressure_ramp_internal
    []
    [coolantPressure]
      boundary = '3'                # Sideset 1 is the inside radius, 2 is the bottom, 3 is the outside clad radius, and 4 is the top of the single element, respectively
      factor = 15.5e6                   # this is PWR type pressure, check to see if approp. for 27(1), is the specified num. for 27(2c)
      function = pressure_ramp_coolant          # use the pressure_ramp function defined above, goes from something to full pressure in ~ 3hr
    []
  []
  [temp_fix]
    type = DirichletBC
    variable = temp
    boundary = '1'
    value = 600
  []
[]

[Materials]
  [elasticity_tensor]
    type = ComputeIsotropicElasticityTensor
    block = 1
    youngs_modulus = 2.e11
    poissons_ratio = 0.3
  []
  [elastic_stress]
    type = ComputeFiniteStrainElasticStress
    block = 1
  []
  [thermal]
    type = HeatConductionMaterial
    block = 1
    specific_heat = 1.0
    thermal_conductivity = 100.
  []
  [density]
    type = StrainAdjustedDensity
    block = 1
    strain_free_density = 6550
  []
[]

[Executioner]
  type = Transient
  solve_type = NEWTON

  petsc_options = '-snes_ksp_ew'
  petsc_options_iname = '-ksp_gmres_restart'
  petsc_options_value = '101'


  line_search = 'none'


  l_max_its = 100
  nl_max_its = 100
  nl_rel_tol = 1e-8
  nl_abs_tol = 1e-5
  l_tol = 1e-5
  start_time = 0.0
  num_steps = 11
  dt = 2.5e6
[]

[Outputs]
  file_base = cumulative_damage_index_test_out
  hide = 'burnup fission_rate'
  [exodus]
    type = Exodus
  []
[]

[Postprocessors]
  [burnup]
    type = ElementAverageValue
    variable = burnup
  []
  [hoop_stress]
    type = ElementAverageValue
    variable = stress_zz
  []
  [CDI]
    type = ElementAverageValue
    variable = cumulative_damage_index
  []
  [Internal_pressure]
    type = FunctionValuePostprocessor
    function = pressure_ramp_internal
    scale_factor = 11e6
  []
  [External_pressure]
    type = FunctionValuePostprocessor
    function = pressure_ramp_coolant
    scale_factor = 15.5e6
  []
[]

(test/tests/cumulative_damage_index/cumulative_damage_index_astm_model.i)

# This tests the AuxKernel CumulativeDamageIndex.
#
# The test consists of a ring with an internal and external pressure.
# The delta pressure (inside - outside) is ramped until the hoop stress
# exceeds a reference stress.  When hoop stress > reference stress and
# bunup > 5000 MWd/MTU, the equation for CDI is evaluated and accumulated.

# Sideset 1 is the inside radius
# Sideset 2 is the bottom
# Sideset 3 - the water side of the clad
# Sideset 4 - the 'top' of the element

[GlobalParams]
  density = 10431.
  displacements = 'disp_x disp_y'
[]

[Mesh]
  coord_type = RZ
  use_displaced_mesh = false
  [mesh]
    type = FileMeshGenerator
    file = single_element.e     # this one has a sideset (3) so we can apply a pressure BC
  []
[]

[Variables]
  [disp_x]
    order = FIRST
    family = LAGRANGE
  []
  [disp_y]
    order = FIRST
    family = LAGRANGE
  []
  [temp]
    order = FIRST
    family = LAGRANGE
    initial_condition = 600
  []

[]

[AuxVariables]
  [cumulative_damage_index]
    order = CONSTANT
    family = MONOMIAL
    block = 1
  []
  [burnup]
    order = FIRST
    family = LAGRANGE
  []
  [fission_rate]
    order = FIRST
    family = LAGRANGE
  []
[]

[Functions]
  [pressure_ramp_coolant]
    type = PiecewiseLinear
    x = '0 1000 86400'
    y = '0 1        1'
  []
  [pressure_ramp_internal]
    type = PiecewiseLinear
    x = '0 1000 86400'
    y = '0 1     3'
    scale_factor = 1
  []
[]

[Physics/SolidMechanics/QuasiStatic]
  [clad]
    strain = SMALL
    incremental = true
    generate_output = 'stress_zz strain_zz'
  []
[]


[Kernels]
  [heat]
    type = HeatConduction
    variable = temp
  []
  [heat_ie]
    type = HeatConductionTimeDerivative
    variable = temp
  []
[]

[AuxKernels]
  [cumulative_damage_index]
    type = CumulativeDamageIndex
    variable = cumulative_damage_index
    temperature = temp
    hoop_stress = stress_zz
    yield_stress = 620e6
    model_option = ASTM
    block = 1
    execute_on = timestep_end
  []
  [burnup]
    type = BurnupAux
    variable = burnup
    fission_rate = fission_rate
    execute_on = timestep_end
  []
  [fission_rate]
    type = ConstantAux
    variable = fission_rate
    value = 1e19
    execute_on = timestep_end
  []
[]


[BCs]
  [no_y_clad_bottom] # pin clad bottom in the axial direction (y)
    type = DirichletBC
    variable = disp_y
    boundary = '2'
    value = 0.0
  []
  [Pressure] #  apply rod plenum pressure on clad inner walls
    [internalPressure]
      boundary = 1
      factor = 20.0e6
      function = pressure_ramp_internal   # use the pressure_ramp function defined above
    []
    [coolantPressure]
       boundary = '3'
       factor = 15.5e6
       function = pressure_ramp_coolant
    []
  []
  [temp_fix]
    type = DirichletBC
    variable = temp
    boundary = '1'
    value = 600
  []
[]

[Materials]
  [elasticity_tensor]
    type = ComputeIsotropicElasticityTensor
    youngs_modulus = 0.9e11
    poissons_ratio = 0.3
  []
  [elastic_stress]
    type = ComputeFiniteStrainElasticStress
  []
  [thermal]
    type = HeatConductionMaterial
    block = 1
    specific_heat = 1.0
    thermal_conductivity = 100.
  []
  [density]
    type = StrainAdjustedDensity
    block = 1
    strain_free_density = 1.0
  []
[]

[Executioner]
  type = Transient
# petsc_options = '-snes_mf_operator -ksp_monitor -snes_ksp_ew'

  solve_type = PJFNK

  petsc_options = '-snes_ksp_ew'
  petsc_options_iname = '-ksp_gmres_restart'
  petsc_options_value = '101'

  line_search = 'none'

  l_max_its = 100
  nl_max_its = 100
  nl_rel_tol = 1e-8
  nl_abs_tol = 1e-5
  l_tol = 1e-5
  start_time = 0.0
  num_steps = 150
  dt = 300
[]

[Outputs]
  hide = 'burnup fission_rate'
  csv = true
  [exodus]
    type = Exodus
  []
[]

[Postprocessors]
  [burnup]
    type = ElementAverageValue
    variable = burnup
  []
  [hoop_stress]
    type = ElementAverageValue
    variable = stress_zz
  []
  [CDI]
    type = ElementAverageValue
    variable = cumulative_damage_index
  []
  [Internal_pressure]
    type = FunctionValuePostprocessor
    function = pressure_ramp_internal
    scale_factor = 11.0e6
  []
  [External_pressure]
    type = FunctionValuePostprocessor
    function = pressure_ramp_coolant
    scale_factor = 15.5e6
  []
[]

(test/tests/cumulative_damage_index/cumulative_damage_index_yield_prop_test.i)

# This tests the AuxKernel CumulativeDamageIndex.  This is a work in progress.
# until we decide what the value of CDI should be.
# The test consists of a ring with an internal and external pressure.
# The delta pressure (inside - outside) is ramped until the hoop stress
# exceeds a reference stress.  When hoop stress > reference stress and
# bunup > 5000 MWd/MTU, the equation for CDI is evaluated and accumulated.
# See the following reference for the CDI equations:
# Y. R. Rashid, A. J. Zangari, C. L. Lin, Water reactor fuel element computer modelling in steady state, transient and accident conditions, in: Proceedings of a technical committee meeting organized by the International Atomic Energy Agency, Preston, 18-22 September, 1988.
# CHANGES
# tkl, 10/17/12 fixed the input mesh so that a pressure boundary condition can be applied to the water side of the clad
# added Sideset 3 - the water side of the clad
# added Sideset 4 - the 'top' of the element
# added a pressure boundary condition to Sideset 3
#
# Updates as of 7/15/2015. Reactivated the test. Re-golded with execute_on = timestep_end for all auxkernels, which gave results consistent with excel calculations (see file).
# Note that it seems wrong to use burnup (a fuel property) to calculate damage in the cladding. That's why a post processor was chosen to take an average of the burnup in all the fuel, and use that as input. The regression test was disabled because of an inconsistency of post processing execution. Specifying execute_on = timestep_end should fix the regression test issue and allow further development of this model. Also added more postprocessors to output.

[GlobalParams]
  displacements = 'disp_x disp_y'
  density = 10431.
  temperature = temp
  volumetric_locking_correction = true
[]

[Mesh]
  coord_type = RZ
  [mesh]
    type = FileMeshGenerator
    file = single_element.e     # this one has a sideset (3) so we can apply a pressure BC
  []
[]

[Variables]
  [temp]
    order = FIRST
    family = LAGRANGE
    initial_condition = 600
  []
[]

[AuxVariables]
  [cumulative_damage_index]
    order = CONSTANT
    family = MONOMIAL
    block = 1
  []
  [burnup]
    order = FIRST
    family = LAGRANGE
  []
  [fission_rate]
    order = FIRST
    family = LAGRANGE
  []
[]

[Functions]
  [pressure_ramp_coolant]
    type = PiecewiseLinear
    x = '0 7.5e6 1.5e7'
    y = '0 1     0.1'
  []
  [pressure_ramp_internal]
    type = PiecewiseLinear
    x = '0 7.5e6 1.5e7'
    y = '0 1     2'
    scale_factor = 1
  []
[]

[Physics/SolidMechanics/QuasiStatic]
  [solid]
    block = 1
    strain = FINITE
    add_variables = true
    decomposition_method = EigenSolution
    generate_output = 'stress_zz strain_zz'
  []
[]

[Kernels]
  [heat]
    type = HeatConduction
    variable = temp
  []
  [heat_ie]
    type = HeatConductionTimeDerivative
    variable = temp
  []
[]

[AuxKernels]
  [cumulative_damage_index]
    type = CumulativeDamageIndex
    variable = cumulative_damage_index
    burnup = burnup
    temperature = temp
    hoop_stress = stress_zz
    yield_stress_property = yield_stress
    clad_type = 'Zr2'
    block = 1
    execute_on = timestep_end
  []
  [burnup]
    type = BurnupAux
    variable = burnup
    fission_rate = fission_rate
    execute_on = timestep_end
  []
  [fission_rate]
    type = ConstantAux
    variable = fission_rate
    value = 1e19
    execute_on = timestep_end
  []
[]


[BCs]
  [no_y_clad_bottom]
    type = DirichletBC
    variable = disp_y
    boundary = '2'
    value = 0.0
  []
  [Pressure]
    [internalPressure]
      boundary = 1
      factor = 11.0e6
      function = pressure_ramp_internal   # use the pressure_ramp function defined above
    []
    [coolantPressure]
       boundary = '3'                # Sideset 1 is the inside radius, 2 is the bottom, 3 is the outside clad radius, and 4 is the top of the single element, respectively
       factor = 15.5e6                   # this is PWR type pressure, check to see if approp. for 27(1), is the specified num. for 27(2c)
       function = pressure_ramp_coolant          # use the pressure_ramp function defined above, goes from something to full pressure in ~ 3hr
    []
  []
  [temp_fix]
    type = DirichletBC
    variable = temp
    boundary = '1'
    value = 600
  []
[]

[Materials]
  [elasticity_tensor]
    type = ComputeIsotropicElasticityTensor
    block = 1
    youngs_modulus = 2.e11
    poissons_ratio = 0.3
  []
  [elastic_stress]
    type = ComputeFiniteStrainElasticStress
    block = 1
  []
  [thermal]
    type = HeatConductionMaterial
    block = 1
    specific_heat = 1.0
    thermal_conductivity = 100.
  []
  [density]
    type = StrainAdjustedDensity
    block = 1
    strain_free_density = 6550
  []
  [yield_stress]
    type = GenericConstantMaterial
    block = 1
    prop_names = 'yield_stress'
    prop_values = '240e6'
  []
[]

[Executioner]
  type = Transient

  solve_type = NEWTON

  petsc_options = '-snes_ksp_ew'
  petsc_options_iname = '-ksp_gmres_restart'
  petsc_options_value = '101'


  line_search = 'none'


  l_max_its = 100
  nl_max_its = 100
  nl_rel_tol = 1e-8
  nl_abs_tol = 1e-5
  l_tol = 1e-5
  start_time = 0.0

  num_steps = 11

  dt = 2.5e6
[]

[Outputs]
  file_base = cumulative_damage_index_yield_prop_test_out
  hide = 'burnup fission_rate'
  [exodus]
    type = Exodus
  []
[]

[Postprocessors]
  [burnup]
    type = ElementAverageValue
    variable = burnup
  []
  [hoop_stress]
    type = ElementAverageValue
    variable = stress_zz
  []
  [CDI]
    type = ElementAverageValue
    variable = cumulative_damage_index
  []
  [Internal_pressure]
    type = FunctionValuePostprocessor
    function = pressure_ramp_internal
    scale_factor = 11e6
  []
  [External_pressure]
    type = FunctionValuePostprocessor
    function = pressure_ramp_coolant
    scale_factor = 15.5e6
  []
[]

(test/tests/cumulative_damage_index/cumulative_damage_index_test.i)

# This tests the AuxKernel CumulativeDamageIndex.  This is a work in progress.
# until we decide what the value of CDI should be.
# The test consists of a ring with an internal and external pressure.
# The delta pressure (inside - outside) is ramped until the hoop stress
# exceeds a reference stress.  When hoop stress > reference stress and
# bunup > 5000 MWd/MTU, the equation for CDI is evaluated and accumulated.
# See the following reference for the CDI equations:
# Y. R. Rashid, A. J. Zangari, C. L. Lin, Water reactor fuel element computer modelling in steady state, transient and accident conditions, in: Proceedings of a technical committee meeting organized by the International Atomic Energy Agency, Preston, 18-22 September, 1988.
# CHANGES
# tkl, 10/17/12 fixed the input mesh so that a pressure boundary condition can be applied to the water side of the clad
# added Sideset 3 - the water side of the clad
# added Sideset 4 - the 'top' of the element
# added a pressure boundary condition to Sideset 3
#
# Updates as of 7/15/2015. Reactivated the test. Re-golded with execute_on = timestep_end for all auxkernels, which gave results consistent with excel calculations (see file).
#
# Note that it seems wrong to use burnup (a fuel property) to calculate damage in the cladding. That's why a post processor was chosen to take an average of the burnup in all the fuel, and use that as input. The regression test was disabled because of an inconsistency of post processing execution. Specifying execute_on = timestep_end should fix the regression test issue and allow further development of this model. Also added more postprocessors to output.
[GlobalParams]
  displacements = 'disp_x disp_y'
  density = 10431.
  temperature = temp
  volumetric_locking_correction = true
[]

[Mesh]
  coord_type = RZ
  [mesh]
    type = FileMeshGenerator
    file = single_element.e     # this one has a sideset (3) so we can apply a pressure BC
  []
[]

[Variables]
  [temp]
    order = FIRST
    family = LAGRANGE
    initial_condition = 600
  []
[]

[AuxVariables]
  [cumulative_damage_index]
    order = CONSTANT
    family = MONOMIAL
    block = 1
  []
  [burnup]
    order = FIRST
    family = LAGRANGE
  []
  [fission_rate]
    order = FIRST
    family = LAGRANGE
  []
[]

[Functions]
  [pressure_ramp_coolant]
    type = PiecewiseLinear
    x = '0 7.5e6 1.5e7'
    y = '0 1     0.1'
  []
  [pressure_ramp_internal]
    type = PiecewiseLinear
    x = '0 7.5e6 1.5e7'
    y = '0 1     2'
    scale_factor = 1
  []
[]

[Physics/SolidMechanics/QuasiStatic]
  [solid]
    block = 1
    strain = FINITE
    add_variables = true
    decomposition_method = EigenSolution
    generate_output = 'stress_zz strain_zz'
  []
[]

[Kernels]
  [heat]
    type = HeatConduction
    variable = temp
  []
  [heat_ie]
    type = HeatConductionTimeDerivative
    variable = temp
  []
[]

[AuxKernels]
  [cumulative_damage_index]
    type = CumulativeDamageIndex
    variable = cumulative_damage_index
    burnup = burnup
    temperature = temp
    hoop_stress = stress_zz
    yield_stress = 240e6
    clad_type = 'Zr2'
    block = 1
    execute_on = timestep_end
  []
  [burnup]
    type = BurnupAux
    variable = burnup
    fission_rate = fission_rate
    execute_on = timestep_end
  []
  [fission_rate]
    type = ConstantAux
    variable = fission_rate
    value = 1e19
    execute_on = timestep_end
  []
[]


[BCs]
  [no_y_clad_bottom]
    type = DirichletBC
    variable = disp_y
    boundary = '2'
    value = 0.0
  []
  [Pressure]
    [internalPressure]
      boundary = 1
      factor = 11.0e6
      function = pressure_ramp_internal
    []
    [coolantPressure]
      boundary = '3'                # Sideset 1 is the inside radius, 2 is the bottom, 3 is the outside clad radius, and 4 is the top of the single element, respectively
      factor = 15.5e6                   # this is PWR type pressure, check to see if approp. for 27(1), is the specified num. for 27(2c)
      function = pressure_ramp_coolant          # use the pressure_ramp function defined above, goes from something to full pressure in ~ 3hr
    []
  []
  [temp_fix]
    type = DirichletBC
    variable = temp
    boundary = '1'
    value = 600
  []
[]

[Materials]
  [elasticity_tensor]
    type = ComputeIsotropicElasticityTensor
    block = 1
    youngs_modulus = 2.e11
    poissons_ratio = 0.3
  []
  [elastic_stress]
    type = ComputeFiniteStrainElasticStress
    block = 1
  []
  [thermal]
    type = HeatConductionMaterial
    block = 1
    specific_heat = 1.0
    thermal_conductivity = 100.
  []
  [density]
    type = StrainAdjustedDensity
    block = 1
    strain_free_density = 6550
  []
[]

[Executioner]
  type = Transient
  solve_type = NEWTON

  petsc_options = '-snes_ksp_ew'
  petsc_options_iname = '-ksp_gmres_restart'
  petsc_options_value = '101'


  line_search = 'none'


  l_max_its = 100
  nl_max_its = 100
  nl_rel_tol = 1e-8
  nl_abs_tol = 1e-5
  l_tol = 1e-5
  start_time = 0.0
  num_steps = 11
  dt = 2.5e6
[]

[Outputs]
  file_base = cumulative_damage_index_test_out
  hide = 'burnup fission_rate'
  [exodus]
    type = Exodus
  []
[]

[Postprocessors]
  [burnup]
    type = ElementAverageValue
    variable = burnup
  []
  [hoop_stress]
    type = ElementAverageValue
    variable = stress_zz
  []
  [CDI]
    type = ElementAverageValue
    variable = cumulative_damage_index
  []
  [Internal_pressure]
    type = FunctionValuePostprocessor
    function = pressure_ramp_internal
    scale_factor = 11e6
  []
  [External_pressure]
    type = FunctionValuePostprocessor
    function = pressure_ramp_coolant
    scale_factor = 15.5e6
  []
[]

Description
Example Input Syntax
Input Parameters
Input Files
References