SS316CreepUpdate

Thermal and irradiation creep for SS AISI 316.

Description

This material, SS316CreepUpdate, computes the steady state creep, the sum of thermal and irradiation creep, using Eq. (1), following Altenbach and Gorash (2013) and Garner and Porter (1988) : $(1)var element = document.getElementById("moose-equation-27368333-ed2c-4540-9f86-f1d137d1061f");katex.render(" \\dot{\\epsilon}_{cr} = a \\exp \\left(\\frac{-Q_{ln}}{RT}\\right) \\sinh \\left[b \\exp \\left(\\frac{-Q_{pw}}{RT}\\right)\\sigma\\right] + (B + D \\dot{S}) \\phi \\sigma", element, {displayMode:true,throwOnError:false});$ where $var element = document.getElementById("moose-equation-73a1fa98-c0d6-4d3a-9ef1-ea9f120ad924");katex.render("\\phi", element, {displayMode:false,throwOnError:false});$ is the neutron flux (n/cm $var element = document.getElementById("moose-equation-383d7f99-9beb-4e8a-aa92-0d058e765280");katex.render("^2", element, {displayMode:false,throwOnError:false});$ -s), $var element = document.getElementById("moose-equation-8b61cef0-08b6-4f8f-9d35-da93e78650e5");katex.render("T", element, {displayMode:false,throwOnError:false});$ is the temperature (K), $var element = document.getElementById("moose-equation-1e922e80-43ac-4dbe-ab14-51f48a55d198");katex.render("\\sigma", element, {displayMode:false,throwOnError:false});$ is the effective stress (MPa), and $var element = document.getElementById("moose-equation-56426c6e-debe-4fa1-998f-b9905c0dbc7e");katex.render("\\dot{\\epsilon}_{cr}", element, {displayMode:false,throwOnError:false});$ is the effective thermal and irradiation creep strain rate (1/s). The parameters for the model are given in Table 1.

Table 1: Parameters for the SS316 Creep Model

Parameters	Values	Unit
$var element = document.getElementById("moose-equation-ce3a88ee-8d0e-410a-aba6-817082e51319");katex.render("a", element, {displayMode:false,throwOnError:false});$	75776	$var element = document.getElementById("moose-equation-ae437dfc-8d79-4484-baee-6db4b6752621");katex.render("h^{-1}", element, {displayMode:false,throwOnError:false});$
$var element = document.getElementById("moose-equation-5ca23e45-074b-42c1-a85f-a47c582b846c");katex.render("b", element, {displayMode:false,throwOnError:false});$	0.7	MPa $var element = document.getElementById("moose-equation-7a571a44-66b1-48d6-8f7b-0c14f39d8847");katex.render("^{-1}", element, {displayMode:false,throwOnError:false});$
$var element = document.getElementById("moose-equation-d3f14cd2-2050-400e-bc33-255e827b7ba8");katex.render("Q_{ln}", element, {displayMode:false,throwOnError:false});$	170000	J/mol
$var element = document.getElementById("moose-equation-2b89269c-d62a-44ba-a334-eadd480860c6");katex.render("Q_{pw}", element, {displayMode:false,throwOnError:false});$	12000	J/mol
$var element = document.getElementById("moose-equation-101457eb-171d-4189-ad3e-476fbec3a2e4");katex.render("B", element, {displayMode:false,throwOnError:false});$	1 $var element = document.getElementById("moose-equation-f554f14a-5233-4844-8204-6f6a901b8460");katex.render("\\times10^{-6}", element, {displayMode:false,throwOnError:false});$	MPa $var element = document.getElementById("moose-equation-6c93ab36-39bc-4431-b0ee-27a131b463ff");katex.render("^{-1}", element, {displayMode:false,throwOnError:false});$ dpa $var element = document.getElementById("moose-equation-d6d24ec2-328d-4ade-ba47-ac0cdbc11b70");katex.render("^{-1}", element, {displayMode:false,throwOnError:false});$
$var element = document.getElementById("moose-equation-79b8299b-d7e0-45e9-bb22-6ef8752c35c6");katex.render("D", element, {displayMode:false,throwOnError:false});$	0.6 $var element = document.getElementById("moose-equation-001eb713-c617-4bf0-aa4a-84b0787160a9");katex.render("\\times10^{-6}", element, {displayMode:false,throwOnError:false});$	MPa $var element = document.getElementById("moose-equation-ff625b74-86bd-45f8-911f-e3cce77de758");katex.render("^{-1}", element, {displayMode:false,throwOnError:false});$
$var element = document.getElementById("moose-equation-d7240ffa-c765-4727-8005-8be0def9e4f1");katex.render("\\dot{S}", element, {displayMode:false,throwOnError:false});$	0.01

Example Input Syntax

[Materials<<<{"href": "../../../syntax/Materials/index.html"}>>>]
  [creep]
    type = SS316CreepUpdate<<<{"description": "Thermal and irradiation creep for SS AISI 316.", "href": "SS316CreepUpdate.html"}>>>
    fast_neutron_flux<<<{"description": "Fast neutron flux in neutrons/m^2s"}>>> = fast_neutron_flux
  []
[]

SS316CreepUpdate must be run in conjunction with the inelastic strain return mapping stress calculator as shown below:

[Materials<<<{"href": "../../../syntax/Materials/index.html"}>>>]
  [radial_return_stress]
    type = ComputeMultipleInelasticStress<<<{"description": "Compute state (stress and internal parameters such as plastic strains and internal parameters) using an iterative process.  Combinations of creep models and plastic models may be used.", "href": "../ComputeMultipleInelasticStress.html"}>>>
    tangent_operator<<<{"description": "Type of tangent operator to return.  'elastic': return the elasticity tensor.  'nonlinear': return the full, general consistent tangent operator."}>>> = elastic
    inelastic_models<<<{"description": "The material objects to use to calculate stress and inelastic strains. Note: specify creep models first and plasticity models second."}>>> = 'creep'
  []
[]

Input Parameters

B1e-06Creep compliance term, considered to be temperature independent for most austenitic stainless steels
Default:1e-06
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:Creep compliance term, considered to be temperature independent for most austenitic stainless steels
D0.006Stress-coupling coefficient, also considered to be relatively constant for austenitic stainless steels
Default:0.006
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:Stress-coupling coefficient, also considered to be relatively constant for austenitic stainless steels
Q_ln170000Activation energy of first thermal exponential
Default:170000
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:Activation energy of first thermal exponential
Q_pw12000Activation energy of second thermal exponential within sinh term
Default:12000
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:Activation energy of second thermal exponential within sinh term
R8.314Universal gas constant
Default:8.314
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:Universal gas constant
S0.01Swelling rate per dpa, constant above approximately 400 C for most austenitic stainless steels
Default:0.01
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:Swelling rate per dpa, constant above approximately 400 C for most austenitic stainless steels
a21.0489Leading term of first thermal exponential
Default:21.0489
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:Leading term of first thermal exponential
absolute_tolerance1e-11Absolute convergence tolerance for Newton iteration
Default:1e-11
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:Absolute convergence tolerance for Newton iteration
acceptable_multiplier10Factor applied to relative and absolute tolerance for acceptable convergence if iterations are no longer making progress
Default:10
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:Factor applied to relative and absolute tolerance for acceptable convergence if iterations are no longer making progress
adaptive_substeppingFalseUse adaptive substepping, where the number of substeps is successively doubled until the return mapping model successfully converges or the maximum number of substeps is reached.
Default:False
C++ Type:bool
Controllable:No
Description:Use adaptive substepping, where the number of substeps is successively doubled until the return mapping model successfully converges or the maximum number of substeps is reached.
automatic_differentiation_return_mappingFalseWhether to use automatic differentiation to compute the derivative.
Default:False
C++ Type:bool
Controllable:No
Description:Whether to use automatic differentiation to compute the derivative.
b0.7Leading term of second thermal exponential within sinh term.
Default:0.7
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:Leading term of second thermal exponential within sinh term.
base_nameOptional parameter that defines a prefix for all material properties related to this stress update model. This allows for multiple models of the same type to be used without naming conflicts.
C++ Type:std::string
Controllable:No
Description:Optional parameter that defines a prefix for all material properties related to this stress update model. This allows for multiple models of the same type to be used without naming conflicts.
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
constant_onNONEWhen ELEMENT, MOOSE will only call computeQpProperties() for the 0th quadrature point, and then copy that value to the other qps.When SUBDOMAIN, MOOSE will only call computeQpProperties() for the 0th quadrature point, and then copy that value to the other qps. Evaluations on element qps will be skipped
Default:NONE
C++ Type:MooseEnum
Options:NONE, ELEMENT, SUBDOMAIN
Controllable:No
Description:When ELEMENT, MOOSE will only call computeQpProperties() for the 0th quadrature point, and then copy that value to the other qps.When SUBDOMAIN, MOOSE will only call computeQpProperties() for the 0th quadrature point, and then copy that value to the other qps. Evaluations on element qps will be skipped
declare_suffixAn optional suffix parameter that can be appended to any declared 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 declared properties. The suffix will be prepended with a '_' character.
fast_neutron_fluxfast_neutron_fluxFast neutron flux in neutrons/m^2s
Default:fast_neutron_flux
C++ Type:MaterialPropertyName
Unit:(no unit assumed)
Controllable:No
Description:Fast neutron flux in neutrons/m^2s
id_wastage_degradation_functionThe optional ID wastage degradation function that takes FCCI effect on cladding into consideration.
C++ Type:FunctionName
Unit:(no unit assumed)
Controllable:No
Description:The optional ID wastage degradation function that takes FCCI effect on cladding into consideration.
max_inelastic_increment0.0001The maximum inelastic strain increment allowed in a time step
Default:0.0001
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:The maximum inelastic strain increment allowed in a time step
maximum_number_substeps25The maximum number of substeps allowed before cutting the time step.
Default:25
C++ Type:unsigned int
Controllable:No
Description:The maximum number of substeps allowed before cutting the time step.
od_wastage_degradation_functionThe optional OD wastage degradation function that takes CCCI effect on cladding into consideration.
C++ Type:FunctionName
Unit:(no unit assumed)
Controllable:No
Description:The optional OD wastage degradation function that takes CCCI effect on cladding into consideration.
relative_tolerance1e-08Relative convergence tolerance for Newton iteration
Default:1e-08
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:Relative convergence tolerance for Newton iteration
temperatureThe coupled temperature
C++ Type:std::vector<VariableName>
Unit:(no unit assumed)
Controllable:No
Description:The coupled temperature
use_substep_integration_errorFalseIf true, it establishes a substep size that will yield, at most,the creep numerical integration error given by substep_strain_tolerance.
Default:False
C++ Type:bool
Controllable:No
Description:If true, it establishes a substep size that will yield, at most,the creep numerical integration error given by substep_strain_tolerance.
use_substeppingNONEWhether and how to use substepping
Default:NONE
C++ Type:MooseEnum
Options:NONE, ERROR_BASED, INCREMENT_BASED
Controllable:No
Description:Whether and how to use substepping

Optional Parameters

apply_strainTrueFlag to apply strain. Used for testing.
Default:True
C++ Type:bool
Controllable:No
Description:Flag to apply strain. Used for testing.
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.
effective_inelastic_strain_nameeffective_creep_strainName of the material property that stores the effective inelastic strain
Default:effective_creep_strain
C++ Type:std::string
Controllable:No
Description:Name of the material property that stores the effective inelastic strain
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
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
substep_strain_tolerance0.1Maximum ratio of the initial elastic strain increment at start of the return mapping solve to the maximum inelastic strain allowable in a single substep. Reduce this value to increase the number of substeps
Default:0.1
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:Maximum ratio of the initial elastic strain increment at start of the return mapping solve to the maximum inelastic strain allowable in a single substep. Reduce this value to increase the number of substeps
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

internal_solve_full_iteration_historyFalseSet true to output full internal Newton iteration history at times determined by `internal_solve_output_on`. If false, only a summary is output.
Default:False
C++ Type:bool
Controllable:No
Description:Set true to output full internal Newton iteration history at times determined by `internal_solve_output_on`. If false, only a summary is output.
internal_solve_output_onon_errorWhen to output internal Newton solve information
Default:on_error
C++ Type:MooseEnum
Options:never, on_error, always
Controllable:No
Description:When to output internal Newton solve information

Debug Parameters

output_propertiesList of material properties, from this material, to output (outputs must also be defined to an output type)
C++ Type:std::vector<std::string>
Controllable:No
Description:List of material properties, from this material, to output (outputs must also be defined to an output type)
outputsnone Vector of output names where you would like to restrict the output of variables(s) associated with this object
Default:none
C++ Type:std::vector<OutputName>
Controllable:No
Description:Vector of output names where you would like to restrict the output of variables(s) associated with this object

Outputs 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

(assessment/metallic_fuel/FBTA/analysis/steady_state/common_inputs/ss316_clad_base.i)
(assessment/MOX/JOYO/MK-I/analysis/MK-I_75MW_master_new_bubble_gb_lim.i)
(test/tests/solid_mechanics/SS316_creep/creep_SS316_rz_therm_only.i)
(assessment/MOX/JOYO/MK-I/analysis/MK-I_50MW_master_old_bubble_gb_lim.i)
(assessment/MOX/JOYO/MK-II/analysis/MK-II_master_old_bubble_gb_lim.i)
(assessment/MOX/JOYO/MK-I/analysis/MK-I_50MW_master_new_bubble_gb_lim.i)
(assessment/MOX/JOYO/B14/PTM003/analysis/b14_ptm003_2DRZ_t.i)
(assessment/nitride/JOYO/JOYO_L4C4/analysis/JOYO_Pin_base.i)
(assessment/metallic_fuel/EBRII/X423/analysis/x423_lm_base.i)
(test/tests/solid_mechanics/SS316_creep/creep_SS316_rz_test.i)
(assessment/MOX/JOYO/B14/PTM002/analysis/b14_ptm002_2DRZ_t.i)
(assessment/MOX/JOYO/B14/PTM010/analysis/b14_ptm010_2DRZ_t.i)
(assessment/MOX/JOYO/B14/PTM001/analysis/b14_ptm001_2DRZ_t.i)
(assessment/MOX/JOYO/MK-I/analysis/MK-I_75MW_master_old_bubble_gb_lim.i)
(test/tests/solid_mechanics/SS316_creep/creep_SS316_3d_test.i)
(assessment/MOX/JOYO/MK-II/analysis/MK-II_master_new_bubble_gb_lim.i)

References

Holm Altenbach and Yevgen Gorash. High-temperature inelastic behavior of the austenitic steel aisi type 316. In Advanced Materials Modelling for Structures, pages 17–30. Springer, 2013.

@incollection{altenbach2013,
    author = "Altenbach, Holm and Gorash, Yevgen",
    title = "High-temperature inelastic behavior of the austenitic steel AISI type 316",
    booktitle = "Advanced Materials Modelling for Structures",
    pages = "17-30",
    year = "2013",
    publisher = "Springer"
}

FA Garner and DL Porter. Irradiation creep and swelling of aisi 316 to exposures of 130 dpa at 385-400 c. Journal of Nuclear Materials, 155:1006–1013, 1988.

@article{garner1988,
    author = "Garner, FA and Porter, DL",
    title = "Irradiation creep and swelling of AISI 316 to exposures of 130 dpa at 385-400 C",
    journal = "Journal of Nuclear Materials",
    volume = "155",
    pages = "1006-1013",
    year = "1988",
    publisher = "Elsevier"
}

(test/tests/solid_mechanics/SS316_creep/creep_SS316_3d_test.i)

#
# The mesh is a 1x1x1 cube with a pressure of 10 MPa on the top face.
# Symmetry boundary conditions on three places provide a uniaxial stress field.
# The temperature is held constant at 973.15.  The yield stress is set at twice
# the load pressure. The solution is advanced through ten time steps of 1e4 for
# a total time of 1e5.
#
# The total strain at time 1e5 can be computed as given by "High Temperature Inelastic Behavior of the AUstenitic Steel AISI Type 316" by Altenbach and Gorash and "Irradiatian Creep and Swelling of AISI 316 to Exposures of 130 dpa at 385-400 degrees C" by Garner and Porter.
#
# e_tot = e_elas +
#         e_thermal_creep +
#         e_irradiation_creep
#
#       = P/E +
#         (a * exp(-Q_ln/(RT)) * sinh(b * exp(-Q_pw/(RT)) * sigma) +
#         [B + DS] * fast_neutron_flux * sigma) * dt
#
#         where P = pressure load
#               E = Youngs modulus
#               a,b = thermal creep coefficients
#               sigma = stress in MPa
#               Q_ln-Q_pw = activation energies
#               B = Irradiation creep coefficient "creep compliance"
#               D = Irradiation creep coefficent "coupling coefficient"
#               S = Swelling Rate
#               dt = total time
#               fast_neutron_flux = 3e17 n/m**2/s (applied) but model expects units of 10**22 n/cm**2/s
#               so, fast_neutron_flux = 3e17 * 0.0001 (convert from cm**2 to m**2) * 1e-22 (convert to 10*22)= 3e-9
#
#               Note that these strains are expressed as %, so they have to be divided by 100, thus  the term 0.01 in the strain calculations found below.
#
#
#
#
#   For this test, the analytical solution at t=1e5 is:
#
#               e_tot = 10e6/1.90e11 +
#                       ((21.04889 * exp(-170000/(8.314*973.15)) *
#                       sinh(0.07 * (10) * exp(-12000/(8.314*973.15)))) +
#                       [1e-6 + 0.6e-2 * 0.01] * 3e-9 * 10 * 0.01) * 1e5
#
#                     = 5.263e-5 + (((4.336e-12 * 1.012) + 1.83e-12 * 0.01) * 1e5)
#                     = 5.263e-5 + 4.386e-7 + 1.83e-9 = 5.307e-5
#   e_elas = 5.263e-5
#   e_thermal_creep = 4.386e-7
#   e_irradiation_creep = 1.83e-9
#   e_tot = 5.307e-5
#
#
#


[Mesh]
  use_displaced_mesh = false
  [mesh]
    type = GeneratedMeshGenerator
    dim = 3
  []
[]

[GlobalParams]
  displacements = 'disp_x disp_y disp_z'
  order = FIRST
  family = LAGRANGE
  temperature = temp
[]

[Variables]
  [temp]
    initial_condition = 973.15
  []
[]

[AuxVariables]
  [stress_yy]
    order = CONSTANT
    family = MONOMIAL
  []
  [elastic_strain_yy]
    order = CONSTANT
    family = MONOMIAL
  []
  [creep_strain_yy]
    order = CONSTANT
    family = MONOMIAL
  []
[]

[Physics]
  [SolidMechanics]
    [QuasiStatic]
      [all]
        strain = SMALL
        incremental = true
        add_variables = true
        generate_output = 'stress_yy elastic_strain_yy creep_strain_yy'
      []
    []
  []
[]

[Functions]
  [top_pull]
    type = PiecewiseLinear
    x = '0 1'
    y = '1 1'
  []
[]

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

[AuxKernels]
  [stress_yy]
    type = RankTwoAux
    rank_two_tensor = stress
    variable = stress_yy
    index_i = 1
    index_j = 1
  []
  [elastic_strain_yy]
    type = RankTwoAux
    rank_two_tensor = elastic_strain
    variable = elastic_strain_yy
    index_i = 1
    index_j = 1
  []
  [creep_strain_yy]
    type = RankTwoAux
    rank_two_tensor = creep_strain
    variable = creep_strain_yy
    index_i = 1
    index_j = 1
    execute_on = timestep_end
  []
[]

[BCs]
  [u_top_pull]
    type = Pressure
    variable = disp_y
    boundary = top
    factor = -10.0e6
    function = top_pull
  []
  [u_bottom_fix]
    type = DirichletBC
    variable = disp_y
    boundary = bottom
    value = 0.0
  []
  [u_back_fix]
    type = DirichletBC
    variable = disp_z
    boundary = back
    value = 0.0
  []
  [u_left_fix]
    type = DirichletBC
    variable = disp_x
    boundary = left
    value = 0.0
  []
  [temp_top_fix]
    type = DirichletBC
    variable = temp
    boundary = top
    value = 973.15
  []
  [temp_bottom_fix]
    type = DirichletBC
    variable = temp
    boundary = bottom
    value = 973.15
  []
[]

[Materials]
  [fast_neutron_flux]
    type = FastNeutronFlux
    calculate_fluence = true
    factor = 3e17
  []
  [elasticity_tensor]
    type = ComputeIsotropicElasticityTensor
    youngs_modulus = 1.90e11
    poissons_ratio = 0.265
  []
  [creep]
    type = SS316CreepUpdate
    fast_neutron_flux = fast_neutron_flux
  []
  [radial_return_stress]
    type = ComputeMultipleInelasticStress
    tangent_operator = elastic
    inelastic_models = 'creep'
  []
  [thermal]
    type = HeatConductionMaterial
    thermal_conductivity = 100.0
  []
  [density]
    type = StrainAdjustedDensity
    strain_free_density = 1.0
  []
[]

[Executioner]
  type = Transient

  solve_type = 'PJFNK'

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

  line_search = 'none'

  l_max_its  = 50
  nl_max_its = 20
  l_tol = 1e-5
  start_time = 0.0
  num_steps = 10
  dt = 10000
[]

[Outputs]
  exodus = true
[]

(test/tests/solid_mechanics/SS316_creep/creep_SS316_3d_test.i)

#
# The mesh is a 1x1x1 cube with a pressure of 10 MPa on the top face.
# Symmetry boundary conditions on three places provide a uniaxial stress field.
# The temperature is held constant at 973.15.  The yield stress is set at twice
# the load pressure. The solution is advanced through ten time steps of 1e4 for
# a total time of 1e5.
#
# The total strain at time 1e5 can be computed as given by "High Temperature Inelastic Behavior of the AUstenitic Steel AISI Type 316" by Altenbach and Gorash and "Irradiatian Creep and Swelling of AISI 316 to Exposures of 130 dpa at 385-400 degrees C" by Garner and Porter.
#
# e_tot = e_elas +
#         e_thermal_creep +
#         e_irradiation_creep
#
#       = P/E +
#         (a * exp(-Q_ln/(RT)) * sinh(b * exp(-Q_pw/(RT)) * sigma) +
#         [B + DS] * fast_neutron_flux * sigma) * dt
#
#         where P = pressure load
#               E = Youngs modulus
#               a,b = thermal creep coefficients
#               sigma = stress in MPa
#               Q_ln-Q_pw = activation energies
#               B = Irradiation creep coefficient "creep compliance"
#               D = Irradiation creep coefficent "coupling coefficient"
#               S = Swelling Rate
#               dt = total time
#               fast_neutron_flux = 3e17 n/m**2/s (applied) but model expects units of 10**22 n/cm**2/s
#               so, fast_neutron_flux = 3e17 * 0.0001 (convert from cm**2 to m**2) * 1e-22 (convert to 10*22)= 3e-9
#
#               Note that these strains are expressed as %, so they have to be divided by 100, thus  the term 0.01 in the strain calculations found below.
#
#
#
#
#   For this test, the analytical solution at t=1e5 is:
#
#               e_tot = 10e6/1.90e11 +
#                       ((21.04889 * exp(-170000/(8.314*973.15)) *
#                       sinh(0.07 * (10) * exp(-12000/(8.314*973.15)))) +
#                       [1e-6 + 0.6e-2 * 0.01] * 3e-9 * 10 * 0.01) * 1e5
#
#                     = 5.263e-5 + (((4.336e-12 * 1.012) + 1.83e-12 * 0.01) * 1e5)
#                     = 5.263e-5 + 4.386e-7 + 1.83e-9 = 5.307e-5
#   e_elas = 5.263e-5
#   e_thermal_creep = 4.386e-7
#   e_irradiation_creep = 1.83e-9
#   e_tot = 5.307e-5
#
#
#


[Mesh]
  use_displaced_mesh = false
  [mesh]
    type = GeneratedMeshGenerator
    dim = 3
  []
[]

[GlobalParams]
  displacements = 'disp_x disp_y disp_z'
  order = FIRST
  family = LAGRANGE
  temperature = temp
[]

[Variables]
  [temp]
    initial_condition = 973.15
  []
[]

[AuxVariables]
  [stress_yy]
    order = CONSTANT
    family = MONOMIAL
  []
  [elastic_strain_yy]
    order = CONSTANT
    family = MONOMIAL
  []
  [creep_strain_yy]
    order = CONSTANT
    family = MONOMIAL
  []
[]

[Physics]
  [SolidMechanics]
    [QuasiStatic]
      [all]
        strain = SMALL
        incremental = true
        add_variables = true
        generate_output = 'stress_yy elastic_strain_yy creep_strain_yy'
      []
    []
  []
[]

[Functions]
  [top_pull]
    type = PiecewiseLinear
    x = '0 1'
    y = '1 1'
  []
[]

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

[AuxKernels]
  [stress_yy]
    type = RankTwoAux
    rank_two_tensor = stress
    variable = stress_yy
    index_i = 1
    index_j = 1
  []
  [elastic_strain_yy]
    type = RankTwoAux
    rank_two_tensor = elastic_strain
    variable = elastic_strain_yy
    index_i = 1
    index_j = 1
  []
  [creep_strain_yy]
    type = RankTwoAux
    rank_two_tensor = creep_strain
    variable = creep_strain_yy
    index_i = 1
    index_j = 1
    execute_on = timestep_end
  []
[]

[BCs]
  [u_top_pull]
    type = Pressure
    variable = disp_y
    boundary = top
    factor = -10.0e6
    function = top_pull
  []
  [u_bottom_fix]
    type = DirichletBC
    variable = disp_y
    boundary = bottom
    value = 0.0
  []
  [u_back_fix]
    type = DirichletBC
    variable = disp_z
    boundary = back
    value = 0.0
  []
  [u_left_fix]
    type = DirichletBC
    variable = disp_x
    boundary = left
    value = 0.0
  []
  [temp_top_fix]
    type = DirichletBC
    variable = temp
    boundary = top
    value = 973.15
  []
  [temp_bottom_fix]
    type = DirichletBC
    variable = temp
    boundary = bottom
    value = 973.15
  []
[]

[Materials]
  [fast_neutron_flux]
    type = FastNeutronFlux
    calculate_fluence = true
    factor = 3e17
  []
  [elasticity_tensor]
    type = ComputeIsotropicElasticityTensor
    youngs_modulus = 1.90e11
    poissons_ratio = 0.265
  []
  [creep]
    type = SS316CreepUpdate
    fast_neutron_flux = fast_neutron_flux
  []
  [radial_return_stress]
    type = ComputeMultipleInelasticStress
    tangent_operator = elastic
    inelastic_models = 'creep'
  []
  [thermal]
    type = HeatConductionMaterial
    thermal_conductivity = 100.0
  []
  [density]
    type = StrainAdjustedDensity
    strain_free_density = 1.0
  []
[]

[Executioner]
  type = Transient

  solve_type = 'PJFNK'

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

  line_search = 'none'

  l_max_its  = 50
  nl_max_its = 20
  l_tol = 1e-5
  start_time = 0.0
  num_steps = 10
  dt = 10000
[]

[Outputs]
  exodus = true
[]

(assessment/metallic_fuel/FBTA/analysis/steady_state/common_inputs/ss316_clad_base.i)

# SS316 Cladding
clad_density = 7980.0 # SS316 RT Density (same as D9)


[Physics/SolidMechanics/QuasiStatic]
  [cladding]
    strain = FINITE
    generate_output = 'stress_xx stress_yy stress_zz vonmises_stress hydrostatic_stress creep_strain_xx creep_strain_yy creep_strain_zz elastic_strain_xx elastic_strain_yy elastic_strain_zz strain_xx strain_yy strain_zz'
    extra_vector_tags = 'ref'
    block = '${cladding_block}'
    eigenstrain_names = 'cladding_thermal_eigenstrain clad_swelling'
    use_automatic_differentiation = true
    volumetric_locking_correction = true
  []
[]

[Materials]
  [clad_elasticity_tensor]
    type = ADSS316ElasticityTensor
    temperature = temp
    id_wastage_degradation_function = id_vpp_func
    od_wastage_degradation_function = od_vpp_func
    block = ${cladding_block}
    elastic_constants_model = legacy_ifr
  []
  [clad_stress]
    type = ADComputeMultipleInelasticStress
    # tangent_operator = nonlinear
    inelastic_models = 'clad_creep'
    block = ${cladding_block}
  []
  [clad_creep]
    type = ADSS316CreepUpdate
    block = ${cladding_block}
    temperature = temp
    fast_neutron_flux = fast_neutron_flux
    id_wastage_degradation_function = id_vpp_func
    od_wastage_degradation_function = od_vpp_func
  []
  [thermal_expansion]
    type = ADSS316ThermalExpansionEigenstrain
    block = ${cladding_block}
    temperature = temp
    stress_free_temperature = 295.0
    eigenstrain_name = cladding_thermal_eigenstrain
    outputs = all
  []
  [clad_swelling]
    type = ADSS316VolumetricSwellingEigenstrain
    eigenstrain_name = clad_swelling
    fast_neutron_fluence = fast_neutron_fluence
    fast_neutron_flux = fast_neutron_flux
    temperature = temp
    outputs = all
    block = ${cladding_block}
  []
  [clad_thermal]
    type = ADSS316Thermal
    block = ${cladding_block}
    temperature = temp
  []
  [clad_density]
    type = ADStrainAdjustedDensity
    block = ${cladding_block}
    strain_free_density = ${clad_density}
  []
  [clad_failure]
    type = D9FailureClad
    block = ${cladding_block}
    method = steady_state
    temperature = temp
    outputs = all
    hoop_stress = stress_zz # Since 2D-RZ
  []
  [wastage_thickness]
    type = ADMetallicFuelWastage
    method = flux_ss316
    temperature = temp
    scale_factor = 1
    boundary = cladding_inside_right
    outputs = all
  []
  [cc_wastage_thickness]
    type = ADMetallicFuelCoolantWastage
    clad_material = SS316
    use_effective_method = true
    temperature = temp
    scale_factor = 1
    boundary = cladding_outside_right
    outputs = all
  []
[]

(assessment/MOX/JOYO/MK-I/analysis/MK-I_75MW_master_new_bubble_gb_lim.i)


initial_fuel_density = 10836.8

[GlobalParams]
  density = ${initial_fuel_density}
  initial_porosity = 0.065
  order = SECOND
  family = LAGRANGE
  energy_per_fission = 3.2e-11 # J/fission
  volumetric_locking_correction = false
  displacements = 'disp_x disp_y'
[]
[Problem]
  type = ReferenceResidualProblem
  extra_tag_vectors = 'ref'
  reference_vector = 'ref'
[]
[Mesh]
  coord_type = RZ
  [smeared_pellet_mesh]
    type = FuelPinMeshGenerator
    pellet_quantity = 1
    pellet_height = 0.6
    pellet_outer_radius = 0.0027
    pellet_mesh_density = customize
    clad_mesh_density = customize
    clad_gap_width = 0.000100
    clad_thickness = 0.00035
    clad_bot_gap_height = 1.0e-3
    bottom_clad_height = 2.24e-3
    top_clad_height = 2.24e-3
    clad_top_gap_height = 0.599
    elem_type = QUAD8
    nx_c = 4
    ny_c = 200
    nx_p = 20
    ny_p = 200
    ny_cu = 3
    ny_cl = 3
  []
  patch_size = 50
  patch_update_strategy = iteration
  partitioner = centroid
  centroid_partitioner_direction = y
[]
[UserObjects]
  [pin_geometry]
    type = FuelPinGeometry
  []
[]
[Variables]
  [temp]
    initial_condition = 295.0
  []
[]
[AuxVariables]
  [pore]
  []
  [fission_rate]
    block = pellet
  []
  [burnup]
    block = pellet
  []
  [gas_gen_3]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [gas_grn_3]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [gas_bdr_3]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [gas_rel_3]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [bbl_bdr_2]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [vcn_bdr_2]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [atm_bbl_bdr]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [vcn_bbl_bdr]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [prs_bbl_bdr]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [prseq_bbl_bdr]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [rad_bbl_bdr]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [vol_bbl_bdr]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [GBCoverage]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [eff_diff_coeff]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [deltav_v0_bd]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [radial_strain]
    order = CONSTANT
    family = MONOMIAL
  []
  [effective_creep_strain]
    block = clad
    order = CONSTANT
    family = MONOMIAL
  []
  [gap_cond]
    order = CONSTANT
    family = MONOMIAL
  []
[]

