Atmospheric Dispersion Modellingof Volcanic Emissions

Department of Earth Sciences Licentiate Thesis

2015

Atmospheric Dispersion Modelling of Volcanic Emissions

Adam Dingwell

Abstract

Gases and particles released by volcanoes pose a serious hazard to humans and society. Emis- sions can be transported over long distances before being reduced to harmless concentrations.

Knowing which areas are, or will be, exposed to volcanic emissions is an important part in reducing the impact on human health or society. In this thesis, the dispersion of volcanic emis- sions is studied using a set of atmospheric models. Two case studies have been performed, one studying potential ash emission from future eruptions on Iceland, and a second covering SO ₂ emissions from Mt. Nyiragongo in D.R. Congo.

The first study covers long range (∼1,000 km) dispersion of fine ash from explosive erup- tions. Three years of meteorological data are used to repeatedly simulate five eruption scenarios.

The resulting concentrations of airborne ash at different times relative the onset of each eruption is compared to current and previous threshold concentrations used by air traffic controllers. The ash hazard showed a seasonal variation, with a higher probability of efficient eastward transport in winter, compared to summer; summer eruptions pose a more persistent hazard.

In the second study, emissions of SO ₂ from passive degassing at Mt. Nyiragongo is studied

over a one–year period. The meteorological impact on the dispersion is studied by assigning

a fixed emission source. Furthermore, flux measurements from the remote sensing data are

used to improve the description of the emission source. Gases are generally transported to the

north-west in June–August and to the south-west in December–January. A diurnal variation due

to land breeze around lake Kivu contributes to high concentrations of SO ₂ along the northern

shore during the night. Daily averaged concentrations in the city of Goma ( ∼15 km SW of the

source) exceeded the European Union’s air quality standard (125 μg/m ³ ) for 120-210 days over

a one year period.

Sammanfattning

Gas- och partikelutsl¨app fr˚an vulkaner utg¨or en fara f¨or m¨anniskor och f¨or v˚art samh¨alle. Ut- sl¨appen kan transporteras ¨over l˚anga avst˚and innan de reduceras till ofarliga halter. Att k¨anna till vilka omr˚aden som uts¨atts, eller kommer uts¨attas, f¨or utsl¨appen ¨ar ett viktigt verktyg f¨or att minska p˚averkan p˚a folkh¨alsa och samh¨allet. I den h¨ar avhandlingen studeras spridningen av ut- sl¨app fr˚an vulkaner med hj¨alp av en upps¨attning atmosf¨arsmodeller. Tv˚a fallstudier har utf¨orts, en fokuserar p˚a vulkanaska fr˚an potentiella framtida utbrott p˚a Island, den andra studerar SO ₂ - ustl¨app fr˚an Nyiragongo i Demokratiska Republiken Kongo.

Den f¨orsta studien beskriver l˚angv¨aga (∼1,000 km) transport av aska fr˚an explosiva utbrott.

Tre ˚ar av meteorologiska data anv¨ands f¨or att modellera spridningen fr˚an fem olika utbrotts- scenarier f¨or varierande v¨adersituationer. Koncentrationen av luftburen aska studeras vid olika tidpunkter relativt utbrottens starttid och j¨amf¨ors med tidigare samt befintliga gr¨ansv¨arden f¨or flygtrafik. Sannolikheten f¨or skadliga halter aska varierar med ˚arstid, med en h¨ogre sannolikhet f¨or effektiv transport ¨osterut under vinterm˚anaderna, j¨amf¨ort med sommarm˚anaderna; sommar- utbrott ¨ar ist¨allet mer ben¨agna att orsaka l˚angvariga problem ¨over specifika omr˚aden.

I den andra studien modelleras utsl¨app av SO ₂ fr˚an passiva utsl¨app vid Nyiragongo ¨over en

ett˚arsperiod. Den meteorologiska effekten p˚a spridningen studeras genom att anv¨anda en kon-

stant utsl¨appsk¨alla. Dessutom studeras spridningen mer i detalj genom att anv¨anda fj¨arranalysdata

f¨or att b¨attre uppskatta utsl¨appen. Gaserna transporteras i regel mot nordv¨ast i juni–augusti och

mot sydv¨ast i december–februari. En sj¨o-/landbriscirkulation runt Kivusj¨on orsakar h¨oga halter

av SO ₂ l¨angs sj¨ons norra strand nattetid. Dygnsmedelkoncentrationer av SO ₂ i provinshuvud-

staden Goma ( ∼15 km sydv¨ast om Nyiragongo) ¨overskred EU-riktlinjer (125 μg/m ³ ) under

120-210 dagar under en ett˚arsperiod.

Science is made up of so many things that appear obvious after they are explained

Pardot Kynes, Muad’Dib

List of papers

This thesis is based on the following papers, which are referred to in the text by their Roman numerals.

I Dingwell, A. and Rutgersson, A. (2014). Estimating volcanic ash hazard in European airspace. Journal of Volcanology and Geothermal Research, 286(0):55–66.

II Dingwell, A.; Rutgersson, A.; Claremar, B.; Arellano, A.; Mapendano, Y. and Galle, B. (2015). Seasonal and diurnal patterns in the dispersion of SO2 from Mt. Nyiragongo. Manuscript prepared for submission Reprints were made with permission from the publishers.

In Paper I, I had the main responsibility for the data analysis and the writing of the paper.

In Paper II, I was responsible for most of the data analysis, except for the processing of the flux-data. I had the main responsibility for writing the paper.

The method for large scale application of the modelling system was mainly de- veloped as part of Paper I. Some adaptations were made in Paper II. During the project, I have also contributed to the development of FLEXPART-WRF.

