Combined effusive-explosive silicic volcanism straddles the multiphase viscous-to-brittle transition

Magma is a viscoelastic fluid that can support fracture propagation when local shear stresses are high, or relax and flow when shear stresses are low. Here we present experiments to confirm this using synthetic and natural magmatic liquids across eruptive conditions and use Maxwell’s linear viscoelasticity to parameterize our results and predict the maximum stresses that can be supported during flow. This model proves universal across a large range of liquid compositions, temperatures, crystallinity and rates of strain relevant to shallow crustal magma ascent. Our results predict that the 2008 Volcán Chaitén eruption resided in the viscous field at the onset of magma ascent, but transitioned to a mixed viscous-brittle regime during degassing, coincident with the observed combined effusive-explosive behaviour during dome extrusion. Taking a realistic maximum effusive ascent rate, we propose that silicic eruptions on Earth may straddle the viscous-to-brittle transition by the time they reach the surface.

A n apparent dichotomy exists in volcanic eruptions on Earth between events that are dominantly effusive-producing lavas-and those that are explosive in characterproducing volcanic ash fragments in plumes or fast-moving density currents [1][2][3][4][5][6] . However, these are end member descriptions and single eruptions can transition from one style to another 4-7 , or exhibit combined effusive-explosive behaviour simultaneously 8,9 . These differences in eruptive behaviour cannot be attributed to variations in primary volatile contents at storage conditions, nor to magma composition alone, and instead are dominantly controlled by the details of ascent dynamics 3,10-13 . Indeed, magmas can undergo local brittle fracturing in otherwise effusive events 11 .
Some of the most spectacular natural field examples from which we can infer this combined behaviour-fracture-formation during otherwise viscous flow of magma-include the preserved Mule Creek volcanic vent 14 (New Mexico, U.S.), and Torfajökull 15,16 or Hrafntinnuhryggur 17 (Iceland). Recently, direct observations of combined explosive-effusive styles have been made at eruptions of Santiaguito volcano (Guatemala) from 2007 to 2015 18 , and at two eruptions in Chile: Volcán Chaitén in 2008 7,19 and Cordón Caulle from 2011 to 2013 8,20 . In these events, detailed documentation could be made of fractures opening through an effusive lava flow or dome lava, and acting as gas-escape pathways and explosive ash vents 8,18 . When such eruptions end, these fractures are often locally welded shut, or filled with ash and pyroclasts that have sintered back together, demonstrating that these fractures formed and then closed all while ambient temperatures remained high 16 .
Local or wholesale fracturing of hot magma occurs because elastic stresses accumulate in the volcanic liquid when the rates of deformation are high 1,11,[21][22][23] . This view of magma rheology is supported by the occurrence of volcanic rocks that can be shown to have fractured and rehealed while hot 14,16,24,25 , experimental evidence demonstrating the feasibility of breaking magmatic liquids or suspensions at eruptive temperatures 21,23,26 , and numerical simulations of fracturing in viscoelastic magma during flow 11,27 . Maxwell proposed the simplest model of how elastic shear stresses are stored or relaxed in fluids; such that relaxation occurs over a characteristic time λ r = μ/G ∞ where μ is the liquid viscosity (Pa.s) and G ∞ is the elastic shear modulus (Pa) 1,28,29 . We can compare this time to the time characteristic of flow λ ¼ 1=_ γ where _ γ is the shear strain rate (s −1 ). This balance between forced flow and relaxation of stress yields a dimensionless strain rate, or a scale-independent Weissenberg number such that Wi ≫ 1 is a first-order constraint of the unrelaxed brittle field, while Wi ≪ 1 represents the viscous field. The value of G ∞ is approximately 10 10 Pa and broadly independent of composition or temperature 28,30 . Viscosity μ, however, is strongly dependent on temperature T and can be predicted using empirical models that additionally account for liquid composition including its dissolved volatile concentrations C H 2 O (in wt.%; where total dissolved water is the dominant volatile phase in most shallow magmas) 31 . Therefore, knowing G ∞ and μðT; C H 2 O Þ, we can predict λ r for any magmatic liquid at any eruptive temperature.
We explore the efficacy of Eq. 1 in scaling transitions from viscous to elastic behaviour by conducting experiments in a largescale press in which we load a cylinder of synthetic or natural volcanic glass, heat it to a set temperature, and compress it in the direction parallel to the axis of rotational symmetry. We control the upper piston velocity at a constant value, which confers a constant axial strain rate _ γ a , and we record the evolution of axial stress σ (computed from axial force). We refer the reader to the methods section for full details of the experimental technique, equipment, the calibration at low Wi, uncertainties, and data treatment details. We characterize the response of the liquid as viscous if the stress stabilizes to a steady state after an initial transient, as brittle if the stress linearly increases with time toward a single large drop to zero stress, or as transitional between these two end-member behaviours if the stabilisation to a steady state is interrupted by drops in stress.