[Functions]
  [power_history] #related to the LHGR at the midplane
    type = PiecewiseLinear
    x = '0 70000 25000000'
    y = '0 38974.7  38974.7'
  []
  [fast_neutron_flux_function]
    type = PiecewiseLinear
    x = '0 70000 25000000'
    y = '0 1.9e+19 1.9e+19'
  []
  [axial_peaking_factors]
    type = PiecewiseBilinear
    x = '0  0.071  0.146  0.221  0.296  0.37  0.443  0.566'
    y = '0 25000000'
    z = '0.889  1.041  1.152  1.173  1.129  0.971  0.782  0.672 0.889  1.041  1.152  1.173  1.129  0.971  0.782  0.672'
    scale_factor = 1
    axis = 1
  []
  [q]
    type = CompositeFunction
    functions = 'power_history axial_peaking_factors'
  []
  [average_power_history]
    type = PiecewiseLinear
    x = '0 70000 25000000'
    y = '0 32000 32000'
  []
  [clad_surface_temp]
    type = PiecewiseBilinear
    x = '0  0.071  0.146  0.221  0.296  0.37  0.443  0.566'
    y = '0 25000000'
    z = '295 295 295 295 295 295 295 295 593.58  606.36  619.13  630.26  640.87  651.76  662.67 673.67'
    scale_factor = 1
    axis = 1
  []
  [pressure_ramp]
    type = PiecewiseLinear
    x = '0 1'
    y = '1 1'
  []
[]
[Physics/SolidMechanics/QuasiStatic]
  [pellets]
    block = pellet
    add_variables = true
    strain = FINITE
    eigenstrain_names = 'fuel_thermal_strain fuel_volumetric_strain'
    generate_output = 'vonmises_stress stress_xx stress_yy stress_zz'
    use_finite_deform_jacobian = true
    extra_vector_tags = 'ref'
  []
  [clad]
    block = clad
    add_variables = true
    strain = FINITE
    eigenstrain_names = 'clad_thermal_eigenstrain'
    generate_output = 'vonmises_stress stress_xx stress_yy stress_zz'
    use_finite_deform_jacobian = true
    extra_vector_tags = 'ref'
  []
[]
[Kernels]
  [gravity]
    type = Gravity
    variable = disp_y
    value = -9.81
    extra_vector_tags = 'ref'
  []
  [heat]
    type = HeatConduction
    variable = temp
    extra_vector_tags = 'ref'
  []
  [heat_ie]
    type = HeatConductionTimeDerivative
    variable = temp
    extra_vector_tags = 'ref'
  []
  [heat_source]
    type = NeutronHeatSource
    variable = temp
    block = pellet
    fission_rate = fission_rate
    extra_vector_tags = 'ref'
  []
[]
[AuxKernels]
  [fission_rate]
    type = FissionRateGeneral
    fission_rate_formulation = MOX
    variable = fission_rate
    block = pellet
    initial_porosity = 0.065
    axial_power_profile = axial_peaking_factors
    rod_ave_lin_pow = power_history
    pellet_diameter = 0.0054
    execute_on = timestep_begin
    porosity = pore
  []
  [burnup]
    type = BurnupAux
    block = pellet
    fission_rate = fission_rate
    variable = burnup
    execute_on = timestep_begin
  []
  [fggen]
    type = MaterialRealAux
    variable = gas_gen_3
    property = gas_concentration_generated_total
    execute_on = timestep_end
  []
  [fggrn]
    type = MaterialRealAux
    variable = gas_grn_3
    property = gas_concentration_intra_total
    execute_on = timestep_end
  []
  [fgbdr]
    type = MaterialRealAux
    variable = gas_bdr_3
    property = gas_concentration_GB_bubble_volume
    execute_on = timestep_end
  []
  [fgrel]
    type = MaterialRealAux
    variable = gas_rel_3
    property = gas_concentration_release_total
    execute_on = timestep_end
  []
  [nbbl2]
    type = MaterialRealAux
    variable = bbl_bdr_2
    property = bubble_GB_surface_density
    execute_on = timestep_end
  []
  [nvcn2]
    type = MaterialRealAux
    variable = vcn_bdr_2
    property = vacancy_concentration_GB_surface
    execute_on = timestep_end
  []
  [atmbbl]
    type = MaterialRealAux
    variable = atm_bbl_bdr
    property = atom_per_bubble_GB
    execute_on = timestep_end
  []
  [vcnbbl]
    type = MaterialRealAux
    variable = vcn_bbl_bdr
    property = vacancy_per_bubble_GB
    execute_on = timestep_end
  []
  [prsbbl]
    type = MaterialRealAux
    variable = prs_bbl_bdr
    property = bubble_GB_pressure
    execute_on = timestep_end
  []
  [prseqbbl]
    type = MaterialRealAux
    variable = prseq_bbl_bdr
    property = bubble_GB_pressure_equilibrium
    execute_on = timestep_end
  []
  [radbbl]
    type = MaterialRealAux
    variable = rad_bbl_bdr
    property = bubble_radius_GB
    execute_on = timestep_end
  []
  [volbbl]
    type = MaterialRealAux
    variable = vol_bbl_bdr
    property = bubble_GB_volume
    execute_on = timestep_end
  []
  [frcvrg]
    type = MaterialRealAux
    variable = GBCoverage
    property = GBCoverage
    execute_on = timestep_end
  []
  [diffc]
    type = MaterialRealAux
    variable = eff_diff_coeff
    property = eff_diff_coeff
    execute_on = timestep_end
  []
  [dvv0bd]
    type = MaterialRealAux
    variable = deltav_v0_bd
    property = deltav_v0_bubble_GB
    execute_on = timestep_end
  []
  [radial_strain]
    type = RankTwoAux
    rank_two_tensor = total_strain
    variable = radial_strain
    index_i = 0
    index_j = 0
    execute_on = timestep_end
  []
  [effective_creep_strain]
    type = MaterialRealAux
    property = effective_creep_strain
    variable = effective_creep_strain
    block = clad
    execute_on = timestep_end
  []
  [conductance]
    type = MaterialRealAux
    property = gap_conductance
    variable = gap_cond
    boundary = 10
    execute_on = 'linear'
  []
[]
[Contact]
  [pellet_clad_mechanical]
    primary = 5
    secondary = 10
    formulation = kinematic
    model = frictionless
    penalty = 1e7
  []
[]
[ThermalContact]
  [thermal_contact]
    type = GasGapHeatTransfer
    variable = temp
    primary = 5
    secondary = 10
    initial_moles = initial_moles
    gas_released = fis_gas_released
    contact_pressure = contact_pressure
    quadrature = true
  []
[]
[BCs]
  [no_x_all]
    type = DirichletBC
    variable = disp_x
    boundary = '12'
    value = 0.0
  []
  [no_y_clad_bottom]
    type = DirichletBC
    variable = disp_y
    boundary = 1
    value = 0.0
  []
  [no_y_fuel_bottom]
    type = DirichletBC
    variable = disp_y
    boundary = 20
    value = 0.0
  []
  [temp_clad_out]
    type = FunctionDirichletBC
    variable = temp
    boundary = '2'
    function = clad_surface_temp
  []
  [Pressure]
    [coolantPressure]
      boundary = '1 2 3'
      factor = 101325
      function = pressure_ramp
    []
  []
  [PlenumPressure]
    [plenumPressure]
      boundary = 9
      initial_pressure = 300000
      startup_time = 0
      R = 8.3143
      output_initial_moles = initial_moles
      temperature = ave_temp_interior
      volume = gas_volume
      material_input = fis_gas_released
      output = plenum_pressure
    []
  []
[]
[Materials]
  [fast_neutron_flux]
    type = FastNeutronFlux
    calculate_fluence = true
    block = clad
    flux_function = fast_neutron_flux_function
  []
  [fuel_thermal]
    type = MAMOXThermal
    block = pellet
    temperature = temp
    Am_content = 0.0
    Np_content = 0.0
    porosity = pore
    output_properties = 'thermal_conductivity'
  []
  [fuel_elasticity_tensor]
    type = MAMOXElasticityTensor
    block = pellet
  []
  [elastic_stress]
    type = ComputeFiniteStrainElasticStress
    block = pellet
  []
  [fuel_thermal_expansion]
    type = MAMOXThermalExpansionEigenstrain
    block = pellet
    temperature = temp
    stress_free_temperature = 295.0
    oxygen_to_metal_ratio = 1.98
    eigenstrain_name = fuel_thermal_strain
  []
  [fuel_volumetric_swelling]
    type = UO2VolumetricSwellingEigenstrain
    gas_swelling_model_type = SIFGRS
    block = pellet
    temperature = temp
    burnup = burnup
    initial_fuel_density = 10836.8
    eigenstrain_name = fuel_volumetric_strain
  []
  [clad_thermal]
    type = SS316Thermal
    block = clad
    temperature = temp
  []
  [clad_density]
    type = StrainAdjustedDensity
    block = clad
    strain_free_density = 8000
  []
  [clad_elasticity_tensor]
    type = SS316ElasticityTensor
    block = clad
    temperature = temp
    elastic_constants_model = legacy_ifr
  []
  [thermal_expansion]
    type = SS316ThermalExpansionEigenstrain
    block = clad
    temperature = temp
    stress_free_temperature = 295.0
    eigenstrain_name = clad_thermal_eigenstrain
  []
  [clad_ss316creep]
    type = SS316CreepUpdate
    block = clad
    temperature = temp
    fast_neutron_flux = fast_neutron_flux
  []
  [clad_stress]
    type = ComputeMultipleInelasticStress
    tangent_operator = elastic
    inelastic_models = 'clad_ss316creep'
    block = clad
  []
  [fission_gas_release]
    type = UO2Sifgrs
    block = pellet
    temperature = temp
    burnup = burnup
    diff_coeff_option = TURNBULL_D1_4D2_4D3
    fission_rate = fission_rate
    grain_radius_const = 8.01e-6 #I'm keeping the grain radius const because the grain growth in MOX is probably different due to high Temp
    bubble_gb_limit = 1.0e+11
  []
  [fuel_density]
    type = StrainAdjustedDensity
    block = pellet
    strain_free_density = ${initial_fuel_density}
  []
[]
[Preconditioning]
  [SMP]
    type = SMP
    full = true
  []
[]
[Executioner]
  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'
  fixed_point_abs_tol = 1e-5
  fixed_point_rel_tol = 1e-5
  fixed_point_max_its = 1
  l_max_its = 70
  l_tol = 8e-3
  nl_max_its = 70
  nl_rel_tol = 1e-5
  nl_abs_tol = 1e-5
  start_time = 0
  n_startup_steps = 1
  end_time = 25000000
  dtmax = 1e6
  dtmin = 0.25
  [TimeStepper]
    type = IterationAdaptiveDT
    dt = 5000
    optimal_iterations = 15
    iteration_window = 2
    linear_iteration_ratio = 100
    growth_factor = 2
    cutback_factor = .5
    force_step_every_function_point = true
    timestep_limiting_function = power_history
  []
[]
[Postprocessors]
  [ave_temp_interior]
    type = SideAverageValue
    boundary = 9
    variable = temp
    execute_on = 'initial linear'
  []
  [average_burnup]
    type = ElementAverageValue
    block = pellet
    variable = burnup
  []
  [clad_inner_vol]
    type = InternalVolume
    boundary = 7
    execute_on = 'initial timestep_end'
  []
  [pellet_volume]
    type = InternalVolume
    boundary = 8
    execute_on = 'initial timestep_end'
  []
  [avg_clad_temp]
    type = SideAverageValue
    boundary = 7
    variable = temp
    execute_on = 'initial timestep_end'
  []
  [fis_gas_produced]
    type = ElementIntegralFisGasGeneratedSifgrs
    block = pellet
    execute_on = 'linear'
  []
  [fis_gas_released]
    type = ElementIntegralFisGasReleasedSifgrs
    block = pellet
    execute_on = 'linear'
  []
  [fis_gas_released_percentage]
    type = FGRPercent
    fission_gas_generated = fis_gas_produced
    fission_gas_released = fis_gas_released
    execute_on = 'linear'
  []
  [fis_gas_grain]
    type = ElementIntegralFisGasGrainSifgrs
    block = pellet
    execute_on = 'linear'
  []
  [fis_gas_boundary]
    type = ElementIntegralFisGasBoundarySifgrs
    block = pellet
    execute_on = 'linear'
  []
  [gas_volume]
    type = InternalVolume
    boundary = 9
    execute_on = 'initial linear'
  []
  [flux_from_fuel]
    type = SideDiffusiveFluxIntegral
    variable = temp
    boundary = 10
    diffusivity = thermal_conductivity
  []
  [rod_total_power]
    type = ElementIntegralPower
    variable = temp
    fission_rate = fission_rate
    block = pellet
  []
  [rod_input_power]
    type = FunctionValuePostprocessor
    function = average_power_history
    scale_factor = 0.6 # rod height
  []
  [average_vonMises_fuel]
    type = ElementAverageValue
    variable = vonmises_stress
    block = pellet
  []
  [average_vonMises_clad]
    type = ElementAverageValue
    variable = vonmises_stress
    block = clad
  []
  [average_strain_rr_fuel]
    type = ElementAverageValue
    variable = radial_strain
    block = pellet
  []
  [average_strain_rr_clad]
    type = ElementAverageValue
    variable = radial_strain
    block = clad
  []
  [average_creep_strain_clad]
    type = ElementAverageValue
    variable = effective_creep_strain
    block = clad
  []
  [ave_pore]
    type = ElementAverageValue
    variable = pore
  []
  [max_pore]
    type = NodalExtremeValue
    value_type = max
    variable = pore
  []
  [min_pore]
    type = NodalExtremeValue
    value_type = min
    variable = pore
  []
[]

[PerformanceMetricOutputs]
[]

[Outputs]
  perf_graph = true
  exodus = true
  color = true
  csv = true
  [console]
    type = Console
    max_rows = 25
  []
  [chkfile]
    type = CSV
    execute_on = FINAL
    show = 'ave_temp_interior fis_gas_released_percentage max_pore'
  []
[]
[MultiApps]
  [sub]
    type = TransientMultiApp
    app_type = BisonApp
    execute_on = timestep_end
    catch_up = true
    max_catch_up_steps = 10
    positions_file = positions.txt
    input_files = MK-I_75MW_sub_new_bubble_gb_lim.i
  []
[]
[Transfers]
  [temp_to_sub]
    type = MultiAppProjectionTransfer
    to_multi_app = sub
    source_variable = temp
    variable = temp
  []
  [pore_from_sub]
    type = MultiAppGeometricInterpolationTransfer
    from_multi_app = sub
    source_variable = pore
    variable = pore
  []
[]
[Debug]
  show_var_residual_norms = true
  show_var_residual = 'temp disp_x disp_y'
[]

(test/tests/solid_mechanics/SS316_creep/creep_SS316_rz_therm_only.i)

#
# This is an axisymmetric RZ model designed to confirm the implementation of
# the thermal creep equations in the paper by Gorash, Altenbach and Lvov. The data in
# Figure 4 of the paper plots the steady state strain rate as a function of
# applied stress. This model uses a pressure BC at one end while holding the
# other end fixed in the y-direction. The strain rate in the y-direction is
# extracted as a postprocessor and compared to the values in Figure 4. The
# irradiation creep term is eliminated by setting the fast neutron flux equal
# to zero in this simulation.
#
# The cli_args option is used in the tests file to create additional tests
# with a variation in the applied stress and temperature. A gnuplot input file
# is also included in the directory to generate a plot of the steady state
# strain rate results.
#
# "Modelling of high-temperature inelastic behaviour of the austenitic steel
# AISI type 316 using a continuum damage mechanics approach", Y. Gorash, H.
# Altenbach and G. Lvov, J. Strain Analysis, 47(4) 229-243, 2012.
#

[GlobalParams]
  displacements = 'disp_x disp_y'
  order = FIRST
  family = LAGRANGE
  temperature = temp
[]

[Mesh]
  coord_type = RZ
  [mesh]
    type = GeneratedMeshGenerator
    dim = 2
    ymax = 0.1
    xmax = 0.01
  []
[]

[Variables]
  [temp]
    initial_condition = 1023.0
  []
[]

[AuxVariables]
  [stress_yy]
    order = CONSTANT
    family = MONOMIAL
  []
  [elastic_strain_yy]
    order = CONSTANT
    family = MONOMIAL
  []
  [creep_strain_yy]
    order = CONSTANT
    family = MONOMIAL
  []
[]

[Physics]
  [SolidMechanics]
    [QuasiStatic]
      [all]
        strain = FINITE
        add_variables = true
        generate_output = 'stress_yy elastic_strain_yy creep_strain_yy'
      []
    []
  []
[]

[Functions]
  [top_pull]
    type = PiecewiseLinear
    x = '0 1e5'
    y = '0 1'
  []
[]

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

[AuxKernels]
  [stress_yy]
    type = RankTwoAux
    rank_two_tensor = stress
    variable = stress_yy
    index_i = 1
    index_j = 1
  []
  [elastic_strain_yy]
    type = RankTwoAux
    rank_two_tensor = elastic_strain
    variable = elastic_strain_yy
    index_i = 1
    index_j = 1
  []
  [creep_strain_yy]
    type = RankTwoAux
    rank_two_tensor = creep_strain
    variable = creep_strain_yy
    index_i = 1
    index_j = 1
    execute_on = timestep_end
  []
[]

[BCs]
  [u_top_pull]
    type = Pressure
    variable = 'disp_y'
    boundary = top
    factor = -30.0e6
    function = top_pull
  []
  [u_bottom_fix]
    type = DirichletBC
    variable = disp_y
    boundary = bottom
    value = 0.0
  []
  [u_left_fix]
    type = DirichletBC
    variable = disp_x
    boundary = left
    value = 0.0
  []
  [temp_fix]
    type = DirichletBC
    variable = temp
    boundary = 'top bottom left right'
    value = 1023.0
  []
[]

[Materials]
  [fast_neutron_flux]
    type = FastNeutronFlux
    calculate_fluence = true
    factor = 0
  []
  [elasticity_tensor]
    type = ComputeIsotropicElasticityTensor
    youngs_modulus = 1.90e11
    poissons_ratio = 0.265
  []
  [creep]
    type = SS316CreepUpdate
    fast_neutron_flux = fast_neutron_flux
    outputs = 'all'
  []
  [radial_return_stress]
    type = ComputeMultipleInelasticStress
    tangent_operator = elastic
    inelastic_models = 'creep'
  []
  [thermal]
    type = HeatConductionMaterial
    thermal_conductivity = 100.0
  []
  [density]
    type = StrainAdjustedDensity
    strain_free_density = 1.0
  []
[]

[Executioner]
  type = Transient
  solve_type = 'PJFNK'

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

  line_search = 'none'

  l_max_its = 50
  nl_max_its = 20
  l_tol = 1e-5
  nl_rel_tol = 1e-5
  nl_abs_tol = 1e-8
  start_time = 0.0
  num_steps = 5000
  end_time = 4e5
  dtmin = 10.0
  dtmax = 1e4

  [TimeStepper]
    type = IterationAdaptiveDT
    dt = 1e3
    optimal_iterations = 20
    iteration_window = 4
    time_t = '0   1.0e5'
    time_dt = '1e3  1e4'
    growth_factor = 2.0
    cutback_factor = 0.5
  []
[]

[Postprocessors]
  [stress_yy]
    type = ElementAverageValue
    variable = stress_yy
  []
  [el_strain_yy]
    type = ElementAverageValue
    variable = elastic_strain_yy
  []
  [crp_strain_yy]
    type = ElementAverageValue
    variable = creep_strain_yy
  []
  [delta_crp_strain_yy]
    type = ChangeOverTimePostprocessor
    postprocessor = crp_strain_yy
  []
  [_dt]
    type = TimestepSize
  []
[]

[Outputs]
  exodus = true
  csv = true
[]

(assessment/MOX/JOYO/MK-I/analysis/MK-I_50MW_master_old_bubble_gb_lim.i)


initial_fuel_density = 10836.8

[GlobalParams]
  density = ${initial_fuel_density}
  initial_porosity = 0.065
  order = SECOND
  family = LAGRANGE
  energy_per_fission = 3.2e-11 # J/fission
  volumetric_locking_correction = false
  displacements = 'disp_x disp_y'
[]
[Problem]
  type = ReferenceResidualProblem
  extra_tag_vectors = 'ref'
  reference_vector = 'ref'
[]
[Mesh]
  coord_type = RZ
  [smeared_pellet_mesh]
    type = FuelPinMeshGenerator
    pellet_quantity = 1
    pellet_height = 0.6
    pellet_outer_radius = 0.0027
    pellet_mesh_density = customize
    clad_mesh_density = customize
    clad_gap_width = 0.000100
    clad_thickness = 0.00035
    clad_bot_gap_height = 1.0e-3
    bottom_clad_height = 2.24e-3
    top_clad_height = 2.24e-3
    clad_top_gap_height = 0.599
    elem_type = QUAD8
    nx_c = 4
    ny_c = 200
    nx_p = 20
    ny_p = 200
    ny_cu = 3
    ny_cl = 3
  []
  patch_size = 50
  patch_update_strategy = iteration
  partitioner = centroid
  centroid_partitioner_direction = y
[]
[UserObjects]
  [pin_geometry]
    type = FuelPinGeometry
  []
[]
[Variables]
  [temp]
    initial_condition = 295.0
  []
[]
[AuxVariables]
  [pore]
  []
  [fission_rate]
    block = pellet
  []
  [burnup]
    block = pellet
  []
  [gas_gen_3]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [gas_grn_3]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [gas_bdr_3]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [gas_rel_3]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [bbl_bdr_2]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [vcn_bdr_2]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [atm_bbl_bdr]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [vcn_bbl_bdr]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [prs_bbl_bdr]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [prseq_bbl_bdr]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [rad_bbl_bdr]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [vol_bbl_bdr]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [GBCoverage]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [eff_diff_coeff]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [deltav_v0_bd]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [radial_strain]
    order = CONSTANT
    family = MONOMIAL
  []
  [effective_creep_strain]
    block = clad
    order = CONSTANT
    family = MONOMIAL
  []
  [gap_cond]
    order = CONSTANT
    family = MONOMIAL
  []
[]

[Functions]
  [power_history] #related to the LHGR at the midplane
    type = PiecewiseLinear
    x = '0 70000 17153028'
    y = '0  25577  25577'
  []
  [fast_neutron_flux_function]
    type = PiecewiseLinear
    x = '0 70000 17153028'
    y = '0 1.2e+19 1.2e+19'
  []
  [axial_peaking_factors]
    type = PiecewiseBilinear
    x = '0  0.071  0.146  0.221  0.296  0.37  0.443  0.566'
    y = '0 17153028'
    z = '0.889  1.041  1.152  1.173  1.129  0.971  0.782  0.672 0.889  1.041  1.152  1.173  1.129  0.971  0.782  0.672'
    scale_factor = 1
    axis = 1
  []
  [q]
    type = CompositeFunction
    functions = 'power_history axial_peaking_factors'
  []
  [average_power_history]
    type = PiecewiseLinear
    x = '0 70000 17153028'
    y = '0 21000 21000'
  []
  [clad_surface_temp]
    type = PiecewiseBilinear
    x = '0  0.075  0.15  0.225  0.3  0.375  0.45  0.525  0.6'
    y = '0 17153028'
    z = '295 295 295 295 295 295 295 295 295 499.9  509.1  517.8  525.42  532.71  540.29  547.7  552.3  554.81'
    scale_factor = 1
    axis = 1
  []
  [pressure_ramp]
    type = PiecewiseLinear
    x = '0 1'
    y = '1 1'
  []
[]
[Physics/SolidMechanics/QuasiStatic]
  [pellets]
    block = pellet
    add_variables = true
    strain = FINITE
    eigenstrain_names = 'fuel_thermal_strain fuel_volumetric_strain'
    generate_output = 'vonmises_stress stress_xx stress_yy stress_zz'
    use_finite_deform_jacobian = true
    extra_vector_tags = 'ref'
  []
  [clad]
    block = clad
    add_variables = true
    strain = FINITE
    eigenstrain_names = 'clad_thermal_eigenstrain'
    generate_output = 'vonmises_stress stress_xx stress_yy stress_zz'
    use_finite_deform_jacobian = true
    extra_vector_tags = 'ref'
  []
[]
[Kernels]
  [gravity]
    type = Gravity
    variable = disp_y
    value = -9.81
    extra_vector_tags = 'ref'
  []
  [heat]
    type = HeatConduction
    variable = temp
    extra_vector_tags = 'ref'
  []
  [heat_ie]
    type = HeatConductionTimeDerivative
    variable = temp
    extra_vector_tags = 'ref'
  []
  [heat_source]
    type = NeutronHeatSource
    variable = temp
    block = pellet
    fission_rate = fission_rate
    extra_vector_tags = 'ref'
  []
[]
[AuxKernels]
  [fission_rate]
    type = FissionRateGeneral
    fission_rate_formulation = MOX
    variable = fission_rate
    block = pellet
    initial_porosity = 0.065
    axial_power_profile = axial_peaking_factors
    rod_ave_lin_pow = power_history
    pellet_diameter = 0.0054
    execute_on = timestep_begin
    porosity = pore
  []
  [burnup]
    type = BurnupAux
    block = pellet
    fission_rate = fission_rate
    variable = burnup
    execute_on = timestep_begin
  []
  [fggen]
    type = MaterialRealAux
    variable = gas_gen_3
    property = gas_concentration_generated_total
    execute_on = timestep_end
  []
  [fggrn]
    type = MaterialRealAux
    variable = gas_grn_3
    property = gas_concentration_intra_total
    execute_on = timestep_end
  []
  [fgbdr]
    type = MaterialRealAux
    variable = gas_bdr_3
    property = gas_concentration_GB_bubble_volume
    execute_on = timestep_end
  []
  [fgrel]
    type = MaterialRealAux
    variable = gas_rel_3
    property = gas_concentration_release_total
    execute_on = timestep_end
  []
  [nbbl2]
    type = MaterialRealAux
    variable = bbl_bdr_2
    property = bubble_GB_surface_density
    execute_on = timestep_end
  []
  [nvcn2]
    type = MaterialRealAux
    variable = vcn_bdr_2
    property = vacancy_concentration_GB_surface
    execute_on = timestep_end
  []
  [atmbbl]
    type = MaterialRealAux
    variable = atm_bbl_bdr
    property = atom_per_bubble_GB
    execute_on = timestep_end
  []
  [vcnbbl]
    type = MaterialRealAux
    variable = vcn_bbl_bdr
    property = vacancy_per_bubble_GB
    execute_on = timestep_end
  []
  [prsbbl]
    type = MaterialRealAux
    variable = prs_bbl_bdr
    property = bubble_GB_pressure
    execute_on = timestep_end
  []
  [prseqbbl]
    type = MaterialRealAux
    variable = prseq_bbl_bdr
    property = bubble_GB_pressure_equilibrium
    execute_on = timestep_end
  []
  [radbbl]
    type = MaterialRealAux
    variable = rad_bbl_bdr
    property = bubble_radius_GB
    execute_on = timestep_end
  []
  [volbbl]
    type = MaterialRealAux
    variable = vol_bbl_bdr
    property = bubble_GB_volume
    execute_on = timestep_end
  []
  [frcvrg]
    type = MaterialRealAux
    variable = GBCoverage
    property = GBCoverage
    execute_on = timestep_end
  []
  [diffc]
    type = MaterialRealAux
    variable = eff_diff_coeff
    property = eff_diff_coeff
    execute_on = timestep_end
  []
  [dvv0bd]
    type = MaterialRealAux
    variable = deltav_v0_bd
    property = deltav_v0_bubble_GB
    execute_on = timestep_end
  []
  [radial_strain]
    type = RankTwoAux
    rank_two_tensor = total_strain
    variable = radial_strain
    index_i = 0
    index_j = 0
    execute_on = timestep_end
  []
  [effective_creep_strain]
    type = MaterialRealAux
    property = effective_creep_strain
    variable = effective_creep_strain
    block = clad
    execute_on = timestep_end
  []
  [conductance]
    type = MaterialRealAux
    property = gap_conductance
    variable = gap_cond
    boundary = 10
    execute_on = 'linear'
  []
[]
[Contact]
  [pellet_clad_mechanical]
    primary = 5
    secondary = 10
    formulation = kinematic
    model = frictionless
    penalty = 1e7
  []
[]
[ThermalContact]
  [thermal_contact]
    type = GasGapHeatTransfer
    variable = temp
    primary = 5
    secondary = 10
    initial_moles = initial_moles
    gas_released = fis_gas_released
    contact_pressure = contact_pressure
    quadrature = true
  []
[]

[BCs]
  [no_x_all]
    type = DirichletBC
    variable = disp_x
    boundary = '12'
    value = 0.0
  []
  [no_y_clad_bottom]
    type = DirichletBC
    variable = disp_y
    boundary = 1
    value = 0.0
  []
  [no_y_fuel_bottom]
    type = DirichletBC
    variable = disp_y
    boundary = 20
    value = 0.0
  []
  [temp_clad_out]
    type = FunctionDirichletBC
    variable = temp
    boundary = '2'
    function = clad_surface_temp
  []
  [Pressure]
    [coolantPressure]
      boundary = '1 2 3'
      factor = 101325
      function = pressure_ramp
    []
  []
  [PlenumPressure]
    [plenumPressure]
      boundary = 9
      initial_pressure = 300000
      startup_time = 0
      R = 8.3143
      output_initial_moles = initial_moles
      temperature = ave_temp_interior
      volume = gas_volume
      material_input = fis_gas_released
      output = plenum_pressure
    []
  []
[]
[Materials]
  [fast_neutron_flux]
    type = FastNeutronFlux
    calculate_fluence = true
    block = clad
    flux_function = fast_neutron_flux_function
  []
  [fuel_thermal]
    type = MAMOXThermal
    block = pellet
    temperature = temp
    Am_content = 0.0
    Np_content = 0.0
    porosity = pore
    output_properties = 'thermal_conductivity'
  []
  [fuel_elasticity_tensor]
    type = MAMOXElasticityTensor
    block = pellet
  []
  [elastic_stress]
    type = ComputeFiniteStrainElasticStress
    block = pellet
  []
  [fuel_thermal_expansion]
    type = MAMOXThermalExpansionEigenstrain
    block = pellet
    temperature = temp
    stress_free_temperature = 295.0
    oxygen_to_metal_ratio = 1.98
    eigenstrain_name = fuel_thermal_strain
  []
  [fuel_volumetric_swelling]
    type = UO2VolumetricSwellingEigenstrain
    gas_swelling_model_type = SIFGRS
    block = pellet
    temperature = temp
    burnup = burnup
    initial_fuel_density = 10836.8
    eigenstrain_name = fuel_volumetric_strain
  []
  [clad_thermal]
    type = SS316Thermal
    block = clad
    temperature = temp
  []
  [clad_density]
    type = StrainAdjustedDensity
    block = clad
    strain_free_density = 8000
  []
  [clad_elasticity_tensor]
    type = SS316ElasticityTensor
    block = clad
    temperature = temp
    elastic_constants_model = legacy_ifr
  []
  [thermal_expansion]
    type = SS316ThermalExpansionEigenstrain
    block = clad
    temperature = temp
    stress_free_temperature = 295.0
    eigenstrain_name = clad_thermal_eigenstrain
  []
  [clad_ss316creep]
    type = SS316CreepUpdate
    block = clad
    temperature = temp
    fast_neutron_flux = fast_neutron_flux
  []
  [clad_stress]
    type = ComputeMultipleInelasticStress
    tangent_operator = elastic
    inelastic_models = 'clad_ss316creep'
    block = clad
  []
  [fission_gas_release]
    type = UO2Sifgrs
    block = pellet
    temperature = temp
    burnup = burnup
    diff_coeff_option = TURNBULL_D1_4D2_4D3
    fission_rate = fission_rate
    grain_radius_const = 8.01e-6 #I'm keeping the grain radius const because the grain growth in MOX is probably different due to high Temp
  []
  [fuel_density]
    type = StrainAdjustedDensity
    block = pellet
    strain_free_density = ${initial_fuel_density}
  []
[]
[Preconditioning]
  [SMP]
    type = SMP
    full = true
  []
[]
[Executioner]
  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'
  fixed_point_abs_tol = 1e-5
  fixed_point_rel_tol = 1e-5
  fixed_point_max_its = 1
  l_max_its = 70
  l_tol = 8e-3
  nl_max_its = 70
  nl_rel_tol = 1e-5
  nl_abs_tol = 1e-5
  start_time = 0
  n_startup_steps = 1
  end_time = 17153028
  dtmax = 1e6
  dtmin = 0.25
  [TimeStepper]
    type = IterationAdaptiveDT
    dt = 5000
    optimal_iterations = 15
    iteration_window = 2
    linear_iteration_ratio = 100
    growth_factor = 2
    cutback_factor = .5
    force_step_every_function_point = true
    timestep_limiting_function = power_history
  []
