Investigation of an Intense Dust Outbreak in the Mediterranean Using XMed-Dry Network, Multiplatform Observations, and Numerical Modeling

Abstract

The Weather Research and Forecasting (WRF) model with online coupled chemistry (WRF-Chem) is applied to study an intense Saharan dust outbreak event affecting the Italian peninsula in 15 and 16 April 2018. According to the MODIS retrievals, this intrusion was characterized by an intense aerosol optical depth (AOD) peak value in the southern Mediterranean. Measurements within the Dry Deposition Network Across the Mediterranean (XMed-Dry) are compared with the output of the WRF-Chem model. XMed-Dry samples from Lecce (Italy), Athens (Greece) and San Lawrenz/Gozo (Malta) were analysed with respect to aerosol particle size distribution, relative dust contribution, and composition. The discrepancy between the model and measured deposition indicate the need to formulate in WRF-Chem more sophisticated deposition schemes, this will need to evaluate the sensitivity of the results to the precise particle size limits chosen for the aerosol model. Moreover, satellite retrievals from MODIS sensors elaborated with the MAIAC algorithm, Aeronet stations, and measurements of ${\mathrm{PM}}_{10}$ at the selected sites were also considered. In a numerical domain that spans the Mediterranean and the northern Saharan desert, two different dust emission schemes, namely Gocart-AFWA and the Shao-2001, were tested and compared with multiplatform observations for simulation period covering the dust outbreak. Actual results indicate that both emission schemes would benefit from replacing the static erodibility map and soil particle distribution with remote sensed and in-situ observational data.

1. Introduction

The Intergovernmental Panel on Climate Change [1] has recognized mineral dust as the main component of atmospheric aerosols, and considered it as an ‘essential climate variable’. It may be defined as one of the largest natural aerosol sources and influencing significantly the weather and climate system [2]. Aerosol particles influence climate by modifying directly the global energy budget by scattering and absorbing long- and short-wave radiation [3] and interacting with clouds modifying the development and occurrence of precipitation (indirect effects).

The Sahara Desert contributes to the total emission of dust to more than half of the global [4]. Strong convective structures generated by the intense surface heating raise the dust particles for thousands of meters into the troposphere, where they are finally transported at the continental and intercontinental scales [5]. It has been estimated that 10% of the Saharan dust amount is transported across the Mediterranean Sea to Europe [6,7] in episodic storms or following periodic patterns [8,9]. In recent studies [10,11], it was reported an average contribution between 30% and 40% to the optical depth over different cities in southern Europe.

The most important parameter describing the level of columnar aerosols in the atmosphere is the Aerosol Optical Depth (AOD). AOD may be measured at ground level with a sun-photometer or from space by Spectroradiometer (e.g., Moderate Resolution Imaging Spectroradiometer, MODIS). At this purpose, a federation of ground-based sun-photometers has been established by NASA and the French CNRS. It includes nearly 1000 sun photometers distributed worldwide within the Aerosol Robotic Network (AERONET) [12] whose data are processed according to the same aerosol retrieval algorithms [13,14]. On the other side, AOD measurements from space started in 1999 by NASA with the launch of the Terra spacecraft that is carrying the MODIS sensor. The characterization of atmospheric aerosol is a core of the MODIS mission [15] and the AOD still represents the most robust aerosol physical parameter derived from space.

To complement measurements, numerical modeling may be considered a useful approach to study the dust transport from the Saharan arid and semi-arid regions toward the Mediterranean area, allowing a comprehensive understanding of the dust intrusions and their effects on air quality [16] and human health [17].

The Weather Research and Forecasting model coupled with chemistry (WRF-Chem; [18]) has been previously used by many research groups to study this argument, among the others we may cite Bossioli et al. [19] that studied the influence of biomass burning during summertime over the eastern Mediterranean, Georgiou et al. [20] who investigated an intense photochemical activity and their impact on air quality over the eastern Mediterranean, Chaibou et al. [21] who evaluated dust extinction and vertical profiles with different dust emission schemes into the WRF-Chem (v4.02) model over North Africa, and Su and Fung [22] and LeGrand et al. [23] that evaluated the model performance in the simulation of dust concentrations over an eastern Asia domain using different dust emission formulations. Previous studies of the WRF-Chem reproduction of the AOD over the Mediterranean area have been performed by Flaunas et al. [24] that made a sensitivity analysis over the broader Mediterranean regions of the WRF-Chem (V3.6.1) model to different dust emission parametrization, by Rizza et al. [9,25], that evaluated the different dust emission schemes in the same version of the WRF-Chem model and by Tsarpalis et al. [26] that implemented a dust wet deposition scheme in a dust module of WRF-Chem (V3.8).

