Characterising the atmospheric conditions leading to large error growth in volcanic ash cloud forecasts

Volcanic ash poses an ongoing risk to the safety of airspace worldwide. The accuracy to which we can forecast volcanic ash dispersion depends on the conditions of the atmosphere into which it is emitted. In this paper we use mete-orological ensemble forecasts to drive a volcanic ash transport and dispersion model for the 2010 Eyjafjallajokull eruption. From analysis of these simulations we determine why the skill of deterministic-meteorological forecasts decrease with increasing ash residence time, and identify the atmospheric conditions in which this drop in skill occurs most rapidly. Large forecast errors are more likely when ash particles encounter regions of large horizontal ﬂow separation in the atmosphere. Nearby ash particle trajectories can rapidly di-verge leading to a reduction in the forecast accuracy of deterministic forecasts which do not represent variability in wind ﬁelds at the synoptic-scale. The ﬂow separation diagnostic identiﬁes where and why large ensemble spread may occur. This diagnostic can be used to alert forecasters to situations in which the ensemble mean is not representative of the individual ensemble member volcanic ash distributions. Knowledge of potential ensemble outliers can be used to assess conﬁdence in the forecast and to avoid potentially dangerous situations in which forecasts fail to predict harmful levels of volcanic ash.


Introduction
Volcanic ash poses a significant hazard to aircraft.It can cause both temporary engine failure and permanent engine damage (Guffanti et al. 2005).Flights are therefore restricted in ash contaminated airspace, which disrupts air traffic leading to the potential for large financial loses.For example the 2010 Eyjafjallajökull eruption grounded over 95,000 flights, costing the airline industry over 1 billion pounds.Analysis of the 1900-2010 Icelandic historical records shows that a volcanic eruption of the size of the 2010 Eyjafjallajökull eruption has a repeat period of between 5 and 10 years (Thordarson and Larson 2007).Worldwide, volcanic eruptions 10 times the size of the 2010 Eyjafjallajökull eruption have repeat periods on a decadal timescale (e.g.Mount St Helens 1980, Hudson 1991, Puyehue 2011).Given the ongoing risk of volcanic eruptions it is important to continually evaluate and improve the accuracy of volcanic ash forecasts to ensure safe and optimised flight operations during future volcanic eruptions.
The volcanic ash advisory centres (VAACs) are responsible for producing volcanic ash cloud analysis and forecasts to assist the aviation community in planning their operations and minimising risks.There are currently 9 VAACs that together provide a comprehensive global modelling and warning system for the aviation community.These 9 VAACs use 6 different volcanic ash transport and dispersion (VATD) models to to produce volcanic ash charts showing the forecast location of volcanic ash in the atmosphere at different flight levels and out to forecast lead-times of 24 hours.VATD models are initialised using data about the location of the eruption, the time at which the eruption started and, if available, information about the plume rise height, vertical profile of volcanic ash and ash size distribution (known collectively as eruption source parameters, ESPs).They also use 3-D winds as input from numerical weather predictions to transport volcanic ash away from the source.Typically the meteorological input used has a horizontal resolution of between 10 and 50km.To represent dispersion on scales smaller than this horizontal diffusion is applied.The diffusion represents the dispersion by unresolved eddies and acts to increase the vertical and lateral spread of volcanic ash clouds (Dacre et al. 2015).This approach assumes that the small-scale dispersion processes are of an eddy viscosity character and thus can be represented using a Gaussian description (Pasquill and Smith 1983).The simulated ash cloud therefore represents the time mean of an ensemble of realisations.
At larger scales however, individual realisations can often display considerable deviations from the ensemble mean (Mylne and Mason 1991).The scale at which this occurs depends on the size of the dispersion processes relative to the width of the time averaged ash cloud.For averaging periods of a few hours, this scale is typically greater than 500 km.Variability on synoptic scales however differs for different atmospheric circulation patterns, meaning that the traditional Gaussian diffusion approach used for small-scale dispersion processes cannot be used.Current operational VATD models do not represent variability at the synoptic scale.They use meteorological input from a single realisation of the flow field to produce a volcanic ash forecast (referred to as deterministic-met volcanic ash forecast in this paper).The aim of this paper is to identify the atmospheric conditions in which there is a higher chance that deterministic-met volcanic ash forecast skill may rapidly decrease and to discuss the potential use of ensemble meteorological input to VATD models as a method to address the missing synoptic-scale variability in volcanic ash forecasts (referred to as ensemble-met volcanic ash forecasts in this paper).
Several studies have investigated the space and time-dependent skill of deterministic-met volcanic ash forecasts.For example, Stunder et al. (2007) analysed the forecast skill for 7 different volcanic eruptions by comparing deterministic-met volcanic ash forecasts with satellite observations.They showed that these forecasts were generally good for short-term (18 hours from start of the eruption) forecasts but that forecast skill appeared to decrease at longer lead-times.This relationship between volcanic ash forecast skill and forecast lead-time is due to (i) increasing errors in the forecast wind fields and ESPs at longer forecast lead-times and (ii) longer lead-time forecasts include particles with longer residence times.These particles experience an accumulation of errors in the wind field leading to larger positional errors on average than particles with shorter residence times.Dacre et al. (2016) examined the second of these sources of error by performing hindcast simulations of the Eyjafjallajökull eruption (using analysis wind fields).They showed that generally skill decreases as the residence time of ash increases but that the rate of skill decrease depends on the meteorological situation.In some situations only the position of ash particles with residence time less than 24 hours are correctly simulated whereas in other situations the position of ash particles with residence times longer than 72 hours can be accurately simulated.
Other studies have shown that the inclusion of buffer zones, to account for positional errors in the deterministic-met volcanic ash clouds, can lead to significant improvement in the agreement with observations (Webster et al. 2012;Grant et al. 2012).These buffer zones are a simplistic attempt to account for uncertainty in the synoptic-scale wind fields.
For some time the use of ensemble-met volcanic ash forecasts has been advocated by the wider volcanic ash community (Bonadonna et al. 2012) as a more rigorous way of accounting for uncertainty in large-scale wind field.Stefanescu et al. (2014) and Madankan et al. (2014) include both ensemble meteorology and an ensemble of ESPs in their study to quantify overall uncertainty in volcanic ash forecasts.They demonstrate that the range of predicted concentrations can be large at forecast lead-times of 48 hours.Similarly Vogel et al. (2014) performed time-lagged ensemble simulations of volcanic ash dispersion from the Eyjafjallajökull plume and found that for some times the spread in ensemble-met forecasts is small but at others it is large.They attribute this to the nonlinear behaviour of the atmosphere.Dare et al. (2016) performed a comparison of both deterministic and ensemble-met volcanic ash forecasts for the 2014 Kelut eruption.They found that both showed good agreement with satellite observations for the first 12 hours from the start of the eruption.However, for longer lead-times (18-24 hours) the ensemble-met forecast showed better agreement with observations than the deterministic-met forecast.
While all these studies demonstrate that ensemble-met forecasts show better agreement with observations than the deterministic-met forecasts, particularly at longer lead-times, the dynamical reasons why they perform better has not been explored.The aim of our study therefore is to illustrate why the skill of deterministic-met forecasts decreases with increasing ash residence time, and furthermore to identify the atmospheric conditions in which this drop in skill occurs most rapidly.These conditions are identified using ECMWF meteorological ensembles as input to the NAME VATD model to simulate an ensemble of particle trajectories.