The following related publication is not included in the thesis:

Brioude, J.; Arnold, D.; Stohl, A.; Cassiani, M.; Morton, D.;

Seibert, P.; Angevine, W.; Evan, S.; Dingwell, A.; Fast, J. D.;

Easter, R. C.; Pisso, I.; Burkhart, J. and Wotawa, G. (2013). The

Lagrangian particle dispersion model FLEXPART-WRF version

3.1. Geoscientific Model Development, 6(6):1889–1904.

Contents

1 Introduction . . . . 11

1.1 Aim of this thesis . . . . 12

2 Modelling tools . . . . 13

2.1 Meteorological data . . . . 13

2.1.1 Downscaling meteorological data . . . . 13

2.1.2 Meteorological model . . . . 14

2.2 Dispersion modelling . . . . 15

2.2.1 Advection of particles . . . . 15

2.2.2 Deposition processes . . . . 16

2.2.3 Aggregation and Chemistry . . . . 18

2.2.4 Plume rise from strong eruptions . . . . 18

3 Model application . . . . 20

4 Results . . . . 22

4.1 Probability of volcanic ash over northern Europe . . . . 22

4.2 Dispersion of SO 2 from Mt. Nyiragongo . . . . 22

5 Concluding remarks . . . . 28

6 Acknowledgments . . . . 29

References . . . . 30

1. Introduction

Volcanoes can emit huge quantities of hazardous ash and gases into the at- mosphere. Both components can have a negative impact on human health or on important components of our society. Ash from large eruptions can re- main in the atmosphere for days or even weeks before removed or dispersed to harmless concentrations. Gas emissions, either during eruptions or as passive degassing, can pose an serious hazard for the environment and human health.

This can occur either during eruptions or as passive degassing over longer pe- riods. Acidic compounds (e.g. SO ₂ , H ₂ S, HCl, HF) affect cloud chemistry and produce acid rain. Exposure to gas plumes can therefore lead to increased wear on man-made structures (e.g. roofs, antennas, machinery), and have an irritating effect on airways when inhaled (Delmelle et al., 2002). Fluoride (F), which is mainly emitted as HF, is toxic in high concentrations and can disturb dental development for children (Baxter et al., 1999; Delmelle et al., 2002) or even prove fatal for grazing animals, which might ingest large quantities of fluoride coated plants (Cronin et al., 2000).

The dispersion of volcanic emissions over Europe has occurred several times in history, the most notable caused by the Lakag´ıgar fissure eruption in Iceland in 1783–1784. The eruption had severe local impact and likely af- fected large parts of Europe as well. Within a year from the eruption, most of Icelands livestock died after ingestion of ash coated grass. During this period, famine and direct exposure to toxic emissions reduced the human pop- ulation of Iceland by 19–22 %, or approximately 10,000 people (Demar´ee and Ogilvie, 2001). A persistent haze was reported over Europe in mid to late June 1783, lasting for about two months and eventually extending as far as Moscow, Baghdad and northern Africa (Stothers, 1996).

Gas emissions from volcanoes is not uncommon, the main contribution is not from eruptions, but rather from passive degassing (?). While individual cases of passive degassing are weaker than events like the Lakag´ıgar eruption, they still pose a regional hazard. Areas downwind of degassing volcanoes often show negative impact on human health and the environment (?).

In April 2010, the eruption of the Icelandic volcano Eyjafjallaj¨okull re-

sulted in almost complete closing of European airspace. Fine grained ash

was transported over 1000 km, reaching continental Europe in a matter of

days. Between th 15th and 22nd of April, around 104,000 flights were can-

celled (EUROCONTROL, 2010). The total emissions of the eruption were

not remarkable, however, the fraction of fine ash produced was unusually high

and unusually persistent northerly winds allowed a large portion of this ash

11

to reach continental Europe (Stohl et al., 2011). Furthermore, the European countries were not prepared for such an events; as a response, new guide-lines for volcanic ash were implemented for aviation and regular exercises were initiated.

The impact of such extreme events is highly dependent on the information available to decision makers. Increased knowledge and better predictions en- able our society to make preparations and better manage the situation when eruptions occur.

1.1 Aim of this thesis

In this thesis, exposure to volcanic ash and gases are studied and techniques for estimating hazards for different types of volcanoes are presented.

• The probability of encountering volcanic ash in European airspace as the result of eruptions on Iceland is covered in Paper I.

• Regional exposure to gases emitted from Mt. Nyiragongo (eastern Demo- cratic Republic of Congo) is covered in Paper II.

12

2. Modelling tools

The dispersion of pollutants in the atmosphere depends on meteorological con- ditions as well as the properties of the emissions in question. The atmosphere is a complex system with interactions on a wide range of scales, from micro- physical processes and chemistry, to global energy fluxes and Rossby waves.

The complexity makes it hard to study the atmosphere as a whole, which is why atmospheric models are important tools for studying the fate of pollu- tants in the atmosphere.

A combination of two models was used throughout this project; a meteoro- logical model was used to produce wind fields and precipitation data for use in a dispersion model. The tools used for large scale application of the models are similar between the studies but has undergone development throughout the project. Most of the work has gone into the dispersion model. The dispersion model has been improved upon as part of the work.

2.1 Meteorological data

When running atmospheric models for research purposes covering past cases, one aims at having the best estimate of the state of the atmosphere to force the model. One type of product intended to provide such information is reanalysis data (Warner, 2011). A reanalysis product is the combined result of a large set of observations and a meteorological model, often run over several decades to provide long time series of consistent data. However, reanalysis data is typically produced by global models— running on relatively coarse resolution.

2.1.1 Downscaling meteorological data