Experimental data of mineral dust deposition are useful to improve the understanding of deposition processes itself but also to support the development and parameterizations of more sophisticated, size dependent deposition schemes in the WRF-Chem model. In this context, a dry deposition measurement network across the Mediterranean (XMed-Dry) has been recently established [27]. In this network, a set of seven newly developed dry deposition-only collectors has been installed at different locations across the Mediterranean (Huelva, Barcelona, Spain; IÎe-Rousse, France; Gozo Island, Malta; Lecce, Italy; Athens, Greece; Nicosia, Cyprus) to capture spatial and temporal variability of atmospheric aerosols including their size and composition [28]. In the context of this study, deposition data from three XMed-Dry stations—namely: Lecce, Athens, and Malta have been utilized to validate the WRF-Chem model. It is utilized to analyze an intense Saharan dust outbreak event affecting the Italian peninsula in 15 and 16 April 2018. This intrusion was characterized by an intense aerosol optical depth (AOD) peak value in the southern Mediterranean, according to the MODIS retrievals.

The WRF-Chem model, that actually implements two different dust emission schemes, namely the GOCART-Air Force Weather Agency (AFWA; [29]) and the University of Cologne (UoC; [30]), will be tested and compared with multiplatform observations, for a simulation period covering the dust outbreak in April, 2018. Generally, analysis shows that the AFWA emission scheme significantly over-predicts the emission of mineral dust, while the UoC parameterization gives overall better results in terms of optical properties and total mass. It is important to remark that rather than discovering new dust-emission parameterizations we are mainly interested in finding the optimal WRF-Chem configuration in term of physics and chemistry setup. We partially did this in a previous work [25] were we tested the UoC dust emission scheme using two different WRF land surface schemes and found very different results only by changing the surface land model.

The paper outline is the following: Section 2 describes materials and methods utilized for this study; Section 3 depicts results and the relative analysis; and finally, in Section 4, the conclusions of the study are drawn.

2. Materials and Methods

2.1. Description of the XMed-Dry Network and the Multiplatform Observations

2.1.1. Dry Deposition Measurements

Aerosol dry deposition was collected at three locations for this study: on the flat roof of the Engineering department building (Ecotekne campus) at Lecce (18.6° E, 40.2° N, 40 m above sea level/15 m above street level), on top of a measurement container in a pine forest at Athens (23.5° E, 37.6° N, 290 m above sea level/10 m above ground level), and on flat roof of a building at the edge of a housing area at San Lawrenz, Gozo (14.1° E, 36.0° N, 130 m above sea level/10 m above street level). A novel dry deposition collector with an active rain protection was used. The geometry of its sampling head follows the ‘flat-plate’ design as described by Waza et al. [31]. The function of the spacers in the given design was replaced by an external mounting frame, thus reducing the disturbance of the flow field between the plates. For additional rain protection, the outer part of the top plate is moved down around the lower plate in case of a positive signal from the rain detector sensor. The sampling substrate is positioned approximately 2 m above the roof. Samples were exposed to the atmosphere for 48 to 72 h, approximately. Particles were collected on a pure carbon adhesive (SpectroTabs, Plano GmbH, Wetzlar, Germany) mounted to a standard electron microscopy aluminum stub and stored in dry boxes at room temperature.

Particles were analyzed by automated scanning electron microscopy (FEI Quanta 400F ESEM, FEI, Eindhoven, The Netherlands) with attached X-ray fluorescence analysis (Oxford X-Max 120, Oxford Instruments, Abingdon, UK) under vacuum condition. Image acquisition was done at 100–150 nm/pixel resolution, and particle size was estimated from their projected cross-sectional area. Particle mass and aerodynamic diameter were estimated from their size assuming preferential flat orientation and a density based on the approximate composition from chemistry (for details refer to Waza et al. [31]). X-Ray fluorescence analysis time was chosen to allow for a minimum of 40,000 X-ray counts per particle. Quantification was carried out by the corresponding control software (Oxford Aztec 4.0). Prior to automated analysis, samples were visually inspected to avoid analysis of damaged samples surface areas. A total of 1500 to 5000 particles per sample were analyzed. More details of the data acquisition, quantification and processing have been previously reported [27].

Particles were classified according to their chemical composition [32]. This detailed classification was then summarized into dust, sea-salt and other particles groups. The dust particle group mainly consists of silicate and carbonate particles, but includes also iron and titanium oxides as well as other minerals, which only contribute negligibly to the total mass. The statistical uncertainty dominating the overall random error was estimated by a bootstrap approach and is reported in terms of the limits of the central 95% confidence interval.

2.1.2. PM₁₀ Data

In this study daily mean ${\mathrm{PM}}_{10}$ (particulate matter with aerodynamic diameter less than 10 µm) concentration values have been obtained by regional or municipal air quality monitoring networks in Lecce, Athens, and Malta by automatic devices. In particular, daily mean ${\mathrm{PM}}_{10}$ values in Lecce have been obtained by the air quality monitoring network of the environmental protection regional agency (ARPA-Puglia) at the Lecce-Garigliano traffic urban monitoring station (18.1° E, 40.3° N). The second ${\mathrm{PM}}_{10}$ measurement point was localized at a rural background site in Għarb, Gozo (14.2° E, 36.1° N). The third measurement point is located at the Demokritos suburban monitoring station in Agia Paraskevi, Attiki, Greece (23.8° E, 38.0° N).

2.1.3. AERONET AOD

