FeCrAlElasticityTensor

Computes Young's modulus and Poisson's ratio as a function of temperature for FeCrAl alloys.

Description

This model computes the isotropic elasticity tensor for FeCrAl alloys by calculating the Young's modulus and Poisson's ratio as a function of temperature. The Young's modulus of various commercial FeCrAl alloys (MA956, PM2000, Kanthal APMT, and Fecralloy) as a function of temperature has been obtained from the datasheets of the manufacturers and are tabulated below. In BISON, the Young's modulus is linearly interpolated between the known values.

Table 1: The Young's modulus as a function of temperature for MA956** (Corporation, 2004)

Temperature (K)	Young's Modulus (GPa)
310.45	172
362.30	170
424.44	167
506.42	161
592.27	153
657.76	146
754.04	139
828.41	131
900.57	121
984.39	114
1066.55	107
1146.44	99.4
1208.19	89.1
1269.02	77.0
1347.65	72.7

Table 2: The Young's modulus as a function of temperature for PM2000 (MatWeb, 2014)

Temperature (K)	Young's Modulus (GPa)
473.15	160
673.15	145
873.15	125
973.15	115
1073.15	110
1223.15	95

Table 3: The Young's modulus as a function of temperature for Kanthal APMT (Sandvik, 2012)

Temperature (K)	Young's Modulus (GPa)
293.15	220
373.15	210
473.15	205
673.15	190
873.15	170
1073.15	150
1273.15	130

For Fecralloy (MatWeb, 2014) the Young's modulus is only given at room temperature (293.15 K) as 180 GPa. The Poisson's ratio of the commercial alloys are independent of temperature and tabulated below.

Table 4: The Poisson's ratio of the commercial FeCrAl alloys

Alloy	Poisson's Ratio
MA956	0.3
PM2000	0.33
Kanthal APMT	0.3
Fecralloy	0.3

For the laboratory optimized FeCrAl alloys from Oak Ridge National Laboratory known as C06M, C35M, and C36M the temperature dependent Young's modulus and Poisson's ratio are given by Field et al. (2018): $var element = document.getElementById("moose-equation-6dc8285c-8494-465a-b343-905440090d2f");katex.render(" E = -5.46\\times 10^{-5}T^2 - 3.85\\times 10^{-2}T + 199", element, {displayMode:true,throwOnError:false});$ $var element = document.getElementById("moose-equation-244ba67e-5bc1-4a01-baa9-9216351c982a");katex.render(" \\nu = 4.46\\times 10^{-5}T + 0.27", element, {displayMode:true,throwOnError:false});$ where $var element = document.getElementById("moose-equation-d66224d7-22ba-4054-b607-db618ffee262");katex.render("E", element, {displayMode:false,throwOnError:false});$ is the Young's modulus (GPa), $var element = document.getElementById("moose-equation-9c675402-4c0b-487e-bd6c-f4f901e06029");katex.render("\\nu", element, {displayMode:false,throwOnError:false});$ is Poisson's ratio, and $var element = document.getElementById("moose-equation-ccaff102-a8c3-464f-91c0-6891647e6425");katex.render("T", element, {displayMode:false,throwOnError:false});$ is the temperature ( $var element = document.getElementById("moose-equation-b022e181-6105-4f24-a750-5ad4da1ce307");katex.render("^\\circ", element, {displayMode:false,throwOnError:false});$ C).

Range of Applicability and Uncertainty

The temperature range over which the models for the elastic constants are valid depends upon the material analyzed. For C06M, C35M, and C36M the correlations for both Young's modulus and Poisson's ratio are valid from 298 K to 1123 K. The piecewise linear data for Young's modulus is valid from 310.45 K to 1347.65 K for MA956, 293.15 K to 1273.15 K for Kanthal APMT, and 473.15 to 1347.65 for PM2000. The Young's modulus of Fecralloy is taken at room temperature (293.15 K).

Young's Modulus Uncertainty

The data for which the models for Young's modulus are based for Kanthal APMT, PM2000, MA956, and Fecralloy did not have associated experimental uncertainty. Therefore, the uncertainty in the models can not be directly computed. For conservatism on the best estimate model it is recommended that an uncertainty of $var element = document.getElementById("moose-equation-66478ae3-5d17-473a-9558-d09ba09873fc");katex.render("\\pm", element, {displayMode:false,throwOnError:false});$ 6% be taken for these alloys at all temperatures.

For the laboratory optimized alloys C06M, C35M, and C36M the experimental data used to develop the model are provided in the handbook. The measurement uncertainty is associated with the different experimental conditions in which the Young's modulus was measured, including heating and cooling. Figure 1 and Figure 2 illustrate the lower and upper bounds of the experimental data provided in the Handbook. A polynomial fit of the same order as the mean equation used for the model was applied to each data set with the equation provided on the figure. These polynomial fits represent the lower and upper uncertainties on the model. Figure 3 shows the Young's modulus model as a solid line with the associated uncertainty bands by dashed lines. Figure 4 illustrates the evolution of the lower and upper uncertainty on the mean as a function of temperature. It is observed that the lower uncertainty band is larger than the upper uncertainty band for all temperatures and the uncertainty is larger at higher temperatures.

Figure 1: Linear fit to lower bound of experimental data from FeCrAl Handbook for Young's modulus model used for C06M, C35M, and C36M

Figure 2: Linear fit to upper bound of experimental data from FeCrAl Handbook for Young's modulus model used for C06M, C35M, and C36M

Figure 3: Young's modulus with associated uncertainty bands for the model used for C06M, C35M, and C36M

Figure 4: Percent difference from the mean for the lower and upper bound uncertainty on Young's modulus

Poisson's Ratio Uncertainty

The data for which the models for Poisson's ratio are based for Kanthal APMT, PM2000, MA956, and Fecralloy did not have associated experimental uncertainty. Therefore, the uncertainty in the models can not be directly computed. For conservatism on the best estimate model it is recommended that an uncertainty of $var element = document.getElementById("moose-equation-085468e9-e7d5-4cfe-a48a-1cf21a6a561a");katex.render("\\pm", element, {displayMode:false,throwOnError:false});$ 10% be taken for these alloys at all temperatures.

