Importance Analyses on Phenomenological Parameters for the Aerosol Dynamics Models in I-COSTA for a Severe Nuclear Power Plant Accident

Abstract

In this study, using in-house code I-COSTA, importance analyses are performed on the phenomenological parameters in the aerosol dynamics using International Standard Problem No. 44. The analyses consider twelve parameters used in multicomponent sectional equations and Mason equations. For the first step of the analysis, Latin hypercube sampling is performed for the aforementioned parameters, and the number of samplings is determined using a comparison of averages and standard deviations between those samplings and the ones gathered from continuous distributions of the parameters. Sensitivity analyses are then performed on the airborne concentrations of the aerosol particles using I-COSTA, and the results are used to obtain the correlation coefficients between the parameters and the airborne concentrations. From the analyses, the dynamic shape factor, which accounts for the drag force of the non-spherical aerosol particles, is found to be one of the most important parameters in the aerosol dynamics. The saturation ratio in the Mason equation is also found to be an important parameter for aerosol particles with high solubility since the mass of the aforementioned particles is sensitive to the hygroscopic growth rate.

1. Introduction

The nuclear accident source term is one of the most important engineering values for safety analyses in nuclear engineering since it is used to determine how a nuclear power plant is designed to manage and mitigate accidents. The source term provides the initial and boundary conditions for the consequence analysis of a nuclear accident so that the licensee of a nuclear power plant can determine the accident management plan. In addition, as the market is increasingly focused on small modular reactors (SMRs), practical determination of the size of the emergency planning zone (EPZ) is one of the most important pending issues for the commercialization of SMRs [1]. In fact, the EPZ size is directly related to the accident source term. The accident source term is the key factor with respect to practical sizing of EPZ, i.e., the size of the EPZ is decreased compared to that of a large-scale nuclear power plant.

The accident source term is determined by the amount of fission products released as a nuclear accident progresses. Except for noble gases (xenon, krypton) and some iodine chemicals, most of the fission products released from the reactor core to the reactor coolant system are in aerosol form; they comprise the condensation of fission products caused by the temperature difference between the reactor core and the reactor coolant system. Therefore, most of the severe accident analysis codes have focused on aerosol dynamics in order to determine the accident source term for a severe accident.

Based on the results of the Phébus FP test [2] and investigations of the Fukushima Daiichi Nuclear Power Plant accident, significant discrepancies have been observed between the measured data and the calculated values obtained from severe accident analysis codes [3,4]. As discussed in the previous paragraph, since the accident source term during a severe accident is governed by aerosol behavior, these discrepancies can be attributed to knowledge gaps in our understanding of aerosol behavior. One of the key phenomena contributing to this gap is the remobilization of potentially volatile fission products (e.g., cesium, iodine, tellurium, ruthenium, etc.) that were deposited on the surfaces of the reactor coolant system and within the containment area. This remobilization has recently been identified as a factor that can lead to significant delayed releases into the environment [5,6,7,8]. According to the observations in the Phébus FP test [2], the deposited fission products may act as a source for delayed release if thermal and/or chemical changes occur in the reactor coolant system; for example, significant remobilization of the cesium deposits from the hot leg of the Phébus experimental circuit was found at the end of the degradation phase, representing about 15~75% of the deposited mass of the cesium. However, most of the severe accident analysis codes cannot consider the aforementioned delayed release. Therefore, to narrow this gap, it is necessary to investigate the mechanism in detail through experimental activities.

Many international cooperative research projects have been conducted to narrow the aforementioned gap. One such example is the ESTER (Experiments on the Source Term for Delayed Release) project, initiated by the OECD/NEA [6,9]. Experimental activities under this project revealed that delayed releases may occur from aerosols deposited on the surfaces of the containment area, producing organic and/or inorganic iodine species [2,10,11,12,13,14]. These iodine compounds can subsequently be transformed into aerosol form through a nucleation process [15,16].

To apply the lessons learned from the experiments to analysis, it is necessary to develop analytical models for the relevant phenomena. As an essential step in model development—particularly concerning aerosol behavior—it is important to consider phenomenological parameters in the aerosol dynamics code. These parameters describe the physical and chemical characteristics and behaviors of aerosol particles formed by various chemical species, such as iodine oxide aerosols, cesium iodide (CsI) aerosols, and others. The aerosol dynamic equations involve numerous phenomenological parameters, some of which are derived from experimental measurements, while others are obtained through calculations. However, it is not feasible to acquire complete information for all parameters associated with the wide range of chemical species that may form aerosol particles.

It is, therefore, important to prioritize parameters related to aerosol behavior based on their significance for the efficient development of analytical models. For effective prioritization, the magnitude of phenomenological uncertainties in aerosol behavior caused by each parameter must be considered. To support this, a number of studies have been conducted to quantify phenomenological uncertainty. Power et al. [17] provided a list of phenomenological parameters and their associated uncertainty ranges based on aerosol dynamics. However, their study does not offer insight into the relative importance of each parameter in terms of its impact on airborne concentrations.

More recently, Malicki and Lind [18] conducted sensitivity studies on aerosol behavior during the Phébus FPT1 test. Drawing on extensive literature reviews of parameters affecting aerosol dynamics, they presented detailed results using the MELCOR code. However, because their analysis was based on data from integral experiments—which involve highly complex phenomena such as core degradation, fission product release, and coupled thermal-hydraulic effects along with their associated uncertainties—it is difficult to isolate uncertainties specifically related to aerosol behavior. Moreover, due to the use of an integral safety analysis code, identifying the detailed mechanisms that cause uncertainties in aerosol behavior from a theoretical perspective remains challenging.

In this study, the authors perform importance analyses on parameters related to the aerosol dynamics. The analyses are conducted by first carrying out random sampling of the parameters via an in-house code named R-SAPhe, which is based on Latin hypercube sampling with the uncertainty ranges. The uncertainty ranges of the parameters are taken from preceding studies on the aerosol dynamics parameters [17,18]. Sensitivity analyses are then performed on the KAEVER-148, KAEVER-186, and KAEVER-187 experiments, which are included in International Standard Problem No. 44 on the aerosol dynamics [19]. The sensitivity analyses are conducted using the in-house code I-COSTA (In-Containment Source Term Analysis), which was developed by one of the authors of this paper [20]. Using the sensitivity analyses results, the authors obtain correlation coefficients on the parameters in the aerosol dynamics and the figure-of-the-merits in the importance analyses via CC-SAPhe. The figure-of-merit considered in this work are the maximum concentrations of the aerosol particles in the atmosphere and the concentrations at the end of the analyses. Since the analyses are performed with the results of separate effect tests on the aerosol dynamics, the effects of the uncertainty in the parameters are shown in terms of the sensitivity of the aerosol concentrations, which can be distinguished from the uncertainty in the measurement during experiments. The parameters identified as important are discussed in relation to the fundamental mechanisms governing specific aerosol behaviors, where these parameters are used to theoretically explain processes such as coagulation, deposition, and other related phenomena.

This paper is organized as follows. Section 2 presents summaries of the models in the I-COSTA code. In this section, the uncertainty ranges of the parameters are discussed as well as the overall procedures for the importance analysis. In Section 3, the importance analyses are performed for the aforementioned three experiments. Finally, Section 4 provides a summary and conclusions of this paper.

2. Summary of the Numerical Methods in the I-COSTA Code

2.1. Summary of Multicomponent Sectional Equations and Mason Equations for Hygroscopic Growth in I-COSTA

The behaviors of the aerosol particles in the containment are analyzed by solving two sets of governing equations. The first are multicomponent sectional equations [21] derived from integration of General Dynamic Equations (GDE) over each discretized section on the size distribution of the aerosol particles. The other are Mason Equations [22] to treat hygroscopic growth of the aerosol particles separately since the multicomponent sectional equations may require a huge computation burden in order to avoid the inherent numerical diffusion.

The multicomponent sectional equations analyze the following four processes: (1) coagulation, (2) interparticle chemical reactions (a topic not covered in this paper), (3) injection of the particles into the volume, and (4) deposition on the surface. If the average mass of the aerosol particle in section l is more than two times greater than that in the previous section, i.e., v_{i + 1} ≥ 2v_i, the equations are expressed as follows:

(1)$\begin{array}{l}{\displaystyle \frac{d{Q}_{l,k}}{dt}}={\displaystyle \frac{1}{2}}\cdot {\displaystyle \sum _{l=1}^{l-2}{}^{1a}\overline{\beta }_{i,l-1,l}\cdot Q_{l-1,k}\cdot {\overline{Q}}_i}+{\displaystyle \frac{1}{2}}\cdot {\displaystyle \sum _{l=1}^{l-2}{}^{1b}\overline{\beta }_{i,l-1,l}\cdot Q_{i,k}\cdot {\overline{Q}}_{l-1}}\\ +{\displaystyle \frac{1}{2}}\cdot {\displaystyle \sum _{j=1}^{l-2}{}^{1a}\overline{\beta }_{l-1,j,l}\cdot Q_{j,k}\cdot {\overline{Q}}_l}+{\displaystyle \frac{1}{2}}\cdot {\displaystyle \sum _{j=1}^{l-2}{}^{1b}\overline{\beta }_{l-1,j,l}\cdot Q_{l-1,k}\cdot {\overline{Q}}_j}\\ +{\displaystyle \frac{1}{2}}\cdot \left[{}^{1a}\overline{\beta }_{l-1,l-1,l}+{}^{1b}\overline{\beta }_{l-1,l-1,l}\right]\cdot Q_{l-1,k}\cdot {\overline{Q}}_{l-1}\\ -{\displaystyle \sum _{i=1}^{l-1}\left[{}^{2a}\overline{\beta }_{i,l}\cdot {\overline{Q}}_i\cdot Q_{l,k}-{}^{2b}\overline{\beta }_{i,l}\cdot {\overline{Q}}_l\cdot Q_{i,k}\right]}-{\displaystyle \frac{1}{2}}\cdot {}^3\overline{\beta }_{l,l}\cdot {\overline{Q}}_l\cdot Q_{l,k}\\ -Q_{l,k}\cdot {\displaystyle \sum _{i=l+1}^m{}^4\overline{\beta }_{i,l}\cdot {\overline{Q}}_i}+{\overline{S}}_{l,k}-{\overline{R}}_{l,k}\cdot Q_{l,k},\end{array}$

where

$Q_{l,k}$: mass concentration of aerosol particles of component k in section l,

${\overline{Q}}_l={\displaystyle \sum _{k=1}^{n_{comp}}{Q}_{l,k}},$

${}^{1a}\overline{\beta }_{i,j,l},{}^{1b}\overline{\beta }_{i,j,l}$: coagulation rate of the aerosol particles in sections lower than l forming an aerosol particle in section l,

n_comp: number of components of the aerosol particles considered in the multicomponent sectional equations,

${}^{2a}\overline{\beta }_{i,j,l}$: coagulation rate of the aerosol particles in section l and sections lower than l forming an aerosol particle larger than those in section l,

${}^{2b}\overline{\beta }_{i,j,l}$: coagulation rate of the aerosol particles in section l and sections lower than l forming an aerosol particle remaining section l,

${}^3\overline{\beta }_{i,l}$: intra-sectional coagulation rate of the aerosol particles in section l forming an aerosol particle larger than those in section l,

${}^4\overline{\beta }_{i,l}$: coagulation rate of the aerosol particles in section l forming an aerosol particle larger than those in section l,

${\overline{S}}_{l,k}$: injection rates of the aerosol particles in section l,

${\overline{R}}_{l,k}$: removal rates of the aerosol particles in section l via deposition.

The coefficients on the coagulation rate, injection rate, and removal rates are section-averaged values calculated via numerical integrations such as an n-Gaussian quadrature rule. For example, in the case of ${}^{1a}\overline{\beta }_{i,j,l}$, and ${\overline{R}}_{l,k}$,

(2)$\begin{array}{l}{}^{1a}\overline{\beta }_{i,i-1,l}={\displaystyle \underset{x_{i-1}}{\overset{x_i}{\int }}{\displaystyle \underset{f_k\left(v_{l-1}-v\right)}{\overset{x_{l-1}}{\int }}\frac{u\cdot \beta \left(u,v\right)}{u\cdot v\cdot \left(x_i-x_{i-1}\right)\cdot \left(x_{l-1}-x_{l-2}\right)}dydx}}\\ ={\displaystyle \sum _{k_1=1}^{n_o}{W}_{k1}}\cdot ⟨{\displaystyle \sum _{k_2=1}^{n_i}\left[X_{k2}\frac{u_{k_2}\cdot \beta \left(u_{k_2},v_{k_1}\right)}{u_{k_2}\cdot v_{k_1}\cdot \left(x_i-x_{i-1}\right)\cdot \left(x_{l-1}-x_{l-2}\right)}\right]}⟩,\end{array}$

(3)${\overline{R}}_{l,k}={\displaystyle \underset{u_{l-1}}{\overset{u_l}{\int }}\frac{R\left(u,t\right)}{\left(x_l-x_{l-1}\right)}dx},$

where

u_l, v_l: mass of a single aerosol particle in section l,

(4)$\beta \left(u,v\right)={\beta }_{brwn}\left(u,v\right)+{\beta }_{grav}\left(u,v\right)+{\epsilon }_0\cdot \sqrt{{\beta }_{turb-sh}\left(u,v\right)+{\beta }_{turb-in}\left(u,v\right)},$

(5)${\beta }_{brwn}\left(u,v\right)={\displaystyle \frac{2\cdot \pi \cdot {\mathrm{\Gamma }}_u\cdot {\mathrm{\Gamma }}_v\cdot \left({\gamma }_u\cdot x_u+{\gamma }_v\cdot x_v\right)}{F_{brwn}}}:\mathrm{Brownian}\mathrm{coagulation}\mathrm{kernel},$

(6)$\begin{array}{ll}{\beta }_{grav}\left(u,v\right)=& {\epsilon }_0\cdot {\displaystyle \frac{\pi }{4}}\cdot {\left({\gamma }_u\cdot x_u+{\gamma }_v\cdot x_v\right)}^2\\ & \times {\left|V_T\left(u\right)-V_T\left(v\right)\right|}^2\end{array}:\mathrm{gravitational}\mathrm{coagulation}\mathrm{kernel},$

(7)${\beta }_{Turb-in}\left(u,v\right)=\sqrt{{\displaystyle \frac{\pi \cdot {\in }_T}{120\cdot \upsilon }}}\cdot {\left({\gamma }_u\cdot x_u+{\gamma }_v\cdot x_v\right)}^3:\mathrm{turbulent}\mathrm{inertial}\mathrm{coagulation}\mathrm{kernel},$

(8)$\begin{array}{l}{\beta }_{Turb-sh}\left(u,v\right)={\displaystyle \frac{0.04029{\rho }_g^{1/4}{\in }_T^{3/4}}{{\eta }^{5/4}}}\\ \times \left({\gamma }_u\cdot x_u+{\gamma }_v\cdot x_v\right)\\ \times \left\{{\displaystyle \frac{{\rho }_u\cdot C_u\cdot x_u^2}{{\chi }_u}}-{\displaystyle \frac{{\rho }_v\cdot C_v\cdot x_v^2}{{\chi }_v}}\right\}\end{array}:\mathrm{turbulent}\mathrm{shear}\mathrm{coagulation}\mathrm{kernel},$

(9)$R\left(u,t\right)=U\left(u,t\right)+V_T\left(u\right)+Diffu\left(u,t\right),$

(10)$U\left(u,t\right)={\displaystyle \frac{2\cdot \eta \cdot C_u\cdot \left({\alpha }_T\cdot K{n}_u+k_f/k_s\right)\cdot \nabla T}{2\cdot \chi \cdot {\rho }_g\cdot T_{\infty }\cdot \left(1+3\cdot C_u\cdot K{n}_u\right)\cdot \left(1+2\cdot {\alpha }_T\cdot K{n}_u+2\cdot k_f/k_s\right)}}:\mathrm{thermophoresis},$

(11)$V_T\left(u\right)={\displaystyle \frac{{\rho }_u\cdot g\cdot x_u^2\cdot C_u}{18\cdot \eta \cdot {\chi }_u}}:\mathrm{gravitational}\mathrm{settling},$

(12)$Diffu\left(u\right)={\displaystyle \frac{k\cdot T_{\infty }\cdot C_u}{3\pi \cdot \eta \cdot {\Delta }_{diff}\cdot x_u}}:\mathrm{diffusiophoresis},$

and the coefficients and nomenclatures used in Equations (2)–(12) are described in the nomenclature section.

Meanwhile, the Mason equations analyze hygroscopic growth that occurs as soluble aerosol particles can adsorb moisture from the atmosphere. Since the hygroscopic growth leads to an increase in the gravitational settling rate, it is an important phenomena in the aerosol dynamics. The Mason equations are expressed with respect to change of the mean radius of aerosol particles of component k, r_k, as shown in the following:

(13)${\displaystyle \frac{d{r}_k}{dt}}={\displaystyle \frac{1}{r_k}}\cdot {\displaystyle \frac{\left(S-S_{r,k}\right)}{a_k+b_k}},$

where

S: saturation ratio, in other words, relative humidity in the containment,