AERONET [33] is a global network of nearly 1000 sun-photometers whose aerosol data are processed according to the same retrieval algorithm and made available online at the following web-portal [12]. For this study, we have utilized Version 3 AOD measurements at level-2 (cloud screened and quality checked) at the stations of Athens (Greece; 23.7° E, 38.0° N, 130 m elevation), Rome Tor Vergata (Italy; 12.6° E, 41.8° N; 130 m elevation) and Gozo (Malta; 14.3° E, 36.0° N, 111 m elevation). The station of Lecce (University) does not have valid level-2 data for the period considered and has been replaced by the above-mentioned Rome Tor Vergata site.

2.1.4. MODIS/MAIAC AOD

The MODIS Multi-Angle Implementation of Atmospheric Correction (MAIAC) is a recent advanced algorithm which uses time series analysis and a combination of pixel- and image-based processing to improve accuracy of cloud detection, aerosol retrievals and atmospheric correction [34]. MAIAC provides a suite of atmospheric and surface products in Hierarchical Data Format for the NASA Earth Observing System (HDF4_EOS), including daily atmospheric properties (MCD19A2). For each Terra/Aqua orbit, the MAIAC daily MCD19A2 file includes, among other things, the: (i) over land AOD and type; and (ii) the over water AOD outside of glint area and Fine Mode Fraction (FMF). Data access is granted through the NASA/USGS Earthdata search portal [35] from which AOD-HDF4 data may be downloaded for the selected domain and time interval of interest. The current MODIS Collection 6 MAIAC algorithm has changed considerably since its first release [36] to adapt to global processing and improve cloud/snow detection, aerosol retrievals and atmospheric correction of MODIS data [37].

2.2. Model Setup and Description

For this study, the WRF-Chem model version 3.6.1 has been utilized in the numerical domain depicted in Figure 1. It is important to remark that all the code alterations reported by Ukhov et al. [38] and LeGrand et al. [23], have been manually incorporated in this version. The numerical domain covers the northern regions of Sahara Desert and part of Europe, with 300 × 300 points with a horizontal grid spacing of 10 km. The simulation started on 10 April 2018 (0000 UTC) and ended on 17 April 2018 (0000 UTC). Boundary and initial conditions are at 1-degree resolution and provided from NCAR/NCEP Final Analysis from Global Forecast System (FNL from GFS).

Based on the WRF setup described in the work by Rizza et al. [39] for the Northern Sahara and the Mediterranean basin, the following physics parameterizations are utilized: the Mellor–Yamada–Nakanishi and Niino (MYNN) 2.5 level turbulent kinetic energy (TKE) parameterization is used to describe the planetary boundary layer and surface layer scheme [40]. The Rapid Update Cycle (RUC) Land Surface Model is chosen to represent the land surface processes [41]. The radiation physics is parameterized using the Rapid Radiative Transfer Model [42] for both short-wave (ra_sw_physics = 4) and long-wave (ra_lw_physics = 4). The Improved Bulk Microphysics Scheme of Thompson et al. [43] is used for the treatment of the microphysics processes (

Table 1. Physics options for the Weather Research and Forecasting model.

Parameterization	Namelist Variable	(Option) Scheme
Surface layer	sf_sfclay_physics	(5) MYNN
Surface model	sf_surface_physics	(3) RUC Land Surface
Boundary layer	bl_pbl_physics	(5) MYNN Level 2.5
Microphysics	mp_physics	(8) Thompson
Longwave radiation	ra_lw_physics	(4) RRTMG Longwave
Shortwave radiation	ra_sw_physics	(4) RRTMG Shortwave
Surface Layer	sf_sfclay_physics	(5) MYNN

).

The WRF-Chem setting related to the aerosol model concerns the hybrid bulk/sectional Georgia Tech/Goddard Global Ozone Chemistry Aerosol Radiation and Transport (GOCART, [44]) model. This is implemented in WRF-Chem model selecting the namelist variable chem_opt = 300 (

Table 2. Weather Research and Forecasting model, coupled with chemistry (WRF-Chem) chemistry setup.

Case	Chem_Opt	DUST_Opt	Dust_Schme	Scheme
run_2	300	4	1	UoC-S01
run_3	300	3	-	AFWA

). As shown in

Table 3. Effective diameter and size bin range for the Air Force Weather Agency AFWA scheme.

Dust-Bin AFWA	Lower-Upper Diameter (µm)	Effective Diameter (µm)
1	0.2–2	1.1
2	2–3.6	2.8
3	3.6–6	4.8
4	6–12	9
5	12–20	16

and

Table 4. Effective diameter and size bin range for the University of Cologne UoC-S01 scheme.

Dust-bin UoC-S01	Lower-Upper Diameter (µm)	Effective Diameter (µm)
1	0–2.5	1.25
2	2.5–5	3.75
3	5–10	7.5
4	10–20	15

it considers five and four dust bins with diameters in the range (0–20) µm for the two dust emission schemes described below.

Table 3 and Table 4 show the effective diameter that is used in the specialized routine that calculates the gravitational settling deposition. In the WRF-Chem releases above 4.0 (June 2018) the number and diameter of dust bins are consistent between the two emission schemes.