[]
[Postprocessors]
  [ave_temp_interior]
    type = SideAverageValue
    boundary = 9
    variable = temp
    execute_on = 'initial linear'
  []
  [average_burnup]
    type = ElementAverageValue
    block = pellet
    variable = burnup
  []
  [clad_inner_vol]
    type = InternalVolume
    boundary = 7
    execute_on = 'initial timestep_end'
  []
  [pellet_volume]
    type = InternalVolume
    boundary = 8
    execute_on = 'initial timestep_end'
  []
  [avg_clad_temp]
    type = SideAverageValue
    boundary = 7
    variable = temp
    execute_on = 'initial timestep_end'
  []
  [fis_gas_produced]
    type = ElementIntegralFisGasGeneratedSifgrs
    block = pellet
    execute_on = 'linear'
  []
  [fis_gas_released]
    type = ElementIntegralFisGasReleasedSifgrs
    block = pellet
    execute_on = 'linear'
  []
  [fis_gas_released_percentage]
    type = FGRPercent
    fission_gas_generated = fis_gas_produced
    fission_gas_released = fis_gas_released
    execute_on = 'linear'
  []
  [fis_gas_grain]
    type = ElementIntegralFisGasGrainSifgrs
    block = pellet
    execute_on = 'linear'
  []
  [fis_gas_boundary]
    type = ElementIntegralFisGasBoundarySifgrs
    block = pellet
    execute_on = 'linear'
  []
  [gas_volume]
    type = InternalVolume
    boundary = 9
    execute_on = 'initial linear'
  []
  [flux_from_fuel]
    type = SideDiffusiveFluxIntegral
    variable = temp
    boundary = 10
    diffusivity = thermal_conductivity
  []
  [rod_total_power]
    type = ElementIntegralPower
    variable = temp
    fission_rate = fission_rate
    block = pellet
  []
  [rod_input_power]
    type = FunctionValuePostprocessor
    function = average_power_history
    scale_factor = 0.6 # rod height
  []
  [average_vonMises_fuel]
    type = ElementAverageValue
    variable = vonmises_stress
    block = pellet
  []
  [average_vonMises_clad]
    type = ElementAverageValue
    variable = vonmises_stress
    block = clad
  []
  [average_strain_rr_fuel]
    type = ElementAverageValue
    variable = radial_strain
    block = pellet
  []
  [average_strain_rr_clad]
    type = ElementAverageValue
    variable = radial_strain
    block = clad
  []
  [average_creep_strain_clad]
    type = ElementAverageValue
    variable = effective_creep_strain
    block = clad
  []
  [ave_pore]
    type = ElementAverageValue
    variable = pore
  []
  [max_pore]
    type = NodalExtremeValue
    value_type = max
    variable = pore
  []
  [min_pore]
    type = NodalExtremeValue
    value_type = min
    variable = pore
  []
[]

[PerformanceMetricOutputs]
[]

[Outputs]
  perf_graph = true
  exodus = true
  color = true
  csv = true
  [console]
    type = Console
    max_rows = 25
  []
  [chkfile]
    type = CSV
    execute_on = FINAL
    show = 'ave_temp_interior fis_gas_released_percentage max_pore'
  []
[]
[MultiApps]
  [sub]
    type = TransientMultiApp
    app_type = BisonApp
    execute_on = timestep_end
    catch_up = true
    max_catch_up_steps = 10
    positions_file = positions.txt
    input_files = MK-I_50MW_sub_old_bubble_gb_lim.i
  []
[]
[Transfers]
  [temp_to_sub]
    type = MultiAppProjectionTransfer
    to_multi_app = sub
    source_variable = temp
    variable = temp
  []
  [pore_from_sub]
    type = MultiAppGeometricInterpolationTransfer
    from_multi_app = sub
    source_variable = pore
    variable = pore
  []
[]
[Debug]
  show_var_residual_norms = true
  show_var_residual = 'temp disp_x disp_y'
[]

(assessment/MOX/JOYO/MK-II/analysis/MK-II_master_old_bubble_gb_lim.i)


initial_fuel_density = 10920.4

[GlobalParams]
  density = ${initial_fuel_density}
  initial_porosity = 0.07
  order = SECOND
  family = LAGRANGE
  energy_per_fission = 3.2e-11 # J/fission
  volumetric_locking_correction = false
  displacements = 'disp_x disp_y'
[]
[Problem]
  type = ReferenceResidualProblem
  extra_tag_vectors = 'ref'
  reference_vector = 'ref'
[]
[Mesh]
  coord_type = RZ
  [smeared_pellet_mesh]
    type = FuelPinMeshGenerator
    pellet_quantity = 1
    pellet_height = 0.55
    pellet_outer_radius = 0.002315
    pellet_mesh_density = customize
    clad_mesh_density = customize
    clad_gap_width = 0.000085
    clad_thickness = 0.00035
    clad_bot_gap_height = 1.0e-3
    bottom_clad_height = 2.24e-3
    top_clad_height = 2.24e-3
    clad_top_gap_height = 0.549
    elem_type = QUAD8
    nx_c = 4
    ny_c = 100
    nx_p = 10
    ny_p = 100
    ny_cu = 3
    ny_cl = 3
  []
  patch_size = 50
  patch_update_strategy = iteration
  partitioner = centroid
  centroid_partitioner_direction = y
[]
[UserObjects]
  [pin_geometry]
    type = FuelPinGeometry
  []
[]
[Variables]
  [temp]
    initial_condition = 295.0
  []
[]
[AuxVariables]
  [pore]
  []
  [fission_rate]
    block = pellet
  []
  [burnup]
    block = pellet
  []
  [gas_gen_3]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [gas_grn_3]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [gas_bdr_3]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [gas_rel_3]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [bbl_bdr_2]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [vcn_bdr_2]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [atm_bbl_bdr]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [vcn_bbl_bdr]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [prs_bbl_bdr]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [prseq_bbl_bdr]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [rad_bbl_bdr]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [vol_bbl_bdr]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [GBCoverage]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [eff_diff_coeff]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [deltav_v0_bd]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [radial_strain]
    order = CONSTANT
    family = MONOMIAL
  []
  [effective_creep_strain]
    block = clad
    order = CONSTANT
    family = MONOMIAL
  []
  [gap_cond]
    order = CONSTANT
    family = MONOMIAL
  []
[]
[Functions]
  [power_history] #related to the LHGR at the midplane
    type = PiecewiseLinear
    x = '0  70000  12697021'
    y = '0  48827.8  48827.8'
  []
  [fast_neutron_flux_function]
    type = PiecewiseLinear
    x = '0  70000  12697021'
    y = '0 2.6e+19 2.6e+19'
  []
  [axial_peaking_factors]
    type = PiecewiseBilinear
    x = '0  0.065  0.134  0.202  0.271  0.339  0.406  0.519'
    y = '0 12697021'
    z = '0.889  1.041  1.152  1.173  1.129  0.971  0.782  0.672 0.889  1.041  1.152  1.173  1.129  0.971  0.782  0.672'
    scale_factor = 1
    axis = 1
  []
  [q]
    type = CompositeFunction
    functions = 'power_history axial_peaking_factors'
  []
  [average_power_history]
    type = PiecewiseLinear
    x = '0  70000  12697021'
    y = '0  40000  40000'
  []
  [clad_surface_temp]
    type = PiecewiseBilinear
    x = '0  0.065  0.134  0.202  0.271  0.339  0.406  0.519'
    y = '0 12697021'
    z = '295  295  295  295  295  295  295  295 416.36  422.49  428.63  434.27  439.36  444.71  450.07  455.48'
    scale_factor = 1
    axis = 1
  []
  [pressure_ramp]
    type = PiecewiseLinear
    x = '0 1'
    y = '1 1'
  []
[]
[Physics/SolidMechanics/QuasiStatic]
  [pellets]
    block = pellet
    add_variables = true
    strain = FINITE
    eigenstrain_names = 'fuel_thermal_strain fuel_volumetric_strain'
    generate_output = 'vonmises_stress stress_xx stress_yy stress_zz'
    use_finite_deform_jacobian = true
    extra_vector_tags = 'ref'
  []
  [clad]
    block = clad
    add_variables = true
    strain = FINITE
    eigenstrain_names = 'clad_thermal_eigenstrain'
    generate_output = 'vonmises_stress stress_xx stress_yy stress_zz'
    use_finite_deform_jacobian = true
    extra_vector_tags = 'ref'
  []
[]
[Kernels]
  [gravity]
    type = Gravity
    variable = disp_y
    value = -9.81
    extra_vector_tags = 'ref'
  []
  [heat]
    type = HeatConduction
    variable = temp
    extra_vector_tags = 'ref'
  []
  [heat_ie]
    type = HeatConductionTimeDerivative
    variable = temp
    extra_vector_tags = 'ref'
  []
  [heat_source]
    type = NeutronHeatSource
    variable = temp
    block = pellet
    fission_rate = fission_rate
    extra_vector_tags = 'ref'
  []
[]
[AuxKernels]
  [fission_rate]
    type = FissionRateGeneral
    fission_rate_formulation = MOX
    variable = fission_rate
    block = pellet
    initial_porosity = 0.07
    axial_power_profile = axial_peaking_factors
    rod_ave_lin_pow = power_history
    pellet_diameter = 0.00463
    execute_on = timestep_begin
    porosity = pore
  []
  [burnup]
    type = BurnupAux
    block = pellet
    fission_rate = fission_rate
    variable = burnup
    execute_on = timestep_begin
  []
  [fggen]
    type = MaterialRealAux
    variable = gas_gen_3
    property = gas_concentration_generated_total
    execute_on = timestep_end
  []
  [fggrn]
    type = MaterialRealAux
    variable = gas_grn_3
    property = gas_concentration_intra_total
    execute_on = timestep_end
  []
  [fgbdr]
    type = MaterialRealAux
    variable = gas_bdr_3
    property = gas_concentration_GB_bubble_volume
    execute_on = timestep_end
  []
  [fgrel]
    type = MaterialRealAux
    variable = gas_rel_3
    property = gas_concentration_release_total
    execute_on = timestep_end
  []
  [nbbl2]
    type = MaterialRealAux
    variable = bbl_bdr_2
    property = bubble_GB_surface_density
    execute_on = timestep_end
  []
  [nvcn2]
    type = MaterialRealAux
    variable = vcn_bdr_2
    property = vacancy_concentration_GB_surface
    execute_on = timestep_end
  []
  [atmbbl]
    type = MaterialRealAux
    variable = atm_bbl_bdr
    property = atom_per_bubble_GB
    execute_on = timestep_end
  []
  [vcnbbl]
    type = MaterialRealAux
    variable = vcn_bbl_bdr
    property = vacancy_per_bubble_GB
    execute_on = timestep_end
  []
  [prsbbl]
    type = MaterialRealAux
    variable = prs_bbl_bdr
    property = bubble_GB_pressure
    execute_on = timestep_end
  []
  [prseqbbl]
    type = MaterialRealAux
    variable = prseq_bbl_bdr
    property = bubble_GB_pressure_equilibrium
    execute_on = timestep_end
  []
  [radbbl]
    type = MaterialRealAux
    variable = rad_bbl_bdr
    property = bubble_radius_GB
    execute_on = timestep_end
  []
  [volbbl]
    type = MaterialRealAux
    variable = vol_bbl_bdr
    property = bubble_GB_volume
    execute_on = timestep_end
  []
  [frcvrg]
    type = MaterialRealAux
    variable = GBCoverage
    property = GBCoverage
    execute_on = timestep_end
  []
  [diffc]
    type = MaterialRealAux
    variable = eff_diff_coeff
    property = eff_diff_coeff
    execute_on = timestep_end
  []
  [dvv0bd]
    type = MaterialRealAux
    variable = deltav_v0_bd
    property = deltav_v0_bubble_GB
    execute_on = timestep_end
  []
  [radial_strain]
    type = RankTwoAux
    rank_two_tensor = total_strain
    variable = radial_strain
    index_i = 0
    index_j = 0
    execute_on = timestep_end
  []
  [effective_creep_strain]
    type = MaterialRealAux
    property = effective_creep_strain
    variable = effective_creep_strain
    block = clad
    execute_on = timestep_end
  []
  [conductance]
    type = MaterialRealAux
    property = gap_conductance
    variable = gap_cond
    boundary = 10
    execute_on = 'linear'
  []
[]
[Contact]
  [pellet_clad_mechanical]
    primary = 5
    secondary = 10
    formulation = kinematic
    model = frictionless
    penalty = 1e7
  []
[]
[ThermalContact]
  [thermal_contact]
    type = GasGapHeatTransfer
    variable = temp
    primary = 5
    secondary = 10
    initial_moles = initial_moles
    gas_released = fis_gas_released
    contact_pressure = contact_pressure
    quadrature = true
  []
[]
[BCs]
  [no_x_all]
    type = DirichletBC
    variable = disp_x
    boundary = '12'
    value = 0.0
  []
  [no_y_clad_bottom]
    type = DirichletBC
    variable = disp_y
    boundary = 1
    value = 0.0
  []
  [no_y_fuel_bottom]
    type = DirichletBC
    variable = disp_y
    boundary = 20
    value = 0.0
  []
  [temp_clad_out]
    type = FunctionDirichletBC
    variable = temp
    boundary = '2'
    function = clad_surface_temp
  []
  [Pressure]
    [coolantPressure]
      boundary = '1 2 3'
      factor = 101325
      function = pressure_ramp
    []
  []
  [PlenumPressure]
    [plenumPressure]
      boundary = 9
      initial_pressure = 300000
      startup_time = 0
      R = 8.3143
      output_initial_moles = initial_moles
      temperature = ave_temp_interior
      volume = gas_volume
      material_input = fis_gas_released
      output = plenum_pressure
    []
  []
[]
[Materials]
  [fast_neutron_flux]
    type = FastNeutronFlux
    calculate_fluence = true
    block = clad
    flux_function = fast_neutron_flux_function
  []
  [fuel_thermal]
    type = MAMOXThermal
    block = pellet
    temperature = temp
    Am_content = 0.0
    Np_content = 0.0
    porosity = pore
    output_properties = 'thermal_conductivity'
  []
  [fuel_elasticity_tensor]
    type = MAMOXElasticityTensor
    block = pellet
  []
  [elastic_stress]
    type = ComputeFiniteStrainElasticStress
    block = pellet
  []
  [fuel_thermal_expansion]
    type = MAMOXThermalExpansionEigenstrain
    block = pellet
    temperature = temp
    stress_free_temperature = 295.0
    oxygen_to_metal_ratio = 1.98
    eigenstrain_name = fuel_thermal_strain
  []
  [fuel_volumetric_swelling]
    type = UO2VolumetricSwellingEigenstrain
    gas_swelling_model_type = SIFGRS
    block = pellet
    temperature = temp
    burnup = burnup
    initial_fuel_density = 10920.4
    eigenstrain_name = fuel_volumetric_strain
  []
  [clad_thermal]
    type = SS316Thermal
    block = clad
    temperature = temp
  []
  [clad_density]
    type = StrainAdjustedDensity
    block = clad
    strain_free_density = 8000
  []
  [clad_elasticity_tensor]
    type = SS316ElasticityTensor
    block = clad
    temperature = temp
    elastic_constants_model = legacy_ifr
  []
  [thermal_expansion]
    type = SS316ThermalExpansionEigenstrain
    block = clad
    temperature = temp
    stress_free_temperature = 295.0
    eigenstrain_name = clad_thermal_eigenstrain
  []
  [clad_ss316creep]
    type = SS316CreepUpdate
    block = clad
    temperature = temp
    fast_neutron_flux = fast_neutron_flux
  []
  [clad_stress]
    type = ComputeMultipleInelasticStress
    tangent_operator = elastic
    inelastic_models = 'clad_ss316creep'
    block = clad
  []
  [fission_gas_release]
    type = UO2Sifgrs
    block = pellet
    temperature = temp
    burnup = burnup
    diff_coeff_option = TURNBULL_D1_4D2_4D3
    fission_rate = fission_rate
    grain_radius_const = 8.01e-6
  []
  [fuel_density]
    type = StrainAdjustedDensity
    block = pellet
    strain_free_density = ${initial_fuel_density}
  []
[]
[Preconditioning]
  [SMP]
    type = SMP
    full = true
  []
[]
[Executioner]
  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'
  fixed_point_abs_tol = 1e-2
  fixed_point_rel_tol = 1e-2
  fixed_point_max_its = 1
  l_max_its = 70
  l_tol = 8e-3
  nl_max_its = 70
  nl_rel_tol = 1e-2
  nl_abs_tol = 1e-2
  start_time = 0
  n_startup_steps = 1
  end_time = 12697021
  dtmax = 2e5
  dtmin = 0.25
  [TimeStepper]
    type = IterationAdaptiveDT
    dt = 5000
    optimal_iterations = 15
    iteration_window = 2
    linear_iteration_ratio = 100
    growth_factor = 2
    cutback_factor = .5
    force_step_every_function_point = true
    timestep_limiting_function = power_history
  []
[]

[Postprocessors]
  [ave_temp_interior]
    type = SideAverageValue
    boundary = 9
    variable = temp
    execute_on = 'initial linear'
  []
  [average_burnup]
    type = ElementAverageValue
    block = pellet
    variable = burnup
  []
  [clad_inner_vol]
    type = InternalVolume
    boundary = 7
    execute_on = 'initial timestep_end'
  []
  [pellet_volume]
    type = InternalVolume
    boundary = 8
    execute_on = 'initial timestep_end'
  []
  [avg_clad_temp]
    type = SideAverageValue
    boundary = 7
    variable = temp
    execute_on = 'initial timestep_end'
  []
  [fis_gas_produced]
    type = ElementIntegralFisGasGeneratedSifgrs
    block = pellet
    execute_on = 'linear'
  []
  [fis_gas_released]
    type = ElementIntegralFisGasReleasedSifgrs
    block = pellet
    execute_on = 'linear'
  []
  [fis_gas_released_percentage]
    type = FGRPercent
    fission_gas_generated = fis_gas_produced
    fission_gas_released = fis_gas_released
    execute_on = 'linear'
  []
  [fis_gas_grain]
    type = ElementIntegralFisGasGrainSifgrs
    block = pellet
    execute_on = 'linear'
  []
  [fis_gas_boundary]
    type = ElementIntegralFisGasBoundarySifgrs
    block = pellet
    execute_on = 'linear'
  []
  [gas_volume]
    type = InternalVolume
    boundary = 9
    execute_on = 'initial linear'
  []
  [flux_from_fuel]
    type = SideDiffusiveFluxIntegral
    variable = temp
    boundary = 10
    diffusivity = thermal_conductivity
  []
  [rod_total_power]
    type = ElementIntegralPower
    variable = temp
    fission_rate = fission_rate
    block = pellet
  []
  [rod_input_power]
    type = FunctionValuePostprocessor
    function = average_power_history
    scale_factor = 0.55 # rod height
  []
  [average_vonMises_fuel]
    type = ElementAverageValue
    variable = vonmises_stress
    block = pellet
  []
  [average_vonMises_clad]
    type = ElementAverageValue
    variable = vonmises_stress
    block = clad
  []
  [average_strain_rr_fuel]
    type = ElementAverageValue
    variable = radial_strain
    block = pellet
  []
  [average_strain_rr_clad]
    type = ElementAverageValue
    variable = radial_strain
    block = clad
  []
  [average_creep_strain_clad]
    type = ElementAverageValue
    variable = effective_creep_strain
    block = clad
  []
  [ave_pore]
    type = ElementAverageValue
    variable = pore
  []
  [max_pore]
    type = NodalExtremeValue
    value_type = max
    variable = pore
  []
  [min_pore]
    type = NodalExtremeValue
    value_type = min
    variable = pore
  []
[]

[PerformanceMetricOutputs]
[]

[Outputs]
  perf_graph = true
  exodus = true
  color = true
  csv = true
  [console]
    type = Console
    max_rows = 25
  []
  [chkfile]
    type = CSV
    execute_on = FINAL
    show = 'ave_temp_interior fis_gas_released_percentage max_pore'
  []
[]
[MultiApps]
  [sub]
    type = TransientMultiApp
    app_type = BisonApp
    execute_on = timestep_end
    catch_up = true
    max_catch_up_steps = 10
    positions_file = positions.txt
    input_files = MK-II_sub_old_bubble_gb_lim.i
  []
[]
[Transfers]
  [temp_to_sub]
    type = MultiAppProjectionTransfer
    to_multi_app = sub
    source_variable = temp
    variable = temp
  []
  [pore_from_sub]
    type = MultiAppGeometricInterpolationTransfer
    from_multi_app = sub
    source_variable = pore
    variable = pore
  []
[]
[Debug]
  show_var_residual_norms = true
  show_var_residual = 'temp disp_x disp_y'
[]

(assessment/MOX/JOYO/MK-I/analysis/MK-I_50MW_master_new_bubble_gb_lim.i)


initial_fuel_density = 10836.8

[GlobalParams]
  density = ${initial_fuel_density}
  initial_porosity = 0.065
  order = SECOND
  family = LAGRANGE
  energy_per_fission = 3.2e-11 # J/fission
  volumetric_locking_correction = false
  displacements = 'disp_x disp_y'
[]
[Problem]
  type = ReferenceResidualProblem
  extra_tag_vectors = 'ref'
  reference_vector = 'ref'
[]
[Mesh]
  coord_type = RZ
  [smeared_pellet_mesh]
    type = FuelPinMeshGenerator
    pellet_quantity = 1
    pellet_height = 0.6
    pellet_outer_radius = 0.0027
    pellet_mesh_density = customize
    clad_mesh_density = customize
    clad_gap_width = 0.000100
    clad_thickness = 0.00035
    clad_bot_gap_height = 1.0e-3
    bottom_clad_height = 2.24e-3
    top_clad_height = 2.24e-3
    clad_top_gap_height = 0.599
    elem_type = QUAD8
    nx_c = 4
    ny_c = 200
    nx_p = 20
    ny_p = 200
    ny_cu = 3
    ny_cl = 3
  []
  patch_size = 50
  patch_update_strategy = iteration
  partitioner = centroid
  centroid_partitioner_direction = y
[]
[UserObjects]
  [pin_geometry]
    type = FuelPinGeometry
  []
[]
[Variables]
  [temp]
    initial_condition = 295.0
  []
[]
[AuxVariables]
  [pore]
  []
  [fission_rate]
    block = pellet
  []
  [burnup]
    block = pellet
  []
  [gas_gen_3]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [gas_grn_3]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [gas_bdr_3]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [gas_rel_3]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [bbl_bdr_2]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [vcn_bdr_2]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [atm_bbl_bdr]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [vcn_bbl_bdr]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [prs_bbl_bdr]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [prseq_bbl_bdr]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [rad_bbl_bdr]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [vol_bbl_bdr]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [GBCoverage]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [eff_diff_coeff]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [deltav_v0_bd]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [radial_strain]
    order = CONSTANT
    family = MONOMIAL
  []
  [effective_creep_strain]
    block = clad
    order = CONSTANT
    family = MONOMIAL
  []
  [gap_cond]
    order = CONSTANT
    family = MONOMIAL
  []
[]

[Functions]
  [power_history] #related to the LHGR at the midplane
    type = PiecewiseLinear
    x = '0 70000 17153028'
    y = '0  25577  25577'
  []
  [fast_neutron_flux_function]
    type = PiecewiseLinear
    x = '0 70000 17153028'
    y = '0 1.2e+19 1.2e+19'
  []
  [axial_peaking_factors]
    type = PiecewiseBilinear
    x = '0  0.071  0.146  0.221  0.296  0.37  0.443  0.566'
    y = '0 17153028'
    z = '0.889  1.041  1.152  1.173  1.129  0.971  0.782  0.672 0.889  1.041  1.152  1.173  1.129  0.971  0.782  0.672'
    scale_factor = 1
    axis = 1
  []
  [q]
    type = CompositeFunction
    functions = 'power_history axial_peaking_factors'
  []
  [average_power_history]
    type = PiecewiseLinear
    x = '0 70000 17153028'
    y = '0 21000 21000'
  []
  [clad_surface_temp]
    type = PiecewiseBilinear
    x = '0  0.075  0.15  0.225  0.3  0.375  0.45  0.525  0.6'
    y = '0 17153028'
    z = '295 295 295 295 295 295 295 295 295 499.9  509.1  517.8  525.42  532.71  540.29  547.7  552.3  554.81'
    scale_factor = 1
    axis = 1
  []
  [pressure_ramp]
    type = PiecewiseLinear
    x = '0 1'
    y = '1 1'
  []
[]
[Physics/SolidMechanics/QuasiStatic]
  [pellets]
    block = pellet
    add_variables = true
    strain = FINITE
    eigenstrain_names = 'fuel_thermal_strain fuel_volumetric_strain'
    generate_output = 'vonmises_stress stress_xx stress_yy stress_zz'
    use_finite_deform_jacobian = true
    extra_vector_tags = 'ref'
  []
  [clad]
    block = clad
    add_variables = true
    strain = FINITE
    eigenstrain_names = 'clad_thermal_eigenstrain'
    generate_output = 'vonmises_stress stress_xx stress_yy stress_zz'
    use_finite_deform_jacobian = true
    extra_vector_tags = 'ref'
  []
[]
[Kernels]
  [gravity]
    type = Gravity
    variable = disp_y
    value = -9.81
    extra_vector_tags = 'ref'
  []
  [heat]
    type = HeatConduction
    variable = temp
    extra_vector_tags = 'ref'
  []
  [heat_ie]
    type = HeatConductionTimeDerivative
    variable = temp
    extra_vector_tags = 'ref'
  []
  [heat_source]
    type = NeutronHeatSource
    variable = temp
    block = pellet
    fission_rate = fission_rate
    extra_vector_tags = 'ref'
  []
[]
[AuxKernels]
  [fission_rate]
    type = FissionRateGeneral
    fission_rate_formulation = MOX
    variable = fission_rate
    block = pellet
    initial_porosity = 0.065
    axial_power_profile = axial_peaking_factors
    rod_ave_lin_pow = power_history
    pellet_diameter = 0.0054
    execute_on = timestep_begin
    porosity = pore
  []
  [burnup]
    type = BurnupAux
    block = pellet
    fission_rate = fission_rate
    variable = burnup
    execute_on = timestep_begin
  []
  [fggen]
    type = MaterialRealAux
    variable = gas_gen_3
    property = gas_concentration_generated_total
    execute_on = timestep_end
  []
  [fggrn]
    type = MaterialRealAux
    variable = gas_grn_3
    property = gas_concentration_intra_total
    execute_on = timestep_end
  []
  [fgbdr]
    type = MaterialRealAux
    variable = gas_bdr_3
    property = gas_concentration_GB_bubble_volume
    execute_on = timestep_end
  []
  [fgrel]
    type = MaterialRealAux
    variable = gas_rel_3
    property = gas_concentration_release_total
    execute_on = timestep_end
  []
  [nbbl2]
    type = MaterialRealAux
    variable = bbl_bdr_2
    property = bubble_GB_surface_density
    execute_on = timestep_end
  []
  [nvcn2]
    type = MaterialRealAux
    variable = vcn_bdr_2
    property = vacancy_concentration_GB_surface
    execute_on = timestep_end
  []
  [atmbbl]
    type = MaterialRealAux
    variable = atm_bbl_bdr
    property = atom_per_bubble_GB
    execute_on = timestep_end
  []
  [vcnbbl]
    type = MaterialRealAux
    variable = vcn_bbl_bdr
    property = vacancy_per_bubble_GB
    execute_on = timestep_end
  []
  [prsbbl]
    type = MaterialRealAux
    variable = prs_bbl_bdr
    property = bubble_GB_pressure
    execute_on = timestep_end
  []
  [prseqbbl]
    type = MaterialRealAux
    variable = prseq_bbl_bdr
    property = bubble_GB_pressure_equilibrium
    execute_on = timestep_end
  []
  [radbbl]
    type = MaterialRealAux
    variable = rad_bbl_bdr
    property = bubble_radius_GB
    execute_on = timestep_end
  []
  [volbbl]
    type = MaterialRealAux
    variable = vol_bbl_bdr
    property = bubble_GB_volume
    execute_on = timestep_end
  []
  [frcvrg]
    type = MaterialRealAux
    variable = GBCoverage
    property = GBCoverage
    execute_on = timestep_end
  []
  [diffc]
    type = MaterialRealAux
    variable = eff_diff_coeff
    property = eff_diff_coeff
    execute_on = timestep_end
  []
  [dvv0bd]
    type = MaterialRealAux
    variable = deltav_v0_bd
    property = deltav_v0_bubble_GB
    execute_on = timestep_end
  []
  [radial_strain]
    type = RankTwoAux
    rank_two_tensor = total_strain
    variable = radial_strain
    index_i = 0
    index_j = 0
    execute_on = timestep_end
  []
  [effective_creep_strain]
    type = MaterialRealAux
    property = effective_creep_strain
    variable = effective_creep_strain
    block = clad
    execute_on = timestep_end
  []
  [conductance]
    type = MaterialRealAux
    property = gap_conductance
    variable = gap_cond
    boundary = 10
    execute_on = 'linear'
  []
[]
[Contact]
  [pellet_clad_mechanical]
    primary = 5
    secondary = 10
    formulation = kinematic
    model = frictionless
    penalty = 1e7
  []
[]
[ThermalContact]
  [thermal_contact]
    type = GasGapHeatTransfer
    variable = temp
    primary = 5
    secondary = 10
    initial_moles = initial_moles
    gas_released = fis_gas_released
    contact_pressure = contact_pressure
    quadrature = true
  []
[]
[BCs]
  [no_x_all]
    type = DirichletBC
    variable = disp_x
    boundary = '12'
    value = 0.0
  []
  [no_y_clad_bottom]
    type = DirichletBC
    variable = disp_y
    boundary = 1
    value = 0.0
  []
  [no_y_fuel_bottom]
    type = DirichletBC
    variable = disp_y
    boundary = 20
    value = 0.0
  []
  [temp_clad_out]
    type = FunctionDirichletBC
    variable = temp
    boundary = '2'
    function = clad_surface_temp
  []
  [Pressure]
    [coolantPressure]
      boundary = '1 2 3'
      factor = 101325
      function = pressure_ramp
    []
  []
  [PlenumPressure]
    [plenumPressure]
      boundary = 9
      initial_pressure = 300000
      startup_time = 0
      R = 8.3143
      output_initial_moles = initial_moles
      temperature = ave_temp_interior
      volume = gas_volume
      material_input = fis_gas_released
      output = plenum_pressure
    []
  []
[]
[Materials]
  [fast_neutron_flux]
    type = FastNeutronFlux
    calculate_fluence = true
    block = clad
    flux_function = fast_neutron_flux_function
  []
  [fuel_thermal]
    type = MAMOXThermal
    block = pellet
    temperature = temp
    Am_content = 0.0
    Np_content = 0.0
    porosity = pore
    output_properties = 'thermal_conductivity'
  []
  [fuel_elasticity_tensor]
    type = MAMOXElasticityTensor
    block = pellet
  []
  [elastic_stress]
    type = ComputeFiniteStrainElasticStress
    block = pellet
  []
  [fuel_thermal_expansion]
    type = MAMOXThermalExpansionEigenstrain
    block = pellet
    temperature = temp
    stress_free_temperature = 295.0
    oxygen_to_metal_ratio = 1.98
    eigenstrain_name = fuel_thermal_strain
  []
  [fuel_volumetric_swelling]
    type = UO2VolumetricSwellingEigenstrain
    gas_swelling_model_type = SIFGRS
    block = pellet
    temperature = temp
    burnup = burnup
    initial_fuel_density = 10836.8
    eigenstrain_name = fuel_volumetric_strain
  []
  [clad_thermal]
    type = SS316Thermal
    block = clad
    temperature = temp
  []
  [clad_density]
    type = StrainAdjustedDensity
    block = clad
    strain_free_density = 8000
  []
  [clad_elasticity_tensor]
    type = SS316ElasticityTensor
    block = clad
    temperature = temp
    elastic_constants_model = legacy_ifr
  []
  [thermal_expansion]
    type = SS316ThermalExpansionEigenstrain
    block = clad
    temperature = temp
    stress_free_temperature = 295.0
    eigenstrain_name = clad_thermal_eigenstrain
  []
  [clad_ss316creep]
    type = SS316CreepUpdate
    block = clad
    temperature = temp
    fast_neutron_flux = fast_neutron_flux
  []
  [clad_stress]
    type = ComputeMultipleInelasticStress
    tangent_operator = elastic
    inelastic_models = 'clad_ss316creep'
    block = clad
  []
  [fission_gas_release]
    type = UO2Sifgrs
    block = pellet
    temperature = temp
    burnup = burnup
    diff_coeff_option = TURNBULL_D1_4D2_4D3
    fission_rate = fission_rate
    grain_radius_const = 8.01e-6 #I'm keeping the grain radius const because the grain growth in MOX is probably different due to high Temp
    bubble_gb_limit = 1.0e+11
  []
  [fuel_density]
    type = StrainAdjustedDensity
    block = pellet
    strain_free_density = ${initial_fuel_density}
  []
