HydrogenPickup

Hydrogen flux BC for various zirconium alloys. Units for H concentration is wt.ppm.

Description

The HydrogenPickup BC is used to model the flux of hydrogen into the clad that is caused by oxide growth. The flux is approximated as a constant fraction of the hydrogen liberated by oxide growth at the interface between the coolant water and the clad.

note:Coupling to Oxide Thickness Variable

This BC must be coupled to an AuxVariable that gives the thickness of the oxide over time when the Zircaloy material type is Zircaloy-4, M5, Zirlo, Optimized Zirlo, or Zircaloy-2 with the model based on linear power. For this to work properly, the AuxKernel must be set to execute on linear.

note:Coupling to Axial Averaged Burnup Variable

This BC must be coupled to an AuxVariable that holds the local axial averaged burnup in the fuel at the corresponding cladding height over time when the Zircaloy material type is Zircaloy-2 and the model based on burnup is being used. For this to work properly, the AuxKernel must be set to execute on linear.

The waterside corrosion reaction of the zirconium alloy cladding generates hydrogen at the water-oxide interface: $var element = document.getElementById("moose-equation-ac404891-8491-496d-9733-010d56fb25b2");katex.render("\\mathrm{Zr} + 2\\mathrm{H}_2\\mathrm{O} \\rightarrow \\mathrm{ZrO}_2 + 2\\mathrm{H}_2", element, {displayMode:true,throwOnError:false});$ Some of the hydrogen migrates through the oxide layer and dissolves in the cladding, where it is slightly soluble. The fraction of the hydrogen produced by the corrosion reaction that ends up in the cladding is termed the hydrogen pickup fraction. Although the mechanism for hydrogen migration through the oxide layer is a topic of active research, the pickup fraction has been well characterized for various zirconium alloys. For example, Couet et al. (2014) found the instantaneous hydrogen pickup fraction to be a strong function of alloying elements and to be a complex function of the exposure time. To enable a wide variety of modeling approaches, there are several available options to specifying the pickup fraction in BISON.

Hydrogen pick up fraction models for Zircaloy-4, M5, Zirlo, and Optimized Zirlo

For Zircaloy-4, M5, Zirlo, and Optimized Zirlo, the average total concentration of hydrogen $var element = document.getElementById("moose-equation-b696caba-d836-4357-9392-6f6e24413dc5");katex.render("C_H", element, {displayMode:false,throwOnError:false});$ (including dissolved hydrogen and hydrogen as ZrH $var element = document.getElementById("moose-equation-0429ef37-09f5-4aa0-b135-cbc6388d9893");katex.render("_x", element, {displayMode:false,throwOnError:false});$ ) in the cladding is roughly proportional to the thickness of the oxide layer (Courty et al., 2014): $(1)var element = document.getElementById("moose-equation-c95ee8e5-9de2-4cda-9cfc-6d2ba5ad0a88");katex.render(" C_H = 1.0 \\times 10^6 \\frac{4f}{R_{PB}} \\frac{V_{oxide}}{V_{tot}-V_{oxide}} \\frac{M_H}{M_{Zr}}", element, {displayMode:true,throwOnError:false});$ where $var element = document.getElementById("moose-equation-d33690a0-bd16-4ec5-9a4c-e2c54ac9f53f");katex.render("f", element, {displayMode:false,throwOnError:false});$ is the hydrogen pickup fraction (unitless), $var element = document.getElementById("moose-equation-f36b82c2-e11b-4e28-a04b-2ac5b28000f1");katex.render("R_{PB}", element, {displayMode:false,throwOnError:false});$ is the Pilling-Bedworth ratio for ZrO $var element = document.getElementById("moose-equation-7994797b-b9b8-4a3e-899d-ccbf6c7b9b7a");katex.render("_2", element, {displayMode:false,throwOnError:false});$ and is equal to 1.56 (unitless), $var element = document.getElementById("moose-equation-04eb84e4-0b28-47d1-a3af-6cf5101a7526");katex.render("V_{oxide}", element, {displayMode:false,throwOnError:false});$ is the volume of the oxide thickness (m $var element = document.getElementById("moose-equation-d9aee349-586f-472a-b16b-16e71851d15b");katex.render("^3", element, {displayMode:false,throwOnError:false});$ ), $var element = document.getElementById("moose-equation-f1d0e9f4-b86b-4ba9-bed2-bb76a349e62b");katex.render("V_{tot}", element, {displayMode:false,throwOnError:false});$ is the initial volume of the cladding (m $var element = document.getElementById("moose-equation-b685e045-f8b4-48d8-b14d-a836c0b45ed1");katex.render("^3", element, {displayMode:false,throwOnError:false});$ ), $var element = document.getElementById("moose-equation-e8691a96-6cd8-4a72-bfc9-bddc72bd4aaf");katex.render("M_H", element, {displayMode:false,throwOnError:false});$ is the molecular weight of hydrogen (g/mol), and $var element = document.getElementById("moose-equation-dabd415b-fec9-40ad-bff4-9693851bf57f");katex.render("M_{Zr}", element, {displayMode:false,throwOnError:false});$ is the molecular weight of zirconium (g/mol). The factor of $var element = document.getElementById("moose-equation-dee989c3-1b37-4bed-9c0c-fac78d576bd6");katex.render("1.0 \\times 10^6", element, {displayMode:false,throwOnError:false});$ is to convert the concentration to wt. parts per million (ppm). For the case of a plate geometry the following holds: $var element = document.getElementById("moose-equation-7da848e7-a393-4827-a4cb-7292d8e4c8e0");katex.render(" \\frac{V_{oxide}}{V_{tot}-V_{oxide}} = \\frac{\\delta}{L - \\delta}", element, {displayMode:true,throwOnError:false});$ where $var element = document.getElementById("moose-equation-2cd729c3-b2a0-48e4-9119-0331ea0a58f9");katex.render("\\delta", element, {displayMode:false,throwOnError:false});$ is the the oxide thickness (m), and $var element = document.getElementById("moose-equation-724fba08-fd89-4b84-b602-17c77e142a60");katex.render("L", element, {displayMode:false,throwOnError:false});$ is the initial cladding thickness (m). For cylindrical geometries, like nuclear fuel cladding, the following is true: $var element = document.getElementById("moose-equation-a7cab0c7-6260-4141-85f3-a430a3dc86f5");katex.render(" \\frac{V_{oxide}}{V_{tot}-V_{oxide}} = \\frac{2 r_{out} \\delta - \\delta^2}{2 r_{out} \\left(L-\\delta\\right) - L^2 + \\delta^2}", element, {displayMode:true,throwOnError:false});$ where $var element = document.getElementById("moose-equation-3f1811f4-a7f7-4b05-8c6a-26b72dd95e5a");katex.render("V_{oxide}=\\pi H \\left(r_{out}^2-\\left(r_{out}-\\delta\\right)^2\\right)", element, {displayMode:false,throwOnError:false});$ and $var element = document.getElementById("moose-equation-dec4bc97-34f1-4152-b374-dc35af5d0e52");katex.render("V_{tot}=\\pi H \\left(r_{out}^2-\\left(r_{out}-L\\right)^2\\right)", element, {displayMode:false,throwOnError:false});$ . This convention corresponds the volume of the oxide consuming some of the initial cladding volume.

To introduce the appropriate amount of hydrogen in the simulated domain, the HydrogenPickup BC introduces a flux on the desired surface of the domain. The flux is defined as $(2)var element = document.getElementById("moose-equation-1fcf1e74-2b97-4f50-adf1-65547cc7f8ac");katex.render(" J = \\frac{V}{S} \\frac{dC_H}{dt},", element, {displayMode:true,throwOnError:false});$ where $var element = document.getElementById("moose-equation-efe6f272-6ce5-4351-861e-63d209e55f78");katex.render("V", element, {displayMode:false,throwOnError:false});$ is the volume of the domain, and $var element = document.getElementById("moose-equation-037a5d7e-c0f8-4e0f-8e9a-99d8cf5e3acb");katex.render("S", element, {displayMode:false,throwOnError:false});$ the oxidized surface of the cladding. For a plate geometry, $var element = document.getElementById("moose-equation-690d9eaf-7fb5-43f2-b6c7-508761c78480");katex.render("\\frac{V}{S}=L", element, {displayMode:false,throwOnError:false});$ with $var element = document.getElementById("moose-equation-f85a79cd-2697-4c97-939b-5b3de9c47e4f");katex.render("L", element, {displayMode:false,throwOnError:false});$ the thickness of the plate. However, for a cylindrical cladding with an outer radius $var element = document.getElementById("moose-equation-6cacedec-4902-4b77-9cbc-2acad3aa763d");katex.render("r_{out}", element, {displayMode:false,throwOnError:false});$ and a thickness $var element = document.getElementById("moose-equation-24e966ad-c4b3-41da-8150-67da0c8614e8");katex.render("L", element, {displayMode:false,throwOnError:false});$ , assuming that oxidation happens only on the outer surface of the cladding, then $(3)var element = document.getElementById("moose-equation-18831881-22a8-4121-97f8-342b755a4eaf");katex.render(" \\frac{V}{S} = L \\left(1.0 - \\frac{L}{2r_{out}}\\right).", element, {displayMode:true,throwOnError:false});$

