Uncertainty propagation in atmospheric dispersion models for radiological emergencies in the pre- and early release phase: summary of case studies

I. Korsakissok, R. Périllat, S. Andronopoulos, P. Bedwell, E. Berge, T. Charnock, G. Geertsema, F. Gering, T. Hamburger, H. Klein, S. Leadbetter, O.C. Lind, T. Pázmándi, Cs. Rudas, B. Salbu, A. Sogachev, N. Syed, J.M. Tomas, M. Ulimoen, H. de Vries and J. Wellings 1 IRSN – Institute for Radiation Protection and Nuclear Safety, Fontenay-aux-Roses, France. 2 PHIMECA engineering, Clermont-Ferrand, France. 3 EEAE/NCSRD Greek Atomic Energy Commission/National Centre for Scientific Research “Demokritos”, Agia Paraskevi, Greece. 4 PHE – Public Health England, Didcot, UK. 5 NMI MET – Norwegian Meteorological Institute, Oslo, Norway. 6 KNMI – Royal Netherlands Meteorological Institute, de Bilt, The Netherlands. 7 BfS Federal Office for Radiation Protection, Neuherberg, Germany. 8 Met Office, Exeter, UK. 9 NMBU/CERAD Norwegian University of Life Sciences, Centre for Environmental Radioactivity, Ås, Norway. 10 EK Centre for Energy Research, Budapest, Hungary. 11 DTU Wind Energy, Roskilde, Denmark. 12 DSA – Norwegian Radiation and Nuclear Safety Authority, Østerås, Norway. 13 RIVM – National Institute for Public Health and the Environment, Bilthoven, The Netherlands.