[]
[Preconditioning]
  [SMP]
    type = SMP
    full = true
  []
[]
[Executioner]
  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'
  fixed_point_abs_tol = 1e-5
  fixed_point_rel_tol = 1e-5
  fixed_point_max_its = 1
  l_max_its = 70
  l_tol = 8e-3
  nl_max_its = 70
  nl_rel_tol = 1e-5
  nl_abs_tol = 1e-5
  start_time = 0
  n_startup_steps = 1
  end_time = 17153028
  dtmax = 1e6
  dtmin = 0.25
  [TimeStepper]
    type = IterationAdaptiveDT
    dt = 5000
    optimal_iterations = 15
    iteration_window = 2
    linear_iteration_ratio = 100
    growth_factor = 2
    cutback_factor = .5
    force_step_every_function_point = true
    timestep_limiting_function = power_history
  []
[]

[Postprocessors]
  [ave_temp_interior]
    type = SideAverageValue
    boundary = 9
    variable = temp
    execute_on = 'initial linear'
  []
  [average_burnup]
    type = ElementAverageValue
    block = pellet
    variable = burnup
  []
  [clad_inner_vol]
    type = InternalVolume
    boundary = 7
    execute_on = 'initial timestep_end'
  []
  [pellet_volume]
    type = InternalVolume
    boundary = 8
    execute_on = 'initial timestep_end'
  []
  [avg_clad_temp]
    type = SideAverageValue
    boundary = 7
    variable = temp
    execute_on = 'initial timestep_end'
  []
  [fis_gas_produced]
    type = ElementIntegralFisGasGeneratedSifgrs
    block = pellet
    execute_on = 'linear'
  []
  [fis_gas_released]
    type = ElementIntegralFisGasReleasedSifgrs
    block = pellet
    execute_on = 'linear'
  []
  [fis_gas_released_percentage]
    type = FGRPercent
    fission_gas_generated = fis_gas_produced
    fission_gas_released = fis_gas_released
    execute_on = 'linear'
  []
  [fis_gas_grain]
    type = ElementIntegralFisGasGrainSifgrs
    block = pellet
    execute_on = 'linear'
  []
  [fis_gas_boundary]
    type = ElementIntegralFisGasBoundarySifgrs
    block = pellet
    execute_on = 'linear'
  []
  [gas_volume]
    type = InternalVolume
    boundary = 9
    execute_on = 'initial linear'
  []
  [flux_from_fuel]
    type = SideDiffusiveFluxIntegral
    variable = temp
    boundary = 10
    diffusivity = thermal_conductivity
  []
  [rod_total_power]
    type = ElementIntegralPower
    variable = temp
    fission_rate = fission_rate
    block = pellet
  []
  [rod_input_power]
    type = FunctionValuePostprocessor
    function = average_power_history
    scale_factor = 0.6 # rod height
  []
  [average_vonMises_fuel]
    type = ElementAverageValue
    variable = vonmises_stress
    block = pellet
  []
  [average_vonMises_clad]
    type = ElementAverageValue
    variable = vonmises_stress
    block = clad
  []
  [average_strain_rr_fuel]
    type = ElementAverageValue
    variable = radial_strain
    block = pellet
  []
  [average_strain_rr_clad]
    type = ElementAverageValue
    variable = radial_strain
    block = clad
  []
  [average_creep_strain_clad]
    type = ElementAverageValue
    variable = effective_creep_strain
    block = clad
  []
  [ave_pore]
    type = ElementAverageValue
    variable = pore
  []
  [max_pore]
    type = NodalExtremeValue
    value_type = max
    variable = pore
  []
  [min_pore]
    type = NodalExtremeValue
    value_type = min
    variable = pore
  []
[]

[PerformanceMetricOutputs]
[]

[Outputs]
  perf_graph = true
  exodus = true
  color = true
  csv = true
  [console]
    type = Console
    max_rows = 25
  []
  [chkfile]
    type = CSV
    execute_on = FINAL
    show = 'ave_temp_interior fis_gas_released_percentage max_pore'
  []
[]
[MultiApps]
  [sub]
    type = TransientMultiApp
    app_type = BisonApp
    execute_on = timestep_end
    catch_up = true
    max_catch_up_steps = 10
    positions_file = positions.txt
    input_files = MK-I_50MW_sub_new_bubble_gb_lim.i
  []
[]
[Transfers]
  [temp_to_sub]
    type = MultiAppProjectionTransfer
    to_multi_app = sub
    source_variable = temp
    variable = temp
  []
  [pore_from_sub]
    type = MultiAppGeometricInterpolationTransfer
    from_multi_app = sub
    source_variable = pore
    variable = pore
  []
[]
[Debug]
  show_var_residual_norms = true
  show_var_residual = 'temp disp_x disp_y'
[]

(assessment/MOX/JOYO/B14/PTM003/analysis/b14_ptm003_2DRZ_t.i)


initial_fuel_density = 11172.82

[GlobalParams]
  density = ${initial_fuel_density}
  order = SECOND
  family = LAGRANGE
  energy_per_fission = 3.2e-11 # J/fission
  volumetric_locking_correction = false
  displacements = 'disp_x disp_y'
[]
[Problem]
  type = ReferenceResidualProblem
  extra_tag_vectors = 'ref'
  reference_vector = 'ref'
[]
[Mesh]
  coord_type = RZ
  [smeared_pellet_mesh]
    type = FuelPinMeshGenerator
    pellet_quantity = 1
    pellet_height = 0.4
    pellet_outer_radius = 0.0027
    pellet_mesh_density = customize
    clad_mesh_density = customize
    clad_gap_width = 0.00008
    clad_thickness = 0.00047
    clad_bot_gap_height = 1.0e-3
    bottom_clad_height = 2.24e-3
    top_clad_height = 2.24e-3
    clad_top_gap_height = 0.685
    elem_type = QUAD8
    nx_c = 4
    ny_c = 100
    nx_p = 20
    ny_p = 100
    ny_cu = 3
    ny_cl = 3
  []
  patch_size = 50
  patch_update_strategy = iteration
  partitioner = centroid
  centroid_partitioner_direction = y
[]
[UserObjects]
  [pin_geometry]
    type = FuelPinGeometry
  []
[]
[Variables]
  [temp]
    initial_condition = 295.0
  []
[]
[AuxVariables]
  [pore]
  []
  [fission_rate]
  []
  [burnup]
  []
  [radial_strain]
    order = CONSTANT
    family = MONOMIAL
  []
  [effective_creep_strain]
    block = clad
    order = CONSTANT
    family = MONOMIAL
  []
[]
[Functions]
  [power_history] #related to the LHGR at the midplane
    type = PiecewiseLinear
    x = ' 0  72000  158040  160200  246600  248400  249000.012  251280'
    y = ' 0  39814.5  39814.5  44289.3  44289.3  53927.4  53927.4  0'
  []
  [fast_neutron_flux_function]
    type = PiecewiseLinear
    x = '0 251280'
    y = '3.3e+15 3.3e+15'
  []
  [f_temp_out_clad]
    type = PiecewiseBilinear
    x = '0.000175  0.0464075  0.0843675  0.1075625  0.152025  0.1994625  0.2464725  0.2947475  0.356915  0.43356  0.49848  0.625  0.700475  0.797485  0.8723425  0.96'
    y = '0 251280'
    z = '295  295  295  295  295  295  295  295  295  295  295  295  295  295  295  295 634.94  662.273  676.998  686.217  706.339  727  743.358  758.311  780.069  799.077  815.576  846.374  860.233  875.494  882.809  889.8'
    scale_factor = 1
    axis = 1
  []
  [axial_peaking_factors]
    type = PiecewiseBilinear
    x = '0.000175  0.0464075  0.0843675  0.1075625  0.152025  0.1994625  0.2464725  0.2947475  0.356915  0.43356  0.49848  0.625  0.700475  0.797485  0.8723425  0.96'
    y = '0 251280'
    z = '0.751  0.752  0.767  0.796  0.82  0.852  0.875  0.915  0.944  0.963  0.988  1  0.985  0.955  0.913  0.846 0.751  0.752  0.767  0.796  0.82  0.852  0.875  0.915  0.944  0.963  0.988  1  0.985  0.955  0.913  0.846'
    scale_factor = 1
    axis = 1
  []
  [q]
    type = CompositeFunction
    functions = 'power_history axial_peaking_factors'
  []
  [average_power_history]
    type = PiecewiseLinear
    x = ' 0  72000  158040  160200  246600  248400  249000.012  251280'
    y = ' 0  34700  34700  38600  38600  47000  47000  0'
  []
[]
[Physics/SolidMechanics/QuasiStatic]
  [pellets]
    block = pellet
    add_variables = true
    strain = FINITE
    eigenstrain_names = 'fuel_thermal_strain fuel_volumetric_strain'
    generate_output = 'vonmises_stress stress_xx stress_yy stress_zz'
    use_finite_deform_jacobian = true
    extra_vector_tags = 'ref'
  []
  [clad]
    block = clad
    add_variables = true
    strain = FINITE
    eigenstrain_names = 'clad_thermal_eigenstrain'
    generate_output = 'vonmises_stress stress_xx stress_yy stress_zz'
    use_finite_deform_jacobian = true
    extra_vector_tags = 'ref'
  []
[]
[Kernels]
  [gravity]
    type = Gravity
    variable = disp_y
    value = -9.81
    extra_vector_tags = 'ref'
  []
  [heat]
    type = HeatConduction
    variable = temp
    extra_vector_tags = 'ref'
  []
  [heat_ie]
    type = HeatConductionTimeDerivative
    variable = temp
    extra_vector_tags = 'ref'
  []
  [heat_source]
    type = NeutronHeatSource
    variable = temp
    fission_rate = fission_rate
    extra_vector_tags = 'ref'
  []
[]
[AuxKernels]
  [fission_rate]
    type = FissionRateGeneral
    fission_rate_formulation = MOX
    variable = fission_rate
    block = pellet
    initial_porosity = 0.1409
    axial_power_profile = axial_peaking_factors
    rod_ave_lin_pow = power_history
    pellet_diameter = 0.0054
    execute_on = timestep_begin
    porosity = pore
    energy_per_fission = 3.2e-11
  []
  [burnup]
    type = BurnupAux
    block = pellet
    fission_rate = fission_rate
    variable = burnup
    execute_on = timestep_begin
  []
  [radial_strain]
    type = RankTwoAux
    rank_two_tensor = total_strain
    variable = radial_strain
    index_i = 0
    index_j = 0
    execute_on = timestep_end
  []
  [effective_creep_strain]
    type = MaterialRealAux
    property = effective_creep_strain
    variable = effective_creep_strain
    block = clad
    execute_on = timestep_end
  []
[]
[Contact]
  [pellet_clad_mechanical]
    primary = 5
    secondary = 10
    formulation = kinematic
    model = frictionless
    penalty = 1e7
  []
[]
[ThermalContact]
  [thermal_contact]
    type = GasGapHeatTransfer
    variable = temp
    primary = 5
    secondary = 10
    initial_moles = initial_moles
    gas_released = fis_gas_released
    contact_pressure = contact_pressure
    quadrature = true
  []
[]
[BCs]
  [temp_clad_outside]
    type = FunctionDirichletBC
    variable = temp
    function = f_temp_out_clad
    boundary = 2
  []
  [no_x_all]
    type = DirichletBC
    variable = disp_x
    boundary = 12
    value = 0.0
  []
  [no_y_clad_bottom]
    type = DirichletBC
    variable = disp_y
    boundary = 1
    value = 0.0
  []
  [no_y_fuel_bottom]
    type = DirichletBC
    variable = disp_y
    boundary = 20
    value = 0.0
  []
  [PlenumPressure]
    [plenumPressure]
      boundary = 9
      initial_pressure = 101325
      startup_time = 0
      R = 8.3143
      output_initial_moles = initial_moles
      temperature = ave_temp_interior
      volume = gas_volume
      material_input = fis_gas_released
      output = plenum_pressure
    []
  []
[]
[Materials]
  [fast_neutron_flux]
    type = FastNeutronFlux
    calculate_fluence = true
    block = clad
    flux_function = fast_neutron_flux_function
  []
  [fuel_thermal]
    type = MAMOXThermal
    temperature = temp
    porosity = pore
    block = pellet
    Am_content = 0.0237
    oxy_to_metal_ratio = 1.961
    output_properties = 'thermal_conductivity'
  []
  [fuel_density]
    type = StrainAdjustedDensity
    block = pellet
    strain_free_density = ${initial_fuel_density}
  []
  [fuel_elasticity_tensor]
    type = MAMOXElasticityTensor
    block = pellet
  []
  [elastic_stress]
    type = ComputeFiniteStrainElasticStress
    block = pellet
    #outputs = exodus
  []
  [fuel_thermal_expansion]
    type = MAMOXThermalExpansionEigenstrain
    block = pellet
    temperature = temp
    stress_free_temperature = 295.0
    oxygen_to_metal_ratio = 1.97
    eigenstrain_name = fuel_thermal_strain
  []
  [fuel_volumetric_swelling]
    type = UO2VolumetricSwellingEigenstrain
    gas_swelling_model_type = SIFGRS
    block = pellet
    temperature = temp
    burnup = burnup
    initial_fuel_density = 11172.82
    eigenstrain_name = fuel_volumetric_strain
  []
  [clad_thermal]
    type = SS316Thermal
    block = clad
    temperature = temp
  []
  [clad_density]
    type = StrainAdjustedDensity
    block = clad
    strain_free_density = 8000
  []
  [clad_elasticity_tensor]
    type = SS316ElasticityTensor
    block = clad
    temperature = temp
    elastic_constants_model = legacy_ifr
  []
  [thermal_expansion]
    type = SS316ThermalExpansionEigenstrain
    block = clad
    temperature = temp
    stress_free_temperature = 295.0
    eigenstrain_name = clad_thermal_eigenstrain
  []
  [clad_ss316creep]
    type = SS316CreepUpdate
    block = clad
    temperature = temp
    fast_neutron_flux = fast_neutron_flux
  []
  [clad_stress]
    type = ComputeMultipleInelasticStress
    tangent_operator = elastic
    inelastic_models = 'clad_ss316creep'
    block = clad
  []
  [fission_gas_release]
    type = UO2Sifgrs
    block = pellet
    temperature = temp
    burnup = burnup
    fission_rate = fission_rate
    grain_radius_const = 11e-06
    bubble_gb_limit = 1.0e+11
  []
[]
[Preconditioning]
  [SMP]
    type = SMP
    full = true
  []
[]
[Executioner]
  type = Transient
  solve_type = 'PJFNK'
  petsc_options = '-snes_ksp_ew'
  petsc_options_iname = '-pc_type -sub_pc_type'
  petsc_options_value = 'asm lu'
  line_search = 'none'
  fixed_point_max_its = 1
  fixed_point_abs_tol = 1e-3
  fixed_point_rel_tol = 1e-3
  l_max_its = 50
  l_tol = 8e-3
  nl_max_its = 50
  nl_rel_tol = 1e-3
  nl_abs_tol = 1e-3
  start_time = 0
  n_startup_steps = 1
  end_time = 251280
  dtmax = 10000
  dtmin = 0.25
  automatic_scaling = true
  compute_scaling_once = false
  [TimeStepper]
    type = IterationAdaptiveDT
    dt = 1e1
    optimal_iterations = 15
    iteration_window = 2
    linear_iteration_ratio = 100
    growth_factor = 2
    cutback_factor = .5
    force_step_every_function_point = true
    timestep_limiting_function = power_history
  []
[]
[Postprocessors]
  [ave_temp_interior]
    type = ElementAverageValue
    variable = temp
    execute_on = 'initial linear'
  []
  [average_burnup]
    type = ElementAverageValue
    variable = burnup
  []
  [ave_pore]
    type = ElementAverageValue
    block = pellet
    variable = pore
  []
  [max_pore]
    type = NodalExtremeValue
    block = pellet
    value_type = max
    variable = pore
  []
  [min_pore]
    type = NodalExtremeValue
    block = pellet
    value_type = min
    variable = pore
  []
  [fis_gas_produced]
    type = ElementIntegralFisGasGeneratedSifgrs
    #  variable = temp
    block = pellet
    execute_on = 'linear'
  []
  [fis_gas_released]
    type = ElementIntegralFisGasReleasedSifgrs
    #  variable = temp
    block = pellet
    execute_on = 'linear'
  []
  [fis_gas_released_percentage]
    type = FGRPercent
    fission_gas_generated = fis_gas_produced
    fission_gas_released = fis_gas_released
    #  variable = temp
    execute_on = 'linear'
  []
  [gas_volume]
    type = InternalVolume
    boundary = 9
    execute_on = 'initial linear'
  []
  [rod_total_power]
    type = ElementIntegralPower
    variable = temp
    fission_rate = fission_rate
    block = pellet
  []
  [rod_input_power]
    type = FunctionValuePostprocessor
    function = average_power_history
    scale_factor = 0.4 # rod height
  []
[]
[VectorPostprocessors]
  [fuel_radial_temperature_Sample1]
    type = LineValueSampler
    variable = temp
    start_point = '0.0 0.283 0.0'
    end_point = '0.0027 0.283 0.0'
    num_points = 200
    execute_on = final
    sort_by = x
    outputs = line_plot
  []
  [radial_porosity_Sample1]
    type = LineValueSampler
    variable = pore
    start_point = '0.0 0.283 0.0'
    end_point = '0.0027 0.283 0.0'
    num_points = 200
    execute_on = final
    sort_by = x
    outputs = line_plot
  []
  [fuel_radial_temperature_Sample2]
    type = LineValueSampler
    variable = temp
    start_point = '0.0 0.347 0.0'
    end_point = '0.0027 0.347 0.0'
    num_points = 200
    execute_on = final
    sort_by = x
    outputs = line_plot
  []
  [radial_porosity_Sample2]
    type = LineValueSampler
    variable = pore
    start_point = '0.0 0.347 0.0'
    end_point = '0.0027 0.347 0.0'
    num_points = 200
    execute_on = final
    sort_by = x
    outputs = line_plot
  []
  [fuel_radial_temperature_Sample3]
    type = LineValueSampler
    variable = temp
    start_point = '0.0 0.2 0.0'
    end_point = '0.0027 0.2 0.0'
    num_points = 200
    execute_on = final
    sort_by = x
    outputs = line_plot
  []
  [radial_porosity_Sample3]
    type = LineValueSampler
    variable = pore
    start_point = '0.0 0.2 0.0'
    end_point = '0.0027 0.2 0.0'
    num_points = 200
    execute_on = final
    sort_by = x
    outputs = line_plot
  []
[]

[PerformanceMetricOutputs]
[]

[Outputs]
  perf_graph = true
  exodus = true
  color = true
  csv = true
  [console]
    type = Console
    max_rows = 25
  []
  [line_plot]
    type = CSV
    execute_on = 'FINAL'
    time_step_interval =  1
    file_base = 1d
    create_final_symlink = true
  []
  [chkfile]
    type = CSV
    execute_on = FINAL
    show = 'ave_temp_interior max_pore'
  []
[]
[MultiApps]
  [sub]
    type = TransientMultiApp
    app_type = BisonApp
    execute_on = timestep_end
    catch_up = true
    max_catch_up_steps = 10
    positions_file = positions.txt
    input_files = b14_ptm003_pore.i
  []
[]
[Transfers]
  [temp_to_sub]
    type = MultiAppProjectionTransfer
    to_multi_app = sub
    source_variable = temp
    variable = temp
  []
  [pore_from_sub]
    type = MultiAppGeometricInterpolationTransfer
    from_multi_app = sub
    source_variable = pore
    variable = pore
  []
[]
[Debug]
  show_var_residual_norms = true
[]

(assessment/nitride/JOYO/JOYO_L4C4/analysis/JOYO_Pin_base.i)

[GlobalParams]
  order = FIRST
  energy_per_fission = ${energy_per_fission}
  displacements = 'disp_x disp_y'
  value_range_behavior = IGNORE
[]

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

[Mesh]
  coord_type = RZ
  [smeared_pellet_mesh]
    type = FuelPinMeshGenerator
    clad_thickness = ${cladding_thickness}
    pellet_outer_radius = ${fuel_radius}
    pellet_height = ${fuel_height}
    clad_top_gap_height = ${plenum_height}
    clad_gap_width = ${cladding_gap_width}
    bottom_clad_height = ${cladding_bottom_top_plug_length}
    top_clad_height = ${cladding_bottom_top_plug_length}
    clad_bot_gap_height = ${gap_bottom_length}
    clad_mesh_density = customize
    pellet_mesh_density = customize
    nx_p = 6
    ny_p = 260
    nx_c = 4
    ny_c = 260
    ny_cu = 3
    ny_cl = 3
    pellet_quantity = 1
    elem_type = QUAD4
  []
  patch_size = 30
  patch_update_strategy = auto
  partitioner = centroid
  centroid_partitioner_direction = y
[]

[Variables]
  [temperature]
    initial_condition = ${initial_temperature}
  []
[]

[Functions]
  [power_history]
    type = PiecewiseLinear
    x = '0 ${time_end_ramp_up} ${time_start_ramp_down} ${time_end_ramp_down}'
    y = '0 ${avg_lin_power} ${avg_lin_power} 0'
  []
  [coolant_wall_temperature]
    type = ParsedFunction
    expression = 'if(t<${time_end_ramp_up}, ${initial_temperature} + t*(${a}*y^2 + ${b}*y + ${c} - ${initial_temperature})/(${time_end_ramp_up}),if(t<${time_start_ramp_down},${a}*y^2 + ${b}*y + ${c}, if(t<${time_end_ramp_down}, ${a}*y^2 + ${b}*y + ${c} + (t-${time_start_ramp_down})*(${initial_temperature}-(${a}*y^2 + ${b}*y + ${c}))/(${time_end_ramp_down}-${time_start_ramp_down}),${initial_temperature})))'
  []
  [coolant_press_ramp]
    type = PiecewiseLinear
    x = '0       ${run_time}'
    y = '${coolant_pressure} ${coolant_pressure}'
  []
  [plenum_pressure]
    type = ConstantFunction
    value = ${initial_plenum_pressure}
  []
  [fast_neutron_flux_function]
    type = PiecewiseLinear
    x = '0 ${time_end_ramp_up} ${run_time}'
    y = '0 ${fast_neutron_flux} ${fast_neutron_flux}'
  []
  [axial_peaking_factors]
    type = PiecewiseBilinear
    x = '${pos_1} ${pos_2} ${pos_3} ${pos_4} ${pos_5} ${pos_6} ${pos_7} ${pos_8} ${pos_9} ${pos_10} ${pos_11} ${pos_12} ${pos_13} ${pos_14} ${pos_15} ${pos_16}'
    y = '0 ${run_time}'
    z = '${pow_1} ${pow_2} ${pow_3} ${pow_4} ${pow_5} ${pow_6} ${pow_7} ${pow_8} ${pow_9} ${pow_10} ${pow_11} ${pow_12} ${pow_13} ${pow_14} ${pow_15} ${pow_16} ${pow_1} ${pow_2} ${pow_3} ${pow_4} ${pow_5} ${pow_6} ${pow_7} ${pow_8} ${pow_9} ${pow_10} ${pow_11} ${pow_12} ${pow_13} ${pow_14} ${pow_15} ${pow_16}'
    scale_factor = 1
    axis = 1
  []
[]

[Physics/SolidMechanics/QuasiStatic]
  [fuel]
    block = pellet
    strain = FINITE
    add_variables = true
    generate_output = 'firstinv_strain stress_xx stress_yy stress_zz vonmises_stress hydrostatic_stress creep_strain_xx creep_strain_yy creep_strain_zz elastic_strain_xx elastic_strain_yy elastic_strain_zz strain_xx strain_yy strain_zz'
    extra_vector_tags = 'ref'
    eigenstrain_names = 'fuel_thermal_expansion solid_swelling_eigenstrain'
    temperature = temperature
  []
  [clad]
    block = clad
    strain = FINITE
    incremental = true
    add_variables = true
    generate_output = 'firstinv_strain 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'
    extra_vector_tags = 'ref'
    eigenstrain_names = 'clad_thermal_expansion'
    temperature = temperature
  []
[]

[Kernels]
  [gravity]
    type = Gravity
    variable = disp_y
    value = -9.81
    extra_vector_tags = 'ref'
    block = 'pellet clad'
  []
  [heat]
    type = HeatConduction
    variable = temperature
    extra_vector_tags = 'ref'
    block = 'pellet clad'
  []
  [heat_ie]
    type = HeatConductionTimeDerivative
    variable = temperature
    extra_vector_tags = 'ref'
    block = 'pellet clad'
  []
  [heat_source]
    type = FissionRateHeatSource
    variable = temperature
    fission_rate = fission_rate
    extra_vector_tags = 'ref'
    energy_deposited_in_fuel = 0.95
    block = 'pellet'
  []
[]

[ThermalContactMortar]
  [thermal_contact]
    secondary_variable = temperature
    primary_boundary = '5' # clad_inside_right
    secondary_boundary = '10' # pellet_outer_radial_surface
    initial_moles = initial_moles # coupling to a postprocessor which supplies the initial plenum/gap gas mass
  []
[]

[Contact]
  [mechanical]
    model = frictionless
    formulation = mortar
    primary = 5 # clad_inside_right
    secondary = 10 # pellet_outer_radial_surface
    c_normal = 1e+11
    correct_edge_dropping = true
  []
[]

[BCs]
  [no_x_all]
    type = DirichletBC
    variable = disp_x
    boundary = '12' # centerline
    value = 0.0
  []
  [no_y_clad]
    type = DirichletBC
    variable = disp_y
    boundary = '1' # clad_outside_bottom
    value = 0.0
  []
  [no_y_fuel]
    type = DirichletBC
    variable = disp_y
    boundary = '20' # bottom_of_bottom_pellet
    value = 0.0
  []
  [Pressure]
    [coolantPressure]
      boundary = '1 2 3' # clad_outside_bottom, clad_outside_right, clad_outside_top
      function = coolant_press_ramp
    []
  []
  [PlenumPressure] # apply plenum pressure on clad inner walls and pellet surfaces
    [plenumPressure]
      boundary = 9 # inside_surfaces
      initial_pressure = ${initial_plenum_pressure}
      startup_time = 0
      R = ${R}
      initial_temperature = ${initial_temperature}
      temperature = ave_temperature_interior # coupling to post processor to get gas temperature approximation
      output_initial_moles = initial_moles # coupling to post processor to get initial fill gas mass
      volume = gas_volume # coupling to post processor to get gas volume
      output = plenum_pressure # coupling to post processor to output plenum/gap pressure
    []
  []
  [clad_outer_temperature]
    type = FunctionDirichletBC
    boundary = '1 2 3' # clad_outside_bottom, clad_outside_right, clad_outside_top
    function = coolant_wall_temperature
    variable = temperature
  []
[]

[Materials]
  [porosity]
    type = GenericConstantMaterial
    block = pellet
    prop_names = porosity
    prop_values = ${initial_porosity}
    outputs = all
  []
  [fission_rate]
    type = FissionRate
    block = pellet
    rod_linear_power = power_history
    axial_power_profile = axial_peaking_factors
    pellet_radius = ${fuel_radius}
    outputs = all
  []
  [fuel_thermal]
    type = MNThermal
    block = pellet
    temperature = temperature
    porosity = porosity
    outputs = all
  []
  [fuel_elasticity_tensor]
    block = pellet
    type = MNElasticityTensor
    temperature = temperature
    porosity = porosity
    output_properties = 'youngs_modulus poissons_ratio'
    outputs = all
  []
  [fuel_thermal_expansion]
    block = pellet
    type = MNThermalExpansionEigenstrain
    eigenstrain_name = fuel_thermal_expansion
    stress_free_temperature = ${initial_temperature}
    temperature = temperature
  []
  [fuel_creep]
    block = pellet
    type = MNCreepUpdate
    max_inelastic_increment = 1e-4
    temperature = temperature
    porosity = porosity
    fission_rate = fission_rate
    outputs = all
    value_range_behavior = EXCEPTION
  []
  [burnup]
    type = Burnup
    block = pellet
    atoms_heavy_metal_per_volume = ${atoms_heavy_metal_per_volume}
    outputs = all
  []
  [burnup_swelling]
    type = BurnupDependentEigenstrain
    block = pellet
    eigenstrain_name = 'solid_swelling_eigenstrain'
    burnup = burnup
    outputs = all
  []
  [fuel_radial_return_stress]
    block = pellet
    type = ComputeMultipleInelasticStress
    inelastic_models = 'fuel_creep'
  []
  [fuel_density]
    block = pellet
    type = StrainAdjustedDensity
    strain_free_density = ${density}
  []

  [clad_elasticity_tensor]
    block = clad
    type = SS316ElasticityTensor
    temperature = temperature
    elastic_constants_model = legacy_ifr
  []
  [clad_thermal_expansion]
    block = clad
    type = SS316ThermalExpansionEigenstrain
    eigenstrain_name = clad_thermal_expansion
    stress_free_temperature = ${initial_temperature}
    temperature = temperature
  []
  [fast_neutron_flux]
    type = FastNeutronFlux
    calculate_fluence = false
    block = clad
    flux_function = fast_neutron_flux_function
  []
  [clad_creep]
    block = clad
    type = SS316CreepUpdate
    fast_neutron_flux = fast_neutron_flux
    temperature = temperature
  []
  [clad_stress]
    block = clad
    type = ComputeMultipleInelasticStress
    inelastic_models = 'clad_creep'
  []
  [clad_thermal]
    type = SS316Thermal
    block = clad
    temperature = temperature
  []
  [clad_density]
    block = clad
    type = DerivativeParsedMaterial
    coupled_variables = 'temperature'
    expression = '-4.454e-5*temperature^2 - 0.4297*temperature + 8089.4'
    property_name = 'density'
  []
[]

[Dampers]
  [disp_x]
    type = MaxIncrement
    variable = disp_x
    max_increment = 1e-4
    min_damping = 1e-4
  []
  [disp_y]
    type = MaxIncrement
    variable = disp_y
    max_increment = 1e-3
    min_damping = 1e-4
  []
  [temperature]
    type = MaxIncrement
    variable = temperature
    max_increment = 50
    min_damping = 1e-4
  []
[]

[Executioner]
  type = Transient
  # With mortar contact
  solve_type = 'PJFNK'
  petsc_options_iname = '-pc_type -pc_factor_mat_solver_package -mat_mffd_err -pc_factor_shift_type -pc_factor_shift_amount'
  petsc_options_value = 'lu    superlu_dist   1e-5          NONZERO               1e-15'
  snesmf_reuse_base = false
  verbose = true
  line_search = 'none'

  l_max_its = 60
  nl_max_its = 20
  nl_rel_tol = 5e-6
  nl_abs_tol = 5e-9

  end_time = ${run_time}
  dtmin = 1
  dtmax = 5e4

  [TimeStepper]
    type = IterationAdaptiveDT
    timestep_limiting_postprocessor = time_step_limit
    force_step_every_function_point = true
    timestep_limiting_function = power_history
    dt = 10
  []
[]

