Abstract

Volcanic ash concentrations in the plume from Sakurajima volcano in Japan are observed from airplanes equipped with optical particle counters and GPS tracking devices. The volcano emits several puffs a day. The puffs are also recorded by the Sakurajima Volcanological Observatory. High concentrations are observed in the puffs and fallout driven by vertical air current, called streak fallout. Puffs dispersion is analyzed by the classical diffusion-advection method and a new gravitational dispersion method. The fluid mechanic of the gravitational dispersion, streak fallout, and classical diffusion-advection theory is described in three separate appendices together with methods to find the time gravitational dispersion constant and the diffusion coefficient from satellite photos. The diffusion-advection equation may be used to scale volcanic eruptions so the same eruption plumes can be scaled to constant flux and wind conditions or two eruptions can be scaled to each other. The dispersion analyses show that dispersion of volcanic plumes does not follow either theories completely. It is most likely diffusion in the interface of the plume and the ambient air, together with gravitational flattening of the plumes core. This means larger boundary concentration gradients and smaller diffusion coefficients than state of the art methods can predict.

1. Introduction

Airborne observations of volcanic ash concentrations with optical particle counters and GPS tracking have recently been taken in use to study volcanic plumes [17].

Three campaigns of airborne observations of volcanic ash and gas content of the plume from the Sakurajima volcano in Japan were performed in 2013, in cooperation between the Universities of Kyoto, University of Iceland, and University of Applied Sciences in Düsseldorf, Germany. The Sakurajima volcano has been in constant eruption since 1955; it does not emit a continuous plume but produces a series of explosion puffs, sometimes many each day [8]. The preliminary results of the first campaign were described in presentations to the IAVCEI conference July 2013 in Kagoshima, Japan [912].

The staggering economic disaster [13, 14], inflicted upon the aviation industry by the Eyjafjallajökull eruption in 2010 and the Grimsvötn eruption in 2011, made clear the importance of ash cloud predictions. These events sparked the research program that produced the Sakurajima campaigns. It is hoped they may help to increase the accuracy of ash cloud predictions.

Ash cloud predictions make use of point source atmospheric dispersion models, for the simulation of continuous, horizontal and neutrally buoyant volcanic plumes. There exist a number of atmospheric dispersion models [15], specially constructed for handling eruption problems both analytically and numerically. In [16] the VAAC Name Model is presented; it is based on the horizontal advection-diffusion equation. The abstract states that they “provide advice on extent of ash” unfortunately the extent was sometimes 40 times too great as shown later in [3]. In this model and most others, vertical diffusion is assumed zero as is discussed in [17, 18]. This is a grave simplification, turbulence is always three-dimensional, but it may be a good approximation, if there is a strong density discontinuity both on the upper side and the underside of the plume. Negative density jumps in the upwards vertical direction make the interface very stable and curb the mixing process, but they also give reason to investigate possible gravitational flattening of the plume. The small vertical mixing approximation seems to be used by the majority of researchers. Suzuki uses it in his much cited publication, [19], and there seems to be a general agreement among researchers that the problem of neutrally buoyant volcanic plumes can be treated as horizontal diffusion and many atmospheric dispersion models that use this approach are now available on the internet [15]. In recent times 3D models are coming into use and they also use the classical diffusion-advection theory.

The most commonly applied method to determine the point source strength in volcanic ash plume simulations is to calculate the eruption output from the plume height. Then an equation originally due to Wilson et al. [20] and Settle [21] is used. It is very simple and most researchers seem to be of the meaning that it can be inaccurate in a variety of situations. Several modifications do exist; one is mentioned in Appendix C on diffusion theory and the Eyjafjallajökull eruption in 2010.

The fluid dynamics in volcanic plumes have been studied by many researchers [2227], but the search for a more accurate method for a plume height-mass output relation than Wilsons is still on. Numerical studies of buoyant plumes by [16, 18, 28] show that simple forms of the diffusion-advection equation can give acceptable results.

In general this paper is about that airborne measurements of volcanic ash concentrations in a plume can be used as a tool in both disaster prevention and volcanology research. It is clear from the Sakurajima results that the fluid dynamics of stratified flows are important as an analysis tool; so two appendices on gravitational effects on dispersion of volcanic plumes are included. Appendix A is on gravitational flattening of plumes, a process that can easily be mistaken for diffusion, Appendix B is on streak fallout, a fallout process many times more effective than fallout by way of terminal velocity of individual grains. These results can be very different from what the classical diffusion-advection theory, summarized in Appendix C, predicts.

2. Airborne Measurement of the Sakurajima Plume in January 2013

2.1. The Measurement Campaign

In 2013 Kyoto University undertook three airborne measurement campaigns in cooperation with the University of Iceland and University of Applied Sciences, Düsseldorf, Germany. The instruments used where two DustMate industrial OPCs, two OPCs from GRIMM for scientific measurements, and a DOAS (differential optical absorption spectrometry) system for detection of sulfur dioxide. The planning, organization, quality assurance, and overview of results from the first campaign on January 15, 16, and 18 are discussed in presentations in the IAVCEI 2013 conference in Kagoshima (IAVCEI, 2013). In this paper we will discuss the observation results measured 15 January 2013 and 27 November 2013.

2.2. Data for the Eruption in Sakurajima 15 January 2013

The Sakurajima Volcanological Observatory (SVO) has provided minute-to-minute data of Sakurajima explosions during the campaign. In Figure 1 the emission data from SVO is shown. Counting consecutive explosions as one event, Figure 1 contains 5 big events, the first 9:01. The puffs shown in Figure 1 can be identified; identification data for the 9:50 event (puff number 1) is shown in Table 1.