S_r,k: effective saturation ratio at the surface of the aerosol particles expressed as

(14)$S_{r,k}=A_{r,k}\cdot \exp \left({\displaystyle \frac{2\cdot M_w\cdot \sigma }{R\cdot T_{\infty }\cdot {\rho }_w\cdot r_k}}\right),$

a_k: thermal conduction of the latent heat associated with condensation from particles to the atmosphere, expressed as

(15)$a_k=\left({\displaystyle \frac{\Delta h_f^2\cdot M_w\cdot {\rho }_w}{R\cdot T_{\infty }\cdot k_a^{*}}}\right),$

b_k: diffusion of water vapor from the atmosphere to the particle surface of component k expressed as

(16)$b_k=\left({\displaystyle \frac{R\cdot T_{\infty }\cdot {\rho }_w}{D_{v,k}^{*}\cdot M_w\cdot p_{sat}\left(T_{\infty }\right)}}\right),$

$D_{v,k}^{*}$: effective vapor diffusion coefficient for component k expressed as

(17)$D_{v,k}^{*}={\displaystyle \frac{D_v}{\left(\frac{r_k}{r_k+{\Delta }_V}\right)+\left(\frac{D_v}{r_k\cdot {\alpha }_c}\cdot {\left(\frac{2\cdot \pi \cdot M_w}{R\cdot T_{\infty }}\right)}^{\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$2$}\right.}\right)}},$

$k_a^{*}$: effective thermal conductivity of atmosphere in the containment expressed as

(18)$k_a^{*}={\displaystyle \frac{k_a}{\left(\frac{r_k}{r_k+{\Delta }_T}\right)+\left(\frac{k_a}{r_k\cdot {\alpha }_T\cdot {\rho }_a\cdot c_{p,a}}\cdot \left(\frac{2\cdot \pi \cdot M_a}{R\cdot T_{\infty }}\right)\right)}}.$

and the coefficients and nomenclatures used in Equations (14)–(18) are described in the nomenclature section.

2.2. Summary of Coupling Scheme of Multicomponent Sectional Equations and Mason Equations for Hygroscopic Growth in I-COSTA

In order to couple the multicomponent sectional equations and the Mason equations, it is necessary to convert the radius change of the aerosol particles obtained from the Mason equations into the change of mass concentration in each section due to hygroscopic growth. This is conducted via interpolation of the distribution of mass concentration within a section and the solution of Mason equations. In I-COSTA, the aforementioned interpolation is completed based on that used in MELCOR [23]. Contrary to the interpolation in the MELCOR code, this interpolation is performed for each component in the section to obtain the mass distribution in each section which depends on the types of the component as well as the size of the particles. For the interpolation, the authors assume that the distribution of mass concentrations for all the sections to be a continuous and smooth function. However, due to the constraint in the multicomponent sectional equation, i.e., v_{i + 1} ≥ 2v_i, the interpolation function would be based on linear interpolation for the logarithmic values of the particle mass for proper construction of the interpolation. Therefore, the distribution of mass concentrations within a section is assumed as follows:

(19)${\displaystyle \frac{d{Q}_{l,k}}{dv}}=c_{l,k}\cdot v^{slo{p}_{l,k}},$

where

v: mass of a single particle,

slop_l,k: slope of the mass concentration function for component k in section l, defined as

(20)$slo{p}_{l,k}={\displaystyle \frac{\ln \left(\frac{Q_{l+1,k}-Q_{l,k}}{v_{l+1,k}-v_{l,k}}\right)-\ln \left(\frac{Q_{l,k}-Q_{l-1,k}}{v_{l,k}-v_{l-1,k}}\right)}{\ln \sqrt{v_{l+1,k}\cdot v_{l,k}}-\ln \sqrt{v_{l,k}\cdot v_{l-1,k}}}}.$

With the results of interpolation, the mass concentrations of aerosol particles of component k remain in section l, and those of particles with a size that increases to another section, e.g., l + 1, can be obtained as follows:

(21)$Q_{l,k}^{*}\left(t_0+\Delta t\right)=f{r}_{l,k}\cdot Q_{l,k}\left(t_0\right),$

(22)$Q_{l+1,k}^{*}\left(t_0+\Delta t\right)=t{r}_{l,k}\cdot Q_{l,k}\left(t_0\right),$

where

$Q_{l,k}^{*}\left(t_0+\Delta t\right)$: mass concentration of aerosol particles of component k that remain in section l after hygroscopic growth for the time step of Δt,

$Q_{l+1,k}^{*}\left(t_0+\Delta t\right)$: mass concentration of aerosol particles of component k that grow in size to section l + 1 after hygroscopic growth for the time step of Δt,

fr_l,k: fraction of the aerosol particles of component k in section l that remain in section l after hygroscopic growth, which is expressed as

(23)$f{r}_{l,k}={\displaystyle \frac{{\left(v_{l,2,k}\right)}^{slo{p}_{l,k}+1}-{\left(v_{l,k}^{*}\right)}^{slo{p}_{l,k}+1}}{{\left(v_{l+1,k}^{*}\right)}^{slo{p}_{l,k}+1}-{\left(v_{l,k}\right)}^{slo{p}_{l,k}+1}}}={\displaystyle \frac{Q_{l,2,k}}{Q_{l,k}}},$

tr_l,k: fraction of the aerosol particles of component k in section l that grow up to section l + 1 after hygroscopic growth, which is expressed as

(24)$t{r}_{l,k}=1-f{r}_{l,k}.$

The relationship among the variables delineated in Equations (19)–(24) is shown in Figure 1.

With the results of fr_l,k and tr_l,k, an elemental transition matrix that satisfies the following can be constructed:

(25)$\exp \left(A_{l,k}\cdot \Delta t\right)\cdot {\overrightarrow{Q}}_{l,k}^{*}\left(t_0\right)=\left[\begin{array}{cc}f{r}_{l,k}& 0\\ t{r}_{l,k}& 0\end{array}\right]\cdot \left[\begin{array}{c}{Q}_{l,k}\left(t_0\right)\\ 0\end{array}\right],$

where

(26)$A_{l,k}=\left[\begin{array}{cc}-{\lambda }_{l,k}& 0\\ {\lambda }_{l,k}& 0\end{array}\right],$

λ_l,k: transition rate of aerosol particles of component k in section l, defined as

(27)${\lambda }_{l,k}=-{\displaystyle \frac{\ln \left(f{r}_{l,k}\right)}{\Delta t}}.$

A global transition rate matrix of hygroscopic growth is obtained with elemental transition rate matrices for all components and sections in the analysis as the follows:

(28)${\displaystyle \sum _{l=1}^{n_{bin}}{\displaystyle \sum _{k=1}^{n_{comp}}{A}_{l,k}}=A},$

where

n_bin: number of sections of the aerosol particles considered in the multicomponent sectional equations.

The process of obtaining the global transition rate matrix can be expressed as illustrated in Figure 2.

In I-COSTA, the following non-homogeneous system of equations is implemented as a way of coupling the multicomponent sectional equations and Mason equations rigorously via the global transition rate matrix:

(29)${\displaystyle \frac{d\overrightarrow{Q}}{dt}}=A\cdot \overrightarrow{Q}+P\left(\overrightarrow{Q}\right),$

where

$P\left(\overrightarrow{Q}\right)$: mass concentration of aerosol particles of component k that remain in section l after hygroscopic growth for the time step of Δt.

2.3. Comparison of the Coupling Scheme in I-COSTA and That in the Conventional Codes

In order to compare the coupling scheme in the I-COSTA and that in conventional code, e.g., MELCOR, comparative calculations are performed on the aerosol behaviors with the hygroscopic growth. In I-COSTA, the coupling scheme is applied to the calculation as discussed in Section 2.2. In the MELCOR code, Mason equations are solved first and then the solutions are used for the initial condition of multicomponent sectional equations for the current time step, since it is assumed that the time scale in hygroscopic growth is much shorter than that of phenomena dealt with in the multicomponent sectional equations in MELCOR [23]. The difference between the two schemes is shown in Figure 3. It should be noted that the main difference is in the use of the transition matrix.

The changes in airborne concentration using a scenario were compared. This scenario involved aerosol particles with intermediate solubility, i.e., with a Van’t Hoff factor of approximately 1. Geometric and thermophysical conditions of the problem are shown in Table 1. In the problem, behaviors of the aerosols are analyzed for 200 s.

The saturation ratio and the effective saturation ratio on the aerosol particles are assumed to be 1.00005 and 1. Temperature-dependent chemical and thermophysical properties used in the calculations are taken from MELCOR 2.2 [23]. The two codes use 1.0 × 10⁻² s as a size of the time step (Δt). Distribution of mass concentrations of the aerosols at the initial condition is shown in Figure 4. Distributions of mass concentrations at t = 0.1, 10, 100, and 200 s are also shown in Figure 4.