For the laboratory optimized alloys C06M, C35M, and C36M the experimental data used to develop the model are provided in the handbook. The measurement uncertainty is associated with the different experimental conditions in which the Poisson's ratio was measured, including heating and cooling. Figure 5 and Figure 6 illustrate the lower and upper bounds of the experimental data provided in the Handbook. A polynomial fit of the same order as the mean equation used for the model was applied to each data set with the equation provided on the figure. These polynomial fits represent the lower and upper uncertainties on the model. Figure 7 shows the Poisson's ratio model as a solid line with the associated uncertainty bands by dashed lines. Figure 8 illustrates the evolution of the lower and upper uncertainty on the mean as a function of temperature. It is observed that the lower uncertainty band is larger than the upper uncertainty band for all temperatures and the uncertainty is larger at higher temperatures.

Figure 5: Linear fit to lower bound of experimental data from FeCrAl Handbook for Poisson's ratio model used for C06M, C35M, and C36M

Figure 6: Linear fit to upper bound of experimental data from FeCrAl Handbook for Poisson's ratio model used for C06M, C35M, and C36M

Figure 7: Poisson's ratio with associated uncertainty bands for the model used for C06M, C35M, and C36M

Figure 8: Percent difference from the mean for the lower and upper bound uncertainty on Poisson's ratio

Example Input Syntax

[Materials<<<{"href": "../../../syntax/Materials/index.html"}>>>]
  [elasticity_tensor]
    type = FeCrAlElasticityTensor<<<{"description": "Computes Young's modulus and Poisson's ratio as a function of temperature for FeCrAl alloys.", "href": "FeCrAlElasticityTensor.html"}>>>
    temperature<<<{"description": "The coupled temperature (K)"}>>> = temperature
    fecral_material_type<<<{"description": "The FeCrAl alloy being used for the cladding. Choices are: APMT C06M C35M C36M MA956 PM2000 FECRALLOY"}>>> = MA956
  []
[]

Input Parameters

temperatureThe coupled temperature (K)
C++ Type:std::vector<VariableName>
Unit:(no unit assumed)
Controllable:No
Description:The coupled temperature (K)

Required Parameters

base_nameOptional parameter that allows the user to define multiple mechanics material systems on the same block, i.e. for multiple phases
C++ Type:std::string
Controllable:No
Description:Optional parameter that allows the user to define multiple mechanics material systems on the same block, i.e. for multiple phases
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
computeTrueWhen false, MOOSE will not call compute methods on this material. The user must call computeProperties() after retrieving the MaterialBase via MaterialBasePropertyInterface::getMaterialBase(). Non-computed MaterialBases are not sorted for dependencies.
Default:True
C++ Type:bool
Controllable:No
Description:When false, MOOSE will not call compute methods on this material. The user must call computeProperties() after retrieving the MaterialBase via MaterialBasePropertyInterface::getMaterialBase(). Non-computed MaterialBases are not sorted for dependencies.
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.
elasticity_tensor_prefactorOptional function to use as a scalar prefactor on the elasticity tensor.
C++ Type:FunctionName
Unit:(no unit assumed)
Controllable:No
Description:Optional function to use as a scalar prefactor on the elasticity tensor.
fecral_material_typeAPMTThe FeCrAl alloy being used for the cladding. Choices are: APMT C06M C35M C36M MA956 PM2000 FECRALLOY
Default:APMT
C++ Type:MooseEnum
Options:APMT, C06M, C35M, C36M, MA956, PM2000, FECRALLOY
Controllable:No
Description:The FeCrAl alloy being used for the cladding. Choices are: APMT C06M C35M C36M MA956 PM2000 FECRALLOY
poissons_ratio_functionPoisson's ratio as a function of temperature.
C++ Type:FunctionName
Unit:(no unit assumed)
Controllable:No
Description:Poisson's ratio as a function of temperature.
youngs_modulus_functionYoung's modulus as a function of temperature.
C++ Type:FunctionName
Unit:(no unit assumed)
Controllable:No
Description:Young's modulus as a function of temperature.

Optional Parameters

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

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

poissons_ratio_scale_factor1Scale factor to be applied to the Poisson's ratio. Used for calibration and sensitivity studies
Default:1
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:Scale factor to be applied to the Poisson's ratio. Used for calibration and sensitivity studies
youngs_modulus_scale_factor1Scale factor to be applied to the Young's modulus. Used for calibration and sensitivity studies
Default:1
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:Scale factor to be applied to the Young's modulus. Used for calibration and sensitivity studies

Advanced: Scaling Factors Parameters

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

Material Property Retrieval Parameters

Input Files

(test/tests/solid_mechanics/fecral_elasticity_tensor/elasticity_PM2000.i)
(test/tests/solid_mechanics/fecral_elasticity_tensor/elasticity_Fecralloy.i)
(test/tests/solid_mechanics/fecral_elasticity_tensor/elasticity_APMT.i)
(examples/accident_tolerant_fuel/uo2_fecral/uo2_fecral.i)
(test/tests/solid_mechanics/fecral_elasticity_tensor/elasticity_C06M.i)
(test/tests/solid_mechanics/fecral_elasticity_tensor/elasticity_MA956.i)

References

Special Metals Corporation. Special Metals Incoloy alloy MA956. www.specialmetals.com/documents/Incoloy+alloy+MA956.pdf, 2004.

@MISC{ma956_2004,
    author = "Corporation, Special Metals",
    title = "{Special Metals Incoloy alloy MA956}",
    howpublished = "www.specialmetals.com/documents/Incoloy+alloy+MA956.pdf",
    year = "2004"
}

K. G. Field, M. A. Snead, Y. Yamamoto, and K. A. Terrani. Handbook on the material properties of fecral alloys for nuclear power production applications. Technical Report ORNL/SPR-2018/905 Rev. 1, Oak Ridge National Laboratory, 2018.

@TECHREPORT{fecral_handbook,
    author = "Field, K. G. and Snead, M. A. and Yamamoto, Y. and Terrani, K. A.",
    title = "Handbook on the Material Properties of FeCrAl Alloys for Nuclear Power Production Applications",
    institution = "Oak Ridge National Laboratory",
    number = "ORNL/SPR-2018/905 Rev. 1",
    year = "2018"
}

MatWeb. Resistalloy International Fecralloy Electrical Resistance Steel. http://www.matweb.com/search/datasheet.aspx?MatGUID=c2427c6297594858bedac2a4e5981d2f, 2014.

@MISC{fecralloy_matweb,
    author = "MatWeb",
    title = "{Resistalloy International Fecralloy Electrical Resistance Steel}",
    howpublished = "http://www.matweb.com/search/datasheet.aspx?MatGUID=c2427c6297594858bedac2a4e5981d2f",
    year = "2014"
}