By downwind tracking, using wind data from Japan Meteorological Agency (JMA), three puffs occurring 15 January 2013 09:50–9:53 are found to be the source of the tops in Figure 2.

The track numbers shown in the left side of Figure 2 refer to the numbered tops in the right side of the picture. The OPC results are filtered with 16 convolutions of a (1/4, 1/2, 1/4) filter kernel to remove turbulent fluctuations in the raw data (Figure 4) here called the F-16 filter. The filter kernel is Gaussian; so the filtered values have concentrations gradients that fit the advection-diffusion model.

Figure 3 shows the big explosion puff at 9:53 o’clock. It is an almost a round ball of ash followed by a tail of gas. The ball is about 1100 meters in diameter.

3. Diffusion Case

3.1. Puff Identification and Diffusion Parameter Estimation

The tracks 1, 2, and 3 in Figure 2 penetrate the big puff in Figure 3 and track 4 penetrates a smaller puff or the tail. The Gaussian distribution Appendix C can be used to find the diffusion coefficient .

Table 2 shows the properties estimated for puff number 1 from the concentrations and point locations in Figure 2. The Gaussian length scale (the constant in Table 2) and the centerline concentration are estimated by maximizing the correlation coefficient (98%) of the theoretical curve (equations (C.1) and (C.13)) and the measurements of the puff number 1 in Figure 2 and minimizing the standard error.

The total radius is estimated from the nearest background concentration value, visible radius is out to the 2000 μgr/m3 value, and effective radius evenly distributed total mass of 363 tons with the center concentration throughout. Here it must be noted that the radius of the original puff in Figure 3 is only 600 m so the center concentration there is 26 times Table 2 value, according to the advection-diffusion model.

The measurements give PM1, PM2.5, PM10, and TSP, and these are used to estimate the grain size distribution given. Variations in grain size distribution inside the puffs are very small so average values only are given. Almost 90% of the dust is aerosol size, that is, grain size below PM10. This is the dust dangerous to jet airliners because it is carried over very large distances due to its low settling velocity. Comparing the total puff mass 363 tons to the SVO estimate we see that aerosol size grains make up 21% of the total erupted mass (1697 tons); this is far higher than in stronger eruptions. If the measured results are scaled up, according to the scaling in Appendix C, to a 4500 m high plume in a 12 m/s wind the results may be seen in Table 4. The scaled puff mass is about 18000 tons erupted in eight minutes. This may be compared to the Eyjafjallajökull 2010 eruption where eruption output was generally 6000–60000 tons/min [29].

Resulting value is shown in Table 2. It is quite lower than in stronger eruptions. When the value is scaled up to the altitude of the Eyjafjallajökull plume discussed earlier, we get a value of  m2/sec. This is about 2/3 of the Eyjafjallajökull value (Table 6 and Appendix C).

The calculated size of the puff is very big. The visible and total radius values show that it reaches all the way to the ground. The distance across the puff from visible concentration to visible concentration in Figure 2 is 3 km.

3.2. The Plume as a Series of Puffs

A volcanic plume may be treated as a series of individual puffs, where each puff is an ash cloud that results from one single explosion. This is the approximation that is used when the differential equation of the advection-diffusion process is derived (Appendix C (C.1)). As the concentration of the volcanic ash does not affect the wind velocity or the mixing properties, the concentration in a given place at any given time, is equal to the sum of the concentrations from all the individual puffs. This method will not be discussed further in this paper.

4. Application of Gravitational Dispersion

Puff identification data is the same as in Figure 2 in chapter 2. However, the filtering needs to be done separately for the ambient air and the ash contaminated air in the puff, as the F-16 filter is a Gaussian filter and smoothens non-Gaussian distributions. We shall apply a F-16X filter that supposes no, or very low, mixing across the concentration discontinuities in Figure 2 by introducing a fixed concentration in a point between the two air masses. The result is shown in Figure 4.

The different filter does not change the overall picture except for the boundary gradients. The center concentration is unchanged. Initial radius is measured from Figure 3, final size from Figure 4 and (A.8) is used for the deformation together with the data in Table 1.

The result is dramatically different, especially the ash mass in the puff. The two last lines in Table 3 show diffusion velocity and puff radius growth according to the measured concentration gradients.

5. Discussion

The scaling of the Eyjafjallajökull plume in Table 2 looks to be quite successful. Using the model laws for continuous plumes gives satisfactory results, when a big 3-minute plume from Sakurajima (as model) is scaled up to Eyjafjallajökull 2010 (as prototype). In disaster prevention research, scaling can be used to scale plumes, measured or simulated, up to scenarios for hazardous events that can be a help for civil protection authorities in hazard assessment.

As a by-product, the self-similarity relations for ash fallout may prove a useful tool. Figure 13 looks promising for a situation with stable weather when no scaling is necessary. When scaling is necessary, as in Figure 14, unscaled values cannot be supposed to show any similarity at all. In the self-similarity relations lies a possibility to improve total estimates for ash output.

However, the advection-diffusion equation (C.1) has to be valid. The diffusion results in Section 3.1 give credible results; the main discrepancy is that the visible puff grows very big and it is not possible to explain the value of the resulting diffusion coefficient .

The gravitational flattening in chapter 4 gives different results. Puff mass is only a fraction of the diffusion estimate. When the boundary is no longer adapted to diffusion theory, the observed sharp boundary and the initial center concentration value are maintained; so this model produces a much smaller puff mass. In the diffusion case; the center concentration is diluted about 26 times during the almost half an hour and 23 kilometers long downwind migration. This difference between the two methods will always exist. To use diffusion models to simulate an observed plume, it simply demands more ash than is the case for gravitational flattening models.