The physics of the emission of dust is described with great details in a recent review by LeGrand et al. [23], in which the authors have analyzed the principal mechanisms responsible for the process of soil particles mobilization. In particular, three basic processes are responsible for the entrainment of dust particles in the atmospheric flows: (1) aerodynamic lift, (2) saltation bombardment, and (3) particle disaggregation. These three processes generate the lifting of a wide size spectrum of soil particles, from sand sized grains at around 70 µm by aerodynamic wind shear, down to smaller dust-sized particles at 0.1 µm (saltation bombardment and disaggregation).

There are actually three different dust emission schemes implemented into the WRF-Chem 3.6.1 model using the GOCART aerosol chemistry core (chem_opt = 300): the GOCART-WRF scheme (dust_opt = 1), the AFWA scheme (dust_opt = 3), and the University of Cologne (UoC-S01) scheme (dust_opt = 4).

Both AFWA and UoC-S01 emission schemes were introduced in WRF-Chem, because the GOCART-WRF original model revealed some issues with the dust emission scheme [23], making this parameterization [45] not appropriate under certain environmental conditions [23].

Specific details of the AFWA dust emission scheme, which may be considered a modification of the Marticorena and Bergametti [46] emission model, is already described in a work by LeGrand et al. [23] so it will not be replicated here. It is important to remark that in this scheme the dust emission is controlled by the saltation of larger particles that are triggered by wind shear at the surface leading to the emission of smaller particles by saltation bombardment. The vertical bulk dust flux is calculated as [23,47]