MatWeb. Schwarzkopf Plansee PM 2000. http://www.matweb.com/search/datasheet.aspx?matguid=21e9ec9a0de24b47bcf69ab11c375567, 2014.

@MISC{pm2000_matweb,
    author = "MatWeb",
    title = "{Schwarzkopf Plansee PM 2000}",
    howpublished = "http://www.matweb.com/search/datasheet.aspx?matguid=21e9ec9a0de24b47bcf69ab11c375567",
    year = "2014"
}

Sandvik. Kanthal APMT Material Database. http://kanthal.com/en/products/material-datasheets/tube/kanthal-apmt/, 2012.

@MISC{kanthal_apmt,
    author = "Sandvik",
    title = "{Kanthal APMT Material Database}",
    howpublished = "http://kanthal.com/en/products/material-datasheets/tube/kanthal-apmt/",
    year = "2012"
}

(test/tests/solid_mechanics/fecral_elasticity_tensor/elasticity_MA956.i)

# This test case is prepared to test the elastic moduli of MA956 FeCrAl alloy.
#
# The geometry represents a ring of MA956 with inner radius of 0.005 m,
# outer radius of 0.0055 m, and a height of 0.01 m.
#
# An axial load is applied to the top of the ring with a value of 100 MPa.
# The strain in the radial direction can be calculated using the Young's modulus
# and the Poisson's ratio from the applied 100.0 MPa pressure.
#
# epsilon_{rr} = nu * sigma_{axial} / E
#
# The Young's modulus (E) for MA956 is calculated via a piecewse linear function.
#
# The Poisson's ratio of MA956 is calculated is taken as 0.3.
#
# The results of BISON compared to the analytical solution for an increase in
# temperature under constant pressure are summarized:
#
#                          BISON             Analytical
#    Temperature (K)    strain_rr (m/m)    strain_rr (m/m)
#        420              1.794104e-4        1.794104e-4
#        540              1.900288e-4        1.900288e-4
#        660              2.057089e-4        2.057089e-4
#        780              2.202522e-4        2.202522e-4
#        900              2.477721e-4        2.477721e-4
#       1020              2.703530e-4        2.703530e-4
#       1140              2.999621e-4        2.999621e-4
#       1260              3.807386e-4        3.807386e-4
#       1380              4.126547e-4        4.126547e-4
#       1500              4.126547e-4        4.126547e-4
#

[Mesh]
  coord_type = RZ
  [mesh]
    type = GeneratedMeshGenerator
    dim = 2
    xmin = 0.005
    xmax = 0.0055
    ymin = 0
    ymax = 0.01
  []
[]

[GlobalParams]
  displacements = 'disp_x disp_y'
[]

[Variables]
  [temperature]
    order = FIRST
    family = LAGRANGE
    initial_condition = 300.0
  []
[]

[Physics]
  [SolidMechanics]
    [QuasiStatic]
      [all]
        strain = FINITE
        block = 0
        add_variables = true
        generate_output = 'elastic_strain_xx'
      []
    []
  []
[]

[Functions]
  [temperature_function]
    type = PiecewiseLinear
    x = '0       10'
    y = '300   1500'
  []
  [pressure_function]
    type = PiecewiseLinear
    x = '0  10'
    y = '1   1'
  []
[]

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

[BCs]
  [top_pressure] # Simplified BC to apply pressure in only disp_y
    type = Pressure
    boundary = top
    variable = disp_y
    function = pressure_function
    factor = 100e6
  []

  [u_bottom_fix]
    type = DirichletBC
    variable = disp_y
    boundary = bottom
    value = 0.0
  []

  [temp_bc]
    type = FunctionDirichletBC
    variable = temperature
    boundary = 'bottom top left right'
    function = temperature_function
  []
[]

[Materials]
  [elasticity_tensor]
    type = FeCrAlElasticityTensor
    temperature = temperature
    fecral_material_type = MA956
  []
  [stress]
    type = ComputeFiniteStrainElasticStress
  []
  [clad_density]
    type = StrainAdjustedDensity
    strain_free_density = 1.0
  []
  [thermal]
    type = FeCrAlThermal
    temperature = temperature
    material = MA956
  []
[]