The temperature difference in an inversion big enough to facilitate the gravitational dispersion is low, indeed hardly measureable. The time constant of the Sakurajima plume spreading in Table 3 compares to the Eyjafjallajökull values in Table 4; so dispersion of the puff in Figure 2 is on the same time scale as the Eyjafjallajökull plume. But assuming diffusion only, the puff also scales to the Eyjafjallajökull plume according to the diffusion-advection scaling rules in Appendix C. Working with horizontal dispersion only, gravitational flattening can therefore easily be mistaken for horizontal diffusion and vice versa.

All this points to that it is necessary to merge the two models; that is, the dispersion is not either diffusion or gravitational but both. To assume zero vertical diffusion together with full horizontal diffusion does not work properly. Assuming zero vertical diffusion, Table 2 mass estimate would go down to 120 tons, but that is still too high.

The boundary data in Table 3 is very interesting. The puff boundary gradients in Figure 4 compare much better to the actual measured gradients (the raw data) than the gradients in Figure 2, as a Gaussian filtering across boundaries flattens the concentrations gradients out in order to produce an optimal fit to the diffusion theory. X-filtering, on the other hand, supposes the measured gradients to represent the average gradient across the boundary. The three last lines in Table 3 show the speed of the diffusion process, if the value was as in Table 2 and the gradients as in Figure 4. This shows the outward speed of the high concentrations for this value of in Table 3 second line from bottom and that the boundary would move 1476 m in all directions in the downwind time and the puff become almost 4 km in diameter. values thus have to be very low, only 10% or so of the diffusion vale, if gravitational dispersion is active.

Some of the difference in puff size estimations in Tables 2 and 3 can possibly be due to fallout. Streak fallout especially like the one observed in July 27, 2013, can explain some of the differences but not all. It cannot be assumed that a puff is sending out streak fallouts until it is under 10% of original size only 20 kilometers from the source. Such process would have been detected and described by the volcanologists. But the diffusion theory gives that the concentration has been 26 times higher than Table 2 value in the newly formed puff Figure 3. This would mean a many times higher density than in the streak fallout in Appendix B and would certainly start large streak fallouts, but not 90% of the plume disappear in unobserved streaks.

If both dispersion methods are active, the plumes core undergoes gravitational deformation, and the concentrations gradients at the boundary produce a diffusion envelope around it. The fluid mechanics of this process demand much more complicated flow model than the simple Bernoulli approximation in Appendix A and will not be attempted in this paper. Besides, the inversion temperature difference in Table 3 is estimated; we would need an observed value in a combined model. But to measure such a small temperature difference is very difficult.

6. Conclusions

Airborne measurements in the plume from the volcano Sakurajima in Japan show very good results that can be used to find the properties of volcanic puffs. To model the dispersion, two methods have to be considered, advection-diffusion method and gravitational dispersion. The fallout can be streak fallout due to vertical gravity currents.

Gravitational Dispersion. This is a new method, only the in Appendix A approximate theory exists. Modeling and dispersion prediction according to this method need very accurate temperature data from the plume center, normally not available. This method explains why there is little or no vertical dispersion but plumes and puffs flatten out horizontally due to density currents instead. It also explains how plume boundaries with large concentration gradients can spread horizontally without any diffusion. Gravitational dispersion is a nonlinear process and scaling is not possible. The approximate theory assumes that plume (or puff) center concentrations are preserved which leads to much smaller estimates for the erupted mass than diffusion models. However, diffusion cannot be totally absent. Over large distances these two models have to be combined.

Streak Fallout. This is a new fallout model; the theory for it is in Appendix B. Streak fallouts carry large quantities of ash to the ground by vertical gravity currents containing all grain size fractions. These currents may be considered as chunks of the ash cloud that fall to the ground with higher velocity than the terminal fallout velocity of the grains in still air. The vertical gravity currents have similarities to microburst and downburst winds but do not reach as high wind velocities as they do. Streak fallouts can deplete the mass in the plume much faster than ordinary fallout can.

Advection-Diffusion Method. Sakurajima eruptions scale to larger eruptions in a convincing manner when this method is used. The scaling also produces self-similarity rules for ordinary fallout that may prove useful in fallout studies (Appendix C). However, the diffusion coefficients that have to be used in the modeling to explain observed dispersion of the plumes are too big. This leads to an overestimation of the erupted mass, especially PM10 and smaller. The commonly used approximation of zero vertical diffusion keeps the overestimation down, but it is still there.

Satellite pictures can also be used to estimate diffusion coefficients by tracking the visible boundary of the plume, as demonstrated by using the data for the Eyjafjallajökull eruption in 2010 in Appendix C. Such estimates of the value are presumably to high as the same data fits very well to the gravitational theory as shown in Appendix A. Center concentrations in the plume will be unaffected by horizontal diffusion for a long time, but fallout, especially streak fallout, will deplete the total mass in the plume.

The dispersion of volcanic plumes is advection-diffusion combined with gravitational flattening. The gravitational flattening is in the plumes core, the diffusion in an outside envelope. Considerable research may be needed in order to establish the new dispersion theory. New data like temperature may be needed, and one second sampling frequency in the OPC measurements is certainly a help.

When diffusion models are used and gravitational flattening is entirely left out, the values have to be very high. Taking the visible limit as 2000 μgr/m3 the ash outside these limits will become a larger and larger part of the total flux as the concentration gradients grow smaller in the downwind direction. This makes the simulated plumes an order of magnitude too wide. The gravitational effects need to be included in ash cloud predictions.

Appendices

A. Gravitational Dispersion of Plumes

A.1. Simplified Model of the Fluid Mechanics of a Plume in a Stable Atmospheric Stratification