The initial meteorological data, as was used in this project, was the ERA In- terim reanalysis product (Dee et al., 2011), provided by the European Centre for Medium range Weather Forecasts (ECMWF). The data was retrieved on a 0.75 ^◦ resolution ( ∼80 km in meridional direction), which is not directly appli- cable for regional studies, where a resolution of several kilometres if typically desired. In order to improve the resolution, a limited area model can be used.

The limited area model is forced by the reanalysis data but runs on a higher

resolution. This process of increasing the resolution by the use of a meteo-

rological model is called downscaling. Downscaling requires large computa-

tional resources compared to simpler methods (interpolation) but allows small

13

scale variations in land-use and topography to be accounted for. An example of this is shown in Figure 2.1; while interpolation produces smooth contours, it lacks the additional information provided by the dynamic downscaling.

Figure 2.1. Comparison of near-surface temperature from basic interpolation (left) versus dynamic downscaling (right) from Dingwell (2012). The left image was pro- duced by interpolating data from ERA-Interim at 0.75 degree resolution (≈ 35 × 85 km) down to 13.5 × 13.5 km. The right image has the same resolution, but was produced by the WRF-model, using the interpolated data from ERA-Interim as forc- ing. Image

If the desired resolution is differs much from the initial resolution, then downscaling can be made in several steps, called nesting to further improve resolution (Jacobson, 2005). There are two reasons why nesting is recom- mended rather than simply setting up a single high resolution domain directly.

First, there could be characteristics in the terrain (e.g. mountains, islands or lakes) which are too small to be resolved by the global model but which have an important influence of the dynamics of the region. With nesting, these in- fluences can be accounted for. Second, numerical problems may arise if there is a too large difference between the input data and the computational grid (Warner, 2011).

2.1.2 Meteorological model

The Advanced Weather Research and Forecasting (WRF) model (Skamarock et al., 2008) was used to downscale meteorological data throughout this project.

The WRF-model is a 3-dimensional, non-hydrostatic, mesoscale meteorologi- cal model. Its modular design allows it to be used in various applications from operational forecasting to idealized scientific cases.

When initiating a run (i.e. cold starting) with the WRF-model, the model

domains are initiated in all grid points using interpolated data from the forcing

data. Some time, usually 6-12 simulation hours are needed for the model

to fully develop the small scale dynamics not present in the coarser data when

14

running from a cold start. This process is referred to as spin-up. Data produced during the spin-up is inaccurate and should be discarded.

After initialization, WRF will continue to apply forcing to the boundaries of the outer domains but no nudging (i.e. adjustment toward forcing data) takes place within the domains. Since the WRF-model makes no attempt at nudging data inside domains during a simulation, a set-up with large outer domains risk deviating too far from the reanalysis data. In order to prevent this, we run the model in segments. Each segment consists of three steps. First, the model is initialized using the interpolated reanalysis data (cold start) and run for 6–12 hours to allow the model to reach a state of quasi–equilibrium under the applied forcing (spin-up). Second, the simulation is paused, and some intermediate processing of model variables is done. The most important part being to overwrite accumulated fields (e.g. precipitation) with values from the end of the previous run segment (used in Papers I and II). Third, the model simulation is continued running for 24–48 hours for the main production of data. We have aimed at keeping each run segment short since this made it easier to fit jobs into the queues at compute centres.

2.2 Dispersion modelling

A dispersion model either works as a stand-alone model using prepared meteo- rological data (offline) or combined with a meteorological model running side- by-side (online). While online models can benefit from allowing emissions to affect the dynamics through various feedback processes, they are generally more computationally demanding. An offline model can run independently of the meteorological model as long as there is available data.

This project uses the (offline) Lagrangian Particle Dispersion Model (LPDM) FLEXPART-WRF (Brioude et al., 2013; Stohl et al., 2005). An LPDM works by releasing a large number of computational particles (e.g. 100 000) from a point, line, area or volume source, located anywhere within the meteorological grid. The computational particles are moved through the domain, according to mean winds and turbulence which is interpolated to the exact position of each par- ticle. Since computations are only made where pollutants are present, LPDM’s are suitable for simulating emission of chemically stable gases or particles from a limited number of sources over large areas. The main downside of Lagrangian models is that they are less suitable for simulating chemistry or aggregation of particles.

2.2.1 Advection of particles

The position of particles in LPDMs is determined by a numerically solving dif-

ferential equations. There are several methods for doing this, with the simplest

15

being the first order Euler solution (as in FLEXPART-WRF):

¯r (t + Δt) = ¯r(t) + ¯V(¯r,t)Δt (2.1) where ¯r (t) is the particle’s position at any given time, t, Δt is the time step of the simulation and ¯r (t + Δt) is the particle’s position at the next step. More elaborate methods can be used to improve precision, such as using a Second order Euler solution (e.g. the NAME mode: Jones et al., 2007) or by using (2.1) to make an initial guess position, which can further be improved by determining the mean wind along the guessed path (e.g. the HYSPLIT model: Draxler and Hess, 1998).

When simulating particulate matter, the dispersion model should take in to account the terminal velocity (i.e. settling due to gravity). For spherical particles larger than ∼10 μm, the terminal velocity, v g , can be calculated as a balance between the downward gravitational force and the upward drag force:

v _g =

 4 3

d p g C _D

 ρ − ρ air

ρ air



≈

 4 3

d p g C _D

 ρ ρ air



(2.2)

where d _p is the particle’s diameter, g is acceleration due gravity, ρ is the particle density and ρ air is the density of air. C D is the drag force which, in FLEXPART and its derivatives, is calculated using the method by N¨aslund and Thaning (1991). In order to accurately describe terminal velocities of particles down to micrometer size, (2.2) can be adjusted by dividingC _D with the Cunningham correction factor as proposed by Cunningham (1910); Knudsen and Weber (1911):