$$\{\begin{array}{c}\begin{array}{cc}{F}_{bulk}=G\times ERO{D}^{\gamma }\times \beta & z_0\le 0.2 \mathrm{m}\end{array}\\ \begin{array}{cc}{F}_{bulk}=0& z_0>0.2 \mathrm{m}\end{array}\end{array}$$

(1)

where G is the total streamwise horizontal saltation flux ([23], eq.13), $\beta ={10}^{0.134\left(\%clay-6\right)}$ is the sandblasting mass efficiency calculated considering only the soil clay fraction; γ is the exponential tuning constant for erodibility. The fraction of erodible grid cell is represented by the EROD parameter which is provided by the WRF-Chem preprocessor (WRF Preprocessing System). Once the total bulk emission ($F_{bulk}$) is established, size-resolved dust emission fluxes (g cm⁻² s⁻¹) are obtained according to a size distribution that is hardcoded in the AFWA dust emission subroutine.

The University of Cologne (UoC-S01) scheme is activated in WRF-Chem by setting dust_opt = 4 and dust_schme = 1, in the chem-section of the WRF input namelist. It was firstly introduced by Shao [30] and further documented in Kang et al. [48], Su and Fung, [22], Rizza et al. [9,25], LeGrand et al. [23], and Zeng et al. [49].

The UoC-S01 scheme may be considered a spectral emission scheme because it is based on a size-resolved dust emission equation by supposing that particles are divided into a finite number of particle size intervals. Then the total dust flux (F) is expressed as an integral of the dust emission rate

$$F={{\displaystyle \int }}_{d_i}^{d_s}{{\displaystyle \int }}_0^{d_i}f\left(d_i,d_s\right)p\left(d_i\right)p\left(d_s\right)\delta d_i\delta d_s$$

(2)

where $f\left(d_i,d_s\right)$ is the emission rate for particles of size $d_i$ generated by the saltation of particles of size $d_s$ [9,30]. Full details of the UoC-S01 formulation may be found in Shao [30], Kang et al., [48] and LeGrand et al. [23].

In this context, UoC-S01 and AFWA schemes in WRF-Chem follow the same general algorithm. The main difference between the two schemes is that the more elaborated UoC-S01 algorithm uses size-resolved saltation particle bins already at the level of the saltation process. Another difference concerns the calculation of the threshold friction velocity, the UoC-S01 scheme uses an additional correction factor caused by surface roughness elements, and commonly referred to as the drag partition correction [23]. Concerning the calculation of the threshold friction velocity, the differences between the two schemes is showed in Figure A1 of Appendix A. It is evident the largest positive difference of the quantity ($u_{\ast }-u_{\ast t})$ for the AFWA scheme (panels b, d) compared to UoC (panels a, c). This quantity is of primarily importance because it defines the emission of dust at each appropriate grid point.

The removal processes consider both dry and wet deposition. A dry deposition scheme [50] accounting for gravitational settling and surface deposition is used to simulate the dry removal of mineral dust. A wet deposition scheme, which considers the total precipitation from both the microphysical and convective parameterizations, is used for natural (sea-spray and mineral dust) aerosols [22]. In particular, the gravitational settling routines actually implemented in WRF-Chem utilizes the Stokes law corrected by the Cunningham slip factor [51] for calculating the terminal fall velocity for each effective diameter (Table 3 and Table 4).

3. Results

3.1. Events Description by ERA5 Reanalysis and MODIS Retrievals

This event was chosen by making a preliminary cross-checking in which we have considered the following criteria: (i) plume intrusion in Mediterranean basin with the highest AOD values; and, (ii) availability of Xmed-dry data. Considering that XMed-dry campaign started in 2017 with a series of Intensive Observation Periods, the period (14–16 April 2018) we have investigated was the only one satisfying the two criteria above.

The synoptic analysis of the dust event is described using the ERA5 reanalysis of the geopotential at 500 hPa. The geopotential is obtained downloading data from the Copernicus Climate data store [52]. In particular, the selected dataset concerns the ERA5 hourly data on pressure levels (1979 to present) from which the geopotential at 500 hPa may be downloaded in GRIdded Binary (GRIB) format. ERA5 represents the fifth generation of European Centre for Medium-Range Weather Forecasts (ECMWF) atmospheric reanalysis of the global climate. It combines model data with observations from across the world into a globally complete and consistent dataset [53]. ERA5 [54] is produced using 4D-Var data assimilation in the ECMWF’s Integrated Forecast System (IFS), with 137 sigma/pressure levels in the vertical and 31 km horizontal resolution.

Figure 2a–f shows the comparison between the geopotential height at an atmospheric pressure of 500 hPa predicted by the WRF-Chem (Figure 2b,d,f) and GRIB data downloaded from ECMWF/ERA5-reanalysis (Figure 2a,c,e). Data are displayed for 14 April at 06:00 UTC (upper panel), 15 April 1at 06:00 UTC (central panel) and 16 April at 18:00 UTC (bottom panel). On 14 April, ERA5-reanalysis (Figure 2a) shows a pressure minimum localized in the internal regions of Algeria fostering a southern circulation in the central Mediterranean basin. The following day (15 April), ERA5-reanalysis depicts that the pressure minimum is dislocated toward Tunisia reinforcing the circulation toward Italy. On 16 April at 18:00 UTC, the weakening of the pressure minimum is evident with the consequent zonal circulation (west–east) on the central and eastern Mediterranean basin.

It is evident the close reproduction of the ERA5 (Figure 2a,c,e) geopotential at 500 hPa by WRF-Chem model (Figure 2b,d,f). This is quite important as the long-range transport of dust plumes is governed by these high-tropospheric winds. These synoptic conditions are extremely favorable for the generation and transport of Saharan dust toward the Mediterranean basin for two reasons: (i) the presence of a low-pressure system and relative high winds speed on the Northern Sahara Desert; and the consequent (ii) southern wind circulation on the Mediterranean basin.

The transport of Saharan dust may be further verified by analyzing Figure 3 that represents the time averaged (14–16 April) map of MODIS/AOD at 550 nm. It is calculated combining the dark target (DT; [15]) and deep blue (DB; [55]) algorithms for land and ocean. The final product depicted at Figure 3 is obtained by merging MODIS-Terra (MOD08_D3 v6) and MODIS-Aqua (MYD08_D3 v6) intermediate maps. The extent of the dust plume represented by the time average map of Figure 3 reports the strong intrusion of Saharan dust (AOD peak ~2.0), from the emission regions encompassing the Sahara Desert portion of Algeria and Tunisia [4] toward southern Italy (15 April) and the eastern Mediterranean regions (16 April).

3.2. Comparison with XMed-Dry Measurements

In

Table 5. Start/end of sampling period for each respective sampling site: namely Lecce (sample code: LES), San Lawrenz-Gozo (sample code: SLS), and Athens (sample code: ATS).

Sampling Site	Start Sampling YYYY-MM-DD UTC	End Sampling YYYY-MM-DD UTC
LES_125	2018-04-11 10:10	2018-04-13 12:10
LES_126	2018-04-13 12:25	2018-04-11 11:35
SLS_148	2018-04-11 18:20	2018-04-13 17:20
SLS_149	2018-04-13 17:30	2018-04-15 10:15
ATS_137	2018-04-11 13:30	2018-04-13 12:20
ATS_138	2018-04-13 12:25	2018-04-16 12:25

the sampling periods are shown in the sampling site of Lecce (sample code: LES), San Lawrenz-Gozo (sample code: SLS) and Athens (sample code: ATS). The sampling periods, from 11 April to 16 April, covered the dust intrusion as can be verified by the presence of dust in the samples.

In

Table 6. Simulated (run_2, run_3) and measured dry deposition rates at the three sampling locations. For the measured dust deposition rates in parentheses the lower and upper limit of the central 95% confidence interval are shown. In the last three columns total deposition rate, non-dust fraction and the average PM₁₀ concentration for the relative period are reported.

Sampling Site	Run_2 (mg m−2day−1)	Run_3 (mg m−2day−1)	Measured Dust Deposition (0.625–10 µm) (mg m−2day−1)	Measured Total Deposition (0.625–10 µm) (mg m−2day−1)	Non-Dust Fraction (0.625–10 µm)	PM10 Average Concentration (µg m−3)
LES_125	3.2	10.8	3.32 (2.75–4.09)	3.5	0.05	33.8
LES_126	4.2	11.3	7.70 (6.82–8.73)	8.1	0.05	74.2
SLS_148	1.7	2.9	10.51 (9.05–12.23)	11.6	0.09	50.5
SLS_149	4.1	16.7	8.83 (7.67–10.23)	17.3	0.49	83.4
ATS_137	1.7	6.7	2.29 (1.83–2.89)	2.3	0.01	31.75
ATS_138	5.5	14.5	2.58 (2.19–3.06)	2.6	0.01	37.95
Average	2.5	10.5	5.9

, for each XMed-Dry investigated sample, it is reported the dust deposition rates by classification with the respective confidence intervals, total deposition rates and non-dust fraction. Dust and total deposition rates as well as the non-dust fraction have been obtained for particles which estimated volume-equivalent diameters range between 0.625 and 10 µm. In the last column of the Table 6 we have reported the average ${\mathrm{PM}}_{10}$ concentration (obtained by daily values shown in

Table 7. Daily mean PM₁₀ at Authority’s station in Gharb (Malta) and predicted with the WRF-Chem configurations run_2 and run_3.

	Malta	Run_2	Run_3
April	${\mathsf{\mu }\mathrm{g} \mathrm{m}}^{-3}$	${\mathsf{\mu }\mathrm{g} \mathrm{m}}^{-3}$	${\mathsf{\mu }\mathrm{g} \mathrm{m}}^{-3}$
12	56.7	73.5	57.5
13	43.4	49.8	161.9
14	48.6	77.0	221.8
15	158.2	115.8	944.3
16	40.4	14.3	29.7

,

Table 8. Daily mean PM₁₀ at Lecce and predicted with the WRF-Chem configurations run_2 and run_3.

	Lecce	Run_2	Run_3
April	${\mathsf{\mu }\mathrm{g} \mathrm{m}}^{-3}$	${\mathsf{\mu }\mathrm{g} \mathrm{m}}^{-3}$	${\mathsf{\mu }\mathrm{g} \mathrm{m}}^{-3}$
12	32.2	41.0	124.8
13	35.5	20.3	25.1
14	35.0	17.2	36.3
15	42.3	23.3	51.5
16	184.0	50.2	214.5

and

Table 9. Daily mean PM₁₀ at Athens and predicted with the WRF-Chem configurations run_2 and run_3.

	Athens	Run_2	Run_3
April	${\mathsf{\mu }\mathrm{g} \mathrm{m}}^{-3}$	${\mathsf{\mu }\mathrm{g} \mathrm{m}}^{-3}$	${\mathsf{\mu }\mathrm{g} \mathrm{m}}^{-3}$
12	27.2	18.6	86.5
13	36.2	29.4	86.2
14	45.2	51.2	129.0
15	35.1	39.7	106.3
16	35.2	181.0	641.8

) for the investigated periods. Details of measurements are already described in paragraph 2.1.1. The deposition rate of the WRF-Chem model is calculated by considering the gravitational settling (GRASET_# variable) and the dust dry deposition for each size bins (DRYDEP_# variable) and summing all bins as depicted in the following Equation (3):

$$\mathrm{dep}\_\mathrm{rate}={{\displaystyle \sum }}_{\mathrm{i}}({\mathrm{GRASET}}_{\mathrm{i}}+{\mathrm{DRYDEP}}_{\mathrm{i}})$$

(3)

This permits a comprehensive comparison between XMed-dry measurements and model outputs. The model results (dep_rate) for the run_2 and run_3 configurations are shown in the second and third columns of Table 6. For each XMed-dry sampling site, the variables of Equation (3) are extracted from the WRF-Chem netCDF output file in correspondence of the appropriate (lat, lon) coordinates of each station. To have a proper comparison between experimental data and model, dep_rate values are normalized using the same units (mg m⁻² day⁻¹).

If we consider the mean values that are calculated by averaging over all periods and stations and reported in the last row of Table 6, it results a positive bias (under-estimation) equal to 2.5 mg m⁻² day⁻¹ for run_2 and a negative bias of −4.6 mg m⁻² day⁻¹ corresponding to a considerable over-estimation for run_3.

The non-dust fraction (sixth column of Table 6) represents the abundance of particle mass not classified as dust. For the presented days, this non-dust fraction was mainly composed of sea-salt, but contained also sulfate particles (sodium/calcium/ammonium/potassium/magnesium sulfate and mixtures thereof) to a lesser extent. As can be seen, for all the samples the dust component is dominant except for the SLS_149 sample (Malta site). In fact, for all samples non-dust fraction ranges between 0.09 and 0.01, while for SLS_149 one non-dust fraction is 0.49, i.e., dust always dominates the total deposition.

Observing average ${\mathrm{PM}}_{10}$ concentration values, obtained by daily values of air quality monitoring networks in Lecce, Athens, and Malta, good site-specific relative agreement with XMed-Dry total deposition rates is found. However, it appears that the same ${\mathrm{PM}}_{10}$ concentration leads to different deposition rates at the different locations, which might indicate local conditions.

3.3. Comparison with AOD/MAIAC and AERONET

The comparison of satellite data with model results is performed by computing the columnar AOD from the extinction coefficient that is calculated by the model optical routines [56] at 0.55 μm ([23], eq.36). Figure 4 shows the binned average of the AOD at 0.55 μm elaborated for 14 April (a–c), 15 April (d–f) and 16 April (g–i). Left column (a,d,g) represents the Modis/MAIAC AOD, central column (b,e,h), the AOD from WRF-Chem/run_2 and finally the right column the AOD from WRF-Chem/run_3 configuration. Both Modis/MAIAC retrievals and model outputs show the generation/emission of mineral dust during 14 April (Figure 4a–c) to the dust source area corresponding to the Chott region in Northern Tunisia [4,9]. The spatial pattern is typical of a low-pressure system with the associated high wind speed at the surface (10 m) as described above. The following days (15 and 16 April) are characterized by the long-range transport toward Italy (Figure 4d–f) and toward Greece on 16 April (Figure 4g–i). More in detail, it may put in evidence: (i) the quite good spatial reproduction; and (ii) the overestimation of model AOD for both configurations, even if the run_2 produces a better comparison for all the investigated days. It should be also notified the consistent presence of clouds that prevented the aerosol retrievals from the MAIAC/AOD algorithm.

To complement the comparison with the satellite retrievals and to better follow the temporal evolution of the AOD field over the central Mediterranean, we used the AOD Level 2.0 (Version 3) AERONET measurements at Athens, Rome (Tor Vergata) and Malta. These data are automatically cloud cleared and quality assured with pre-field and post-field calibration applied.

Figure 5a–c depicts the hourly resolved AOD at 550 nm from the run_2 simulation (red line), run_3 (blue line) and the corresponding AERONET measurements (black crosses).

The measurements show AOD peaks (<1.0) on 16 April at Athens, on 15 and 16 April at Rome Tor Vergata (~1.0) and the highest peak (~1.5) on 15 April at Malta. These measurements agree with the synoptic transport of the dust plume that is described in Section 3.1. The comparison between AERONET measurements with the two WRF-Chem runs reveals clearly the optimal reproduction of this columnar optical property by the run_2 configuration, with a good performance in terms of peak value and timing. On the other side, the AFWA setup (run_3) produces a relevant overestimation of the AOD, with several unrealistic peak values registered during the simulation period.

3.4. Comparison with MERRA2 Reanalysis and Surface PM₁₀

The final set of comparison is realized considering the dust column mass density (DUCM) and in situ particulate mass concentration (${\mathrm{PM}}_{10}$) data. The variable DUCM represents the columnar content of dust aerosol particles in the atmosphere. This quantity is calculated in the WRF-Chem output system by multiplying for each grid point the dust mass mixing ratio (${\mathsf{\mu }\mathrm{g} \mathrm{kg}}^{-1})$ with the dry-air density (${\mathrm{kg} \mathrm{m}}^{-3}$) and summing over each vertical layer depth (m). The resulting quantity is expressed in ${\mathsf{\mu }\mathrm{g} \mathrm{m}}^{-2}$ and gives an indication of the amount of dust in the atmosphere.

In the panels of Figure 6 we have reported the comparison of the daily averaged DUCM from WRF-Chem run_2 and run_3 with the corresponding M2T1NXAER-DUCMASS variable of the MERRA2 Time-averaged, Single Level, Assimilation, Aerosol Diagnostics (V5.12.4) (GMAO, 2015). More in detail, in Figure 6 is reported the daily averaged maps for dust column mass density for 14 April (a–c), 15 April (d–f), and 16 April (g–i). Left column (a,d,g) represents run_2 outputs, central column (b,e,h) the run_3 outputs and right column (c,f,i) represents the MERRA2 reanalysis maps. There is in general an agreement for the spatial pattern of the daily averaged field between MERRA2 reanalysis and the WRF-Chem predictions. Model runs and reanalysis show that the dust plume emitted on 14 April is transported toward Southern/Central Italy on 15 April and then dislocated eastward in the following day. This spatial pattern is quite well reproduced by both schemes, we have to notify once more the considerable overprediction in term of intensity of the run_3 configuration and the good matching between run_2 and MERRA2, especially during 15 April.

If we consider particulate mass concentration (${\mathrm{PM}}_{10}$) data, in the second column of Table 7, Table 8 and Table 9, we have reported the daily mean ${\mathrm{PM}}_{10}$ concentrations measured at Gharb (Malta; Table 7), Lecce (Table 8), and Athens (Table 9) monitoring stations. In third and fourth columns are reported the corresponding ${\mathrm{PM}}_{10}$ daily average from the WRF-Chem run_2 and run_3, respectively. Being surface values, all the quantities are expressed in μg m⁻³. If we consider the dust-grain distribution showed in Table 3 and Table 4, the dust contribution to ${\mathrm{PM}}_{10}$ is calculated, in the WRF-Chem model, using the formulas [38]

$${\mathrm{PM}}_{10} = \mathrm{dust}\_1 + \mathrm{dust}\_2 + \mathrm{dust}\_3 \mathrm{run}\_2\text{-}\mathrm{scheme}$$

$${\mathrm{PM}}_{10} = \mathrm{dust}\_1 + \mathrm{dust}\_2 + \mathrm{dust}\_3 + 0.737 * \mathrm{dust}\_4 \mathrm{run}\_3\text{-}\mathrm{scheme}$$

Due to the proximity to North African coasts (Figure 1), Malta is impacted by the dust plume on 15 April with a quite high daily mean value (158 μg m⁻³) reported by the authority’s measurement station in Gharb (Malta; Table 7). The XMed-Dry sample SLS_149 (Table 6) revealed that, during the sampling period from 13 (17:30 UTC) to 15 (10:15 UTC) April 2018, the sea salt deposition rate (not showed) is nearly as high as the dust rate for the same size range of 0.625–10 μm. On the other side, on 15 April the run_2 output underestimated the daily mean ${\mathrm{PM}}_{10}$ concentration by about 40 ${\mathsf{\mu }\mathrm{g} \mathrm{m}}^{-3}$. This bias was caused by neglecting the contribution of sea spray to ${\mathrm{PM}}_{10}$ concentrations in the present study. In conformity with the above results, run_3 configuration provides a considerable over-prediction for this and the other ${\mathrm{PM}}_{10}$ stations described below.

Considering the ${\mathrm{PM}}_{10}$ concentration values of Lecce station (Table 8), it may be noted that the dust plume invested the city on 16 April, with a daily average value that is quite high (184 ${\mathsf{\mu }\mathrm{g} \mathrm{m}}^{-3});$ on the other side the XMed-Dry sample LES_126 revealed only 5% of non-dust composition, in the range of 0.625–10 μm. Model output from run_2 does not reproduce the intensity of the daily mean ${\mathrm{PM}}_{10}$ value, even if the concentration for the day 16 (50.2 ${\mathsf{\mu }\mathrm{g} \mathrm{m}}^{-3}$) is doubled with respect to day 15 (23.3 ${\mathsf{\mu }\mathrm{g} \mathrm{m}}^{-3}$).

Finally, the Demokritos ${\mathrm{PM}}_{10}$ station in Athens is analyzed and the relative measurements are reported in Table 9. The dust plume reached the capital of Greece on the afternoon of 17 April with a peak value of 175 ${\mathsf{\mu }\mathrm{g} \mathrm{m}}^{-3}$ (not showed) at 16:00 UTC. The model configuration run_2 anticipated the peak in the late afternoon of 16 April. The XMed-Dry sample ATS_138 reported a 1% of non-dust composition, in the range of 0.625–10 μm, indicating a ‘pure’ dust event.

This analysis makes clear the importance of XMed-Dry measurements as it reveals the composition and the size distribution of the collected aerosol samples, this is very important when comparing model output and in situ measurements of PM data.

4. Discussion and Conclusions

The ERA5-reanalysis and the Modis/AOD retrievals highlighted a plume intrusion of Saharan dust on 15 April in central/southern Italy. This plume was dislocated toward the Eastern Mediterranean regions on the following day (16 April).

To simulate this event, two different dust emission schemes that are actually implemented in WRF-Chem under the GOCART aerosol package, are tested and compared with experimental data. These include the dry-deposition XMed experimental data, satellite retrievals, reanalysis products, and in situ measurements. The comparison concerns both optical AOD and mass properties (deposition, DUCM, and ${\mathrm{PM}}_{10}$).

The above analysis clearly demonstrates that the AFWA dust emission scheme (run_3) significantly over-predicts the amount of mineral dust in atmosphere, while the UoC-S01 parameterization (run_2) gives overall better results expressed in term of optical properties, dry deposition, total mass, and the surface ${\mathrm{PM}}_{10}$ loading.

This may be explained by considering that the UoC-S01 scheme (run_2), compared to the AFWA scheme (run_3), produces a spatial distribution of emissions that are more localized and restricted to intense emission sites in the semi-arid region [4,25] of the Algerian desert. The modeled AOD for the UoC-S01 scheme corresponds to an intense pulse transported northward on 15 April and eastward the following day. In this context, the spatial extent of the dust plume and the modeled AOD intensity is considerably lower compared to the AFWA scheme. This outcome confirms a sensitivity study performed by Flaunas et al. [24] in the Sahara Desert showing that the AFWA scheme produced exceptionally high AOD values over North Africa compared to the UoC-S01 scheme. It may also be mentioned recent studies in southwest Asia by Khang et al. [48] and LeGrand et al. [23], with analogous results.

However, also on the deposition end of the dust cycle there are probably systematic differences. As it became visible from the comparison of measured ${\mathrm{PM}}_{10}$ with deposition, varying effective deposition velocities exist at the different locations. This and the discrepancy between the model and measured deposition indicate the need to apply more sophisticated size depending deposition schemes. These could for example include ultra-Stokesian effects [57], which might affect deposition velocities by 10% for particles with d >10 µm, or could be non-Stokesian deposition models (a list is given e.g., by Waza et al. [31]), which might modify effective deposition velocities by more than an order of magnitude. While these models are available in principle, the choice and implementation might be difficult due to their diversity [31] and missing information on boundary conditions (e.g., relevant friction velocities). Future investigations need to yield a more robust statistical assessment here. This will need to assess then also the sensitivity of the results to the precise particle size limits chosen for the model comparison. Finally, it is important to remark that when comparing model output and in situ measurements of PM data, the combined analysis of the Xmed-dry samples may provide additional information in the size distribution and composition of aerosol particles.

It is important to point out that, in getting these results, we have performed all the code alterations as suggested by LeGrand et al. [23] and Ukhov et al. [38] that have been incorporated in WRF-Chem starting from version 4.0 (June 2018). However, both emission schemes would benefit from replacing the static erodibility map and soil particle distribution with remote sensed and in-situ observational data. Future works will go in this direction considering the WRF-Chem release above its version 4.0.