Methodology a. Meteorological fields
In order to determine the uncertainty associated with the synoptic scale meteorological flow field an ensemble of meteorological forecasts are used.Each forecast is produced from perturbed initial conditions that represent the likely initial analysis error distribution.In this paper the European Centre for Medium Range Weather Forecasting (ECMWF) Integrated Forecasting System (cycle 41r1) has been used to create bespoke ensemble forecasts of the meteorological conditions during the 2010 eruption of Eyjafjallajökull.Global forecasts are initialised every 12 hours between 00 UTC on 1 May -12 UTC on 8 May 2010.Each forecast is 42 hours long and has 20 ensemble members.Data is archived every 6 hours on 26 levels and at T639 spectral truncation (approximately 32km horizontal grid spacing).Initial perturbations are constructed using the singular-vector approach (Buizza and Palmer 1995) and model uncertainty is taken into account through the use of a simple stochastic physics scheme (Buizza et al. 1995).Data is extracted from the ECMWF archive at 0.25 • × 0.25 • on a regular lat/lon grid and several fields (surface stresses, sensible heat flux and precipitation fields) are post-processed as data extracted from the ECMWF archive cannot be used directly as input to the VATD model described in section b.

b. NAME dispersion simulations
The VATD model used in this study is the Numerical Atmospheric-dispersion Modelling Environment (NAME).NAME is used by the London Volcanic Ash Advisory centre to forecast the spatial distribution of volcanic ash following an eruption.In this study we use NAME III (version 6.3) and ECMWF numerical weather prediction meteorological data to disperse particles released into the atmosphere at the position of the Eyjafjallajökull volcano in Iceland.The dispersion of volcanic ash by small-scale three-dimensional atmospheric turbulence and unresolved mesoscale motions are parametrized within NAME using random-walk techniques.The aim of the randomwalk dispersion is to compute an ensemble of random trajectories of Lagrangian particles through a flow field whose statistics are based on observations of vertical and horizontal velocity variances and diffusivities (Thomson and Wilson 2013).The position of the particles is tracked for 42 hours to create particle trajectories.The volcanic ash density is assumed to be 2300 kg m −3 based on the value used in the operational version of NAME (Leadbetter and Hort 2011) and the particle size is assumed to be 2 µm based on in-situ observations of the ash cloud by the FAAM aircraft over and around the UK in the Eyjafjallajokull ash cloud (Johnson et al. 2012).Particles are subject to removal processes including sedimentation, wet and dry deposition (Jones et al. 2007).Note that the choice of particle size does not affect the conclusions reached in the paper.