[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  = 100
  nl_max_its = 100
  nl_rel_tol = 1e-8
  nl_abs_tol = 1e-8
  l_tol      = 1e-8
  start_time = 0.0
  end_time   = 10
  dt         = 1
[]

[Postprocessors]
  [elastic_strain_rr]
    type = ElementAverageValue
    variable = elastic_strain_xx
  []
  [temp]
    type = AverageNodalVariableValue
    variable = temperature
    execute_on = 'initial timestep_end'
  []
[]

[Outputs]
  exodus = true
[]

(test/tests/solid_mechanics/fecral_elasticity_tensor/elasticity_PM2000.i)

# This test case is prepared to test the elastic moduli of PM2000 FeCrAl alloy.
#
# The geometry represents a ring of PM2000 with inner radius of 0.005 m,
# outer radius of 0.0055 m, and a height of 0.01 m.
#
# An axial load is applied to the top of the ring with a value of 100 MPa.
# The strain in the radial direction can be calculated using the Young's modulus
# and the Poisson's ratio from the applied 100.0 MPa pressure.
#
# epsilon_{rr} = nu * sigma_{axial} / E
#
# The Young's modulus (E) for PM2000 is calculated via a piecewse linear function.
#
# The Poisson's ratio of PM2000 is calculated is taken as 0.33.
#
# The results of BISON compared to the analytical solution for an increase in
# temperature under constant pressure are summarized:
#
#                          BISON             Analytical
#    Temperature (K)    strain_rr (m/m)    strain_rr (m/m)
#        420              2.062500e-4        2.062500e-4
#        540              2.129221e-4        2.129221e-4
#        660              2.260487e-4        2.260487e-4
#        780              2.456911e-4        2.456911e-4
#        900              2.697952e-4        2.697952e-4
#       1020              2.929232e-4        2.929232e-4
#       1140              3.194115e-4        3.194115e-4
#       1260              3.473684e-4        3.473684e-4
#       1380              3.473684e-4        3.473684e-4
#       1500              3.473684e-4        3.473684e-4
#

[Mesh]
  coord_type = RZ
  [mesh]
    type = GeneratedMeshGenerator
    dim = 2
    xmin = 0.005
    xmax = 0.0055
    ymin = 0
    ymax = 0.01
  []
[]

[GlobalParams]
  displacements = 'disp_x disp_y'
[]

[Variables]
  [temperature]
    order = FIRST
    family = LAGRANGE
    initial_condition = 300.0
  []
[]

[Physics]
  [SolidMechanics]
    [QuasiStatic]
      [all]
        strain = FINITE
        block = 0
        add_variables = true
        generate_output = 'elastic_strain_xx'
      []
    []
  []
[]

[Functions]
  [temperature_function]
    type = PiecewiseLinear
    x = '0       10'
    y = '300   1500'
  []
  [pressure_function]
    type = PiecewiseLinear
    x = '0  10'
    y = '1   1'
  []
[]

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

[BCs]
  [top_pressure] # Simplified BC to apply pressure in only disp_y
    type = Pressure
    boundary = top
    variable = disp_y
    function = pressure_function
    factor = 100e6
  []

  [u_bottom_fix]
    type = DirichletBC
    variable = disp_y
    boundary = bottom
    value = 0.0
  []

  [temp_bc]
    type = FunctionDirichletBC
    variable = temperature
    boundary = 'bottom top left right'
    function = temperature_function
  []
[]

[Materials]
  [elasticity_tensor]
    type = FeCrAlElasticityTensor
    temperature = temperature
    fecral_material_type = PM2000
  []
  [stress]
    type = ComputeFiniteStrainElasticStress
  []
  [clad_density]
    type = StrainAdjustedDensity
    strain_free_density = 1.0
  []
  [thermal]
    type = FeCrAlThermal
    temperature = temperature
    material = PM2000
  []
[]

[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  = 100
  nl_max_its = 100
  nl_rel_tol = 1e-8
  nl_abs_tol = 1e-8
  l_tol      = 1e-8
  start_time = 0.0
  end_time   = 10
  dt         = 1
[]

[Postprocessors]
  [elastic_strain_rr]
    type = ElementAverageValue
    variable = elastic_strain_xx
  []
  [temp]
    type = AverageNodalVariableValue
    variable = temperature
    execute_on = 'initial timestep_end'
  []
[]

[Outputs]
  exodus = true
[]

(test/tests/solid_mechanics/fecral_elasticity_tensor/elasticity_Fecralloy.i)

# This test case is prepared to test the elastic moduli of Fecralloy FeCrAl alloy.
#
# The geometry represents a ring of Fecralloy with inner radius of 0.005 m,
# outer radius of 0.0055 m, and a height of 0.01 m.
#
# An axial load is applied to the top of the ring with a value of 100 MPa.
# The strain in the radial direction can be calculated using the Young's modulus
# and the Poisson's ratio from the applied 100.0 MPa pressure.
#
# epsilon_{rr} = nu * sigma_{axial} / E
#
# The Young's modulus (E) for Fecralloy is taken as 180 GPa.
#
# The Poisson's ratio of Fecralloy is calculated is taken as 0.3.
#
# The results of BISON compared to the analytical solution for an increase in
# temperature under constant pressure are summarized:
#
#                          BISON             Analytical
#    Temperature (K)    strain_rr (m/m)    strain_rr (m/m)
#        420              1.666667e-4        1.666667e-4
#        540              1.666667e-4        1.666667e-4
#        660              1.666667e-4        1.666667e-4
#        780              1.666667e-4        1.666667e-4
#        900              1.666667e-4        1.666667e-4
#       1020              1.666667e-4        1.666667e-4
#       1140              1.666667e-4        1.666667e-4
#       1260              1.666667e-4        1.666667e-4
#       1380              1.666667e-4        1.666667e-4
#       1500              1.666667e-4        1.666667e-4
#

[Mesh]
  coord_type = RZ
  [mesh]
    type = GeneratedMeshGenerator
    dim = 2
    xmin = 0.005
    xmax = 0.0055
    ymin = 0
    ymax = 0.01
  []
[]

[GlobalParams]
  displacements = 'disp_x disp_y'
[]

[Variables]
  [temperature]
    order = FIRST
    family = LAGRANGE
    initial_condition = 300.0
  []
[]

[Physics]
  [SolidMechanics]
    [QuasiStatic]
      [all]
        strain = FINITE
        block = 0
        add_variables = true
        generate_output = 'elastic_strain_xx'
      []
    []
  []
[]

[Functions]
  [temperature_function]
    type = PiecewiseLinear
    x = '0       10'
    y = '300   1500'
  []
  [pressure_function]
    type = PiecewiseLinear
    x = '0  10'
    y = '1   1'
  []
[]

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

[BCs]
  [top_pressure] # Simplified BC to apply pressure in only disp_y
    type = Pressure
    boundary = top
    variable = disp_y
    function = pressure_function
    factor = 100e6
  []

  [u_bottom_fix]
    type = DirichletBC
    variable = disp_y
    boundary = bottom
    value = 0.0
  []

  [temp_bc]
    type = FunctionDirichletBC
    variable = temperature
    boundary = 'bottom top left right'
    function = temperature_function
  []
[]

[Materials]
  [elasticity_tensor]
    type = FeCrAlElasticityTensor
    temperature = temperature
    fecral_material_type = FECRALLOY
  []
  [stress]
    type = ComputeFiniteStrainElasticStress
  []
  [clad_density]
    type = StrainAdjustedDensity
    strain_free_density = 1.0
  []
  [thermal]
    type = FeCrAlThermal
    temperature = temperature
    material = FECRALLOY
  []
[]

[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  = 100
  nl_max_its = 100
  nl_rel_tol = 1e-8
  nl_abs_tol = 1e-8
  l_tol      = 1e-8
  start_time = 0.0
  end_time   = 10
  dt         = 1
[]

[Postprocessors]
  [elastic_strain_rr]
    type = ElementAverageValue
    variable = elastic_strain_xx
  []
  [temp]
    type = AverageNodalVariableValue
    variable = temperature
    execute_on = 'initial timestep_end'
  []
[]

[Outputs]
  exodus = true
[]

(test/tests/solid_mechanics/fecral_elasticity_tensor/elasticity_APMT.i)

# This test case is prepared to test the elastic moduli of APMT FeCrAl alloy.
#
# The geometry represents a ring of APMT with inner radius of 0.005 m,
# outer radius of 0.0055 m, and a height of 0.01 m.
#
# An axial load is applied to the top of the ring with a value of 100 MPa.
# The strain in the radial direction can be calculated using the Young's modulus
# and the Poisson's ratio from the applied 100.0 MPa pressure.
#
# epsilon_{rr} = nu * sigma_{axial} / E
#
# The Young's modulus (E) for APMT is calculated via a piecewse linear function.
#
# The Poisson's ratio of APMT is calculated is taken as 0.3.
#
# The results of BISON compared to the analytical solution for an increase in
# temperature under constant pressure are summarized:
#
#                          BISON             Analytical
#    Temperature (K)    strain_rr (m/m)    strain_rr (m/m)
#        420              1.444687e-4        1.444687e-4
#        540              1.500103e-4        1.500103e-4
#        660              1.570794e-4        1.570794e-4
#        780              1.673033e-4        1.673033e-4
#        900              1.793025e-4        1.793025e-4
#       1020              1.931558e-4        1.931558e-4
#       1140              2.093291e-4        2.093291e-4
#       1260              2.284583e-4        2.284583e-4
#       1380              2.307692e-4        2.307692e-4
#       1500              2.307692e-4        2.307692e-4
#

[Mesh]
  coord_type = RZ
  [mesh]
    type = GeneratedMeshGenerator
    dim = 2
    xmin = 0.005
    xmax = 0.0055
    ymin = 0
    ymax = 0.01
  []
[]

[GlobalParams]
  displacements = 'disp_x disp_y'
[]

[Variables]
  [temperature]
    order = FIRST
    family = LAGRANGE
    initial_condition = 300.0
  []
[]

[Physics]
  [SolidMechanics]
    [QuasiStatic]
      [all]
        strain = FINITE
        block = 0
        add_variables = true
        generate_output = 'elastic_strain_xx'
      []
    []
  []
[]

[Functions]
  [temperature_function]
    type = PiecewiseLinear
    x = '0       10'
    y = '300   1500'
  []
  [pressure_function]
    type = PiecewiseLinear
    x = '0  10'
    y = '1   1'
  []
[]

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

[BCs]
  [top_pressure] # Simplified BC to apply pressure in only disp_y
    type = Pressure
    boundary = top
    variable = disp_y
    function = pressure_function
    factor = 100e6
  []

  [u_bottom_fix]
    type = DirichletBC
    variable = disp_y
    boundary = bottom
    value = 0.0
  []

  [temp_bc]
    type = FunctionDirichletBC
    variable = temperature
    boundary = 'bottom top left right'
    function = temperature_function
  []
[]

[Materials]
  [elasticity_tensor]
    type = FeCrAlElasticityTensor
    temperature = temperature
    fecral_material_type = APMT
  []
  [stress]
    type = ComputeFiniteStrainElasticStress
  []
  [clad_density]
    type = StrainAdjustedDensity
    strain_free_density = 1.0
  []
  [thermal]
    type = FeCrAlThermal
    temperature = temperature
    material = APMT
  []
[]

[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  = 100
  nl_max_its = 100
  nl_rel_tol = 1e-8
  nl_abs_tol = 1e-8
  l_tol      = 1e-8
  start_time = 0.0
  end_time   = 10
  dt         = 1
[]

[Postprocessors]
  [elastic_strain_rr]
    type = ElementAverageValue
    variable = elastic_strain_xx
  []
  [temp]
    type = AverageNodalVariableValue
    variable = temperature
    execute_on = 'initial timestep_end'
  []
[]

[Outputs]
  exodus = true
[]

(examples/accident_tolerant_fuel/uo2_fecral/uo2_fecral.i)


initial_fuel_density = 10431.0

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

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

[Mesh]
  coord_type = RZ
  displacements = 'disp_x disp_y'
  patch_size = 10 # For contact algorithm
  patch_update_strategy = auto
  partitioner = centroid
  centroid_partitioner_direction = y
  [mesh]
    type = FileMeshGenerator
    file = uo2_fecral_smeared.e
  []
[]

[Variables]
  [disp_x]
  []
  [disp_y]
  []
  [temp]
    initial_condition = 293.0
  []
[]

[UserObjects]
  [pin_geometry]
    type = FuelPinGeometry
    clad_inner_wall = 5
    clad_outer_wall = 2
    clad_top = 3
    clad_bottom = 1
    pellet_exteriors = 8
  []
[]

[AuxVariables]
  [fast_neutron_flux]
    block = clad
  []
  [fast_neutron_fluence]
    block = clad
  []
  [grain_radius]
    block = pellet_type_1
    initial_condition = 10e-6
  []
  [gap_cond]
    order = CONSTANT
    family = MONOMIAL
  []
  [coolant_htc]
    order = CONSTANT
    family = MONOMIAL
  []
  [oxide_thickness]
    order = CONSTANT
    family = MONOMIAL
  []
  [total_hoop_strain]
    order = CONSTANT
    family = MONOMIAL
  []
  [creep_strain_hoop]
    order = CONSTANT
    family = MONOMIAL
  []
  [hoop_stress]
    order = CONSTANT
    family = MONOMIAL
  []
  [creep_rate]
    order = CONSTANT
    family = MONOMIAL
  []
  [mass_gain]
    order = CONSTANT
    family = MONOMIAL
  []
[]

[Functions]
  [power_history]
    type = PiecewiseLinear
    x = '0 1e4 1e8'
    y = '0 2.5e4 2.5e4'
    scale_factor = 1
  []

  [axial_peaking_factors]
    type = ParsedFunction
    expression = 1
  []

  [pressure_ramp]
    type = PiecewiseLinear
    x = '-200 0 1e8'
    y = '6.537e-3 1 1'
    scale_factor = 15.5e6
  []

  [mass_flux_func]
    type = PiecewiseLinear
    x = '-200 0  1e8'
    y = '3800. 3800. 3800.'
  []
[]

[Physics/SolidMechanics/QuasiStatic]
  [pellets]
    block = pellet_type_1
    strain = FINITE
    eigenstrain_names = 'fuel_relocation_strain fuel_thermal_strain fuel_volumetric_strain'
    generate_output = 'vonmises_stress stress_xx stress_yy stress_zz'
    extra_vector_tags = 'ref'
    temperature = temp
  []
  [clad]
    block = clad
    strain = FINITE
    eigenstrain_names = 'clad_thermal_eigenstrain clad_volumetric_strain'
    generate_output = 'vonmises_stress stress_xx stress_yy stress_zz'
    extra_vector_tags = 'ref'
    temperature = temp
  []
[]

[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_type_1
    burnup_function = burnup
    extra_vector_tags = 'ref'
  []
[]

[Burnup]
  [burnup]
    block = pellet_type_1
    rod_ave_lin_pow = power_history
    axial_power_profile = axial_peaking_factors
    num_radial = 81
    num_axial = 11
    fuel_pin_geometry = pin_geometry
    fuel_volume_ratio = 1.0
    RPF = RPF
  []
[]

[AuxKernels]
  [fast_neutron_flux]
    type = FastNeutronFluxAux
    variable = fast_neutron_flux
    block = clad
    rod_ave_lin_pow = power_history
    axial_power_profile = axial_peaking_factors
    factor = 3e13
    execute_on = timestep_begin
  []
  [fast_neutron_fluence]
    type = FastNeutronFluenceAux
    variable = fast_neutron_fluence
    block = clad
    fast_neutron_flux = fast_neutron_flux
    execute_on = timestep_begin
  []
  [grain_radius]
    type = GrainRadiusAux
    block = pellet_type_1
    variable = grain_radius
    temperature = temp
    execute_on = linear
  []
  [hoop_stress]
    type = RankTwoScalarAux
    rank_two_tensor = stress
    variable = hoop_stress
    scalar_type = HoopStress
    execute_on = timestep_end
  []
  [total_hoop_strain]
    type = RankTwoScalarAux
    rank_two_tensor = total_strain
    variable = total_hoop_strain
    scalar_type = HoopStress
    execute_on = timestep_end
  []
  [creep_strain_hoop]
    type = RankTwoScalarAux
    rank_two_tensor = creep_strain
    variable = creep_strain_hoop
    scalar_type = HoopStress
    block = clad
    execute_on = timestep_end
  []
  [conductance]
    type = MaterialRealAux
    property = gap_conductance
    variable = gap_cond
    boundary = 10
  []
  [coolant_htc]
    type = MaterialRealAux
    property = coolant_channel_htc
    variable = coolant_htc
    boundary = 2
  []
  [creep_rate]
    type = MaterialRealAux
    variable = creep_rate
    property = creep_rate
    execute_on = timestep_end
    block = clad
  []
  [oxide]
    type = MaterialRealAux
    variable = oxide_thickness
    property = scale_thickness
    boundary = 2
  []
  [mass_gain]
    type = MaterialRealAux
    variable = mass_gain
    property = oxide_mass_gain
    boundary = 2
  []
[]

[Contact]
  [pellet_clad_mechanical]
    primary = 5
    secondary = 10
    normal_smoothing_distance = 0.1
    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 = 1020
    value = 0.0
  []

  [Pressure]
    [coolantPressure]
      boundary = '1 2 3'
      function = pressure_ramp
    []
  []

  [PlenumPressure]
    [plenumPressure]
      boundary = 9
      initial_pressure = 2.0e6
      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
    []
  []

[]

[CoolantChannel]
  [convective_clad_surface]
    boundary = '1 2 3'
    variable = temp
    inlet_temperature = 580 # K
    inlet_pressure = pressure_ramp # Pa
    inlet_massflux = mass_flux_func # kg/m^2-sec
    rod_diameter = 9.5e-3 # m
    rod_pitch = 1.26e-2 # m
    linear_heat_rate = power_history
    axial_power_profile = axial_peaking_factors
    oxide_thickness = oxide_thickness
  []
[]

[Materials]
  [fuel_thermal]
    type = UO2Thermal
    block = pellet_type_1
    thermal_conductivity_model = NFIR
    temperature = temp
    burnup_function = burnup
  []
  [fuel_elasticity_tensor]
    type = ComputeIsotropicElasticityTensor
    block = pellet_type_1
    youngs_modulus = 2.0e11
    poissons_ratio = 0.345
  []
  [elastic_stress]
    type = ComputeSmearedCrackingStress
    block = pellet_type_1
    cracking_stress = 1.68e8
    inelastic_models = 'fuel_creep'
    softening_models = exponential_softening
    shear_retention_factor = 0.1
    max_stress_correction = 0
    cracked_elasticity_type = DIAGONAL
    output_properties = crack_damage
    outputs = exodus
  []
  [exponential_softening]
    type = ExponentialSoftening
  []
  [fuel_creep]
    type = UO2CreepUpdate
    block = pellet_type_1
    burnup_function = burnup
    temperature = temp
  []
  [fuel_relocation]
    type = UO2RelocationEigenstrain
    block = pellet_type_1
    burnup_function = burnup
    rod_ave_lin_pow = power_history
    axial_power_profile = axial_peaking_factors
    fuel_pin_geometry = 'pin_geometry'
    relocation_activation1 = 5000
    relocation_model = ESCORE_modified
    eigenstrain_name = fuel_relocation_strain
  []
  [fuel_thermal_expansion]
    type = UO2ThermalExpansionMATPROEigenstrain
    block = pellet_type_1
    temperature = temp
    stress_free_temperature = 295.0
    eigenstrain_name = fuel_thermal_strain
  []
  [fuel_volumetric_swelling]
    type = UO2VolumetricSwellingEigenstrain
    gas_swelling_model_type = SIFGRS
    block = pellet_type_1
    temperature = temp
    burnup_function = burnup
    initial_fuel_density = 10431.0
    eigenstrain_name = fuel_volumetric_strain
  []
  [clad_thermal]
    type = FeCrAlThermal
    material = C35M
    block = clad
    temperature = temp
  []
  [clad_elasticity_tensor] # isotropic elasticity tensor for Zry cladding
    type = FeCrAlElasticityTensor
    temperature = temp
    fecral_material_type = C35M
    block = clad
  []
  [clad_stress] # stress update class to govern the return mapping algorithm for creep
    type = ComputeMultipleInelasticStress
    tangent_operator = elastic
    inelastic_models = 'clad_creep clad_plasticity'
    block = clad
  []
  [clad_creep]
    type = FeCrAlCreepUpdate
    block = clad
    temperature = temp
    fecral_material_type = C35M
    fast_neutron_flux = fast_neutron_flux
    model_irradiation_creep = true
    model_thermal_creep = true
    max_inelastic_increment = 1e-4
  []
  [thermal_expansion]
    type = FeCrAlThermalExpansionEigenstrain
    block = clad
    temperature = temp
    fecral_material_type = C35M
    stress_free_temperature = 295.0
    eigenstrain_name = clad_thermal_eigenstrain
  []
  [irradiation_swelling]
    type = FeCrAlVolumetricSwellingEigenstrain
    block = clad
    temperature = temp
    fast_neutron_fluence = fast_neutron_fluence
    eigenstrain_name = clad_volumetric_strain
  []
  [clad_plasticity]
    type = FeCrAlPlasticityUpdate
    block = clad
    hardening_constant = 2.5e9
    temperature = temp
    yield_stress = 500.0
  []
  [fission_gas_release]
    type = UO2Sifgrs
    block = pellet_type_1
    temperature = temp
    burnup_function = burnup
    grain_radius = grain_radius
    gbs_model = true
  []
  [clad_density]
    type = StrainAdjustedDensity
    block = clad
    strain_free_density = 7250.0
  []
  [fuel_density]
    type = StrainAdjustedDensity
    block = pellet_type_1
    strain_free_density = ${initial_fuel_density}
  []

  [failure_criterion]
    type = FeCrAlCladdingFailure
    boundary = '2 5'
    hoop_stress = hoop_stress
    failure_criterion = UTS
    temperature = temp
  []
  [oxidation]
    type = FeCrAlOxidation
    reactor_type = PWR
    boundary = 2
  []
[]

[Dampers]
  [limitT]
    type = BoundingValueNodalDamper
    max_value = 3200.0
    min_value = 293.0
    variable = temp
  []
  [limitX]
    type = MaxIncrement
    max_increment = 1e-5
    variable = disp_x
  []
[]

[Executioner]
  type = Transient
  solve_type = 'PJFNK'

  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 = 8e-3

  nl_max_its = 25
  nl_rel_tol = 1e-5
  nl_abs_tol = 1e-10

  start_time = -200
  n_startup_steps = 1
  end_time = 1e8

  dtmax = 1e5
  dtmin = 1

  [TimeStepper]
    type = IterationAdaptiveDT
    dt = 2.0e2
    force_step_every_function_point = true
    timestep_limiting_function = power_history
    max_function_change = 5e5
    optimal_iterations = 10
    iteration_window = 2
    linear_iteration_ratio = 100
    growth_factor = 2.0
    timestep_limiting_postprocessor = material_timestep
  []

  [Quadrature]
    order = FIFTH
    side_order = SEVENTH
  []
[]

[Postprocessors]
  [ave_temp_interior]
    type = SideAverageValue
    boundary = 9
    variable = temp
    execute_on = 'initial linear'
  []

  [avg_clad_temp]
    type = SideAverageValue
    boundary = 7
    variable = temp
  []

  [fis_gas_produced]
    type = ElementIntegralFisGasGeneratedSifgrs
    block = pellet_type_1
  []

  [fis_gas_released]
    type = ElementIntegralFisGasReleasedSifgrs
    block = pellet_type_1
  []
  [fis_gas_grain]
    type = ElementIntegralFisGasGrainSifgrs
    block = pellet_type_1
  []
  [fis_gas_boundary]
    type = ElementIntegralFisGasBoundarySifgrs
    block = pellet_type_1
  []

  [gas_volume]
    type = InternalVolume
    boundary = 9
    execute_on = 'initial linear'
  []

  [_dt]
    type = TimestepSize
  []

  [num_lin_it]
    type = NumLinearIterations
  []
  [num_nonlin_it]
    type = NumNonlinearIterations
  []
  [tot_lin_it]
    type = CumulativeValuePostprocessor
    postprocessor = num_lin_it
  []
  [tot_nonlin_it]
    type = CumulativeValuePostprocessor
    postprocessor = num_nonlin_it
  []
  [alive_time]
    type = PerfGraphData
    section_name = Root
    data_type = TOTAL
  []

  [rod_total_power]
    type = ElementIntegralPower
    variable = temp
    burnup_function = burnup
    block = pellet_type_1
  []

  [alhr_input]
    type = FunctionValuePostprocessor
    function = power_history
  []
  [average_burnup]
    type = ElementAverageValue
    block = pellet_type_1
    variable = burnup
  []
  [oxide_thickness]
    type = ElementExtremeValue
    block = clad
    variable = oxide_thickness
  []
  [mass_gain]
    type = ElementExtremeValue
    block = clad
    variable = mass_gain
  []
  [fis_gas_percent]
    type = FGRPercent
    fission_gas_released = fis_gas_released
    fission_gas_generated = fis_gas_produced
  []
  [material_timestep]
    type = MaterialTimeStepPostprocessor
    block = clad
  []
[]

[Outputs]
  perf_graph = true
  time_step_interval =  1
  exodus = true
  csv = true
  print_linear_residuals = true
  color = false

  [console]
    type = Console
    max_rows = 25
  []
[]

(test/tests/solid_mechanics/fecral_elasticity_tensor/elasticity_C06M.i)

# This test case is prepared to test the elastic moduli of C06M FeCrAl alloy.
#
# The geometry represents a ring of C06M with inner radius of 0.005 m,
# outer radius of 0.0055 m, and a height of 0.01 m.
#
# An axial load is applied to the top of the ring with a value of 100 MPa.
# The strain in the radial direction can be calculated using the Young's modulus
# and the Poisson's ratio from the applied 100.0 MPa pressure.
#
# epsilon_{rr} = nu * sigma_{axial} / E
#
# The Young's modulus (E) for C06M is calculated via:
#
# E = -5.46e-5 * (T - 273.15)^2 - 3.85e-2 * (T - 273.15) + 1.99e2) * 1e9
#
# where E is in Pa and T is the temperature in K.
#
# The Poisson's ratio of C06M is calculated via:
#
# nu = 4.46e-5 * (t - 273.15) + 0.27
#
# The results of BISON compared to the analytical solution for an increase in
# temperature under constant pressure are summarized:
#
#                          BISON             Analytical
#    Temperature (K)    strain_rr (m/m)    strain_rr (m/m)
#        420              1.439097e-4        1.439097e-4
#        540              1.525125e-4        1.525125e-4
#        660              1.632723e-4        1.632723e-4
#        780              1.768440e-4        1.768440e-4
#        900              1.942209e-4        1.942209e-4
#       1020              2.169732e-4        2.169732e-4
#       1140              2.477254e-4        2.477254e-4
#       1260              2.912040e-4        2.912040e-4
#       1380              3.568533e-4        3.568533e-4
#       1500              4.666521e-4        4.666521e-4
#
# The answers are the same for the FeCrAlloy alloys C35M and C36M.

[Mesh]
  coord_type = RZ
  [mesh]
    type = GeneratedMeshGenerator
    dim = 2
    xmin = 0.005
    xmax = 0.0055
    ymin = 0
    ymax = 0.01
  []
[]

[GlobalParams]
  displacements = 'disp_x disp_y'
[]

[Variables]
  [temperature]
    order = FIRST
    family = LAGRANGE
    initial_condition = 300.0
  []
[]

[Physics]
  [SolidMechanics]
    [QuasiStatic]
      [all]
        strain = FINITE
        block = 0
        add_variables = true
        generate_output = 'elastic_strain_xx'
      []
    []
  []
[]

[Functions]
  [temperature_function]
    type = PiecewiseLinear
    x = '0       10'
    y = '300   1500'
  []
  [pressure_function]
    type = PiecewiseLinear
    x = '0  10'
    y = '1   1'
  []
[]

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

[BCs]
  [top_pressure] # Simplified BC to apply pressure in only disp_y
    type = Pressure
    boundary = top
    variable = disp_y
    function = pressure_function
    factor = 100e6
  []

  [u_bottom_fix]
    type = DirichletBC
    variable = disp_y
    boundary = bottom
    value = 0.0
  []

  [temp_bc]
    type = FunctionDirichletBC
    variable = temperature
    boundary = 'bottom top left right'
    function = temperature_function
  []
[]

[Materials]
  [elasticity_tensor]
    type = FeCrAlElasticityTensor
    temperature = temperature
    fecral_material_type = C06M
  []
  [stress]
    type = ComputeFiniteStrainElasticStress
  []
  [clad_density]
    type = StrainAdjustedDensity
    strain_free_density = 1.0
  []
  [thermal]
    type = FeCrAlThermal
    temperature = temperature
    material = C06M
  []
[]

[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  = 100
  nl_max_its = 100
  nl_rel_tol = 1e-8
  nl_abs_tol = 1e-8
  l_tol      = 1e-8
  start_time = 0.0
  end_time   = 10
  dt         = 1
[]

[Postprocessors]
  [elastic_strain_rr]
    type = ElementAverageValue
    variable = elastic_strain_xx
  []
  [temp]
    type = AverageNodalVariableValue
    variable = temperature
    execute_on = 'initial timestep_end'
  []
[]

[Outputs]
  exodus = true
[]

(test/tests/solid_mechanics/fecral_elasticity_tensor/elasticity_MA956.i)

# This test case is prepared to test the elastic moduli of MA956 FeCrAl alloy.
#
# The geometry represents a ring of MA956 with inner radius of 0.005 m,
# outer radius of 0.0055 m, and a height of 0.01 m.
#
# An axial load is applied to the top of the ring with a value of 100 MPa.
# The strain in the radial direction can be calculated using the Young's modulus
# and the Poisson's ratio from the applied 100.0 MPa pressure.
#
# epsilon_{rr} = nu * sigma_{axial} / E
#
# The Young's modulus (E) for MA956 is calculated via a piecewse linear function.
#
# The Poisson's ratio of MA956 is calculated is taken as 0.3.
#
# The results of BISON compared to the analytical solution for an increase in
# temperature under constant pressure are summarized:
#
#                          BISON             Analytical
#    Temperature (K)    strain_rr (m/m)    strain_rr (m/m)
#        420              1.794104e-4        1.794104e-4
#        540              1.900288e-4        1.900288e-4
#        660              2.057089e-4        2.057089e-4
#        780              2.202522e-4        2.202522e-4
#        900              2.477721e-4        2.477721e-4
#       1020              2.703530e-4        2.703530e-4
#       1140              2.999621e-4        2.999621e-4
#       1260              3.807386e-4        3.807386e-4
#       1380              4.126547e-4        4.126547e-4
#       1500              4.126547e-4        4.126547e-4
#

[Mesh]
  coord_type = RZ
  [mesh]
    type = GeneratedMeshGenerator
    dim = 2
    xmin = 0.005
    xmax = 0.0055
    ymin = 0
    ymax = 0.01
  []
[]

[GlobalParams]
  displacements = 'disp_x disp_y'
[]

[Variables]
  [temperature]
    order = FIRST
    family = LAGRANGE
    initial_condition = 300.0
  []
[]

[Physics]
  [SolidMechanics]
    [QuasiStatic]
      [all]
        strain = FINITE
        block = 0
        add_variables = true
        generate_output = 'elastic_strain_xx'
      []
    []
  []
[]

[Functions]
  [temperature_function]
    type = PiecewiseLinear
    x = '0       10'
    y = '300   1500'
  []
  [pressure_function]
    type = PiecewiseLinear
    x = '0  10'
    y = '1   1'
  []
[]

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

[BCs]
  [top_pressure] # Simplified BC to apply pressure in only disp_y
    type = Pressure
    boundary = top
    variable = disp_y
    function = pressure_function
    factor = 100e6
  []

  [u_bottom_fix]
    type = DirichletBC
    variable = disp_y
    boundary = bottom
    value = 0.0
  []

  [temp_bc]
    type = FunctionDirichletBC
    variable = temperature
    boundary = 'bottom top left right'
    function = temperature_function
  []
[]

[Materials]
  [elasticity_tensor]
    type = FeCrAlElasticityTensor
    temperature = temperature
    fecral_material_type = MA956
  []
  [stress]
    type = ComputeFiniteStrainElasticStress
  []
  [clad_density]
    type = StrainAdjustedDensity
    strain_free_density = 1.0
  []
  [thermal]
    type = FeCrAlThermal
    temperature = temperature
    material = MA956
  []
[]

[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  = 100
  nl_max_its = 100
  nl_rel_tol = 1e-8
  nl_abs_tol = 1e-8
  l_tol      = 1e-8
  start_time = 0.0
  end_time   = 10
  dt         = 1
[]

[Postprocessors]
  [elastic_strain_rr]
    type = ElementAverageValue
    variable = elastic_strain_xx
  []
  [temp]
    type = AverageNodalVariableValue
    variable = temperature
    execute_on = 'initial timestep_end'
  []
[]

[Outputs]
  exodus = true
[]

Description
Range of Applicability and Uncertainty
Example Input Syntax
Input Parameters
Input Files
References