Figure 5 shows a density stratification in the atmosphere with a volcanic plume drifting along with the wind velocity (coming out of the layer of the paper). For simplification, the density of the plume is assumed to be the average of the densities of the lower layer and the upper layer . This will keep the plume buoyant, floating half-submerged in the heavier air. Then the plume can be assumed to be symmetric. The plume does not have to be this perfect in shape or composition, but it simplifies the mathematical problem without too much loss of generality. In treating this problem we can let the density difference out everywhere, except in the gravity term (Boussinesq approximation).

Figure 5 shows the static pressures, , inside the plume, that is, pressure as it would be when horizontal velocities are zero on the average in diffusion-advection theory. This static pressure distribution means total hydrostatic balance in any vertical while there will be a negative pressure gradient in the horizontal direction, and this means flow away from the horizontal symmetry line in Figure 5. This means that the plume expands in the horizontal direction.

The horizontal expansion velocity will increase from zero in the centerline to full value in the ends. There are some flow resistance terms due to entrainment of the cold outside air into the plume, diffusion of ash through the interface and there can also be turbulent shear stress and pressure resistance in the interface, at least in theory, but these will be let out for a moment.

The horizontal outwards flow flow must satisfy the continuity equation. Therefore it can be modeled by a stream function: The velocity will increase monotonically from the center in point 0 to full value in point 2. The mathematical representation for the stream function (A.1) in its simplest form is as follows: Equation (A.2) gives the velocities and in the boundary points and , respectively. The is a local time scale, , and vary with time but not independently and cannot vary with or .

If diffusion is excluded for a moment, the cross-sectional area of the plume in Figure 5 is constant as it flattens out. Using the ellipse as an approximation for the plume it gives us the condition = constant as long as there is no entrainment. This is a reasonable assumption if the plume is buoyant in the wind and migrating with the wind velocity , without any velocity gradients acting on it.

In point 0 (Figure 5) there is no velocity, vertical or horizontal; so the easiest way to find is the Bernoullis equation along the streamline 1–0–2 Here is the average density and later will be used. If , there would be a local static overpressure of in point 0 ( acceleration of gravity). In (A.3) and (A.4), is this pressure somewhat modified by the flow but has the same value in both equations and can be eliminated. The outside pressure difference so the following simple differential equation system can be found to determine The is a convenient length scale and is time. may be interpreted as the radius of the plume in the beginning. Equation (A.5) may be solved for and inserted for , and is a downwind coordinate. The result will be a nonlinear ordinary differential equation for is the time constant of the plume spreading, it is different from . Here, , where would to be located if the plume does reach that far back. Another initial condition may be used if it is introduced in (A.5) and (A.6). If is assumed, (A.6) may be approximated by Equation (A.7) has the surprising property that the relation has only two parameters, and .

Equation (A.7) is derived using Bernoulli’s equation that assumes no flow resistance. It is therefore necessary to introduce a correction factor into the equation in order to be able to compare it with field data.

If the ash cloud is not a continuous plume, but an isolated puff, (A.7) takes the following form:

A.2. Comparison with Field Data

In [30] Andradóttir et al. use the diffusion equation to analyze the spreading of the Eyjafjallajökull plume 2010 on 5 different dates. Table 4 is compiled from their wind and temperature data using (A.7) on the boundaries of the visible plume as seen in satellite photos, instead of the diffusion theory.

Both correlation and the root mean square error (Rms) are very satisfactory. However, the results are biased, positive errors in the middle but negative in the ends indicating a slow reduction in the , or a smaller value in the last 40 km of the path than in the first. In using (A.7), and can be assumed piecewise constant along the plume, this would increase the correlation coefficient and reduce the bias, but the only available temperature profile data is from Keflavik International Airport 200 km to the west. The data shows that inversions do exist in the approximate level of the plume, but there is no data about changes in the properties of the inversion as no other temperature profile data is available.

Wind shear produces diffusion; Figure 6 shows an example.

The transparent plume is diffusion from the over- and underside of the plume in the crosswind. In the middle is the main plume. Fitting (A.7) to it suggests a temperature inversion of 2°C if a value of 0,1 is used.

A.3. Discussion

Plumes riding in stable temperature stratification will have a tendency to spread out like an oil slick on water because of the gravitational effects in the stratified flow. This complicates the dispersion process and makes it more difficult to model the dispersion. The gravitational effect will by time make the plumes very thin (i.e., in transport over long distances). This has been observed in plumes over Europe from the 2010 Eyjafjallajökull and 2011 Grímsvötn eruptions; see Figure 14 in [4].

The treatment in the previous chapter shows that in density stratification a continuous plume will flatten out to the sides under influence of gravity. The fluid mechanics of the real problem are presumably more complicated than the simplified theory of the symmetrical plume, but the necessary temperature data for better analysis do not exist.

However, instead of in (A.7), there should be where is a correction coefficient of order of magnitude 0,1. In order to find , the temperature data has to be very accurate. The most accessible data for is from radiosondes but they are rare. Data on the horizontal spreading of plumes is accessible from satellite photos.

The correction factor has most likely the form . This is due to the following physical processes: (1) conversion of pressure energy into turbulent energy; (2) local cooling and mixing at the plume boundary; (3) wind shear and internal waves; and (4) acceleration terms in the initial phase when . Then it must be noted that is a correction factor on dissipated energy. The correction factor on velocities calculated from (A.7) and (A.8) is . The expected value of the velocity correction factors in Table 4 that corresponds to is around .

These are preliminary results that have to be verified with in situ measurements and numerical modeling. Small temperature inversions can easily explain the lateral spreading of volcanic plumes seen in satellite photos. Consequently, we have to rethink the diffusion problem. It may be the most interesting part that spreading of volcanic plumes by diffusion comes on top of the gravitational effect. Diffusion coefficients estimated from the total lateral spreading without regard to the gravitational effects will therefore be orders of magnitude to high.

