Inverting Sediment Bedforms for Exploring the Hazard of Volcanic Density Currents Directly in the Field

7 Pyroclastic density currents are ground hugging gas-particle flows moving at high speed down the 8 volcano slope. They are among the most hazardous events of explosive volcanism, causing 9 devastation and deaths 1,2 . Because of the hostile nature they cannot be analyzed directly and most of 10 their fluid dynamic behavior is reconstructed by the deposits left in the geological record, which 11 frequently show peculiar structures such as bedforms of the types of ripples and dunes 3,4 . In this paper, 12 we simplify a set of equations that link flow behavior to particle motion and deposition. This allows, 13 for the first time, the build up of a phase diagram by which the hazard of dilute pyroclastic density 14 currents can be explored easily and quickly by inverting bedforms wavelength and grain size.

related to flow parameters such as the densiometric Froude number ′10 and the critical Shields 27 number 11 . For symbols see table1. 28 ′ = √ ′ is a balance between inertial and gravitational effects, with ′ = ( − ) 29 representing the reduced gravity, g the gravity acceleration, V the current velocity, H the current 30 depth. 31 is the density of the fluid-particle mixture with s particle density, f fluid density and C particle

108
We never found antidunes, in fact their interpretation has always been questioned in volcanic 109 deposits 11,27 . 110 The software PYFLOW 2.0 by Dioguardi and Mele 31 has been used to plot data in Fig 1. It was  111 implemented here so to obtain both the impact parameters of the current together with Srw and Qb. 112 The software employs sediment data that result from time-consuming laboratory analyses, which 113 involve technologies and calculation resources not available to all scientists (see the Method section). 114 The aim of this paper is to rearrange and simplify the dataset in order to construct a phase diagram 115 by which to invert Wl and D of PDCs' deposits bedforms and obtain the impact parameters directly 116 in the field, without the need of the extra terms that require extensive work in the laboratory. 117 By means of regression analysis we obtained three fitting laws (Fig. 3a, b and c) that correlate just 118 few of the many terms of the formulas of (5), (6) and (7)

127
A fourth fitting law between Srw/Qb vs Wl with a good correlation was obtained by selecting 32 128 deposits characterized by well exposed bedforms ranging between ripples and dunes (Fig. 3d). 129 The fitting laws allow interpreting Qb and Srw in terms of the deposit formation processes, and also to 130 relate them to the current's flow parameters. In Fig. 3a a relationship between Qb and C * 2 is shown. 131 Since C is directly proportional to mix and mix * 2 is the turbulent shear stress of the current 32 , it 132 means that Qb is proportional to the shear stress, which confirms the finding of sedimentary currents 29 . 133 The relationship between D * and C * 2 of Fig. 3b implies that shear stress is proportional to bedforms 134 that can be solved numerically, once D and Wl are specified, to obtain * , C, Srw and Qb. Important 152 information on the hazard of dilute PDCs can be obtained from the first three parameters. C and * 153 serve for the calculation of Pdyn (2), since C is used for obtaining mix in (2) by means of (1) and * 154 is used for the calculation of V by means of the law of the wall of a turbulent boundary layer 32 155 which is the physical model of PDCs that we employ (see the method section), where V(y) is the 157 velocity profile of the stratified flow 32 , and ks is the substrate roughness. 158 When comparing results obtained by (8) with those resulting from PYFLOW 2, the average absolute 159 error of * is 28% and that of C is 30%. This means that a good approximation can be achieved for 160 exploring the range of impact parameters by means of the simplified formulas, without the terms that 161 involve extensive laboratory analysis. 162 The absolute error of Srw is about 45%. While it is larger than that of * and C, we discuss also the 163 role of Srw because it allows the calculation of flow duration, t 15 , which is an important factor of 164 hazard. The total time of aggradation is a proxy of flow duration, t, which is equal to deposit thickness, 165 Hdep, divided by Ar, the aggradation rate. Sedimentation occurs by continuous aggradation during the 166 passage of the current, and is equal to Srw divided by one meter, which is the reference width of 167 the sedimentation rate per unit width, (see Dellino et al. 28 ). Therefore, flow duration, which 168 approximates the time in which harmful concentrations of ash are suspended in the current to which 169 a human being can be exposed, can be calculated by means of Srw. With our model a reasonable 170 approximation can be achieved also on such a relevant parameter of PDCs. 171 In Fig. 4, which was constructed by means of (8), the main flow variables and impact parameters are 172 shown as a function of D and Wl. The Wl range was set between 10 and 300 cm. Bedforms with larger 173 Wl can be found in the geologic record of volcanic deposits, but this scenario is out of the range of 174 applicability of our model. We are, in fact, considering bedforms that develop on an almost flat 175 surface. Much larger bedforms, instead, typically develop as an interplay between the current's flow dynamics and large ground morphology elements 27 (e.g. ridges, big obstacles). The range of D of Fig.  177 4 was set between 4 and -2 phi (0.0064 mm and 4 mm respectively). We do not include coarser values 178 because, in volcanic sediments, larger sizes (coarse lapilli and bombs) do not form dunes, but 179 lenticular beds representing highly concentrated traction-carpets at the base of PDCs 34,35 , to which 180 our model does not apply.