As shown in Figure 4, aerosol particles grow more slowly in I-COSTA than in MELCOR. For example, at t = 200 s, relative mass concentration in sections higher than 10 (d ≥ $3.05\times {10}^{-7}\mathrm{m}$) is 0.23 in MELCOR, whereas it is 0.19 in I-COSTA. The reason for this difference is that, as Figure 3 shows, MELCOR uses the solutions of the Mason equations at each time step as the initial conditions for the multicomponent sectional equations. This approach is unrealistic because aerosol particles grow through simultaneous processes of hygroscopic growth and coagulation. In contrast, I-COSTA explicitly couples the two equations and solves them simultaneously. Such differences in numerical methods can lead to variations in the predicted airborne concentrations of aerosol particles in the containment atmosphere during a severe accident—particularly for cesium-related aerosols (e.g., CsOH, CsI, and others), whose behavior strongly depends on hygroscopic growth due to their high solubility in water. It is important to note that cesium is one of the most significant nuclides, as it has a relatively long half-life and a high dose conversion factor, both of which contribute to a high exposure dose rate in the event of a large release caused by containment failure during a severe accident.

The aforementioned differences in the coupling scheme result in variations in the airborne concentrations of aerosols, as demonstrated in the KAEVER experiments referenced in the work by one of the authors [21]. In that study, I-COSTA produced more conservative results than other conventional codes when aerosol particles exhibited low solubility. However, when the particles had high solubility, the coupling scheme yielded results similar to those of the other codes. Currently, information on the solubility of iodine-related aerosol particles is not available, and thus, the coupling scheme used in I-COSTA should be taken into consideration.

2.4. Selected Parameters and Their Uncertainty Range for Importance Analysis

As discussed in the previous section, I-COSTA considers multicomponent sectional equations and Mason equations. For the derivation of the two aforementioned equations from the motion of particles in the gaseous media, numerous phenomenological parameters are used to describe characteristics of the aerosols, such as correction for the deviation from a perfect spherical shape of the particles, non-continuum regime, etc. A database for those parameters is obtained through experiments and/or theoretical derivation, and the parameters are known to have uncertain values, which induces phenomenological uncertainty in the aerosol dynamics [17,18]. In this paper, the authors have considered 12 parameters for the importance analysis. The 12 parameters are incorporated into the coagulation and deposition kernels, as shown in Equations (2)–(12), and are essential for determining the airborne concentrations of aerosol particles in the containment area. It is important to note that the accident source term for radiological consequence analysis is directly determined by the aerosol concentration at the time of containment failure. The characteristics of each parameter are discussed below:

(1) Sticking efficiency (${\epsilon }_0$) is considered in the coagulation of two particles, as shown in Equations (4), (6) and (A2) in Appendix A. It considers the effect via Van der Walls forces, changes in surface free-energies, and/or chemical reaction [16]. The computer codes on the aerosol dynamics, e.g., MAEROS [21,24], MELCOR [23], etc., have used unity as the value sticking efficiency. According to a previous study [17], the sticking efficiency has an uncertainty range of 0.1 to 1.

(2) Slip correction factor ($f_{slip}$) in Equation (A1), used in Equation (8), arises in the various expressions on the coagulation kernels; it accounts for deviation from continuum mechanics and is derived from the experimental data on spherical parameters. Most of the computer codes for the aerosol dynamics have used 1.37 as a default value. The uncertainty range of the slip correction factor is from 1.1 to 1.4, according to previous studies [17,18,19].

(3) Diffusion boundary thickness (${\Delta }_{diff}$) is considered in the Brownian diffusion used for the deposition kernel, as shown in Equation (12). According to the computer codes for aerosol dynamics, such as MAEROS [25], MELCOR [24], etc., the default value of 1.0 × 10⁻⁵ m is widely used, and its uncertainty range is from 10⁻⁶ to 10⁻⁴ m [14,15,16].

When two particles coagulate, the resulting combined particle does not form a perfect sphere, which affects the behavior of aerosols during a severe accident. In addition, most of the equations for aerosol coagulation are derived theoretically, assuming that the resulting particles are perfectly dense spherical particles. However, in reality, the particles do not have the aforementioned geometry. Therefore, in the aerosol dynamics, two correction factors are introduced to consider the deviation between theoretical derivations and the actual particle characteristics: (4) a dynamic shape factor (${\chi }_u$) and (5) a collision shape factor (${\gamma }_u$), as shown in Equations (5)–(8), (11), etc. The dynamic shape factor describes the difference in the drag force between a perfectly spherical particle and a non-spherical, porous particle. Meanwhile, the collision shape factor accounts for the difference in terms of spatial extent between a perfectly spherical particle and a non-spherical, porous particle. In the current aerosol dynamics codes, the values of the two shape factors are assumed to be unity. In this work, based on the literature, the authors assumed that the two parameters have an uncertainty range from 1 to 4 in the case of severe accident conditions [17,18,19].

(6) The thermal accommodation coefficient (α_T) used in Equations (10) and (18) accounts for the interaction between gas and aerosol particles in terms of temperature, and it is known to be strongly dependent on the structure of gas molecules [17,18,19]. The conventional codes have used unity as the value of the coefficient. From a previous study [18], the authors assume an uncertainty range from 0.5 to 1.5 for the importance analysis.

During the motion of particles along the streamline of the fluid, turbulent coagulation may occur as particles collide with each other due to differences in inertia, impaction by structure change in the medium, etc. In turbulent coagulation, the turbulent motion of the fluid is one of the key factors in analyzing the phenomenon. (7) The turbulent energy dissipation rate (${\in }_T$) is considered in the turbulent inertial and shear coagulations, as shown in Equations (7) and (8). It is the rate at which turbulent kinetic energy is converted into thermal energy affecting turbulent motion. In the conventional computer codes, a default value of 10⁻³ m²/sec³ is used. In the importance analysis, based on the literature review, a range of 5 × 10⁻⁴–1.5 × 10⁻³ is considered [17,18,19].

For the hygroscopic growth of the aerosol particles, the condensation of the vapor to the surface of the aerosol particles is affected by (8) the saturation ratio (S), (9) the effective diffusion coefficient of vapor ($D_{v,k}^{*}$), and (10) the effective thermal conductivity of the atmosphere in the containment ($k_a^{*}$), as shown in Equation (14). According to the experimental activities on the aerosol particles [26], the three parameters have an uncertainty range of $\pm$7% during the experiment. The authors have considered the aforementioned range for the importance analysis with the input data for ISP No44 [20].

The aerosol particles in the containment area during a severe accident are not composed of a single substance. They are a mixture of the vapor and various kinds of fission products. Therefore, the (11) density of aerosol particles (${\rho }_u$) has uncertainty. In the conventional computer codes [24,25], a value of 1000 kg/m³ is used as a default value, and it is known to be underestimated [17,18] since the value considers the aerosol particles to be wet particles. In the importance analysis, an uncertainty range of 1000–5000 kg/m³ is considered, based on the literature review [17,18].

For the deposition of the aerosol particles, thermophoresis, which is induced by a thermal gradient in the particle, is considered an important mechanism. In the deposition kernel, (12) the ratio of thermal conductivity of the atmosphere to that of particles ($k_f/k_s$) is used to consider the aforementioned effect, as shown in Equation (10). The conventional computer codes [24,25] have used a value of 0.05 as a default value, and an uncertainty range of 2 × 10⁻⁴~5 × 10⁻³ is considered in the importance analysis, based on the literature review [17,18]. The aforementioned phenomenological parameters are summarized in Table 2.

2.5. Summaries of the Framework for Importance Analysis

With the parameters discussed in the previous section, the authors performed an importance analysis on the aerosol behaviors in the containment area. The first step of the analysis is to sample the values of the parameters randomly within the range of uncertainty discussed in Section 2.3. The sampling is based on Latin hypercube sampling, in which the uncertainty range is divided into n non-overlapping intervals on the basis of equal probability for each parameter. The division is based on the probability density function (PDF) given by users. A cumulative distribution function (CDF) is then calculated for each parameter to sample a single value for each parameter in each interval randomly. With n different values of each parameter generated from previous procedures, n vectors are generated by mixing the order of n values. The mixing is completed by picking up the parameters randomly. Note that the dimension of the aforementioned vector is equal to the total number of parameters considered in the importance analysis. The discussed process is conducted using a stand-alone in-house code named R-SAPhe.

