Health hazard prospecting by modeling wind transfer of metal-bearing dust from mining waste dumps: application to Jebel Ressas Pb–Zn–Cd abandoned mining site (Tunisia)

Abstract

This work presents a modeling approach to simulate spatial distribution of metal contamination in aerosols with evaluation of health hazard. This approach offers the advantage to be non-intrusive, less expensive than sampling and laboratory analyses. It was applied to assess the impact of metal-bearing dust from mining wastes on air quality for a nearby community and agricultural lands in Jebel Ressas (Tunisia) locality. Dust emission rates were calculated using existing parameterization adapted to the contamination source composed of mining wastes. Metal concentrations were predicted using a Gaussian model (fugitive dust model) with, as input: emission rates, dump physical parameters and meteorological data measured in situ for 30 days in summertime. Metal concentration maps were built from calculated PM10 particle concentrations. They evidence the areas where Pb and Cd concentrations exceeded WHO guidelines (0.5 and 0.005 µg/m³, respectively). Maximum concentrations of Pb and Cd in PM10 are, respectively, of 5.74 and 0.0768 µg/m³ for measured wind speed values up to 22 m/s. Preferential areas of contamination were determined in agricultural lands to the NW from the source dump where Pb and Cd exceeded guidelines up to a distance of 1,200 m. The secondary spreading directions were SW and E, toward the village. Health hazard prospecting shown that a major part of the village was exposed to contaminated dust and that daily hazard quotient (HQ) values reached locally 118 and 158, respectively, for Pb and Cd during the study period. However, HQ variations in the village are high, both temporally and geographically.

Supplementary data: methodology of dust emission and aerosol concentration modeling

Dust emission calculation

The dust emission rate is the sum of the direct soil particle aerodynamic entrainment [f (µg/m²/s)] and of the vertical flux [F (µg/m²/s)] resulting from saltation and sand-blasting processes.

The aerodynamic entrainment function (f) gives the emission rate of mobilized particles on the source surface in the absence of saltation at low wind speed (Loosmore and Hunt 2000 in Shao 2008).

$$ f = 3.6U^{*3} $$

where $U^{*}$ (m/s) is the friction velocity

$$ U* = \frac{KU(z)}{{\ln \frac{z}{{z_{0} }}}} $$

K: von Karman’s constant, z: height above the surface, m, z ₀: roughness height, m, U: wind speed, m/s, F is the emission flux generated by saltation when wind speed is high enough.

The vertical flux F of particles emitted is first based on the calculation of the horizontal flux of soil saltating particles G (g/m²/s) according to Marticorena and Bergametti (1995):

$$ G = EC\frac{{\rho_{a} }}{g}U^{*3} \int\limits_{dp} {\left[ {1 + \frac{{U*_{t} }}{U*}} \right]} \left[ {1 - \left( {\frac{{U^{*}_{t} }}{{U^{*} }}} \right)^{2} } \right]\delta S_{rel} (dp)\delta dp $$

E: ratio of erodible to total surface (taken to 1 here), C: constant of proportionality with a value of 2.61 determined from wind tunnel experiments (White 1979), ρ _a: air density, g: gravity acceleration (m/s²), d _p: soil particles diameter (µm), S _rel: relative surface occupied by a soil aggregate bin of diameter d _p (m²), $U^{*}$: friction velocity of wind (m/s). It depends on wind velocity and on aerodynamic roughness z ₀ of the eroded surface. $U_{t}^{*}$: threshold friction velocity of wind (m/s): minimal wind friction velocity able to mobilize sand particles by saltation (cf Marticorena and Bergametti 1995 for more details).

Initially, this parameterization is designed for desert dust emissions, and the soil aggregate distribution is related to a given textural class using standard distribution (e.g., Zakey et al. 2006). Here, the soil aggregate distribution has been defined directly from in situ measurements according to the following method.

Samples from DIII surface have been taken in 5 different points by scraping off the surface on 1–2 cm over several hundreds of cm². The samples were homogenized to make one representative sample from which 500 g were quartered, got rid of the fraction more than 2 mm and then analyzed for grain size.

First, dry grain size separation was performed with an AFNOR column: 1,600, 1,250, 630, 400, 355, 250, 200, 150, 125, 80 and 63 µm.

Then, grain size distribution of fraction below 63 µm was determined using a Coulter LS 200 laser granulometer with a small water volume. Measurement range is between 0.393 and 905 µm.