In simulation models, plumes are normally assumed to disperse because of horizontal diffusion; vertical diffusion is normally left out entirely. This will only be true if there is strong density stratification that prevents vertical mixing and the wind is piecewise constant in direction and velocity along the path. Volcanic plumes usually find neutral buoyancy in a stable stratification like that in the troposphere. But then, gravitational flattening cannot be left out.

B. Streak Fallout of Volcanic Ash

B.1. Fluid Mechanics of Streak Fallouts

Normal fallout of particle grains from a volcanic plume is when the ash particles fall through the air with the terminal fallout velocity. Streak fallout of volcanic ash was only known from visual observation of fallout from volcanic plumes like Figure 7. Even the name is new, it is chosen because the streak fallout leaves in the air almost vertical line for a short moment, these lines are actually streaklines (not streamlines or pathlines) in the fluid mechanical meaning of the words. In a flow field, a streakline is made visible by constant injection of a dye in a fixed point in the flow.

To understand streak fallout it has to be recalled that in the terminal velocity situation the grains are affecting the air with a force equal to their weight. If this force is large enough, we have a vertical down flow of air going, or a downward vertical current powered by the density difference between the air in the current and the ambient air in the main plume. In high concentrations, the velocity in the vertical current can be many times the terminal fall velocity in still air of the grains in the streak.

There exists a meteorological phenomenon that resembles streak fallout; this is called a microburst. It is a much localized column of sinking air, producing divergent, and sometimes damaging straight-line winds at the surface. These downfalls of air are associated with single convective storms [31].

There are many laboratory investigations and analysis of the fluid mechanics of dense jets and plumes. Unfortunately many of them suffer from a mix up in Lagrangian and Eulerian parameters of the flow. In [32] these traps are bypassed; so it is used as basis for this analysis. Treating the streak as a dense plume or jet, momentum and mass exchange with the ambient air is by entrainment only.

The empirical rule in treating dense and buoyant jets and plumes is to put the ratio of the entrainment velocity and the average plume velocity as a constant. From [32, eq. (12.8)], we have and that dense plumes run on constant densimetric Froude number . The is a velocity distribution constant. Here we use a velocity profile that gives for a round plume; this is a little higher than the used for a planar plume [32, Eqs. (12.8) and (12.6)]. With as the radius of the streak and the acceleration of gravity this gives The volcano Sakurajima in Japan, Figure 9, emits ash in several explosions almost every day. There have been three airborne measurement campaigns, and in one of them big streak fallout was detected as the plane hit some streaks on July 27, 2013. Figure 8 shows a picture of big streak fallout and two smaller streaks encountered in the same trip. Figure 9 shows a measurement of the TSP (total suspended particles) in the streak, by a DustMate OPC optical particle counter.

Figure 10 shows the ash on the airplane when it came home. Figure 11 shows the results from the SkyOPC, optical particle counter in the plane. The two meters do not agree completely because the DustMate is saturated above 6000 micrograms/m3 and there are random fluctuations in the concentration values. Consequently the SkyOPC data are used in the analysis and Figure 11 shows the observed values, both the raw data and F-16 filtered values; the filtering is necessary to eliminate the random fluctuations clearly seen in Figure 11 as before in Figure 2.

The radius of the streak fallout is 666 m in Figure 11. If the temperature difference between the ambient air and the underside of the plume is small, the relative density difference between the streak fallout and the ambient air will be in kg/kg. Figure 11 shows in μgr/m3. Using the red data line in Figure 11, numerical integration of the concentration and velocity profiles gives the mass flow in the streak fallout  kg/s or 1,4 kg/s/km2 at the measurement level using (B.1); this corresponds to an average current velocity the streak fallout of  m/s. The sides of the plume slope together upwards at the rate (); so the outflow radius from the main plume is about 50 m smaller than the measured radius. This means 17% higher average concentration up there, or 17% dilution by entrainment, but this does not necessarily mean higher density difference up in the volcanic plume’s underside, as there is some small temperature difference between the ambient air and the underside of the plume cloud.

In the measurement results (Figure 11), there are two smaller streaks. Comparing different streaks the scale between the flows will be , where is the scale for the total ash flow, is the length scale (radius ratio), and is the concentration scale. The flow velocity in streak fallouts, scales in the ratio . For the two small streaks in Figure 11, the data gives and 1/3, respectively. The velocity  m/s equals the terminal fallout velocity of a 66 micron grain () while in these smaller streaks is approximately and their total load in kg/sec under 20% of the big streak. They also have a lower portion of fine particles (<PM10), 10–15% instead of 25% in the big streak.

This means that the coarse size grains are falling faster in these small streaks than the downwards air velocity , and in doing so they outrun the fine particle load (<33 micron), which is left in midair. Then the fine particles are no longer streak fallout, but ordinary fallout governed by terminal fall velocity. Small streaks can thus disintegrate or “die” in midair when the large grains in the streak outrun the density current and diminish the relative density difference and the air velocity falls below the terminal fall velocity of the large grains in the streak. This will eventually happen to all streak fallout columns, if it does not happen in midair; they must fan out when the current closes in on the ground and loses the downwards velocity. They are therefore difficult to observe except in airborne measurements. Results like Figure 11 cannot be obtained on the ground.

B.2. Discussion

The physical effect of streak fallout is to convey to the ground large flow of ash without the aero dynamical sorting in different grain size diameters as in normal fallout. Whole chunks can fall from the volcanic plume; one such may be seen in the left side of Figure 7. If streak fallout activity dominates over the effect of the normal fallout process, the total content of ash in volcanic plumes will be reduced much faster than ordinary fallout processes can explain, because the downwards air current takes all grains downwards at almost the same speed, the small diameter grains too.

