Ensemble Dispersion Simulation of a Point-Source Radioactive Aerosol Using Perturbed Meteorological Fields over Eastern Japan

: We conducted single-model initial-perturbed ensemble simulations to quantify uncertainty in aerosol dispersion modeling, focusing on a point-source radioactive aerosol emitted from the Fukushima Daiichi Nuclear Power Plant (FDNPP) in March 2011. The ensembles of the meteorological variables were prepared using a data assimilation system that consisted of a non-hydrostatic weather-forecast model with a 3-km horizontal resolution and a four-dimensional local ensemble transform Kalman ﬁlter (4D-LETKF) with 20 ensemble members. The emission of radioactive aerosol was not perturbed. The weather and aerosol simulations were validated with in-situ measurements at Hitachi and Tokai, respectively, approximately 100 km south of the FDNPP. The ensemble simulations provided probabilistic information and multiple case scenarios for the radioactive aerosol plumes. Some of the ensemble members successfully reproduced the arrival time and intensity of the radioactive aerosol plumes, even when the deterministic simulation failed to reproduce them. We found that a small ensemble spread of wind speed produced large uncertainties in aerosol concentrations. M.K. (Mizuo Kajino), and M.K. (Masaru Kunii); Software, M.K. (Mizuo Kajino) and M.K. (Masaru Kunii); Validation, T.T.S.; Formal Analysis, T.T.S.; Investigation, T.T.S.; Data Curation, T.T.S. and M.K. (Mizuo Kajino); Writ-ing—Original Draft Preparation, T.T.S.; Writing—Review & Editing, T.T.S.; Visualization, T.T.S.; Funding Acquisition,


Introduction
Ensemble simulation is a set of multiple numerical simulations that have slightly different initial conditions, boundary conditions, parameters, or models that are all geophysically plausible. Such a simulation enables the estimation of the predictability or reliability of the model simulation by providing a spread of ensemble forecasts. The simulation is most certain if the ensemble members are close to each other; otherwise, the ensemble provides a possible range of different events. Thus, probabilistic model information can be obtained from an ensemble simulation. Additionally, an ensemble simulation provides an ensemble average that is often more accurate than a deterministic single simulation because the model errors tend to be averaged out.
From a scientific viewpoint, the model predictability indicates not only the imperfection of simulation models but also the Lorenz's deterministic chaos of Earth systems. The error growth and propagation in the model simulation depend on the chaotic advection, diffusion, precipitation, thermodynamics, and chemistry, which all should be explored in detail. From a practical perspective, probabilistic model information complements deterministic model information, especially for atmospheric forecasts. Therefore, ensemble prediction systems (EPSs) have been developed worldwide by operational weather forecast centers. These systems have adopted initial-condition ensemble simulations that are suitable for error growth evaluation.
However, it is difficult to generate the ensemble perturbations of initial conditions because randomly chosen (Monte Carlo) perturbations are likely to fade away or fail to