188
In Fig. 4a the value of velocity V, as averaged in the first 1000 cm of the current, which was obtained 189 by integrating (9) over flow height, with ks = 10 cm, is plotted against Wl. We chose this depth- shear stress is needed to move larger grain sizes. The volumetric concentration, C, ranges from less 201 than 0.001 to about 0.017 (Fig. 4b). It decreases as Wl increases, and it does so for all grain sizes, 202 although at a rate that decreases at decreasing D, because a higher concentration favors a higher Srw 203 (see (7)) and a larger Srw/Qb ratio, hence a smaller Wl. The evident change in trend with decreasing 204 grain size can be explained by the fact that the finer the particles, the lower the concentration required which is a value that does not cause severe damages to buildings 1,13 , to almost 30 kPa with larger Wl 212 and coarser D, which can destroy even the more resistant, modern buildings of reinforced concrete 1,13 . 213 The sedimentation rate (Fig. 4e) increases as grain size coarsens, meaning that with finer sizes flow 214 duration is longer, as it is expected since finer sizes result in a smaller settling velocity. As far as the 215 flow duration with longer bedforms. With the coarsest sizes, instead, the sedimentation rate decreases 217 as Wl increases, meaning a longer flow duration with longer bedforms. 218 The ranges of Wl and D used in (8) for obtaining the trends of Fig. 4 replicate the ranges in our dataset, 219 and result in parameters that span from currents that do not impact severely on structures, to values 220 of devastating effects. Such a range well represents the situation of large-scale PDCs whose strength 221 decreases along runout 15 , and change from totally destructive flows around the volcano to residual 222 currents that in the distal outreach do not possess a high strength but can still be rich of ash. Such fine 223 glassy material can be highly dangerous to breath even at concentrations lower than 0.001 38 if flow 224 duration t, which can be calculated by means of Srw, lasts more than a couple of minutes. Thus, with 225 our model it is possible to invert bedforms of past eruptions, and follow the different aspects of PDCs 226 hazard as they evolve along flow runout. 227 In order to help scientists not availing of numerical resources to take advantage of our results, we 228 solved (8) at discrete intervals of D and Wl and constructed a phase diagram where the stability fields 229 of Pdyn, C and Srw are represented inside a grid (Fig. 5). The values are averaged among the four 230 neighboring grid points and the uncertainty is expressed in terms of one standard deviation.

238
In the supplementary Information, additional diagrams with ks = 10 cm and s = 1000 kg/m 3 ; ks = 30 239 cm and s = 2000 kg/m 3 ; and ks = 10 cm and s = 1000 kg/m 3 are included (Supp. Fig. 1,2 and 3   240 respectively), and a table is also provided (Supp. Tab.1) where the values of * , C and Sr are set at 241 half phi intervals of D in relation to Wl. By means of these data, and specifying in (4) and (9) the 242 value of ks, s, and H at which to integrate V, more precise data of the impact parameters can be 243 obtained. 244 With our diagrams and tables at hand it is thus possible for every scientist working on hazardous 245 volcanoes to make an exploratory hazard assessment by means simply of the wavelength and grain 246 size of bedforms. It is true that bedforms are not always well exposed in their complete longitudinal 247 profile, because of truncations due to erosion. Sometimes they are also difficult to measure precisely, 248 because a direct access to the deposit is hard. Anyway, our experience tells that dilute PDCs most 249 always leave well-preserved bedforms as a trace of their passage. Scientists working on active 250 volcanoes are encouraged to look for good outcrops where bedforms can be measured. By means of 251 our phase diagrams, now they have a tool for exploring the behavior of hazardous pyroclastic density 252 currents directly in the field. 253 254

Method 255
The reconstruction of the impact parameters of PDCs is based on a flow mechanical model that starts 256 with the assumption that the turbulent current is velocity and density stratified 12,40 . In the stratified 257 multiphase gas-particle current, the basal part is a shear flow that moves attached to the ground and 258 has a density higher than atmosphere (Fig. 6). The upper part is buoyant, because particle 259 concentration decreases with height down to a value that, combined with the effect of gas temperature, 260 makes the mixture density lower than the surrounding atmosphere.

287
The shear-flow height and density are obtained by solving the system of (11) and (12), which is valid 288 for a turbulent current 289 where is the shear-driving stress of the flow moving down an inclined slope of angle . 292 The density profile, which is a function of concentration, particle density and gas density, is: 293 The gas density and Rouse number are obtained by solving numerically the following system: 295 Equation (14)