From a value for the fallout, 1,4 kg/s/km2 with 25% in the aerosol range, it could look like the main plume is quickly drained for the bulk of the ash load. This is probably not so, each streak lives for short time only, there does not have to be a great number of large streaks in the air at any given time.

Near the ground the current fans out and becomes horizontal. Then the streaks are turned into ordinary fallout as the coarse grains fall down to earth but the aerosol size grains (particles < 10 micron) get mixed into the ambient air. This makes the streaks disappear before they reach the ground. Aero dynamical sorting in different grain size diameters is therefore partially active, but the sediments on the ground will be of very mixed grain size fractions. But while the current is still vertical the boundary to ambient air is rather sharp; so the streak fallout lines in Figure 7 have just the appearance that is to be expected for a dense current on its way down.

The density difference is on the one hand because the hot gas in the plume is of somewhat lower density than the ambient air; on the other hand there is a high ash concentration . The onset of streaks is a complicated instability phenomenon created by the temperature difference between the main plume and ambient air below and more or less impossible to predict in the time domain. Accumulation due to gravity increases the density of the underside of the cloud so the Rayleigh number of the interface is brought up to a critical value, a downward flow will start, but the critical value of the Rayleigh number is completely in the dark. Particle aggregation can play a role here, but streak fallout is so effective because of the velocity of the vertical current, not the higher terminal velocity of aggregated particles. Streak fallouts may be responsible for a large part of the volcanic ash fallout, and if that is so, a high portion of the fallout is in the aerosol range and this has to be accounted for in fallout estimations.

However, detailed fluid mechanical description of the mechanics of a streak fallout demands data that is not available for the moment but will hopefully be so in the future.

C. Conventional Dispersion Theory for Volcanic Plumes

C.1. The Diffusion-Advection Equation

In [19] Suzuki uses the following equation for the dispersion of the plume (diffusion-advection equation): Here, is the concentration of ash, is the wind velocity, is the eddy diffusivity or diffusion coefficient, and is the gradient/divergence operator in the horizontal coordinates (downwind) and (sideways). Other versions of (C.1) exist in the papers cited in the introduction, but the results will essentially be the same. There is no coordinate, as most researchers assume little or no vertical diffusion, that is, the (horizontal) ≫ (vertical) that can be counted as zero. Physically this means that the plume preserves the vertical thickness while the horizontal width increases due to diffusion. There is no physical reason for that vertical diffusion should always be very small in the same time as horizontal diffusion is large, and this is actually in contrast to the established fact that turbulence is three-dimensional. However, if the atmospheric stability is large, that is, if the plume is riding in a stable temperature inversion, this would result in a small vertical diffusion, but then we have the situation treated in Appendix A.

Most researchers of horizontal dispersion of volcanic ash use the point source approximation. This results in simple boundary and initial value problem that has the Gaussian plume as a solution to (C.1). The scales of this plume are derived in the following.

C.2. Scales of the Diffusion-Advection Equation

Concentration enters the equation in a linear manner; its true value has to come from the boundary conditions which is the top of the eruption column or the first horizontal part of the plume. Calling this boundary value we can insert and divide out of the equation. represents the source strength; it may be scaled separately according to a height/output formula if one exists.

Imagine two eruptions, denoted for model and for prototype, and we ask the question if the eruption can be a model of the prototype in the sense that there is geometrical and dynamical similarity between the two. We can insert dimensionless variables in (C.1) for the coordinates and , where and are scaling constants. ξ and are nondimensional variables for the downwind coordinate and the time. Omitting the subscripts and , we now get two equations, one for model, another for prototype, both looking like (C.2), The boundary and initial conditions for (C.2) are in ; in for . The condition for dynamic similarity is now that the two constants in (C.2) are equal in the model and the prototype. When this is the case, the model and the prototype share the same differential equation with the same boundary condition so the distribution of the dimensionless concentration will be the same function, , in the model and in the prototype. There are two constants in the equation and two scales to find. There are two natural scales, first the column height and then the wind velocity so from the first constant we get the scales for length, velocity, and time Equation (C.3) assures geometric similarity as this choice of scales ensures that . Now dynamic similarity requires that the second constant, , is the same for model and prototype. This requires that the scale for the diffusion coefficient, , is as follows: There are arguments both for and against if (C.4) holds or not. The critical point is the structure of the turbulence in the plume. Few researchers do report data on the turbulent structure of buoyant volcanic plumes, but in [33] Mikkelsen et al. discuss eddy diffusivity scaling in surface plumes in a great detail and come up with two candidates for diffusivity scaling: Richardson scaling and Bachelor scaling. Without getting too deep into turbulence theory, it can be stated that both scaling types require that turbulent shear is the source of the diffusion, but volcanic plumes are usually high above the planetary boundary layer and turbulent shear so high up is locally generated, usually connected to temperature stratification and cloud formation.

The validity of (C.4) has to be tested in each individual case. The true meaning of may be seen from the following simple treatment: is the downwind ( component) of the mass flux of the ash flow through a small unit area perpendicular to the instantaneous velocity. Let be the time dependent velocity in the -direction and let be the time dependent fluctuation, now : as the time average of the fluctuations is zero. Here is the correlation coefficient between the concentrations and the velocity, and and are their standard deviations. It is left to the reader to obtain the similar expression for the component of the velocity. In the diffusion-advection relation (C.1), the vector (subscript means vertical velocity fluctuations) is replaced by , according to the definition of in [34], which is the closure equation for the turbulent diffusion problem. This works most of the time in laboratory experiments in a confined space but has never been proved to work for buoyant plumes high up in the free atmosphere. However, the use of (C.1), for this problem, is state of the art and here it is suggested that the diffusivity scales according to (C.4).

