14.3. Fuel and Fission-gas Ejection from the Pins

14.3.1. Physical Model and Assumptions

In the previous section, the sinks for the in-pin motion due to fuel and fission-gas ejection appeared in the fuel, free fission gas, and dissolved fission-gas mass conservation equations (see Eq. (14.2-9), 14.2-15, and 14.2-20). These sink terms are, of course, source terms in the channel thermal-hydraulics treatment.

The approach taken for the calculation of the fuel and fission-gas ejection from the pins is based on the assumption that there will always be pressure equilibrium established between a cavity node and the adjacent channel node if fuel and fission-gas ejection occurs. This assumption is justified for cladding ruptures with dimensions larger than one-pin diameter because of the very short distance between cavity nodes and corresponding coolant channel nodes. For pin-hole type failures leading to a pure gas release, which is not currently treated in PLUTO2, an orifice equation for calculating this gas release may be more appropriate. In the older PLUTO code [14-1, 14-3, 14-4], both an asymptotic orifice equation and an ejection calculation, based on pressure equilibration, were computed for the fuel and gas ejections, and the smaller of the two predictions was used. The asymptotic orifice equation usually predicted the higher ejection rates. These were causing a higher pressure in the channel node than in the adjacent cavity node which was considered non-physical.

To achieve the pressure equilibrium between the cavity and channel nodes in PLUTO2, an estimated amount of fuel and gas (with the same volume ratio as present in the cavity node) is ejected into the channel and the resulting cavity and channel pressures are calculated. This calculation is repeated several times with updated estimates for the ejected fuel and fission-gas masses until a pressure equilibrium is achieved. This approach does not consider (or require) the cladding rupture area, which is an advantage because the evolution of the cladding rupture is not known. The cladding rupture sizes found in the post-test examinations of TREAT tests H5 [14-24] and J1 [14-25], which were terminated at a time when the pins were still largely intact, show large enough rupture sizes to justify the assumption of pressure equilibration between the cavity node behind the rupture and the adjacent channel node. Moreover, these rupture areas are considerably larger than the cross section of the molten cavity (50% areal melt fraction in an FFTF type pin means a molten cross sectional area of about 0.1 cm2; cladding rupture areas in H5 and J1 are about 1 cm2). This means that the controlling aspect of the fuel and gas ejection (at least in mild TOP accidents) is really the cross sectional area of the cavity, which is taken into account in the in-pin motion calculation in PLUTO2. In the SAS3D model SAS/FCI, which does not consider the in-pin motion mechanistically, the cross sectional area for the orifice equation, used in the calculation of the SAS/FCI fuel and gas ejection, has to be set to about twice the cross sectional area of the cavity near the failure location [14-3]. Under LOF-driven-TOP conditions, when axial cladding rupture propagation is likely (see Section 14.3.3), the in-pin motion is not dominant in controlling the ejection because several contiguous cavity nodes may eject fuel and gas simultaneously. In addition, the time scales are very short under these conditions, leading to less axial fuel and gas convection into an ejection node than can actually be ejected.

Fuel and fission gas are being ejected from the fuel pins with the same volume ratio as exists in the ejecting cavity node. For conditions involving high void fractions in the ejecting cavity node, it is possible that fission gas or fuel vapor will be ejected preferentially because the remaining fuel may be in an annular flow configuration.

Another assumption about the fuel and gas ejection is that no backflow of materials from the channel to the pin is allowed when the pressure in the coolant channel node adjacent to the rupture becomes higher than the pressure in the ejection node. To model this backflow would be difficult and is not warranted because only a small amount of liquid sodium flowing back into the pin would lead to enough vaporization and pressurization to inhibit further backflow of materials. However, the potential pressure increase in the cavity node may diminish the in-pin fuel motion. This could be beneficial if midplane failures are assumed. However, most channel pressurizations that have been observed experimentally and that may be due to FCIs are not sufficient in magnitude or duration to cause concern about this assumption.

14.3.2. Numerical Solution of the Fuel and Gas Ejection Calculation