[Postprocessors]
  # elemental temperatures
  [temperature_fuel_avg]
    type = ElementAverageValue
    variable = temperature
    block = pellet
    execute_on = 'initial timestep_end'
  []
  [temperature_fuel_max]
    type = ElementExtremeValue
    variable = temperature
    block = pellet
  []
  [temperature_fuel_min]
    type = ElementExtremeValue
    variable = temperature
    block = pellet
    value_type = min
  []

  [temperature_clad_avg]
    type = ElementAverageValue
    variable = temperature
    block = clad
    execute_on = 'initial timestep_end'
  []
  [temperature_clad_max]
    type = ElementExtremeValue
    variable = temperature
    block = clad
  []
  [temperature_clad_min]
    type = ElementExtremeValue
    variable = temperature
    block = clad
    value_type = min
  []

  [ave_temperature_interior] # average temperature of the cladding interior and all pellet exteriors
    type = SideAverageValue
    boundary = '9' # inside_surfaces
    variable = temperature
    execute_on = 'initial linear'
  []
  [temperature_fuel_centerline_avg]
    type = AxisymmetricCenterlineAverageValue
    boundary = '12' # centerline
    variable = temperature
  []
  [temperature_fuel_centerline_max]
    type = NodalExtremeValue
    boundary = '12' # centerline
    variable = temperature
  []
  [temperature_fuel_centerline_min]
    type = NodalExtremeValue
    boundary = '12' # centerline
    variable = temperature
    value_type = min
  []
  [temperature_fuel_surface_avg]
    type = SideAverageValue
    boundary = '10' # pellet_outer_radial_surface
    variable = temperature
  []
  [temperature_fuel_surface_max]
    type = NodalExtremeValue
    boundary = '10' # pellet_outer_radial_surface
    variable = temperature
  []
  [temperature_fuel_surface_min]
    type = NodalExtremeValue
    boundary = '10' # pellet_outer_radial_surface
    variable = temperature
    value_type = min
  []

  # stresses
  [stress_vonmises_fuel_avg]
    type = ElementAverageValue
    variable = vonmises_stress
    block = pellet
  []
  [stress_vonmises_fuel_max]
    type = ElementExtremeValue
    variable = vonmises_stress
    block = pellet
  []
  [stress_vonmises_fuel_min]
    type = ElementExtremeValue
    variable = vonmises_stress
    value_type = min
    block = pellet
  []
  [stress_hydro_fuel_avg]
    type = ElementAverageValue
    variable = hydrostatic_stress
    block = pellet
  []
  [stress_hydro_fuel_max]
    type = ElementExtremeValue
    variable = hydrostatic_stress
    block = pellet
  []
  [stress_hydro_fuel_min]
    type = ElementExtremeValue
    variable = hydrostatic_stress
    value_type = min
    block = pellet
  []

  # strain information
  [strain_solid_swelling_fuel_avg]
    type = ElementAverageValue
    variable = solid_swelling
    block = pellet
  []
  [strain_volumetric_fuel_avg]
    type = ElementAverageValue
    variable = firstinv_strain
    block = pellet
  []
  [strain_axial_fuel_avg]
    type = ParsedPostprocessor
    pp_names = 'disp_y_fuel_top_surface_avg disp_y_fuel_bottom_surface_avg'
    expression = '(disp_y_fuel_top_surface_avg - disp_y_fuel_bottom_surface_avg) / ${fuel_height}'
  []
  [disp_y_fuel_top_surface_avg]
    type = SideAverageValue
    variable = disp_y
    boundary = '11' # top_of_top_pellet
  []
  [disp_y_fuel_top_surface_max]
    type = NodalExtremeValue
    variable = disp_y
    boundary = '11' # top_of_top_pellet
  []
  [disp_y_fuel_bottom_surface_avg]
    type = SideAverageValue
    variable = disp_y
    boundary = '20' # bottom_of_bottom_pellet
  []
  [disp_y_fuel_bottom_surface_max]
    type = NodalExtremeValue
    variable = disp_y
    boundary = '20' # bottom_of_bottom_pellet
  []
  [disp_x_fuel_radial_surface_max]
    type = NodalExtremeValue
    variable = disp_x
    boundary = '10' # pellet_outer_radial_surface
  []
  [disp_x_fuel_radial_surface_avg]
    type = SideAverageValue
    variable = disp_x
    boundary = '10' # pellet_outer_radial_surface
  []

  [burnup_avg]
    type = ElementAverageValue
    block = pellet
    variable = burnup
  []
  [burnup_max]
    type = ElementExtremeValue
    block = pellet
    variable = burnup
  []
  [fission_rate_avg]
    type = ElementAverageValue
    variable = fission_rate
    block = pellet
  []
  [porosity_fuel_avg]
    type = ElementAverageValue
    variable = porosity
    block = pellet
  []
  [porosity_fuel_max]
    type = ElementExtremeValue
    variable = porosity
    block = pellet
  []
  [porosity_fuel_min]
    type = ElementExtremeValue
    variable = porosity
    value_type = min
    block = pellet
  []

  [swelling_dia_percent]
    type = ParsedPostprocessor
    pp_names = 'disp_x_fuel_radial_surface_avg'
    expression = 'disp_x_fuel_radial_surface_avg / 3.14159 / ${fuel_radius} * 100'
  []
  [fuel_volume]
    type = VolumePostprocessor
    block = pellet
    execute_on = 'TIMESTEP_END INITIAL'
    use_displaced_mesh = true
  []
  [swelling_vol_percent]
    type = ParsedPostprocessor
    pp_names = 'fuel_volume'
    expression = '(fuel_volume - ${fuel_volume}) / ${fuel_volume} * 100'
  []
  [gas_volume]
    type = InternalVolume
    boundary = '9' # inside_surfaces
    execute_on = 'initial linear'
  []

  # extras
  [actual_time_step_limit]
    type = MaterialTimeStepPostprocessor
    outputs = none
    block = 'pellet'
  []
  [time_step_limit]
    type = ParsedPostprocessor
    expression = 'if(actual_time_step_limit > 1e6, 1e6, actual_time_step_limit)'
    pp_names = 'actual_time_step_limit'
  []
[]

[PerformanceMetricOutputs]
  outputs = 'performance_metrics performance_metrics_sync exodus console base_out'
[]

[Outputs]
  perf_graph = true
  csv = true
  sync_times = '1e2 1e3 5e3 1e4 5e4 1e5 5e5 1e6 2.5e6 5e6 1e7 2e7 ${fparse run_time -1e3} ${run_time}'
  file_base = '${group_name}_nominal_Pin'
  [base_out]
    type = CSV
    file_base = '${group_name}_Pin_base_out'
    show = 'burnup_max temperature_fuel_avg temperature_clad_avg ave_temperature_interior swelling_vol_percent disp_x_fuel_radial_surface_max'
    sync_only = true
  []
  [checkpoint]
    type = Checkpoint
    time_step_interval =  10
  []
  [exodus]
    type = Exodus
  []
  [sync]
    type = CSV
    file_base = '${group_name}_nominal_Pin_sync'
    sync_only = true
  []
  [console]
    type = Console
    show = 'time_step_size temperature_fuel_avg temperature_fuel_centerline_max stress_vonmises_fuel_max stress_hydro_fuel_max stress_hydro_fuel_min strain_axial_fuel_avg burnup_avg fission_rate_avg porosity_fuel_avg time_step_limit disp_x_fuel_radial_surface_max disp_x_fuel_radial_surface_avg swelling_dia_percent swelling_vol_percent fuel_volume'
  []
  [performance_metrics]
    type = CSV
    file_base = '${group_name}_nominal_performance_metrics_Pin'
    show = 'simulation_alive_time number_linear_iterations number_nonlinear_iterations time_step_size total_linear_iterations total_nonlinear_iterations physical_memory_use number_dofs number_nonlinear_variables residual_compute_time jacobian_compute_time'
  []
  [performance_metrics_sync]
    type = CSV
    sync_only = true
    file_base = '${group_name}_nominal_performance_metrics_Pin_sync'
    show = 'simulation_alive_time number_linear_iterations number_nonlinear_iterations time_step_size total_linear_iterations total_nonlinear_iterations physical_memory_use number_dofs number_nonlinear_variables residual_compute_time jacobian_compute_time'
  []
[]

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

(assessment/metallic_fuel/EBRII/X423/analysis/x423_lm_base.i)

[GlobalParams]
  density = ${fuel_density}
  order = FIRST
  family = LAGRANGE
  energy_per_fission = 3.2e-11  # J/fission
  volumetric_locking_correction = true
  displacements = 'disp_x disp_y'
  X_Pu = ${fuel_pu}
  X_Zr = ${fuel_zr}
[]

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

[Mesh]
  coord_type = RZ
  [gen]
    type = FIPDRodletMeshGenerator
    fipd_geom_file = ${raw '../../../../../../fipd-bison-integration-data/X423/ ${pin_id} / ${pin_id} _design.csv'}
    fipd_as_fabricated_file = ${raw '../../../../../../fipd-bison-integration-data/X423/ ${pin_id} / ${pin_id} _as_fabricated.csv'}

    gap_bottom_length = 0.31e-3 # arbitrary
    cladding_bottom_plug_length = 2.24e-3 # arbitrary
    cladding_top_plug_length = 2.24e-3 # arbitrary

    cladding_sidewall_radial_elements = 10
    cladding_sidewall_axial_element_numbers = '2 150 150'
    cladding_top_plug_radial_elements = 10
    cladding_top_plug_axial_elements = 5
    cladding_bottom_plug_radial_elements = 10
    cladding_bottom_plug_axial_elements = 5
    fuel_radial_elements = 6
    fuel_axial_element_intervals = '0 1'
    fuel_axial_element_numbers = '150'
    use_default_cladding_sidewall_axial_element_intervals = true
    elem_type = QUAD4
  []
[]

[Variables]
  [temp]
    initial_condition = 298
  []
[]

[AuxVariables]
  # Aux variables for output
  [creep_strain_mag]
    order = CONSTANT
    family = MONOMIAL
  []
  [gap_cond]
    order = CONSTANT
    family = MONOMIAL
  []
  [coolant_htc]
    order = CONSTANT
    family = MONOMIAL
  []
  [cumulative_damage_index]
    order = CONSTANT
    family = MONOMIAL
  []
  [element_failed]
    order = CONSTANT
    family = MONOMIAL
  []
  [solid_swell]
    block = fuel
    order = CONSTANT
    family = MONOMIAL
  []
  [gas_swell]
    block = fuel
    order = CONSTANT
    family = MONOMIAL
  []
  [volumetric_strain]
    block = fuel
    order = CONSTANT
    family = MONOMIAL
  []
  [hoop_stress]
    order = CONSTANT
    family = MONOMIAL
  []
  [hoop_creep_strain]
    order = CONSTANT
    family = MONOMIAL
  []
  [hoop_elastic_strain]
    order = CONSTANT
    family = MONOMIAL
  []
  [total_hoop_strain]
    order = CONSTANT
    family = MONOMIAL
  []
  [func_val1]
  []
  [func_val2]
  []
  # AuxVariables used for thermal expansion correction
  [fuel_thermal_strain_xx]
    order = CONSTANT
    family = MONOMIAL
    block = fuel
  []
  [fuel_thermal_strain_yy]
    order = CONSTANT
    family = MONOMIAL
    block = fuel
  []
  [fuel_thm_exp]
    order = CONSTANT
    family = MONOMIAL
    block = fuel
  []
  [clad_thermal_eigenstrain_xx]
    order = CONSTANT
    family = MONOMIAL
    block = cladding
  []
  [clad_thm_exp]
    order = CONSTANT
    family = MONOMIAL
    block = cladding
  []
[]

[Functions]  #copied from fipd-tdep
  [clad_od_temp]
    type =  FIPDAxialProfileFunction
    data_file = ${raw '../../../../../../fipd-bison-integration-data/X423/ ${pin_id} /clad_od_temp_history_ ${pin_id} .csv'}
    use_metadata = true
    mesh_generator = gen
  []
  [power_history]
    type = PiecewiseLinear
    data_file = ${raw '../../../../../../fipd-bison-integration-data/X423/ ${pin_id} /power_history_ ${pin_id} .csv'}
  []
  [pwr_axial_peaking_factors]
    type = FIPDAxialProfileFunction
    data_file = ${raw '../../../../../../fipd-bison-integration-data/X423/ ${pin_id} /peakingfactor_power_relative_ ${pin_id} .csv'}
    use_metadata = true
    mesh_generator = gen
    zero_ends = true
    data_shift_type = peaking
  []
  [pwr_axial_peaking_factors_elongate]
    type = FIPDAxialProfileFunction
    data_file = ${raw '../../../../../../fipd-bison-integration-data/X423/ ${pin_id} /peakingfactor_power_relative_ ${pin_id} .csv'}
    use_metadata = true
    mesh_generator = gen
    zero_ends = true
    data_shift_type = peaking
    fuel_elongation_pp = max_fuel_elongation
  []
  [fflux_axial_peaking_factors]
    type = FIPDAxialProfileFunction
    data_file = ${raw '../../../../../../fipd-bison-integration-data/X423/ ${pin_id} /peakingfactor_flux_relative_ ${pin_id} .csv'}
    use_metadata = true
    mesh_generator = gen
    zero_ends = true
    data_shift_type = peaking
    extrapolate_to_zero = true
  []
  [fflux_axial_peaking_factors_elongate]
    type = FIPDAxialProfileFunction
    data_file = ${raw '../../../../../../fipd-bison-integration-data/X423/ ${pin_id} /peakingfactor_flux_relative_ ${pin_id} .csv'}
    use_metadata = true
    mesh_generator = gen
    zero_ends = true
    data_shift_type = peaking
    extrapolate_to_zero = true
    fuel_elongation_pp = max_fuel_elongation
  []
  [flux_history]
    type = PiecewiseLinear
    data_file = ${raw '../../../../../../fipd-bison-integration-data/X423/ ${pin_id} /flux_history_ ${pin_id} .csv'}
  []
  [coolant_press_ramp]
    type = PiecewiseLinear
    x = '0 56465640'
    y = '0.151e6 0.151e6'
  []
  [id_vpp_func]
    type = MetallicFuelWastageDegradationFunction
    vectorpostprocessor_name = id_wastage
    argument_column = y
    wastage_type = ID
    value_column = wastage_thickness
    use_metadata = true
    degradation_factor = 0.001
    mesh_generator = 'gen'
    transition_width = 1E-4
  []
  [od_vpp_func]
    type = MetallicFuelWastageDegradationFunction
    vectorpostprocessor_name = od_wastage
    argument_column = y
    wastage_type = OD
    value_column = cc_wastage_thickness
    use_metadata = true
    degradation_factor = 0.001
    mesh_generator = 'gen'
    transition_width = 1E-4
  []
  [ci_temp]
    type = PiecewiseLinearFromVectorPostprocessor
    argument_column = y
    component = y
    value_column = temp
    vectorpostprocessor_name = clad_inn_temp
  []
  [na_vol]
    type = MeshPropertyFunction
    mesh_generator = gen
    mesh_property_name = sodium_volume
    scale_factor = -1.0
  []
[]

[Physics/SolidMechanics/QuasiStatic]
  temperature = temp
  [fuel]
    strain = FINITE
    add_variables = true
    generate_output = 'stress_xx stress_yy stress_zz vonmises_stress hydrostatic_stress creep_strain_xx creep_strain_yy creep_strain_zz elastic_strain_xx elastic_strain_yy elastic_strain_zz strain_xx strain_yy strain_zz'
    extra_vector_tags = 'ref'
    block = fuel
    eigenstrain_names = 'fuel_thermal_strain fuel_volumetric_strain'
  []
  [clad]
    strain = FINITE
    add_variables = true
    generate_output = 'stress_xx stress_yy stress_zz vonmises_stress hydrostatic_stress creep_strain_xx creep_strain_yy creep_strain_zz elastic_strain_xx elastic_strain_yy elastic_strain_zz strain_xx strain_yy strain_zz'
    extra_vector_tags = 'ref'
    block = cladding
    eigenstrain_names = 'clad_swelling clad_thermal_eigenstrain'
  []
[]

[Kernels]
  # Define kernels for the various terms in the PDE system
  [gravity]
    type = Gravity
    variable = disp_y
    value = -9.81
    extra_vector_tags = 'ref'
  []
  [heat]
    type = HeatConduction
    variable = temp
    extra_vector_tags = 'ref'
  []
  [heat_ie_f]
    type = HeatConductionTimeDerivative
    variable = temp
    extra_vector_tags = 'ref'
    block = fuel
    density_name = density
  []
  [heat_ie_c]
    type = HeatConductionTimeDerivative
    variable = temp
    extra_vector_tags = 'ref'
    block = cladding
    density_name = density
  []
  [heat_source]
    type = FissionRateHeatSource
    variable = temp
    fission_rate = fission_rate
    block = fuel
    extra_vector_tags = 'ref'
  []
[]

[AuxKernels]
  [conductance]
    type = MaterialRealAux
    property = gap_conductance
    variable = gap_cond
    boundary = fuel_outer_radial_surface
  []
  [cdf_amount]
    block = cladding
    type = MaterialRealAux
    property = cdf_failure
    variable = cumulative_damage_index
  []
  [failed_element]
    boundary = 'cladding_outside_right'
    type = MaterialRealAux
    property = failed
    variable = element_failed
  []
  [gas_swell]
    type = MaterialRealAux
    variable = gas_swell
    property = gas_swelling
    execute_on = timestep_end
  []
  [solid_swell]
    type = MaterialRealAux
    variable = solid_swell
    property = solid_swelling
    execute_on = timestep_end
  []
  [volumetric_strain]
    type = RankTwoScalarAux
    rank_two_tensor = total_strain
    variable = volumetric_strain
    scalar_type = VolumetricStrain
    execute_on = timestep_end
    block = fuel
  []
  [hoop_stress]
    type = RankTwoAux
    rank_two_tensor = stress
    variable = hoop_stress
    index_j = 2
    index_i = 2
    execute_on = timestep_end
  []
  [hoop_creep_strain]
    type = RankTwoAux
    rank_two_tensor = creep_strain
    variable = hoop_creep_strain
    index_j = 2
    index_i = 2
    execute_on = timestep_end
    block = cladding
  []
  [hoop_elastic_strain]
    type = RankTwoAux
    rank_two_tensor = elastic_strain
    variable = hoop_elastic_strain
    index_j = 2
    index_i = 2
    execute_on = timestep_end
    block = cladding
  []
  [total_hoop_strain]
    type = RankTwoAux
    rank_two_tensor = total_strain
    variable = total_hoop_strain
    index_j = 2
    index_i = 2
    execute_on = timestep_end
    block = cladding
  []
  [func_val1]
    type = FunctionAux
    function = id_vpp_func
    variable = func_val1
    block = cladding
  []
  [func_val2]
    type = FunctionAux
    function = od_vpp_func
    variable = func_val2
    block = cladding
  []
  # AuxKernels used to correct thermal expansion
  [fuel_thermal_strain_xx]
    type = RankTwoAux
    rank_two_tensor = fuel_thermal_strain
    variable = fuel_thermal_strain_xx
    index_j = 0
    index_i = 0
    execute_on = 'initial timestep_end'
    block = fuel
  []
  [fuel_thermal_strain_yy]
    type = RankTwoAux
    rank_two_tensor = fuel_thermal_strain
    variable = fuel_thermal_strain_yy
    index_j = 1
    index_i = 1
    execute_on = 'initial timestep_end'
    block = fuel
  []
  [fuel_thm_exp]
    type = SpatialUserObjectAux
    variable = fuel_thm_exp
    execute_on = 'initial timestep_end'
    user_object = fuel_thm_exp
    block = fuel
  []
  [clad_thermal_eigenstrain_xx]
    type = RankTwoAux
    rank_two_tensor = clad_thermal_eigenstrain
    variable = clad_thermal_eigenstrain_xx
    index_j = 0
    index_i = 0
    execute_on = 'initial timestep_end'
    block = cladding
  []
  [clad_thm_exp]
    type = SpatialUserObjectAux
    variable = clad_thm_exp
    execute_on = 'initial timestep_end'
    user_object = clad_thm_exp
    block = cladding
  []
[]

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

[ThermalContact]
  [thermal_contact]
    type = GapHeatTransfer
    variable = temp
    primary = cladding_inside_right
    secondary = fuel_outer_radial_surface
    quadrature = true
    gap_conductivity = 61.0
    min_gap = 0.5e-03
  []
[]

[BCs]
  [no_x_all]
    type = DirichletBC
    variable = disp_x
    boundary = centerline
    value = 0.0
  []
  [no_y_fuel]
    type = DirichletBC
    variable = disp_y
    boundary = fuel_bottom
    value = 0.0
  []
  [no_y_clad]
    type = DirichletBC
    variable = disp_y
    boundary = cladding_outside_bottom
    value = 0.0
  []
  [fuel_top_temp]
    type = FunctionDirichletBC
    boundary = fuel_top
    variable = temp
    function = ci_temp
  []
  [surf] #copied from fipd-tdep
    type = FunctionDirichletBC
    variable = temp
    boundary = 'cladding_outside_bottom cladding_outside_right cladding_outside_top'
    function = clad_od_temp
  []
  [Pressure]
    [coolantPressure]
      boundary = 'cladding_outside_bottom cladding_outside_right cladding_outside_top'
      function = coolant_press_ramp
    []
  []
  [PlenumPressure]
    [plenumPressure]
      boundary = 'fuel_outside_all cladding_inside_all'
      initial_pressure = 0.084e6 # Pa
      startup_time = 0
      R = 8.3143
      temperature = ave_temp_interior
      volume = gas_volume
      output = plenum_pressure
      material_input = fis_gas_released
    []
  []
[]

[Materials]
  [fission_rate]
    type = UPuZrFissionRate
    block = fuel
    rod_linear_power = power_history
    axial_power_profile = pwr_axial_peaking_factors
    use_metadata = true
    mesh_generator = gen
    outputs = all
  []
  [fission_rate_elongate]
    type = UPuZrFissionRate
    block = cladding
    fission_rate_name = fission_rate
    rod_linear_power = power_history
    axial_power_profile = pwr_axial_peaking_factors_elongate
    use_metadata = true
    mesh_generator = gen
    outputs = all
  []
  [burnup]
    type = UPuZrBurnup
    initial_X_Pu = ${fuel_pu}
    initial_X_Zr = ${fuel_zr}
    density = ${fuel_density}
    outputs = all
    block = fuel
  []
  [burnup_elongate]
    type = UPuZrBurnup
    initial_X_Pu = ${fuel_pu}
    initial_X_Zr = ${fuel_zr}
    density = ${fuel_density}
    outputs = all
    block = cladding
    burnup_name = burnup
  []
  [fast_neutron_flux]
    type = FastNeutronFlux
    calculate_fluence = true
    rod_ave_lin_pow = flux_history
    axial_power_profile = fflux_axial_peaking_factors
    block = fuel
    factor = 1.0
    outputs = all
  []
  [fast_neutron_flux_elongate]
    type = FastNeutronFlux
    calculate_fluence = true
    rod_ave_lin_pow = flux_history
    axial_power_profile = fflux_axial_peaking_factors_elongate
    block = cladding
    factor = 1.0
    outputs = all
  []
  [fuel_elasticity_tensor]
    type = UPuZrElasticityTensor
    block = fuel
    temperature = temp
  []
  [fuel_elastic_stress]
    type = ComputeMultipleInelasticStress
    tangent_operator = nonlinear
    inelastic_models = 'fuel_upuzrcreep'
    block = fuel
  []
  [fuel_upuzrcreep]
    type = UPuZrCreepUpdate
    block = fuel
    temperature = temp
    porosity = porosity
    max_inelastic_increment = 2e-3
  []
  [fuel_thermal_expansion]
    type = ComputeThermalExpansionEigenstrain
    block = fuel
    thermal_expansion_coeff = 1.18e-5
    temperature = temp
    stress_free_temperature = 295.0
    eigenstrain_name = fuel_thermal_strain
  []
  [fuel_volumetric_swelling]
    type = UPuZrVolumetricSwellingEigenstrainLM
    block = fuel
    use_preset_bubble_size = true
    anisotropic_factor = 1.26
    temperature = temp
    fission_rate = fission_rate
    burnup = burnup
    fis_gas_ret = fis_gas_ret
    hydrostatic_stress = hydrostatic_stress
    eigenstrain_name = fuel_volumetric_strain
    gas_swelling_scale_factor = 1.0
    outputs = all
  []
  [metal_fuel_thermal]
    type = UPuZrThermal
    block = fuel
    spheat_model = savage
    thcond_model = lanl
    porosity = porosity
    temperature = temp
  []
  [fuel_density]
    type = StrainAdjustedDensity
    displacements = 'disp_x disp_y'
    block = fuel
    strain_free_density = ${fuel_density}
  []
  [fission_gas_behavior]
    type = FgrUPuZrLM
    block = fuel
    temperature = temp
    fission_rate = fission_rate
    epsilon_c = 0.36
  []
  [clad_elasticity_tensor]
    type = SS316ElasticityTensor
    temperature = temp
    id_wastage_degradation_function = id_vpp_func
    od_wastage_degradation_function = od_vpp_func
    block = cladding
    elastic_constants_model = legacy_ifr
  []
  [clad_stress]
    type = ComputeMultipleInelasticStress
    tangent_operator = nonlinear
    inelastic_models = 'clad_ss316creep'
    block = cladding
  []
  [clad_ss316creep]
    type = SS316CreepUpdate
    block = cladding
    temperature = temp
    fast_neutron_flux = fast_neutron_flux
    id_wastage_degradation_function = id_vpp_func
    od_wastage_degradation_function = od_vpp_func
  []
  [clad_swelling]
    type = SS316VolumetricSwellingEigenstrain
    eigenstrain_name = clad_swelling
    fast_neutron_fluence = fast_neutron_fluence
    fast_neutron_flux = fast_neutron_flux
    temperature = temp
    outputs = all
  []
  [thermal_expansion]
    type = SS316ThermalExpansionEigenstrain
    block = cladding
    temperature = temp
    stress_free_temperature = 295.0
    eigenstrain_name = clad_thermal_eigenstrain
  []
  [clad_thermal]
    type = SS316Thermal
    block = cladding
    temperature = temp
  []
  [clad_density]
    type = StrainAdjustedDensity
    block = cladding
    strain_free_density = 7874.0
  []
  [longSS316_failure]
    type = D9FailureClad
    block = cladding
    method = steady_state
    temperature = temp
    outputs = all
    hoop_stress = stress_zz # Since 2D-RZ
  []
  [wastage_thickness]
    type = MetallicFuelWastage
    method = flux_ss316
    temperature = temp
    scale_factor = 1
    boundary = cladding_inside_right
    outputs = all
  []
  [cc_wastage_thickness]
    type = MetallicFuelCoolantWastage
    clad_material = SS316
    use_effective_method = true
    temperature = temp
    scale_factor = 1
    boundary = cladding_outside_right
    outputs = all
  []
[]

[Dampers]
  [limitT]
    type = MaxIncrement
    max_increment = 50
    variable = temp
  []
[]

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

[Executioner]
  type = Transient
  solve_type = 'PJFNK'
  petsc_options = '-snes_ksp_ew'
  petsc_options_iname = '-pc_type -pc_factor_mat_solver_package -ksp_gmres_restart'
  petsc_options_value = 'lu       superlu_dist                  51'
  line_search = 'none'

  l_max_its = 100
  l_tol = 1e-3
  nl_max_its = 50
  nl_rel_tol = 1e-5
  nl_abs_tol = 1e-7

  end_time = ${time_last}
  dtmin = 1
  dtmax = ${max_dt}

  [Quadrature]
    order = fifth
    side_order = seventh
  []
  [TimeStepper]
    type = IterationAdaptiveDT
    timestep_limiting_function = power_history
    #max_function_change = 300 # Removed to decrease run time
    timestep_limiting_postprocessor = creep_timestep
    dt = 100
    iteration_window = 2
    optimal_iterations = 10
    force_step_every_function_point = true
  []
[]

[Postprocessors]
  [_dt]
    type = TimestepSize
    outputs = 'csv_general console'
  []
  [num_lin_it]
    type = NumLinearIterations
    outputs = csv_general
  []
  [num_nonlin_it]
    type = NumNonlinearIterations
    outputs = csv_general
  []
  [tot_lin_it]
    type = CumulativeValuePostprocessor
    postprocessor = num_lin_it
    outputs = csv_general
  []
  [tot_nonlin_it]
    type = CumulativeValuePostprocessor
    postprocessor = num_nonlin_it
    outputs = csv_general
  []
  [alive_time]
    type = PerfGraphData
    section_name = Root
    data_type = TOTAL
    outputs = csv_general
  []
  [ave_temp_interior]
    type = SideAverageValue
    boundary = cladding_inside_top
    variable = temp
    execute_on = 'initial linear'
    outputs = csv_general
  []
  [approx_FCT]
    type = AverageNodalVariableValue
    boundary = centerline
    variable = temp
    outputs = csv_general
  []
  [max_approx_FCT]
    type = TimeExtremeValue
    value_type = max
    postprocessor = approx_FCT
    outputs = csv_general
  []
  [ave_FST]
    type = SideAverageValue
    boundary = fuel_outer_radial_surface
    variable = temp
    outputs = csv_general
  []
  [max_ave_FST]
    type = TimeExtremeValue
    value_type = max
    postprocessor = ave_FST
    outputs = csv_general
  []
  [ave_CIT]
    type = SideAverageValue
    boundary = cladding_inside_right
    variable = temp
    outputs = csv_general
  []
  [max_ave_CIT]
    type = TimeExtremeValue
    value_type = max
    postprocessor = ave_CIT
    outputs = csv_general
  []
  [avg_clad_temp]
    type = ElementAverageValue
    variable = temp
    block = cladding
    outputs = csv_general
  []
  [peak_clad_temp]
    type = ElementExtremeValue
    variable = temp
    value_type = max
    block = cladding
    outputs = csv_general
  []
  [peak_fuel_temp]
    type = ElementExtremeValue
    variable = temp
    value_type = max
    block = fuel
    outputs = csv_general
  []
  [max_hydro]
    type = ElementExtremeValue
    variable = hydrostatic_stress
    value_type = max
    block = fuel
    outputs = csv_general
  []
  [min_hydro]
    type = ElementExtremeValue
    variable = hydrostatic_stress
    value_type = min
    block = fuel
    outputs = csv_general
  []
  [avg_hydro]
    type = ElementAverageValue
    variable = hydrostatic_stress
    block = fuel
    outputs = csv_general
  []
  [peak_porosity]
    type = ElementExtremeValue
    variable = porosity
    value_type = max
    block = fuel
    outputs = csv_general
  []
  [clad_inner_vol]
    type = InternalVolume
    boundary = cladding_inside_all
    outputs = csv_general
  []
  [pellet_volume]
    type = InternalVolume
    boundary = fuel_outside_all
    outputs = csv_general
  []
  [gas_volume]
    type = InternalVolume
    boundary = 'fuel_outside_all cladding_inside_all'
    execute_on = 'initial timestep_end'
    addition = na_vol
    outputs = csv_general
  []
  [clad_fuel_gap]
    type = NodalExtremeValue
    variable = penetration
    boundary = fuel_outer_radial_surface
    outputs = csv_general
  []
  [max_cont_press]
    type = NodalExtremeValue
    variable = contact_pressure
    boundary = fuel_outer_radial_surface
    outputs = csv_general
  []
  [flux_from_clad]
    type = SideDiffusiveFluxIntegral
    variable = temp
    boundary = cladding_inside_right
    diffusivity = thermal_conductivity
    outputs = csv_general
  []
  [flux_from_fuel]
    type = SideDiffusiveFluxIntegral
    variable = temp
    boundary = fuel_outer_radial_surface
    diffusivity = thermal_conductivity
    outputs = csv_general
  []
  [rod_total_power]
    type = ElementIntegralPower
    variable = temp
    use_material_fission_rate = true
    fission_rate_material = fission_rate
    block = fuel
    outputs = csv_general
  []
  [LHGR_W_per_cm]
    type = FunctionValuePostprocessor
    function = power_history
    scale_factor = 0.01
    outputs = csv_general
  []
  [average_burnup]
    type = ElementAverageValue
    block = fuel
    variable = burnup
    outputs = csv_general
  []
  [max_cdf]
    type = ElementExtremeValue
    value_type = max
    variable = cumulative_damage_index
    outputs = csv_general
  []
  [fis_gas_produced]
    type = ElementIntegralMaterialProperty
    mat_prop = fis_gas_prod
    block = fuel
    outputs = csv_general
  []
  [fis_gas_released]
    type = ElementIntegralMaterialProperty
    mat_prop = fis_gas_rel
    block = fuel
    execute_on = 'initial timestep_end'
    outputs = csv_general
  []
  [creep_timestep]
    type = MaterialTimeStepPostprocessor
    block = fuel
    outputs = 'csv_general console'
  []
  [hydrostatic_stress]
    type = ElementAverageValue
    variable = hydrostatic_stress
    execute_on = 'initial timestep_end'
    block = fuel
    outputs = csv_general
  []
  [solid_swelling]
    type = ElementAverageValue
    variable = solid_swell
    block = fuel
    outputs = csv_general
  []
  [gas_swelling]
    type = ElementAverageValue
    variable = gas_swell
    block = fuel
    outputs = csv_general
  []
  [volumetric_strain]
    type = ElementAverageValue
    variable = volumetric_strain
    block = fuel
    outputs = csv_general
  []
  [fission_rate]
    type = ElementAverageValue
    variable = fission_rate
    block = fuel
    outputs = csv_general
  []
  [porosity]
    type = ElementAverageValue
    variable = porosity
    block = fuel
    outputs = csv_general
  []
  [gaseous_porosity]
    type = ElementAverageValue
    variable = gaseous_porosity
    block = fuel
    outputs = csv_general
  []
  [fis_gas_percent]
    type = FGRPercent
    fission_gas_released = fis_gas_released
    fission_gas_generated = fis_gas_produced
    outputs = csv_general
  []
  [max_clad_hoop_creep]
    type = ElementExtremeValue
    value_type = max
    block = cladding
    variable = hoop_creep_strain
    outputs = csv_general
  []
  [max_clad_creep_strain_mag]
    type = ElementExtremeValue
    value_type = max
    block = cladding
    variable = creep_strain_mag
    outputs = csv_general
  []
  [max_total_hoop_strain]
    type = ElementExtremeValue
    value_type = max
    block = cladding
    variable = total_hoop_strain
    outputs = csv_general
  []
  [max_fuel_elongation]
    type = NodalExtremeValue
    variable = disp_y
    boundary = fuel_top
    outputs = csv_general
  []
  [avg_fuel_ax_thm_str]
    type = AxisymmetricCenterlineAverageValue
    variable = fuel_thermal_strain_yy
    boundary = fuel_inner_radial_surface
    outputs = csv_general
  []
  [max_clad_elongation]
    type = NodalExtremeValue
    variable = disp_y
    boundary = 'cladding_outside_top cladding_outside_right'
    outputs = csv_general
  []
  [max_wastagethickness]
    type = ElementExtremeValue
    value_type = max
    variable = wastage_thickness
    outputs = 'console'
  []
  [avg_fuel_temp]
    type = ElementAverageValue
    variable = temp
    block = fuel
    execute_on = 'initial timestep_end'
    allow_duplicate_execution_on_initial = true
    outputs = csv_general
  []