Sensitivity analyses on the change of aerosol concentrations are carried out via I-COSTA with the aforementioned n vectors. The correlation coefficients are then calculated using CC-SAPhe in order to show the relationship between the parameters and the figure-of-merit in the importance analyses. The correlation coefficients lie between −1 and 1; a value of −1 means that a parameter has a perfect negative linear relationship with the figure-of-merit, i.e., the parameter increases, the figure-of-merit decreases. A value of 1 means, in contrast, means that a parameter has perfect a positive linear relationship with the figure-of-merit, i.e., the parameter increases, the figure-of-merit increases. Meanwhile, a value of 0 means there is no linear relationship. According to previous studies, if the absolute value of the correlation coefficient is higher than 0.7, it has strong a linear relationship with the parameters [24]. In this analysis, four types of correlation coefficients are considered: Pearson coefficient, Spearmann coefficient, standardized regression coefficient, and partial correlation coefficient [19,25,26,27]. The procedure for the analyses is summarized in Figure 5.

3. Numerical Results

3.1. Computation Conditions on the Importance Analysis

In this section, importance analyses are conducted on the phenomenological parameters in I-COSTA. For this analysis, International Standard Problem No. 44 (ISP-44) [20] is considered. ISP-44 consists of five experiments performed using three types of aerosols—CsI, Ag, and CsOH—to examine aerosol behavior with varying solubility. The characteristics of these aerosols are presented in Table 3.

Among the five ISP-44 experiments, three—KAEVER-148, KAEVER-186, and KAEVER-187—are considered in this study. The KAEVER-148 experiment used Ag aerosols, which have the lowest solubility and exhibit the slowest hygroscopic growth. Importance analysis of this experiment highlights the sensitivity of aerosol behavior to parameters used in the general dynamic equations.

The KAEVER-186 experiment involved a combination of Ag and CsOH aerosols, with CsOH showing the highest solubility and thus the fastest hygroscopic growth. Compared to KAEVER-148, an importance analysis of KAEVER-186 better illustrates the influence of parameters related to hygroscopic growth on aerosol behavior.

The KAEVER-187 experiment used a mixture of all three aerosol types. In comparison to the previous two experiments, the importance analysis of KAEVER-187 reflects the behavior of a broader range of aerosols under conditions closely resembling those in a containment during a severe accident.

All experiments were conducted in the KAEVER facility. The geometric conditions of the KAEVER facility are summarized in Table 4.

For the sensitivity analyses on the change of airborne concentrations of the aerosol particles, I-COSTA is used as discussed in the previous section. The implicit Euler method is used as a numerical method for the time discretization. The Newton method is then used for solving the discretized equations obtained via the implicit method. In the analysis, the step doubling method is used to determine the size of the time step automatically.

In the discretized equations, section-averaged coefficients are required for the calculation of the multicomponent sectional equations, and they are obtained via numerical integrations with Gaussian quadrature. The computation conditions are summarized in Table 5.

For the importance analysis, it is necessary to determine the appropriate number of samples for the Latin hypercube sampling of phenomenological parameters in aerosol dynamics. In this study, sensitivity analyses on the number of samples were conducted by comparing the averages and standard deviations obtained from various sample sizes with those derived from the continuous distributions of the parameters. Based on this comparison, a sample size of 2000 was selected, as both the average and standard deviation from 2000 samples exhibited a relative difference of less than 10⁻⁵ compared to those from the continuous distributions, as shown in Table 6.

3.2. Importance Analyses on the Phenomenological Parameters for the Aerosol Dynamics

With the aforementioned conditions, the importance analyses are performed on the phenomenological parameters for the aerosol dynamics in I-COSTA. The results of the sensitivity analysis on the KAEVER-148 experiment are shown in Figure 6. The results for the KAEVER-186 and the KAVER-187 experiments are shown in Figure 7 and Figure 8, respectively.

As shown in Figure 6, the change of airborne concentrations of the aerosol particles obtained from the KAEVER-148 experiment is within the range of the sensitivity analysis results via I-COSTA. For the KAEVER-186 and KAEVER-187 experiments, the experimental data on the change of the airborne concentrations are also within the range of the sensitivity analysis results. In terms of the maximum airborne concentrations, the range of sensitivity for the KAEVER-148 experiment is $\pm$30% of the values obtained with the reference phenomenological parameters shown in Table 1. Similarly, the range of sensitivity analyses obtained from the other two experiments are $\pm$40%, respectively. According to the experimental conditions [20], the range of uncertainty in the measurement of aerosol concentration is known to be $\pm$7% of the measured values for most of the experiments in ISP44. This means that the range of sensitivity analyses is much larger than the uncertainty range of the measurement during the experiment. Therefore, we can conclude that the uncertainty in the phenomenological parameters in the aerosol dynamics cause larger uncertainties than those from measurement in terms of the airborne aerosol concentrations. It is thus necessary to analyze the aforementioned parameters in order to enhance the accuracy in the change of aerosol concentrations. Figure 9 shows the correlation coefficients of the phenomenological parameters for the maximum airborne concentrations, and for the airborne concentrations at the end of the experiments obtained from the KAEVER-148 experiment. The correlation coefficients obtained from the other experiments are shown in Figure 10 and Figure 11 for the maximum airborne concentrations and the concentrations at the end of each experiment, respectively.

As shown in Figure 9, Figure 10 and Figure 11, the dynamic shape factor shows a strong linear relationship with the airborne concentrations. The airborne concentrations of the aerosol particles decrease as the dynamic shape factor increases. In the aerosol dynamics, the dynamic shape factor is considered in the rate of coagulation of two particles and deposition of the aerosol particles. In the calculation of the coagulation rate, as shown in Equation (4), the dynamic shape factor is used in Brownian, gravitational, and turbulent coagulation kernels, and it is inversely proportional to the rates since the drag force on the non-spherical particle is greater than that on the perfectly spherical particle. In other words, the higher the value of the dynamic shape factor, the smaller the values are in the coagulation kernels. This leads to slower coagulation rates of the aerosol particles, which causes the mass of the aerosol particles in the atmosphere to increase more slowly; this is a key factor in determining the deposition rate.

The dynamic shape factor is also used in the deposition kernels, as shown in Equation (9). As with the coagulation kernels, the dynamic shape factor is inversely proportional to the kernels for the gravitational settling and turbulent deposition, which means that the higher dynamic shape factor, the slower the deposition rate of the aerosol particles. The combined effects of the dynamic shape factor on the coagulation and deposition kernels implies that the dynamic shape factor should have a strong linear relationship with the airborne concentrations of the aerosol particles.

The saturation ratio also shows a strong linear relationship with the airborne concentrations of the aerosol particles. In I-COSTA, the saturation ratio is used in the Mason equations, which determine the hygroscopic rate shown in Equation (13). A higher saturation ratio leads to a faster growth rate of the aerosol particles. A fast growth rate is related to a fast increase of the aerosol particle mass, which leads to a faster particle deposition rate. The saturation ratio, therefore, must show a negative linear relationship with the airborne mass of the aerosol particles.

4. Summary and Conclusions

In this study, importance analyses on the twelve phenomenological parameters for the aerosol dynamics were performed with I-COSTA. The twelve parameters were selected from a review of the literature on the uncertainty of the phenomenological parameters [17,18,19] and uncertainty during measurement [26]. In the importance analyses, the parameters were sampled via Latin hypercube sampling. Sensitivity analyses were then performed on the change of the airborne concentrations of KAEVER-148, KAEVER-186, and KAEVER-187 experiments, which are included in International Standard Problem No. 44 [20]. With the results of the sensitivity analyses, the correlation coefficients between the phenomenological parameters and the maximum airborne concentrations of the aerosol particles were obtained. The correlation coefficients between the aforementioned parameters and the airborne concentrations were also obtained. The main findings from the importance analyses on the phenomenological parameters are as follows:

The range of the phenomenological parameters taken from previous studies leads to $\pm$50% of the airborne concentration of the aerosol particles obtained from experimental values.

The range of sensitivity of the maximum airborne aerosol concentrations is seven times the uncertainty range of the measurement during experiment, i.e., the uncertainty range of the measurement is known to be $\pm$7% of the concentrations [17].

The larger range of sensitivity compared to the range of measurement appears to be caused by the uncertainty range of the phenomenological parameters. In order to make the accuracy of the aerosol dynamics model comparable to that of the experimental measurement, the uncertainty of the parameters should be reduced.

From the correlation coefficients between the phenomenological parameters and the airborne concentrations, the dynamic shape factor and the saturation ratio were found to be important parameters with respect to the airborne concentrations of the aerosol particles.

From the aforementioned main findings, the dynamic shape factor should be refined in terms of various severe accident conditions, such as pressure, temperature, dose, etc., to gain a better understanding of the aerosol concentrations. Since the dynamic shape factor is known to depend heavily on the pressure of the gas, which would be linked to the difference in the drag force of non-spherical particles and that of perfectly spherical particles [28,29,30,31,32], it is necessary to construct a database for severe accident conditions so that the uncertainty range of the source term during a severe accident analysis can be reduced.