The rate of change of total hydrogen concentration in plate geometries was derived by (Courty et al., 2014), $(4)var element = document.getElementById("moose-equation-09843db5-d0d3-4ff4-b9cf-9d34f52d8f03");katex.render(" \\frac{dC_H}{dt} = \\frac{4f}{R_{PB}} \\frac{M_H}{M_{Zr}} \\frac{L}{(L-\\delta)^2} \\frac{d \\delta}{dt},", element, {displayMode:true,throwOnError:false});$ which results in a hydrogen flux of $(5)var element = document.getElementById("moose-equation-6a0168b4-6f05-4c70-9877-dcacc2825f05");katex.render(" J = 1.0 \\times 10^6 \\frac{4f}{R_{PB}} \\frac{M_H}{M_{Zr}} \\frac{L^2}{(L-\\delta)^2} \\frac{d \\delta}{dt}", element, {displayMode:true,throwOnError:false});$ in wt.ppm-m / s. For a cylindrical geometry the hydrogen flux is given by $(6)var element = document.getElementById("moose-equation-d9bd33ce-7459-4797-82a4-abcb57fc5d4d");katex.render(" J = 1.0 \\times 10^6 \\frac{4f}{R_{PB}} L \\left(2.0 - \\frac{L}{r_{out}}\\right) \\frac{M_H}{M_{Zr}} \\left(-\\frac{2L\\left(L-2r_{out}\\right)\\left(r_{out} - \\delta\\right)}{\\left(L-\\delta\\right)^2\\left(-2r_{out}+\\delta+L\\right)^2}\\right) \\frac{d \\delta}{dt}", element, {displayMode:true,throwOnError:false});$

To minimize timestep dependence of the oxide thickness on the calculated flux of hydrogen, $var element = document.getElementById("moose-equation-01577cf7-5b59-4ff6-9f17-bfc3811397b5");katex.render("\\delta_{avg}", element, {displayMode:false,throwOnError:false});$ is used in the preceding equations in BISON. Where: $var element = document.getElementById("moose-equation-6c1d7ca6-2a4c-4aa4-81fd-b15aa737b303");katex.render(" \\delta_{avg} = \\frac{\\delta + \\delta_{old}}{2.0}", element, {displayMode:true,throwOnError:false});$ and $var element = document.getElementById("moose-equation-73a8082e-24a6-450c-8270-4c0f5fdf8f23");katex.render("d \\delta = \\delta - \\delta_{old}", element, {displayMode:false,throwOnError:false});$ .

The first approach to specifying the pickup fraction is through a user specified constant value. The second allows the user to select the zirconium alloy of interest and the hydrogen pickup fraction is then computed based upon the recommendations of Clifford (2020). For

Zircaloy-4(zircaloy_material_type = STRESS_RELIEF_ANNEALED) $var element = document.getElementById("moose-equation-df141620-49c9-422e-996f-633e20368fe8");katex.render("f=0.2", element, {displayMode:false,throwOnError:false});$
ZIRLO® (zircaloy_material_type = ZIRLO) $var element = document.getElementById("moose-equation-1628b18a-ec0b-4db8-a0a0-93b940bc798e");katex.render("f=0.25", element, {displayMode:false,throwOnError:false});$
Optimized ZIRLO™ (zircaloy_material_type = PARTIAL_RECRYSTALLIZATION_ANNEALED) $var element = document.getElementById("moose-equation-322fd45b-9018-410d-9359-d798201e641c");katex.render("f=0.25", element, {displayMode:false,throwOnError:false});$
M5® (zircaloy_material_type = M5) $var element = document.getElementById("moose-equation-654f23a5-b2da-43c6-a92c-6431c24a7b66");katex.render("f=0.15", element, {displayMode:false,throwOnError:false});$ .

Hydrogen pick up fraction models for Zircaloy-2

For Zircaloy-2 (zircaloy_material_type = RECRYSTALLIZATION_ANNEALED), primarily used in a boiling water reactor (BWR), there are two distinct options.

Hydrogen pick up model for Zircaloy-2 based on fuel burnup

In one option (BURNUP_MODEL), the hydrogen concentration is a function of local axial burnup of the fuel adjacent to the cladding (Clifford, 2020): $var element = document.getElementById("moose-equation-613e3d4c-e97e-4f0c-a92f-5820259260cf");katex.render(" C_H = \\begin{cases} 1.4 \\left(47.8 \\exp\\left(\\frac{-1.3}{1.0 + Bu}\\right)+0.316Bu\\right), & \\text{for } Bu \\leq 50~\\text{GWd/tU}\\\\ 1.4 \\left(28.9 + \\exp\\left(0.117\\left(Bu-20\\right)\\right) \\right), & \\text{for } Bu > 50~\\text{GWd/tU} \\end{cases}", element, {displayMode:true,throwOnError:false});$ where $var element = document.getElementById("moose-equation-cd03783f-7bf3-466b-8945-41a4cd1d1f23");katex.render("Bu", element, {displayMode:false,throwOnError:false});$ is the local axial burnup (GWd/tU). Since $var element = document.getElementById("moose-equation-723fa3a5-db84-4a97-b62f-0129762fffc3");katex.render("C_H", element, {displayMode:false,throwOnError:false});$ is not a function of oxide thickness in the case of Zircaloy-2, the flux is obtained by taking the time derivative of $var element = document.getElementById("moose-equation-59d8f3a0-dd6a-4542-b940-b0f799b6b4b5");katex.render("C_H", element, {displayMode:false,throwOnError:false});$ where the local axial burnup is time dependent and using the same $var element = document.getElementById("moose-equation-5f66d759-3ba3-47cc-ba7a-a41c9381b613");katex.render("\\frac{V}{S}", element, {displayMode:false,throwOnError:false});$ ratios for plate and cylindrical geometries as the for the other alloys. No conversion factor is necessary for this correlation as the correlation returns a value in ppm. Thus, for the Zircaloy-2 hydrogen model the flux boundary condition for cylindrical geometries is given by: $var element = document.getElementById("moose-equation-71509510-c9e0-49c3-9d0c-57543ea00d10");katex.render(" J = \\begin{cases} L \\left(1.0 - \\frac{L}{2r_{out}}\\right) \\left(\\frac{86.996}{\\left(1+Bu\\right)^2} \\exp\\left(\\frac{-1.3}{1+Bu}\\right)+0.4424\\right)\\frac{dBu}{dt}, & \\text{for } Bu \\leq 50~\\text{GWd/tU}\\\\ L \\left(1.0 - \\frac{L}{2r_{out}}\\right) 0.1638 \\exp\\left(0.117\\left(Bu-20.0\\right)\\right)\\frac{dBu}{dt}, & \\text{for } Bu > 50~\\text{GWd/tU} \\end{cases}", element, {displayMode:true,throwOnError:false});$ and for plate specimens: $var element = document.getElementById("moose-equation-0bebf371-596b-4838-8e95-99f8a4807ea7");katex.render(" J = \\begin{cases} L \\left(\\frac{86.996}{\\left(1+Bu\\right)^2} \\exp\\left(\\frac{-1.3}{1+Bu}\\right)+0.4424\\right)\\frac{dBu}{dt}, & \\text{for } Bu \\leq 50~\\text{GWd/tU}\\\\ 0.1638 L \\exp\\left(0.117\\left(Bu-20.0\\right)\\right)\\frac{dBu}{dt}, & \\text{for } Bu > 50~\\text{GWd/tU} \\end{cases}", element, {displayMode:true,throwOnError:false});$

As before, with the oxide thickness, the average burnup is used to minimize timestep dependence. Therefore: $var element = document.getElementById("moose-equation-0d96ecff-b5cf-47b1-809a-8fdb2af9d3ca");katex.render(" Bu_{avg} = \\frac{Bu + Bu_{old}}{2.0}", element, {displayMode:true,throwOnError:false});$ and $var element = document.getElementById("moose-equation-f697e2a9-59b8-4afb-a74e-6228d26a6778");katex.render("d Bu = Bu - Bu_{old}", element, {displayMode:false,throwOnError:false});$ .

note:Initial hydrogen content for Zircaloy-2 in a BWR environment