[]

[VectorPostprocessors]
  [clad_x_disp]
    type = NodalValueSampler
    variable = disp_x
    boundary = cladding_outside_right
    sort_by = y
    outputs = none
  []
  [fuel_cl_temp]
    type = NodalValueSampler
    variable = temp
    boundary = centerline
    sort_by = y
    outputs = none
  []
  [fuel_surf_temp]
    type = NodalValueSampler
    variable = temp
    boundary = fuel_outer_radial_surface
    sort_by = y
    outputs = none
  []
  [clad_inn_temp]
    type = NodalValueSampler
    variable = temp
    boundary = cladding_inside_right
    sort_by = y
    execute_on = 'initial timestep_end'
    allow_duplicate_execution_on_initial = true
    outputs = none
  []
  [clad_out_temp]
    type = NodalValueSampler
    variable = temp
    boundary = cladding_outside_right
    sort_by = y
    outputs = none
  []
  [id_wastage]
    type = FuelRodLineValueSampler
    variable = wastage_thickness
    material = 'clad'
    fraction = 0.0
    num_points = 600
    orientation = 'vertical'
    fuel_pin_geometry = 'pin_geometry'
    execute_on = 'initial timestep_end'
    allow_duplicate_execution_on_initial = true
    outputs = none
  []
  [od_wastage]
    type = FuelRodLineValueSampler
    variable = cc_wastage_thickness
    material = 'clad'
    fraction = 1.0
    num_points = 600
    orientation = 'vertical'
    fuel_pin_geometry = 'pin_geometry'
    execute_on = 'initial timestep_end'
    allow_duplicate_execution_on_initial = true
    outputs = none
  []
  # PIE Comparison VPPs
  [nrad_comparison_0]
    type = FIPDAxialPIEComparison
    boundary = fuel_outer_radial_surface
    sort_by = y
    csv_file = ${raw '../../../../../../fipd-bison-integration-data/X423/ ${pin_id} /X423_ ${pin_id} _NRAD.csv'}
    variable = disp_x
    thermal_strain_variable = fuel_thm_exp
    involved_component = fuel
    mesh_generator = gen
    series_type_to_read = 'Fuel O.D. (mils)'
    outputs = csv_vpp_0
    enable = ${enable_0}
    execute_on = 'initial timestep_end'
    extra_pp_output_name = 'extra_csv extra_csv_0'
    max_fuel_elongation_pp_name = max_fuel_elongation
    avg_fuel_ax_thm_str_pp_name = avg_fuel_ax_thm_str
    gas_swelling_pp_name = gas_swelling
    solid_swelling_pp_name = solid_swelling
  []
  [nrad_comparison_a]
    type = FIPDAxialPIEComparison
    boundary = fuel_outer_radial_surface
    sort_by = y
    csv_file = ${raw '../../../../../../fipd-bison-integration-data/X423/ ${pin_id} /X423A_ ${pin_id} _NRAD.csv'}
    variable = disp_x
    thermal_strain_variable = fuel_thm_exp
    involved_component = fuel
    mesh_generator = gen
    series_type_to_read = 'Fuel O.D. (mils)'
    outputs = csv_vpp_a
    enable = ${enable_a}
    execute_on = 'initial timestep_end'
    extra_pp_output_name = 'extra_csv extra_csv_a'
    max_fuel_elongation_pp_name = max_fuel_elongation
    avg_fuel_ax_thm_str_pp_name = avg_fuel_ax_thm_str
    gas_swelling_pp_name = gas_swelling
    solid_swelling_pp_name = solid_swelling
  []
  [nrad_comparison_b]
    type = FIPDAxialPIEComparison
    boundary = fuel_outer_radial_surface
    sort_by = y
    csv_file = ${raw '../../../../../../fipd-bison-integration-data/X423/ ${pin_id} /X423B_ ${pin_id} _NRAD.csv'}
    variable = disp_x
    thermal_strain_variable = fuel_thm_exp
    involved_component = fuel
    mesh_generator = gen
    series_type_to_read = 'Fuel O.D. (mils)'
    outputs = csv_vpp_b
    enable = ${enable_b}
    execute_on = 'initial timestep_end'
    extra_pp_output_name = 'extra_csv extra_csv_b'
    max_fuel_elongation_pp_name = max_fuel_elongation
    avg_fuel_ax_thm_str_pp_name = avg_fuel_ax_thm_str
    gas_swelling_pp_name = gas_swelling
    solid_swelling_pp_name = solid_swelling
  []
  [nrad_comparison_c]
    type = FIPDAxialPIEComparison
    boundary = fuel_outer_radial_surface
    sort_by = y
    csv_file = ${raw '../../../../../../fipd-bison-integration-data/X423/ ${pin_id} /X423C_ ${pin_id} _NRAD.csv'}
    variable = disp_x
    thermal_strain_variable = fuel_thm_exp
    involved_component = fuel
    mesh_generator = gen
    series_type_to_read = 'Fuel O.D. (mils)'
    outputs = csv_vpp_c
    enable = ${enable_c}
    execute_on = 'initial timestep_end'
    extra_pp_output_name = 'extra_csv extra_csv_c'
    max_fuel_elongation_pp_name = max_fuel_elongation
    avg_fuel_ax_thm_str_pp_name = avg_fuel_ax_thm_str
    gas_swelling_pp_name = gas_swelling
    solid_swelling_pp_name = solid_swelling
  []
[]

[UserObjects]
  [pin_geometry]
    type = FuelPinGeometry
    clad_bottom = cladding_outside_bottom
    clad_inner_wall = cladding_inside_right
    clad_outer_wall = cladding_outside_right
    clad_top = cladding_outside_top
    pellet_exteriors = fuel_outside_all
  []
  [fuel_thm_exp]
    type = LayeredAverage
    variable = fuel_thermal_strain_xx
    direction = y
    num_layers = 1000
    block = fuel
  []
  [clad_thm_exp]
    type = LayeredAverage
    variable = clad_thermal_eigenstrain_xx
    direction = y
    num_layers = 1000
    block = cladding
  []
[]

[Outputs]
  perf_graph = true
  color = false
  [console]
    type = Console
    max_rows = 25
    time_step_interval =  1
    output_linear = true
    sync_times = ${time_spots}
  []
  [csv_vpp_0]
    type = CSV
    sync_only = true
    sync_times = ${time_spots_0}
    enable = ${enable_0}
    create_latest_symlink = true
  []
  [csv_vpp_a]
    type = CSV
    sync_only = true
    sync_times = ${time_spots_a}
    enable = ${enable_a}
    create_latest_symlink = true
  []
  [csv_vpp_b]
    type = CSV
    sync_only = true
    sync_times = ${time_spots_b}
    enable = ${enable_b}
    create_latest_symlink = true
  []
  [csv_vpp_c]
    type = CSV
    sync_only = true
    sync_times = ${time_spots_c}
    enable = ${enable_c}
    create_latest_symlink = true
  []
  [csv_general]
    type = CSV
    sync_times = ${time_spots}
  []
  [extra_csv]
    type = CSV
    sync_only = true
    sync_times = ${time_spots_extra}
  []
  [extra_csv_0]
    type = CSV
    sync_only = true
    sync_times = ${time_spots_0}
    enable = ${enable_0}
  []
  [extra_csv_a]
    type = CSV
    sync_only = true
    sync_times = ${time_spots_a}
    enable = ${enable_a}
  []
  [extra_csv_b]
    type = CSV
    sync_only = true
    sync_times = ${time_spots_b}
    enable = ${enable_b}
  []
  [extra_csv_c]
    type = CSV
    sync_only = true
    sync_times = ${time_spots_c}
    enable = ${enable_c}
  []
  [exodus]
    type = Exodus
    sync_only = true
    sync_times = ${time_spots}
  []
[]

(test/tests/solid_mechanics/SS316_creep/creep_SS316_rz_test.i)

#
# This is an axisymmetric version of the 3D test. Since the boundary
# conditions are similar (i.e., uniaxial applied stress of 10 MPa) the
# analytical solution is the same.
#
#
# e_tot = e_elas +
#         e_thermal_creep +
#         e_irradiation_creep
#
#       = P/E +
#         (a * exp(-Q_ln/(RT)) * sinh(b * exp(-Q_pw/(RT)) * sigma) +
#         [B + DS] * fast_neutron_flux * sigma) * dt
#
#         where P = pressure load
#               E = Youngs modulus
#               a,b = thermal creep coefficients
#               sigma = stress in MPa
#               Q_ln-Q_pw = activation energies
#               B = Irradiation creep coefficient "creep compliance"
#               D = Irradiation creep coefficent "coupling coefficient"
#               S = Swelling Rate
#               dt = total time
#               fast_neutron_flux = 3e17 n/m**2/s (applied) but model expects units of 10**22 n/cm**2/s
#               so, fast_neutron_flux = 3e17 * 0.0001 (convert from cm**2 to m**2) * 1e-22 (convert to 10*22)= 3e-9
#
#               Note that these strains are expressed as %, so they have to be divided by 100, thus  the term 0.01 in the strain calculations found below.
#
#
#
#
#   For this test, the analytical solution at t=1e5 is:
#
#               e_tot = 10e6/1.90e11 +
#                       ((21.04889 * exp(-170000/(8.314*973.15)) *
#                       sinh(0.07 * (10) * exp(-12000/(8.314*973.15)))) +
#                       [1e-6 + 0.6e-2 * 0.01] * 3e-9 * 10 * 0.01) * 1e5
#
#                     = 5.263e-5 + (((4.336e-12 * 1.012) + 1.83e-12 * 0.01) * 1e5)
#                     = 5.263e-5 + 4.386e-7 + 1.83e-9 = 5.307e-5
#   e_elas = 5.263e-5
#   e_thermal_creep = 4.386e-7
#   e_irradiation_creep = 1.83e-9
#   e_tot = 5.307e-5
#
#
#


[Mesh]
  coord_type = RZ
  use_displaced_mesh = false
  [mesh]
    type = GeneratedMeshGenerator
    dim = 2
  []
[]

[GlobalParams]
  displacements = 'disp_x disp_y'
  order = FIRST
  family = LAGRANGE
  temperature = temp
[]

[Variables]
  [temp]
    initial_condition = 973.15
  []
[]

[AuxVariables]
  [stress_yy]
    order = CONSTANT
    family = MONOMIAL
  []
  [elastic_strain_yy]
    order = CONSTANT
    family = MONOMIAL
  []
  [creep_strain_yy]
    order = CONSTANT
    family = MONOMIAL
  []
[]

[Physics]
  [SolidMechanics]
    [QuasiStatic]
      [all]
        strain = SMALL
        incremental = true
        add_variables = true
        generate_output = 'stress_yy elastic_strain_yy creep_strain_yy'
      []
    []
  []
[]

[Functions]
  [top_pull]
    type = PiecewiseLinear
    x = '0 1'
    y = '1 1'
  []
[]

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

[AuxKernels]
  [stress_yy]
    type = RankTwoAux
    rank_two_tensor = stress
    variable = stress_yy
    index_i = 1
    index_j = 1
  []
  [elastic_strain_yy]
    type = RankTwoAux
    rank_two_tensor = elastic_strain
    variable = elastic_strain_yy
    index_i = 1
    index_j = 1
  []
  [creep_strain_yy]
    type = RankTwoAux
    rank_two_tensor = creep_strain
    variable = creep_strain_yy
    index_i = 1
    index_j = 1
    execute_on = timestep_end
  []
[]

[BCs]
  [u_top_pull]
    type = Pressure
    variable = 'disp_y'
    boundary = top
    factor = -10.0e6
    function = top_pull
  []
  [u_bottom_fix]
    type = DirichletBC
    variable = disp_y
    boundary = bottom
    value = 0.0
  []
  [temp_fix]
    type = DirichletBC
    variable = temp
    boundary = top
    value = 973.15
  []
[]

[Materials]
  [fast_neutron_flux]
    type = FastNeutronFlux
    calculate_fluence = true
    factor = 3e17
  []
  [elasticity_tensor]
    type = ComputeIsotropicElasticityTensor
    youngs_modulus = 1.90e11
    poissons_ratio = 0.265
  []
  [creep]
    type = SS316CreepUpdate
    fast_neutron_flux = fast_neutron_flux
  []
  [radial_return_stress]
    type = ComputeMultipleInelasticStress
    tangent_operator = elastic
    inelastic_models = 'creep'
  []
  [thermal]
    type = HeatConductionMaterial
    thermal_conductivity = 100.0
  []
  [density]
    type = StrainAdjustedDensity
    strain_free_density = 1.0
  []
[]

[Executioner]
  type = Transient

  solve_type = 'PJFNK'

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

  line_search = 'none'

  l_max_its  = 50
  nl_max_its = 20
  l_tol = 1e-5
  start_time = 0.0
  num_steps = 10
  dt = 10000
[]

[Outputs]
  exodus = true
[]

(assessment/MOX/JOYO/B14/PTM002/analysis/b14_ptm002_2DRZ_t.i)


initial_fuel_density = 11057.75

[GlobalParams]
  density = ${initial_fuel_density}
  order = SECOND
  family = LAGRANGE
  energy_per_fission = 3.2e-11 # J/fission
  volumetric_locking_correction = false
  displacements = 'disp_x disp_y'
[]
[Problem]
  type = ReferenceResidualProblem
  extra_tag_vectors = 'ref'
  reference_vector = 'ref'
[]
[Mesh]
  coord_type = RZ
  [smeared_pellet_mesh]
    type = FuelPinMeshGenerator
    pellet_quantity = 1
    pellet_height = 0.4
    pellet_outer_radius = 0.002675
    pellet_mesh_density = customize
    clad_mesh_density = customize
    clad_gap_width = 0.000105
    clad_thickness = 0.00047
    clad_bot_gap_height = 1.0e-3
    bottom_clad_height = 2.24e-3
    top_clad_height = 2.24e-3
    clad_top_gap_height = 0.685
    elem_type = QUAD8
    nx_c = 4
    ny_c = 100
    nx_p = 20
    ny_p = 100
    ny_cu = 3
    ny_cl = 3
  []
  patch_size = 50
  patch_update_strategy = iteration
  partitioner = centroid
  centroid_partitioner_direction = y
[]
[UserObjects]
  [pin_geometry]
    type = FuelPinGeometry
  []
[]
[Variables]
  [temp]
    initial_condition = 295.0
  []
[]
[AuxVariables]
  [pore]
  []
  [fission_rate]
  []
  [burnup]
  []
  [radial_strain]
    order = CONSTANT
    family = MONOMIAL
  []
  [effective_creep_strain]
    block = clad
    order = CONSTANT
    family = MONOMIAL
  []
[]
[Functions]
  [power_history] #related to the LHGR at the midplane
    type = PiecewiseLinear
    x = ' 0  72000  158040  160200  246600  248400  249000.012  251280'
    y = ' 0  39814.5  39814.5  44289.3  44289.3  53927.4  53927.4  0'
  []
  [fast_neutron_flux_function]
    type = PiecewiseLinear
    x = '0 251280'
    y = '3.3e+15 3.3e+15'
  []
  [f_temp_out_clad]
    type = PiecewiseBilinear
    x = '0.000175  0.0464075  0.0843675  0.1075625  0.152025  0.1994625  0.2464725  0.2947475  0.356915  0.43356  0.49848  0.625  0.700475  0.797485  0.8723425  0.96'
    y = '0 251280'
    z = '295  295  295  295  295  295  295  295  295  295  295  295  295  295  295  295 634.94  662.273  676.998  686.217  706.339  727  743.358  758.311  780.069  799.077  815.576  846.374  860.233  875.494  882.809  889.8'
    scale_factor = 1
    axis = 1
  []
  [axial_peaking_factors]
    type = PiecewiseBilinear
    x = '0.000175  0.0464075  0.0843675  0.1075625  0.152025  0.1994625  0.2464725  0.2947475  0.356915  0.43356  0.49848  0.625  0.700475  0.797485  0.8723425  0.96'
    y = '0 251280'
    z = '0.751  0.752  0.767  0.796  0.82  0.852  0.875  0.915  0.944  0.963  0.988  1  0.985  0.955  0.913  0.846 0.751  0.752  0.767  0.796  0.82  0.852  0.875  0.915  0.944  0.963  0.988  1  0.985  0.955  0.913  0.846'
    scale_factor = 1
    axis = 1
  []
  [q]
    type = CompositeFunction
    functions = 'power_history axial_peaking_factors'
  []
  [average_power_history]
    type = PiecewiseLinear
    x = ' 0  72000  158040  160200  246600  248400  249000.012  251280'
    y = ' 0  34700  34700  38600  38600  47000  47000  0'
  []
[]
[Physics/SolidMechanics/QuasiStatic]
  [pellets]
    block = pellet
    add_variables = true
    strain = FINITE
    eigenstrain_names = 'fuel_thermal_strain fuel_volumetric_strain'
    generate_output = 'vonmises_stress stress_xx stress_yy stress_zz'
    use_finite_deform_jacobian = true
    extra_vector_tags = 'ref'
  []
  [clad]
    block = clad
    add_variables = true
    strain = FINITE
    eigenstrain_names = 'clad_thermal_eigenstrain'
    generate_output = 'vonmises_stress stress_xx stress_yy stress_zz'
    use_finite_deform_jacobian = true
    extra_vector_tags = 'ref'
  []
[]
[Kernels]
  [gravity]
    type = Gravity
    variable = disp_y
    value = -9.81
    extra_vector_tags = 'ref'
  []
  [heat]
    type = HeatConduction
    variable = temp
    extra_vector_tags = 'ref'
  []
  [heat_ie]
    type = HeatConductionTimeDerivative
    variable = temp
    extra_vector_tags = 'ref'
  []
  [heat_source]
    type = NeutronHeatSource
    variable = temp
    fission_rate = fission_rate
    extra_vector_tags = 'ref'
  []
[]
[AuxKernels]
  [fission_rate]
    type = FissionRateGeneral
    fission_rate_formulation = MOX
    variable = fission_rate
    block = pellet
    initial_porosity = 0.1372
    axial_power_profile = axial_peaking_factors
    rod_ave_lin_pow = power_history
    pellet_diameter = 0.00535
    execute_on = timestep_begin
    porosity = pore
    energy_per_fission = 3.2e-11
  []
  [burnup]
    type = BurnupAux
    block = pellet
    fission_rate = fission_rate
    variable = burnup
    execute_on = timestep_begin
  []
  [radial_strain]
    type = RankTwoAux
    rank_two_tensor = total_strain
    variable = radial_strain
    index_i = 0
    index_j = 0
    execute_on = timestep_end
  []
  [effective_creep_strain]
    type = MaterialRealAux
    property = effective_creep_strain
    variable = effective_creep_strain
    block = clad
    execute_on = timestep_end
  []
[]
[Contact]
  [pellet_clad_mechanical]
    primary = 5
    secondary = 10
    formulation = kinematic
    model = frictionless
    penalty = 1e7
  []
[]
[ThermalContact]
  [thermal_contact]
    type = GasGapHeatTransfer
    variable = temp
    primary = 5
    secondary = 10
    initial_moles = initial_moles
    gas_released = fis_gas_released
    contact_pressure = contact_pressure
    quadrature = true
  []
[]
[BCs]
  [temp_clad_outside]
    type = FunctionDirichletBC
    variable = temp
    function = f_temp_out_clad
    boundary = 2
  []
  [no_x_all]
    type = DirichletBC
    variable = disp_x
    boundary = 12
    value = 0.0
  []
  [no_y_clad_bottom]
    type = DirichletBC
    variable = disp_y
    boundary = 1
    value = 0.0
  []
  [no_y_fuel_bottom]
    type = DirichletBC
    variable = disp_y
    boundary = 20
    value = 0.0
  []
  [PlenumPressure]
    [plenumPressure]
      boundary = 9
      initial_pressure = 101325
      startup_time = 0
      R = 8.3143
      output_initial_moles = initial_moles
      temperature = ave_temp_interior
      volume = gas_volume
      material_input = fis_gas_released
      output = plenum_pressure
    []
  []
[]
[Materials]
  [fast_neutron_flux]
    type = FastNeutronFlux
    calculate_fluence = true
    block = clad
    flux_function = fast_neutron_flux_function
  []
  [fuel_thermal]
    type = MAMOXThermal
    temperature = temp
    porosity = pore
    block = pellet
    Am_content = 0.0237
    oxy_to_metal_ratio = 1.982
    output_properties = 'thermal_conductivity'
  []
  [fuel_density]
    type = StrainAdjustedDensity
    block = pellet
    strain_free_density = ${initial_fuel_density}
  []
  [fuel_elasticity_tensor]
    type = MAMOXElasticityTensor
    block = pellet
  []
  [elastic_stress]
    type = ComputeFiniteStrainElasticStress
    block = pellet
  []
  [fuel_thermal_expansion]
    type = MAMOXThermalExpansionEigenstrain
    block = pellet
    temperature = temp
    stress_free_temperature = 295.0
    oxygen_to_metal_ratio = 1.98
    eigenstrain_name = fuel_thermal_strain
  []
  [fuel_volumetric_swelling]
    type = UO2VolumetricSwellingEigenstrain
    gas_swelling_model_type = SIFGRS
    block = pellet
    temperature = temp
    burnup = burnup
    initial_fuel_density = 11057.75
    eigenstrain_name = fuel_volumetric_strain
  []
  [clad_thermal]
    type = SS316Thermal
    block = clad
    temperature = temp
  []
  [clad_density]
    type = StrainAdjustedDensity
    block = clad
    strain_free_density = 8000
  []
  [clad_elasticity_tensor]
    type = SS316ElasticityTensor
    block = clad
    temperature = temp
    elastic_constants_model = legacy_ifr
  []
  [thermal_expansion]
    type = SS316ThermalExpansionEigenstrain
    block = clad
    temperature = temp
    stress_free_temperature = 295.0
    eigenstrain_name = clad_thermal_eigenstrain
  []
  [clad_ss316creep]
    type = SS316CreepUpdate
    block = clad
    temperature = temp
    fast_neutron_flux = fast_neutron_flux
  []
  [clad_stress]
    type = ComputeMultipleInelasticStress
    tangent_operator = elastic
    inelastic_models = 'clad_ss316creep'
    block = clad
  []
  [fission_gas_release]
    type = UO2Sifgrs
    block = pellet
    temperature = temp
    burnup = burnup
    fission_rate = fission_rate
    grain_radius_const = 12e-06
    bubble_gb_limit = 1.0e+11
  []
[]
[Preconditioning]
  [SMP]
    type = SMP
    full = true
  []
[]
[Executioner]
  type = Transient
  solve_type = 'PJFNK'
  petsc_options = '-snes_ksp_ew'
  petsc_options_iname = '-pc_type -sub_pc_type'
  petsc_options_value = 'asm lu'
  line_search = 'none'
  fixed_point_max_its = 1
  fixed_point_abs_tol = 1e-3
  fixed_point_rel_tol = 1e-3
  l_max_its = 50
  l_tol = 8e-3
  nl_max_its = 50
  nl_rel_tol = 1e-3
  nl_abs_tol = 1e-3
  start_time = 0
  n_startup_steps = 1
  end_time = 251280
  dtmax = 10000
  dtmin = 0.25
  [TimeStepper]
    type = IterationAdaptiveDT
    dt = 1e1
    optimal_iterations = 15
    iteration_window = 2
    linear_iteration_ratio = 100
    growth_factor = 2
    cutback_factor = .5
    force_step_every_function_point = true
    timestep_limiting_function = power_history
  []
[]
[Postprocessors]
  [ave_temp_interior]
    type = ElementAverageValue
    variable = temp
    execute_on = 'initial linear'
  []
  [average_burnup]
    type = ElementAverageValue
    variable = burnup
  []
  [ave_pore]
    type = ElementAverageValue
    block = pellet
    variable = pore
  []
  [max_pore]
    type = NodalExtremeValue
    block = pellet
    value_type = max
    variable = pore
  []
  [min_pore]
    type = NodalExtremeValue
    block = pellet
    value_type = min
    variable = pore
  []
  [fis_gas_produced]
    type = ElementIntegralFisGasGeneratedSifgrs
    block = pellet
    execute_on = 'linear'
  []
  [fis_gas_released]
    type = ElementIntegralFisGasReleasedSifgrs
    block = pellet
    execute_on = 'linear'
  []
  [fis_gas_released_percentage]
    type = FGRPercent
    fission_gas_generated = fis_gas_produced
    fission_gas_released = fis_gas_released
    execute_on = 'linear'
  []
  [gas_volume]
    type = InternalVolume
    boundary = 9
    execute_on = 'initial linear'
  []
  [rod_total_power]
    type = ElementIntegralPower
    variable = temp
    fission_rate = fission_rate
    block = pellet
  []
  [rod_input_power]
    type = FunctionValuePostprocessor
    function = average_power_history
    scale_factor = 0.4 # rod height
  []
[]
[VectorPostprocessors]
  [fuel_radial_temperature_Sample1]
    type = LineValueSampler
    variable = temp
    start_point = '0.0 0.283 0.0'
    end_point = '0.002675 0.283 0.0'
    num_points = 200
    execute_on = final
    sort_by = x
    outputs = line_plot
  []
  [radial_porosity_Sample1]
    type = LineValueSampler
    variable = pore
    start_point = '0.0 0.283 0.0'
    end_point = '0.002675 0.283 0.0'
    num_points = 200
    execute_on = final
    sort_by = x
    outputs = line_plot
  []
  [fuel_radial_temperature_Sample2]
    type = LineValueSampler
    variable = temp
    start_point = '0.0 0.347 0.0'
    end_point = '0.002675 0.347 0.0'
    num_points = 200
    execute_on = final
    sort_by = x
    outputs = line_plot
  []
  [radial_porosity_Sample2]
    type = LineValueSampler
    variable = pore
    start_point = '0.0 0.347 0.0'
    end_point = '0.002675 0.347 0.0'
    num_points = 200
    execute_on = final
    sort_by = x
    outputs = line_plot
  []
  [fuel_radial_temperature_Sample3]
    type = LineValueSampler
    variable = temp
    start_point = '0.0 0.2 0.0'
    end_point = '0.002675 0.2 0.0'
    num_points = 200
    execute_on = final
    sort_by = x
    outputs = line_plot
  []
  [radial_porosity_Sample3]
    type = LineValueSampler
    variable = pore
    start_point = '0.0 0.2 0.0'
    end_point = '0.002675 0.2 0.0'
    num_points = 200
    execute_on = final
    sort_by = x
    outputs = line_plot
  []
[]

[PerformanceMetricOutputs]
[]

[Outputs]
  perf_graph = true
  exodus = true
  color = true
  csv = true
  [console]
    type = Console
    max_rows = 25
  []
  [line_plot]
    type = CSV
    execute_on = 'FINAL'
    time_step_interval =  1
    file_base = 1d
    create_final_symlink = true
  []
  [chkfile]
    type = CSV
    execute_on = FINAL
    show = 'ave_temp_interior max_pore'
  []
[]
[MultiApps]
  [sub]
    type = TransientMultiApp
    app_type = BisonApp
    execute_on = timestep_end
    catch_up = true
    max_catch_up_steps = 10
    positions_file = positions.txt
    input_files = b14_ptm002_pore.i
  []
[]
[Transfers]
  [temp_to_sub]
    type = MultiAppProjectionTransfer
    to_multi_app = sub
    source_variable = temp
    variable = temp
  []
  [pore_from_sub]
    type = MultiAppGeometricInterpolationTransfer
    from_multi_app = sub
    source_variable = pore
    variable = pore
  []
[]
[Debug]
  show_var_residual_norms = true
[]

(assessment/MOX/JOYO/B14/PTM010/analysis/b14_ptm010_2DRZ_t.i)


initial_fuel_density = 10964.6

[GlobalParams]
  density = ${initial_fuel_density}
  order = SECOND
  family = LAGRANGE
  energy_per_fission = 3.2e-11 # J/fission
  volumetric_locking_correction = false
  displacements = 'disp_x disp_y'
[]
[Problem]
  type = ReferenceResidualProblem
  extra_tag_vectors = 'ref'
  reference_vector = 'ref'
[]
[Mesh]
  coord_type = RZ
  [smeared_pellet_mesh]
    type = FuelPinMeshGenerator
    pellet_quantity = 1
    pellet_height = 0.4
    pellet_outer_radius = 0.002675
    pellet_mesh_density = customize
    clad_mesh_density = customize
    clad_gap_width = 0.000105
    clad_thickness = 0.00047
    clad_bot_gap_height = 1.0e-3
    bottom_clad_height = 2.24e-3
    top_clad_height = 2.24e-3
    clad_top_gap_height = 0.685
    elem_type = QUAD8
    nx_c = 4
    ny_c = 100
    nx_p = 20
    ny_p = 100
    ny_cu = 3
    ny_cl = 3
  []
  patch_size = 50
  patch_update_strategy = iteration
  partitioner = centroid
  centroid_partitioner_direction = y
[]
[UserObjects]
  [pin_geometry]
    type = FuelPinGeometry
  []
[]
[Variables]
  [temp]
    initial_condition = 295.0
  []
[]
[AuxVariables]
  [pore]
  []
  [fission_rate]
  []
  [burnup]
  []
  [radial_strain]
    order = CONSTANT
    family = MONOMIAL
  []
  [effective_creep_strain]
    block = clad
    order = CONSTANT
    family = MONOMIAL
  []
[]
[Functions]
  [power_history] #related to the LHGR at the midplane
    type = PiecewiseLinear
    x = ' 0  72000  158040  160200  246600  248400  249000.012  251280'
    y = ' 0  39814.5  39814.5  44289.3  44289.3  53927.4  53927.4  0'
  []
  [fast_neutron_flux_function]
    type = PiecewiseLinear
    x = '0 251280'
    y = '3.3e+15 3.3e+15'
  []
  [f_temp_out_clad]
    type = PiecewiseBilinear
    x = '0.000175  0.0464075  0.0843675  0.1075625  0.152025  0.1994625  0.2464725  0.2947475  0.356915  0.43356  0.49848  0.625  0.700475  0.797485  0.8723425  0.96'
    y = '0 251280'
    z = '295  295  295  295  295  295  295  295  295  295  295  295  295  295  295  295 634.94  662.273  676.998  686.217  706.339  727  743.358  758.311  780.069  799.077  815.576  846.374  860.233  875.494  882.809  889.8'
    scale_factor = 1
    axis = 1
  []
  [axial_peaking_factors]
    type = PiecewiseBilinear
    x = '0.000175  0.0464075  0.0843675  0.1075625  0.152025  0.1994625  0.2464725  0.2947475  0.356915  0.43356  0.49848  0.625  0.700475  0.797485  0.8723425  0.96'
    y = '0 251280'
    z = '0.751  0.752  0.767  0.796  0.82  0.852  0.875  0.915  0.944  0.963  0.988  1  0.985  0.955  0.913  0.846 0.751  0.752  0.767  0.796  0.82  0.852  0.875  0.915  0.944  0.963  0.988  1  0.985  0.955  0.913  0.846'
    scale_factor = 1
    axis = 1
  []
  [q]
    type = CompositeFunction
    functions = 'power_history axial_peaking_factors'
  []
  [average_power_history]
    type = PiecewiseLinear
    x = ' 0  72000  158040  160200  246600  248400  249000.012  251280'
    y = ' 0  34700  34700  38600  38600  47000  47000  0'
  []
[]
[Physics/SolidMechanics/QuasiStatic]
  [pellets]
    block = pellet
    add_variables = true
    strain = FINITE
    eigenstrain_names = 'fuel_thermal_strain fuel_volumetric_strain'
    generate_output = 'vonmises_stress stress_xx stress_yy stress_zz'
    use_finite_deform_jacobian = true
    extra_vector_tags = 'ref'
  []
  [clad]
    block = clad
    add_variables = true
    strain = FINITE
    eigenstrain_names = 'clad_thermal_eigenstrain'
    generate_output = 'vonmises_stress stress_xx stress_yy stress_zz'
    use_finite_deform_jacobian = true
    extra_vector_tags = 'ref'
  []