Results
Healed brittle fractures in effusive or combined explosiveeffusive volcanic events. There is a large body of textural evidence for brittle fractures in silicic magma that have healed, all while temperatures remain high 14,16,24,25,32,33 . In Fig. 1, we provide an example from Hrafntinnuhryggur, in Krafla, Iceland (photo taken by Jenny Schauroth). Other examples have been documented in the products of both dominantly effusive 11,17 , and mixed explosive-effusive 24 silicic eruptions. These healed, or partially healed, fractures often contain welded fragments 9,34 . We take these observations as a motivating starting point in the design of our experimental campaign to understand the conditions required for fracture formation in otherwise viscous magma.
Experimental results. We map our experimental data according to the response of the liquid in the regime space between λ and λ r , such that curves of constant Wi can be plotted. We additionally reanalyze published results for single-phase liquids for which sufficient detail was available to compute λ r and λ (refs. 22,26,35 ), and map them in this same regime space (we use a constant G ∞ = 10 10 Pa after refs. 28,30 ). Since the Weissenberg number accounts for the temperature-and compositional dependence via the incorporation of μ, and in the absence of any other timedependent process, the transition from viscous to brittle behaviour will occur at a single Wi. Our experiments imply an important threshold at a critical Wi c~1 0 −2 , and thus delineate macroscopically viscous (Wi < 10 −2 ) from macroscopically brittle results (Wi > 10 −2 ). If Wi is increased further, then brittle behaviour is localized in a single explosive fracture event in the 1 cm Fig. 1 Photomicrographs of obsidian from Hrafntinnuhryggur, Krafla volcano (Iceland) collected in situ toward the margins of an effusive dyke-geometry conduit 17 with textural evidence for fracture formation and subsequent healing while at high temperature (marked with arrows). These images were taken by Jenny Schauroth system, leading to a single large stress drop for Wi ≥ 0.04, which represents the dimensionless strain rates required for strain localization in a confined sample, or thorough shattering of the glass in an unconfined sample.
Viscous dissipation as heat. There exists an additional process in which heat is produced by internal friction in flowing viscous liquids, which is parameterized here by the power density Φ g ¼ 2μ_ γ 2 . Such heat can be conductively lost out of a system of a given size with a power density Φ l ¼ ∇q ¼ κ∇ 2 T where q is the heat flux (Wm −2 ) out of the system of a given area, and κ is the thermal conductivity (W m −1 K −1 ). The ratio of these powers is a Brinkman number Br = Φ g /Φ l such that Br>1 yields net heat gains due to viscous dissipation as heat 22,36 . We can constrain this for our system, showing that for high λ r (high μ), the samples will fracture before the accumulation of heat will occur. By contrast, at low λ r , heat storage prior to brittle behaviour will occur, potentially reducing λ r further and pushing the system away from the critical Wi = 10 −2 . As Br is scale dependent, these regimes (of heat retention vs. dissipation) shift with system size. Here for our experimental data, we can plot the threshold Br = 1, calibrating Φ l ≈ 10 4.5 W m −3 for the scale of our samples and setup 22 (Fig. 2). Data that showed a smooth, monotonic, timedependent decrease in the peak stresses reached, can be characterized as heating by viscous dissipation, and are consistent with the region of Fig. 2d in which Br > 1 (below the line Br = 1 in this space).
Extensions to account for multiphase magma. Magmas are rarely single-phase liquids. Two-phase crystal-bearing or bubbly magmas generally break at lower applied stresses or strain rates than the pure liquid phase 21 . To test this effect, we compile data from published sources 21 collected using an analogue for rhyolitic magmas that contained crystals, and sheared in a simple-shear, torsional Paterson press. We find that in the case of crystalbearing magmas, we can account for this effect by scaling the critical Weissenberg number Wi c for the effect of crystals, by defining the two-phase critical Weissenberg number as Wi ′ c ¼ Wi c ð1 À ϕ x =ϕ m Þ where ϕ x and ϕ m are the crystal volume fraction, and the maximum packing crystal volume fraction, respectively. Following Cordonnier et al. 21 we find that we must exchange μ for η in the definition of Wi c , where η is the measured suspension viscosity. The efficacy of this result is shown in Fig. 3 where ϕ m = 0.74 for the moderately polydisperse crystal size distribution is given in ref. 21 . Figure 3 demonstrates that Wi ′ c is shifted to higher λ as the fraction of crystals is increased (e.g. compare Fig. 3a with Fig. 3d).   Results of this study are compared with a reanalysis of previous experimental results for magmatically-relevant pure liquids 22,26,35 . All data are consistent with viscous behaviour at Wi < 0.01 and brittle behaviour at Wi ≳ 0:04. Where G ∞ was not known, we have taken 10 10 Pa for silicate liquids 28,30 An equivalent way to scale for the presence of crystals is by making substitutions in the definitions of λ r and λ, redefining Here γ′ is the strain rate in the liquid amplified locally as liquid has to move more rapidly around crystals for a given bulk strain rate on the system _ γ as described above. This results in a critical Weissenberg number as Wi ¼ λ ′ r =λ′, which is equivalent to the Wi for singlephase liquids, and permits all data to be plotted together in a master regime diagram for multiphase and single-phase magmas where the breaking point Wi c is at a constant position (Fig. 4).
The same scaling approach as above has been attempted for analogues of bubbly magma that are rapidly decompressed 37 . Here, the characteristic decompression time appears to be a firstorder constraint of λ, separating results in which the bubbly samples expanded viscously from those which fragmented by bubble-overpressure at Wi~10 −2 , as predicted by our experimental constraints. However, in detail, λ has a more complex relationship with decompression rate, and would require further work to scale across all Wi.
A universal map for shear fractures in magma. In light of the multiphase constraint presented above, we can collapse the data for crystal-bearing, and single-phase liquid samples to the same regime map given in Fig. 2, when λ r is replaced by λ ′ r and λ by λ′ in the case of multiphase magmas. This gives a unified description in which all magma begins to break as Wi increases above 10 −2 (Fig. 4).
A framework for predicting shear-fracture in silicic eruptions. In this framework, some observations at silicic volcanoes can be tied to specific physical quantities via Eq. 1 with predictive power. To give a range of examples, we compile the range of λ r we can calculate for individual eruptions in which some transitional explosive-effusive behaviour has been observed. We do this for three 'families' of silicic eruptions. Using our analysis, the critical deformation time λ c or critical scaled deformation time λ ′ c can be predicted for all eruptions of these types, thresholding the minimum deformation rates associated with the onset of brittle behaviour, or associated with mixed behaviour that is locally brittle. Moreover, as all magmas degas during ascent and eruption, producing spatially heterogeneous, and generally increasing μ (and therefore an increasing λ r ), Fig. 4 demonstrates how the resultant range of relaxation times and evolution thereof would confer drastically different deformation behaviour and almost inevitably results in brittle behaviour during ascent and degassing.
A detailed look at the 2008 eruption of Volcán Chaitén can provide further insights. Using the temperature (825°C during initial ascent 7 ), water contents 24 and average strain rates constrained for the magma storage conditions and initial propagation to the surface 7 , we can provide a range of possible λ r and λ for the early phase of the eruption. This range plots in the viscous field (Fig. 4). This indicates that conduit-scale shear fracturing was not likely to be a dominant component of the initially Plinian explosive phase of the eruption. In this initial phase, the classical magma fragmentation by bubble overpressure process is likely dominant. As the eruption progressed, the magma degassed, and combined explosive-effusive behaviour resulted in dome extrusion punctuated by ash-venting 9,19,41 . Similar behaviour has been inferred to have occurred during silicic eruptions elsewhere 42 . Using the water contents measured in the dome lava 24 and the dome extrusion rates 41 , we can compute a new range of λ r and λ for the mid-stages to late-stages of the eruption, which straddles the critical Weissenberg number threshold (Fig. 4). This analysis implies that even as contiguous magma reached the surface, and the bulk of the explosive behaviour ceased, conduit-scale brittle fractures could open up, which is consistent with the observation of combined explosiveeffusive eruption styles.
To generalize beyond a single eruption, we use the observation that constraints of silicic (rhyolitic and dacitic) magma ascent rates at most eruptions delineate dominantly explosive eruptions from effusive eruptions 43 . We note that the ascent rate that separates one eruptive style from another occurs around 0:01 < u < 0:5 m:s À1 (ref. 43 ), which therefore represent the maxima for an effusive eruption. As an approximate scaling, we note that the average conduit strain rate h _ γi will be h _ γi % u=L, where L is the conduit radius. For bracketing values 5 < L < 50 m, we can predict a wide range of possible h _ γi for silicic eruptions, which in turn, give a range of 10 < hλi < 5000. Using Fig. 4 as a guide, we see that for silicic magmas ascending with this range of deformation timescales, will dominantly be in the viscous field Wi < 10 À2 when close to their shallow storage conditions. However, using the range of λ r plotted for silicic eruptions at degassed or partially crystallized conditions equivalent to the shallowest portions of final ascent, the same deformation timescales clearly make brittle behaviour very likely, even for the upper end of possible conduit radii and using the conservative estimates of h _ γi. We therefore hypothesize that ascent-driven locally brittle behaviour at silicic volcanoes is common, and may be a key process in the recently observed combined explosive-effusive behaviour.
In our experiments, the bulk axial stress σ is approximately related to the mean shear stress τ by τ = σ/3 (by using the Trouton ratio 28 ). If we constrain our analysis to the regime space Br < 1, we can demonstrate that in the viscous regime, the peak shear stress value τ m achieved in any experiment at steady state, scales with Wi via τ m = G ∞ Wi predicted for viscoelastic liquids (via Eq. 1 assuming τ m ¼ μ _ γ in the low-Wi regime). In the brittle regime, τ m is broadly independent of Wi and is given by τ m ¼ 10 À2 G 1 , which for magmatic liquids is τ m ¼ Oð10 8 Þ Pa, consistent with the order of magnitude for the strength of glassy materials 28,44 . Using Eq. 1, we can additionally reanalyze published experimental data collected using a large range of silicate glass-forming synthetic and natural magmatic  Andesite glass (Webb & Dingwell, 1990a) Nephelinite glass (Webb & Dingwell, 1990a) Rhyolite glass (Webb & Dingwell, 1990a) Tholeiite glass (Webb & Dingwell, 1990a) 26,28,35,[44][45][46] , all collected over a wide range of temperatures, demonstrating consistency with our analysis (Fig. 5a). Our hypotheses are supported by observing the signatures of acoustic emission events-representing episodes of fracture nucleation and rapid growth. These increase above background noise at Wi ≳ 3 10 À3 , and become dominant at the threshold beyond which we observe macroscopic brittleness at Wi > 10 −2 (Fig. 5b). Time-dependent drops in stress indicating brittleness are not directly observed in our macroscopic mechanical data as low as Wi = 10 −3 (Fig. 2), possibly because the damage rate is low and small nucleating fractures are partially healing on the timescale of deformation. While our AE observations potentially lower the critical Wi c , we note that our use of cylindrical samples, for which the Trouton ratio for shear stresses and strain rates is only an approximate correction, make this difficult to interpret at this resolution.
We have found that, in general 0.01 < Wi < 0.04 represents a transitional window from viscous to brittle behaviour. To understand the Wi thresholds quantitatively, we compile data for viscous volcanic liquids deformed in the frequency domain where the amplitude of strain was sufficiently low to avoid fracturing 47,48 . In these experiments an oscillation of strain at forcing frequency ω results in a measured evolution of shear stress. In turn, this can be decomposed into the real G′ and imaginary G″ components of the elastic modulus. The normalized frequency ωλ r is a Deborah number (the frequency domain counterpart to the Weissenberg number 49 ), collapsing all data to a temperature-independent description. Such data demonstrates that the onset of effective non-zero accumulation of elastic stress also occurs at De~10 −2 (Fig. 5c), consistent with changes in τ m (Wi), acoustic emission data, and the brittle onset (Figs. 2 and  4).
Herein, we provide a non-dimensional framework for understanding the onset of brittle shear fractures in silicic magmas, accounting for multiphase heterogeneities of suspended crystals. Through the application of this framework to a few type-localities worldwide, including two recently observed rhyolitic domeforming eruptions, we show that it is likely that shear fracturing is an important process in the shallow portions of effusive magma ascent with concomitant degassing. At storage conditions and the early parts of ascent, most silicic magmas appear to reside at conditions amenable to viscous behaviour without brittleness, however, we show that it is almost an inevitability that during ascent and degassing of volatiles, magmas will straddle the viscous-to-brittle transition at the Earth's surface. Therefore, models of magma ascent that include the acknowledgment that ascent can be accommodated by sliding on fracture networks, are preferred over those that omit this. Finally, an implication of this framework is that most surficial silicic magmas will have undergone the transition to fracture-hosted open-system dynamics in which parts of the gas phase is strongly decoupled from the residual magma.