Because the equations in Ref. (Clifford, 2020) states that the amount of hydrogen is equal to $var element = document.getElementById("moose-equation-a51c228b-3c94-4c36-9d3d-dacfbcab729c");katex.render("C_H = 18.238", element, {displayMode:false,throwOnError:false});$ wt.ppm at 0 burnup (and not 0 wt.ppm), 18.238 wt.ppm of hydrogen should be initially added to the cladding for the simulation to match the equation provided by Ref. (Clifford, 2020).

Hydrogen pick up model for Zircaloy-2 based on linear heat rate

For Zircaloy-2 in BWR, it is also possible to use a model based on LHR, or linear power (Nissen et al., 2014). When using option LINEAR_POWER_MODEL the hydrogen flux is computed just as in Eq. (5) and Eq. (6) depending on the sample geometry, but with the hydrogen pick up fraction defined as $(7)var element = document.getElementById("moose-equation-9f748a91-f3fc-4980-8b71-50c0e06d5775");katex.render(" f = f_{Zr2}^0 \\times \\begin{cases} \\alpha, & \\text{for } \\text{LHR} \\leq \\text{LHR}_{lim}~\\text{W/m}\\\\ 1, & \\text{for } \\text{LHR} > \\text{LHR}_{lim}~\\text{W/m}, \\end{cases}", element, {displayMode:true,throwOnError:false});$ where $var element = document.getElementById("moose-equation-eb826fa2-c0e6-4557-accf-10d5ec72d690");katex.render("f_{Zr2}^0", element, {displayMode:false,throwOnError:false});$ is the hydrogen pick up fraction for Zircaloy-2 in BWR at high LHR, $var element = document.getElementById("moose-equation-d3f0ba85-7591-4ce7-a4b6-9d8045018018");katex.render("\\alpha>1", element, {displayMode:false,throwOnError:false});$ is the enhancement multiplier for the hydrogen pick up fraction low LHR values, and LHR $var element = document.getElementById("moose-equation-fbef3e6e-f98b-4f52-bfda-4abca9fdc908");katex.render("_{lim}", element, {displayMode:false,throwOnError:false});$ is the LHR threshold value. Since Ref. (Nissen et al., 2014) does not provide quantitative values for this model, $var element = document.getElementById("moose-equation-128ffae7-08c5-4356-89d4-f273247b3f67");katex.render("f_{Zr2}", element, {displayMode:false,throwOnError:false});$ , $var element = document.getElementById("moose-equation-70e0c335-9e65-47e5-b8e0-7f333ebddd35");katex.render("\\alpha", element, {displayMode:false,throwOnError:false});$ , and LHGR $var element = document.getElementById("moose-equation-c3280f94-0aa4-4994-ad01-3a852fa766d9");katex.render("_{lim}", element, {displayMode:false,throwOnError:false});$ are user defined parameters and the current default values should not be granted meaning. This model is an alternative to the empirical bunrup-dependent approach. Since LHR is proportional to the heat flux on the cladding outer surface, it relates to temperature gradient within the cladding and can be used a proxy variable for the corrosion environment. Since hydrogen tends to diffuse down temperature gradients, a lower LHR might encourage hydrogen diffusion within the cladding thickness and further hydrogen pick up (Nissen et al., 2014).

note:Coupling with oxyde growth models

The hydrogen pickup model for options others than using the burnup-based model for Zircaloy-2 in a BWR environment can be coupled to any of the oxide growth models described in the Cladding and Corrosion Page.

Example Input Syntax

Using the default zirconium alloy of Zircaloy-4:

[BCs<<<{"href": "../../syntax/NuclearMaterials/BCs/index.html"}>>>]
  [hydrogen_pickup]
    type = HydrogenPickup<<<{"description": "Hydrogen flux BC for various zirconium alloys. Units for H concentration is wt.ppm.", "href": "HydrogenPickup.html"}>>>
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = hydrogen_in_solution_ppm
    boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = right
    oxide_thickness<<<{"description": "Coupled oxide thickness."}>>> = oxide_thickness
    clad_thickness<<<{"description": "Initial thickness of clad."}>>> = 660e-6
    clad_outer_radius<<<{"description": "Initial outer radius of the clad."}>>> = 660e-6
  []
[]

Using a user supplied hydrogen pickup fraction:

[BCs<<<{"href": "../../syntax/NuclearMaterials/BCs/index.html"}>>>]
  [hydrogen_pickup]
    type = HydrogenPickup<<<{"description": "Hydrogen flux BC for various zirconium alloys. Units for H concentration is wt.ppm.", "href": "HydrogenPickup.html"}>>>
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = hydrogen_in_solution_ppm
    boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = left
    oxide_thickness<<<{"description": "Coupled oxide thickness."}>>> = oxide_thickness
    input_pickup_fraction<<<{"description": "The hydrogen pickup fraction if specified by the user."}>>> = 0.15
    clad_thickness<<<{"description": "Initial thickness of clad."}>>> = 660e-6
  []
[]

Input Parameters

boundaryThe list of boundary IDs from the mesh where this object applies
C++ Type:std::vector<BoundaryName>
Controllable:No
Description:The list of boundary IDs from the mesh where this object applies
variableThe name of the variable that this residual object operates on
C++ Type:NonlinearVariableName
Unit:(no unit assumed)
Controllable:No
Description:The name of the variable that this residual object operates on

Required Parameters

axial_power_profileaxial power peaking function. Relevant to Zr-2 with linear power model.
C++ Type:FunctionName
Unit:(no unit assumed)
Controllable:No
Description:axial power peaking function. Relevant to Zr-2 with linear power model.
clad_outer_radius0.00418Initial outer radius of the clad.
Default:0.00418
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:Initial outer radius of the clad.
clad_thickness0.00066Initial thickness of clad.
Default:0.00066
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:Initial thickness of clad.
displacementsThe displacements
C++ Type:std::vector<VariableName>
Unit:(no unit assumed)
Controllable:No
Description:The displacements
fuel_pin_geometryName of the UserObject that reads the pin geometry from the mesh.
C++ Type:UserObjectName
Controllable:No
Description:Name of the UserObject that reads the pin geometry from the mesh.
hydrogen_pickup_model_Zr2The hydrogen pickup model for Zr-2, either based on burnup or linear power. Relevant to Zr-2, i.e., RECRYSTALLIZATION_ANNEALED
C++ Type:MooseEnum
Options:BURNUP_MODEL, LINEAR_POWER_MODEL
Controllable:No
Description:The hydrogen pickup model for Zr-2, either based on burnup or linear power. Relevant to Zr-2, i.e., RECRYSTALLIZATION_ANNEALED
input_linear_power_threshold10000The linear power threshold value if specified by the user. Relevant to Zr-2 with linear power model.
Default:10000
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:The linear power threshold value if specified by the user. Relevant to Zr-2 with linear power model.
input_pickup_fractionThe hydrogen pickup fraction if specified by the user.
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:The hydrogen pickup fraction if specified by the user.
input_pickup_fraction_multiplier1.5The hydrogen pickup fraction multiplier if specified by the user. Relevant to Zr-2 with linear power model.
Default:1.5
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:The hydrogen pickup fraction multiplier if specified by the user. Relevant to Zr-2 with linear power model.
linear_heat_rate_variableVariable holding the average linear heating rate in pellet in W/m
C++ Type:std::vector<VariableName>
Unit:(no unit assumed)
Controllable:No
Description:Variable holding the average linear heating rate in pellet in W/m
local_axial_burnupCoupled local axial burnup.
C++ Type:std::vector<VariableName>
Unit:(no unit assumed)
Controllable:No
Description:Coupled local axial burnup.
oxide_thicknessCoupled oxide thickness.
C++ Type:std::vector<VariableName>
Unit:(no unit assumed)
Controllable:No
Description:Coupled oxide thickness.
rod_ave_lin_powlinear power function. Relevant to Zr-2 with linear power model.
C++ Type:FunctionName
Unit:(no unit assumed)
Controllable:No
Description:linear power function. Relevant to Zr-2 with linear power model.
zircaloy_material_typeSTRESS_RELIEF_ANNEALEDType of zircaloy material properties to use in calculating creep. Note: ESCORE_IRRADIATIONGROWTHZR4 is not valid.
Default:STRESS_RELIEF_ANNEALED
C++ Type:MooseEnum
Options:STRESS_RELIEF_ANNEALED, RECRYSTALLIZATION_ANNEALED, PARTIAL_RECRYSTALLIZATION_ANNEALED, ZIRLO, M5, ESCORE_IRRADIATIONGROWTHZR4
Controllable:No
Description:Type of zircaloy material properties to use in calculating creep. Note: ESCORE_IRRADIATIONGROWTHZR4 is not valid.

Optional Parameters

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

Contribution To Tagged Field Data Parameters