C cun = 1 + 2λ d _p



A 0 + Q · e ^{−C 2 dp λ}



(2.3) where λ is the mean free path of gas molecules in air. A 0 , Q and C are dimen- sionless constants, which in FLEXPART-WRF are set to 1.257, 0.400 and 1.10, respectively.

Furthermore, (2.2) assumes spherical particles, particles of different shape will have a different terminal velocity. This is usually corrected by assigning an equivalent spherical diameter (i.e. aerodynamic diameter) for particles when setting up simulations.

2.2.2 Deposition processes

Dry deposition in models are usually defined in terms of the deposition veloc- ity, V _d , defined as the fraction between the flux, F, and the concentration, C, of given species. FLEXPART-WRF calculates the flux using the resistance method (Wesely, 1989):

V _d = 1

r _A + r B + r C

(2.4)

16

where r _A is the aerodynamic resistance of the surface layer and r _B is the viscous resistance of the quasi-laminar sublayer, the implementation of these is explained by Stohl et al. (2005). The final term, r _C is the bulk surface resistance, which is dependent on thickness and type of vegetation, presence soil, water, or buildings, etc., see Wesely (1989); Wesely and Hicks (2000) for details. An illustration of this is shown in Figure 2.2.

A

B

C

net flux

Figure 2.2. Illustration of the different sub-processes of dry deposition. (A) repre- sents the surface layer, where the resistance is determined by turbulent eddies. (B) is the laminar sublayer, a thin layer of air in contact with surfaces, where turbulence is negligible and viscous forces dominate. Transport through this layer depend on Brow- nian motion and inertia of heavier particles. (C) represents the bulk surface properties, which combines numerous processes for absorption of gases by the surface.

In FLEXPART-WRF r C is estimated using land-use categories following We- sely (1989); Wesely and Hicks (2000). A problem with the land-use data was en- countered in Paper II, when the spatial resolution was set to a much higher value than the typical resolution in FLEXPART. The land-use data was on a coarse 0.3 ^◦ grid, resulting in patchy deposition in high resolution grids. Therefore, FLEXPART-WRF was modified to read land-use data from the WRF model. At a horizontal resolution of 2 km, the model results differ significantly depending on which land-use data is used, an example of this is shown in Figure 2.3.

Particles behave differently than gases, mainly by the exclusion of surface resistances and the added effect of gravitational settling. In FLEXPART-WRF the deposition velocity of particles is given by:

v d = 1

r _A + r B + r A r _B v _g + v g (2.5) where r _B is the same as for gases but including an impaction factor.

Wet deposition in FLEXPART-WRF has recently undergone major changes.

Wet deposition was first officially added by Brioude et al. (2013), following the same model as was used in FLEXPART 9.02. Currently, the model documenta- tion describes a slightly different approach than what the model i programmed to do. An up-to-date description of wet deposition of particles was therefore included in Paper I.

17

0.75 1.25

0.95 1.05

0.85 1.15 1.5

1.3

1.1 0.9 0.7 0.5 (A) (B)