c. SEVIRI satellite observations
To qualitatively evaluate the performance of the NAME forecasts we compare simulated ash cloud distributions with data from the Spinning Enhanced Visible and Infrared imager (SEVIRI).SEVIRI volcanic ash retrievals are calculated using brightness temperature difference measurements (see Francis et al. (2012) for more details).The advantage of using volcanic ash retrievals from an instrument onboard a geostationary satellite is that they are available at high temporal resolution, every hour, allowing us to track the evolution of the volcanic ash cloud and to interpolate between timesteps when water or ice clouds obscure the volcanic ash.Following the method of Harvey and Dacre (2016) we composite satellite observations over a 5 hour time window.This has been shown to be sufficient to create a continuous time series while remaining highly correlated with the noncompostied fields.The satellite volcanic ash retrievals are averaged onto a 0.5 • × 0.5 • latitude/longitude grid to allow direct comparison with the NAME output.

d. Ensemble spread and flow separation diagnostics
One measure of the uncertainty in meteorological flow conditions is the time evolution of spatial spread in particle trajectories.In this paper the ensemble spread is calculated using the rootmean-square (RMS) distance between individual ensemble particle positions (1 particle from each ensemble simulation), (x i ), and the mean position of the particles, (x i ), summed over all N particles (thus N equals 20 as there are 20 ensemble simulations).The distance is measured perpendicular to the mean direction travelled by the particles during the previous 10 minutes to capture lateral spreading of the trajectories only.

RMS
The diagnostic used to characterise the synoptic-scale flow conditions is the 2-D horizontal flow separation diagnostic introduced in Dacre et al. (2016).This flow separation is calculated as the velocity gradient perpendicular to the flow.
where v is the velocity vector, q is the wind speed, n is distance in the direction perpendicular to the flow, and x and y are distances in longitude and latitude directions, respectively.Where this diagnostic is positive, the atmospheric flow separates, and where it is negative, the flow contracts.
Thus it is a good diagnostic for identifying where particle trajectories will spread apart.The flow separation diagnostic is related to the 3-D Lyapunov exponents used by Legras et al. (2005) and Pisso and Legras (2008) to characterise the rate of separation of infinitesimally close trajectories in phase space.and (f) the majority of the volcanic ash is transported cyclonically and is advected towards Europe.

a. Satellite-detected ash clouds
In contrast in figures 1(c) and (d) the majority of the volcanic ash cloud continues to travel anticyclonically and is advected into the North Atlantic.For this example, the deterministic-met forecast shown in figures 1(c) and (d) would be considered a good forecast as it closely matches the evolution of the ash cloud seen in the satellite observations.However the deterministic-met forecast shown in figures 1(e) and (f) would be considered a poor forecast as it does not forecast the observed ash in the North Atlantic.This is despite both forecasts using equally plausible realisations of the flow field.This example highlights the danger of using a single deterministic-met flow field as input to a VATD model to forecast the ash cloud distribution.These 2 ensemble members are chosen because they exhibit very different volcanic ash cloud evolutions, the other 18 ensemble members result in ash distributions which resemble a mixture of the two extremes.

c. Flow separation
In this section we explain why the ensemble member forecasts differ so much from each other.
In order to do this we examine the flow pattern at approximately 50 • N and 15 • W, the location at which the ash particle trajectories show an increase in spread.In order to isolate transport by the resolved-scale flow it is not subject to perturbations representing unresolved eddies, hence its smooth trajectory.
The black star indicates the location of the particle at the time of the flow separation field.12 hours after the particle is released into the atmosphere (figure 2(a)) the particle is at 57 • N, 13 • W where the streamlines are roughly parallel to one another and hence flow separation is small.24 hours after the particle is released into the atmosphere (figure 2(b)) the particle is at 51 • N, 17 • W and is in a region of strong positive flow separation.The streamlines spread apart as they approach the point of intersection between the trough and ridge region (known as a col or saddle point).
It is difficult to analyse the along-trajectory flow separation in this Eulerian framework, therefore figure 2(c) shows the flow separation extracted at the relevant time along the Lagrangian particle trajectory.This Lagrangian analysis demonstrates that the particle advected in this deterministicmet forecast enters a region of strong flow separation at 52 • N, 17 • W. In order to determine whether this is specific to a single ensemble-met member forecast or to all of the meteorological ensemble forecasts initialised at 06UTC on 6 May the along-trajectory flow separation has been calculated for each of the meteorological ensemble forecast members.