Methodology
The EnKF is an approximate treatment of the Kalman filter for application to highdimensional systems such as the atmosphere cf. [17,18]. The Kalman filter defines an analysis as an arithmetic weighted mean of forecasts and observations, imposing a minimum variance estimation of the analysis error on the weight optimization [19]. We simultaneously obtain data assimilation products (i.e., analysis) and perturbations (i.e., ensemble members) by repeating the EnKF procedure. The analysis is the mean of the ensemble members. The perturbations generated by the EnKF are qualitatively superior to a random perturbation because they reflect the model uncertainty distribution and are flow-dependent, similar to the SV method. Furthermore, Wang et al. [20] reported that when the ensemble size is small, the EnKF method has a statistical advantage because the other methods consistently generate symmetric positive-negative paired ensemble mem-bers to keep the average value and thus cannot make statistical full use of the ensemble dimension. In contrast, the EnKF does not generate symmetric pairs but keeps the overall average. While the EnKF has several numerical implementation methods, the square root EnKF is implemented in this study. In the square root EnKF, each ensemble member retains its identity through the data assimilation cycle because its relative position in the state space among the ensemble members is invariant.
In this study, the model domain covered eastern Japan as shown in Figure 1 and its horizontal resolution was set to 3 km, which represented a typical grid scale for the regional simulations implemented for the FDNPP accident cf. [37,39]. The model settings, such as map projection, vertical coordinate, turbulence scheme, convective scheme, and terrain features, were the same as those of the 3-km grid simulation performed by Sekiyama et al. [37,39]. The domain consists of 215 × 259 horizontal grid points in the Lambert conformal projection and 60 vertical levels including 11 levels below 1 km above ground level. The terrain features were generated from the global digital elevation data with a horizontal grid spacing of 30 arc seconds (GTOPO30) provided from the U.S. Geological Survey. The turbulence scheme was based on the improved Mellor-Yamada level 3 closure model [42,43]. A cumulus parameterization was not used in this study. The data assimilation system was initiated at 06:00 UTC on 10 March 2011 with 20 ensemble members and a 3-h time window. The assimilation settings, such as time slots (3 h), prognostic variables (three wind components, temperature, pressure, water vapor mixing ratio, and water/ice microphysics variables), inflation scheme (adaptively multiplicative factors at each grid point), and covariance localization (1/e 0.5 within 150 km in the horizontal and 0.2 natural-logarithmical p-coordinate in the vertical), were the same as those of the 3-km grid simulation performed by Sekiyama et al. [37,39]. We obtained the initial condition and the boundary conditions from the JMA operational global 15-km grid analysis.
We assimilated JMA's operational observation dataset, which was integrated and quality-controlled for the JMA mesoscale Non-hydrostatic-model four-dimensional Variational data Assimilation system (JNoVA) [44], similarly to Kunii [21] and Sekiyama et al. [37]. Additionally, we assimilated surface wind observations acquired by the Automated Meteorological Data Acquisition System (AMeDAS) similarly to Sekiyama et al. [39]. AMeDAS is a nationwide meteorological observation network managed by JMA. The data assimilation system generated 20 ensemble members every 3 h (hereafter called ensemble analysis members) and simultaneously calculated the mean value of the ensemble members as a deterministic analysis. The forecasts were calculated by the identical JMA-NHM using these 20 ensemble members and a deterministic analysis as the initial conditions. Hereafter, these forecasts are called ensemble forecast members and a deterministic forecast, respectively.
We used the emission scenario of the radioactive Cs-137, which was released from the FDNPP, estimated by the Japan Atomic Energy Agency (JAEA) [55][56][57]. Cs-containing sulfate-organics-mixed aerosol particles [58] were injected at every time step into a grid cell above the FDNPP. The emission scenario has been revised by JAEA several times after 2012, e.g., [59]. However, the difference between the previous ones and the revised ones is not very large in comparison with the dispersion model uncertainty (i.e., ensemble spread) of Cs-137 concentrations in this study. Since we focus on the quantification of the dispersion model uncertainty, the revision of the Cs-137 emission scenario scarcely affects the results of the uncertainty evaluation in this study.
The data assimilation cycle and the dispersion simulations were performed continuously from 11 March to 1 April 2011. The "analysis" run contained 20 ensemble simulations and one deterministic simulation. Although the analysis is only provided every 3 h, the dispersion simulations require the meteorological variables at much smaller time intervals. Therefore, the meteorological variables were generated between the 3-h analysis points by the 3-h forecast runs using the identical JMA-NHM. The variables were stored at every 10 min of simulation time and inputted into the RAQM2, and the variables were linearly interpolated during each 10-min interval. Additionally, we performed a "forecast" run for 24 h in two specific periods (15 March and 21 March, see Section 3). Each forecast was started at 21:00 local time (JST) of the previous day using the "analysis" run as the initial condition. The "forecast" run also contained 20 ensemble simulations and one deterministic simulation.
Atmosphere 2021, 12, x FOR PEER REVIEW 4 of 14 are described in work of Kajino et al. [46] and Sekiyama et al. [37,39]. Note that the RAQM2 used in this study implements simplified aerosol dynamics compared with those of Kajino et al. [46] by assuming perpetual particle size distribution similarly to Sekiyama et al. [37,39]. The combination of the JMA-NHM, the LETKF, the JNoVA+AMeDAS observations, and the RAQM2 has been successfully used for the Fukushima radioactive pollution simulation [37,[39][40][41][51][52][53][54]. We used the emission scenario of the radioactive Cs-137, which was released from the FDNPP, estimated by the Japan Atomic Energy Agency (JAEA) [55][56][57]. Cs-containing sulfate-organics-mixed aerosol particles [58] were injected at every time step into a grid cell above the FDNPP. The emission scenario has been revised by JAEA several times after 2012, e.g., [59]. However, the difference between the previous ones and the revised ones is not very large in comparison with the dispersion model uncertainty (i.e., ensemble spread) of Cs-137 concentrations in this study. Since we focus on the quantification of the dispersion model uncertainty, the revision of the Cs-137 emission scenario scarcely affects the results of the uncertainty evaluation in this study.
The data assimilation cycle and the dispersion simulations were performed continuously from 11 March to 1 April 2011. The "analysis" run contained 20 ensemble simulations and one deterministic simulation. Although the analysis is only provided every 3 h, the dispersion simulations require the meteorological variables at much smaller time intervals. Therefore, the meteorological variables were generated between the 3-h analysis points by the 3-h forecast runs using the identical JMA-NHM. The variables were stored at every 10 min of simulation time and inputted into the RAQM2, and the variables were linearly interpolated during each 10-min interval. Additionally, we performed a "forecast" run for 24 h in two specific periods (15 March and 21 March, see Section 3). Each forecast was started at 21:00 local time (JST) of the previous day using the "analysis" run as the initial condition. The "forecast" run also contained 20 ensemble simulations and one deterministic simulation.

Results and Discussion
Here, we focused on March 15 (Period 1) and 21 (Period 2), 2011 local time to investigate the radioactive plumes that were carried landward. Nakajima et al. [60] pointed out that the plume intrusion inland occurred twice on a large scale; i.e., 15 March and 20-21

Results and Discussion
Here, we focused on March 15 (Period 1) and 21 (Period 2), 2011 local time to investigate the radioactive plumes that were carried landward. Nakajima et al. [60] pointed out that the plume intrusion inland occurred twice on a large scale; i.e., 15 March and 20-21 March. The comparison of the radioactive aerosol concentration was performed at the model grid corresponding to the location of Tokai, where JAEA has been operationally monitoring radionuclide concentrations and clearly detected highly radioactive plumes on both 15 and 21 March, 2011 [61]. The JAEA Tokai facilities are located approximately 100 km south of the FDNPP (Figure 1). The meteorological components were compared The wind speed (u and v) ensemble had a small spread and was almost synchronized with the JMA operational analysis winds. The wind errors of the analysis members (i.e., the distance from filled circles to a bunch of black lines) were larger than the difference between the wind ensemble members (i.e., the spread width of a bunch of black lines). Figure 2b shows the forecast (initiated at 21:00 14 March 2011 JST) in which the wind ensemble had a large spread. However, the wind forecast spread was still smaller than the wind forecast error (i.e., the averaged distance from filled circles to a bunch of black lines).
March. The comparison of the radioactive aerosol concentration was performed at the model grid corresponding to the location of Tokai, where JAEA has been operationally monitoring radionuclide concentrations and clearly detected highly radioactive plumes on both 15 and 21 March, 2011 [61]. The JAEA Tokai facilities are located approximately 100 km south of the FDNPP (Figure 1). The meteorological components were compared at the model grid corresponding to the location of Hitachi (Figure 1), where the nearest AMeDAS station (10 km north) to Tokai was located. AMeDAS observations, and the JMA operational 5-km gridded analysis. The wind speed (u and v) ensemble had a small spread and was almost synchronized with the JMA operational analysis winds. The wind errors of the analysis members (i.e., the distance from filled circles to a bunch of black lines) were larger than the difference between the wind ensemble members (i.e., the spread width of a bunch of black lines). Figure 2b shows the forecast (initiated at 21:00 14 March 2011 JST) in which the wind ensemble had a large spread. However, the wind forecast spread was still smaller than the wind forecast error (i.e., the averaged distance from filled circles to a bunch of black lines). In contrast to the wind speed ensemble, the Cs-137 ensemble concentrations had a large spread. Some of the ensemble members successfully represented the real peak concentration but failed to represent the peak timing, appearing two hours early. The deterministic analysis underestimated the real concentration. In Figure 3b, the Cs-137 forecast members presented a larger ensemble spread for the concentration. Some of the Cs-137 forecast members presented large overestimations before and after the peak. In contrast to the wind speed ensemble, the Cs-137 ensemble concentrations had a large spread. Some of the ensemble members successfully represented the real peak concentration but failed to represent the peak timing, appearing two hours early. The deterministic analysis underestimated the real concentration. In Figure 3b, the Cs-137 forecast members presented a larger ensemble spread for the concentration. Some of the Cs-137 forecast members presented large overestimations before and after the peak.

Period 1 (15 March 2011)
In this period, the radioactive aerosol plume projected from the FDNPP to the south, coastwise, and then swept across Tokai (Figure 4). The percentile distribution of the 20 ensemble members was narrow in the analysis (Figure 4a) but relatively broad in the forecast (Figure 4b), which was in agreement with the time series of the ensemble analysis/forecast concentrations at Tokai that are shown in Figure 3. The threshold (15 Bq/m 3 ) used here was defined as the half value of the air quality standard of Japan's radioisotope regulations. The percentile distribution of aerosol concentrations can be usefully applied to probabilistic forecasts such as the chance-of-rain forecast. The probabilistic forecasts provide multiple scenarios for environmental pollution or disasters. Generally, the accuracy of the forecast decreases with time, and consequently the percentile distribution tends to diffuse with time. Figure 4 shows that the forecast percentile distributions were very similar to the analysis percentile distributions in the short forecasts (02:00 15 March and 06:00 15 March), indicating the high accuracy of the forecast. However, the distributions were less similar in the longer forecast (10:00 15 March), in which the forecast percentile distribution was diffused.
The ensemble analysis of the surface wind speed (u and v) exhibited a relatively small spread, which was usually less than 1 m/s even though the analysis errors (i.e., analysis minus observation) were generally 1 or 2 m/s. For example, the relative standard deviation (RSD) of the ensemble analysis for the specific 4 h during Period 1 (02:00-08:00 15 March) was 5% on average (Table 1). By contrast, the ensemble of the surface Cs-137 concentration had a large spread, in which some members occasionally presented almost zero concentrations, whereas others presented very high concentrations. The RSD of the Cs-137 concentration analysis for the same time periods was 93 % on average (Table 1). This result indicates that a small ensemble spread in meteorology produces a large ensemble spread for aerosol concentration. Thus, the uncertainty on the concentration is amplified in comparison with that on the wind field.
The ensemble forecast exhibited the same behavior as the analysis. The ensemble spread was relatively small in the meteorological simulation and very large in the dispersion simulation. The forecast RSD of the wind speed for the same time periods mentioned above (02:00-08:00 15 March) was 7% on average (Table 1). In contrast, the forecast RSD of the Cs-137 concentration was 82% (Table 1). A comparison of the analysis and forecast RSDs indicates that the errors of the dispersion models are not linearly correlated with the errors of the meteorological models. In this period, the radioactive aerosol plume projected from the FDNPP to the south, coastwise, and then swept across Tokai (Figure 4). The percentile distribution of the 20 ensemble members was narrow in the analysis (Figure 4a) but relatively broad in the forecast (Figure 4b), which was in agreement with the time series of the ensemble analysis/forecast concentrations at Tokai that are shown in Figure 3. The threshold (15 Bq/m3) used here was defined as the half value of the air quality standard of Japan's radioisotope regulations. The percentile distribution of aerosol concentrations can be usefully applied to probabilistic forecasts such as the chance-of-rain forecast. The probabilistic forecasts provide multiple scenarios for environmental pollution or disasters. Generally, the accuracy of the forecast decreases with time, and consequently the percentile distribution tends to diffuse with time. Figure 4 shows that the forecast percentile distributions were very similar to the analysis percentile distributions in the short forecasts (02:00 15 March and 06:00 15 March), indicating the high accuracy of the forecast. However, the distributions were less similar in the longer forecast (10:00 15 March), in which the forecast percentile distribution was diffused.   The ensemble analysis of the surface wind speed (u and v) exhibited a relatively small spread, which was usually less than 1 m/s even though the analysis errors (i.e., analysis minus observation) were generally 1 or 2 m/s. For example, the relative standard deviation (RSD) of the ensemble analysis for the specific 4 h during Period 1 (02:00-08:00 15 March) was 5% on average (Table 1). By contrast, the ensemble of the surface Cs-137 concentration had a large spread, in which some members occasionally presented almost zero concentrations, whereas others presented very high concentrations. The RSD of the Cs-137 concentration analysis for the same time periods was 93 % on average (Table 1). This result indicates that a small ensemble spread in meteorology produces a large ensemble spread for aerosol concentration. Thus, the uncertainty on the concentration is amplified in comparison with that on the wind field.
The ensemble forecast exhibited the same behavior as the analysis. The ensemble spread was relatively small in the meteorological simulation and very large in the dispersion simulation. The forecast RSD of the wind speed for the same time periods mentioned above (02:00-08:00 15 March) was 7% on average (Table 1). In contrast, the forecast RSD of the Cs-137 concentration was 82% (Table 1). A comparison of the analysis and forecast RSDs indicates that the errors of the dispersion models are not linearly correlated with the errors of the meteorological models.

Period 2 (21 March 2011)
In comparison with Period 1, the ensemble spreads tended to be larger in Period 2. Figure 5 shows time series of the meteorological ensembles at the Hitachi AMeDAS station similar to those in Figure 2 but for the time period from 21:00 20 March to 21:00 21 March, 2011 JST. Note that it was raining or snowing on this day in a wide area of eastern Japan. As seen in Figure 5b, the ensemble forecast (initiated at 21:00 20 March 2011 JST) presented chaotic motions of the wind speed (u and v). Furthermore, the ensemble spread of the precipitation forecast was extremely large.
The time series of the Cs-137 concentrations at the Tokai JAEA station are shown in Figure 6 similar to those in Figure 3, except that they are in the time period from 21:00 20 March to 21:00 21 March 2011 JST. The ensemble spreads were very large, similar to those in Period 1. As shown in Figure 6a, the deterministic analysis failed to represent the peak timing with a two-hour delay. However, some of the ensemble analysis members successfully represented the real peak timing and concentration. In contrast to the deterministic analysis, the deterministic forecast completely failed to represent the plume arrival  Figure 6b). Many of the ensemble forecast members behaved similarly to the deterministic forecast. This failure was caused by the wet deposition (below-cloud scavenging) that occurred before the plume arrived at Tokai in the forecast simulations.     In Period 2, the RSD of the wind speed analysis for the specified 4 h (04:00-10:00 21 March) was 10% on average (Table 1). However, the RSD of the Cs-137 concentration analysis for the same time periods was 77% on average (Table 1). These scores also indicate that the uncertainty on the concentration is amplified as shown in Period 1. This finding implies that if we deterministically pursue an accurate dispersion, it is necessary to unrealistically increase the accuracy of the meteorological simulation. As expected, the percentile distributions of the Cs-137 concentration (Figure 7) tended to be broader than those in Period 1 in both the analysis and forecast before raining (approximately 08:00 JST). The percentile distributions rapidly shrank after the rain because the radioactive aerosols were deposited through precipitation. In the forecast, the precipitation area or timing was slightly shifted from the actual values, and consequently, the distribution of wet deposition was inadequately distorted. The forecast RSD of the wind speed for the same time periods was 23% on average and that of the Cs-137 concentration was 235% ( Table 1). The extremely large forecast errors of the Cs-137 concentration were caused by the forecast errors of the precipitation (not only the strength but also the timing). Specifically, the deterministic forecast almost completely failed to reproduce the Cs-137 plume arriving at Tokai (Figure 6b). In comparison with the extreme error of the Cs-137 concentration, the precipitation error was moderate as shown in the right panel of Figure 5b.
On this day, light precipitation was widespread over Japan (Figure 8a; derived from the JMA Radar/rain-gauge Analyzed Precipitation data [62]) and reasonably reproduced by the deterministic analysis (Figure 8b). The deterministic forecast (Figure 8c) produced a different distribution from the observed values ( Figure 8a) and the analyzed values (Figure 8b). This difference caused the Cs-137 concentration to have a high error, indicating that the error of the dispersion models is not only amplified in comparison with the error of the wind speed but also crucially magnified by the error of the precipitation because of the high sensitivity of aerosol deposition on precipitation. However, even so, some of the ensemble members were successful in reproducing the high concentration at Tokai. This is the advantage of ensemble simulations. The forecast RSD of the wind speed for the same time periods was 23% on average and that of the Cs-137 concentration was 235% ( Table 1). The extremely large forecast errors of the Cs-137 concentration were caused by the forecast errors of the precipitation (not only the strength but also the timing). Specifically, the deterministic forecast almost completely failed to reproduce the Cs-137 plume arriving at Tokai (Figure 6b). In comparison with the extreme error of the Cs-137 concentration, the precipitation error was moderate as shown in the right panel of Figure 5b.
On this day, light precipitation was widespread over Japan (Figure 8a; derived from the JMA Radar/rain-gauge Analyzed Precipitation data [62]) and reasonably reproduced by the deterministic analysis (Figure 8b). The deterministic forecast (Figure 8c) produced a different distribution from the observed values ( Figure 8a) and the analyzed values (Figure 8b). This difference caused the Cs-137 concentration to have a high error, indicating that the error of the dispersion models is not only amplified in comparison with the error of the wind speed but also crucially magnified by the error of the precipitation because of the high sensitivity of aerosol deposition on precipitation. However, even so, some of the ensemble members were successful in reproducing the high concentration at Tokai. This is the advantage of ensemble simulations.

