Experimental constraints on degassing and permeability in volcanic conduit flow

Abstract

This study assesses the effect of decompression rate on two processes that directly influence the behavior of volcanic eruptions: degassing and permeability in magmas. We studied the degassing of magma with experiments on hydrated natural rhyolitic glass at high pressure and temperature. From the data collected, we defined and characterized one degassing regime in equilibrium and two regimes in disequilibrium. Equilibrium bubble growth occurs when the decompression rate is slower than 0.1 MPa s⁻¹, while higher rates cause porosity to deviate rapidly from equilibrium, defining the first disequilibrium regime of degassing. If the deviation is large enough, a critical threshold of super-saturation is reached and bubble growth accelerates, defining the second disequilibrium regime. We studied permeability and bubble coalescence in magma with experiments using the same rhyolitic melt in open degassing conditions. Under these open conditions, we observed that bubbles start to coalesce at ~43 vol% porosity, regardless of decompression rate. Coalescence profoundly affects bubble texture and size distributions, and induces the melt to become permeable. We determined coalescence to occur on a time scale (~180 s) independent of decompression rate. We parameterized and incorporated our experimental results into a 1D conduit flow model to explore the implications of our findings on eruptive behavior of rhyolitic melts with low crystal contents stored in the upper crust. Compared to previous models that assume equilibrium degassing of the melt during ascent, the introduction of disequilibrium degassing reduces the deviation from lithostatic pressure by ~25%, the acceleration at high porosities (>50 vol%) by a factor 5, and the associated decompression rate by an order of magnitude. The integration of the time scale of coalescence to the model shows that the transition between explosive and effusive eruptive regimes is sensitive to small variations of the initial magma ascent speed, and that flow conditions near fragmentation may significantly be affected by bubble coalescence and gas escape.

Appendix

Our simple model of magma ascent through a conduit before fragmentation is based on the approach first developed by Wilson (1980). We closely follow Mastin and Ghiorso (2000)’s solution, so that when the decompression rates are small enough to allow degassing to occur in equilibrium, the output of our model is undistinguishable from the one of Mastin and Ghiorso (2000). The novelty we introduce is the effect of disequilibrium degassing, for which we develop empirical relations between pressure, decompression rate, and porosity for the regimes of bubble growth, and their respective range of applicability [Eqs. (A10)–(A15)]. The model is one-dimensional and thus flow properties are averaged across the cross-sectional area of the conduit at any given depth. The flow is considered isothermal and homogeneous (i.e., bubbles rise at the same speed as the magma), the conduit is a vertical cylinder with impermeable rigid walls, and the gas phase is H₂O. Although those assumptions are not strictly applicable to the natural system, they have been widely used to provide a reasonable approximation of magma velocities in a conduit during an eruption. We consider that the mass flux Q remains constant, so that:

$$ Q = \rho \pi r^{2} v $$

(A1)

where ρ is the mixture density, v its speed, and r the conduit radius. The flow below fragmentation is laminar, thus the steady-state conservation of mass and momentum read (e.g., Woods 1995):

$$\frac{{d{\left( {\rho v} \right)}}} {{dz}} = 0$$

(A2)

$$\rho v\frac{{dv}} {{dz}} = - \frac{{dP}} {{dz}} - \rho g - f$$

(A3)

where P is the magma pressure, g is the acceleration of gravity (9.81 m s^-2), z is the height within the conduit, and f a friction factor given by (Mastin and Ghiorso 2000):

$$f = \frac{{8\mu v}} {{r^{2} }} + 0.0025\frac{{\rho v^{2} }} {r}$$

(A4)

where μ is the magma viscosity. The presence of bubbles in the magma causes viscosity to change as a function of gas fraction and shear (e.g., Pal 2003), but these complex changes cannot be fully evaluated by a model that uses depth-averaged quantities. Thus, following Mastin and Ghiorso (2000), we use the pure melt viscosity in Equ. (A4) (Hess and Dingwell 1996). When bubbles are in equilibrium with the liquid, the solubility law gives the weight fraction x of water remaining in the liquid (e.g., Woods 1995):

$$ x = k{\sqrt P } $$

(A5)

with Henry’s constant κ=3.44×10⁻⁶ kg^–0.5 m^0.5 s determined from the best-fit parameters of the second-order polynomial regression of our experimental set in equilibrium (Table 1). The difference with the value determined by Mangan and Sisson (2000) for the same PCD rhyolite (4.15×10⁻⁶) is mainly due to the regression equation they used (non-zero intercept second order polynomial regression). The mixture density ρ is calculated using the perfect gas law:

$$\rho = \alpha \rho _{g} + {\left( {1 - \alpha } \right)}\;\rho _{l} = \alpha \frac{M} {{RT}}P + {\left( {1 - \alpha } \right)}\;\rho _{l} $$

(A6)

