Preprint

Article

Airspace Contamination by Volcanic Ash from Sequences of Etna Paroxysms: Coupling the WRF‐Chem Dispersion Model with Near‐Source L‐Band Radar Observations

Altmetrics

Downloads

111

Views

Comments

A peer-reviewed article of this preprint also exists.

This version is not peer-reviewed

Submitted:

27 June 2023

Posted:

29 June 2023

You are already at the latest version

Alerts

Abstract

Volcanic emissions (ash, gas, aerosols) dispersed in the atmosphere during explosive eruptions generate hazards affecting aviation, human health, air quality and the environment. We document for the first time the contamination of airspace by very fine volcanic ash due to sequences of tran-sient ash plumes from Mount Etna. The atmospheric dispersal of sub-10 μm (PM10) ash is mod-elled using the WRF-Chem model coupled online with meteorology and aerosols and offline with Mass Eruption Rates (MER) derived from near-vent Doppler radar measurements and inferred plume altitudes. We analyse two sequences of paroxysms with widely varied volcanological conditions and contrasted meteorological synoptic patterns in October–December 2013 and on 3-5 December 2015. We analyse the PM10 ash dispersal simulation maps in terms of time-averaged columnar ash density, concentration at specified flight levels averaged over the entire sequence interval, and daily average concentration during selected paroxysm days at these flight levels. The very fine ash from such eruption sequences is shown to easily contaminate the airspace around the volcano within a radius of about 1000 km in a matter of a few days. Synoptic patterns with relatively weak tropospheric currents lead to the accumulation of PM10 ash at a re-gional scale all around Etna. In this context, closely interspersed paroxysms tend to accumulate very fine ash more diffusively in the lower troposphere and in stretched ash clouds higher up in the troposphere. Low-pressure, high-winds weather systems tend to stretch ash clouds into ~100 km-wide clouds forming large-scale vortices 800-1600 km in diameter. Daily average PM10 ash concentrations commonly exceed the aviation hazard threshold up to 1000 km downwind from the volcano and up to the upper troposphere for intense paroxysms. Vertical distributions show ash cloud thicknesses in the range 0.7–3 km, and PM10 sometimes stagnating at ground level, represent a potential health hazard.

Keywords:

Subject: Environmental and Earth Sciences - Atmospheric Science and Meteorology

1. Introduction

1.1. Volcanic ash hazards to aviation

During explosive eruptions, large quantities of volcanic material (gas, ash, aerosols) are injected into the atmosphere as an ash plume, reaching tropospheric and sometimes even stratospheric heights. The subsequent ash cloud drifts according to the dominant winds in altitude, generating hazards for aviation, human health, and the environment (pollution) both in the short and long terms. Whereas lapilli and coarse ash (> 1 mm in diameter) fall out within an hour, fine ash (< 1 mm, atmospheric residence > 30 min) and especially very fine ash (< 30 μm, atmospheric residence from three hours to several days) are agents for volcanic cloud hazards to aircraft and aviation services [1]. Impacts can be engine power loss, blockage of sensors, impaired vision owing to windscreen damage, contamination of aircraft ventilation and pressurisation systems, health risks to crew and passengers, abrasive damage of aircraft exterior, malfunction of electric/electronic elements with potential short-circuits, manoeuvre for cloud avoidance, flight disruptions, deposits of volcanic ash on runways, limited ground operations, airspace and airport closure, and all subsequent implications for travelling passengers [2]. Because very fine ash is often dispersed thousands of kilometres from the source volcano, it may remain for weeks in the atmosphere in amounts hardly detectable by aircraft, Air Traffic Control and weather radars, or even visually. Ash products from satellite imagery used for real-time monitoring (e.g., Jiménez-Escalona et al. [3]) sometimes also fail to detect ash due to (i) overcast weather, (ii) mask effects from very fine ash acting as cloud condensation and ice nuclei resulting in a mixed plume/cloud containing liquid water, SO₂ and ice, as reported by Marchese et al. [4] for the 3 December 2015 paroxysm of Etna, and (iii) dilute particle concentrations below detection threshold in distal areas. Even when ash is detected, communication delays to aircraft from the hazard occurrence time remain another limitation. In our study, we model the dispersion of very fine ash with particles less than 10 μm in size (PM10), as these are also well-known to affect air quality and human health (sub-10 μm inhalable volcanic ash, including the subset of respirable volcanic ash below 4 μm) (e.g., Eychenne et al. [5]). Very fine volcanic ash, SO2 and aerosols further impact the environment and climate through radiative forcing (e.g., Sellito et al. [6]).

The modelling of volcanic plume transport has attracted increasing attention to the natural hazard research community because of the impact of volcanic emissions dispersal on aviation and the economy. The global economic impact of the eruption of Eyjafjallajökull in 2010 has been estimated by Oxford Economics to be ~US$ 5 billion, with about US$ 2 billions of direct financial loss to the aviation industry. According to EUROCONTROL [7], it could be more than halved owing to improvements achieved in the air sector regarding crisis management and optimised operational performances in about a decade. For example, cost-based trajectory optimization based on dispersion models and flight planning software decreases cost and disruption of air traffic deriving from volcanic events [8].

The number of reported aircraft encounters with volcanic ash clouds causing various degrees of airframe or engine damage amounted to 129 between 1953 and 2009 [9] and to 113 between 2010 and 2014 [10]; yet, not all encounters are reported and therefore their number remains underestimated. While the most damaging encounters have occurred within 24 h of eruption onset and/or within 1000 km of the source, less safety-significant but still economically damaging encounters have occurred at greater distances and extended times [9]. In parallel, the total number of passengers transported each year also increases, including by a factor of 2.8 in the two decades preceding 2019 when it reached over 4.5 billion [2].