C.3. Scaling of the Source Strength

The source strength in diffusion scaling is the ash concentration in the beginning of the plume layer meters above sea level. This concentration is orders of magnitude lower than the concentration in the ejected jet, not only because of the airborne magma that falls to the ground just besides the crater, but also because of the large entrainment of ambient air into the buoyant part of the vertical plume. This fallout from and entrainment in the vertical plume is a process quite different from dispersion by advection-diffusion.

The most popular formula for the total mass output is that the erupted mass in kg/s is proportional to . This is a scaling relation by Wilson et al. [20] and Settle 1978 [21]. Woodhouse et al. 2011 show very nicely how it is derived [27]. Woods solves the fluid dynamical problem in [22] and his treatment shows that the results of the Wilson and Settle relation can easily be wrong by an order of magnitude. In spite of these shortcomings the formula is popular among scientists. A way to start up a point source model is to assume that a fixed portion of the erupted material survives the transport to the elevation and thus becomes the boundary value for the mass transport in the plume. This seems to be the assumption of the many researchers using the Wilson and Settle relation. It can therefore be stated that the scaling (the is the cross section area of the plume; the symbol is used for “proportional to”) is state of the art, when no other information is available on the concentration value.

The area has to scale according to the length scale squared so the final scaling relation for the source strength renders the concentration scale that corresponds to the height/output relation : When and are known, the scale is and this can be used as valid for individual grain size fractions.

C.4. Ash Fallout Scaling

Ash fallout has been modeled successfully by Suzuki [19]. Several authors have worked on this problem. In [16] the diffusion-advection equation is used and the fallout evaluated in the limit , the result shown in their [19, eq. (5)]. A further study of the fallout process shows that the horizontal parts of all fallout trajectories scale according to the length scale for all particles independent of size. Therefore the vertical parts of the trajectories have to scale in the length scale to, and the fallout velocity of ash, , has to be scaled from the model to the prototype; this adds up to The stands for the settling velocity, or terminal fallout velocity, it depends on the aerodynamic diameter and density of the ash. The dependence is quite well known, but the formulas for the settling velocity may be very different, although the render the same result. Bonadonna and Phillips publish in [35] a set of very simple formulas for all particle sizes with a clear range of validity. The finest fraction  mm follows Stokes law; the coarsest fraction  mm has a constant drag coefficient. In between these limits . The settling velocity of Bonadonna and Phillips [35] gives us the scale for the diameters of the various fractions that follows the trajectories that scale in between and : This means that if we select a place in the prototype and find grains in the ash sediment of category 1, we should find grains in the corresponding place in model of diameter .

There is a difficulty in this approach if the grain size scaling places the model grains in a category different from the prototype grains. This is especially true for category 3; here the scale is that can be a large number. However, this is not as bad as it seems; this scale is only valid when both and are in the same category. As an example we can take ; then the prototype grains have to be 5 mm, but grains that large fall to the ground in few minutes, so they hardly leave the area where fallout from the vertical column heaps up and diffusion theories do not apply anyway. It can also be inferred from [35] results that grains  mm do not get far away from the vent in eruptions of moderate size. In practice, accidental grains of large sizes can be found far away from the crater; this is the effect of the turbulent fluctuations, which are suppressed in the closure equation in [34].

This leaves us with a scaling problem if category 1 grains in the model, correspond to category 2 grains in the prototype when scaled. However, as the settling velocity is to scale according to , the true grain size scale can easily be found in the following manner: Here the constants (category 2)/ and (category 1)/ according to [35, eq. (A4)] for the fallout velocities. The constant is 0,080–0,1 millimeters, depending on location and type of eruption. An evaluation of may be needed in the grain size range micrometers. An educated guess for how much of the total ash will be in this size range is about 5%. Smaller sizes will scale according to category 1 ().

C.5. Scaling by Use of Self-Similarity of Fallout Trajectories

Close to the ground the wind varies with elevation and the trajectories of ash particles falling down to earth become very complicated paths. In estimating ash fallout they have to be modeled, directly or indirectly. Using the settling velocity alone, with or without a random component due to diffusion included in the estimation, has proved to give good results [16, 19, 35].

The settling is an aerodynamic process that sorts the fallout in finer and finer grain sizes with increasing distance from the vent. Fallout samples taken will be well sorted; that is, the grain size distribution curve will show that the majority of the sample is in a relatively narrow diameter range around the 50% diameter, . In order to produce this sampling result, the total grain size distribution [36] of all the airborne ash has to be well graded; that is, it has to contain all grain sizes found in all the fallout samples. Such grain size distributions are more or less straight lines in a lin-log plot or where is the weight ration of grains of smaller diameter than () and and are constants. Equation (C.10) can very well be approximated by a log-normal probability density distribution and this seems to be a popular approach; see [35] and many of the papers in there.

Using the simple fact that in our model all particles in the prototype have to fall the same vertical distance, , down to earth, we can deduce the following self-similarity relation for two particle samples, 1 and 2, in different downwind distances from the vent, and , and having two different due to the aerodynamic sorting The last part assumes that scales according to . Here denotes the time it takes the particles to travel the distance down to earth along its trajectory and end up in the distance from the vent. Equation (C.11) does not include the mixing effects from the diffusion and the velocity profile in the shear layer so it is approximate only and has to be verified a posteori. The last part of (C.11) suggests that the 50% diameter might be used instead of the settling velocity of the particular trajectory. This is the median diameter; it is a straight forward quantity to measure in each fallout sample. We have to use such a measure when we use the self-similarity scaling relation (C.11) because the plume has a finite thickness that creates more than one trajectory for each from different elevations in the plume. There will therefore be some mixing of different grain sizes in each sample.