d. Trajectory spread
To establish if trajectory spreading always rapidly increases after trajectories encounter regions of positive flow separation similar experiments were performed for meteorological ensemble forecasts initialised at 06UTC on the 15 April -7 May 2010.For each of these ensemble forecasts a single particle were released at a height corresponding to the observed plume top from the Eyjafjallajökull volcano.It is well known that in low wind-speed conditions wind direction can vary significantly in a short period of time causing particle trajectories can rapidly diverge (Venkatram et al. 2004).In this paper we choose to focus on the less well studied uncertainty occurring in moderate-strong wind conditions and thus only analyse the situations in which the wind speed at the release height was greater than 10m s −1 .Figure 3 shows the ensemble-met member forecasts

Discussion and Conclusions
In this paper we examine the atmospheric flow characteristics that lead to volcanic ash cloud bifurcation and a reduction in forecast skill.We performed multiple forecasts using the UK Met Office volcanic ash transport and dispersion model (NAME) and input from ensemble meteorological flow fields from the ECMWF ensemble prediction system.
In moderate to strong wind situations the atmospheric conditions leading to large variability in volcanic ash particle positions are associated with large flow separation.When ash particles encounter regions of large horizontal flow separation their future trajectories are very sensitive to their position at that time.Nearby ash particle trajectories can rapidly diverge leading to a reduction in forecast accuracy for deterministic-met volcanic ash forecasts.Potentially leading to predictions of ash-free airspace in regions that are in-reality contaminated with ash or vice versa.
In order to fully represent the synoptic-scale meteorological uncertainty ensemble-met volcanic ash forecasts are needed.When volcanic ash clouds encounter regions of large flow separation the individual ensemble-met members may display considerable deviations from the ensemble mean.
2-D fields of positive flow separation could be used as a flag to alert forecasters to this potential risk and the individual ensemble-met member forecasts analysed.A combination of the flow separation diagnostic and ensemble volcanic ash forecasts will help to identify where and why large uncertainty in the forecast occurs and provide an estimate of the confidence of the forecast.
For example, a forecaster could reduce the size of the hazardous area whenever high confidence in the ash cloud forecast was indicated.Reductions in the hazard area would avoid unnecessary disruption to airspace.
In this paper we have only considered the uncertainty in the horizontal wind fields.Uncertainty also exists in the magnitude and location of precipitation which leads to wet-deposition of volcanic ash.This uncertainty may also cause large errors in the magnitude of volcanic ash forecasts as precipitation is a very efficient removal mechanism.We have also not considered the uncertainty associated with the volcanic eruption source parameters (ESPs).The best way to combine the meteorological and ESP uncertainty and effective ways of communicating this uncertainty with users is the subject of future work.

Figure 1
Figure 1(a) and (b) show the ash cloud from the Eyjafjallajökull eruption, as detected by the

Figures 2
Figures 2(a) and (b) shows the streamlines and flow separation at 12UTC and 18UTC on 6 May Figure 2(d) shows the evolution of flow separation along 20 particle trajectories released at the same time, a single particle trajectory in each ensemble-met forecast.The flow separation in each ensemble-met forecast is very similar up until the point at which the trajectories start to diverge.This is expected since the regions of positive and negative flow separation are spatially coherent.It also illustrates how the trajectory separation rapidly increases after the point at which the trajectories encounter the region of positive flow separation.Performing an ensemble-met volcanic ash forecast for this case accounts for the variability in the synoptic flow field and is necessary to fully encompass the ash cloud distribution uncertainty due to the flow field.
with the 4 highest (figures 3(a)-(d)) and 4 lowest (figures 3(e)-(h)) trajectory spreads.Individual particle trajectories correspond to a single particle released at the same time in each ensemble-met member forecast.It can be seen that on some days, figures 3(a)-(d), the trajectories diverge after encountering regions of positive flow separation, consistent with the analysis for the 6 May 2010 (figures 2(d)).By comparison on other days, figures 3(e)-(h) the trajectories remain close together for 42 hours.

Figure 4
Figure 4 quantitatively describes the relationship between residence time and ensemble-met FIG. 2. Flow separation averaged from 325-225hPa (filled contours) overlaid with 275hPa streamlines (grey) at (a) 18TUC on 6 May 2010, and (b) 06 UTC on 7 May 2010.42hr particle trajectory initialised at 06 UTC on 6 May 2010 (thick black line) and position of particle at time of flow separation and streamline fields (black star).(c) Flow separation along the particle trajectory shown in (a) and (b) from 12 hrs residence time onwards.(d) Flow separation along 20 particle trajectories advected by 20 different forecast wind fields.