In the G formulation, the roughness length z ₀ (µm) is a very sensitive parameter. We have thus determined the roughness from direct measurement of the surface of DIII mining waste dump.

The surface of dump DIII contains erodible and non-erodible or rough elements. There are two types of rough elements:

• Plants and stones that are consistently exposed to the wind. The roughness caused by these elements is uniform for all wind directions. Then, it is called isotropic roughness.

• Ripples with curvilinear forms have different fronts of exposure to wind according to their directions (Fig. 9). Thus, the number of fronts (n) in front of the wind as well as their height (h) and their widths (b) are variable depending on the direction of the wind. The effect of such anisotropic roughness varies with the orientation of the ripple with respect to the direction of the wind. Therefore, wind erosion will not be homogeneous and amounts of emitted dust will vary with wind direction. This effect of anisotropic roughness has been studied to assess the consequences of the presence of grooves of labor in the agricultural land on wind erosion (Saleh et al. 1997; Armbrust et al. 1964).

  Fig. 9

  

  Schematic drawing of a surface waste ripple exposed to wind from different directions with parameters of shape used for aerodynamic roughness calculation

On DIII, surface rough elements are spread in a heterogeneous way. We distinguished between four different parcels in their aspects of surface roughness which was measured for three wind directions:

• in the N direction which is the dominant wind direction measured in situ during the period between July 13th and August 12th, 2009.

• in the SE direction which is the second dominant direction and for which strong winds were common during the period between July 13th and August 12th, 2009.

• in the NW direction which is the annual prevailing direction known for the region.

We performed measurements of roughness on representative surfaces between 0.25 and 4 m² according to the dimensions of the rough elements.

For each parcel, we proceeded as follows:

• 1. We determined the number of isotropic elements (plants, rocks), their widths and their heights.

• 2. For each wind direction, we determined the number of ripple-fronts exposed to the wind, the width (b) and the height (h) of each front. λ represents the totality (n) of fronts (bh) exposed to the wind on a surface (s) with λ = nbh/s.

Pending on the λ value, z ₀ can be calculated using one of the following equations (Marshall 1971; Jarvis et al. 1976; Garrat 1977; Raupach et al. 1980; Raupach 1991) where n, b and h are, respectively, the number, width and the height of the rough element and s is the measurement surface. So that, z ₀ = 0.005 h, if λ < 0.1 and z ₀ = (0.479 λ − 0.001)h, if λ > 0.1.

The aerodynamic roughness was calculated for the other wind directions in the same way.

An overall roughness z _0g which takes into account the different types of surfaces was calculated by summing the different roughness observed on the dump and weighted by the fractions of surface they occupy.

From the saltation flux, the aerosol vertical emission flux is finally calculated according to Alfaro et al. (1997) approach as:

$$ F = G\left[ {\frac{\pi }{6}\rho_{\text{p}} \beta p_{{i{\text{dp}}}} \frac{{d_{i}^{3} }}{{e_{i} }}} \right] $$

where i represents three modes of a predefined log-normal emission distribution. G: the saltation flux (µg/m²/s), ρ _p: density of particles, β: proportionality constant, p _idp : fraction of kinetic energy of particles of each diameter class dp able to generate aerosols of mode i, d _i : mass of aerosols of mode i (g), e _i : kinetic energy of aerosols of mode i (g cm²/s²).

The total vertical flux is obtained by integrating all flux resulting from each diameter class of surface particles and then by summing the quantities of emitted aerosol.

The model we use did not consider the local turbulence related to the topography of the site. We alleviate this point by calculating emission rate with wind gusts measured every minute. Also, the surface moisture was not considered as we worked in summer season.

Calculation of airborne dust concentrations

PM10 concentrations (µg/m³) have been calculated with the fugitive dust model (FDM, Winges 1992). FDM is a Gaussian plume model designed for computing particle concentrations and deposition rates from fugitive dust emissions. The model is well accepted by USEPA for this purpose. It incorporates transport, dispersion and deposition of pollutants in the atmosphere, using input data for particulate matter (particle diameter, density) and air parameters (wind velocity, wind direction, temperature) measured in the site together with flux of emitted PM10 from the emission model. For the 30 days, 720 data for hourly mean emission rate of PM10 were integrated on FDM.

Pb and Cd concentrations in PM10 were calculated using metal content of dump DIII surface.

Source and receptors

The source, DIII waste dump, is assimilated to a rectangle of 235 m × 150 m.

We considered 2,000 receptor points obtained with a 100 × 100 m gridding over a surface of 20 km². We also added 30 points over the village.