[]
[Kernels]
  [gravity]
    type = Gravity
    variable = disp_y
    value = -9.81
    extra_vector_tags = 'ref'
  []
  [heat]
    type = HeatConduction
    variable = temp
    extra_vector_tags = 'ref'
  []
  [heat_ie]
    type = HeatConductionTimeDerivative
    variable = temp
    extra_vector_tags = 'ref'
  []
  [heat_source]
    type = NeutronHeatSource
    variable = temp
    fission_rate = fission_rate
    extra_vector_tags = 'ref'
  []
[]
[AuxKernels]
  [fission_rate]
    type = FissionRateGeneral
    fission_rate_formulation = MOX
    variable = fission_rate
    block = pellet
    initial_porosity = 0.143
    axial_power_profile = axial_peaking_factors
    rod_ave_lin_pow = power_history
    pellet_diameter = 0.00535
    execute_on = timestep_begin
    porosity = pore
    energy_per_fission = 3.2e-11
  []
  [burnup]
    type = BurnupAux
    block = pellet
    fission_rate = fission_rate
    variable = burnup
    execute_on = timestep_begin
  []
  [radial_strain]
    type = RankTwoAux
    rank_two_tensor = total_strain
    variable = radial_strain
    index_i = 0
    index_j = 0
    execute_on = timestep_end
  []
  [effective_creep_strain]
    type = MaterialRealAux
    property = effective_creep_strain
    variable = effective_creep_strain
    block = clad
    execute_on = timestep_end
  []
[]
[Contact]
  [pellet_clad_mechanical]
    primary = 5
    secondary = 10
    formulation = kinematic
    model = frictionless
    penalty = 1e7
  []
[]
[ThermalContact]
  [thermal_contact]
    type = GasGapHeatTransfer
    variable = temp
    primary = 5
    secondary = 10
    initial_moles = initial_moles
    gas_released = fis_gas_released
    contact_pressure = contact_pressure
    quadrature = true
  []
[]
[BCs]
  [temp_clad_outside]
    type = FunctionDirichletBC
    variable = temp
    function = f_temp_out_clad
    boundary = 2
  []
  [no_x_all]
    type = DirichletBC
    variable = disp_x
    boundary = 12
    value = 0.0
  []
  [no_y_clad_bottom]
    type = DirichletBC
    variable = disp_y
    boundary = 1
    value = 0.0
  []
  [no_y_fuel_bottom]
    type = DirichletBC
    variable = disp_y
    boundary = 20
    value = 0.0
  []
  [PlenumPressure]
    [plenumPressure]
      boundary = 9
      initial_pressure = 101325
      startup_time = 0
      R = 8.3143
      output_initial_moles = initial_moles
      temperature = ave_temp_interior
      volume = gas_volume
      material_input = fis_gas_released
      output = plenum_pressure
    []
  []
[]
[Materials]
  [fast_neutron_flux]
    type = FastNeutronFlux
    calculate_fluence = true
    block = clad
    flux_function = fast_neutron_flux_function
  []
  [fuel_thermal]
    type = MAMOXThermal
    temperature = temp
    porosity = pore
    block = pellet
    Am_content = 0.0237
    oxy_to_metal_ratio = 1.99
    output_properties = 'thermal_conductivity'
  []
  [fuel_density]
    type = StrainAdjustedDensity
    block = pellet
    strain_free_density = ${initial_fuel_density}
  []
  [fuel_elasticity_tensor]
    type = MAMOXElasticityTensor
    block = pellet
  []
  [elastic_stress]
    type = ComputeFiniteStrainElasticStress
    block = pellet
  []
  [fuel_thermal_expansion]
    type = MAMOXThermalExpansionEigenstrain
    block = pellet
    temperature = temp
    stress_free_temperature = 295.0
    oxygen_to_metal_ratio = 2.0
    eigenstrain_name = fuel_thermal_strain
  []
  [fuel_volumetric_swelling]
    type = UO2VolumetricSwellingEigenstrain
    gas_swelling_model_type = SIFGRS
    block = pellet
    temperature = temp
    burnup = burnup
    initial_fuel_density = 10964.6
    eigenstrain_name = fuel_volumetric_strain
  []
  [clad_thermal]
    type = SS316Thermal
    block = clad
    temperature = temp
  []
  [clad_density]
    type = StrainAdjustedDensity
    block = clad
    strain_free_density = 8000
  []
  [clad_elasticity_tensor]
    type = SS316ElasticityTensor
    block = clad
    temperature = temp
    elastic_constants_model = legacy_ifr
  []
  [thermal_expansion]
    type = SS316ThermalExpansionEigenstrain
    block = clad
    temperature = temp
    stress_free_temperature = 295.0
    eigenstrain_name = clad_thermal_eigenstrain
  []
  [clad_ss316creep]
    type = SS316CreepUpdate
    block = clad
    temperature = temp
    fast_neutron_flux = fast_neutron_flux
  []
  [clad_stress]
    type = ComputeMultipleInelasticStress
    tangent_operator = elastic
    inelastic_models = 'clad_ss316creep'
    block = clad
  []
  [fission_gas_release]
    type = UO2Sifgrs
    block = pellet
    temperature = temp
    burnup = burnup
    fission_rate = fission_rate
    grain_radius_const = 10e-06
    bubble_gb_limit = 1.0e+11
  []
[]
[Preconditioning]
  [SMP]
    type = SMP
    full = true
  []
[]
[Executioner]
  type = Transient
  solve_type = 'PJFNK'
  petsc_options = '-snes_ksp_ew'
  petsc_options_iname = '-pc_type -sub_pc_type'
  petsc_options_value = 'asm lu'
  line_search = 'none'
  fixed_point_max_its = 1
  fixed_point_abs_tol = 1e-3
  fixed_point_rel_tol = 1e-3
  l_max_its = 50
  l_tol = 8e-3
  nl_max_its = 50
  nl_rel_tol = 1e-3
  nl_abs_tol = 1e-3
  start_time = 0
  n_startup_steps = 1
  end_time = 251280
  dtmax = 10000
  dtmin = 0.25
  automatic_scaling = true
  compute_scaling_once = false
  [TimeStepper]
    type = IterationAdaptiveDT
    dt = 1e1
    optimal_iterations = 15
    iteration_window = 2
    linear_iteration_ratio = 100
    growth_factor = 2
    cutback_factor = .5
    force_step_every_function_point = true
    timestep_limiting_function = power_history
  []
[]
[Postprocessors]
  [ave_temp_interior]
    type = ElementAverageValue
    variable = temp
    execute_on = 'initial linear'
  []
  [average_burnup]
    type = ElementAverageValue
    variable = burnup
  []
  [ave_pore]
    type = ElementAverageValue
    block = pellet
    variable = pore
  []
  [max_pore]
    type = NodalExtremeValue
    block = pellet
    value_type = max
    variable = pore
  []
  [min_pore]
    type = NodalExtremeValue
    block = pellet
    value_type = min
    variable = pore
  []
  [fis_gas_produced]
    type = ElementIntegralFisGasGeneratedSifgrs
    block = pellet
    execute_on = 'linear'
  []
  [fis_gas_released]
    type = ElementIntegralFisGasReleasedSifgrs
    block = pellet
    execute_on = 'linear'
  []
  [fis_gas_released_percentage]
    type = FGRPercent
    fission_gas_generated = fis_gas_produced
    fission_gas_released = fis_gas_released
    execute_on = 'linear'
  []
  [gas_volume]
    type = InternalVolume
    boundary = 9
    execute_on = 'initial linear'
  []
  [rod_total_power]
    type = ElementIntegralPower
    variable = temp
    fission_rate = fission_rate
    block = pellet
  []
  [rod_input_power]
    type = FunctionValuePostprocessor
    function = average_power_history
    scale_factor = 0.4 # rod height
  []
[]
[VectorPostprocessors]
  [fuel_radial_temperature_Sample1]
    type = LineValueSampler
    variable = temp
    start_point = '0.0 0.283 0.0'
    end_point = '0.002675 0.283 0.0'
    num_points = 200
    execute_on = final
    sort_by = x
    outputs = line_plot
  []
  [radial_porosity_Sample1]
    type = LineValueSampler
    variable = pore
    start_point = '0.0 0.283 0.0'
    end_point = '0.002675 0.283 0.0'
    num_points = 200
    execute_on = final
    sort_by = x
    outputs = line_plot
  []
  [fuel_radial_temperature_Sample2]
    type = LineValueSampler
    variable = temp
    start_point = '0.0 0.347 0.0'
    end_point = '0.002675 0.347 0.0'
    num_points = 200
    execute_on = final
    sort_by = x
    outputs = line_plot
  []
  [radial_porosity_Sample2]
    type = LineValueSampler
    variable = pore
    start_point = '0.0 0.347 0.0'
    end_point = '0.002675 0.347 0.0'
    num_points = 200
    execute_on = final
    sort_by = x
    outputs = line_plot
  []
  [fuel_radial_temperature_Sample3]
    type = LineValueSampler
    variable = temp
    start_point = '0.0 0.2 0.0'
    end_point = '0.002675 0.2 0.0'
    num_points = 200
    execute_on = final
    sort_by = x
    outputs = line_plot
  []
  [radial_porosity_Sample3]
    type = LineValueSampler
    variable = pore
    start_point = '0.0 0.2 0.0'
    end_point = '0.002675 0.2 0.0'
    num_points = 200
    execute_on = final
    sort_by = x
    outputs = line_plot
  []
[]

[PerformanceMetricOutputs]
[]

[Outputs]
  perf_graph = true
  exodus = true
  color = true
  csv = true
  [console]
    type = Console
    max_rows = 25
  []
  [line_plot]
    type = CSV
    execute_on = 'FINAL'
    time_step_interval =  1
    file_base = 1d
    create_final_symlink = true
  []
  [chkfile]
    type = CSV
    execute_on = FINAL
    show = 'ave_temp_interior max_pore'
  []
[]
[MultiApps]
  [sub]
    type = TransientMultiApp
    app_type = BisonApp
    execute_on = timestep_end
    catch_up = true
    max_catch_up_steps = 10
    positions_file = positions.txt
    input_files = b14_ptm010_pore.i
  []
[]
[Transfers]
  [temp_to_sub]
    type = MultiAppProjectionTransfer
    to_multi_app = sub
    source_variable = temp
    variable = temp
  []
  [pore_from_sub]
    type = MultiAppGeometricInterpolationTransfer
    from_multi_app = sub
    source_variable = pore
    variable = pore
  []
[]
[Debug]
  show_var_residual_norms = true
[]

(assessment/MOX/JOYO/B14/PTM001/analysis/b14_ptm001_2DRZ_t.i)

initial_fuel_density = 11057.75

[GlobalParams]
  density = ${initial_fuel_density}
  order = SECOND
  family = LAGRANGE
  energy_per_fission = 3.2e-11 # J/fission
  volumetric_locking_correction = false
  displacements = 'disp_x disp_y'
[]
[Problem]
  type = ReferenceResidualProblem
  extra_tag_vectors = 'ref'
  reference_vector = 'ref'
[]
[Mesh]
  coord_type = RZ
  [smeared_pellet_mesh]
    type = FuelPinMeshGenerator
    pellet_quantity = 1
    pellet_height = 0.4
    pellet_outer_radius = 0.0027
    pellet_mesh_density = customize
    clad_mesh_density = customize
    clad_gap_width = 0.00008
    clad_thickness = 0.00047
    clad_bot_gap_height = 1.0e-3
    bottom_clad_height = 2.24e-3
    top_clad_height = 2.24e-3
    clad_top_gap_height = 0.685
    elem_type = QUAD8
    nx_c = 4
    ny_c = 100
    nx_p = 20
    ny_p = 100
    ny_cu = 3
    ny_cl = 3
  []
  patch_size = 50
  patch_update_strategy = iteration
  partitioner = centroid
  centroid_partitioner_direction = y
[]
[UserObjects]
  [pin_geometry]
    type = FuelPinGeometry
  []
[]
[Variables]
  [temp]
    initial_condition = 295.0
  []
[]
[AuxVariables]
  [pore]
  []
  [fission_rate]
  []
  [burnup]
  []
  [radial_strain]
    order = CONSTANT
    family = MONOMIAL
  []
  [effective_creep_strain]
    block = clad
    order = CONSTANT
    family = MONOMIAL
  []
[]
[Functions]
  [power_history] #related to the LHGR at the midplane
    type = PiecewiseLinear
    x = ' 0  72000  158040  160200  246600  248400  249000.012  251280'
    y = ' 0  39814.5  39814.5  44289.3  44289.3  53927.4  53927.4  0'
  []
  [fast_neutron_flux_function]
    type = PiecewiseLinear
    x = '0 251280'
    y = '3.3e+15 3.3e+15'
  []
  [f_temp_out_clad]
    type = PiecewiseBilinear
    x = '0.000175  0.0464075  0.0843675  0.1075625  0.152025  0.1994625  0.2464725  0.2947475  0.356915  0.43356  0.49848  0.625  0.700475  0.797485  0.8723425  0.96'
    y = '0 251280'
    z = '295  295  295  295  295  295  295  295  295  295  295  295  295  295  295  295 634.94  662.273  676.998  686.217  706.339  727  743.358  758.311  780.069  799.077  815.576  846.374  860.233  875.494  882.809  889.8'
    scale_factor = 1
    axis = 1
  []
  [axial_peaking_factors]
    type = PiecewiseBilinear
    x = '0.000175  0.0464075  0.0843675  0.1075625  0.152025  0.1994625  0.2464725  0.2947475  0.356915  0.43356  0.49848  0.625  0.700475  0.797485  0.8723425  0.96'
    y = '0 251280'
    z = '0.751  0.752  0.767  0.796  0.82  0.852  0.875  0.915  0.944  0.963  0.988  1  0.985  0.955  0.913  0.846 0.751  0.752  0.767  0.796  0.82  0.852  0.875  0.915  0.944  0.963  0.988  1  0.985  0.955  0.913  0.846'
    scale_factor = 1
    axis = 1
  []
  [q]
    type = CompositeFunction
    functions = 'power_history axial_peaking_factors'
  []
  [average_power_history]
    type = PiecewiseLinear
    x = ' 0  72000  158040  160200  246600  248400  249000.012  251280'
    y = ' 0  34700  34700  38600  38600  47000  47000  0'
  []
[]
[Physics/SolidMechanics/QuasiStatic]
  [pellets]
    block = pellet
    add_variables = true
    strain = FINITE
    eigenstrain_names = 'fuel_thermal_strain fuel_volumetric_strain'
    generate_output = 'vonmises_stress stress_xx stress_yy stress_zz'
    use_finite_deform_jacobian = true
    extra_vector_tags = 'ref'
  []
  [clad]
    block = clad
    add_variables = true
    strain = FINITE
    eigenstrain_names = 'clad_thermal_eigenstrain'
    generate_output = 'vonmises_stress stress_xx stress_yy stress_zz'
    use_finite_deform_jacobian = true
    extra_vector_tags = 'ref'
  []
[]
[Kernels]
  [gravity]
    type = Gravity
    variable = disp_y
    value = -9.81
    extra_vector_tags = 'ref'
  []
  [heat]
    type = HeatConduction
    variable = temp
    extra_vector_tags = 'ref'
  []
  [heat_ie]
    type = HeatConductionTimeDerivative
    variable = temp
    extra_vector_tags = 'ref'
  []
  [heat_source]
    type = NeutronHeatSource
    variable = temp
    fission_rate = fission_rate
    extra_vector_tags = 'ref'
  []
[]
[AuxKernels]
  [fission_rate]
    type = FissionRateGeneral
    fission_rate_formulation = MOX
    variable = fission_rate
    block = pellet
    initial_porosity = 0.1372
    axial_power_profile = axial_peaking_factors
    rod_ave_lin_pow = power_history
    pellet_diameter = 0.0054
    execute_on = timestep_begin
    porosity = pore
    energy_per_fission = 3.2e-11
  []
  [burnup]
    type = BurnupAux
    block = pellet
    fission_rate = fission_rate
    variable = burnup
    execute_on = timestep_begin
  []
  [radial_strain]
    type = RankTwoAux
    rank_two_tensor = total_strain
    variable = radial_strain
    index_i = 0
    index_j = 0
    execute_on = timestep_end
  []
  [effective_creep_strain]
    type = MaterialRealAux
    property = effective_creep_strain
    variable = effective_creep_strain
    block = clad
    execute_on = timestep_end
  []
[]
[Contact]
  [pellet_clad_mechanical]
    primary = 5
    secondary = 10
    formulation = kinematic
    model = frictionless
    penalty = 1e7
  []
[]
[ThermalContact]
  [thermal_contact]
    type = GasGapHeatTransfer
    variable = temp
    primary = 5
    secondary = 10
    initial_moles = initial_moles
    gas_released = fis_gas_released
    contact_pressure = contact_pressure
    quadrature = true
  []
[]
[BCs]
  [temp_clad_outside]
    type = FunctionDirichletBC
    variable = temp
    function = f_temp_out_clad
    boundary = 2
  []
  [no_x_all]
    type = DirichletBC
    variable = disp_x
    boundary = 12
    value = 0.0
  []
  [no_y_clad_bottom]
    type = DirichletBC
    variable = disp_y
    boundary = 1
    value = 0.0
  []
  [no_y_fuel_bottom]
    type = DirichletBC
    variable = disp_y
    boundary = 20
    value = 0.0
  []
  [PlenumPressure]
    [plenumPressure]
      boundary = 9
      initial_pressure = 101325
      startup_time = 0
      R = 8.3145
      output_initial_moles = initial_moles
      temperature = ave_temp_interior
      volume = gas_volume
      material_input = fis_gas_released
      output = plenum_pressure
    []
  []
[]
[Materials]
  [fast_neutron_flux]
    type = FastNeutronFlux
    calculate_fluence = true
    block = clad
    flux_function = fast_neutron_flux_function
  []
  [fuel_thermal]
    type = MAMOXThermal
    temperature = temp
    porosity = pore
    block = pellet
    Am_content = 0.0237
    oxy_to_metal_ratio = 1.982
    output_properties = 'thermal_conductivity'
  []
  [fuel_density]
    type = StrainAdjustedDensity
    block = pellet
    strain_free_density = ${initial_fuel_density}
  []
  [fuel_elasticity_tensor]
    type = MAMOXElasticityTensor
    block = pellet
  []
  [elastic_stress]
    type = ComputeFiniteStrainElasticStress
    block = pellet
  []
  [fuel_thermal_expansion]
    type = MAMOXThermalExpansionEigenstrain
    block = pellet
    temperature = temp
    stress_free_temperature = 295.0
    oxygen_to_metal_ratio = 1.98
    eigenstrain_name = fuel_thermal_strain
  []
  [fuel_volumetric_swelling]
    type = UO2VolumetricSwellingEigenstrain
    gas_swelling_model_type = SIFGRS
    block = pellet
    temperature = temp
    burnup = burnup
    initial_fuel_density = ${initial_fuel_density}
    eigenstrain_name = fuel_volumetric_strain
  []
  [clad_thermal]
    type = SS316Thermal
    block = clad
    temperature = temp
  []
  [clad_density]
    type = StrainAdjustedDensity
    block = clad
    strain_free_density = 8000
  []
  [clad_elasticity_tensor]
    type = SS316ElasticityTensor
    block = clad
    temperature = temp
    elastic_constants_model = legacy_ifr
  []
  [thermal_expansion]
    type = SS316ThermalExpansionEigenstrain
    block = clad
    temperature = temp
    stress_free_temperature = 295.0
    eigenstrain_name = clad_thermal_eigenstrain
  []
  [clad_ss316creep]
    type = SS316CreepUpdate
    block = clad
    temperature = temp
    fast_neutron_flux = fast_neutron_flux
  []
  [clad_stress]
    type = ComputeMultipleInelasticStress
    tangent_operator = elastic
    inelastic_models = 'clad_ss316creep'
    block = clad
  []
  [fission_gas_release]
    type = UO2Sifgrs
    block = pellet
    temperature = temp
    burnup = burnup
    fission_rate = fission_rate
    grain_radius_const = 14e-06 #I'm keeping the grain radius const because the grain growth in MOX is probably different due to high Temp
    bubble_gb_limit = 1.0e+11
  []
[]
[Preconditioning]
  [SMP]
    type = SMP
    full = true
  []
[]
[Executioner]
  type = Transient
  solve_type = 'PJFNK'
  petsc_options = '-snes_ksp_ew'
  petsc_options_iname = '-pc_type -sub_pc_type'
  petsc_options_value = 'asm lu'
  line_search = 'none'
  fixed_point_max_its = 1
  fixed_point_abs_tol = 1e-3
  fixed_point_rel_tol = 1e-3
  l_max_its = 50
  l_tol = 8e-3
  nl_max_its = 50
  nl_rel_tol = 1e-3
  nl_abs_tol = 1e-3
  start_time = 0
  n_startup_steps = 1
  end_time = 251280
  dtmax = 10000
  dtmin = 0.25
  [TimeStepper]
    type = IterationAdaptiveDT
    dt = 1e1
    optimal_iterations = 15
    iteration_window = 2
    linear_iteration_ratio = 100
    growth_factor = 2
    cutback_factor = .5
    force_step_every_function_point = true
    timestep_limiting_function = power_history
  []
[]
[Postprocessors]
  [ave_temp_interior]
    type = ElementAverageValue
    variable = temp
    execute_on = 'initial linear'
  []
  [average_burnup]
    type = ElementAverageValue
    variable = burnup
  []
  [ave_pore]
    type = ElementAverageValue
    block = pellet
    variable = pore
  []
  [max_pore]
    type = NodalExtremeValue
    block = pellet
    value_type = max
    variable = pore
  []
  [min_pore]
    type = NodalExtremeValue
    block = pellet
    value_type = min
    variable = pore
  []
  [fis_gas_produced]
    type = ElementIntegralFisGasGeneratedSifgrs
    #  variable = temp
    block = pellet
    execute_on = 'linear'
  []
  [fis_gas_released]
    type = ElementIntegralFisGasReleasedSifgrs
    #  variable = temp
    block = pellet
    execute_on = 'linear'
  []
  [fis_gas_released_percentage]
    type = FGRPercent
    fission_gas_generated = fis_gas_produced
    fission_gas_released = fis_gas_released
    #  variable = temp
    execute_on = 'linear'
  []
  [gas_volume]
    type = InternalVolume
    boundary = 9
    execute_on = 'initial linear'
  []
  [rod_total_power]
    type = ElementIntegralPower
    variable = temp
    fission_rate = fission_rate
    block = pellet
  []
  [rod_input_power]
    type = FunctionValuePostprocessor
    function = average_power_history
    scale_factor = 0.4 # rod height
  []
[]
[VectorPostprocessors]
  [fuel_radial_temperature_Sample1]
    type = LineValueSampler
    variable = temp
    start_point = '0.0 0.283 0.0'
    end_point = '0.0027 0.283 0.0'
    num_points = 200
    execute_on = final
    sort_by = x
    outputs = line_plot
  []
  [radial_porosity_Sample1]
    type = LineValueSampler
    variable = pore
    start_point = '0.0 0.283 0.0'
    end_point = '0.0027 0.283 0.0'
    num_points = 200
    execute_on = final
    sort_by = x
    outputs = line_plot
  []
  [fuel_radial_temperature_Sample2]
    type = LineValueSampler
    variable = temp
    start_point = '0.0 0.347 0.0'
    end_point = '0.0027 0.347 0.0'
    num_points = 200
    execute_on = final
    sort_by = x
    outputs = line_plot
  []
  [radial_porosity_Sample2]
    type = LineValueSampler
    variable = pore
    start_point = '0.0 0.347 0.0'
    end_point = '0.0027 0.347 0.0'
    num_points = 200
    execute_on = final
    sort_by = x
    outputs = line_plot
  []
  [fuel_radial_temperature_Sample3]
    type = LineValueSampler
    variable = temp
    start_point = '0.0 0.2 0.0'
    end_point = '0.0027 0.2 0.0'
    num_points = 200
    execute_on = final
    sort_by = x
    outputs = line_plot
  []
  [radial_porosity_Sample3]
    type = LineValueSampler
    variable = pore
    start_point = '0.0 0.2 0.0'
    end_point = '0.0027 0.2 0.0'
    num_points = 200
    execute_on = final
    sort_by = x
    outputs = line_plot
  []
[]

[PerformanceMetricOutputs]
[]

[Outputs]
  perf_graph = true
  exodus = true
  color = true
  csv = true
  [console]
    type = Console
    max_rows = 25
  []
  [line_plot]
    type = CSV
    execute_on = 'FINAL'
    time_step_interval =  1
    file_base = 1d
    create_final_symlink = true
  []
  [chkfile]
    type = CSV
    execute_on = FINAL
    show = 'ave_temp_interior fis_gas_released_percentage max_pore'
  []
[]
[MultiApps]
  [sub]
    type = TransientMultiApp
    app_type = BisonApp
    execute_on = timestep_end
    catch_up = true
    max_catch_up_steps = 10
    positions_file = positions.txt
    input_files = b14_ptm001_pore.i
  []
[]
[Transfers]
  [temp_to_sub]
    type = MultiAppProjectionTransfer
    to_multi_app = sub
    source_variable = temp
    variable = temp
  []
  [pore_from_sub]
    type = MultiAppGeometricInterpolationTransfer
    from_multi_app = sub
    source_variable = pore
    variable = pore
  []
[]
[Debug]
  show_var_residual_norms = true
[]

(assessment/MOX/JOYO/MK-I/analysis/MK-I_75MW_master_old_bubble_gb_lim.i)


initial_fuel_density = 10836.8

[GlobalParams]
  density = ${initial_fuel_density}
  initial_porosity = 0.065
  order = SECOND
  family = LAGRANGE
  energy_per_fission = 3.2e-11 # J/fission
  volumetric_locking_correction = false
  displacements = 'disp_x disp_y'
[]
[Problem]
  type = ReferenceResidualProblem
  extra_tag_vectors = 'ref'
  reference_vector = 'ref'
[]
[Mesh]
  coord_type = RZ
  [smeared_pellet_mesh]
    type = FuelPinMeshGenerator
    pellet_quantity = 1
    pellet_height = 0.6
    pellet_outer_radius = 0.0027
    pellet_mesh_density = customize
    clad_mesh_density = customize
    clad_gap_width = 0.000100
    clad_thickness = 0.00035
    clad_bot_gap_height = 1.0e-3
    bottom_clad_height = 2.24e-3
    top_clad_height = 2.24e-3
    clad_top_gap_height = 0.599
    elem_type = QUAD8
    nx_c = 4
    ny_c = 200
    nx_p = 20
    ny_p = 200
    ny_cu = 3
    ny_cl = 3
  []
  patch_size = 50
  patch_update_strategy = iteration
  partitioner = centroid
  centroid_partitioner_direction = y
[]
[UserObjects]
  [pin_geometry]
    type = FuelPinGeometry
  []
[]
[Variables]
  [temp]
    initial_condition = 295.0
  []
[]
[AuxVariables]
  [pore]
  []
  [fission_rate]
    block = pellet
  []
  [burnup]
    block = pellet
  []
  [gas_gen_3]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [gas_grn_3]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [gas_bdr_3]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [gas_rel_3]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [bbl_bdr_2]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [vcn_bdr_2]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [atm_bbl_bdr]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [vcn_bbl_bdr]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [prs_bbl_bdr]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [prseq_bbl_bdr]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [rad_bbl_bdr]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [vol_bbl_bdr]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [GBCoverage]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [eff_diff_coeff]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [deltav_v0_bd]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [radial_strain]
    order = CONSTANT
    family = MONOMIAL
  []
  [effective_creep_strain]
    block = clad
    order = CONSTANT
    family = MONOMIAL
  []
  [gap_cond]
    order = CONSTANT
    family = MONOMIAL
  []
[]

[Functions]
  [power_history] #related to the LHGR at the midplane
    type = PiecewiseLinear
    x = '0 70000 25000000'
    y = '0 38974.7  38974.7'
  []
  [fast_neutron_flux_function]
    type = PiecewiseLinear
    x = '0 70000 25000000'
    y = '0 1.9e+19 1.9e+19'
  []
  [axial_peaking_factors]
    type = PiecewiseBilinear
    x = '0  0.071  0.146  0.221  0.296  0.37  0.443  0.566'
    y = '0 25000000'
    z = '0.889  1.041  1.152  1.173  1.129  0.971  0.782  0.672 0.889  1.041  1.152  1.173  1.129  0.971  0.782  0.672'
    scale_factor = 1
    axis = 1
  []
  [q]
    type = CompositeFunction
    functions = 'power_history axial_peaking_factors'
  []
  [average_power_history]
    type = PiecewiseLinear
    x = '0 70000 25000000'
    y = '0 32000 32000'
  []
  [clad_surface_temp]
    type = PiecewiseBilinear
    x = '0  0.071  0.146  0.221  0.296  0.37  0.443  0.566'
    y = '0 25000000'
    z = '295 295 295 295 295 295 295 295 593.58  606.36  619.13  630.26  640.87  651.76  662.67 673.67'
    scale_factor = 1
    axis = 1
  []
  [pressure_ramp]
    type = PiecewiseLinear
    x = '0 1'
    y = '1 1'
  []
[]
[Physics/SolidMechanics/QuasiStatic]
  [pellets]
    block = pellet
    add_variables = true
    strain = FINITE
    eigenstrain_names = 'fuel_thermal_strain fuel_volumetric_strain'
    generate_output = 'vonmises_stress stress_xx stress_yy stress_zz'
    use_finite_deform_jacobian = true
    extra_vector_tags = 'ref'
  []
  [clad]
    block = clad
    add_variables = true
    strain = FINITE
    eigenstrain_names = 'clad_thermal_eigenstrain'
    generate_output = 'vonmises_stress stress_xx stress_yy stress_zz'
    use_finite_deform_jacobian = true
    extra_vector_tags = 'ref'
  []
[]
[Kernels]
  [gravity]
    type = Gravity
    variable = disp_y
    value = -9.81
    extra_vector_tags = 'ref'
  []
  [heat]
    type = HeatConduction
    variable = temp
    extra_vector_tags = 'ref'
  []
  [heat_ie]
    type = HeatConductionTimeDerivative
    variable = temp
    extra_vector_tags = 'ref'
  []
  [heat_source]
    type = NeutronHeatSource
    variable = temp
    block = pellet
    fission_rate = fission_rate
    extra_vector_tags = 'ref'
  []
[]
[AuxKernels]
  [fission_rate]
    type = FissionRateGeneral
    fission_rate_formulation = MOX
    variable = fission_rate
    block = pellet
    initial_porosity = 0.065
    axial_power_profile = axial_peaking_factors
    rod_ave_lin_pow = power_history
    pellet_diameter = 0.0054
    execute_on = timestep_begin
    porosity = pore
  []
  [burnup]
    type = BurnupAux
    block = pellet
    fission_rate = fission_rate
    variable = burnup
    execute_on = timestep_begin
  []
  [fggen]
    type = MaterialRealAux
    variable = gas_gen_3
    property = gas_concentration_generated_total
    execute_on = timestep_end
  []
  [fggrn]
    type = MaterialRealAux
    variable = gas_grn_3
    property = gas_concentration_intra_total
    execute_on = timestep_end
  []
  [fgbdr]
    type = MaterialRealAux
    variable = gas_bdr_3
    property = gas_concentration_GB_bubble_volume
    execute_on = timestep_end
  []
  [fgrel]
    type = MaterialRealAux
    variable = gas_rel_3
    property = gas_concentration_release_total
    execute_on = timestep_end
  []
  [nbbl2]
    type = MaterialRealAux
    variable = bbl_bdr_2
    property = bubble_GB_surface_density
    execute_on = timestep_end
  []
  [nvcn2]
    type = MaterialRealAux
    variable = vcn_bdr_2
    property = vacancy_concentration_GB_surface
    execute_on = timestep_end
  []
  [atmbbl]
    type = MaterialRealAux
    variable = atm_bbl_bdr
    property = atom_per_bubble_GB
    execute_on = timestep_end
  []
  [vcnbbl]
    type = MaterialRealAux
    variable = vcn_bbl_bdr
    property = vacancy_per_bubble_GB
    execute_on = timestep_end
  []
  [prsbbl]
    type = MaterialRealAux
    variable = prs_bbl_bdr
    property = bubble_GB_pressure
    execute_on = timestep_end
  []
  [prseqbbl]
    type = MaterialRealAux
    variable = prseq_bbl_bdr
    property = bubble_GB_pressure_equilibrium
    execute_on = timestep_end
  []
  [radbbl]
    type = MaterialRealAux
    variable = rad_bbl_bdr
    property = bubble_radius_GB
    execute_on = timestep_end
  []
  [volbbl]
    type = MaterialRealAux
    variable = vol_bbl_bdr
    property = bubble_GB_volume
    execute_on = timestep_end
  []
  [frcvrg]
    type = MaterialRealAux
    variable = GBCoverage
    property = GBCoverage
    execute_on = timestep_end
  []
  [diffc]
    type = MaterialRealAux
    variable = eff_diff_coeff
    property = eff_diff_coeff
    execute_on = timestep_end
  []
  [dvv0bd]
    type = MaterialRealAux
    variable = deltav_v0_bd
    property = deltav_v0_bubble_GB
    execute_on = timestep_end
  []
  [radial_strain]
    type = RankTwoAux
    rank_two_tensor = total_strain
    variable = radial_strain
    index_i = 0
    index_j = 0
    execute_on = timestep_end
  []
  [effective_creep_strain]
    type = MaterialRealAux
    property = effective_creep_strain
    variable = effective_creep_strain
    block = clad
    execute_on = timestep_end
  []
  [conductance]
    type = MaterialRealAux
    property = gap_conductance
    variable = gap_cond
    boundary = 10
    execute_on = 'linear'
  []