Figure 2.3. Difference between FLEXPART-WRF simulations using different land- use data (figure from Paper II. Both figures show average values from a simulation with land-use data from WRF divided by results from a simulation using the old built- in data. (A) shows differences in accumulated dry-deposition of SO 2 over one year.

(B) shows difference in average SO 2 -concentration 0-500 m a.g.l. during September–

November 2010.

2.2.3 Aggregation and Chemistry

Chemistry and particle aggregation affect the aerodynamic properties of emis- sions as well as the removal rate. These processes are not easily included in LPDMs, some attempts exist where gas-phase chemistry is calculated on fixed Eulerian grids (e.g. Oettl and Uhrner, 2011) or with simple redistribution of mass within computational particles (e.g. Businger et al., 2015). FLEXPART-WRF can only handle chemistry as an additional removal process, conversion between modelled species is not possible. An attempt to estimate the model error arising from the lack of sulphur chemistry is made in Paper II.

2.2.4 Plume rise from strong eruptions

In Paper I, the dispersion of volcanic ash from a range of eruption scenarios was studied. An important part in setting up the study was assigning the emission altitude. In reality, the emissions take place at the vent of the volcano, a mixture of ash and gases are released into the atmosphere at high velocity. Ambient air entrained into the plume is heated; convection carries the emissions to higher al- titudes (Sparks and Wilson, 1976; Sparks, 1986). This process can be studied in advanced Eulerian models (Textor et al., 2006), however, the high computational demand of such models, currently make them unsuitable for long range simula- tions. Instead, simplifications need to be made regarding the emission source.

A common method is to set up a source volume were the ash is released. This volume typically represent part of the plume where convection is still important.

Emission sources in earlier studies and operational models were assigned as

vertical columns with uniform mass released throughout the column (Witham

et al., 2007). However, since volcanic plumes often are convective (Sparks,

18

1986), the majority of fine ash should be expected at the top of the eruption col- umn due to their smaller settling velocity (Carey and Sparks, 1986; Koyaguchi and Ohno, 2001). Other studies, (e.g. Peterson and Dean, 2008; Steensen et al., 2013), use a top-weighted function to determine the particle release in the source column. However, this introduces another error, since concentration in the source column increases with altitude. If the plume is convective, concentration should decrease with altitude due to volume expansion. This process becomes more important with taller the eruption columns. In Paper I, the source volume was designed to account for an increase in fine ash at higher altitudes by distributing mass over ten stacked source segments according to a Poisson distribution. Sim- ilar methods have previously been used by Peterson and Dean (2008); Steensen et al. (2013), but without accounting for varying sizes of the source volumes. The width of each segment, d, varied between layers determined by the mass release:

d = M ˙ (z i )

U (z i ) ˜CΔz (2.6)

where ˙ M(z i ) is the mass release rate in a given source segment (z i ), U (z i ) is the average wind speed within a segment, Δz is the segment thickness, and ˜C is a constant proportional to the concentration within each segment. This method reduces the overestimation of concentration in the upper portions of the plume.

The method can be further improved by including atmospheric density profiles.

In Paper II, the plume rise was typically much lower (i.e. several 100 metres), making plume expansion less important. Therefore, the classic approach of a uniform source was used instead.

19

3. Model application

The studies covered in this thesis apply the same modelling system to two highly different cases. Different methods for building statistics over exposure to emis- sions were used. In Paper I, five different eruption scenarios were studied, rep- resenting potential future eruptions on Iceland. Each scenario was repeatedly simulated, for varying meteorological conditions, to determine which areas are most likely to be exposed to hazardous concentrations of ash. The area covered by the study is shown in Figure 3.1A. In order to minimize systematic errors from any diurnal variation, eruption scenarios were initiated every six hours, produc- ing four cases for each scenario and day. Each of these cases was simulated for 96 hours. In total, the meteorological model runs covered 3 years, corresponding to an equivalent of 240 years of dispersion simulations. The data was analyzed relative to the onset of each eruption case, and hourly average concentrations were compared to thresholds used by air traffic control (ATC) for declaration of no-flight zones.

~

Provincial capital Goma

(population: ~1 000 000) Lake Kivu Mt. Nyamuragira (3058 m)

Mt. Nyiragongo (3470 m)

d04

60°N

50°N

40°N

30°N

(A) (B)

Figure 3.1. Geographical coverage of the two studies included in this thesis. (A) shows three nested domain, all of which were used in the dispersion simulations in Paper I. (B) shows the single domain (third order nest) used in the dispersion simula- tions in Paper II

In Paper II, a case of passive degassing was studied. Flux data from the Net-

work for Observation of Volcanic and Atmospheric Change (NOVAC) by Galle

et al. (2010) was used to study the dispersion of emissions from Mt. Nyiragongo,

an active volcano in eastern Democratic Republic of Congo (DRC), see Fig-

ure 3.1B. The flux data was used as emission strength when available, however,

since the observations could only be made in daylight, there were regular gaps in

20

the data. In order to cover night-time cases as well as other periods lacking data, any gaps in the data were filled by randomly sampling flux data from the avail- able measurements. This procedure was repeated 30 times, creating 30 different time series but with common emissions whenever observations were available.

A dispersion simulation was conducted for each of the emissions cases, forming an ensemble with 30 members.

21

4. Results

4.1 Probability of volcanic ash over northern Europe

The probability of exposure to ash as a result of future volcanic eruptions was studied in Paper I. The hazard was found to vary with season, with the strongest difference seen between summer (June–August) and winter (December–January).

An example of this is shown in Figure 4.1; the hazard from a summer time erup- tion is more persistent but initially lower for areas east of Iceland. Wintertime eruptions have a higher probability of affecting most of the Scandinavian penin- sula and areas around the Baltic sea. Since the polar front is weaker during sum- mer, the westerly winds are weaker and less persistent, resulting in a higher prob- ability of ash being transported in other directions. In winter, however, the strong polar front, results in strong westerly winds most of the time, transporting ash eastward in most of the cases. This transport is both more frequent as well as more efficient, resulting in high concentrations of ash at greater distances from the source.

A comparison was made between different eruption scenarios, an example of this is shown in Figure 4.2A-E, using a threshold value of 0.2 mg/m ³ ; this threshold corresponds to the no-flight condition in use prior to the eruption in 2010 (Webster et al., 2012). Note that the Eyafjallaj¨okull eruption case (covering the most intense phases of the eruption) has a probability below 10 % for most land areas beyond Iceland. The weakest scenario (Mt. Spurr) shows no exceedances for this period, partly due to the short duration of the eruption. In general, ash at higher altitude is more consistently transported eastwards compared to lower levels. However, the longest transport is seen in the mid-level ( ∼6–11km)where the jet stream is expected.

A comparison is also made with earlier results by Leadbetter and Hort (2011) (Figure 4.2F), who used a similar approach. The scenario set up by Leadbetter and Hort (2011) should, in terms of emission strength, correspond to the weakest scenario in Paper I (i.e. Figure 4.2A). However, it turns out to be more similar to the Askja-1875 scenario (Figure 4.2C), suggesting that ash hazard is much lower than previously estimated.

4.2 Dispersion of SO ₂ from Mt. Nyiragongo

The geographical distribution of gases emitted from Mt. Nyiragongo was stud-

ied in Paper II. At first only the meteorological influence on the dispersion was

22



  



















































































    

    

    

    

    

    

    

    

    

    

    

    

















































Figure 4.1. Probability of ash concentration exceeding the no-flight threshold (2.0 mg/m ³ ) for an eruption similar to Askja-1875 (from Paper I). This scenario is de- signed to be roughly 10 times stronger than the Eyjafjallaj¨okull eruption in 2010. The probability is calculated for three different periods: yearly (A–D), winter (E–H) and summer (I–L). Different time periods are given, specifically (top to bottom) 0–24, 24–

48, 48–72 and 72–96 hours relative the onset of the eruption. All cases are for flight levels FL200–FL350 (∼ 6–11 km a.s.l.).

23

  













  


  

















  



  























  





 











  !

"#$ ! % & '( '')%!*+

Figure 4.2. Areas where the probability to of exceeding hourly average concentrations of 0.2 mg/m ³ is at least 10 % within 24–48 h after the onset of an eruption (from Paper I). Five different eruption scenarios, based on historic events, are presented (A-E), as well as results from Leadbetter and Hort (2011) (F). Contours are given for three different flight levels (given in units of 100 ft.) as indicated by the grey-scale.

studied by using a constant emission source. The mean dispersion direction over the whole year was to the west (as expected). However, a seasonal shift was seen with the higher portion of the plume being transported further north in December–

January compared to the annual mean. This can be seen in the plume cross section shown in Figure 4.3, which stretches along the line marked in Figure 3.1B. Lower portions of the plume was instead transported to the south. In June–August, the pattern was reversed with northward transport of the lower and southward trans- port of the higher parts of the plume. This skewness corresponds with the migra- tion of the Inter Tropical Convergence Zone (ITCZ), which is located to the north in July, and to the southwest in January. Surface winds converge at the ITCZ, which is why lower portions of the plume vary with season. At higher altitudes, air flows away from the ITCZ which causes the skewness seen in Figure 4.3.

Figure 4.3 also shows a diurnal variation. In daytime, emissions are more likely to remain over land (i.e. north of 1.6 ^◦ S in Figure 4.3). In nighttime, emis- sions tend to form a shallow layer with high concentration of SO ₂ over lake Kivu.

This is related to a lake-/land breeze cycle forming over lake Kivu, especially

evident around the equinoxes, when the ITCZ has the least meridional influence

on the flow. Lake Kivu lies on a rift valley with steep slopes (escarpments) on its

east and west sides, with a valley located about halfway on the eastern side. The

escarpments create a channel effect, preventing air exchange other than at the

24

DJF S ON

8 6 4 2

Height above sea level [km] JJA MAM

Local Time

06:00 13:00

S eason

10⁵ 10³

10² 10⁴ ng/m3

10 8 6 4 2

10 8 6 4 2

10 8 6 4 2

-1.8° -1.6° -1.4° -1.2° -1.8° -1.6° -1.4° -1.2°

Position along cross section (degress latitude) 10

8 6 4 2

10 8 6 4 2

10 8 6 4 2

10 8 6 4 2

Daily

-1.8° -1.6° -1.4° -1.2°

Figure 4.3. Modelled seasonal average diurnal variation of SO ₂ along the cross section marked in Figure 3.1. The leftmost column shows average concentrations for each 3- month period of a full year. The middle middle column shows average concentrations between 05:30–06:30, local time (i.e. around sunrise) while the rightmost columns shows average concentrations between 12:30–13:30 (i.e. at maximum insolation).

The two leftmost columns represent the diurnal extremes regarding the concentrations over lake Kivu. The data covers one year (April 2010 through March 2011) and uses a constant emission source. Figure from Paper II

25

Fraction of days above 125 μ g/m

Month Kitchanga

Kahe Masisi Sake Goma E Goma W -1°

-2°

27° 28° 29°

(A) (B)

Kitchanga

Kahe Masisi

Sake Goma E

Goma W

Figure 4.4. Monthly variation of SO ₂ exposure in populated areas near Mt. Nyi- ragongo (from Paper II). (A) Map showing the total number of days during which the 24-hour average SO 2 concentration exceeded 125 μg m ⁻³ (reference value used by the European Union). (B) Fraction of day when the reference value was exceeded each month for locations marked in (A). Values are presented for 12 months, from April 2010 to February 2011. The lines represent the average result from an ensemble of 30 dispersion simulations, the shaded areas show the 25–75 percentiles of ensemble members.

northern and southern shores (or through the eastern valley). This concentrates the lake-/land breeze flow at the northern and southern shores.

Using the flux data from the NOVAC network, an estimate of past SO ₂ concen- trations in populated areas was made. Results from this are shown in Figure 4.4, where average concentrations are compared to reference values in use by the Eu- ropean Union. Exceeding these values over longer periods could be connected to an increase in mortality (WHO, 2006). In the city of Goma ( ∼1,000,000 in- habitants), the concentrations exceeded the reference value a total of 120–210 days, with the highest exposure in the western parts. Over the region, the highest exposure was seen along the shore of Lake Kivu (Goma and Sake) and to the NW of Mt. Nyiragongo (Kahe). The area with at least 10 days above the reference value extends to about 170 km NW of Mt. Nyiragongo.

The exposure strongly depends on the season. In April–August (2010), areas to the NW saw the highest exposure. Kitchanga ( ∼85kmNW)andKahe(∼35km NW) exceeded the reference value ∼45 % and >80 % of the time, respectively.

During the same period, concentrations in Sake (W) and Western Goma (S) where above the the reference value ∼20 % of the days, while eastern Goma was below most of the time. Most of the exceedances in Sake and Goma were related to the land-breeze forming during the night. In September, the conditions changed.

Goma and Sake saw a steep increase in SO 2 while concentrations in Kitchanga

26

and Kahe decreased. The highest exposure at the norhtern shore of lake Kivu was reached in November 2010 through March 2011. Western Goma exceeded the reference value 90 % of the days, with slightly better conditions in Sake.

27

5. Concluding remarks

In Paper I we showed seasonal variations in dispersion for different eruption scenarios on Iceland. The more persistent westerly winds in winter, makes it more likely for ash to be transported westward, typically reaching further than if the same eruption would occur in summer. Summertime eruption are less pre- dictable over long term, since the otherwise dominating westerly winds are less frequent. However, we also found that winter conditions vary more between the years studied compared to summer. In order to further improve the accuracy of the probability analysis, it might be more important to study additional winters than full years.

Paper II demonstrated the first attempt at a regional dispersion study around Mt. Nyiragongo. The combination of observational data from DOAS retrieval with a regional modelling system is an important step in studying health prob- lems in the region. The study shows the importance of studying both daytime and nighttime conditions as well as seasonal variations. The modelled air concentra- tions of SO ₂ suggest suggest that a large number of people are being exposed to hazardous amounts of SO 2 , which can have both long- and short term negative impact on human health.

While it has been shown that the modelling system can be applied to differ- ent cases of volcanic emissions, continued development is still important. A diagnostic plume rise scheme could be included to better describe weaker vol- canic plumes (or industrial sources), instead of using external tools to adjust it. A built-in scheme would benefit from continuously updated meteorological fields, allowing development over time. HCl and HF are two important gases in vol- canic plumes, which are currently not included in FLEXPART-WRF. Adding them would require new wet scavenging and dry deposition coefficients. In addi- tion, SO 2 is only removed through dry deposition, despite being a highly soluble gas; adding wet deposition for SO ₂ should be considered.

28

6. Acknowledgments

I would like to thank my supervisors Anna Rutgersson and Bj¨orn Claremar, for their support throughout this project. Special thanks to Anna for giving me this opportunity in the first place.

I thank Jerome Brioude (NOAA) for being so supportive with addressing is- sues with the dispersion model and helping me join the development team. San- tiago Arellano helped a lot with the NOVAC-data and had great patience with my many questions.

I would also like to thank the many PhD-students in CNDS and at the depart- ment, especially my office buddies, Marc and Jean-Marc for plenty of disastrous discussions.

Finally, a big thank you to friends and family, for putting up with my weird working hours and all modelling jargon which I keep blabbing on about.

29

References

Baxter, P. J., Baubron, J.-C., and Coutinho, R. (1999). Health hazards and disaster potential of ground gas emissions at furnas volcano, s˜ao miguel, azores. Journal of Volcanology and Geothermal Research, 92(1–2):95 – 106.

Brioude, J., Arnold, D., Stohl, A., Cassiani, M., Morton, D., Seibert, P., Angevine, W., Evan, S., Dingwell, A., Fast, J. D., Easter, R. C., Pisso, I., Burkhart, J., and Wotawa, G. (2013). The Lagrangian particle dispersion model FLEXPART-WRF version 3.1. Geoscientific Model Development, 6(6):1889–1904.

Businger, S., Huff, R., Pattantyus, A., Horton, K., Sutton, A. J., Elias, T., and Cherubini, T. (2015). Observing and forecasting vog dispersion from Kilauea volcano, Hawaii. Bulletin of the American Meteorological Society.

Carey, S. and Sparks, R. (1986). Quantitative models of the fallout and dispersal of tephra from volcanic eruption columns. Bulletin of Volcanology, 48:109–125.

Cronin, S. J., Manoharan, V., Hedley, M. J., and Loganathan, P. (2000). Fluoride: A review of its fate, bioavailability, and risks of fluorosis in grazed-pasture systems in new zealand. New Zealand Journal of Agricultural Research, 43(3):295–321.

Cunningham, E. (1910). On the velocity of steady fall of spherical particles through fluid medium. Proceedings of the Royal Society of London. Series A, Containing Papers of a Mathematical and Physical Character, 83(563):357–365.

Dee, D. P., Uppala, S. M., Simmons, A. J., Berrisford, P., Poli, P., Kobayashi, S., Andrae, U., Balmaseda, M. A., Balsamo, G., Bauer, P., Bechtold, P., Beljaars, A.

C. M., van de Berg, L., Bidlot, J., Bormann, N., Delsol, C., Dragani, R., Fuentes, M., Geer, A. J., Haimberger, L., Healy, S. B., Hersbach, H., H¨olm, E. V., Isaksen, L., K˚allberg, P., K¨ohler, M., Matricardi, M., McNally, A. P., Monge-Sanz, B. M., Morcrette, J.-J., Park, B.-K., Peubey, C., de Rosnay, P., Tavolato, C., Th´epaut, J.-N., and Vitart, F. (2011). The ERA-Interim reanalysis: configuration and performance of the data assimilation system. Quarterly Journal of the Royal Meteorological Society, 137(656):553–597.

Delmelle, P., Stix, J., Baxter, P., Garcia-Alvarez, J., and Barquero, J. (2002).

Atmospheric dispersion, environmental effects and potential health hazard associated with the low-altitude gas plume of masaya volcano, nicaragua. Bulletin of Volcanology, 64:423–434. 10.1007/s00445-002-0221-6.

Demar´ee, G. R. and Ogilvie, A. E. J. (2001). History and Climate: Memories of the Future?, chapter Bons baisers d’Islande: Climatic, environmental, and human dimensions impacts of the Lakag´ıgar eruption (1783–1784) in Iceland, pages 219–246. Plenum Publishing Co.,N.Y.

Dingwell, A. (2012). Dispersion of ash and gases from volcanic eruptions.

Introduction essay for PhD-studies at the Department of Earth Sciences, Uppsala University and the Centre for Natural Disaster Science.

Draxler, R. R. and Hess, G. D. (1998). An overview of the HYSPLIT 4 modelling system for trajectories, dispersion and deposition. Aust. Met. Mag., 47:295–308.

30

EUROCONTROL (2010). Ash-cloud of April and May 2010: Impact on Air Traffic.

Galle, B., Johansson, M., Rivera, C., Zhang, Y., Kihlman, M., Kern, C., Lehmann, T., Platt, U., Arellano, S., and Hidalgo, S. (2010). Network for Observation of Volcanic and Atmospheric Change (NOVAC)–A global network for volcanic gas monitoring: Network layout and instrument description. Journal of Geophysical Research (Atmospheres), 115:D5304.

Jacobson, M. Z. (2005). Fundamentals of Atmospheric Modeling. Cambridge University Press, 2 edition.

Jones, A., Thomson, D., Hort, M., and Devenish, B. (2007). The U.K. met office’s next-generation atmospheric dispersion model, NAME III. In Borrego, C. and Norman, A.-L., editors, Air Pollution Modeling and Its Application XVII, pages 580–589. Springer US. 10.1007/978-0-387-68854-1 62.

Knudsen, M. and Weber, S. (1911). Luftwiderstand gegen die langsame bewegung kleiner kugeln. Annalen der Physik, 341(15):981–994.

Koyaguchi, T. and Ohno, M. (2001). Reconstruction of eruption column dynamics on the basis of grain size of tephra fall deposits: 1. methods. Journal of

Geophysical Research: Solid Earth, 106(B4):6499–6512.

Leadbetter, S. J. and Hort, M. C. (2011). Volcanic ash hazard climatology for an eruption of hekla volcano, iceland. Journal of Volcanology and Geothermal Research, 199(3–4):230 – 241.

N¨aslund, E. and Thaning, L. (1991). On the Settling Velocity in a Nonstationary Atmosphere. Aerosol Science and Technology, 14(2):247–256.

Oettl, D. and Uhrner, U. (2011). Development and evaluation of GRAL-C dispersion model, a hybrid Eulerian–Lagrangian approach capturing NO–NO2–O3 chemistry.

Atmospheric Environment, 45(4):839 – 847.

Peterson, R. A. and Dean, K. G. (2008). Forecasting exposure to volcanic ash based on ash dispersion modeling. Journal of Volcanology and Geothermal Research, 170(3–4):230 – 246.

Skamarock, W. C., Klemp, J. B., Dudhia, J., Gill, D. O., Barker, D. M., Duda, M. G., Huang, X.-Y., Wang, W., and Powers, J. G. (2008). A Description of the Advanced Research WRF Version 3. NCAR technical note, National Center for Atmospheric Research, Boulder, Colorado, USA.

Sparks, R. S. J. (1986). The dimensions and dynamics of volcanic eruption columns.

Bulletin of Volcanology, 48:3–15.

Sparks, R. S. J. and Wilson, L. (1976). A model for the formation of ignimbrite by gravitational column collapse. Journal of the Geological Society, 132(4):441–451.

Steensen, T., Stuefer, M., Webley, P., Grell, G., and Freitas, S. (2013). Qualitative comparison of Mount Redoubt 2009 volcanic clouds using the PUFF and WRF-Chem dispersion models and satellite remote sensing data. Journal of Volcanology and Geothermal Research, 259(0):235 – 247. The 2009 Eruption of Redoubt Volcano, Alaska.

Stohl, A., Forster, C., Frank, A., Seibert, P., and Wotawa, G. (2005). Technical note:

the Lagrangian particle dispersion model FLEXPART version 6.2. Atmospheric Chemistry and Physics, 5(9):2461–2474.

Stohl, A., Prata, A. J., Eckhardt, S., Clarisse, L., Durant, A., Henne, S., Kristiansen, N. I., Minikin, A., Schumann, U., Seibert, P., Stebel, K., Thomas, H. E.,

Thorsteinsson, T., Tørseth, K., and Weinzierl, B. (2011). Determination of time-

31

and height-resolved volcanic ash emissions and their use for quantitative ash dispersion modeling: the 2010 Eyjafjallaj¨okull eruption. Atmospheric Chemistry and Physics, 11(9):4333–4351.

Stothers, R. B. (1996). The great dry fog of 1783. Climatic Change, 32:79–89.

Textor, C., Graf, H., Herzog, M., Oberhuber, J., Rose, W. I., and Ernst, G. (2006).

Volcanic particle aggregation in explosive eruption columns. part i:

Parameterization of the microphysics of hydrometeors and ash. Journal of Volcanology and Geothermal Research, 150(4):359 – 377.

Warner, T. T. (2011). Numerical Weather and Climate Prediction. Cambridge University Press.

Webster, H. N., Thomson, D. J., Johnson, B. T., Heard, I. P. C., Turnbull, K., Marenco, F., Kristiansen, N. I., Dorsey, J., Minikin, A., Weinzierl, B., Schumann, U., Sparks, R. S. J., Loughlin, S. C., Hort, M. C., Leadbetter, S. J., Devenish, B. J., Manning, A. J., Witham, C. S., Haywood, J. M., and Golding, B. W. (2012).

Operational prediction of ash concentrations in the distal volcanic cloud from the 2010 Eyjafjallaj¨okull eruption. Journal of Geophysical Research (Atmospheres), 117.

Wesely, M. (1989). Parameterization of surface resistances to gaseous dry deposition in regional-scale numerical models. Atmospheric Environment (1967),

23(6):1293–1304.

Wesely, M. and Hicks, B. (2000). A review of the current status of knowledge on dry deposition. Atmospheric Environment, 34(12–14):2261–2282.

WHO (2006). Air quality guidelines: global update 2005: particulate matter, ozone, nitrogen dioxide, and sulfur dioxide. World Health Organization, 3rd edition.

Witham, C. S., Hort, M. C., Potts, R., Servranckx, R., Husson, P., and Bonnardot, F.

(2007). Comparison of VAAC atmospheric dispersion models using the 1 November 2004 Grimsv¨otn eruption. Meteorological Applications, 14(1):27–38.

32