Conclusions
We conducted ensemble simulations for the dispersion of a point-source aerosol using perturbed meteorological fields. The ensemble simulations provided probabilistic information and multiple case scenarios for the aerosol dispersion. We found that a small ensemble spread of wind speed resulted in a large uncertainty in aerosol concentrations, i.e., the uncertainty on the aerosol dispersion was amplified in comparison with that on the wind simulation. This finding implies that a high accuracy of dispersion modeling requires much higher accuracy of meteorological modeling, thus representing a limitation of deterministic dispersion simulations for analyzing/predicting the location and intensity of aerosol plumes. Therefore, the probabilistic information of ensemble simulations exhibits great potential for aerosol analysis and prediction. The deterministic simulation did not provide the best analysis/prediction in this study. However, some of the ensemble members successfully reproduced the arrival time and intensity of the aerosol plumes. With only a deterministic simulation, it is not possible to account for another event. Regrettably, in the field of atmospheric chemistry modeling, too much emphasis has been placed on deterministic simulations uncritically. The usefulness of ensemble simulations should be recognized to a greater extent.
The errors in the aerosol simulation were not only cumulative with the errors in wind speed simulation but also crucially magnified by the errors in the precipitation simulation because of the dependence of aerosol deposition on precipitation. Although the singlemodel initial-perturbed ensemble simulation as used in this study is a powerful tool to explore probabilistic analysis/prediction, the limitations of the single-model simulation should also be considered because the single-model simulation implements only a single module for the precipitation.