The numerical evaluation of the amount of fuel and fission gas ejected during a PLUTO2 time step is performed as a second step after the fuel mass conservation Eq. (14.2-9) has been solved and advanced to the end of the time step without accounting for the ejection term \(S'_{\text{fuca,ej}}\). This section describes how the fuel smear density value obtained at the end of the time step without accounting for ejection \({\rho'}_{\text{fuca}}^{\text{old}}\) is corrected for fuel ejection, thus obtaining the actual fuel smear density at the end of time step the \({\rho'}_{\text{fuca}}^{\text{new}}\) and the resulting equilibrated pressure of the ejecting cavity cell.

If \(\text{FN}_{\text{k}}\) is the fraction of the fuel in a cavity node that is ejected during a time step, the new generalized fuel smear density in the cavity node will be:

(14.3-1)\[ {\rho'}_{\text{fuca}}^{\text{new}} = {\rho'}_{\text{fuca}}^{\text{old}} - {\rho'}_{\text{fuca}}^{\text{old}} \cdot \text{FN}_{\text{k}} = {\rho'}_{\text{fuca}}^{\text{old}} - \text{FF}\]

where \(\text{FF} = {\rho'}_{\text{fuca}}^{\text{old}} \cdot \text{FN}_{\text{K}}\) and is actually the same as the sink term for fuel ejection, \({S'}_{\text{fuca, ej}} \cdot \Delta t_{\text{PL}}\) (see Eq. (14.2-15)).

Since the fission gas and fuel are assumed to be ejected with the same volume ratio as present in the ejecting cell, the new generalized free-fission-gas smear density is

(14.3-2)\[ {\rho'}_{\text{fica}}^{\text{new}} = {\rho'}_{\text{fica}}^{\text{old}} - {\rho'}_{\text{fica}}^{\text{old}} \cdot \text{FN}_{\text{k}} = {\rho'}_{\text{fica}}^{\text{old}} - \text{FF} \cdot \frac{{\rho'}_{\text{fica}}^{\text{old}}}{{\rho'}_{\text{fuca}}^{\text{old}}}\]

The new cavity fuel fraction will be

(14.3-3)\[ \theta_{\text{fica}}^{\text{new}} = \theta_{\text{fica}}^{\text{old}} - \frac{\text{FF}}{\rho_{\text{fuca}}}\]

where \(\rho_{\text{fuca}} = \rho_{\text{fuca}} \left( T_{\text{fuca}} \right)\) is the theoretical fuel density.

By using Eq. (14.3-1), Eq. (14.3-2), and Eq. (14.3-3) in the cavity pressure calculation of Eq. (14.2-39) and Eq. (14.2-40) and by dropping the superscripts “old”, one obtains

(14.3-4)\[\begin{split} P_{\text{ca}}^{\text{new}} = P_{\text{fvca}} + R_{\text{fi}} \cdot T_{\text{fuca}} \cdot \left( {\rho'}_{\text{fica}} - \text{FF} \cdot \frac{{\rho'}_{\text{fica}}}{{\rho'}_{\text{fuca}}} \right) \\ \big/ \left\{ \theta_{\text{ca}} - \theta_{\text{fuca}} + \frac{\text{FF}}{\rho_{\text{fuca}}\left( T \right)} + \left\lbrack \theta_{\text{fuca}} - \frac{\text{FF}}{\rho_{\text{fuca}}\left( T \right)} \right\rbrack K_{\text{fu}} \cdot \left( P_{\text{ca}}^{\text{new}} - P_{\text{fuca}} \right) \right\} \\\end{split}\]

where

\(\theta_{\text{ca}} - \theta_{\text{fuca}} = \theta_{\text{fica}}\)

By collecting the terms with FF on the left-hand side one obtains

(14.3-5)\[\begin{split} \text{FF} \cdot \left\lbrack \frac{\left( P_{\text{ca}}^{\text{new}} - P_{\text{fvca}} \right)}{\rho_{\text{fuca}}} + \frac{{\rho'}_{\text{fica}}}{{\rho'}_{\text{fuca}}} \cdot R_{\text{fi}} \cdot T_{\text{fuca}} \right. \\ \left. - \frac{K_{\text{fu}}}{\rho_{\text{fu}}} \cdot \left( P_{\text{ca}}^{\text{new}} - P_{\text{fvca}} \right)^{2} \right\rbrack = R_{\text{fi}} \cdot T_{\text{fuca}} \cdot {\rho'}_{\text{fica}} \\ - \theta_{\text{fica}} \left( P_{\text{ca}}^{\text{new}} - P_{\text{fvca}} \right) - \theta_{\text{fica}} \cdot K_{\text{fu}} \cdot \left( P_{\text{ca}}^{\text{new}} - P_{\text{fvca}} \right)^{2}\end{split}\]

This equation is used to calculate the fuel ejection FF for an estimated \(P_{\text{ca}}^{\text{new}}\). FF is used, in turn, to update the fission-gas and fuel densities in the channel node adjacent to the ejecting cavity node. A new channel pressure is then calculated based on the updated fuel and gas densities. The method of calculating the channel pressures is discussed in Section 14.4.5. If the new channel pressure is not within 1% of the estimated new cavity pressure, a better estimate for the new cavity pressure is made and the calculational sequence for the fuel and gas ejection is repeated.

Because the conditions in the coolant channel node can vary from normal coolant flow to fully voided, the procedure of estimating the new cavity pressure is somewhat complex. It has evolved through trial and error and it usually leads to convergence in a few iterations. If the iteration does not converge, it usually indicates that a non-physical condition has developed in the pin or channel. When the axial extent of the cladding rupture includes more numerical nodes as a result of axial pin failure propagation (which is likely under LOF′d′TOP conditions), the ejection calculation will be performed for each of the axial failure nodes.

The ejection calculation from the pins was originally formulated to allow a simultaneous ejection from the three pin failure groups at any given axial location. As mentioned earlier, the second and third pin failure groups are not yet operational. Because of the complexity of the current ejection calculation from only one pin failure group, it is expected that a simultaneous ejection from three pin failure groups into the same coolant channel node may become too complicated. Therefore, a simpler approach, which would only allow the ejection from the cavity node with the highest pressure during a PLUTO2 time step, may have to be adopted. In this approach it would, of course, still be possible for the different pin failure groups to eject fuel simultaneously at different axial locations.

The fuel and gas ejection calculation is done at the end of subroutine PL1PIN which was discussed in Section 14.1.2. The ejection calculation is done in the calculational sequence after all in-pin as well as channel mass and energy conservation have been solved and before the momentum equations for both the in-pin and channel flows are solved. The ejection calculation thus provides updated pressures for the momentum equations.

14.3.3. Axial Pin Failure Propagation

Rapid axial pin failure propagation due to the stress concentrations at the edges of the initial rupture (zipper-type failure propagation) has not been accepted as a mechanism for rapid axial failure propagation (i.e., on a millisecond or submillisecond time scale) [14-26]. However, if the conditions in a fuel-pin node are such that failure of the cladding should occur irrespective of cladding stress concentrations due to an adjacent cladding rupture and whether the potential new failure node location is adjacent to the initial failure or not, there is no reason to disallow this type of failure propagation.

The failure of nodes other than the initial one is determined by PLUTO2 because the DEFORM calculation is not done in a calculational channel once PLUTO2 (or LEVITATE) has been initiated. There are two options for pin failure propagation in PLUTO2. The first option is a mechanistic one based on the burst strength of the cladding. The second option allows the pin failure propagation to proceed according to a combination of input failure criteria based on fuel melt fraction, cladding temperature, and pressure differential between cavity and channel.

The mechanistic option (KFAILP=0) assumes that the solid fuel surrounding the molten cavity is completely cracked and that the radial stress at the fuel-cladding interface can therefore be calculated by reducing the cavity pressure by a geometrical reduction factor.

(14.3-6)\[ \sigma_{\text{radial,fc,K}} = \frac{P_{\text{ca,K}} \cdot D_{\text{ca,K}}}{\left( 2R_{\text{cl,fc,K}} \right)}\]

where

\(R_{\text{cl,fc,K}}\) is the inner radius of the cladding at axial node \(K\).

\(D_{\text{ca,K}}\) is the diameter of the molten cavity at axial node \(K\).

The hoop stress in the cladding can be calculated from

(14.3-7)\[ \sigma_{\text{hoop,K}} = \left( \sigma_{\text{radial,fc,K}} - P_{\text{ch,i}} \right) \frac{R_{\text{cl,fc,K}}}{\Delta R_{\text{cl,K}}}\]

where \(\Delta R_{\text{cl,K}}\) is the cladding thickness at axial node \(K\), and subscript \(i\) refers to the channel node which is adjacent to cavity node \(K\); \(i = K + \text{IDIFF}\) where IDIFF is the number of channel nodes below cavity node \(K = 1\).

By substituting Eq. (14.3-6) for the radial stress in Eq. (14.3-7),

(14.3-8)\[ \sigma_{\text{hoop,K}} = \frac{\left( P_{\text{ca,K}} \cdot D_{\text{ca,K}} - 2 \cdot P_{\text{ch,i}} \cdot R_{\text{cl,fc,K}} \right)}{\left( 2\Delta R_{\text{cl,K}} \right)}\]

Failure is assumed to occur when the following two conditions are satisfied:

(14.3-9)\[ \sigma_{\text{hoop,K}} > \sigma_{\text{UTS}} \left( T_{\text{cl,in,K}} \right)\]
(14.3-10)\[ \text{AREA}_{\text{ca,K}} > \text{FNARME} \cdot \pi \cdot R_{\text{cl,fc,K}}^{2}\]

where

FNARME = an input melt fraction (or more precisely cavity cross sectional area divided by cross sectional area of the fuel) which should be set to \(\geq 0.2\) because the PLUTO2 in-point motion model requires a moderately sized cavity diameter to run well.

\(\sigma_{\text{UTS}}\) = the burst strength of the cladding that is a function of cladding temperature. This is calculated from the middle cladding node temperature in function subroutine UTS, which is described in Chapter 8 on fuel behavior.

A problem with this burst failure criterion is that it is not compatible with a melt fraction failure criterion, cladding strain failure criterion or any other criteria that are not directly related to the cladding strength. Since these latter criteria may be used in DEFORM to predict the initial failure, care has to be taken that the cladding will not rip open along sizable length as soon as PLUTO2 takes over the calculation. For example, this could happen if the initial pin failure was assumed to occur due to a melt fraction criterion at a time when the cavity pressure was already high enough to burst the cladding at several axial elevations. If a melt fraction criterion was chosen for the initial failure prediction, the problem that several nodes will instantaneously rupture once PLUTO2 initiates, can be avoided by setting the input parameter FNARME (see Eq. (14.3-10)) equal to or greater than the input melt fraction criterion for the initial failure.

The non-mechanistic input failure propagation criterion (KFAILP=1) involves several separate criteria. Failure occurs when the cladding middle node and outer surface temperatures, the cross sectional area of the cavity, and the pressure difference between cavity and channel exceed input criteria.

(14.3-11)\[ \begin{aligned} T_{\text{cl,in,K}} && \text{and} && T_{\text{cl,os,K}} > \text{TEFAIL} \end{aligned}\]

and

(14.3-12)\[ \text{AREA}_{\text{ca,K}} > \text{FNARME} \cdot \pi \cdot R_{\text{cl,fc,K}}^{2}\]

and

(14.3-13)\[ P_{\text{ca,K}} > P_{\text{ch,i}} + \text{PRFAIL}\]

where

TEFAIL, FNARME, and PRFAIL are input parameters.

Subscripts in and os refer to the middle and outer cladding nodes, respectively. This combined criterion can of course be simplified by setting some of the input values to extreme numbers. For example, if TEFAIL is set to zero and PRFAIL to \(-10^{20}\), this criterion will reduce to a pure melt fraction criterion. This input failure propagation criterion may be useful for comparison calculations with codes that allow failure propagation only according to one of the above non-mechanistic criteria. Moreover, if the initial failure criterion is a melt fraction criterion, one may also want to do the failure propagation based on a melt fraction criterion or on a melt fraction criterion combined with the condition that some overpressure has to exist in the molten cavity nodes.

14.3.4. Complete Pin Disruption and Switch to LEVITATE

Once the cladding has melted at a certain axial location and extensive fuel melting has also occurred at the same location, complete fuel-pin disruption will eventually occur. This cannot be treated with PLUTO2, and therefore, a switch to LEVITATE is made under these conditions. For each pin node, PLUTO2 checks whether the middle and outer cladding temperatures have exceeded the liquidus and whether the fuel-surface temperature has exceeded the solidus.

\(T_{\text{cl, in, K}} > T_{\text{cl, liq}}\)

and

\(T_{\text{cl, os, K}} > T_{\text{cl, liq}}\)

and

\(T_{\text{fu, os, K}} > T_{\text{fu,sol}}\)

where

cl, in refers to the middle cladding node,

cl, os refers to the outer cladding node,

fu, os refers to the outermost fuel node.

If the above conditions are met, control will be transferred from PLUTO2 to LEVITATE. These checks are performed at the end of the PLUTO2 driver routine. At that location in the code, it is also checked whether cladding melting has occurred in more than NCPLEV nodes (NCPLEV is input). If this criterion is met, a switch to LEVITATE will also occur. Moreover, a switch to LEVITATE is also made if the fuel temperature in the cavity or in the channel exceeds 4000 K because significant fuel vapor pressures will develop beyond this temperature.