Some researchers suppose that the trajectories are straight lines. They are more like a parabola, but this does not affect the scaling if the wind profile is self-similar for all . Large surface structures change the wind profile and the wind velocity changes from time to time so (C.11) can easily fail. Never the less, (C.11) can be quiet useful as it is universal in character. The Mount St. Helens 1980 data follows it nicely [37].

C.6. Case Study, Dispersion of the Eyjafjallajökull Plume in Iceland 2010
C.6.1. Self-Similarity

There are many studies of the Eyjafjallajökull eruption in Iceland in 2010. Here two are mainly used [29, 38].

Figure 12 is grain size data from a day-by-day record of the 2010 eruption shown in [29, Figure 2]. The authors of [29] use a modified form of Sparks equation as plume height-mass output relation. April 14–16th had relatively steady wind and output. Figure 12 shows 4 samples and Figure 13 shows their self-similarity relation. The aerodynamic sorting from the well graded sample (2 km) to the other well sorted samples is clearly seen in Figure 12. The power law trend line in Figure 13 has the exponent −0,926 instead of −1,000. This is a small deviation from the self-similarity scaling equation (C.11).

To use this scaling as in Figure 13, has to be constant. For variable , the scaling relations equations (C.3)–(C.5) are used to scale the eruption to standard values in output and wind speed before the self-similarity relations are applied. The results are in Table 1.

The data listed in Table 5 is from and data in [29, 38]. The part of their data, where (C.10) gives stable values for estimates of , is used in Table 5. Figure 14 shows that the scaled data supports the suggested self-similarity relation nicely (green line), but the unscaled data shows no self-similarity rule. Few values, shown by bold in Table 5, do not scale. The reason for that is unknown but there are several possibilities. To read from Table 5 if the scaling is successful, , , and (theoretical value of ) have to be compared.

C.6.2. Finding the Value from Satellite Photo Data

Close to the crater there is no better estimate for the concentration distribution at the zero buoyancy level than a constant , inside the plume’s effective cross-sectional area, the top hat distribution used in [22, 27]. The analytical solution to (C.1) is the Gaussian normal distribution: This corresponds to that is the effective width () and height () of the plume with respect to the centerline distribution . is the crosswind effective area as before. When the Gaussian distribution is used as a solution to (C.1), one finds that  Kt with .

In a continuous plume in a steady wind the effective width of the mixing layers on each side will increase downwind if , due to the diffusion. In the center region, will be constant until the mixing layers on both sides meet in the center of the plume. This is called the point . There is no fallout term in (C.13) so the effective width and height of the plume are constant. Downwind of that point will diminish and increase. We also follow the “state-of-the-art” method and assume no vertical diffusion. The no fallout assumption will be discussed later.

The visible plume is expanding to both sides with the diffusion velocity that carries visible concentrations sideways in all directions. Satellite photos like Figure 15 show the plume between the two visible sides, but the contamination beyond this limit is not zero. Figure 15(a) shows a close-up of the crater region. The total visible plume is in Figure 15(c) and the mixing layer with the assumed concentration profile in a cross section is in Figure 15(b). The red lines mark the widening of the visible plume and will be used to find .

Tiesi et al. [28] use the following equation: accepting his model and modifying (C.14) to fit the mixing layer we end up with in (C.15) is chosen instead of the area in order to make the horizontal effective width equal to in the point . In this model, the mixing layer (Figure 15(b)) follows (C.15) with a and effective width (denoted in Figure 15(b)), on each side. As the mixing layer is upstream of , the value of in (C.15) is negative. In (C.15), decreases backwards from the point until it becomes zero where the mixing layer begins. We want to find this length of the mixing layer and start by assuming it is the length in Figure 15(a) assuming to be in between and .

The red lines in Figure 15(a) are the visible boundary of the plume; it is the locus (visible ). According to (C.12) the in Figure 15(b) will increase as along these red lines as increases according to (C.15). Thus the visible plume widens in the direction according to On May 11, 2010 the visible concentration and the width of the plume were measured from an airplane carrying an OPC and a GPS tracking device [4]; the original data is displayed in [1]. On this day the maximum plume height and wind are approximately stable. can now be evaluated by inserting (C.15) into (C.16), differentiating with respect to , insert the basic data and evaluate in the point   ().

The basic data is Eruption mass output: 1,74 108 g/s in [29].  m/s;  m (4000 ft.);  g/m3 [1]. The data measured from the satellite photograph Figure 15(a) is displayed in Table 6.

In Table 6   is used for and the position is assumed in . Afterwards, is found 753 m, meaning that is in this distance to the right of . It should be moved to there, a new boundary drawn, a new found, the process repeated and a new value of found. However, this resulted in an insignificant change in the value. seems a credible value. It also compares well with [28]. is a measure for the diffusion velocity, it is approximately 1,5 m/s, or the same order of magnitude as can be expected in the planetary boundary layer just above the ground in this wind conditions.

Fallout will affect , but not the determination of if it does not affect the visible boundary of the plume significantly. From the value we can estimate the (or ) where centerline concentration is thinned down below visible concentration; this place can be seen in Figure 15(c), but we need to know a fallout coefficient to do so. To use a fallout coefficient and assume is a rather crude, but popular approximation and it gives an idea of the average fallout in percentage of the mass flux in the plume. In our case the original satellite picture shows the plume disappearing about 800 km downstream (Figure 15(c)). The result of this is  km−1 so the average fallout is a little under 1% per km. The true value is presumably much higher close to the vent but likely low enough so the constant flux approximation behind (C.15) is justified for the 13 km from to in Figure 15(a).

Conflict of Interests

The authors declare that there is no conflict of interests regarding the publication of this paper.