Conclusions
We conducted ensemble simulations for the dispersion of a point-source aerosol using perturbed meteorological fields. The ensemble simulations provided probabilistic information and multiple case scenarios for the aerosol dispersion. We found that a small ensemble spread of wind speed resulted in a large uncertainty in aerosol concentrations, i.e., the uncertainty on the aerosol dispersion was amplified in comparison with that on the wind simulation. This finding implies that a high accuracy of dispersion modeling requires much higher accuracy of meteorological modeling, thus representing a limitation of deterministic dispersion simulations for analyzing/predicting the location and intensity of aerosol plumes. Therefore, the probabilistic information of ensemble simulations exhibits great potential for aerosol analysis and prediction. The deterministic simulation did not provide the best analysis/prediction in this study. However, some of the ensemble members successfully reproduced the arrival time and intensity of the aerosol plumes. With only a deterministic simulation, it is not possible to account for another event. Regrettably, in the field of atmospheric chemistry modeling, too much emphasis has been placed on deterministic simulations uncritically. The usefulness of ensemble simulations should be recognized to a greater extent.
The errors in the aerosol simulation were not only cumulative with the errors in wind speed simulation but also crucially magnified by the errors in the precipitation simulation because of the dependence of aerosol deposition on precipitation. Although the singlemodel initial-perturbed ensemble simulation as used in this study is a powerful tool to explore probabilistic analysis/prediction, the limitations of the single-model simulation should also be considered because the single-model simulation implements only a single module for the precipitation.