control_tagsAdds user-defined labels for accessing object parameters via control logic.
C++ Type:std::vector<std::string>
Controllable:No
Description:Adds user-defined labels for accessing object parameters via control logic.
diag_save_inThe name of auxiliary variables to save this BC's diagonal jacobian contributions to. Everything about that variable must match everything about this variable (the type, what blocks it's on, etc.)
C++ Type:std::vector<AuxVariableName>
Unit:(no unit assumed)
Controllable:No
Description:The name of auxiliary variables to save this BC's diagonal jacobian contributions to. Everything about that variable must match everything about this variable (the type, what blocks it's on, etc.)
enableTrueSet the enabled status of the MooseObject.
Default:True
C++ Type:bool
Controllable:Yes
Description:Set the enabled status of the MooseObject.
implicitTrueDetermines whether this object is calculated using an implicit or explicit form
Default:True
C++ Type:bool
Controllable:No
Description:Determines whether this object is calculated using an implicit or explicit form
save_inThe name of auxiliary variables to save this BC's residual contributions to. Everything about that variable must match everything about this variable (the type, what blocks it's on, etc.)
C++ Type:std::vector<AuxVariableName>
Unit:(no unit assumed)
Controllable:No
Description:The name of auxiliary variables to save this BC's residual contributions to. Everything about that variable must match everything about this variable (the type, what blocks it's on, etc.)
seed0The seed for the master random number generator
Default:0
C++ Type:unsigned int
Controllable:No
Description:The seed for the master random number generator
skip_execution_outside_variable_domainFalseWhether to skip execution of this boundary condition when the variable it applies to is not defined on the boundary. This can facilitate setups with moving variable domains and fixed boundaries. Note that the FEProblem boundary-restricted integrity checks will also need to be turned off if using this option
Default:False
C++ Type:bool
Controllable:No
Description:Whether to skip execution of this boundary condition when the variable it applies to is not defined on the boundary. This can facilitate setups with moving variable domains and fixed boundaries. Note that the FEProblem boundary-restricted integrity checks will also need to be turned off if using this option
use_displaced_meshFalseWhether or not this object should use the displaced mesh for computation. Note that in the case this is true but no displacements are provided in the Mesh block the undisplaced mesh will still be used.
Default:False
C++ Type:bool
Controllable:No
Description:Whether or not this object should use the displaced mesh for computation. Note that in the case this is true but no displacements are provided in the Mesh block the undisplaced mesh will still be used.

Advanced Parameters

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

Material Property Retrieval Parameters

Input Files

(examples/spent_fuel/full_life_cycle_coarse/newHydrideModel.i)
(test/tests/hydrogen/hydrogen_pickup.i)
(test/tests/hydrogen/hydrogen_pickup_zr2_linear_power.i)
(examples/hydride_rim/hydride_rim.i)
(examples/non-cylindrical_fuel/2D/non-cyl_hydride.i)
(test/tests/hydrogen/hydrogen_pickup_zr2_burnup.i)

References

Paul Clifford. Pressurized-water reactor control rod ejection and boiling-water reactor control rod drop accidents. Technical Report Regulatory Guide RG 1.236, United States Nuclear Regulatory Commission, 2020.

@techreport{nrc_reg_guide_1236,
    author = "Clifford, Paul",
    institution = "United States Nuclear Regulatory Commission",
    pages = "1-29",
    title = "Pressurized-water reactor control rod ejection and boiling-water reactor control rod drop accidents",
    number = "Regulatory Guide RG 1.236",
    year = "2020"
}

Adrien Couet, Arthur T. Motta, and Robert J. Comstock. Hydrogen pickup measurements in zirconium alloys: relation to oxidation kinetics. Journal of Nuclear Materials, 451(1 - 3):1 – 13, 2014. URL: http://www.sciencedirect.com/science/article/pii/S0022311514001081, doi:http://dx.doi.org/10.1016/j.jnucmat.2014.03.001.

@article{couet2014,
    author = "Couet, Adrien and Motta, Arthur T. and Comstock, Robert J.",
    title = "Hydrogen pickup measurements in zirconium alloys: Relation to oxidation kinetics",
    journal = "Journal of Nuclear Materials",
    volume = "451",
    number = "1 - 3",
    pages = "1 - 13",
    year = "2014",
    doi = "http://dx.doi.org/10.1016/j.jnucmat.2014.03.001",
    url = "http://www.sciencedirect.com/science/article/pii/S0022311514001081"
}

Olivier Courty, Arthur T. Motta, and Jason D. Hales. Modeling and simulation of hydrogen behavior in Zircaloy-4 fuel cladding. Journal of Nuclear Materials, 452(1):311–320, 2014. doi:10.1016/j.jnucmat.2014.05.013.

@article{COURTY2014311,
    author = "Courty, Olivier and Motta, Arthur T. and Hales, Jason D.",
    title = "Modeling and simulation of hydrogen behavior in {Z}ircaloy-4 fuel cladding",
    journal = "Journal of Nuclear Materials",
    volume = "452",
    number = "1",
    pages = "311-320",
    year = "2014",
    issn = "0022-3115",
    doi = "10.1016/j.jnucmat.2014.05.013"
}

K. L. Nissen, W. Goll, V. I. Arimescu, G. Ledergerber, and C. Hellwig. Hydrogen uptake of BWR fuel rods. Power history effects at long irradiation times. In Proceedings of the Water Reactor Fuel Performance Meeting (WRFPM) 2014, Tokyo, Japan. Framatome, 9 2014.

@INPROCEEDINGS{Nissen2014,
    author = "Nissen, K. L. and Goll, W. and Arimescu, V. I. and Ledergerber, G. and Hellwig, C.",
    booktitle = "Proceedings of the Water Reactor Fuel Performance Meeting (WRFPM) 2014, Tokyo, Japan",
    month = "9",
    organization = "Framatome",
    title = "Hydrogen uptake of {BWR} fuel rods. {P}ower history effects at long irradiation times",
    year = "2014"
}

(test/tests/hydrogen/hydrogen_pickup.i)

# This input is designed to test the HydrogenPickup boundary condition. It
# computes the average total concentration of hydrogen C_H due to oxidation
# within cladding. This is computed via:
#
# C_H = 1e6 * 4f / R_{PB} * delta /(L-delta) * M_H / M_{Zr}
#
# where f is the pickup fraction, R_{PB} is the Pilling-Bedworth ratio (1.56),
# delta is the oxide thickness, L is the initial thickness of the cladding,
# M_H is the molecular weight of hydrogen (1.00794 g/mol), and M_{Zr} is the
# molecular weight of zirconium (91.224 g/mol).
#
# There are 5 different options for the hydrogen pickup fraction. The user can
# can supply the zircaloy_material_type (default stress relief annealed, i.e.,
# Zircaloy-4). Other options included ZIRLO, M5, Zircaloy-2, and user input.
#
# This input tests the case for a slab of Zircaloy-4 (stress relief annealed). Additional
# tests use cli_args in the tests spec file in the directoy of this input to
# test the ZIRLO, M5, and user options. Cases for each in a ring (RZ)
# configuration are also tested via cli_args.
#
# The oxide thickness is ramped over 2e8 seconds from 0 to 100 microns. 1e6 s
# timesteps are used. The value for the oxide thickness used in the calculation
# of the hydrogen flux is the average over the timestep to minimize timestep
# dependence of the solution.
#
# The hydrogen pickup value for Zircaloy-4 is 0.2 (20%).
#
# The calculated average hydrogen content is:
#
# C_H = 1e6 * 4 * 0.2 / 1.56 * 100e-6 / (660e-6 - 100e-6) * 1.00794 / 91.224.
#
# C_H = 1011.8192 wt. ppm.
#
# At 2e8 seconds BISON calculates a value of 1011.823 wt. ppm. One can
# approach more closely to the analytical solution by refining the timestep size
# and mesh.

[Mesh]
  [generated_mesh]
    type = GeneratedMeshGenerator
    dim = 1
    xmin = 0
    xmax = 660e-6
    nx = 10
  []
[]

[Variables]
  [hydrogen_in_solution_ppm]
  []
[]

[AuxVariables]
  [oxide_thickness]
    order = CONSTANT
    family = MONOMIAL
  []
  [grad_hydrogen_in_solution_ppm_x]
    order = CONSTANT
    family = MONOMIAL
  []
[]

[Functions]
  [oxide_thickness_function]
    type = PiecewiseLinear
    x = '0 2e8'
    y = '0 100e-6'
  []
[]

[Kernels]
  [dCprecdt]
    type = TimeDerivative
    variable = hydrogen_in_solution_ppm
  []
  [diffusion]
    type = MatDiffusion
    diffusivity = diffusivity
    variable = hydrogen_in_solution_ppm
  []
[]