From the perspective of safety analysis—specifically nuclear safety regulation—the main findings of this study can be used to evaluate uncertainties in the accident source term [33] arising from limited information on the shape as well as the physical and chemical characteristics of aerosol particles released during a severe accident. In the licensing review process for radiological consequences, these uncertainties can help define the uncertainty range of the nuclear accident source term due to aerosol behavior. This, in turn, supports more informed and reasonable regulatory decisions regarding radiological consequences.

Throughout the long history of source term research, the effective dose of Cs type aerosol particles such as CsI and CsOH has been highlighted in terms of the public [17,34], as these particles are long-lived isotopes with high affinity to human organs [34]. As the delayed release of the fission product becomes an important area of research, given the findings from the Fukushima accident [35] and the international cooperative research project [9], it is also necessary to specify the shape and the density of the aerosol particles formed from chemical species, such as iodine oxide in the gaseous phase, etc., especially for coupling the aerosol dynamics code and the iodine chemistry code to evaluate the delayed release of the fission products. Note that, as discussed, the shape and density are known to be the main parameters that determine the dynamic shape factor.

As discussed in Section 1, this aspect may become one of the most critical issues during the licensing review process, as SMR licensees seek to minimize the size of the EPZ to enhance the value of SMRs by enabling construction closer to high population density areas. Minimizing the EPZ requires a best-estimate evaluation of the accident source term, for which detailed analyses of aerosol behavior are essential. In such analyses, it is necessary to develop a database on the characteristics of aerosol particles that could be released from the reactor. The main findings of this study will support the prioritization of parameters needed to construct that database.

Data Availability Statement

Data can be provided upon request.

Conflicts of Interest

The authors declare no conflicts of interest.

Nomenclature

The following nomenclatures are used in this manuscript:

A	Global transition rate matrix
Al,k	Elemental transition rate matrix for aerosol particles of component k in section l
ak	Thermal conduction of the latent heat associated with condensation from particles to the atmosphere
Ar,k	Activity of the aerosol particles
bk	Diffusion of water vapor from the atmosphere to the particle surface of component k
$D_{v,k}$	Vapor diffusion coefficient for component k
$D_{v,k}^{*}$	Effective vapor diffusion coefficient for component k
$f_{slip}$	Slip correction factor
$f_k(u)$	Radius function of the aerosol particles for component k with respect to the mass of particles u
frl,k	Fraction of the aerosol particles of component k in section l that remain in section l after hygroscopic growth
$k_a$	Thermal conductivity of atmosphere in the containment
$k_a^{*}$	Effective thermal conductivity of atmosphere in the containment
Mw	Molecular weight of the water
nbin	Number of sections of the aerosol particles considered in the multicomponent sectional equations
ncomp	Number of components of the aerosol particles considered in the multicomponent sectional equations
no, ni,	Number of Gaussian quadrature sets
$P\left(\overrightarrow{Q}\right)$	Vector form of the multicomponent sectional equations of the aerosol particles
$Q_{l,k}$	Mass concentration of aerosol particles of component k in section l
$Q_{l,k}^{*}\left(t_0+\Delta t\right)$	Mass concentration of aerosol particles of component k that remain section l after hygroscopic growth for the time step of Δt
$Q_{l+1,k}^{*}\left(t_0+\Delta t\right)$	Mass concentration of aerosol particles of component k that grow up to section l + 1 after hygroscopic growth for the time step of Δt
R	Gas constant
${\overline{R}}_{l,k}$	Removal rates of the aerosol particles in section l via deposition,
S	Saturation ratio, in other words, relative humidity in the containment
${\overline{S}}_{l,k}$	Injection rates of the aerosol particles in section l
Sr,k	Effective saturation ratio at the surface of the aerosol particles
slopl,k	Slope of the mass concentration function for component k in section l
T∞	Temperature of the atmosphere in the containment
trl,k	Fraction of the aerosol particles of component k in section l that grow up to section l + 1 after hygroscopic growth
U(u,t)	Deposition via thermophoresis,
VT(u)	Deposition via gravitational settling
Diffu(u)	Deposition via diffusiophoresis,
$u_{k_2}$	Abscissas for Gaussian quadrature
ul	Mass of a single aerosol particle in section l
v	Mass of a single particle
vl	Mass of a single aerosol particle in section l
$W_{k_1},X_{k2}$	Weights for Gaussian quadrature
xl	Radius of a single aerosol particle in section l
αc	Condensation coefficient
αT	Thermal accommodation coefficient
$\beta$(u,v)	Kernel for Brownian, gravitational, and turbulent shear coagulations for particles having masses of u and v
${\beta }_{brwn}$(u,v)	Brownian coagulation kernel,
${\beta }_{grav}$(u,v)	Gravitational coagulation kernel
${\beta }_{Turb-in}$(u,v)	Turbulent inertial coagulation kernel
${\beta }_{Turb-shear}$(u,v)	Turbulent shear coagulation kernel
${}^{1a}\overline{\beta }_{i,j,l},{}^{1b}\overline{\beta }_{i,j,l}$	Coagulation rate of the aerosol particles in sections lower than l forming an aerosol particle in section l
${}^{2b}\overline{\beta }_{i,j,l}$	Coagulation rate of the aerosol particles in section l and sections lower than l forming an aerosol particle larger than those in section l
${}^3\overline{\beta }_{i,l}$	Intra-sectional coagulation rate of the aerosol particles in section l forming an aerosol particle larger than those in section l
${}^4\overline{\beta }_{i,l}$	Coagulation rate of the aerosol particles in section l forming an aerosol particle larger than those in section l
${\Delta }_{diff}$	Diffusion boundary thickness
${\gamma }_u$	Collision shape factor of an aerosol particle with mass u
${\epsilon }_0$	Sticking efficiency of the two aerosol particles
$\eta$	Gas viscosity
$\lambda$	Mean free path of an aerosol particle,
λl,k	Transition rate of aerosol particles of component k in section l
ρw	Density of the water
σ	Surface tension of water
${\chi }_u$	Dynamic shape factor of an aerosol particle with mass u
${\in }_T$	Turbulent energy dissipation rate

Footnotes

Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Figures and Tables

Figure 1 Concept of change of mass concentration due to hygroscopic growth.

Figure 2 Construction of a global transition rate matrix of hygroscopic growth [19].

Figure 3 Comparison of the coupling scheme in MELCOR and I-COSTA used in the aerosol dynamics. (1) Refer to Equations (21) and (22). (2) Refer to Equation (1). (3) Refer to Equation (29).

Figure 4 Distribution of mass concentrations for the test problem; (a) at t = 0; (b) at t = 0.1; (c) at t = 10; (d) at t = 100; (e) at t = 200.

Figure 5 Importance analysis scheme for the phenomenological parameters with I-COSTA.

Figure 6 Sensitivity analysis on the KAEVER-148 experiment with sampling of phenomenological parameters.

Figure 7 Sensitivity analysis on the KAEVER-186 experiment with sampling of phenomenological parameters.

Figure 8 Sensitivity analysis on the KAEVER-187 experiment with sampling of phenomenological parameters.

Figure 9 Correlation coefficients for the KAEVER-148 experiments. (a) With the maximum airborne concentrations, (b) at the end of experiment.

Figure 10 Correlation coefficients for the KAEVER-186 experiments. (a) With the maximum airborne concentrations, (b) at the end of experiment.

Figure 11 Correlation coefficients for the KAEVER-187 experiments. (a) With the maximum airborne concentrations, (b) at the end of experiment.

Table 1

Geometric and thermophysical conditions for the test problems.

Parameter	Value
Van’t Hoff factor	1.00
Volume of chamber [m3]	2.86
Temperature [K]	298
Saturation ratio	1.00005
Effective saturation ratio	1.0
Pressure [Pa]	$1.01\times {10}^5$
Density of aerosol particle [kg/m3]	1000
Number of sections (nsection)	20

Table 2

Phenomenological parameters in the analysis.