Methods
A high-temperature, high-load uniaxial press for magma deformation. We use a uniaxial compression press built by Voggenreiter GmbH. A linear variable differential transformer (LVDT) with a resolution of 10 −6 m and a total 150 mm travel range is used with a hydraulic motor system to control the motion of the top piston. Within the working range 8.3 × 10 −7 ms −1 to 1.0 × 10 −2 ms −1 , the displacement rates are well controlled using an Instron® electronic control interface and their WaveMatrix™ dynamic testing software. The force on the bottom piston is monitored using a Lorenz Messtechnik GmbH K11 load cell with a range of 300 kN and an approximate accuracy of 0.05% of the measured force. Surrounding a length of the pistons up to 15 times the length of the sample is a three-zone split furnace from Gero GmbH with a limiting temperature of 1100°C accurate to within 2°C. Thermal gradients were minimized by thorough insulation and temperature was recorded at several points and inside the samples with K-type thermocouple arrays. At equilibrium, the hot zone is 0.12 m long. Outside of two cooling rings-and therefore protected from temperature effects-are two piezoelectric AE broadband transducers with 125 kHz central frequency. The AE signals are fed through a 40 dB buffered pre-amplifier, before being recorded at a Richter data acquisition system from Applied Seismology Consultants, which has a 20 MHz sample rate.
Syn-deformation acoustic emission testing. AE event onsets were triggered and recorded automatically from the continuous acoustic data using an adaptation of an autoregressive-Akaike-Information-Criterion (AR-AIC) picker 50 . The AR-AIC picker follows a defined workflow: (i) detection of the onset of a waveform above the baseline using an STA-LTA detector, (ii) de-noising of the acoustic signal, and (iii) AIC computation where the minimum indicates the arrival time. The STA-LTA window was set to 1 and 20 ms, respectively and the STA/LTA threshold was 2. The amplitude in dB of each single event was also computed (based on a resistance reference standard value of 10 kΩ). The amplitude threshold was set at 68 dB to eliminate noise further. Finally, the acoustic emission data was cut at time of contact and time of failure (or end of the experiment), consistent with the mechanical data, allowing excellent synchronization.
Materials and samples. Samples were cored from obsidian blocks from Hrafntinnuhryggur, Krafla, Iceland, or pre-prepared from extruded soda-lime-silica glass rods from Schott GmbH, Germany. All samples were cut to 20 mm diameter and an initial 40 mm height, termed L 0 . The press was set to use a constant velocity v of the upper piston, which was tuned so that the axial strain rate in the sample _ γ a ¼ v=L 0 covered a wide range. Using Eq. 1, and the details in the main text, we predicted λ r , which allowed us to predict what range of _ γ was required to conduct experiments that straddle Wi = 1 by a few orders of magnitude either side for each temperature. We used temperature 550-800°C for the soda-lime-silica glass, and 750-1000°C for the obsidian. In the Wi ( 1 regime, we used the Gent 51 method   Fig. 5 The consequence of viscoelasticity for the emplacement behaviour of magmas. a At Wi < 10 −2 , the maximum shear stress is τ m = G ∞ Wi, and at Wi > 10 −2 , τ m = 10 −2 G ∞ , which in the case of magmatic liquids is τ m ¼ Oð10 2 Þ MPa, consistent with a reanalysis of data for synthetic 28,[44][45][46] and natural silicate liquids 35 . Here τ m is computed from σ m using the Trouton ratio for cylindrical geometry where appropriate. b The acoustic event number confirming fracture-formation during deformation of magmatic liquids at Wi ≳ 3 10 À3 and dominant acoustic activity at Wi = 0.01. c The viscoelasticity of a rhyolitic liquid oscillated at forcing frequency ω to low strain amplitudes, yielding a real G′ and imaginary G″ component of the elastic modulus 47, 48 as a function of the Deborah number De = ωλ r and showing the onset of elasticity at De = 10 −2 , consistent with all data (Figs. 2 and 4) to confirm our calculated μ value, through a direct measurement, validating the method and providing confidence when moving to experiments at Wi > 1. For the soda-lime-silica glass, we found that our measured μ agreed well with the temperature-dependent value provided by Schott GmbH, and for the obsidian, we found that μ agreed well with the Hess and Dingwell 31 model for C H 2 O ¼ 0:14 wt. %, consistent with measurements presented for the same material from the same location 17 . For the obsidian, this confirmed that degassing by diffusion out of the sample during the test was negligible, and we confirmed that bubbles did not nucleate during the test. Mechanical and acoustic data was analysed and visualized using custom-made Python™ algorithms or the ObsPy package 50 . For data compiled from the literature, we carefully applied the Trouton ratio correction for any value of σ, to extract an estimate for τ (used in Fig. 5). Otherwise, we used published values of μ and the general approximation for G ∞ of 10 10 Pa quoted in the main text.

Data availability
All processed and raw data is available from the corresponding author on request.