Where ρ _l is the liquid density, ρ _g is the gas density, α its volume fraction, M is the water molecular weight (18×10⁻³ kg mol⁻¹), T is the mixture temperature, and R is the universal gas constant (8.3144 J mole⁻¹ °K⁻¹). Combining the derivative of Eq. (A6) with Eqs. (A1)–(A3) leads to:

$$\frac{{dP}} {{dz}} = {\left[ {\frac{{d\alpha }} {{dz}}v^{2} {\left( {\frac{{MP}} {{RT}} - \rho _{l} } \right)} - \rho g - f} \right]}{\left[ {1 - \alpha \frac{{v^{2} M}} {{RT}}} \right]}^{{ - 1}} $$

(A7)

The evolution of pressure in the conduit is thus controlled by the degassing rate dα/dz, which has to be determined for each experimentally determined regime: equilibrium, slow growth, and fast growth. In the equilibrium regime, the gas volume fraction α _equ (Fig. 6) is calculated from the difference between the maximum amount of water that can be dissolved in the magma at a given pressure P and the original water content at the initial pressure P ₀ :

$$\alpha _{{equ}} = \frac{\beta } {{MP + \beta }}\;{\text{with}}\;\beta = RT\rho _{l} \kappa {\left( {{\sqrt {P_{0} } } - {\sqrt P }} \right)}$$

(A8)

The derivative of Eq. (A8) gives the equilibrium-degassing rate:

$${\left( {\frac{{d\alpha }} {{dz}}} \right)}_{{equ}} = \frac{{dP}} {{dz}}M{\left( { - \beta - \frac{{{\sqrt P }}} {2}} \right)}{\left( {MP + \beta } \right)}^{{ - 2}} $$

(A9)

The slow growth regime is a linear function of α and P (Fig. 5; see Table 3 for empirical constant value):

$${\left( {\frac{{d\alpha }} {{dz}}} \right)}_{{slow}} = A_{z} \frac{{dP}} {{dz}}$$

(A10)

In the fast growth regime, Δα is given by the linear regression for each decompression rate (Fig. 5):

$$\frac{{d(\Delta \alpha )}} {{dz}} = A_{s} \frac{{dP}} {{dz}}$$

(A11)

where A _s varies as a function of the decompression rate according to:

$$A_{s} = a_{2} \frac{{dP}} {{dt}} + b_{2} $$

(A12)

Hence the derivative of α with depth becomes:

$${\left( {\frac{{d\alpha }} {{dz}}} \right)}_{{fast}} = {\left( {\frac{{d\alpha }} {{dz}}} \right)}_{{equ}} - A_{s} \frac{{dP}} {{dz}}$$

(A13)

Note the differential form of Eq. (A10) allows a melt to reach atmospheric pressure super-saturated. We also need to determine the boundaries between each regime of degassing. Experimental data (Fig. 3) show that the critical decompression rate between the equilibrium and the slow growth regimes varies little between high and low pressure, so we used a unique median value of 0.15 MPa s⁻¹. The boundary between the slow growth and the fast growth regimes is reached when the difference between the disequilibrium (α) and the equilibrium (α _equ ) porosities Δα is greater than a critical value Δα _cr (Fig. 6). The relation between Δα _cr and the decompression rate is given by the regression:

$$\ln (\Delta \alpha _{{cr}} ) = a_{1} \ln {\left( {\frac{{dP}} {{dt}}} \right)} + b_{1} $$

(A14)

Experimental data were used for the regression when available and interpolated between the slow and fast growth curves when no data point fell on Δα _cr (Figs. 5 and 13). The nonlinear form of Eq. (A14) has been preferred over a linear form to ensure a realistic (asymptotic) behavior of the critical porosity at large decompression rates. In disequilibrium degassing, the melt viscosity is calculated from the water content at P _right (Eq. (A5), Fig. 6), because the amount of water in bubbles at a given pressure corresponds to the equivalent equilibrium value α _right . P _right is obtained using Eq. (A8):

Fig. 13

Nonlinear regression of Δα _cr (maximum porosity deviation from equilibrium degassing) as a function of decompression rate. Circles are data from Gardner et al. (1999) and squares are data from this study

$$MP_{{right}} \alpha _{{right}} + TR\kappa \rho _{l} {\left( {\alpha _{{right}} + 1} \right)}{\left( {{\sqrt {P_{{right}} } } - {\sqrt {P_{0} } }} \right)} = 0$$

(A15)

The degassing behavior of Eq. (A7) is thus given by (dα/dz)_equ [Eq. (A9)] if dP/dt<0.15 MPa s⁻¹; (dα/dz)_slow [Eq. (A10)] if dP/dt>0.15 MPa s⁻¹ and α<α _cr; and (dα/dz)_fast [Eq. (A13)] if dP/dt>0.15 MPa s⁻¹ and α>α _cr with the porosity α being calculated from the previous distance step. Equation (A7) is solved using a fourth-order Runge-Kutta algorithm with a constant distance step of 0.5 m below 50 vol% porosity and 0.05 m above.