Parameters[Units]	Description	Mean	Min.	Max.
Sticking efficiency (${\epsilon }_0$)(P_stick)	Effect of Van der Walls force	1.0	0.1	1.0
Slip correction factor ($f_{slip}$)	Deviation from continuum mechanics	1.37	1.1	1.3
Diffusion boundary thickness (${\Delta }_{diff}$)	Deposition via Brownian diffusion	10−5	10−6	10−4
Dynamic shape factor (${\chi }_u$)	Drag force of non-spherical particles	1.0	1.0	4.0
Collision shape factor (${\gamma }_u$)	Spatial extent of non-spherical particle	1.0	1.0	4.0
Thermal accommodation coefficient (αT)	Interaction between gas and particle in terms of temperature	1.0	0.5	1.5
Turbulent energy dissipation rate (${\in }_T$)	Rate at which turbulent kinetic energy is converted into thermal energy	10−3	5 × 10−4	1.5 × 10−3
Saturation ratio (S)	Humidity	Experimental condition	95%	105%
Effective diffusion coefficient of vapor ($D_{v,k}^{*}$)	Diffusion coefficient for vapor	Experimental condition	95%	105%
Effective thermal conductivity of the atmosphere in the containment ($k_a^{*}$)	Thermal conductivity on the particle	Experimental condition	93%	107%
Density of aerosol particles (${\rho }_u$)	Mixed particle density	1000.0	1000.0	5000.0
Ratio of thermal conductivity of the atmosphere to that of particles ($k_f/k_s$)	Ratio of thermal conductivity of gas to that of particle	0.037	0.0002	0.055

Table 3

Types of aerosol particles considered in KAEVER experiments.

Experiments	Aerosols	Solubility(Van’t Hoff Factor)
KAEVER-148	Ag	0.08
KAEVER-186	Ag	0.08
	CsOH	2.0
KAEVER-187	Ag	0.08
	CsI	1.68
	CsOH	2.0

Table 4

Geometric conditions of KAEVER experiments.

Parameter		Value
Volume [m3]		1.00
Deposition surface [m3]	Ceiling	5.677
	Vertical wall	26.2665
	Floor	5.6667

Table 5

Computation conditions of I-COSTTA.

Parameter		Value
Numerical method for solving Equation (29)	Time step size (Δt, s)	Determined by step doubling method (0.05~0.2)
	Time discretization	Implicit Euler method
	Method for solving time discretized equations	Newton method
Section-averaged coefficients	Quadrature	Gaussian quadrature
	Number of sets of weights and abscissas	20
Number of sections		20

Table 6

Sensitivity studies on the average and standard deviation for various number of samplings in Latin hypercube sampling.

Case	Average (1)	Standard Deviation (2)
10	2.45 × 10−3	5.79 × 10−2
50	1.54 × 10−3	4.68 × 10−3
100	5.10 × 10−4	1.94 × 10−3
500	6 × 10−5	1.19 × 10−4
1000	2.46 × 10−5	7.96 × 10−5
2000	1.27 × 10−6	1.42 × 10−6

(1) $\underset{i=param}{\max }\left\{abs\left({\displaystyle \frac{{\overline{\xi }}_i-{\widehat{\xi }}_i}{{\widehat{\xi }}_i}}\right)\right\},$ (2) $\underset{i=param}{\max }\left\{abs\left({\displaystyle \frac{{\overline{\xi }}_i-{\widehat{\xi }}_i}{{\widehat{\xi }}_i}}\right)\right\}$, where $\overline{\xi }$: average of the phenomenological parameter obtained from Latin Hypercube sampling, ${\widehat{\xi }}_i$: average of the phenomenological parameter obtained from continuous distribution, ${\overline{\sigma }}_i$: standard deviation of the phenomenological parameter obtained from Latin Hypercube sampling, ${\widehat{\sigma }}_i$: standard deviation of the phenomenological parameter obtained from continuous distribution.

Appendix A. Description on the Coefficients Used in Equations (1)–(29)

The followings are the coefficients used in Equations (1)–(29) of the paper(A1)$$C_u=1+K{n}_u\cdot \left(f_{slip}+0.4\cdot \exp \left({\displaystyle \raisebox{1ex}{$-1.1$}\!\left/ \!\raisebox{-1ex}{$K{n}_u$}\right.}\right)\right):\mathrm{Cunningham}\mathrm{slip}\mathrm{correction}\mathrm{factor},$$ (A2)$$F_{brwn}={\displaystyle \frac{x_u+x_v}{x_u+x_v+2\cdot g_{u,v}}}+{\displaystyle \frac{8\cdot \left({\mathrm{\Gamma }}_u+{\mathrm{\Gamma }}_v\right)}{{\overline{V}}_{u,v}\cdot \left(x_u+x_v\right)\cdot {\epsilon }_0}},$$(A3)$$g_{u,v}=\sqrt{\left(g_u^2+g_v^2\right)},$$(A4)$${\overline{V}}_{u,v}=\sqrt{\left(V_u^2+V_v^2\right)},$$(A5)$$g_u={\displaystyle \frac{1}{3\cdot x_u\cdot l_u}}\left\{{\left(x_u+l_u\right)}^3-{\left(x_u^2+l_u^2\right)}^{3/2}\right\}-x_u,$$(A6)$$l_u={\displaystyle \frac{8\cdot {\mathrm{\Gamma }}_u}{\pi \cdot V_u}},$$(A7)$$V_u=\sqrt{{\displaystyle \frac{8\cdot k\cdot T_{\infty }}{\pi \cdot u}}},$$ (A8)$${\mathrm{\Gamma }}_u={\displaystyle \frac{k\cdot T_{\infty }}{3\cdot \pi \cdot x_u\cdot \eta \cdot {\chi }_u}}\cdot C_u,$$ (A9)$$K{n}_u={\displaystyle \frac{2\lambda }{x_u}}:\mathrm{Cunningham}\mathrm{slip}\mathrm{correction}\mathrm{factor},$$$f_k(u)$: radius function of the aerosol particles for component k with respect to the mass of particles expressed as(A10)$$f_k\left(u\right)=x,$$$W_{k_1},X_{k2}$: weights for Gaussian quadrature, satisfying(A11)$${\displaystyle \sum _{k_1=1}^{n_o}{W}_{k_1}}=x_i-x_{i-1},$$(A12)$${\displaystyle \sum _{k_2=1}^{n_i}{X}_{k_2}}=x_{l-1}-x_{i-2},$$$u_{k_2},v_{k_1}$: weights for Gaussian quadrature, satisfying(A13)$$u_{k_2}\in \left[x_{i-1},x_i\right],$$(A14)$$v_{k_1}\in \left[f_k\left(v_{i-1}-v_{k_1}\right),x_{l-1}\right].$$

References

1. Hart, M.; Herviou, K.; Wu, J.; Taehee, K.; Shapovalov, A.; Thomas, K.; Wallace, D. SMR Regulators’ Forum, Pilot Project Report: Report from Working Group on Emergency Planning Zone; International Atomic Energy Agency (IAEA): Wien, Austria, 2018; Available online: https://www.iaea.org/sites/default/files/18/01/smr-rf-report-appendix-iv-29012018.pdf (accessed on 31 May 2025).