Introduction
In the framework of the European project CONFI-DENCE, Work Package 1 (WP1) focused on the uncertainties in the pre-and early phase of a radiological emergency, when environmental observations are not available and the evaluation of the environmental and health impact of the accident heavily relies on atmospheric dispersion modelling. The results of these simulations, including dose calculations, are used to infer recommendations for the protection of the population. However, the model outputs are subject to uncertainties (Gering et al., 2016;Sørensen et al., 2016Sørensen et al., , 2019, stemming from meteorological data (stochastic and modelling uncertainties) and lack of knowledge of the source term. The physical and numerical approximations made within the dispersion and dose models introduce additional uncertainties in the output. The first task of WP1 was to identify and characterize these input uncertainties . Then, several case studies were designed, for which different participants propagated the input uncertainties through their atmospheric dispersion models. For each case study, this resulted in several "ensembles" of results, consisting of a few tens to several hundreds of simulations. These ensembles were compared, to evaluate whether the inter-ensemble variability was significant or not by comparison to the intra-ensemble variability (given by the ensembles' spread). In the Fukushima case, the ensembles were also compared to radiological observations in the environment ( 137 Cs air concentration and gamma dose rates), in order to evaluate their "quality" by using statistical indicators. This also enables the identification of instances where measurements do not fall within the range of model uncertainty, which is indicative that the full model uncertainty has not been taken into account. The results of these case studies were also used to discuss various ways of presenting uncertainties, including different thresholds, maps and indicators. The ensembles outputs were also used as input by other WPs and for interactive discussion during panels and workshops. This paper describes the case studies, with a summary of the scenario, input data and participants (Sect. 2); then, some illustrative results are given for each case in Section 3; finally, a summary of the main findings and discussions is given in Section 4.

Description of the case studies
The summary of case studies scenarios (meteorology, source term) and participants is given in the Table 1.
Two hypothetical accident scenarios in Europe were designed: the "Western Norway" (WN) case, a release from a nuclear vessel west of the Norwegian coast, and the "Radiological Ensemble Modelling" (REM) case, a release from the Borssele Nuclear Power Plant in the Netherlands, but with a source term scaled for a 900 MWe reactor (instead of 485 MWe). Finally, three days of the Fukushima accident, from March 14th to March 16th, 2011, were also simulated. For each scenario, several participants carried out ensemble simulations using their respective atmospheric dispersion model (Tab. 1).

REM case study
For the REM case study, two release scenarios were considered: a 4-hour "short release" and a 72-hour "long release". These are described in Korsakissok et al. (2019aKorsakissok et al. ( , 2019c. For the short release, the overall emitted quantities of radionuclides are shown in Table 2. The starting time of release was given with an uncertainty of ±6 hours; the effective release height was 50 m ± 50 m and the released quantity varied between a third and a factor of 3 of the values given in Table 2. Each participant was free to choose how to take these uncertainties into account, depending on their computational capabilities. This resulted in ensembles with numbers of members ranging from 10 to 650 simulations. For the long release, an ensemble of 10 source terms, with complex kinetics, was extracted from a source term database built in the European project FASTNET using the severe accident code ASTEC (Chevalier-Jabet, 2019a, 2019b. The spread of the released quantities as a function of time for each radionuclide was designed to be representative of uncertainties stemming from the modelling of reactor physics, iodine chemistry, and other sources of uncertainty such as the unknown status of safety devices. Two meteorological scenarios were also considered: a case with a well-established wind direction on January 11th, 2017 (called REM1), representative of a small meteorological variability, and a "warm front passage" on January 12th, 2017 (called REM2), with high precipitation and turning winds resulting in a larger variability (Geertsema et al., 2019). The combination of the release and meteorological scenarios lead to three case studies: REM1 and REM2 with the short release (24 hours simulation each), and REM-L with the long release, covering 72 hours of simulation between January 11th and 13th, 2017, therefore, comprising both REM1 and REM2 meteorological situations. The meteorological data for these 72 hours was provided by KNMI, using the HARMONIE-AROME high-resolution meteorological model. A hybridlagged ensemble was constructed, combining different versions of the model and different forecast lead times, to provide 10 meteorological members representing the meteorological uncertainty.

Western Norway case study
The WN case, described in Berge et al. (2019), assumes a total release of 5.41 · 10 17 Bq over a 7-hour period from a hypothetical fire in a floating power plant 100 km off the west coast of Norway. The fire starts at 09 UTC on March 16th, 2017 and it is located at 52.5°N and 4°E (see Fig. 6). A lognormal particle size distribution is assumed for UO 2 , U 3 O 8 and RuO 2 (RuOxid) (see Tab. 3), while a single particle size is assumed for 137 Cs, 134 Cs, 144 Ce, fly ash, and 131 I. The particles have different densities. Five different emission scenarios are considered in which the fractional release in the first hour varies from 14.3% in scenario 1 to 90% of the total release in scenario 5. Constant release is assumed for the last 6 hours of the accident. During the first hour, the emissions are evenly distributed between 20 m and 500 m height above sea level; after this, the emissions are evenly distributed between 20 m and 100 m. The 10 meteorological ensembles are taken from the HARMONIE-AROME model operated with 2.5 km horizontal resolution (Müller et al., 2017). The HARMO-NIE-AROME model is run for a 66-hour period starting at 06 UTC on March 16th, 2017. The dispersion analysis focuses on the 24-hour period starting at the time of the release. The two dispersion models used are SNAP (Severe Nuclear Accident Program, Bartnicki et al., 2011) and DIPCOT (Dispersion over Complex Terrain, Andronopoulos et al., 2009) which is a part of more complex decision support system JRODOS used for the creating and running worldwide accident scenarios (Ievdin et al., 2010;Landman et al., 2016).

Fukushima case study
The Fukushima disaster occurred in March 2011, triggered by an earthquake followed by a tsunami. Radioactive materials S58 I. Korsakissok et al.: Radioprotection 2020, 55(HS1), S57-S68  were released into the atmosphere by the damaged reactors of the Fukushima Daiichi Nuclear Power Plant (FDNPP), mostly during the first three weeks following the earthquake. It resulted in significant contamination of the air and ground in Japan. According to the observations, two release periods were responsible for the highest contamination on Honshu Island, the first on March 14-16th and the second on March 20-22nd.
The first 3-day period was retained for study in CONFIDENCE WP1 (Korsakissok et al., 2019b). The ECMWF Integrated Forecast System was used to create a 72-hour ensemble forecast starting at 00:00 UTC on 14 March 2011 with the specification given in Table 1. Nine source terms from the literature were used to describe the release uncertainty. All source terms were inferred from the combination of facility events, radiological observations, meteorological and dispersion modelling. The total release considered over the three days varied by a factor of 7 (3.1-21.4 PBq) for 137 Cs and a factor of 9 (46-400 PBq) for 131 I between the different source terms.
In this study, two observation datasets are used. The first dataset consists of the hourly concentrations of 137 Cs retrieved by Tsuruta et al. (2014). The data were obtained from the automated air quality monitoring network and give information on the temporal variation of 137 Cs concentration close to ground level (108 stations). The second observation dataset is gamma dose rate measurements, provided by automated  Activity (Bq) 1.95 · 10 17 1.95 · 10 17 9.5 · 10 15 1 · 10 16 6 · 10 16 6.95 · 10 16 5.41 · 10 17 stations, with a 10-minute time step. There are 88 stations spread over Japan, although the spatial coverage is heterogeneous. The dose rate readings include the contribution from all radionuclides, including the short-lived species that could not be detected by other monitoring surveys due to the delay between contamination and measurement. They are composed of two parts: the direct plume contribution ("cloud-shine") and the gamma-rays emitted by radionuclides deposited on the ground ("ground-shine"). Cloud-shine is usually responsible for peak values observed during plume passage, whereas ground-shine constitutes a lasting contribution that continues after the plume has left the area and decreases due to radioactive decay and soil migration.

Illustrative results
The aim of this section is to give an overview of the various ways of presenting ensemble results and uncertainties, for various scenarios and different participants.

REM case studies
The REM case study with the short release is used here to illustrate the difference between the participants' ensembles. Figure 1 shows the median of the 137 Cs deposition 24 hours after the release, for the seven participants who contributed to the REM1 case study. Although the pattern is globally similar, there are clear differences that may be due to the different types of models, wet deposition schemes, diffusion schemes and interpolation methods. It should be emphasised that, in this particular case, no source term uncertainties were taken into account. Therefore, all participants performed 10 simulations (using the 10 meteorological ensembles) and only the dispersion and dose models differ between the participants. Figure 2 shows frequency maps of threshold exceedance for 137 Cs deposition, for a reference threshold of 10 kBq/m 2 . In each location, the frequency of threshold exceedance is given by the ratio between the number of simulations in the ensemble that give a deposition above the threshold and the total number of simulations: a probability of 100% means that all simulations (or "members") of the ensemble are above the threshold, 0% means that no member exceeds the value in the cell. The probability of threshold exceedance reflects the "level of agreement" between the different simulations for a given ensemble provided by a participant. If there is a small uncertainty, the frequency of threshold exceedance is very high within a particular area whilst it is zero outside that area, as can be seen in the REM1 case study without source perturbations (Fig. 2, left). If uncertainties related to the source term are included, there is a larger surface area of non-zero frequencies of threshold exceedance, with globally lower frequencies (Fig. 2, right). In the REM1 case, there is a small meteorological variability (Fig. 2, left). When the uncertainty  in the release time is included, some simulations feature a plume travelling North-East and therefore a non-zero probability of threshold exceedance appears in this area (Fig. 2, right). The last example for the REM case study is given for the REM2 "warm front passage" scenario. It consists of a comparison of the maximum distance of threshold exceedance (given by each ensemble member) for all participants. This is done by using box plots for three different variables and thresholds shown in Figure 3. The ensemble median is given in red, the blue boxes correspond to the 25th-75th percentiles, and the full ensemble spread is represented by the dashed lines (with outliers as blue crosses). Here, the inter-ensemble variability is illustrated by the range of variation of the medians. It is larger for 137 Cs deposition than for inhalation thyroid dose. The ensemble spread also varies more in the case of 137 Cs deposition. However, in the case of inhalation thyroid dose, some ensembles feature outliers with a maximum distance of more than 100 km (to be compared with the global median values ranging from 15 to 40 km). This raises the  question of how to use such ensembles in an operational context and how to take into account these "worst cases".

Western Norway case study
In the WN case, a low-pressure system in the Norwegian Sea gave rise to west and southwest winds toward the Norwegian west coast. The radioactive plume from the source was dispersed toward east and northeast and considerable wet deposition occurred when the plume hit the mountains of western Norway. Figures 4 and 5 present the 10th, 50th and 90th percentiles of accumulated concentration in air and accumulated deposition (dry plus wet) for the period 09 UTC 16.03.2017 to 09 UTC 17.03.2017. The results were obtained by combining 10 meteorological ensemble members and the 5 emission scenarios, giving in total 50 members in each of the two dispersion ensembles. Qualitatively, the two models give similar results, but the SNAP model yields the largest accumulated concentrations, while the DIPCOT model yields the largest total depositions. For example, the DIPCOT model shows a large area where the 50th percentile of total deposition is above 3 · 10 7 Bq/m 2 (Fig. 5), while the SNAP model gives values below 3 · 10 7 Bq/m 2 and even 1 · 10 7 Bq/m 2 in the same area. This could be due to differences in parameterization of wet and dry deposition, the treatment of the vertical and horizontal transport and the particle size distribution.
More detailed studies have been carried out at the source and for a focus area, Vikedal, situated in a fjord at the foot of a mountainous area approximately 45 km inland and about 150 km northeast of the source (see Fig. 6).
A summary of the spread in the dispersion calculations at the source and at Vikedal is presented in Figure 7 for accumulated air concentration and deposition. The ratio (90th percentile/10th percentile), representative of the ensemble's spread, is shown Figure 7a for the ten meteorological members combined with emission scenario 3 (50% release of the total emission the first hour). The same ratios for five emissions scenarios combined with ensemble member 0 are presented in Figure 7b. In addition, the 50th percentiles of the two dispersion models are compared in Figure 7c. The ratio (90th percentile/10th percentile) of the ten ensemble runs ranges from about 1 to 3 in the source area and from 2 to 4 at Vikedal. DIPCOT gives the largest spread at the source for both variables and at Vikedal for accumulated concentrations. The ratio (90th percentile/ 10th percentile) of the five emissions scenarios is of the same magnitude as the spread due to the meteorological ensembles. The ratio between the 50th percentile value of

S64
I. Korsakissok et al.: Radioprotection 2020, 55(HS1), S57-S68 the two models, representing the inter-model variability, is rather high (5.3 and 6.9) for total deposition at the source. The other ratios are ranging from 1 to 3, i.e. they are of similar magnitudes as the spread due to the meteorological ensembles and the emissions scenarios.

Fukushima case study
This section presents some comparisons of ensemble results to radiological observations, to illustrate how to assess the quality of an ensemble and thus, of the uncertainty propagation. This is done by plotting time series of the observed quantities (air concentrations or gamma dose rates) measured by the stations, along with the simulated time series obtained for each member within the ensemble (for a given participant). This provides a view of the ensemble's ability to encompass the observations. This is illustrated here for stations at Fukushima city, located 58 km north-west of FDNPP, in an area that was highly contaminated due to wet deposition during the episode of 14-16th March. Air concentrations are shown in Figure 8 and gamma dose rates in Figure 9, for four (out of six) participants.
The gamma dose rate observations shown in Figure 9 are typical of a wet deposition pattern: there is first an increase that corresponds to the beginning of the plume passage and/or rain, but no significant decrease after the plume has left the area, due to the dominating contribution of radionuclides deposited on the ground to the total measured gamma dose rate. The initial increase of the gamma dose rate corresponds to the beginning of the rain (between 03:00 and 07:00 UTC on March 15th) and is likely due to the scavenging of a plume in altitude, since concentrations at ground level are not high enough to explain the deposition values and subsequent gamma dose rate measured (as shown in Fig. 8). Overall, results for both air concentration and gamma dose rates show that all ensembles are able to encompass the observations, although some members highly overestimate the contamination. For instance, the maximum observed air concentration on the station is around 30 Bq/m 3 while maximum values for some members are of the order of a few hundreds, up to a thousand, Bq/m 3 (Fig. 8). When looking at the dark blue lines, representing the 25-75th percentile, globally the air concentration observations are encompassed, while the gamma dose rates are underestimated. The timing of the increase of the gamma dose rate is well encompassed by the models, thanks to the meteorological ensemble, which has a sufficient variability in the rain timing. The prediction of this timing is usually difficult for deterministic simulations, as light rains are not well forecasted (Mathieu et al., 2018a).

Summary and findings
This paper describes the different case studies of uncertainty propagation carried out within WP1 of CONFI-DENCE. Various scenarios were explored, including different meteorological scenarios (well established wind direction with small uncertainties and warm front for REM cases, storm for WN case, and stable situation with snow for the Fukushima accident) and different meteorological resolutions and models to describe the related uncertainties. The uncertainties related to the source term were also taken into account, especially those related to the release time (for REM case study) and kinetics (for Western Norway). The range of variation for the Fukushima source terms found in the literature showed that, even several years after the accident and with the help of numerous observations, source term uncertainties are still very high (up to a factor of seven on the released quantities for the 3day period).
The results of the REM hypothetical accident scenario highlight the importance of taking into account perturbations in source parameters and not only meteorological uncertainties. Indeed, perturbing the release time introduced a much larger variability, which triggered lower probabilities but over a larger overall contaminated area. It also showed that intermodel variability is not negligible, although secondary (in this case) to meteorological and source term uncertainties.
The WN case study showed that the spread due to meteorological uncertainties, emission and dispersion model formulation is of similar magnitude near the source and at an inland site located ca. 150 km northeast of the source. The maps of accumulated concentrations and depositions are qualitatively similar for the two models. However, it is clearly seen, that the SNAP model yields higher concentrations, but lower depositions than the DIPCOT model. A more in-depth analysis, such as analyses of the 3-D transport, wet-and dry deposition and particle size distributions would be needed to better understand the different behavior of the two dispersion models.
Finally, the comparisons to environmental observations over Japan in the case of the Fukushima accident showed that the ensembles encompass the observations reasonably well. However, there are large variations between the participants when looking at medians or 25-75th percentiles. In addition, some ensemble members hugely overestimate the observations (by several decades), which raises the issue of the use of such ensembles for decision making.
A range of different representations of the uncertainties has been used to analyze the case studies; the need for various, complementary outputs, including graphical representations and statistical indicators , has been highlighted.