[AuxKernels]
  [oxide]
    type = FunctionAux
    boundary = right
    variable = oxide_thickness
    function = oxide_thickness_function
    execute_on = 'linear timestep_end'
  []
  [grad_hydrogen_in_solution_ppm_x]
    type = VariableGradientComponent
    variable = grad_hydrogen_in_solution_ppm_x
    component = x
    gradient_variable = hydrogen_in_solution_ppm
    execute_on = 'linear timestep_end'
  []
[]

[Materials]
  [diffusivity]
    type = GenericConstantMaterial
    prop_values = 2.36e-14
    prop_names = diffusivity
  []
[]

[BCs]
  [hydrogen_pickup]
    type = HydrogenPickup
    variable = hydrogen_in_solution_ppm
    boundary = right
    oxide_thickness = oxide_thickness
    clad_thickness = 660e-6
    clad_outer_radius = 660e-6
  []
[]

[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'

  l_max_its = 60
  l_tol = 8e-3

  nl_max_its = 30
  nl_rel_tol = 1e-8
  nl_abs_tol = 1e-17

  end_time = 2e8

  [TimeStepper]
    type = ConstantDT
    dt = 1e6
  []

[]

[Postprocessors]
  [hydrogen_gradient]
    type = ElementExtremeValue
    value_type = max
    variable = 'grad_hydrogen_in_solution_ppm_x'
  []
  [oxide_thickness]
    type = ElementExtremeValue
    variable = oxide_thickness
  []
  [hydrogen_content]
    type = ElementAverageValue
    variable = hydrogen_in_solution_ppm
  []
[]

[Outputs]
  csv = true
[]

(examples/hydride_rim/hydride_rim.i)

# This input file shows how to use the hydride action
# along with hydrogen pick-up and oxidation

[Mesh]
  [generated_mesh]
    type = GeneratedMeshGenerator
    dim = 2
    xmin = 0
    xmax = 660e-6
    ymin = 0
    ymax = 660e-8
    nx = 40
    ny = 1
    bias_x = 1.15
  []
[]

[Variables]
  [hydrogen_in_solution_ppm]
    initial_condition = 0
    scaling = 1e6
  []
  [temperature]
    initial_condition = 615
  []
  [hydrogen_as_hydride_ppm]
    order = FIRST
    family = LAGRANGE
    initial_condition = 127
    scaling = 1e6
  []
[]

[AuxVariables]
  [oxide_thickness]
    # This needs to be monomial because it uses material properties (thermal_conductivity).
    order = CONSTANT
    family = MONOMIAL
  []
  # Used to prevent negative concentrations
  [bounds_dummy]
    order = FIRST
    family = LAGRANGE
  []
[]

[Bounds]
  # To prevent negative concentrations
  [Cprec_lower_bound]
    type = ConstantBounds
    variable = bounds_dummy
    bounded_variable = hydrogen_as_hydride_ppm
    bound_type = lower
    bound_value = 0.
  []
[]

[CladdingHydrides]
  [hydrides]
    block = 0
    temperature = temperature
    hydrogen_in_solution_ppm = hydrogen_in_solution_ppm
    hydrogen_as_hydride_ppm = hydrogen_as_hydride_ppm
  []
[]

[Kernels]
  # This gives a linear profile in temperature within a few seconds.
  [heatflow]
    type = HeatConduction
    variable = temperature
  []
 [dTdt]
    type = HeatConductionTimeDerivative
    variable = temperature
  []
[]

[AuxKernels]
  [oxide]
    type = MaterialRealAux
    boundary = left
    variable = oxide_thickness
    property = oxide_scale_thickness
  []
[]

[BCs]
  [leftt]
    type = DirichletBC
    variable = temperature
    boundary = left
    value = 603
  []
  [rightt]
    type = DirichletBC
    variable = temperature
    boundary = right
    value = 633
  []

  [hydrogen_pickup]
    type = HydrogenPickup
    variable = hydrogen_in_solution_ppm
    boundary = left
    oxide_thickness = oxide_thickness
    input_pickup_fraction = 0.15
    clad_thickness = 660e-6
  []
[]

[Materials]
  # Thermal conductivity for oxide aux kernel and for the temperature profile.
  [clad_thermal]
    type = HeatConductionMaterial
    block = 0
    thermal_conductivity = 16.0  # W/m-K
    specific_heat = 330.0        # J/kg-K
  []
  [density]
    type = ParsedMaterial
    block = 0
    property_name = density
    expression = 6551.0             # kg/m^3
  []
  [ZryOxidation]
    type = ZryOxidation
    boundary = left
    clad_inner_radius = 660e-6
    clad_outer_radius = 0
    temperature = temperature
  []
[]

[Executioner]
  type = Transient
  solve_type = 'PJFNK'

  petsc_options = '-snes_ksp_ew'
  petsc_options_iname = '-pc_type -pc_factor_mat_solver_package -snes_type'
  petsc_options_value = 'lu superlu_dist vinewtonrsls' # This petsc option helps prevent negative concentrations

  line_search = 'none'

  l_max_its = 60
  l_tol = 8e-3

  nl_max_its = 30
  nl_rel_tol = 1e-5
  nl_abs_tol = 1e-6

  end_time = 2e6

  dt = 1
  dtmax = 1e4

  [TimeStepper]
    type = IterationAdaptiveDT
    dt = 1
    growth_factor = 1.2
    optimal_iterations = 5
    iteration_window = 1
    linear_iteration_ratio = 100
  []

[]

[Postprocessors]
  [peak_value_hydride]
    type = ElementExtremeValue
    variable = 'hydrogen_as_hydride_ppm'
    block = 0
    value_type = 'max'
  []
[]

[PerformanceMetricOutputs]
[]

[VectorPostprocessors]
  [content]
    type = LineValueSampler
    variable = 'temperature hydrogen_in_solution_ppm hydrogen_as_hydride_ppm oxide_thickness'
    sort_by = x
    start_point = '0 300e-8 0'
    end_point = '660e-6 300e-8 0'
    num_points = 200
    outputs = profile
  []
[]

[Outputs]
  perf_graph = true
  exodus = true
  [console]
    type = Console
    max_rows = 25
  []
  [profile]
    type = CSV
    execute_on = 'final'
  []
[]

(examples/spent_fuel/full_life_cycle_coarse/newHydrideModel.i)

# This model is a linear element, 10 discrete fuel pellet stack (pellet_type_1) with a fine mesh.
# Modifying the base model to simulate the complete fuel life cycle from
# irradiation through dry storage
# Irradiation Time 3 years (6% burnup, ~ 60 MWd/kgU)
# Spent Fuel Pool 3 years
# Vacuum Drying 8 hours
# Dry Cask Storage (DCSS) 2 years
#
store_end = 4.732416e8 # 2 yrs storage

[GlobalParams]
  order = FIRST
  family = LAGRANGE
[]

[Mesh]
  coord_type = RZ
  displacements = 'disp_x disp_y'
  [mesh]
    type = FileMeshGenerator
    file = 'coarse10_rz_med2.e'
  []
[]

[Variables]
  [hydrogen_in_solution_ppm]
    initial_condition = 20 # Initial H from manufacture of rod.
    scaling = 1e6
  []
  [hydrogen_as_hydride_ppm]
    order = FIRST
    family = LAGRANGE
    initial_condition = 0
    scaling = 1e6
  []
[]

[AuxVariables]
  [temp]
  []
  [disp_x]
  []
  [disp_y]
  []
  [oxide_thickness]
    # This needs to be monomial because it uses material properties (thermal_conductivity) on a boundary.
    order = CONSTANT
    family = MONOMIAL
  []
  # Used to prevent negative concentrations
  [bounds_dummy]
    order = FIRST
    family = LAGRANGE
  []
[]

[Bounds]
  # To prevent negative concentrations
  [Cprec_lower_bound]
    type = ConstantBounds
    variable = bounds_dummy
    bounded_variable = hydrogen_as_hydride_ppm
    bound_type = lower
    bound_value = 0.
  []
[]

[CladdingHydrides]
  [hydrides]
    block = 'clad 3'
    temperature = temp
    hydrogen_in_solution_ppm = hydrogen_in_solution_ppm
    hydrogen_as_hydride_ppm = hydrogen_as_hydride_ppm
  []
[]

[BCs]
  [hydrogen_pickup]
    type = HydrogenPickup
    variable = hydrogen_in_solution_ppm
    boundary = '2'
    oxide_thickness = oxide_thickness
    input_pickup_fraction = 0.15
    clad_thickness = 560e-6
  []
[]

[Functions]
  [interpolated_temp]
    type = SolutionFunction
    from_variable = temp
    solution = read_thermo_solution
  []
  [interpolated_x]
    type = SolutionFunction
    from_variable = disp_x
    solution = read_thermo_solution
  []
  [interpolated_y]
    type = SolutionFunction
    from_variable = disp_y
    solution = read_thermo_solution
  []
  [interpolated_oxide]
    type = SolutionFunction
    from_variable = oxide_thickness
    solution = read_thermo_solution
  []
[]

[UserObjects]
  [read_thermo_solution]
    type = SolutionUserObject
    mesh = 'discrete_out.e'
    execute_on = 'nonlinear'
    system_variables = 'temp disp_x disp_y oxide_thickness'
  []
[]

[AuxKernels]
  [interp_temp]
    type = FunctionAux
    variable = temp
    function = interpolated_temp
  []
  [interp_x]
    type = FunctionAux
    variable = disp_x
    function = interpolated_x
  []
  [interp_y]
    type = FunctionAux
    variable = disp_y
    function = interpolated_y
  []
  [interp_oxide]
    type = FunctionAux
    variable = oxide_thickness
    boundary = '2'
    function = interpolated_oxide
  []
[]

[Materials]
  [clad_thermal]
    type = HeatConductionMaterial
    block = 'clad 3'
    thermal_conductivity = 16.0
    specific_heat = 330.0
  []
[]

[Executioner]
  type = Transient
  solve_type = 'PJFNK'

  petsc_options = '-snes_ksp_ew'
  petsc_options_iname = '-pc_type -pc_factor_mat_solver_package -snes_type'
  petsc_options_value = 'lu superlu_dist vinewtonrsls' # This petsc option helps prevent negative concentrations

  line_search = 'none'

  l_max_its = 50
  l_tol = 1e-3
  nl_max_its = 10
  nl_rel_tol = 1e-4
  nl_abs_tol = 1e-8

  start_time = 0
  end_time = ${store_end}

  dtmax = 1e6
  dtmin = 100

  [TimeStepper]
    type = IterationAdaptiveDT
    dt = 2e2
    optimal_iterations = 20
    iteration_window = 6
  []
  [Quadrature]
    order = fifth
    side_order = seventh
  []
[]

[Postprocessors]
  [max_hydrogen]
    type = NodalExtremeValue
    variable = 'hydrogen_in_solution_ppm'
    block = 'clad'
  []
  [max_hydride]
    type = ElementExtremeValue
    variable = 'hydrogen_as_hydride_ppm'
    block = 'clad'
  []
  [avg_hydrogen]
    type = ElementAverageValue
    variable = 'hydrogen_in_solution_ppm'
    block = 'clad'
  []
  [avg_hydride]
    type = ElementAverageValue
    variable = 'hydrogen_as_hydride_ppm'
    block = 'clad'
  []
[]

[PerformanceMetricOutputs]
[]

[Outputs]
  perf_graph = true
  exodus = true
  csv = true
  [console]
    type = Console
    output_file = 1
    output_linear = true
    max_rows = 25
    all_variable_norms = true
  []
  [chkfile]
    type = CSV
    show = 'max_hydrogen max_hydride'
    execute_on = final
  []
[]

(test/tests/hydrogen/hydrogen_pickup.i)

# This input is designed to test the HydrogenPickup boundary condition. It
# computes the average total concentration of hydrogen C_H due to oxidation
# within cladding. This is computed via:
#
# C_H = 1e6 * 4f / R_{PB} * delta /(L-delta) * M_H / M_{Zr}
#
# where f is the pickup fraction, R_{PB} is the Pilling-Bedworth ratio (1.56),
# delta is the oxide thickness, L is the initial thickness of the cladding,
# M_H is the molecular weight of hydrogen (1.00794 g/mol), and M_{Zr} is the
# molecular weight of zirconium (91.224 g/mol).
#
# There are 5 different options for the hydrogen pickup fraction. The user can
# can supply the zircaloy_material_type (default stress relief annealed, i.e.,
# Zircaloy-4). Other options included ZIRLO, M5, Zircaloy-2, and user input.
#
# This input tests the case for a slab of Zircaloy-4 (stress relief annealed). Additional
# tests use cli_args in the tests spec file in the directoy of this input to
# test the ZIRLO, M5, and user options. Cases for each in a ring (RZ)
# configuration are also tested via cli_args.
#
# The oxide thickness is ramped over 2e8 seconds from 0 to 100 microns. 1e6 s
# timesteps are used. The value for the oxide thickness used in the calculation
# of the hydrogen flux is the average over the timestep to minimize timestep
# dependence of the solution.
#
# The hydrogen pickup value for Zircaloy-4 is 0.2 (20%).
#
# The calculated average hydrogen content is:
#
# C_H = 1e6 * 4 * 0.2 / 1.56 * 100e-6 / (660e-6 - 100e-6) * 1.00794 / 91.224.
#
# C_H = 1011.8192 wt. ppm.
#
# At 2e8 seconds BISON calculates a value of 1011.823 wt. ppm. One can
# approach more closely to the analytical solution by refining the timestep size
# and mesh.

[Mesh]
  [generated_mesh]
    type = GeneratedMeshGenerator
    dim = 1
    xmin = 0
    xmax = 660e-6
    nx = 10
  []
[]

[Variables]
  [hydrogen_in_solution_ppm]
  []
[]

[AuxVariables]
  [oxide_thickness]
    order = CONSTANT
    family = MONOMIAL
  []
  [grad_hydrogen_in_solution_ppm_x]
    order = CONSTANT
    family = MONOMIAL
  []
[]

[Functions]
  [oxide_thickness_function]
    type = PiecewiseLinear
    x = '0 2e8'
    y = '0 100e-6'
  []
[]

[Kernels]
  [dCprecdt]
    type = TimeDerivative
    variable = hydrogen_in_solution_ppm
  []
  [diffusion]
    type = MatDiffusion
    diffusivity = diffusivity
    variable = hydrogen_in_solution_ppm
  []
[]

[AuxKernels]
  [oxide]
    type = FunctionAux
    boundary = right
    variable = oxide_thickness
    function = oxide_thickness_function
    execute_on = 'linear timestep_end'
  []
  [grad_hydrogen_in_solution_ppm_x]
    type = VariableGradientComponent
    variable = grad_hydrogen_in_solution_ppm_x
    component = x
    gradient_variable = hydrogen_in_solution_ppm
    execute_on = 'linear timestep_end'
  []
[]

[Materials]
  [diffusivity]
    type = GenericConstantMaterial
    prop_values = 2.36e-14
    prop_names = diffusivity
  []
[]

[BCs]
  [hydrogen_pickup]
    type = HydrogenPickup
    variable = hydrogen_in_solution_ppm
    boundary = right
    oxide_thickness = oxide_thickness
    clad_thickness = 660e-6
    clad_outer_radius = 660e-6
  []
[]

[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'

  l_max_its = 60
  l_tol = 8e-3

  nl_max_its = 30
  nl_rel_tol = 1e-8
  nl_abs_tol = 1e-17

  end_time = 2e8

  [TimeStepper]
    type = ConstantDT
    dt = 1e6
  []

[]

[Postprocessors]
  [hydrogen_gradient]
    type = ElementExtremeValue
    value_type = max
    variable = 'grad_hydrogen_in_solution_ppm_x'
  []
  [oxide_thickness]
    type = ElementExtremeValue
    variable = oxide_thickness
  []
  [hydrogen_content]
    type = ElementAverageValue
    variable = hydrogen_in_solution_ppm
  []
[]

[Outputs]
  csv = true
[]

(test/tests/hydrogen/hydrogen_pickup_zr2_linear_power.i)

# This input is designed to test the HydrogenPickup boundary condition. It
# computes the average total concentration of hydrogen C_H due to oxidation
# within cladding for Zircaloy-2. C_H is derived as:
#
# C_H = 1e6 * 4f / R_{PB} * delta /(L-delta) * M_H / M_{Zr}
#
# where f is the pickup fraction, R_{PB} is the Pilling-Bedworth ratio (1.56),
# delta is the oxide thickness, L is the initial thickness of the cladding,
# M_H is the molecular weight of hydrogen (1.00794 g/mol), and M_{Zr} is the
# molecular weight of zirconium (91.224 g/mol).
#
# For Zircaloy-2, the hydrogen pickup fraction is not a constant, and
# is expressed as:
#
# f = f_{Zr2}^0 * \alpha,
#
# where \alpha = 1 where the linear heat rate is higher than a threshold
# value LHR_{lim}, and \alpha>1 otherwise, leading to enhanced pick up
# at lower linear power.
#
# The calculated average hydrogen content for high linear power is:
#
# C_H = 1e6 * 4 * 0.2 / 1.56 * 100e-6 / (660e-6 - 100e-6) * 1.00794 / 91.224.
#
# C_H = 1011.8192 wt. ppm. BISON currently predicts 1011.82339 wt. ppm at
# 2e8 seconds
#
# The calculated average hydrogen content for low linear power is:
#
# C_H = 1e6 * 4 * 0.3 / 1.56 * 100e-6 / (660e-6 - 100e-6) * 1.00794 / 91.224.
#
# C_H = 1517.7289 wt. ppm. BISON currently predicts 1517.735086 wt. ppm at
# 2e8 seconds
#

[Mesh]
  [generated_mesh]
    type = GeneratedMeshGenerator
    dim = 1
    xmin = 0
    xmax = 660e-6
    nx = 10
  []
[]

[Variables]
  [hydrogen_in_solution_ppm]
    initial_condition = 0
  []
[]

[AuxVariables]
  [oxide_thickness]
    order = CONSTANT
    family = MONOMIAL
  []
  [grad_hydrogen_in_solution_ppm_x]
    order = CONSTANT
    family = MONOMIAL
  []
[]

[Functions]
  [oxide_thickness_function]
    type = PiecewiseLinear
    x = '0 2e8'
    y = '0 100e-6'
  []
  [linear_power_history]
    type = ConstantFunction
    value = 9999
  []
  [axial_peaking_factors]
    type = ConstantFunction
    value = 1
  []
[]

[Kernels]
  [dCprecdt]
    type = TimeDerivative
    variable = hydrogen_in_solution_ppm
  []
  [diffusion]
    type = MatDiffusion
    diffusivity = diffusivity
    variable = hydrogen_in_solution_ppm
  []
[]

[AuxKernels]
  [oxide]
    type = FunctionAux
    boundary = right
    variable = oxide_thickness
    function = oxide_thickness_function
    execute_on = 'linear timestep_end'
  []
  [grad_hydrogen_in_solution_ppm_x]
    type = VariableGradientComponent
    variable = grad_hydrogen_in_solution_ppm_x
    component = x
    gradient_variable = hydrogen_in_solution_ppm
    execute_on = 'initial linear timestep_end'
  []
[]

[BCs]
  [hydrogen_pickup]
    type = HydrogenPickup
    variable = hydrogen_in_solution_ppm
    boundary = right
    zircaloy_material_type = RECRYSTALLIZATION_ANNEALED
    hydrogen_pickup_model_Zr2 = LINEAR_POWER_MODEL
    rod_ave_lin_pow = linear_power_history
    axial_power_profile = axial_peaking_factors
    oxide_thickness = oxide_thickness
    clad_thickness = 660e-6
    clad_outer_radius = 660e-6
  []
[]

[Materials]
  [diffusivity]
    type = GenericConstantMaterial
    prop_values = 2.36e-14
    prop_names = diffusivity
  []
[]

[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'

  l_max_its = 60
  l_tol = 8e-3

  nl_max_its = 30
  nl_rel_tol = 1e-8
  nl_abs_tol = 1e-17

  end_time = 2e8

  [TimeStepper]
    type = ConstantDT
    dt = 1e6
  []

[]

[Postprocessors]
  [hydrogen_gradient]
    type = ElementExtremeValue
    value_type = max
    variable = 'grad_hydrogen_in_solution_ppm_x'
  []
  [oxide_thickness]
    type = ElementExtremeValue
    variable = oxide_thickness
  []
  [hydrogen_content]
    type = ElementAverageValue
    variable = hydrogen_in_solution_ppm
    execute_on = 'initial linear timestep_end'
  []
[]

[Outputs]
  csv = true
[]

(examples/hydride_rim/hydride_rim.i)

# This input file shows how to use the hydride action
# along with hydrogen pick-up and oxidation

[Mesh]
  [generated_mesh]
    type = GeneratedMeshGenerator
    dim = 2
    xmin = 0
    xmax = 660e-6
    ymin = 0
    ymax = 660e-8
    nx = 40
    ny = 1
    bias_x = 1.15
  []
[]

[Variables]
  [hydrogen_in_solution_ppm]
    initial_condition = 0
    scaling = 1e6
  []
  [temperature]
    initial_condition = 615
  []
  [hydrogen_as_hydride_ppm]
    order = FIRST
    family = LAGRANGE
    initial_condition = 127
    scaling = 1e6
  []
[]

[AuxVariables]
  [oxide_thickness]
    # This needs to be monomial because it uses material properties (thermal_conductivity).
    order = CONSTANT
    family = MONOMIAL
  []
  # Used to prevent negative concentrations
  [bounds_dummy]
    order = FIRST
    family = LAGRANGE
  []
[]

[Bounds]
  # To prevent negative concentrations
  [Cprec_lower_bound]
    type = ConstantBounds
    variable = bounds_dummy
    bounded_variable = hydrogen_as_hydride_ppm
    bound_type = lower
    bound_value = 0.
  []
[]

[CladdingHydrides]
  [hydrides]
    block = 0
    temperature = temperature
    hydrogen_in_solution_ppm = hydrogen_in_solution_ppm
    hydrogen_as_hydride_ppm = hydrogen_as_hydride_ppm
  []
[]

[Kernels]
  # This gives a linear profile in temperature within a few seconds.
  [heatflow]
    type = HeatConduction
    variable = temperature
  []
 [dTdt]
    type = HeatConductionTimeDerivative
    variable = temperature
  []
[]

[AuxKernels]
  [oxide]
    type = MaterialRealAux
    boundary = left
    variable = oxide_thickness
    property = oxide_scale_thickness
  []
[]

[BCs]
  [leftt]
    type = DirichletBC
    variable = temperature
    boundary = left
    value = 603
  []
  [rightt]
    type = DirichletBC
    variable = temperature
    boundary = right
    value = 633
  []

  [hydrogen_pickup]
    type = HydrogenPickup
    variable = hydrogen_in_solution_ppm
    boundary = left
    oxide_thickness = oxide_thickness
    input_pickup_fraction = 0.15
    clad_thickness = 660e-6
  []
[]

[Materials]
  # Thermal conductivity for oxide aux kernel and for the temperature profile.
  [clad_thermal]
    type = HeatConductionMaterial
    block = 0
    thermal_conductivity = 16.0  # W/m-K
    specific_heat = 330.0        # J/kg-K
  []
  [density]
    type = ParsedMaterial
    block = 0
    property_name = density
    expression = 6551.0             # kg/m^3
  []
  [ZryOxidation]
    type = ZryOxidation
    boundary = left
    clad_inner_radius = 660e-6
    clad_outer_radius = 0
    temperature = temperature
  []
[]

[Executioner]
  type = Transient
  solve_type = 'PJFNK'

  petsc_options = '-snes_ksp_ew'
  petsc_options_iname = '-pc_type -pc_factor_mat_solver_package -snes_type'
  petsc_options_value = 'lu superlu_dist vinewtonrsls' # This petsc option helps prevent negative concentrations

  line_search = 'none'

  l_max_its = 60
  l_tol = 8e-3

  nl_max_its = 30
  nl_rel_tol = 1e-5
  nl_abs_tol = 1e-6

  end_time = 2e6

  dt = 1
  dtmax = 1e4

  [TimeStepper]
    type = IterationAdaptiveDT
    dt = 1
    growth_factor = 1.2
    optimal_iterations = 5
    iteration_window = 1
    linear_iteration_ratio = 100
  []

[]

[Postprocessors]
  [peak_value_hydride]
    type = ElementExtremeValue
    variable = 'hydrogen_as_hydride_ppm'
    block = 0
    value_type = 'max'
  []
[]

[PerformanceMetricOutputs]
[]

[VectorPostprocessors]
  [content]
    type = LineValueSampler
    variable = 'temperature hydrogen_in_solution_ppm hydrogen_as_hydride_ppm oxide_thickness'
    sort_by = x
    start_point = '0 300e-8 0'
    end_point = '660e-6 300e-8 0'
    num_points = 200
    outputs = profile
  []
[]

[Outputs]
  perf_graph = true
  exodus = true
  [console]
    type = Console
    max_rows = 25
  []
  [profile]
    type = CSV
    execute_on = 'final'
  []
[]

(examples/non-cylindrical_fuel/2D/non-cyl_hydride.i)

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

[Problem]
 kernel_coverage_check = false
[]

[Mesh]
  displacements = 'disp_x disp_y'
  [mesh]
    type = FileMeshGenerator
    file = non-cyl_mesh_2d_fine_clad.e
  []
[]

[Variables]
  [hydrogen_in_solution_ppm]
    initial_condition = 20 # Initial H from manufacture of rod.
    scaling = 1e6
    block = 'cladding'
  []
  [hydrogen_as_hydride_ppm]
    order = FIRST
    family = LAGRANGE
    initial_condition = 0
    scaling = 1e6
    block = 'cladding'
  []
[]

[AuxVariables]
  [temp]
  []
  [disp_x]
  []
  [disp_y]
  []
  [oxide_thickness]
    order = CONSTANT
    family = MONOMIAL
  []
  # Used to prevent negative concentrations
  [bounds_dummy]
    order = FIRST
    family = LAGRANGE
  []
[]

[Functions]
  [interpolated_temp]
    type = SolutionFunction
    from_variable = temp
    solution = read_thermo_solution
  []
  [interpolated_x]
    type = SolutionFunction
    from_variable = disp_x
    solution = read_thermo_solution
  []
  [interpolated_y]
    type = SolutionFunction
    from_variable = disp_y
    solution = read_thermo_solution
  []
  [interpolated_oxide]
    type = SolutionFunction
    from_variable = oxide_thickness
    solution = read_thermo_solution
  []
[]

[CladdingHydrides]
  [hydrides]
    block='cladding'
    temperature = temp
    hydrogen_in_solution_ppm = hydrogen_in_solution_ppm
    hydrogen_as_hydride_ppm = hydrogen_as_hydride_ppm
  []
[]

[UserObjects]
  [read_thermo_solution]
    type = SolutionUserObject
    mesh = 'non-cyl_base_irrad_out.e'
    execute_on = 'nonlinear'
    system_variables = 'temp disp_x disp_y oxide_thickness'
  []
[]

[AuxKernels]
  [interp_temp]
    type = FunctionAux
    variable = temp
    function = interpolated_temp
  []
  [interp_x]
    type = FunctionAux
    variable = disp_x
    function = interpolated_x
  []
  [interp_y]
    type = FunctionAux
    variable = disp_y
    function = interpolated_y
  []
  [interp_oxide]
    type = FunctionAux
    variable = oxide_thickness
    boundary = 'side'
    function = interpolated_oxide
  []
[]

[BCs]
  [hydrogen_pickup]
    type = HydrogenPickup
    variable = hydrogen_in_solution_ppm
    boundary = 'side'
    oxide_thickness = oxide_thickness
    input_pickup_fraction = 0.15
    clad_thickness = 677e-6
  []
[]

[Materials]
  [zirc_thermal]
    type = HeatConductionMaterial
    block = 'cladding displacer fuel'
    thermal_conductivity = 16.0
    specific_heat = 330.0
  []
  [zirc_density]
    type = StrainAdjustedDensity
    block = 'cladding displacer fuel'
    strain_free_density = 6551.0
  []
[]

[Executioner]
  type = Transient
  solve_type = PJFNK
  petsc_options_iname = '-pc_type -pc_factor_mat_solver_package -ksp_gmres_restart -snes_ksp_ew_rtol0 -snes_ksp_ew_rtolmax -snes_ksp_ew_gamma -snes_ksp_ew_alpha -snes_ksp_ew_alpha2 -snes_ksp_ew_threshold -snes_type'
  petsc_options_value = ' lu                       superlu_dist               101                 0.5                  0.9                  1                  2                   2                    0.1 vinewtonrsls'

  line_search = 'none'

  start_time = 0
  end_time = 8e7

  dtmax = 1e6 # Value used for speed of example case, 1e4 would be more appropriate
  dtmin = 1.0

  nl_rel_tol = 1e-4
  nl_abs_tol = 1e-8
  l_max_its = 50
  l_tol = 8e-3
  nl_max_its = 30

  [TimeStepper]
    type = IterationAdaptiveDT
    dt = 2e2
    optimal_iterations = 15
    iteration_window = 3
    linear_iteration_ratio = 100
    growth_factor = 2
    cutback_factor = .5
  []
  [Quadrature]
    order = fifth
    side_order = seventh
  []
[]

[Postprocessors]
  [max_hydrogen]
    type = NodalExtremeValue
    variable = 'hydrogen_in_solution_ppm'
    block = 'cladding'
  []
  [max_hydride]
    type = ElementExtremeValue
    variable = 'hydrogen_as_hydride_ppm'
    block = 'cladding'
  []
  [avg_hydrogen]
    type = ElementAverageValue
    variable = 'hydrogen_in_solution_ppm'
    block = 'cladding'
  []
  [avg_hydride]
    type = ElementAverageValue
    variable = 'hydrogen_as_hydride_ppm'
    block = 'cladding'
  []
[]

[PerformanceMetricOutputs]
[]

[Outputs]
  perf_graph = true
  exodus = true
  csv = true
  [console]
    type = Console
    max_rows = 15
  []
  [chkfile]
    type = CSV
    show = 'max_hydrogen max_hydride'
    execute_on = final
  []
[]

(test/tests/hydrogen/hydrogen_pickup_zr2_burnup.i)

# This input is designed to test the HydrogenPickup boundary condition. It
# computes the average total concentration of hydrogen C_H due to oxidation
# within cladding for Zircaloy-2. This is computed via:
#
# C_H = 1.4 * (47.8 * exp(-1.3 / (1+Bu)) + 0.316 * Bu) below 50 GWd/tU
# C_H = 1.4 * (28.9 + exp(0.117 * (Bu-20)))            above 50 GWd/tU
#
# where Bu is the axial average burnup in the fuel next to the cladding (GWd/tU).
# This tests the low burnup case with the burnup linearly ramped over 2e8 s
# to 0.05 FIMA (47.5 GWd/tU) using 1e6 timesteps.
#
# The calculated average hydrogen content is:
#
# C_H = 1.4 * (47.8 * exp(-1.3/(1+47.5)) + 0.316 * 47.5)
# C_H = 86.164 wt. ppm.
#

[Mesh]
  [generated_mesh]
    type = GeneratedMeshGenerator
    dim = 1
    xmin = 0
    xmax = 660e-6
    nx = 10
  []
[]

[Variables]
  [hydrogen_in_solution_ppm]
    initial_condition = 18.238
  []
[]

[AuxVariables]
  [grad_hydrogen_in_solution_ppm_x]
    order = CONSTANT
    family = MONOMIAL
  []
  [burnup]
  []
  [local_axial_burnup]
    order = CONSTANT
    family = MONOMIAL
  []
[]

[Functions]
  [burnup_function]
    type = PiecewiseLinear
    x = '0 2e8'
    y = '0 0.05'
  []
[]

[Kernels]
  [dCprecdt]
    type = TimeDerivative
    variable = hydrogen_in_solution_ppm
  []
  [diffusion]
    type = MatDiffusion
    diffusivity = diffusivity
    variable = hydrogen_in_solution_ppm
  []
[]

[AuxKernels]
  [grad_hydrogen_in_solution_ppm_x]
    type = VariableGradientComponent
    variable = grad_hydrogen_in_solution_ppm_x
    component = x
    gradient_variable = hydrogen_in_solution_ppm
    execute_on = 'initial linear timestep_end'
  []
  [burnup]
    type = FunctionAux
    variable = burnup
    function = burnup_function
    execute_on = 'initial linear timestep_end'
  []
  [local_axial_burnup]
    type = SpatialUserObjectAux
    variable = local_axial_burnup
    user_object = local_axial_burnup
    execute_on = 'initial linear timestep_end'
  []
[]

[BCs]
  [hydrogen_pickup]
    type = HydrogenPickup
    variable = hydrogen_in_solution_ppm
    boundary = right
    zircaloy_material_type = RECRYSTALLIZATION_ANNEALED
    hydrogen_pickup_model_Zr2 = BURNUP_MODEL
    local_axial_burnup = local_axial_burnup
    clad_thickness = 660e-6
    clad_outer_radius = 660e-6
  []
[]

[Materials]
  [diffusivity]
    type = GenericConstantMaterial
    prop_values = 2.36e-14
    prop_names = diffusivity
  []
[]

[UserObjects]
  [local_axial_burnup]
    type = LayeredAverage
    variable = burnup
    block = 0
    direction = y
    num_layers = 1
    execute_on = 'initial linear timestep_end'
  []
[]

[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'

  l_max_its = 60
  l_tol = 8e-3

  nl_max_its = 30
  nl_rel_tol = 1e-8
  nl_abs_tol = 1e-17

  end_time = 2e8

  [TimeStepper]
    type = ConstantDT
    dt = 1e6
  []

[]

[Postprocessors]
  [hydrogen_gradient]
    type = ElementExtremeValue
    value_type = max
    variable = 'grad_hydrogen_in_solution_ppm_x'
  []
  [hydrogen_content]
    type = ElementAverageValue
    variable = hydrogen_in_solution_ppm
    execute_on = 'initial linear timestep_end'
  []
[]

[Outputs]
  csv = true
[]

Description
Hydrogen pick up fraction models for Zircaloy - 4 , M 5 , Zirlo , and Optimized Zirlo
Hydrogen pick up fraction models for Zircaloy - 2
Example Input Syntax
Input Parameters
Input Files
References