[]
[Contact]
  [pellet_clad_mechanical]
    primary = 5
    secondary = 10
    formulation = kinematic
    model = frictionless
    penalty = 1e7
  []
[]
[ThermalContact]
  [thermal_contact]
    type = GasGapHeatTransfer
    variable = temp
    primary = 5
    secondary = 10
    initial_moles = initial_moles
    gas_released = fis_gas_released
    contact_pressure = contact_pressure
    quadrature = true
  []
[]
[BCs]
  [no_x_all]
    type = DirichletBC
    variable = disp_x
    boundary = '12'
    value = 0.0
  []
  [no_y_clad_bottom]
    type = DirichletBC
    variable = disp_y
    boundary = 1
    value = 0.0
  []
  [no_y_fuel_bottom]
    type = DirichletBC
    variable = disp_y
    boundary = 20
    value = 0.0
  []
  [temp_clad_out]
    type = FunctionDirichletBC
    variable = temp
    boundary = '2'
    function = clad_surface_temp
  []
  [Pressure]
    [coolantPressure]
      boundary = '1 2 3'
      factor = 101325
      function = pressure_ramp
    []
  []
  [PlenumPressure]
    [plenumPressure]
      boundary = 9
      initial_pressure = 300000
      startup_time = 0
      R = 8.3143
      output_initial_moles = initial_moles
      temperature = ave_temp_interior
      volume = gas_volume
      material_input = fis_gas_released
      output = plenum_pressure
    []
  []
[]
[Materials]
  [fast_neutron_flux]
    type = FastNeutronFlux
    calculate_fluence = true
    block = clad
    flux_function = fast_neutron_flux_function
  []
  [fuel_thermal]
    type = MAMOXThermal
    block = pellet
    temperature = temp
    Am_content = 0.0
    Np_content = 0.0
    porosity = pore
    output_properties = 'thermal_conductivity'
  []
  [fuel_elasticity_tensor]
    type = MAMOXElasticityTensor
    block = pellet
  []
  [elastic_stress]
    type = ComputeFiniteStrainElasticStress
    block = pellet
  []
  [fuel_thermal_expansion]
    type = MAMOXThermalExpansionEigenstrain
    block = pellet
    temperature = temp
    stress_free_temperature = 295.0
    oxygen_to_metal_ratio = 1.98
    eigenstrain_name = fuel_thermal_strain
  []
  [fuel_volumetric_swelling]
    type = UO2VolumetricSwellingEigenstrain
    gas_swelling_model_type = SIFGRS
    block = pellet
    temperature = temp
    burnup = burnup
    initial_fuel_density = 10836.8
    eigenstrain_name = fuel_volumetric_strain
  []
  [clad_thermal]
    type = SS316Thermal
    block = clad
    temperature = temp
  []
  [clad_density]
    type = StrainAdjustedDensity
    block = clad
    strain_free_density = 8000
  []
  [clad_elasticity_tensor]
    type = SS316ElasticityTensor
    block = clad
    temperature = temp
    elastic_constants_model = legacy_ifr
  []
  [thermal_expansion]
    type = SS316ThermalExpansionEigenstrain
    block = clad
    temperature = temp
    stress_free_temperature = 295.0
    eigenstrain_name = clad_thermal_eigenstrain
  []
  [clad_ss316creep]
    type = SS316CreepUpdate
    block = clad
    temperature = temp
    fast_neutron_flux = fast_neutron_flux
  []
  [clad_stress]
    type = ComputeMultipleInelasticStress
    tangent_operator = elastic
    inelastic_models = 'clad_ss316creep'
    block = clad
  []
  [fission_gas_release]
    type = UO2Sifgrs
    block = pellet
    temperature = temp
    burnup = burnup
    diff_coeff_option = TURNBULL_D1_4D2_4D3
    fission_rate = fission_rate
    grain_radius_const = 8.01e-6 #I'm keeping the grain radius const because the grain growth in MOX is probably different due to high Temp
  []
  [fuel_density]
    type = StrainAdjustedDensity
    block = pellet
    strain_free_density = ${initial_fuel_density}
  []
[]
[Preconditioning]
  [SMP]
    type = SMP
    full = true
  []
[]
[Executioner]
  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'
  fixed_point_abs_tol = 1e-5
  fixed_point_rel_tol = 1e-5
  fixed_point_max_its = 1
  l_max_its = 70
  l_tol = 8e-3
  nl_max_its = 70
  nl_rel_tol = 1e-5
  nl_abs_tol = 1e-5
  start_time = 0
  n_startup_steps = 1
  end_time = 25000000
  dtmax = 1e6
  dtmin = 0.25
  [TimeStepper]
    type = IterationAdaptiveDT
    dt = 5000
    optimal_iterations = 15
    iteration_window = 2
    linear_iteration_ratio = 100
    growth_factor = 2
    cutback_factor = .5
    force_step_every_function_point = true
    timestep_limiting_function = power_history
  []
[]
[Postprocessors]
  [ave_temp_interior]
    type = SideAverageValue
    boundary = 9
    variable = temp
    execute_on = 'initial linear'
  []
  [average_burnup]
    type = ElementAverageValue
    block = pellet
    variable = burnup
  []
  [clad_inner_vol]
    type = InternalVolume
    boundary = 7
    execute_on = 'initial timestep_end'
  []
  [pellet_volume]
    type = InternalVolume
    boundary = 8
    execute_on = 'initial timestep_end'
  []
  [avg_clad_temp]
    type = SideAverageValue
    boundary = 7
    variable = temp
    execute_on = 'initial timestep_end'
  []
  [fis_gas_produced]
    type = ElementIntegralFisGasGeneratedSifgrs
    block = pellet
    execute_on = 'linear'
  []
  [fis_gas_released]
    type = ElementIntegralFisGasReleasedSifgrs
    block = pellet
    execute_on = 'linear'
  []
  [fis_gas_released_percentage]
    type = FGRPercent
    fission_gas_generated = fis_gas_produced
    fission_gas_released = fis_gas_released
    execute_on = 'linear'
  []
  [fis_gas_grain]
    type = ElementIntegralFisGasGrainSifgrs
    block = pellet
    execute_on = 'linear'
  []
  [fis_gas_boundary]
    type = ElementIntegralFisGasBoundarySifgrs
    block = pellet
    execute_on = 'linear'
  []
  [gas_volume]
    type = InternalVolume
    boundary = 9
    execute_on = 'initial linear'
  []
  [flux_from_fuel]
    type = SideDiffusiveFluxIntegral
    variable = temp
    boundary = 10
    diffusivity = thermal_conductivity
  []
  [rod_total_power]
    type = ElementIntegralPower
    variable = temp
    fission_rate = fission_rate
    block = pellet
  []
  [rod_input_power]
    type = FunctionValuePostprocessor
    function = average_power_history
    scale_factor = 0.6 # rod height
  []
  [average_vonMises_fuel]
    type = ElementAverageValue
    variable = vonmises_stress
    block = pellet
  []
  [average_vonMises_clad]
    type = ElementAverageValue
    variable = vonmises_stress
    block = clad
  []
  [average_strain_rr_fuel]
    type = ElementAverageValue
    variable = radial_strain
    block = pellet
  []
  [average_strain_rr_clad]
    type = ElementAverageValue
    variable = radial_strain
    block = clad
  []
  [average_creep_strain_clad]
    type = ElementAverageValue
    variable = effective_creep_strain
    block = clad
  []
  [ave_pore]
    type = ElementAverageValue
    variable = pore
  []
  [max_pore]
    type = NodalExtremeValue
    value_type = max
    variable = pore
  []
  [min_pore]
    type = NodalExtremeValue
    value_type = min
    variable = pore
  []
[]

[PerformanceMetricOutputs]
[]

[Outputs]
  perf_graph = true
  exodus = true
  color = true
  csv = true
  [console]
    type = Console
    max_rows = 25
  []
  [chkfile]
    type = CSV
    execute_on = FINAL
    show = 'ave_temp_interior fis_gas_released_percentage max_pore'
  []
[]
[MultiApps]
  [sub]
    type = TransientMultiApp
    app_type = BisonApp
    execute_on = timestep_end
    catch_up = true
    max_catch_up_steps = 10
    positions_file = positions.txt
    input_files = MK-I_75MW_sub_old_bubble_gb_lim.i
  []
[]
[Transfers]
  [temp_to_sub]
    type = MultiAppProjectionTransfer
    to_multi_app = sub
    source_variable = temp
    variable = temp
  []
  [pore_from_sub]
    type = MultiAppGeometricInterpolationTransfer
    from_multi_app = sub
    source_variable = pore
    variable = pore
  []
[]
[Debug]
  show_var_residual_norms = true
  show_var_residual = 'temp disp_x disp_y'
[]

(test/tests/solid_mechanics/SS316_creep/creep_SS316_3d_test.i)

#
# The mesh is a 1x1x1 cube with a pressure of 10 MPa on the top face.
# Symmetry boundary conditions on three places provide a uniaxial stress field.
# The temperature is held constant at 973.15.  The yield stress is set at twice
# the load pressure. The solution is advanced through ten time steps of 1e4 for
# a total time of 1e5.
#
# The total strain at time 1e5 can be computed as given by "High Temperature Inelastic Behavior of the AUstenitic Steel AISI Type 316" by Altenbach and Gorash and "Irradiatian Creep and Swelling of AISI 316 to Exposures of 130 dpa at 385-400 degrees C" by Garner and Porter.
#
# e_tot = e_elas +
#         e_thermal_creep +
#         e_irradiation_creep
#
#       = P/E +
#         (a * exp(-Q_ln/(RT)) * sinh(b * exp(-Q_pw/(RT)) * sigma) +
#         [B + DS] * fast_neutron_flux * sigma) * dt
#
#         where P = pressure load
#               E = Youngs modulus
#               a,b = thermal creep coefficients
#               sigma = stress in MPa
#               Q_ln-Q_pw = activation energies
#               B = Irradiation creep coefficient "creep compliance"
#               D = Irradiation creep coefficent "coupling coefficient"
#               S = Swelling Rate
#               dt = total time
#               fast_neutron_flux = 3e17 n/m**2/s (applied) but model expects units of 10**22 n/cm**2/s
#               so, fast_neutron_flux = 3e17 * 0.0001 (convert from cm**2 to m**2) * 1e-22 (convert to 10*22)= 3e-9
#
#               Note that these strains are expressed as %, so they have to be divided by 100, thus  the term 0.01 in the strain calculations found below.
#
#
#
#
#   For this test, the analytical solution at t=1e5 is:
#
#               e_tot = 10e6/1.90e11 +
#                       ((21.04889 * exp(-170000/(8.314*973.15)) *
#                       sinh(0.07 * (10) * exp(-12000/(8.314*973.15)))) +
#                       [1e-6 + 0.6e-2 * 0.01] * 3e-9 * 10 * 0.01) * 1e5
#
#                     = 5.263e-5 + (((4.336e-12 * 1.012) + 1.83e-12 * 0.01) * 1e5)
#                     = 5.263e-5 + 4.386e-7 + 1.83e-9 = 5.307e-5
#   e_elas = 5.263e-5
#   e_thermal_creep = 4.386e-7
#   e_irradiation_creep = 1.83e-9
#   e_tot = 5.307e-5
#
#
#


[Mesh]
  use_displaced_mesh = false
  [mesh]
    type = GeneratedMeshGenerator
    dim = 3
  []
[]

[GlobalParams]
  displacements = 'disp_x disp_y disp_z'
  order = FIRST
  family = LAGRANGE
  temperature = temp
[]

[Variables]
  [temp]
    initial_condition = 973.15
  []
[]

[AuxVariables]
  [stress_yy]
    order = CONSTANT
    family = MONOMIAL
  []
  [elastic_strain_yy]
    order = CONSTANT
    family = MONOMIAL
  []
  [creep_strain_yy]
    order = CONSTANT
    family = MONOMIAL
  []
[]

[Physics]
  [SolidMechanics]
    [QuasiStatic]
      [all]
        strain = SMALL
        incremental = true
        add_variables = true
        generate_output = 'stress_yy elastic_strain_yy creep_strain_yy'
      []
    []
  []
[]

[Functions]
  [top_pull]
    type = PiecewiseLinear
    x = '0 1'
    y = '1 1'
  []
[]

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

[AuxKernels]
  [stress_yy]
    type = RankTwoAux
    rank_two_tensor = stress
    variable = stress_yy
    index_i = 1
    index_j = 1
  []
  [elastic_strain_yy]
    type = RankTwoAux
    rank_two_tensor = elastic_strain
    variable = elastic_strain_yy
    index_i = 1
    index_j = 1
  []
  [creep_strain_yy]
    type = RankTwoAux
    rank_two_tensor = creep_strain
    variable = creep_strain_yy
    index_i = 1
    index_j = 1
    execute_on = timestep_end
  []
[]

[BCs]
  [u_top_pull]
    type = Pressure
    variable = disp_y
    boundary = top
    factor = -10.0e6
    function = top_pull
  []
  [u_bottom_fix]
    type = DirichletBC
    variable = disp_y
    boundary = bottom
    value = 0.0
  []
  [u_back_fix]
    type = DirichletBC
    variable = disp_z
    boundary = back
    value = 0.0
  []
  [u_left_fix]
    type = DirichletBC
    variable = disp_x
    boundary = left
    value = 0.0
  []
  [temp_top_fix]
    type = DirichletBC
    variable = temp
    boundary = top
    value = 973.15
  []
  [temp_bottom_fix]
    type = DirichletBC
    variable = temp
    boundary = bottom
    value = 973.15
  []
[]

[Materials]
  [fast_neutron_flux]
    type = FastNeutronFlux
    calculate_fluence = true
    factor = 3e17
  []
  [elasticity_tensor]
    type = ComputeIsotropicElasticityTensor
    youngs_modulus = 1.90e11
    poissons_ratio = 0.265
  []
  [creep]
    type = SS316CreepUpdate
    fast_neutron_flux = fast_neutron_flux
  []
  [radial_return_stress]
    type = ComputeMultipleInelasticStress
    tangent_operator = elastic
    inelastic_models = 'creep'
  []
  [thermal]
    type = HeatConductionMaterial
    thermal_conductivity = 100.0
  []
  [density]
    type = StrainAdjustedDensity
    strain_free_density = 1.0
  []
[]

[Executioner]
  type = Transient

  solve_type = 'PJFNK'

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

  line_search = 'none'

  l_max_its  = 50
  nl_max_its = 20
  l_tol = 1e-5
  start_time = 0.0
  num_steps = 10
  dt = 10000
[]

[Outputs]
  exodus = true
[]

(assessment/MOX/JOYO/MK-II/analysis/MK-II_master_new_bubble_gb_lim.i)


initial_fuel_density = 10920.4

[GlobalParams]
  density = ${initial_fuel_density}
  initial_porosity = 0.07
  order = SECOND
  family = LAGRANGE
  energy_per_fission = 3.2e-11 # J/fission
  volumetric_locking_correction = false
  displacements = 'disp_x disp_y'
[]
[Problem]
  type = ReferenceResidualProblem
  extra_tag_vectors = 'ref'
  reference_vector = 'ref'
[]
[Mesh]
  coord_type = RZ
  [smeared_pellet_mesh]
    type = FuelPinMeshGenerator
    pellet_quantity = 1
    pellet_height = 0.55
    pellet_outer_radius = 0.002315
    pellet_mesh_density = customize
    clad_mesh_density = customize
    clad_gap_width = 0.000085
    clad_thickness = 0.00035
    clad_bot_gap_height = 1.0e-3
    bottom_clad_height = 2.24e-3
    top_clad_height = 2.24e-3
    clad_top_gap_height = 0.549
    elem_type = QUAD8
    nx_c = 4
    ny_c = 100
    nx_p = 10
    ny_p = 100
    ny_cu = 3
    ny_cl = 3
  []
  patch_size = 50
  patch_update_strategy = iteration
  partitioner = centroid
  centroid_partitioner_direction = y
[]
[UserObjects]
  [pin_geometry]
    type = FuelPinGeometry
  []
[]
[Variables]
  [temp]
    initial_condition = 295.0
  []
[]
[AuxVariables]
  [pore]
  []
  [fission_rate]
    block = pellet
  []
  [burnup]
    block = pellet
  []
  [gas_gen_3]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [gas_grn_3]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [gas_bdr_3]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [gas_rel_3]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [bbl_bdr_2]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [vcn_bdr_2]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [atm_bbl_bdr]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [vcn_bbl_bdr]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [prs_bbl_bdr]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [prseq_bbl_bdr]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [rad_bbl_bdr]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [vol_bbl_bdr]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [GBCoverage]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [eff_diff_coeff]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [deltav_v0_bd]
    order = CONSTANT
    family = MONOMIAL
    block = pellet
  []
  [radial_strain]
    order = CONSTANT
    family = MONOMIAL
  []
  [effective_creep_strain]
    block = clad
    order = CONSTANT
    family = MONOMIAL
  []
  [gap_cond]
    order = CONSTANT
    family = MONOMIAL
  []
[]
[Functions]
  [power_history] #related to the LHGR at the midplane
    type = PiecewiseLinear
    x = '0  70000  12970000'
    y = '0  48827.8  48827.8'
  []
  [fast_neutron_flux_function]
    type = PiecewiseLinear
    x = '0  70000  12970000'
    y = '0 2.6e+19 2.6e+19'
  []
  [axial_peaking_factors]
    type = PiecewiseBilinear
    x = '0  0.065  0.134  0.202  0.271  0.339  0.406  0.519'
    y = '0 12970000'
    z = '0.889  1.041  1.152  1.173  1.129  0.971  0.782  0.672 0.889  1.041  1.152  1.173  1.129  0.971  0.782  0.672'
    scale_factor = 1
    axis = 1
  []
  [q]
    type = CompositeFunction
    functions = 'power_history axial_peaking_factors'
  []
  [average_power_history]
    type = PiecewiseLinear
    x = '0  70000  12970000'
    y = '0  40000  40000'
  []
  [clad_surface_temp]
    type = PiecewiseBilinear
    x = '0  0.065  0.134  0.202  0.271  0.339  0.406  0.519'
    y = '0 12970000'
    z = '295  295  295  295  295  295  295  295 416.36  422.49  428.63  434.27  439.36  444.71  450.07  455.48'
    scale_factor = 1
    axis = 1
  []
  [pressure_ramp]
    type = PiecewiseLinear
    x = '0 1'
    y = '1 1'
  []
[]

[Physics/SolidMechanics/QuasiStatic]
  [pellets]
    block = pellet
    add_variables = true
    strain = FINITE
    eigenstrain_names = 'fuel_thermal_strain fuel_volumetric_strain'
    generate_output = 'vonmises_stress stress_xx stress_yy stress_zz'
    use_finite_deform_jacobian = true
    extra_vector_tags = 'ref'
  []
  [clad]
    block = clad
    add_variables = true
    strain = FINITE
    eigenstrain_names = 'clad_thermal_eigenstrain'
    generate_output = 'vonmises_stress stress_xx stress_yy stress_zz'
    use_finite_deform_jacobian = true
    extra_vector_tags = 'ref'
  []
[]
[Kernels]
  [gravity]
    type = Gravity
    variable = disp_y
    value = -9.81
    extra_vector_tags = 'ref'
  []
  [heat]
    type = HeatConduction
    variable = temp
    extra_vector_tags = 'ref'
  []
  [heat_ie]
    type = HeatConductionTimeDerivative
    variable = temp
    extra_vector_tags = 'ref'
  []
  [heat_source]
    type = NeutronHeatSource
    variable = temp
    block = pellet
    fission_rate = fission_rate
    extra_vector_tags = 'ref'
  []
[]
[AuxKernels]
  [fission_rate]
    type = FissionRateGeneral
    fission_rate_formulation = MOX
    variable = fission_rate
    block = pellet
    initial_porosity = 0.07
    axial_power_profile = axial_peaking_factors
    rod_ave_lin_pow = power_history
    pellet_diameter = 0.00463
    execute_on = timestep_begin
    porosity = pore
  []
  [burnup]
    type = BurnupAux
    block = pellet
    fission_rate = fission_rate
    variable = burnup
    execute_on = timestep_begin
  []
  [fggen]
    type = MaterialRealAux
    variable = gas_gen_3
    property = gas_concentration_generated_total
    execute_on = timestep_end
  []
  [fggrn]
    type = MaterialRealAux
    variable = gas_grn_3
    property = gas_concentration_intra_total
    execute_on = timestep_end
  []
  [fgbdr]
    type = MaterialRealAux
    variable = gas_bdr_3
    property = gas_concentration_GB_bubble_volume
    execute_on = timestep_end
  []
  [fgrel]
    type = MaterialRealAux
    variable = gas_rel_3
    property = gas_concentration_release_total
    execute_on = timestep_end
  []
  [nbbl2]
    type = MaterialRealAux
    variable = bbl_bdr_2
    property = bubble_GB_surface_density
    execute_on = timestep_end
  []
  [nvcn2]
    type = MaterialRealAux
    variable = vcn_bdr_2
    property = vacancy_concentration_GB_surface
    execute_on = timestep_end
  []
  [atmbbl]
    type = MaterialRealAux
    variable = atm_bbl_bdr
    property = atom_per_bubble_GB
    execute_on = timestep_end
  []
  [vcnbbl]
    type = MaterialRealAux
    variable = vcn_bbl_bdr
    property = vacancy_per_bubble_GB
    execute_on = timestep_end
  []
  [prsbbl]
    type = MaterialRealAux
    variable = prs_bbl_bdr
    property = bubble_GB_pressure
    execute_on = timestep_end
  []
  [prseqbbl]
    type = MaterialRealAux
    variable = prseq_bbl_bdr
    property = bubble_GB_pressure_equilibrium
    execute_on = timestep_end
  []
  [radbbl]
    type = MaterialRealAux
    variable = rad_bbl_bdr
    property = bubble_radius_GB
    execute_on = timestep_end
  []
  [volbbl]
    type = MaterialRealAux
    variable = vol_bbl_bdr
    property = bubble_GB_volume
    execute_on = timestep_end
  []
  [frcvrg]
    type = MaterialRealAux
    variable = GBCoverage
    property = GBCoverage
    execute_on = timestep_end
  []
  [diffc]
    type = MaterialRealAux
    variable = eff_diff_coeff
    property = eff_diff_coeff
    execute_on = timestep_end
  []
  [dvv0bd]
    type = MaterialRealAux
    variable = deltav_v0_bd
    property = deltav_v0_bubble_GB
    execute_on = timestep_end
  []
  [radial_strain]
    type = RankTwoAux
    rank_two_tensor = total_strain
    variable = radial_strain
    index_i = 0
    index_j = 0
    execute_on = timestep_end
  []
  [effective_creep_strain]
    type = MaterialRealAux
    property = effective_creep_strain
    variable = effective_creep_strain
    block = clad
    execute_on = timestep_end
  []
  [conductance]
    type = MaterialRealAux
    property = gap_conductance
    variable = gap_cond
    boundary = 10
    execute_on = 'linear'
  []
[]
[Contact]
  [pellet_clad_mechanical]
    primary = 5
    secondary = 10
    formulation = kinematic
    model = frictionless
    penalty = 1e7
  []
[]
[ThermalContact]
  [thermal_contact]
    type = GasGapHeatTransfer
    variable = temp
    primary = 5
    secondary = 10
    initial_moles = initial_moles
    gas_released = fis_gas_released
    contact_pressure = contact_pressure
    quadrature = true
  []
[]
[BCs]
  [no_x_all]
    type = DirichletBC
    variable = disp_x
    boundary = '12'
    value = 0.0
  []
  [no_y_clad_bottom]
    type = DirichletBC
    variable = disp_y
    boundary = 1
    value = 0.0
  []
  [no_y_fuel_bottom]
    type = DirichletBC
    variable = disp_y
    boundary = 20
    value = 0.0
  []
  [temp_clad_out]
    type = FunctionDirichletBC
    variable = temp
    boundary = '2'
    function = clad_surface_temp
  []
  [Pressure]
    [coolantPressure]
      boundary = '1 2 3'
      factor = 101325
      function = pressure_ramp
    []
  []
  [PlenumPressure]
    [plenumPressure]
      boundary = 9
      initial_pressure = 300000
      startup_time = 0
      R = 8.3143
      output_initial_moles = initial_moles
      temperature = ave_temp_interior
      volume = gas_volume
      material_input = fis_gas_released
      output = plenum_pressure
    []
  []
[]

[Materials]
  [fast_neutron_flux]
    type = FastNeutronFlux
    calculate_fluence = true
    block = clad
    flux_function = fast_neutron_flux_function
  []
  [fuel_thermal]
    type = MAMOXThermal
    block = pellet
    temperature = temp
    Am_content = 0.0
    Np_content = 0.0
    porosity = pore
    output_properties = 'thermal_conductivity'
  []
  [fuel_elasticity_tensor]
    type = MAMOXElasticityTensor
    block = pellet
  []
  [elastic_stress]
    type = ComputeFiniteStrainElasticStress
    block = pellet
  []
  [fuel_thermal_expansion]
    type = MAMOXThermalExpansionEigenstrain
    block = pellet
    temperature = temp
    stress_free_temperature = 295.0
    oxygen_to_metal_ratio = 1.98
    eigenstrain_name = fuel_thermal_strain
  []
  [fuel_volumetric_swelling]
    type = UO2VolumetricSwellingEigenstrain
    gas_swelling_model_type = SIFGRS
    block = pellet
    temperature = temp
    burnup = burnup
    initial_fuel_density = 10920.4
    eigenstrain_name = fuel_volumetric_strain
  []
  [clad_thermal]
    type = SS316Thermal
    block = clad
    temperature = temp
  []
  [clad_density]
    type = StrainAdjustedDensity
    block = clad
    strain_free_density = 8000
  []
  [clad_elasticity_tensor]
    type = SS316ElasticityTensor
    block = clad
    temperature = temp
    elastic_constants_model = legacy_ifr
  []
  [thermal_expansion]
    type = SS316ThermalExpansionEigenstrain
    block = clad
    temperature = temp
    stress_free_temperature = 295.0
    eigenstrain_name = clad_thermal_eigenstrain
  []
  [clad_ss316creep]
    type = SS316CreepUpdate
    block = clad
    temperature = temp
    fast_neutron_flux = fast_neutron_flux
  []
  [clad_stress]
    type = ComputeMultipleInelasticStress
    tangent_operator = elastic
    inelastic_models = 'clad_ss316creep'
    block = clad
  []
  [fission_gas_release]
    type = UO2Sifgrs
    block = pellet
    temperature = temp
    burnup = burnup
    diff_coeff_option = TURNBULL_D1_4D2_4D3
    fission_rate = fission_rate
    grain_radius_const = 8.01e-6
    bubble_gb_limit = 1.0e+11
  []
  [fuel_density]
    type = StrainAdjustedDensity
    block = pellet
    strain_free_density = ${initial_fuel_density}
  []
[]

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

[Executioner]
  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'
  fixed_point_abs_tol = 1e-3
  fixed_point_rel_tol = 1e-3
  fixed_point_max_its = 1
  l_max_its = 70
  l_tol = 8e-3
  nl_max_its = 70
  nl_rel_tol = 1e-3
  nl_abs_tol = 1e-3
  start_time = 0
  n_startup_steps = 1
  end_time = 12970000
  dtmax = 5e5
  dtmin = 0.25
  [TimeStepper]
    type = IterationAdaptiveDT
    dt = 5000
    optimal_iterations = 15
    iteration_window = 2
    linear_iteration_ratio = 100
    growth_factor = 2
    cutback_factor = .5
    force_step_every_function_point = true
    timestep_limiting_function = power_history
  []
[]
[Postprocessors]
  [ave_temp_interior]
    type = SideAverageValue
    boundary = 9
    variable = temp
    execute_on = 'initial linear'
  []
  [average_burnup]
    type = ElementAverageValue
    block = pellet
    variable = burnup
  []
  [clad_inner_vol]
    type = InternalVolume
    boundary = 7
    execute_on = 'initial timestep_end'
  []
  [pellet_volume]
    type = InternalVolume
    boundary = 8
    execute_on = 'initial timestep_end'
  []
  [avg_clad_temp]
    type = SideAverageValue
    boundary = 7
    variable = temp
    execute_on = 'initial timestep_end'
  []
  [fis_gas_produced]
    type = ElementIntegralFisGasGeneratedSifgrs
    block = pellet
    execute_on = 'linear'
  []
  [fis_gas_released]
    type = ElementIntegralFisGasReleasedSifgrs
    block = pellet
    execute_on = 'linear'
  []
  [fis_gas_released_percentage]
    type = FGRPercent
    fission_gas_generated = fis_gas_produced
    fission_gas_released = fis_gas_released
    execute_on = 'linear'
  []
  [fis_gas_grain]
    type = ElementIntegralFisGasGrainSifgrs
    block = pellet
    execute_on = 'linear'
  []
  [fis_gas_boundary]
    type = ElementIntegralFisGasBoundarySifgrs
    block = pellet
    execute_on = 'linear'
  []
  [gas_volume]
    type = InternalVolume
    boundary = 9
    execute_on = 'initial linear'
  []
  [flux_from_fuel]
    type = SideDiffusiveFluxIntegral
    variable = temp
    boundary = 10
    diffusivity = thermal_conductivity
  []
  [rod_total_power]
    type = ElementIntegralPower
    variable = temp
    fission_rate = fission_rate
    block = pellet
  []
  [rod_input_power]
    type = FunctionValuePostprocessor
    function = average_power_history
    scale_factor = 0.55 # rod height
  []
  [average_vonMises_fuel]
    type = ElementAverageValue
    variable = vonmises_stress
    block = pellet
  []
  [average_vonMises_clad]
    type = ElementAverageValue
    variable = vonmises_stress
    block = clad
  []
  [average_strain_rr_fuel]
    type = ElementAverageValue
    variable = radial_strain
    block = pellet
  []
  [average_strain_rr_clad]
    type = ElementAverageValue
    variable = radial_strain
    block = clad
  []
  [average_creep_strain_clad]
    type = ElementAverageValue
    variable = effective_creep_strain
    block = clad
  []
  [ave_pore]
    type = ElementAverageValue
    variable = pore
  []
  [max_pore]
    type = NodalExtremeValue
    value_type = max
    variable = pore
  []
  [min_pore]
    type = NodalExtremeValue
    value_type = min
    variable = pore
  []
[]

[PerformanceMetricOutputs]
[]

[Outputs]
  perf_graph = true
  exodus = true
  color = true
  csv = true
  [console]
    type = Console
    max_rows = 25
  []
  [chkfile]
    type = CSV
    execute_on = FINAL
    show = 'ave_temp_interior fis_gas_released_percentage max_pore'
  []
[]
[MultiApps]
  [sub]
    type = TransientMultiApp
    app_type = BisonApp
    execute_on = timestep_end
    catch_up = true
    max_catch_up_steps = 10
    positions_file = positions.txt
    input_files = MK-II_sub_new_bubble_gb_lim.i
  []
[]
[Transfers]
  [temp_to_sub]
    type = MultiAppProjectionTransfer
    to_multi_app = sub
    source_variable = temp
    variable = temp
  []
  [pore_from_sub]
    type = MultiAppGeometricInterpolationTransfer
    from_multi_app = sub
    source_variable = pore
    variable = pore
  []
[]
[Debug]
  show_var_residual_norms = true
  show_var_residual = 'temp disp_x disp_y'
[]

Description
Example Input Syntax
Input Parameters
Input Files
References