2. Haste, T.; Payot, F.; Bottomley, P.D.W. Transport and deposition in the Phébus FP circuit. Ann. Nucl. Energy; 2013; 61, pp. 102-121. [DOI: https://dx.doi.org/10.1016/j.anucene.2012.10.032]

3. Mascari, F.; Bersano, A.; Massone, M.; Agnello, G.; Coindreau, O.; Beck, S.; Tiborcz, L.; Paci, S.; Angelucci, M.; Herranz, L.E. . Main outcomes of the Phebus FPT1 uncertainty and sensitivity analysis in the EU-MUSA project. Ann. Nucl. Energy; 2024; 196, 110205. [DOI: https://dx.doi.org/10.1016/j.anucene.2023.110205]

4. Yu, M. Main outputs from the OECD/NEA ARC-F Project. Proceedings of the 20th International Topical Meeting on Nuclear Reactor Thermal Hydraulics (NURETH-20); Washington, DC, USA, 20–25 August 2023.

5. Doi, T.; Masumoto, K.; Toyoda, A.; Tanaka, A.; Shibata, Y.; Hirose, K. Anthropogenic radionuclides in the atmosphere observed at Tsukuba: Characteristics of the radionuclides derived from Fukushima. J. Environ. Radioact.; 2013; 122, pp. 55-62. [DOI: https://dx.doi.org/10.1016/j.jenvrad.2013.02.001] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/23542231]

6. Grégoire, A.C.; Le Fessant, E.; Cantrel, L. ESTER Program—State of the Art on Revaporization from RCS Surfaces; IRSN No 2021-00282 IRSN: Fontenay-aux-Roses, France, 2021.

7. Katata, G.; Chino, M.; Kobayashi, T.; Terada, H.; Ota, M.; Nagai, H.; Kajino, M.; Draxler, R.; Hort, M.C.; Malo, A. . Detailed source term estimation of the atmospheric release for the Fukushima Daiichi Nuclear Power Station accident by coupling simulations of an atmospheric dispersion model with an improved deposition scheme and oceanic dispersion model. Atmos. Chem. Phys.; 2015; 15, pp. 1029-1070. [DOI: https://dx.doi.org/10.5194/acp-15-1029-2015]

8. Terada, H.; Katata, G.; Chino, M.; Nagai, H. Atmospheric discharge and dispersion of radionuclides during the Fukushima Dai-ichi Nuclear Power Plant accident. Part II: Verification of the source term and analysis of regional-scale atmospheric dispersion. J. Environ. Radioact.; 2012; 112, pp. 141-154. [DOI: https://dx.doi.org/10.1016/j.jenvrad.2012.05.023] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/22721917]

9. Experiments on Source Term for Delayed Releases (ESTER) Project. 2020; Available online: http://www.oecd-nea.org/jcms/pl_58926/experiments-on-source-term-for-delayed-releases-ester-project (accessed on 31 May 2025).

10. Colombani, J.; Pascal, C.; Martinet, L.; Duffieux, C.; Boy, G. AER-11 STEM/EPICUR Test Report; Rapport no PLEIADE PSN-RES/SEREX/2014-00365 Institut de Radioprotection et de Surete Nucleaire: Caradache, France, 2014.

11. Colombani, J.; Planteur, S. AER-4, AER-5 and AER-6 STEM/EPICUR Test Report Formation and Release of Volatile Iodine Species from CsI Aerosols Under Radiation; Rapport no PLEIADE PSN-RES/SEREX/2013-00616 Institut de Radioprotection et de Surete Nucleaire: Caradache, France, 2013.

12. Leroy, O. ESTER-EPICUR LOCOV Test Report; Rapport no 2022-00313 Institut de Radioprotection et de Surete Nucleaire: Caradache, France, 2023.

13. Leroy, O.; Bosland, L. Study of the stability of iodine oxides (IxOy) aerosols in severe accident conditions. Ann. Nucl. Energy; 2023; 181, 109526. [DOI: https://dx.doi.org/10.1016/j.anucene.2022.109526]

14. Dickinson, S.; Auvinen, A.; Ammar, Y.; Bosland, L.; Clément, B.; Funke, F.; Glowa, G.; Kärkelä, T.; Powers, D.A.; Tietze, S. . Experimental and modelling studies of iodine oxide formation and aerosol behaviour relevant to nuclear reactor accidents. Ann. Nucl. Energy; 2014; 74, pp. 200-207. [DOI: https://dx.doi.org/10.1016/j.anucene.2014.05.012]

15. Funke, F.; Langrock, G.; Kanzleiter, T.; Poss, G.; Fischer, K.; Kühnel, A.; Weber, G.; Allelein, H.-J. Iodine oxide in large-scale THAI tests. Nucl. Eng. Des.; 2012; 245, pp. 206-222. [DOI: https://dx.doi.org/10.1016/j.nucengdes.2012.01.005]

16. Allelein, H.-J.; Auvinen, A.; Ball, J.; Gűntay, S.; Herranz, L.E.; Hidaka, A.; Jones, A.V.; Kissane, M.; Powers, D.; Weber, G. State-of-the-Art Report on Nuclear Aerosols; NEA/CSNI/R(2009)5 Organisation for Economic Co-operation and Development (OECD) Nuclear Energy Agency (NEA): Boulogne-Billancourt, France, 2009.

17. Powers, D.A.; Washington, K.E.; Burson, S.B.; Sprung, J.L. A Simplified Model of Aerosol Removal by Natural Processes in Reactor Containments; NUREG/CR-6189 U.S. NRC: Washington, DC, USA, 1996.

18. Malicki, M.; Lind, T. Parametric MELCOR 2.2 sensitivity and uncertainty study with a focus on aerosols based on Phebus test FPT1. Prog. Nucl. Energy; 2023; 158, 104609. [DOI: https://dx.doi.org/10.1016/j.pnucene.2023.104609]

19. Firnhaber, M.; Fischer, K.; Schwarz, S.; Weber, G. International Standard Problem ISP44 KAEVER Experiments on the Behavior of Core-melt Aerosols in a LWR Containment–Comparison Report; NEA/CSNI/R (2003) 5 Nuclear Energy Agency (NEA): Boulogne-Billancourt, France, 2003.

20. Lee, Y.; Cho, Y.J.; Lim, K. Coupling scheme of multicomponent sectional equations and Mason equations via transition rate matrix of hygroscopic growth applied to international standard problem No. 44. Ann. Nucl. Energy; 2019; 127, pp. 437-449. [DOI: https://dx.doi.org/10.1016/j.anucene.2018.12.028]

21. Gelbard, F.; Seinfeld, J.H. Simulation of multicomponent aerosol dynamics. J. Colloid Interface Sci.; 1980; 78, pp. 485-501. [DOI: https://dx.doi.org/10.1016/0021-9797(80)90587-1]

22. Pruppacher, H.R.; Klett, J.D. Microphysics of Clouds and Precipitation; D. Reidel Publishing Co.: Dordrecht, The Netherlands, 1980.

23. Gauntt, R.O.; Cole, R.K.; Erickson, C.M.; Gido, R.I.; Gasser, R.D.; Rodriguez, S.B.; Young, M.F. MELCOR Computer Code Manuals Vol. 2: Reference Manual, Version 1.8.6; NUREG/CR-6119 Sandia National Laboratories: Albuquerque, NM, USA, 2001.

24. Gelbard, F. MAEROS User Manual; NUREG/CR-1391, SAND80-0822 Sandia National Laboratories: Albuquerque, NM, USA, 1982.

25. Rovelli, G.; Miles, R.E.H.; Reid, J.P.; Clegg, S.L. Accurate Measurements of Aerosol Hygroscopic Growth over a Wide Range in Relative Humidity. J. Phys. Chem. A.; 2016; 120, pp. 4376-4388. [DOI: https://dx.doi.org/10.1021/acs.jpca.6b04194] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/27285052]

26. Ratner, B. The correlation coefficient: Its values range between +1/−1, or do they?. J. Target. Meas. Anal. Market.; 2009; 17, pp. 139-142. [DOI: https://dx.doi.org/10.1057/jt.2009.5]

27. Cohen, J.; Cohen, P.; West, S.G.; Aiken, L.S. Applied Multiple Regression/Correlation Analysis for the Behavioral Sciences; 3rd ed. Lawrence Erlbaum Associate Inc.: Mahwah, NJ, USA, 2003.

28. Iman, R.L.; Shortencarier, M.J. A FORTRAN 77 Program and User’s Guide for the Generation of Latin Hypercube and Random Samples for Use with Computer Models; NUREC/CR-3624, SAND83-2365 Sandia National Laboratories: Albuquerque, NM, USA, 1984.

29. Lee, Y.; Cho, Y.J. Preliminary importance analyses on model for pH in the presence of organic impurities in the aqueous phase for a severe accident of a nuclear power plant. Nucl. Eng. Technol.; 2024; 56, pp. 2079-2091. [DOI: https://dx.doi.org/10.1016/j.net.2024.01.016]

30. Hansson, H.C.; Ahlberg, M.S. Dynamic shape factors of sphere aggregates in an electric-field and their dependence on the Knudsen number. J. Aero. Sci.; 1985; 16, pp. 69-79. [DOI: https://dx.doi.org/10.1016/0021-8502(85)90021-7]

31. Hochrainer, D.; Hänel, G. Der dynamische formfaktor nicht-kugelf¨ormiger teilchen als funktion des luftdrucks. J. Aero. Sci.; 1975; 6, pp. 97-103. [DOI: https://dx.doi.org/10.1016/0021-8502(75)90002-6]

32. Wang, S.; Zhou, K.; Lu, X.; Chen, H.; Yang, F.; Li, Q.; Chen, J.; Prather, K.A.; Yang, X.; Wang, X. Online shape and density measurement of single aerosol particles. J. Aero. Sci.; 2022; 159, 105880. [DOI: https://dx.doi.org/10.1016/j.jaerosci.2021.105880]

33. Soffer, L.; Burson, S.B.; Ferrell, C.M.; Lee, R.Y.; Ridgely, J.N. Accident Source Terms for Light-Water Nuclear Power Plants; NUREG-1465 U.S. NRC: Washington, DC, USA, 1995.

34. ICRP. Annuals of the ICRP. ICRP Publication 119. Compendium of Dose Coefficients Based on ICRP Publication 60; Elsevier: Amsterdam, The Netherlands, 2012; Volume 41 (Suppl. 1).

35. Fukushima Daiichi Nuclear Power Station Accident Information Collection and Evaluation (FACE) Project. 2022; Available online: https://www.oecd-nea.org/jcms/pl_70741/fukushima-daiichi-nuclear-power-station-accident-information-collection-and-evaluation-face-project (accessed on 31 May 2025).