As demonstrated in this paper, Mount Etna emissions present significant hazards to regional air traffic, in addition to significant economic losses. The closest international aeronautic infrastructure, the Catania airport, must frequently limit its services because of Etna ash emissions, leading to numerous rerouted or cancelled flights and airport closures. About 11 million annual commercial passengers have transited through Catania airport in 2022 (6.4 and 7.1 million in 2013 and 2015, respectively), for about 76000 flights in 2022 involving 86 routes to domestic and European destinations (https://cataniaairport.com/statistics/, accessed on 3 May 2023).

1.2. Volcanic Ash Advisories for Etna

With the purpose of keeping aviators informed about volcanic hazards, the International Civil Aviation Organization (ICAO) established nine Volcanic Ash Advisory Centres (VAAC) covering the world as part of the International Airways Volcano Watch (IAVW). Each of the nine centres under the IAVW responds to reports of volcanic ash within their region and provides the aviation community with information on the ash cloud's current location and forecasts of future extent in case of significant volcanic eruptions. Ash warnings issued by the VAACs are in the form of standardised Volcanic Ash Advisory messages (VAA, sometimes including Graphics: VAG) based on Volcano Observatory Notices for Aviation (VONA), aircraft in flight (AIREP) and governmental agencies utilising meteorological ground stations and satellite imagery (SIGMET, NOTAM, ASHTAM) (https://www.ssd.noaa.gov/VAAC/vaac.html).

The Toulouse VAAC of Météo-France (https://vaac.meteo.fr/), with London VAAC as a backup, is in charge of the nowcasting of volcanic emissions dispersal in a very large zone of the globe including the Mediterranean basin hosting very active volcanoes like Etna and Stromboli. Between 2006 and 2022, the Toulouse VAAC has released 792 Volcanic Ash Advisories (VAA, most including Graphics - VAG) for Etna, including 22 and 29 during/following the eruptions from the two sequences analysed in this paper (Tables 1, 2, A1). The reaction time of the VAACs depends mainly on the on-site volcano observatory that must timely provide a VONA with critical information about the on-going activity and input parameters. At Etna, the all-year 24/7 monitoring is mainly carried out by INGV-Osservatorio Etneo, transmitting VONA, but also by information to the Italian Department of Civil Protection (https://www.protezionecivile.gov.it/en/) through specific communication procedures during crises.

As for volcanic hazards to aviation, present challenges include the detection of volcanic ash clouds, the accurate forecasting of their dispersion, and the timely and targeted communication of this information. Dispersion model simulations may contribute to evaluate the space-time distribution of volcanic constituents (solid particles, gases, aerosols), including for ash dispersal nowcasting as implemented by the VAACs during eruptions. Understanding engine and airframe tolerances to ash ingestion and gas effects would further inform the operational risk management of airlines.

1.3. Improving ash dispersal modelling: coupling WRF-Chem and time-varying Eruption Source Parameters

Despite the knowledge gained with remote sensing data and the numerical modelling of volcanic eruptions, the identification of hazard areas around active volcanoes still remains a considerable scientific challenge [11]. The parameters that most influence the forecast of ash column loading distributions for mapping flight hazard areas are the free-troposphere turbulence levels, precipitation threshold for wet deposition, and eruption source parameters such as plume height and mass eruption rate (MER), and the onset and end times of the paroxysm [12,13,14]. Egan et al. [15] have implemented a simplified volcanic ash aggregation scheme into the Weather Research Forecasting model with Chemistry (WRF-Chem; Grell et al. [16]) for computing the aggregation rate coupled to the atmospheric environment at each model step. Harvey et al. [13] have listed the observation of plume height through mobile radars among the top research priorities for reducing uncertainties in forecasts of long-range transport of volcanic ash. WRF-Chem model has been recently utilised [17] to model the volcanic plumes associated with Etna’s paroxysmal episodes on 3–7 December 2015 that led to the closure of the nearby Catania International Airport. It emerged that the WRF–Chem model reasonably reproduced the distribution of SO₂ and volcanic ash, as compared with multi-platform data from satellite sensors and Doppler radars. More recently, Rizza et al. [18] have demonstrated the need for time-dependent Eruption Source Parameters (ESP) together with the correct specification of injection height to properly describe the transport of Etna volcanic ash in the Mediterranean basin. The online coupling between chemistry/aerosols with meteorology represents an innovative development of current volcanic ash transport and dispersion models (VATDM) as it avoids any space-time interpolation by considering the two-way feedback for all the processes (meteorology, chemistry, and physics) at the time step level, in a context of time-varying ESP [18].

The main aim of this paper is to document the large-scale contamination of airspace by a succession of closely interspersed explosive eruptions, providing insight into hazards to air traffic and atmospheric dispersion of volcanic pollutants. It also demonstrates the usefulness and high potential for nowcasting of ground-based near-vent radar observations to provide the most realistic source term estimates of volcanic ash that will undergo long-range transport. To these purposes, we use MER retrievals from a L-band Doppler radar (carrier frequency 1.274 GHz) monitoring Etna (VOLDORAD-2B, here noted V2B; [19, 20]) in input to WRF-Chem atmospheric model simulations to evaluate the regional-scale dispersion and accumulation of very fine volcanic ash released by successive paroxysms. This full coupling methodology, using measured time-evolutive source parameters, accurately reproduces the atmospheric dispersal from volcanic transient emissions, as attested by satellite imagery from various gas and ash sensors for the 23 November 2013 paroxysm [18]. The higher time resolution in the definition of the meteorological conditions and the specification of eruption conditions allows the reduction of related uncertainties present in the modelling of transport and dispersal of tephra in the atmosphere.

In the present study, we analyse two test-case eruptive sequences of Etna, the first one including 6 successive paroxysms from the New SouthEast crater (NSE) that occurred in 38 days between October 26 and December 2, 2013, the second sequence including 4 paroxysms from the Voragine crater (VOR) that occurred in less than 3 days between December 3 and 5, 2015. The methods and data used in this study are described in Section 2, which also describes the details of the two eruptive sequences, the volcanic emission estimates determined using the V2B radar data and the model setup. In section 3, we discuss the modelling results, in terms of meteorological and ash transport characteristics in the numerical model domain. Conclusions and aims for future works are presented in Section 4.

2. Materials and Methods

2.1. Description of the two eruptive sequences

2.1.1. The 2013 sequence: October 26 to December 2

The first simulation reproduces the sequence of paroxysms that occurred at the New South East (NSE) crater from 26 October to 2 December 2013. After a 6-month quiescence, a new paroxysm began on 26 October at NSE crater. For this event, the SO2 mass detected by the OMI (Ozone Monitoring Instrument) sensor onboard Aura satellite peaked at 0.77 kilotons (https://so2.gsfc.nasa.gov/index.html, accessed on 05 May 2023). The October-December sequence includes six paroxysms featuring usual lava fountain feeding many kilometres high tephra plumes, ranging in duration from 3.2 to 10 hours and from 2.9x10⁴ to 3.9x10⁵ kg s^-1 in average MER (peaks: 2.1x10⁵ to 8.9x10⁶ kg s^-1) as retrieved from VOLDORAD-2B radar (Table 1). Only the last three paroxysms punctually exceeded 10⁶ kg s^-1, the shortest and most violent of the series occurred under strong westerly wind on November 23 [18, 21, 22]. This paroxysm averaged nearly 4x10⁵ kg s^-1 during 193 min including 3.5x10⁶ kg s^-1 sustained during its 16 min climax phase which released 80% of the total erupted mass of pyroclasts at maximum velocities over 200 m s^-1 [20]. The six successive paroxysms released in total about 18.9 x 10⁹ kg of tephra with respective mass contributions of 18, 40, 24, 5, 8, and 5 %.

Plume top altitudes for this sequence (8.1–11.1 km a.s.l.) can be found in Corradini et al. [23]. The Toulouse VAAC Météo-France released 23 VAA for Etna during the period of the six paroxysms for flight levels up to FL380 (11.58 km a.s.l.) with wind speed estimates up to 35 m s^-1.

2.1.2. The 2015 sequence: December 3 to 5

An exceptionally rapid sequence of eruptive episodes occurred in December 2015 generating powerful eruptions columns during paroxysms, tall lava fountains, Strombolian activity, lava flow emission and persistent ash emission, involving all four summit craters at different times [24, 25]. This activity culminated with four paroxysmal events on December 3–5, 2015 from the Voragine (VOR) crater constituting the second sequence of our analysis (Table 2, Figure 2). As most paroxysms at Etna, they were preceded by Strombolian activity and displayed lava fountains and subsequent high tephra plumes. These paroxysms were short (less than 1.5 h) and not accompanied by lava flows and, as compared with the 2013 sequence, they were characterised by exceptional eruptive power and rapid occurrence rate. They rank among the most violent in the last two decades, in part due to the depth of their magma source region 1.5 km b.s.l. [26,27].

The first episode on December 3 was the most powerful, with lava fountain velocities inferred from V2B averaging 200 m s^-1 during the climax and peaking at 250–290 m s^-1. An ash plume ascent velocity of 5.6 m s^-1 (10 km in height above crater level in about 30 min) is inferred from the timing of the V2B near-source power and velocity increase and the plume top altitudes from MSG-SEVIRI HOTVOLC data and vertical temperature profiles [28], the maximum plume top altitude of 13.5 ± 1.5 km a.s.l. being reached around 03:00 UTC. The following paroxysms gradually decreased in MERs, total mass of erupted pyroclasts, and lava fountain average height [26].

These paroxysms caused the atmospheric dispersion and widespread ash fallout to more than 100 km away, with disruptive consequences to aviation and aeronautical infrastructures and services. In particular, the Catania international airport reported repeated closures of airspace sectors between 3 and 9 December 2015, and at least 57 cancelled flights, 56 diverted flights and tens of modified flights on 4–5 December. The Toulouse VAAC Météo-France released 29 VAA in 7 days (December 3–9) for flight levels up to FL350 (10660 m a.s.l.) and wind speeds up to 15 m s^-1, while INGV released 16 VONA in four days (December 2–5) reporting on these paroxysms.

The regional environmental impact of volcanic emissions from this short series of paroxysms can be apprehended from estimates of released products amounts. According to Bonaccorso and Calvari [26], the total production of pyroclasts was 10x10⁶ m³ DRE (Dense Rock Equivalent, i.e. about 25x10⁹ kg) with respective mass contributions of 41, 27, 19, and 13 % for the successive events. The V2B estimates point to 11x10⁹ kg with successive mass contributions of 84.1, 6.5, 4.8, and 4.6 %. Differences among estimates can stem from both the numerous hypotheses used by the former to derive magma volumes from strainmeter and camera measurements and from the possible underestimation by V2B of the last three weaker paroxysms due to the marginal position of the Voragine crater with respect to the radar beam. The SO₂ mass detected by the OMI sensor peaked at 20 kilotons in early December (https://so2.gsfc.nasa.gov/index.html, accessed on 16 May 2023).

The definition of eruption source parameters (ESP), namely the volcano location in the grid domain, start and end of each single paroxysm (UTC time), the total grain size distribution (TGSD), MER, and the injection height was pre-processed as described in the following paragraphs.

2.2. The VOLDORAD-2B (V2B) Doppler Radar System

2.2.1. Description

Between 2009 and 2022, INGV-OE monitored Etna tephra emissions using a volcano Doppler radar of OPGC (VOLDORAD-2B), jointly operated from the Montagnola station about 3 km south of the New SE Crater. This L-band (23.5 cm-wavelength) pulse radar probed 13 atmospheric volumes, 150 m-deep, right above the summit craters using a fixed-pointing antenna beam. The echo power and Doppler velocity parameters (level 1 data) in each range bin are computed at intervals of 0.22 s from 3 incoherent integrations of the Doppler spectra (raw data). These near-source measurements of emitted tephra are unique as they allow maximum ejection velocities at the base of the eruptive column and MER (level 2 products) to be determined, potentially in real-time [19,20] and to a satisfying level of uncertainty as revealed by comparisons with other sensors [29,30,31]. During the 2011–2015 period, the V2B radar successfully recorded about 50 paroxysms of vigorous lava fountains feeding tephra plumes up to over 13 km high. Level 1 data are available from an open-access database of the Observatoire de Physique du Globe de Clermont-Ferrand (http://voldorad.opgc.fr/bddtr.php, accessed on 16 June 2023; [32]), as well as more elaborated (level 2) products upon request, such as MER, onset and end of fountaining phase, and eruptive crater.

2.2.2. Data elaboration

We used the methodology of Freret-Lorgeril et al. [20] to compute the MER from the V2B data. They showed that the product of the measured echo power and radial velocities of rising tephra in the beam volumes located immediately above the eruptive crater is proportional to the MER (i.e., a MER proxy) and correlates to the observed plume heights for many paroxysms at Etna. Their calibration constant was then found by comparing the radar proxies with MER derived from a plume ascent model matching observed plume heights and including wind vertical profiles [33].

In this study, the MER was computed from range bins displaying the maximum power, at 3135 and 3285 m in slant distance from radar for paroxysms originating from the New SE Crater during the first sequence (26 October–2 December 2013, Table 1) and from 3 range gates for the second sequence of paroxysms originated from the Voragine Crater (Table 2). The V2B MER data were elaborated to be compatible with the WRF-Chem input system; in particular, they were averaged on a 30-second interval to be integrated in the numerical model at each time step. The corresponding time series of the MER are depicted in Figure 1 and Figure 2 for the 2013 and 2015 eruptive sequences respectively.

Figure 1. Mass Eruption Rate from V2B radar system ingested in WRF-Chem volcanic package for the 2013 sequence: (a) full sequence, (b) NSE2 paroxysm (c) NSE3 paroxysm (d) NSE4 paroxysm (e)NSE5 paroxysm, and (f) NSE6 paroxysm; Julian day 290 corresponds to October 17, 2013.

Figure 2. Mass Eruption Rate from V2B radar system ingested in WRF-Chem volcanic package for the 2015 sequence of paroxysms from Voragine crater: (a) full sequence; (b) VOR1; (c) VOR2 and VOR3; (d) VOR4. Julian day 335 corresponds to December 2, 2015.

2.3. Determination of the eruption height

In case of missing experimental data for the model injection height, the zero-order semi-empirical relationship of Mastin et al. [34] linking column height with MER observations can be utilised. It is expressed by the following eq.1 [35], and reported in Figure 3:

h = 5 \times 10^{- 4} {(\frac{M E R}{ρ})}^{0.241}

(1)

where h is expressed in metres (above the vent), ρ denotes the magma density (2600 kg m^-3) and MER is obtained by the above-described elaborations of the V2B radar data and reported in Figure 1 and Figure 2.

For the range of mass eruption rates considered in the two series of eruptions [0–6] x 10⁶ kg s^-1, the corresponding eruption height has an upper limit of 12–13 km (above the vent), in agreement with Mastin et al. [34] (their Figure 1).

The emission profile is specified according to the vertical profile derived from the umbrella cloud model in which 75% of the erupted mass is constrained in the umbrella cloud and the remaining 25% is linearly distributed from the umbrella base down to the vent [35].

2.4. Model description and setup

For this study, the WRF-Chem model version 4.3.1 has been utilised in a numerical domain covering the central Mediterranean region, with 310 × 320 grid points and a horizontal grid spacing of 6 km (1854 × 1914 km², Table 3). Boundary and initial conditions are at 1-degree resolution and provided from NCAR/NCEP Final Analysis from Global Forecast System (FNL from GFS). The first sequence simulation started at 00:00 UTC on October 25 and ended at 00:00 UTC on December 04, 2013. The second sequence run started at 00:00 UTC on December 2 and ended at 00:00 UTC on December 08, 2015.

The full set of parameterisations of the WRF-Chem model are reported in Table 3.

Based on the WRF setup recommended by Shi et al. [36] and Rizza et al. [37]the following physics schemes are utilised: the Mellor–Yamada–Janjic (MYJ; [38]) and the Eta similarity theory [39] schemes are used to describe the planetary boundary-layer (bl_pbl_physics = 2) and the surface layer (sf_sfclay_physics = 2) processes, respectively. The land surface exchange coefficients and the associated fluxes of momentum, heat, and humidity are represented by the Noah-MP Land Surface Model (sf_surface_physics = 4; [40]).

The radiative schemes are parameterized using the Goddard radiation model [41] for shortwave (ra_sw_physics = 5) and longwave (ra_lw_physics = 5) components. The microphysics parameterization considers the new one-moment Goddard four-class ice (4ICE) scheme developed by Lang et al. [42] and implemented in the early versions (3.4.1) of the WRF package (mp_physics = 7). This scheme has prognostic variables for cloud ice, snow, graupel, and hail.

In WRF-Chem, the total grain size distribution (TGSD) consists of ten bins of volcanic ash (vash_#) particles with a grain diameter range starting from 2 mm down to less than 3.9 μm [35] and expressed in mixing ratio units (

μ

g kg^-1). The corresponding mass fraction percentage is described by the E1 distribution and reported in Table 4. It has been recently introduced inside the WRF-Chem volcanic package by Rizza et al. [18] and derived from the Poret et al. [22] analysis of the experimental data from the November 23, 2013 Etna eruption.

3. Results and discussion: hazard maps

3.1. Atmospheric circulation during the events: main mid-tropospheric flows

In this section we describe the general synoptic situation in terms of prevailing mid-troposphere flows, i.e., at 500 hPa, which also represents an important flight level (see next section). In Figure 4 and Figure 5 are reported the geopotential height and the wind vectors derived by WRF-Chem simulations. Data are shown in correspondence of the eruption timeline for both the 2013 and the 2015 sequences; for each event, the chosen times are close to the paroxysms and are reported in Table 5.

Although studied as a unique sequence, the 2013 eruptions lasted 38 days and the atmospheric configurations continually changed; it nevertheless generated some recurrent flows. For some of the six 2013 eruptive events several similarities are in fact noticeable in terms of synoptic configuration and prevailing upper-level winds. For NSE2, NSE4 e NSE5 (Figure 4b, 4d, 4e, respectively) some upper-level troughs, although of different amplitude and depth, dominated the central Tyrrhenian favouring south-westerly winds on southern Italy and on the Sicily region in particular.

For NSE6 (Figure 4f) a trough is also visible, but located south of the Strait of Sicily, between the southern coasts of the region and north Africa. This configuration favoured the advection of south-easterly upper-level winds over eastern Sicily.

An upper-level ridge is visible on southern Mediterranean for the NSE3 case (Figure 4c); this structure is connected to a trough on the Balearics, and the main mid-tropospheric flows over north-eastern Sicily followed the curvature of the ridge, i.e., with currents toward east-southeast.

Finally, a different (from the previous ones) configuration is noticeable for NSE1 (Figure 4a), characterised by a wide high-pressure area over the whole Italian Peninsula. In this situation, the upper-level currents over Sicily were weaker and mainly directed towards the southwest, also driven by a lower pressure area over North Africa.

For the 2015 sequence (Figure 5), that lasted about 3 days, the main synoptic configuration remained quite stationary and has been characterised by a very wide high-pressure system persistently dominating the whole Mediterranean Basin, only slowly moving westward during the sequence. Some local (and weak) high- and low-pressure areas, visible on Figure 5, developed over the Mediterranean during this sequence, maintaining relatively weak mid-tropospheric currents; these winds flowed over north-eastern Sicily firstly toward east (Figure 5a, 5b) and later, after the formation of a local low over north Africa, toward north (Figure 5a, 5b).

To highlight some general differences between the two sequences it can be seen how the 2013 eruptive events were characterised by a succession of different baric configurations, most of the time characterised by the presence of low-pressure areas over the central Mediterranean. On the contrary the 2015 events, that occurred in a few days, were characterised by a prevailing and quite stationary wide high-pressure system, with smaller pressure gradients in the Mediterranean. The resulting circulations at 500 hPa were also quite different, with mid-tropospheric currents on average in the range 20-40 m s^-1 for the 2013 (Figure A2) and less than 10 m s^-1 for the 2015 (Figure A3) sequence, respectively. These results, based on WRF-Chem simulations, represent a robust meteorological basis for the evaluation of the regional-scale ash transport after volcanic eruptions (see next section).

3.2. Analysis of the volcanic ash dispersion patterns

The objective of this work is to evaluate the contamination of airspace by a succession of Mount Etna explosive paroxysms and the implications for hazards, in particular with respect to air traffic in the Mediterranean area. For this purpose, we define the vertically integrated (

\bar{p_{10}}

) mass, that is calculated considering the ash bins of the E1 distribution with diameter less than 10 μm (Table 4), as follows:

\bar{p_{10}} = \int_{0}^{t o p} < P M 10 > ρ_{a i r} d z

(2)

P M 10 = 0.5 \times v a s h_10 + v a s h_9 + 0.5 \times v a s h_8

(3)

where PM10 is the sum of the finest ash bins (units μg kg^-1) and

ρ_{a i r}

is the density of the air (kg m^-3). Brackets denote the time average over the simulation period and the vertical integration (represented by the overbar) is performed considering the whole vertical domain (top = 18 km). The resulting

\bar{p_{10}}

units are expressed in μg m^-2 which can be considered as a columnar ash concentration.

The model simulations produce hazard maps of

\bar{p_{10}}

to inform on the atmospheric ash dispersal of PM10 but it should be reminded that time- and spatially-averaged representations presented below are intended to guide aviation in terms of hazard appraisal and risk mitigation of Etna paroxysms. These representations should not be confused with instantaneous ash concentrations that could be met by aircrafts or measured by instruments at a given moment. To circumvent the lack of vertical information in columnar ash concentration representations, we derive ash concentrations at specific flight levels in Section 3.2.2 by performing a WRF-Chem data interpolation (3D) for a given pressure level (specialised function of NCAR Command Language named wrf_user_intrp3d). To further image the ash dispersal with a better time resolution, in Section 3.2.3 we also compute the daily average ash concentration at specific flight levels for four selected paroxysms.

3.2.1. Time-averaged columnar ash concentration maps

Figure 6 and Figure 7 show the results for the 2013 and the 2015 sequences, respectively. Both figures were built considering the time average on the full simulated sequences. For the 2013 sequence, a successive spatial filter is applied to leave out days without ash in the domain.

These two figures give a quantitative estimation of the global contamination of airspace by a succession of eruptions. The average spatial pattern outlined by the 2013 sequence (Figure 6) shows basically three affected regions in terms of vertically-integrated columnar ash concentrations: a red-shaded region around 10 g m^-2 within a radius of approximately 100 km, a yellow-shaded region of the order of 0.5 g m^-2 within a ~ 500 km radius, and a green/blue-shaded area with columnar ash content lower than 0.1 g m^-2 within a ~ 800 km radius. This average spatial pattern is overprinted by the dispersal from each event oriented by its specific wind conditions, conferring the starfish-like overall pattern.

The distribution of columnar ash concentrations from the 2015 short sequence (Figure 7) presents a similar starfish-like spatial pattern, with remarkable differences consisting in: (i) a more extended red-shaded region centred on Sicily and south Italy and also more intense in term of magnitude, within a radius of approximately 150 km from the eruptive crater and (ii) a strong eastward upper tropospheric transport caused by the VOR1 eruption.

As highlighted in the previous section, these differences are mainly related to the different mid-tropospheric transport conditions during the two sequences and also by the very rapid succession of eruptive episodes in December 2015, that generated powerful eruptions columns up to the stratosphere.

This analysis may be considered a preliminary indication of potential hazard zones around Mount Etna, which can be re-evaluated regularly as meteorological conditions change.

3.2.2. Time-averaged ash concentration maps at specified flight levels

To evaluate the contamination of airspace and the related hazards to aviation, it is necessary to point out that the forecasts by VAACs (https://www.ssd.noaa.gov/VAAC/vaac.html) use the following concentration thresholds to define the ensuing levels of ash contamination: low (200–2000 μg m⁻³), medium (2000–4000 μg m⁻³), and high (> 4000 μg m⁻³).

The flight level (FL), represents the altitude of an aircraft at the standard air pressure, expressed in hundreds of feet (1 FL = 100 ft). It may be considered as a proxy for elevation, and it may be associated with a given pressure level as reported in Table 6. These levels correspond to the upper airspace characterised by ICAO, as the zone where most air transport operations take place.

Figure 8 displays the time-averaged PM10 concentration for the 2013 (upper row) and 2015 (bottom row) sequences, at flight levels FL180 (panels a,d), FL240 (panels b,e) and FL300 (panels c,f) corresponding to the pressure levels of 500, 400 and 300 hPa (Table 6). The label-scale of Figure 8 is saturated at 4000 μg m^-3 (corresponding to 10^3.6) which represents the threshold defined by ICAO for the high level of contamination and the consequent prohibited airspace.

The transport of volcanic ash at these flight levels has different patterns for each respective sequence. The effects of the low-pressure weather conditions that occurred during the 2013 sequence is particularly evident in Figures 8a-c from the circular to elliptical stretched vortices (800 to over 1600 km in diameter) that transported volcanic ash counter-clockwise at large distances. This is particularly evident for the 2013 sequence for FL180 (Figure 8a) and FL240 (Figure 8b) flight levels at which ash remains concentrated mainly along these vortices in narrow bands generally less than 100 km wide. Ash concentrations in these bands range from about 1 to over 4000 μg m^-3. Being well localised under these high wind conditions, ash can remain concentrated over long distances, up to 25 μg m^-3 (10^1.4) still at 1000 km eastward for the most violent paroxysm of the sequence that occurred on 23 November 2013 (NSE4). At FL300 (300 hPa) in the upper troposphere, vortices are still present, but ash is more widely dispersed in the north direction owing to tropospheric wind at 300 hPa (Figure 8c).

Contrastingly, the 2015 sequence, occurring in prevalently high-pressure conditions and weak winds (see Figure A3), caused high ash concentration at all flight levels, as outlined by the red shaded area in Figures 8d-f. At FL300 (Figure 8f), mostly the intense eastward transport of very fine ash from VOR1 by tropospheric wind at 300 hPa is visible, due to the combination of the eruption intensity and high-pressure low wind conditions. Considering this is a time-averaged map, this potentially indicates high ash hazards at all flight levels during - and for days following - this eruptive sequence. The concentration of very fine ash appears more localised spatially going upward into the upper troposphere, although still reaching no fly zone level.

What can be inferred from this analysis on the PM10 maps at these flight levels is that the long-range transport is strictly correlated with the synoptic circulation. In particular, synoptic patterns with relatively weak tropospheric currents are conditions leading to the accumulation of volcanic ash in the atmosphere after an eruption. This seems a rather obvious conclusion, but it is important to highlight the fundamental role played by the meteorological component in VATDM models and the need for their online coupling with the numerical aerosol packages.

This result is in part confirmed by analysing the VAA issued by Toulouse-VAAC that are reported in Table A1: a total of 19 VAA were issued in the period December 3 to 6, of which 12 were classified as “Red”, 6 “Orange”, and 1 “Yellow”.

3.2.3. Daily average ash concentration maps of selected paroxysms

In the following Figure 9 and Figure 10 we have analysed the daily average (from 00:00 UTC until 23:59) ash concentrations of two powerful paroxysms of the first sequence having the most contrasted durations on November 23 (NSE4, 3.2h) and 28 (NSE5, 8.3h), and the most contrasted in MER for the second sequence, on December 3 (VOR1) and 5 (VOR3-4). It is calculated using the following quantity:

p_{10 D} = [P M 10]

(4)

where square brackets denote daily time average on the selected day. The left column of both figures (panels a-d) refers to FL180, the central to FL240 (panels b-e) and the right column to FL300 (panels c-f). The contour shading matches the above-mentioned four hazard levels used by aviation for volcanic ash concentration reported in Table 7.

For the NSE4 paroxysm, the geopotential height at 500 hPa (Figure 4d) shows a wide low-pressure area on the central Mediterranean basin favouring south-westerly strong upper-tropospheric winds (> 30 ms^-1) on southern Italy. The very fine ash was rapidly transported toward Albania (in about 6h) and the Balkan regions at relatively low concentrations (Figure 9a).

The very fine ash after the NSE5 paroxysm (Figure 9d,e,f,) was transported eastward in a straight line by the intense zonal currents at 500 hPa (Figure 4e).

These maps indicate low potential risks for aviation at FL180 to FL300. The shading contours (blue-green) indicate in fact the low hazard level. On the other side all the 21 VAA by Toulouse VAAC reported in Table A1 indicate “NIL” alert level (NIL = no ash cloud is produced) for the whole period.

The corresponding analysis of the 2015 sequence is depicted in Figure 10. Panels a-b-c in the upper row report the daily p_10D concentration at FL180, FL240, and FL300 for December 3, 2015 (VOR1), while the lower panels (d-e-f) describe the equivalent for December 5, 2015 (VOR3, and VOR4).

For both cases, the p_10D daily patterns, in agreement with the high pressure meteorological conditions (Figure 5a,d), indicate a red-shaded region in proximity of Mount Etna for FL180 and FL240 (panels a-b-c-d). FL300 (Figure 10c) relative to VOR1 again reveals the presence of p_10D ash fraction at upper tropospheric level (about 9200 m) in the form of wave-like pattern which is the artificial effect of the averaging (daily) of pulsed ash plumes. To better visualise these maps, a high-resolution Figure A1 is provided in the appendix with a zoom on the Mount Etna area. This is consistent with the 19 VAA issued by the Toulouse VAAC during the period 3 to 6 December 2015, 12 of which were classified as “Red” and 6 as “Orange”.

3.2.4. Vertical distribution of very fine ash

Animated videos (gif files) showing the fast volcanic PM10 space-time evolution of the successive emissions (often 'pulse-like') due to the transient nature of their source can be found at https://doi.org/10.25519/GG1K-2H34 for the 2013 and 2015 sequences. Figure 11 are extracted snapshots for different paroxysms informing on the vertical distribution of PM10 along cross-sections indicated in Figure 6 and Figure 7. The ash cloud thickness is a particularly valuable information for aviation hazard assessment and countermeasures in case of ash encounter, and for satellite mass retrievals as generally a non-measured average thickness is assumed. For instance, the ash cloud inherited from the 26 Oct. 2013 tephra plume (NSE1) is seen slowly drifting to the SW after about 1.5 day between 1.9 and 3.7 km a.s.l., i.e., with an inferred vertical thickness of about 1.8 km (Figure 11a). After the top part of the VOR1 plume had raised above 12 km a.s.l. in the minutes following the eruption onset, some very fine ash from the middle region of the eruptive column is observed after 7 h rising to the NE above 11.8 km a.s.l. with an inferred vertical thickness of about 2.5–3 km (Figure 11c), while the lower part of the eruptive column feeds another PM10 cloud, 700 m thick, drifting at about 2 km a.s.l. laterally to the SW 1.5 day after the onset (Figure 11d). The modelled ash cloud from the 4 Dec 2015 evening paroxysm has a vertical thickness of 1.4 km at an altitude of 5 km a.s.l. when the following paroxysm occurs (Figure 11f). Note that an ongoing vertical plume emission coexisting with ash from a previous plume is visible in Figure 11e–f.

Another interesting feature seen in Figure 11b,e is ash stagnating for several days at ground level. PM10 concentrations are commonly above 500 μg m^-3, as during four days in late October to early November 2013 following NSE1. PM10 concentrations at ground level can even punctually exceed the aviation threshold as observed following VOR2 around the occurrence time of VOR3. The repeated ash emissions of Etna recurrently exposing surrounding inhabitants, animals, and plants to volcanic PM10 and gas may therefore potentially represent health hazards over the long term.

4. Conclusions

We evaluate for the first time the contamination of airspace by very fine ash from a succession of ash plumes produced by Mount Etna and the related hazards to aviation in the Mediterranean Basin. To this aim, the atmospheric dispersal of sub-10 μm (PM10) ash was simulated using the WRF-Chem model coupled online with meteorology and aerosols and offline with Mass Eruption Rates (MER) derived from near-vent Doppler radar measurements and inferred plume altitudes. The pioneering direct input of systematic observations, best resolved in time and thus improving the model accuracy, is all the more so valuable as paroxysms are transients involving time-varying eruption source parameters. PM10 ash concentration hazard maps are provided for two sequences of paroxysms with widely varied volcanological conditions and contrasted occurrence interval and meteorological synoptic patterns: 6 ash plumes produced in 38 days during low-pressure, high-wind conditions in October - December 2013, and 4 ash plumes released in 3 days during high-pressure, low-wind conditions on 3-5 December 2015. Average MERs covered two orders of magnitude (2.9 x 10⁴ – 1.7 x 10⁶ kg s^-1), durations spanned from less than 1 hour to nearly half a day, and plume top altitudes ranged from 8.1 to 13.5 km a.s.l. The synoptic conditions at 500 hPa were quite different, with mid-tropospheric current speeds on average in the range 20–40 ms^-1 and less than 10 ms^-1 for the 2013 and 2015 sequences, respectively.

We analyse the ash PM10 dispersal in terms of time-averaged columnar ash concentration (

\bar{p_{10}}

), of concentration at specified flight levels (FL180, FL240, and FL300) averaged over the entire sequence interval (PM10), and of daily average concentration (p_10D) during selected paroxysm days at these flight levels. Consistently with the Volcanic Ash Advisories issued by the Toulouse VAAC, they may be considered as first-approach model hazard maps to aircrafts and, although certainly not covering all possible volcanological and meteorological configurations for hazard appraisal at Etna, the following conclusions can be drawn:

The very fine ash from sequences of Mount Etna paroxysms is shown to easily contaminate the airspace around the volcano within a radius of about 1000 km in a matter of days. The airspace from many countries around the Mediterranean Basin is impacted, including most of southern Europe from the Balearic Islands westward, south of France, the whole Italy, Greece and western coast of Turkey and the Balkans eastward, to beyond the Alps northward and to Malta and the African northern coast (from Algeria to Libya) southward.
Low-pressure weather systems favour the trapping and circulation of very fine ash in the whole troposphere within this area, yet at low concentration, generally below 1 μg m^-3. In this meteorological context, high winds tend to stretch ash clouds into ~100 km-wide clouds forming large-scale vortices 800–1600 km in diameter where PM10 ash concentrations can still exceed the aviation hazard threshold up to 1000 km downwind from the volcano, a distance reached in about 10h (e.g., NSE5).
High-pressure, low-wind conditions tend to favour the accumulation of PM10 ash in a wide atmospheric region surrounding Etna. In this context, closely interspersed paroxysms tend to accumulate very fine ash, more diffusively in the lower troposphere and in stretched ash clouds higher up in the troposphere.
In all volcanological and meteorological configurations simulated, the lower troposphere appears particularly prone to accumulate diffuse PM10 ash during sequences of eruption; this is likely to particularly affect take-off and landing aircrafts in regional airports.
High MER paroxysms propel ash up to the upper troposphere, where most of the air traffic occurs, and sometimes also into the lower stratosphere, whichever the weather conditions. High-troposphere ash clouds from Etna appear as a pulsed feature resulting mostly from the short-lived climax phase.
Daily average PM10 tropospheric ash concentrations commonly exceed the aviation hazard threshold, up to 1000 km downwind from the volcano and up to the upper troposphere for intense paroxysms.
The thickness of modelled PM10 ash clouds generated from different parts of the eruptive columns ranges from 0.7 to 3 km.
Potential health hazards may stem from the stagnation for PM10 ash at ground level for several days, commonly above 500 μg m^-3, and sometimes punctually exceeding the aviation threshold.

Improving the accuracy of ash dispersal models is necessary both for nowcasting and forecasting purposes and to compensate for ash detection gaps from monitoring systems. Future WRF-Chem model developments will include refining injection heights and the initial columnar ash distribution. The semi-empirical law from Mastin et al. [34]’s dataset used in our study averages all observed ash plume heights from mainly large silicic explosive eruptions that have occurred in a variety of wind conditions, tending to overestimate height from a given MER (Freret-Lorgeril et al. [20]; their Figure 7). Other models relating MER and plume height advantageously consider the vertical wind profile (e.g., Degruyter and Bonadonna [33]). Wind profiles may come from the meteorological data online-coupled to WRF-Chem and cross-checked from atmospheric balloons measurements achieved several times a day in the region around Etna. Additional datasets for wind profiles at high resolution (30 km) could be considered, such as the ERA5 data available from the Climate Data Store in netCDF format [44]. Plume heights controlled by wind conditions specific to the eruption could thus be systematically computed from the radar-derived MER. This would lead to a more accurate definition of the injection height in the simulations.

As for future research on volcanic ash forecasts, multi-model multi-source term ensemble approaches seem to represent a promising development, as explored for instance by Plu et al. [45]. They suggest that quantiles of ash concentrations can be relevant products for air traffic management, that can be used for route optimization in the areas where ash does not pose a direct and urgent threat to aviation. Cost and disruption of air traffic could be eliminated to a great extent by including the results of dispersion models into flight planning software to apply cost-based trajectory optimizations.

Supplementary Materials

The following supporting information can be downloaded at: https://doi.org/10.25519/GG1K-2H34, Video S1: Ash dispersion modelling from successive paroxysms of Etna using WRF-Chem and near-source VOLDORAD radar measurements – Dataset.

Author Contributions

Conceptualization, F.D. and U.R.; methodology, F.D.; software, M.M., C.B. and E.M. ; validation, U.R., G.C., A.S. and E.A.; formal analysis, U.R. and F.D.; investigation, X.X.; resources, G.P. and S.M.; data curation, F.D.; writing—original draft preparation, U.R.; writing—review and editing, F.D.; visualization, M.M.; supervision, U.R.; project administration, S.M. and G.P.; funding acquisition, F.D. All authors have read and agreed to the published version of the manuscript.”

Funding

This research received no external funding

Data Availability Statement

Not Applicable.

Acknowledgments

Measurements at Etna using the OPGC VOLDORAD-2B radar were carried out in the frame of a collaborative research agreement between the Observatoire de Physique du Globe de Clermont-Ferrand (OPGC, Université Clermont Auvergne, Clermont-Ferrand, France), the French CNRS, and the Instituto Nazionale di Geofisica e Vulcanologia, Osservatorio Etneo, sezione di Catania (INGV-OE). We are deeply indebted to Claude Hervier, Patrick Freville, Yannick Guéhenneux, and Philippe Cacault at OPGC as well as INGV-OE staff including M. Coltelli and M. Prestifilippo, for their support in radar monitoring and data management. This study used the VOLDORAD open-access database (http://voldorad.opgc.fr/ (accessed on 22 Dec 2022)) of OPGC—Université Clermont Auvergne, with support from the EPOS, EUROVOLC (#731070, including research funds for F.D. and contract funding for C.B.) and MED-SUV (#308665, F.D.) European programs, and from the Service National d’Observation en Volcanologie of the French CNRS-INSU (F.D.). We thank V. Freret-Lorgeril for providing useful comments on an early version of the manuscript.

Appendix A

Table A1. VAA issued by the Toulouse VAAC MeteoFrance for 2013 and 2015 eruptive sequences; NIL = no ash cloud is produced; green =normal; yellow=advisory; orange=watch, red=warning.

N	VA ADVISORY (UTC)	AVIATION COLOUR CODE	ERUPTION DETAILS
1	ETNA - 2013-10-26 06:00	NIL	STARTED AT 0200Z
2	ETNA - 2013-10-26 06:00	NIL	STARTED AT 0200Z
3	ETNA - 2013-10-26 11:30	NIL	STOPPED AROUND 1100Z
4	ETNA - 2013-10-26 17:30	NIL	STOPPED AROUND 1100Z
5	ETNA - 2013-10-27 08:13	NIL	CLOUD SEEMS TO BE COMPOSED OF WATER VAPOUR, NO SIGNAL OF ASH NEITHER VOLCANIC GAS ON SAT IMAGERY
6	ETNA - 2013-10-28 16:00	NIL	VA NOT IDENTIFIABLE
7	ETNA - 2013-10-29 10:50	NIL	UNKNOWN
8	ETNA - 2013-11-11 03:43	NIL	UNKNOWN
9	ETNA - 2013-11-11 10:11	NIL	CLOUD IDENTIFIABLE ON WEBCAM MAY CONTAIN VA
10	ETNA - 2013-11-16 23:03	NIL	IN PROGRESS LOW INTENSITY
11	ETNA - 2013-11-17 02:17	NIL	GOING ON
12	ETNA - 2013-11-17 05:09	NIL	ERUPTION STOPPED ABOUT 0500Z
13	ETNA - 2013-11-17 11:29	NIL	ERUPTION ENDED
14	ETNA - 2013-11-23 10:07	NIL	ASH CLOUD OF SEVERE INTENSITY STARTS AT 0930Z
15	ETNA - 2013-11-23 11:17	NIL	ERUPTION ENDED AT 1030Z
16	ETNA - 2013-11-23 14:26	NIL	ERUPTION ENDED AT 1030Z
17	ETNA - 2013-11-23 20:19	NIL	ERUPTION ENDED
18	ETNA - 2013-11-28 17:37	NIL	ERUPTION HAS STARTED AT 1730Z, GOING ON
19	ETNA - 2013-11-28 19:50	NIL	ERUPTION HAS RESTARTED AT 1930Z, GOING ON
20	ETNA - 2013-11-28 23:45	NIL	ERUPTION HAS STOPPED AT 2330Z
21	ETNA - 2013-12-02 19:00	NIL	ERUPTION HAS STARTED AT 1820Z, GOING ON
22	ETNA - 2013-12-02 23:40	NIL	ERUPTION ENDED AT 2300Z

1	ETNA - 2015-12-02 23:00	ORANGE	STROMBOLIAN ACTIVITY
2	ETNA - 2015-12-03 02:41	RED	EXPLOSIVE ERUPTION OCCURRED AT 02000Z
3	ETNA - 2015-12-03 04:00	ORANGE	ERUPTION AND ASH EMISSION DECREASING. VA IDENTIFIABLE FM SATELLITE IMAGERY
4	ETNA - 2015-12-03 10:00	ORANGE	SMALL ACTIVITY IN VICINITY OF VOLCANO
5	ETNA - 2015-12-03 14:00	YELLOW	NO SIGNIFICANT ASH EMISSION
6	ETNA - 2015-12-04 09:40	RED	ASH EMISSION VISIBLE ON WEBCAM
7	ETNA - 2015-12-04 10:45	RED	ERUPTION ONGOING.
8	ETNA - 2015-12-04 15:45	RED	ASH PLUME NEAR SUMMIT VOLCANO
9	ETNA - 2015-12-04 21:00	RED	ERUPTION STARTED AT 2045Z
10	ETNA - 2015-12-05 03:00	RED	STROMBOLIAN EXPLOSIONS
11	ETNA - 2015-12-05 08:45	RED	STROMBOLIAN EXPLOSIONS
12	ETNA - 2015-12-05 15:00	ORANGE	SOME VOLCANIC ASH NEAR THE SUMMIT.
13	ETNA - 2015-12-05 15:05	RED	SIGNIFICANT EMISSION OF ASH OVER THE VOLCANO
14	ETNA - 2015-12-05 20:31	ORANGE	EXPLOSIVE ACTIVITY AND SIGNIFICANT ASH EMISSION STOPPED
15	ETNA - 2015-12-06 11:45	ORANGE	VOLCANIC ASH NEAR THE SUMMIT
16	ETNA - 2015-12-06 12:00	ORANGE	VOLCANIC ASH NEAR THE SUMMIT
17	ETNA - 2015-12-06 13:00	RED	INCREASING ACTIVITY
18	ETNA - 2015-12-06 17:55	RED	ONGOING MODERATE INTENSITY
19	ETNA - 2015-12-06 23:00	RED	EXPLOSIVE ACTIVITY
20	ETNA - 2015-12-07 03:00	RED	ACTIVITY STILL ONGOING
21	ETNA - 2015-12-07 09:00	ORANGE	ACTIVITY STILL ONGOING
22	ETNA - 2015-12-07 15:00	RED	VA NOT IDENTIFIABLE FM SAT DATA, WINDS FL100 280/10KT FL300 290/15KT
23	ETNA - 2015-12-07 20:50	RED	SPORADIC ERUPTIONS STILL GO ON
24	ETNA - 2015-12-08 03:00	RED	SPORADIC ERUPTIONS STILL GO ON
25	ETNA - 2015-12-08 09:00	ORANGE	STROMBOLIAN ACTIVITY HAS DECREASED
26	ETNA - 2015-12-09 04:00	ORANGE	ERUPTION STILL GOING ON
27	ETNA - 2015-12-09 09:17	RED	ERUPTION STILL GOING ON
28	ETNA - 2015-12-09 14:50	RED	SEEMS DECREASING
29	ETNA - 2015-12-09 20:45	ORANGE	WEAK ERUPTIVE ACTIVITY IS ONGOING

Figure A1. Daily average of p_10D particle distribution for the December 2015 sequence at flight levels FL180 (5500 m a.s.l.; left panels) and FL240 (7300 m a.s.l.; right panels) on 3 December (top row a,b), 4 December (middle row c,d), 5 December (bottom row e,f). Units are μg m^-3, label bar is reported in log10 scale with intervals matching the hazard levels for aviation.

Figure A2. WRF wind speeds (m s^-1) and wind vectors at 500 hPa for the different 2013 paroxysms, namely (a) NSE1, (b) NSE2, (c) NSE3, (d) NSE4, (e) NSE5, and (f) NSE6.

Figure A3. WRF wind speeds (m s^-1) and wind vectors at 500 hPa for the different 2015 paroxysms, namely (a) VOR1, (b) VOR2, (c) VOR3, and (d) VOR4.

References

Prata, F.; Rose, B. Volcanic ash hazards to aviation. In The Encyclopedia of Volcanoes, 2nd ed.; Sigurdsson H., Ed.; Academic Press: 2015, pp. 911–934, Elsevier. [CrossRef]
International Civil Aviation Organization (ICAO). Flight safety and volcanic ash, 2012. Doc 9974, ANB/487. https://www.icao.int/publications/Documents/9974_en.pdf.
Jiménez-Escalona, J.C.; Poom-Medina, J.L.; Roberge, J.; Aparicio-García, R.S.; Avila-Razo, J.E.; Huerta-Chavez, O.M.; Da Silva, R.F. Recognition of the Airspace Affected by the Presence of Volcanic Ash from Popocatepetl Volcano Using Historical Satellite Images. Aerospace 2022, 9, 308. [Google Scholar] [CrossRef]
Marchese, F.; Falconieri, A.; Filizzola, C.; Pergola, N.; Tramutoli, V. Investigating Volcanic Plumes from Mt. Etna Eruptions of December 2015 by Means of AVHRR and SEVIRI Data. Sensors 2019, 19, 1174. [Google Scholar] [CrossRef]
Eychenne, J.; Gurioli, L.; Damby, D.; Belville, C.; Schiavi, F.; Marceau, G.; Szczepaniaks, C.; Blavignacs, C.; Laumonier, M.; Gardes, E.; Le Pennec, J.-L.; Nedelec, J.-M.; Blanchon, L.; Sapin, V. Spatial distribution and physicochemical properties of respirable volcanic ash from the 16-17 August 2006 Tungurahua eruption (Ecuador), and alveolar epithelium response in-vitro. GeoHealth 2022, 6, p.e2022GH000680. [CrossRef]
Sellitto, P.; Salerno, G.; Corradini, S.; Xueref-Remy, I.; Riandet, A.; Bellon, C.; et al. Volcanic emissions, plume dispersion, and downwind radiative impacts following Mount Etna series of eruptions of February 21–26, 2021. J. Geophys. Res. Atmos. 2023, 128, e2021JD035974. [Google Scholar] [CrossRef]
EUROCONTROL - 11 years after the eruption of Icelandic volcano Eyjafjallajökull. Available online: https://www.eurocontrol.int/news/11-years-after-eruption-icelandic-volcano-eyjafjallajokull (accessed on 23 June 2023).
Hirtl, M.; Arnold, D.; Baro, R.; Brenot, H.; Coltelli, M.; Eschbacher, K.; Hard-Stremayer, H.; Lipok, F.; Maurer, C.; Meinhard, D. and Mona, L. A volcanic-hazard demonstration exercise to assess and mitigate the impacts of volcanic ash clouds on civil and military aviation. Nat. Hazards Earth Syst. Sci. 2020, 20(6), 1719–1739. [CrossRef]
Guffanti, M.; Casadevall, T.J.; Budding, K. Encounters of aircraft with volcanic ash clouds: A compilation of known incidents, 1953–2009: U.S. Geological Survey Data Series, 2010, 545, ver. 1.0, 12 p., plus 4 appendixes including the compilation database. https://pubs.usgs.gov/ds/545/.
Christmann, C.; Nunes, R.R.; Schmitt, A.R. Recent encounters of aircraft with volcanic ash clouds. Deutscher Luft- und Raumfahrtkongress (DLR), 2015, DocumentID: 370124. https://www.dglr.de/publikationen/2015/370124.pdf.
Thouret J.-C.; Charbonnier S. Assessment, Delineation of Hazard Zones and Modeling of Volcanic Hazards. In Hazards and Monitoring of Volcanic Activity 1, J.-F. Lénat Ed. ISTE-Géosciences, 2022, ISBN 978-1-78945-043-9, 151–184. [CrossRef]
Folch, A. A review of tephra transport and dispersal models: evolution, current status, and future perspectives. J. Volcanol. Geotherm. Res. 2012, 235-236, 96–115. [CrossRef]
Harvey, N.J.; Huntley, N.; Dacre, H.F.; Goldstein, M.; Thomson, D. and Webster, H. Multi-level emulation of a volcanic ash transport and dispersion model to quantify sensitivity to uncertain parameters. Nat. Hazards Earth Syst. Sci. 2018, 18(1), 41–63. [CrossRef]
Aubry, T.J.; Engwell, S.; Bonadonna, C.; Carazzo, G.; Scollo, S.; Van Eaton, A.R.; Taylor, I.A.; Jessop, D.; Eychenne, J.; Gouhier, M.; Mastin, L.G.; Wallace, K.L.; Biass, S.; Bursik, M.; Grainger, R.G.; Jellinek, A.M.; Schmidt, A. The Independent Volcanic, Eruption Source Parameter Archive (IVESPA, versio 1.0): a new observational database to support explosive eruptive column model validation and development. J. Volcanol. Geotherm. Res. 2021, 417, 107295. [CrossRef]
Egan, S. D.; Stuefer, M.; Webley, P. W.; Lopez, T.; Cahill, C. F. and Hirtl, M. Modeling volcanic ash aggregation processes and related impacts on the April–May 2010 eruptions of Eyjafjallajökull volcano with WRF-Chem. Nat. Hazards Earth Syst. Sci. 2020, 20, 2721–2737. https://doi.org/10.5194/nhess-20-2721-2020. [CrossRef]
Grell, G.A.; Peckham, S.E.; Schmitz, R.; McKeen, S.A.; Frost, G.; Skamarock, W.C.; Eder, B. Fully coupled “online” chemistry within the WRF model. Atmos. Environ. 2005, 39, 6957–6976. [Google Scholar] [CrossRef]
Rizza, U.; Brega, E.; Caccamo, M.T.; Castorina, G.; Morichetti, M.; Munaò, G.; Passerini, G.; Magazù, S. Analysis of the ETNA 2015 Eruption Using WRF–Chem Model and Satellite Observations. Atmosphere 2020, 11, 1168. [Google Scholar] [CrossRef]
Rizza, U.; Donnadieu, F.; Magazu, S.; Passerini, G.; Castorina, G.; Semprebello, A.; Morichetti, M.; Virgili, S.; Mancinelli, E. Effects of Variable Eruption Source Parameters on Volcanic Plume Transport: Example of the 23 November 2013 Paroxysm of Etna. Remote Sens. 2021, 13, 4037. [Google Scholar] [CrossRef]
Donnadieu, F.; Freville, P.; Hervier, C.; Coltelli, M.; Scollo, S.; Prestifilippo, M.; Valade, S.; Rivet, S.; Cacault, P. Near-source Doppler radar monitoring of tephra plumes at Etna. J. Volcanol. Geotherm. Res. 2016, 312, 26–39. [Google Scholar] [CrossRef]
Freret-Lorgeril, V.; Donnadieu, F.; Scollo, S.; Provost, A.; Fréville, P.; Guéhenneux, Y.; Hervier, C.; Prestifilippo, M.; Coltelli, M. Mass Eruption Rates of Tephra Plumes during the 2011–2015 Lava Fountain Paroxysms at Mt. Etna From Doppler Radar Retrievals. Front. Earth Sci. 2018, 6, 73. [Google Scholar] [CrossRef]
Bonaccorso, A.; Calvari, S.; Linde, A.; Sacks, S. Eruptive processes leading to the most explosive lava fountain at Etna volcano: The 23 November 2013 episode. Geophys. Res. Lett. 2014, 41, 4912–4919. [Google Scholar] [CrossRef]
Poret, M.; Corradini, S.; Merucci, L.; Costa, A.; Andronico, D.; Montopoli, M.; Vulpiani, G.; Freret-Lorgeril, V. Reconstructing volcanic plume evolution integrating satellite and ground-based data: Application to the 23 November 2013 Etna eruption. Atmos. Chem. Phys. Discuss. 2018, 18, 4695–4714. [Google Scholar] [CrossRef]
Corradini, S.; Guerrieri, L.; Lombardo, V.; Merucci, L.; Musacchio, M.; Prestifilippo, M.; Scollo, S.; Silvestri, M.; Spata, G.; Stelitano, D. Proximal Monitoring of the 2011–2015 Etna Lava Fountains Using MSG-SEVIRI Data. Geosciences 2018, 8, 140. [Google Scholar] [CrossRef]
Corsaro, R.A.; Andronico, D.; Behncke, B.; Branca, S.; Caltabiano, T.; Ciancitto, F.; Cristaldi, A.; De Beni, E.; La Spina, A.; Lodato, L.; Miraglia, L.; Neri, M.; Salerno, G.; Scollo, S.; Spata, G. Monitoring the December 2015 summit eruptions of Mt. Etna (Italy): Implications on eruptive dynamics. J. Volcanol. Geotherm. Res. 2017, 341, 53-69, ISSN 0377-0273. [CrossRef]
Global Volcanism Program. Report on Etna (Italy) (Crafford, A.E., and Venzke, E., eds.). Bulletin of the Global Volcanism Network, 2017, 42:5. Smithsonian Institution. [CrossRef]
Bonaccorso, A.; Calvari, S. A new approach to investigate an eruptive paroxysmal sequence using camera and strainmeter networks: Lessons from the 3–5 December 2015 activity at Etna volcano. Earth Planet. Sc. Lett. 2017, 475, 231-241, ISSN 0012-821X. [CrossRef]
Bonforte, A.; Cannavo, F.; Gambino, S.; Guglielmo, F. Combining High- and Low-Rate Geodetic Data Analysis for Unveiling Rapid Magma Transfer Feeding a Sequence of Violent Summit Paroxysms at Etna in Late 2015. Appl. Sci. 2021, 11(10), 4630. [Google Scholar] [CrossRef]
Donnadieu, F.; Freret-Lorgeril, V.; Gouhier, M.; Coltelli, M.; Scollo, S.; Fréville, P.; Hervier, C.; Prestifilippo, M. The 3 December 2015 paroxysm of Voragine crater at Etna: insights from Doppler radar measurements. Geophys. Res. Abstracts, 3 December 2016; -1. [Google Scholar]
Marzano, F.S.; Mereu, L.; Scollo, S.; Donnadieu, F.; Bonadonna, C. Tephra Mass Eruption Rate from X-Band and L-Band Microwave Radars during the 2013 Etna Explosive Lava Fountain. IEEE Trans. Geosci. Remote Sens. 2020, 58, 3314–3327. [Google Scholar] [CrossRef]
Freret-Lorgeril, V.; Bonadonna, C.; Corradini, S.; Donnadieu, F.; Guerrieri, L.; Lacanna, G.; Marzano, F.; Mereu, L.; Merucci, L.; Ripepe, M.; et al. Examples of Multi-Sensor Determination of Eruptive Source Parameters of Explosive Events at Mount Etna. Remote Sens. 2021, 13, 2097. [Google Scholar] [CrossRef]
Mereu, L.; Scollo, S.; Bonadonna, C.; Donnadieu, F.; Freret-Lorgeril, V.; Marzano, F. S. Ground-Based Remote Sensing and Uncertainty Analysis of the Mass Eruption Rate Associated With the 3–5 December 2015 Paroxysms of Mt. Etna. IEEE J. Sel. Top. Appl. Earth Obs. Remote Sens. 2022, 15, 504–518. [Google Scholar] [CrossRef]
Donnadieu, F.; Freville, P.; Rivet, S.; Hervier, C.; Cacault, P. The Volcano Doppler Radar Data Base of Etna (VOLDORAD-2B); Université Clermont Auvergne CNRS: Clermont-Ferrand, France, 2015. [Google Scholar] [CrossRef]
Degruyter, W.; Bonadonna, C. Improving on mass flow rate estimates of volcanic eruptions. Geophys. Res. Lett. 2012, 39:L16308. [CrossRef]
Mastin, L.G.; Guffanti, M.; Servranckx, R.; Webley, P.; Barsotti, S.; Dean, K.; et al. A multidisciplinary effort to assign realistic source parameters to models of volcanic ash-cloud transport and dispersion during eruptions. J. Volcanol. Geotherm. Res. 2009, 186(1-2), 10-21. [CrossRef]
Stuefer, M.; Freitas, S.R.; Grell, G.; Webley, P.; Peckham, S.; McKeen, S.A.; Egan, S.D. Inclusion of ash and SO2 emissions from volcanic eruptions in WRF-Chem: Development and some applications. Geosci. Model Dev. 2013, 6, 457–468. [Google Scholar] [CrossRef]
Shi, J.J.; Matsui, T.; Tao, W.K.; Tan, Q.; Peters-Lidard, C.; Chin, M.; Pickering, K.; Guy, N.; Lang, S.; Kemp, E.M. Implementation of an aerosol–cloud-microphysics–radiation coupling into the NASA unified WRF: Simulation results for the 6–7 August 2006 AMMA special observing period. Q. J. R. Meteorol. Soc. 2014, 140, 2158–2175. [Google Scholar] [CrossRef]
Rizza, U.; Avolio, E.; Morichetti, M.; Di Liberto, L.; Bellini, A.; Barnaba, F.; Virgili, S.; Passerini, G.; Mancinelli, E. On the Interplay between Desert Dust and Meteorology Based on WRF-Chem Simulations and Remote Sensing Observations in the Mediterranean Basin. Remote Sens. 2023, 15, 435. [Google Scholar] [CrossRef]
Janjic, Z.I. The Step–Mountain Eta Coordinate Model: Further developments of the convection, viscous sublayer, and turbulence closure schemes. Mon. Wea. Rev. 1994, 122, 927–945. [Google Scholar] [CrossRef]
Janjic, Z.I. The surface layer in the NCEP Eta Model. Eleventh conference on numerical weather prediction, Norfolk, VA, 19–23 August 1996. Amer. Meteor. Soc. 1996, 354–355, Boston, MA.
Niu, G.-Y.; Yang, Z.-L.; Mitchell, K.E.; Chen, F.; Ek, M.B.; Barlage, M.; Kumar, A.; Manning, K.; Niyogi, D.; Rosero, E.; Tewari, M.; Xia, Y. The community Noah land surface model with multiparameterization options (Noah–MP): 1. Model description and evaluation with local–scale measurements. J. Geophys. Res. 2011, 116, D12109. [CrossRef]
Chou, M.D.; Suarez, M.J. A solar radiation parameterization for atmospheric studies. NASA Tech. Memo. 1999, 15, 40. [Google Scholar]
Lang, S.E.; Tao, W.K.; Chern, J.D.; Wu, D.; Li, X. Benefits of a fourth ice class in the simulated radar reflectivities of convective systems using a bulk microphysics scheme. J. Atmos. Sci. 2014, 71, 3583–3612. [Google Scholar] [CrossRef]
UK Civil Aviation Authority: CAP1236: Guidance Regarding Flight Operations in the Vicinity of Volcanic Ash, 33 pp., Available online:. Available online: https://publicapps.caa.co.uk/docs/33/CAP%201236%20FEB17.pdf (accessed on 21 May 2023).
Climate Data Store. Available online: https://cds.climate.copernicus.eu/#!/home (accessed on 4 September 2022).
Plu, M.; Scherllin-Pirscher, B.; Arnold Arias, D.; Baro, R.; Bigeard, G.; Bugliaro, L. et al. An ensemble of state-of-the-art ash dispersion models: towards probabilistic forecasts to increase the resilience of air traffic against volcanic eruptions. Nat. Hazards Earth Syst. Sci. 2021, 21(10), 2973-2992. [CrossRef]

Figure 3. The injection height from Mastin et al. [34]’s formulation for plume height (Eq. 1); violet points: 2013 sequence; green points: 2015 sequence.

Figure 4. WRF Geopotential height [m] and wind barbs at 500 hPa for the 2013 paroxysms, namely (a) NSE1, (b) NSE2, (c) NSE3, (d) NSE4, (e) NSE5, and (f) NSE6.

Figure 5. WRF Geopotential height [m] and wind vectors at 500 hPa for the 2015 paroxysms, namely (a) VOR1, (b) VOR2, (c) VOR3, and (d) VOR4.

Figure 6. Average distribution of the total columnar

\bar{p_{10}}

ash concentration related to the 2013 paroxysms. Data without ash in the domain are filtered out. Units are μg m-2 and reported in log10 scale. The dotted line denotes the cross-section transect.

Figure 6. Average distribution of the total columnar

\bar{p_{10}}

ash concentration related to the 2013 paroxysms. Data without ash in the domain are filtered out. Units are μg m-2 and reported in log10 scale. The dotted line denotes the cross-section transect.

Figure 7. Average distribution of the total columnar

\bar{p_{10}}

ash concentration related to the 2015 paroxysms. Units are μg m-2 and reported in log10 scale. The dotted line denotes the cross-section transect.

Figure 7. Average distribution of the total columnar

\bar{p_{10}}

ash concentration related to the 2015 paroxysms. Units are μg m-2 and reported in log10 scale. The dotted line denotes the cross-section transect.

Figure 8. Top (a-b-c) Time-averaged PM10 concentration maps for the Oct-Dec 2013 sequence at different flight levels, namely (a) FL180, (b) FL240, and (c) FL300. Bottom (d-e-f): time-averaged PM10 concentration maps for the December 2015 sequence at the same flight levels as in a,b,c. Units are μg m-3 and reported in log10 scale. The upper scale limit of 4000 μg m-3 (3.6 in log10 scale) corresponds to prohibited airspace.

Figure 9. Daily average p10D ash concentrations at different flight levels for November 23, 2013 (first row); November 28, 2013 (second row). Panels (a, d) for FL180; panels (b, e) for FL240 and panels (c, f,) for FL300. Units are μg m^-3, label bar is reported in log10 scale with intervals matching the hazard levels for aviation.

Figure 10. Daily average p10D concentrations at different flight levels for December 3, 2015 (first row); December 5 (second row). Panels (a, d) for FL180; panels (b, e) for FL240 and panels (c, f) for FL300. Units are μg m^-3, label bar is reported in log10 scale with intervals matching the hazard levels for aviation.

Figure 11. Vertical distribution of PM10 hourly concentrations (in μg m^-3) as a function of air pressure (in hPa, cf. corresponding heights in Table 6) in altitude (i.e., from 0 to 11.8 km a.s.l.) and position along the SW-NE transects shown in Figures 6-7. (a) 27 Oct. 2013 at 20:00 UTC, (b) 31 Oct 2013 at 03:00 UTC, (c) 3 Dec. 2015 at 11:00 UTC, (d) 4 Dec 2015 at 15:00 UTC, (e) 4 Dec. 2015 at 21:00 UTC, (f) 5 Dec. 2015 at 15:00 UTC. The 4000 μg m^-3 concentration corresponds to prohibited airspace. Horizontal scales are 750 km (a-b) and 900 km (c-f). The topography of Etna and Aspromonte in Calabria appears in white.

Table 1. Onset/end of each paroxysm of the 2013 sequence and Mass Eruption Rates (MER) obtained by the VOLDORAD-2B radar system. Volcanic Explosivity Index (VEI) from MSVOLSO2L4 database (nd = no data or undetermined). Number of Volcanic Ash Advisories (VAA) released from Toulouse VAAC-MeteoFrance (https://vaac.meteo.fr/, accessed on 16 May 2023).

	2013	Julian day	Start (UTC)	End (UTC)	Duration	Bins	V2B MER (kg s^-1)	VEI	VAA
Paroxysm	Date	Start	HH:MM	HH:MM	HH:MM	Range in m (bin number)	Average (max)
NSE1	Oct. 26	299	01:35	10:27	08:52	3135–3285 (3; 4)	2.9 10⁴ (2.1 10⁵)	2	7
NSE2	Nov. 11	315	00:01	11:52	11:51	3135–3285 (3; 4)	3.6 10⁴ (2.3 10⁵)	2	2
NSE3	Nov. 16-17	320-321	22:14 (16/11)	04:35 (17/11)	05:21	3135–3285 (3; 4)	4.3 10⁴ (2.8 10⁵)	nd	4
NSE4	Nov. 23	327	07:13	10:26	03:13	3135–3285 (3; 4)	3.9 10⁵ (8.9 10⁶)	2	4
NSE5	Nov. 28	332	15:15	23:35	08:20	3135–3285 (3; 4)	2.5 10⁵ (1.3 10⁶)	nd	3
NSE6	Dec. 2	336	19:08	22:42	03:34	3135–3285 (3; 4)	2.5 10⁵ (1.8 10⁶)	2	2

Table 2. Onset/end of each paroxysm of the 2015 sequence and Mass Eruption Rates (MER) obtained by the V2B radar system. Volcanic Explosivity Index (VEI) from MSVOLSO2L4 database (nd = no data or undetermined). Number of Volcanic Ash Advisories (VAA) released from Toulouse VAAC-MeteoFrance (https://vaac.meteo.fr/, accessed on 16 May 2023), totalling 29 for 3–9 December.

	2015	Julian day	Start (UTC)	End (UTC)	Duration	Bins	V2B MER (kg s^-1)	VEI	VAA
Paroxysm	Date	start	HH:MM	HH:MM	HH:MM	Range in m (bin number)	Aver (max)
VOR1	Dec. 3	336	02:00	03:31	01:31	3885–4185 (8; 9; 10)	1.7 10⁶ (4.5 10⁶)	2	5
VOR2	Dec. 4	337	09:03	10:14	01:11	3735–4035 (7; 8; 9)	1.7 10⁵ (5.0 10⁵)	nd	3
VOR3	Dec. 4	337	20:26	21:15	00:49	3735–4035 (7; 8; 9)	1.8 10⁵ (1.2 10⁶)	nd	2
VOR4	Dec. 5	338	14:45	16:10	01:25	3735–4035 (7; 8; 9)	1.0 10⁵ (4.1 10⁵)	nd	4

Table 3. Physical and chemical options of the WRF-Chem model, in parenthesis the namelist value.

Microphysics	(7) – 4Ice Goddard Scheme
LW/SW radiation	(5,5) – New Goddard Shortwave and Longwave Schemes
Surface Layer	(2) – Eta Similarity Scheme
PBL	(2) – Mellor–Yamada–Janjic Scheme
Land surface	(4) – Noah–MP Land Surface Model
Initial and Boundary Conditions	FNL-GFS
chem_opt	(300) – GOCART aerosol model
TGSD	E1 distribution
ESP	MER from V2B radar data Injection heights from Mastin et al. [34]

Table 4. Total grain size distribution (E1) from Poret et al. [22] and Rizza et al. [18].

vash_#	1	2	3	4	5	6	7	8	9	10
µm	1000 2000	500 1000	250 500	125 250	62.5 125	31.25 62.5	15.625 31.25	7.8125 15.625	3.90625 7.8125	0 3.90625
wt%	16.7	8.3	10.4	12.5	6.4	12.5	14.6	8.3	6.2	4.2

Table 5. Times of the meteorological snapshots for the synoptic analysis.

Sequence	Paroxysm	Date	UTC time HH:MM
2013	NSE1	October 26	12:00
2013	NSE2	November 11	12:00
2013	NSE3	November 17	00:00
2013	NSE4	November 23	12:00
2013	NSE5	November 28	12:00
2013	NSE6	December 3	00:00
2015	VOR1	December 3	00:00
2015	VOR2	December 4	12:00
2015	VOR3	December 5	00:00
2015	VOR4	December 5	12:00

Table 6. Flight levels in the airspace, their corresponding pressure level, and the associated height above sea level (a.s.l.).

Pressure level (hPa)	Flight level (FL)	Height (feet, a.s.l.)	Height (m, a.s.l.)
300	FL300	30000	9200
400	FL240	24000	7300
500	FL180	18000	5500

Table 7. Shading intervals reported in Figure 9 and Figure 10 and their connection with the hazard levels defined by UK Civil Aviation Authority [43] and used by ICAO.

p10D interval (μg m^-3)	log10 scale	Shading	Hazard level
0–200	0–2.2	blue	none
200–2000	2.2–3.3	green	low
2000–4000	3.3–3.6	red	medium
> 4000	> 3.6	dark red	high

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.

© 2023 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).

Copyright: This open access article is published under a Creative Commons CC BY 4.0 license, which permit the free download, distribution, and reuse, provided that the author and preprint are cited in any reuse.

MDPI Initiatives

Important Links

Choose an area of interest and we will send you notifications of new preprints at your preferred frequency.

Disclaimer