## Abstract

The growth of lava domes is often associated with cyclic variations of ground deformation, seismicity and mass flux of gas and magma. We present a model of cyclic volcanic activity which is controlled by degassing of supersaturated magma, magma flow into the conduit, gas escape from the permeable magma, deformation of the conduit walls and the friction between the walls and the plug at the top of the conduit. When the difference between magma pressure and ambient pressure exceeds the static friction, motion begins, bubbles expand and overpressure relaxes. Bubble expansion builds permeability, allows gas escape and faster depressurization. Depressurization and crystallization of the magma build supersaturation and gas diffusion from melt to bubbles. Gas flux into bubbles and magma flux from the chamber act to increase pressure. The rate of extrusion is controlled by the gas pressure, driving the motion, and by the rate- and state-dependent friction along shear zones between the plug and the host rock. When the magma overpressure drops to the dynamic strength of the slip surfaces, the plug sticks and blocks the vent. As bubble volume is now constant, exsolution of gas from the supersaturated melt leads to pressurization and begins a new cycle.

The growth of lava domes is commonly characterized by repetitive cyclic patterns of ground deformation, degassing, seismicity, dome extrusion and explosive eruptions, with timescales typically of hours to days. Such cyclic patterns have been observed at the Soufrière Hills Volcano (SHV), Montserrat (Voight *et al.* 1999), Mount Pinatubo, Phillipines (Denlinger & Hoblitt 1999) and Merapi, Indonesia (Voight *et al.* 2000). At SHV, the cycles have been observed from tilt, seismic and gas emission data (Green & Neuberg 2006; Voight *et al.* 1998; Voight *et al.* 1999; Watson *et al.* 2000). They are also closely associated with episodic lava extrusion and repetitive vulcanian explosions (Connor *et al.* 2003; Druitt *et al.* 2002; Formenti *et al.* 2003).

## Soufrière Hills Volcano (SHV) cycles

At the SHV, two tiltmeters were stationed in 1997 within 400 m of the growing andesitic dome and several months of observation yielded a detailed time series of data (Denlinger & Hoblitt 1999; Voight *et al.* 1998, 1999). The main observations were of ground tilt, seismicity and visual phenomena. The geophysical data were also supplemented with scans of gas emission (SO_{2}) using the COSPEC instrument (Edmonds *et al.* 2003; Watson *et al.* 2000) and recording of vulcanian explosions (Druitt *et al.* 2002). Since the data have already been described in several publications, here we synthesize the key observations (Fig. 1).

Each cycle begins with an inflation of ground tilt often accompanied by a swarm of shallow earthquakes, which have been classified as hybrid events on SHV (Miller *et al.* 1998). Hybrid earthquakes are characterized by emergent signals, an absence of detectable S wave arrivals and energy predominantly in the 1–5 Hz range (Neuberg *et al.* 1998; Neuberg & O'Gorman 2002). Their sources have been located at shallow depths (<2 km) and there appears to be a continuum of hybrid events from mixed high and low frequency components to pure, almost monochromatic long period events. Hybrid earthquakes begin when the tilt goes through its inflection point, build-up to a maximum as the peak in the tilt cycle is approached and decline as the peak is passed with hybrid earthquakes being few or absent in the deflation phase of a tilt cycle (Green & Neuberg 2006). Tilt cycles can occur without hybrid swarms. Green & Neuberg (2006) also showed that the seismicity consists of several families of distinctive earthquakes, each family having similar spectral characteristics.

The peak in tilt and the following deflationary phase are associated with a variety of volcanic phenomena. The most commonly observed phenomenon is that the deflation is accompanied by onset of elevated rockfall activity of the dome and sometimes generation of pyroclastic flows. The rockfall activity in the deflation periods and the infrequent occurrence of rockfall seismicity during the inflation period were taken by Voight *et al.* (1999) to indicate much faster lava extrusion rates during deflation than during inflation. They further proposed stick-slip behaviour on the basis of these observations and this was followed by some heuristic models by Denlinger & Hoblitt (1999) and Wylie *et al.* (1999) and seismic observations consistent with stick-slip behaviour (Green *et al.* 2006). In some cases, the maximum tilt and onset of deflation were marked by intense ash and gas emissions. Ash venting commonly occurs for up to an hour after the tilt peak. Watson *et al.* (2000) found that in June 1997 SO_{2} fluxes fluctuated over a cycle with the peak SO_{2} emission occurring about an hour after the peak in ground tilt. Edmonds *et al.* (2003) noted that during December 1999 and January 2000, SO_{2} emissions appeared to correlate with relative seismic amplitude measurements (RSAM).

The most dramatic phenomena, however, were intense vulcanian explosions which occurred at the peak in ground tilt in early August 1997 (Voight *et al.* 1998). Unfortunately, the destruction of the tiltmeter by one of these explosions prevented further observations, but it is inferred that 74 repetitive vulcanian explosions occurred in the period between 21 September to 22 October 1997, with an average cycle of ten hours (Druitt *et al.* 2002). Detailed analysis of the explosions (Druitt *et al.* 2002; Formenti *et al.* 2003) indicated that each explosion evacuated several hundred metres of the upper conduit in a few tens of seconds. Statistical distribution of repose periods between vulcanian explosions (Connor *et al.* 2003) can be well described by a log logistic model, which requires the influence of two counteracting processes. Connor *et al.* (2003) proposed that these processes were the build-up of gas pressure in supersaturated magma and the escape of gas from magma by permeable flow. The former increases gas pressure and the latter decreases gas pressure. Jacquet *et al.* (2006) found that the same data treated as a time series showed memory effects. They interpreted the 50–60 hour memory in the time series as the repeated decompression of magma batches during their ascent in the conduit by successive vulcanian explosions. Clarke *et al.* (2007) studied the petrological characteristics of the ejecta from vulcanian explosions and found that they originated at pressures of up to 40 MPa, consistent with evacuation of the conduit to depths of up to 2 km in the largest vulcanian explosions. They observed that dense clasts gave low pressure estimates while the density of vesicular clasts increased with estimated pressure. Their results were interpreted in terms of a dense degassed cap of lava overlying a gas-rich and pressurized column of magma prior to an explosion.

### The time period of the cycles

The time period of the tilt cycles varied between four and 36 hours (Voight *et al.* 1999), as measured by the intervals between peaks. The period and amplitude could remain stable for several days, but eventually changed or broke down. Data from June to August 1997 showed that the period and amplitude of the tilt could change systematically over periods of several weeks. After an abrupt onset of high amplitude and short period tilt cycles, the period gradually increased and amplitude decreased (Voight *et al.* 1999). Concurrently, the seismicity would decline and eventually almost cease. Such an episode would then be interrupted by another abrupt decrease of period and increase in amplitude, accompanied by a marked increase in seismicity. Watts *et al.* (2002) showed that these abrupt changes were associated with the extrusion of a new lobe of lava. Costa *et al.* (2007) interpreted these longer-term cycles in terms of dynamic flow in a dyke with the timescale of several weeks being controlled by the elastic deformation of the dyke walls. This inference is supported by analysis of the trends in the tilt data over the six to seven-week cycles, which are consistent with a model of a pressurized dyke with a NW–SE orientation from a pressure source at 880–1230 m depth (Hautmann *et al.* 2008).

In detail, tilt cycles are characterized by roughly similar periods of inflation and deflation. In most cases, the inflation deformation was largely recovered during the deflation, indicating almost fully recoverable strain. Voight *et al.* (1999) estimated a depth of more than 400 m below the dome as the top of the shallow pressure source. For the majority of tilt cycles, the main phenomenon of the deflation period was the onset of numerous rockfalls. This led Voight *et al.* (1999) to propose that the deflation was associated with active dome extrusion. Alternatively, Green *et al.* (2006) proposed that the tilt cycles were better explained by shear stresses acting along a conduit over a distance of *c*. 1000 m, rather than a single pressure source. Their model avoids having to invoke unreasonably high pressures and a large geologically implausible volume for the pressure source.

### Models

It has been suggested that pressure build-up in the magma within the conduit is caused by: (i) variations of viscosity along the conduit (Connor *et al.* 2003; Jaquet *et al.* 2006; Melnik & Sparks 1999; Sparks 1997; Voight *et al.* 1999), (ii) influx of magma from the chamber to a conduit that is partly blocked from above by the resistance of magma (Denlinger & Hoblitt 1999; Wylie *et al.* 1999) and by (iii) volatile diffusion from the supersaturated melt into the bubbles (Sparks & Young 2002; Watson *et al.* 2000). The volatile diffusion is augmented by growth of microlite crystals from degassed melt which can achieve pressure increases of 1 MPa per 1% of crystallization in a volume-constrained, sealed system (Sparks 1997; Stix *et al.* 1997). All these cases invoke the presence of a low permeability plug of highly viscous and partially crystallized magma at the upper part of the conduit, which confines the volume of the conduit. Relaxation of the overpressurized magma has been attributed to magma extrusion and sliding of the plug along the conduit walls when the yield strength between the plug and the wall rocks is exceeded (Connor *et al.* 2003; Jaquet *et al.* 2006; Sparks 1997).

In this paper, we present a model that accounts for the cyclic activity of a volcano using the data from the SHV as an example. We estimate the role of different processes controlling the rates of pressure build-up and relaxation in a supersaturated magma capped with a plug that sticks when pressure falls and slips when pressure is high.

## Observations indicating gas pressurization and stick-slip motion in lava domes

### Gas emission

Strong gas emission during the deflation periods implicates gas escape and depressurization (Watson *et al.* 2000). The association of vulcanian explosions with the inflation maxima is consistent with gas pressurization as the main cause of inflation. In general, the repetitive vulcanian explosions are thought to be caused by the same basic cyclic process as the tilt cycles (Druitt *et al.* 2002; Jaquet *et al.* 2006).

### Observations from SHV and other lava domes

At SHV observations of lava dome behaviour, morphology and structure are consistent with a stick-slip style of extrusion during exogenous growth (Sparks *et al.* 2000). Common morphologies of the lava are spines and shear lobes bounded by cylindrical, relatively smooth surfaces (Cashman *et al.* 2007; Iverson *et al.* 2006; Nakada & Motomura 1999; Sparks *et al.* 2000; Watts *et al.* 2002), which can be interpreted as faulted surfaces rooted in the upper conduit (Tuffen & Dingwell 2005). They indicate that the upper conduit wall is a zone of localized deformation and that the extrusion process can be seen as faulting process acting on a rigid plug of degassed magma which occupies the upper conduit. Blocks ejected during the 16 September 1996 subplinian explosive eruption at SHV gives some insights into the conditions in the upper conduit. Robertson *et al.* (1998) documented the density and petrology of ejecta. Blocks from the lava dome and possibly uppermost conduit have low water content and high density (2200–2400 kg/m^{3}) and are consistent with an uppermost region of degassed lava. Rather dense pumiceous blocks (1200–2000 kg/m^{3}) with low water content (Harford & Sparks 2001) may represent deeper zones of bubbly magma. Pressure estimations based on the shooting distances of ballistic blocks indicate pressures up to 25 MPa for these ejecta (Robertson *et al.* 1998) consistent with the interpretation that they represent cleaning out of the uppermost conduit. Rare, but important, rock types in the ejecta are highly sheared vesicular blocks with characteristic cataclastic textures of ground hornblende and plagioclase phenocrysts (Fig. 2). These samples have not been described previously, but are important in the context of this study as they are thought to represent fault rock formed at the conduit wall by stick-slip behaviour.

Direct evidence for plug ascent along cylindrical shear surfaces was recently provided by Bluth & Rose (2004), who observed gas and ash venting from a ring-shaped set of fractures at the summit of Santiaguito volcano. They interpret the correlation between ash bursts and measured extrusion rates as an indication for incremental plug flow. During an eruption, the magma ascends a few centimetres, with most of the shear localized along the conduit walls, about 50 m in diameter. Time interval between venting events is of the order of tens of minutes.

Plug ascent along localized shear surfaces was also observed in the recent dome extrusion at Mount St. Helens (Cashman *et al.* 2007; Iverson *et al.* 2006) and at Unzen volcano (Nakada & Motomura 1999). Seismicity associated with these extrusions was related to stick-slip motion of the plug. The period of the seismic drumbeats is of the order of minutes (Iverson *et al.* 2006), much less than the period of tilt cycles discussed in the present paper (see discussion).

Tuffen *et al.* (2003) and Tuffen & Dingwell (2005) described fault textures in rhyolites from Iceland and interpreted them in terms of multiple seismogenic cycles. They suggested that each cycle includes stick-slip motion, which opens pathways for gas escape.

The final matter is to consider decompression of the ascending magma during a cycle. In the case of SHV, the time-averaged extrusion rates are well constrained at 4–7 m^{3}/s in the May–August 1997 period when the tilt cycles were documented (Sparks *et al.* 1998). At a characteristic average extrusion rate of 5 m^{3}/s, conduit area of 700 m^{2} (Melnik & Sparks 2002) and average tilt cycle of 12 hours, the vertical displacement of the extruded magma is estimated at 300 m. Taking magma density at 2300 kg/m^{3}, then the change in magmastatic pressure for a parcel of ascending magma is estimated at *c*. 7 MPa. Displacement may be larger if the cross-sectional area of the conduit narrows with depth (e.g. Couch *et al.* 2003*b*; Jaquet *et al.* 2006) suggesting that decompression may actually be larger. Following Voight *et al.* (1999), it is surmised that ascent and decompression occurred mainly during deflation and that flow rates were approximately twice the average. We conclude that decompression during pulsatory ascent associated with a single cycle may lead to supersaturation of several MPa, provided that ascent is sufficiently fast, so that the system cannot maintain equilibrium. Crystallization of the melt also contributes to the build-up of supersaturation.

## The model of cyclic activity in domes

### Overview

We consider a cylindrical conduit filled with magma and topped by a rigid plug (Fig. 3). The two major pressurization processes considered are volatile diffusion from melt to bubbles and magma flow from the feeding source region into the blocked conduit. The drive for volatile diffusion is the supersaturation that develops due to decompression of volatile-bearing magma and due to crystallization which leads to enrichment of volatiles in the residual melt. As the magma pressurizes, the surrounding conduit walls deform, allowing magma expansion and reduction in the rate of pressurization. With efficient pressurization, magma pressure beneath the plug-base surpasses the strength of the plug and initiates motion. As the magma flows, bubbles expand and pressure decreases. In addition, bubbles merge and coalesce and new pathways are opened for gas escape and further depressurization. When magma pressure falls below the dynamic friction of the plug, motion stops and a new cycle begins. For low depressurization rate, a steady-state plug motion is expected.

### Pressures and forces driving the plug

The largest pressure component of the magma in the conduit is associated with gravity and is termed the magma static pressure, *P*_{ms}. This pressure at the base of the plug (at *Z*_{pb}) is:
1

Using the average density, the equation simplifies to *P*_{ms}=ρ*gZ*_{pb}.

The ambient magma pressure just below the depth of the plug base, *P*_{a}, is a sum of the weight of the plug, *P*_{ms}, and additional pressure driving the plug motion, τ_{dr}:
2

In slowly ascending bubbly magma, the gas pressure, *P*_{g}, closely approaches the ambient pressure (Lensky *et al.* 2004). Hence, the driving pressure may be well approximated by the difference between the gas pressure and the magma static pressure:
3

Volatile diffusion and magma flow from the deep feeding source lead to increase of the driving pressure. The pressurization stage comes to its end when the stress conditions meet the criterion for the onset of slip and plug motion, τ_{dr}≥*P*_{mo}. Following the Byerlee law (Byerlee 1967), slip starts when the ratio between the driving stress, *P*_{mo}, and normal stress, σ_{N}, overcomes the internal friction, *f*_{mo}. Typical values of *f*_{mo} range between 0.6 and 0.8 for granite samples (e.g. Marone 1998).

The criterion for plug arrest, *P*_{ar}, is the decrease of τ_{dr}/σ_{N} below the dynamic friction, *f*_{ar}, which in seismological modelling of the tectonic earthquake source is taken as 0.2–0.3 (BenZion & Rice 1993). In this study, we adopted the values of *f*_{mo}=0.6 and *f*_{ar}=0.3. For a first order approximation, the horizontal normal stress (σ_{N}) component in an elastic plug is proportional to the vertical stress, ρ · *g* · *z*, multiplied by the ratio ν/(1−ν) (Jaeger & Cook 1969) where ν is the Poisson's ratio.
4

The Poisson ratio varies from 0 to 0.5 with typical value of *c*. 0.3. Hence, we assume the horizontal stress (normal stress) is *c*. 0.5 · ρ · *g* · *z*, and the average value along the plug is 0.25 · ρ · *g* · *Z*_{pb}. Finally, using the static and dynamic friction values, and a plug that is 500 m long, the onset of motion occurs at driving pressure, *P*_{mo}=13.0 MPa, and the plug is arrested at *P*_{ar}=12.1 MPa. These values are used for all model realizations presented in this study.

### Pressurization

#### Rigid conduit (including diffusion and magma flow)

We consider a batch of magma in the conduit below the plug. The magma is considered as a suspension of equally sized and evenly distributed gas bubbles of radius *R*, each surrounded by a spherical shell of melt of outer radius *S* (Proussevitch *et al.* 1993). The shells partially overlap so that the volume fraction of gas in the magma, α, equals the volume of fraction of the bubble in a single cell (bubble+melt shell). No mass is transferred between cells, and the pressure at the cell boundaries is the ambient pressure of the suspension (Lyakhovsky *et al.* 1996). The gas pressure within the bubble, *P*_{g}, is in equilibrium with the dissolved gas at the bubble wall, *C*_{R}=*C*(*r*=*R*), approximated by Henry's solubility law:
5
where *K*_{H} is Henry's (for silicic magmas *n*=2 and *K*_{H}=4.11 · 10^{−6} Pa^{−0.5} (Burnham 1975)).

When the plug sticks and a new cycle begins, ambient and gas pressure are low, but the magma is supersaturated. In this model, the conduit walls are considered rigid and the plug is blocked, hence, magma is confined to a limited volume and bubbles cannot expand. Thus, volatile exsolution into the bubbles leads to pressurization of the bubbles and the magma. The time for pressure build-up due to exsolution is controlled by the kinetics of water diffusion, as described by the equation of mass balance of the gas in the bubble:
6
where ρ_{g} and ρ_{m} are the densities of gas and melt and *D* is the diffusivity of water in the melt. Under the relevant conditions (*T*=820–880 °C, *P*=10^{5}–10^{7} Pa), the gas is approximated as an ideal gas:
7
where *P*_{g} is the gas pressure, *M* is the molecular weight of water, *G* the gas constant and *T* is the temperature of the gas and melt (see values in Table 1). In the case of lava domes, with their high viscous melt and low ascent rate, the diffusive flux of volatiles into bubbles can be regarded as quasi-static, (Lensky *et al.* 2004; Lyakhovsky *et al.* 1996; Navon *et al.* 1998):
8
where *C*_{R} and *C*_{S} are the volatile concentration at the bubble-melt interface and at the shell boundary. Initially, the concentration of volatiles away from the bubble is equal to the initial concentration, *C*_{S}=*C*_{i}. Later, as water diffuses into the bubbles, *C*_{S} declines. During the numerical procedure *C*_{S} is adjusted so that the total volatile mass in the finite cell (bubble+melt shell) is conserved and (∂*C*/∂*r*)_{r=S}=0, while the concentration gradient is approximated by equation 8 (e.g. Lyakhovsky *et al.* 1996). Substituting equations 7–8 into equation 6 yields the equation for the rate of pressure build-up:
9

Initially, *C*_{R}<*C*_{S}, and with time *C*_{R} increases while *C*_{S} decreases until equilibration is achieved when *C*_{R}=*C*_{S} and diffusive flux ceases. The timescale can be estimated if *C*_{S} is assumed to remain constant. In this case, gas pressure exponentially approaches the equilibrium gas pressure, *P*_{0}, (*P*_{g}=*P*_{0}−(*P*_{0}−*P*_{ar})*e*^{−t/τ}) with a timescale .

The gas pressure, *P*_{0}, and density, ρ_{g,0}, at equilibrium are calculated following Lensky *et al.* (2006). Onset of plug motion occurs when *P*_{g} exceeds *P*_{mo} after time .

The above approximation underestimates the diffusive flux at the initial stage of pressurization, since *C*_{S} actually decreases with time, and hence, the onset of plug motion is achieved earlier than this analytical approximation. For typical values for SHV diffusion coefficient *D*∼1 µm^{2}/s and bubbles radius of 100 µm, *C*_{0}ρ_{m}/ρ_{g,0}∼1, leading to τ-values of the order of hours, which is the upper limit for the duration of diffusive controlled pressurization.

Since the magma in the upper conduit is evolved and partially degassed, we estimate bubble size based on observations on clasts ejected during vulcanian eruptions. The bubble size chosen (130–220 µm) are between the size of observed vesicles in SHV dome and pumice samples (Clarke *et al.* 2007; Formenti & Druitt 2003; Murphy *et al.* 2000) and the size of bubbles in the magma chamber, for which a radius of *c*. 45 µm was estimated by Voight *et al.* (2006) based on deformation data. Vesicularity varies between 10–30%, leading to bubble number density of *c*. 10^{10} m^{−3}, consistent with these observations. The temperature is kept constant at 850 °C (Murphy *et al.* 2000) and the water content of the melt varies in the range 1.5–1.6 wt%. Since the melt (but not the bulk rock) is rhyolitic, we calculate diffusion coefficients using the equation of Zhang & Behrens (2000). The diffusivity in these conditions is around 4×10^{−12} m^{2}/s. Hence the pressurization timescale (*R*^{2}/*D*) is of the order of a few hours; smaller bubbles of 50 µm (number density of 10^{12} m^{−3}) will reduce the timescale to an hour.

Pressure build-up is very sensitive to the initial volatile content (*C*_{i}), which reflects the potential energy for diffusive pressurization. Figure 4a presents the pressure build-up due to diffusion of volatiles in a rigid conduit (i.e. – constant volume) with water content at the range of *C*_{i}=1.5–1.6 wt%. Higher *C*_{i} means that the last equilibration was at greater depth and that decompression rate was higher.

Another process that influences the pressure in the conduit is the flow of magma into, or out of, the conduit. Magma is supplied from the feeding source region, but when pressure in the conduit increases, backflow into the chamber is also possible. This process is controlled by the difference between the magmatic overpressure in the source, *P*_{op}, and the ambient pressure. Assuming that the flux of the magma, *Q*, is proportional to the pressure difference with a flow rate coefficient *K*_{mf}, and that *P*_{a}∼*P*_{g},
10

For incompressible melt, the rate of gas volume change is equal to the flux of magma. Combining equation 10 with equation of state of the ideal gas equation 7, the rate of pressurization governed by the magma flow is:
11
where *H*α stands for the total volume of gas bubbles per unit area of the conduit cross-section. While *K*_{mf} is an important parameter that controls the dynamic of a volcano, it is not well confined by observations. High *K*_{mf} values mean that the magma can easily flow into the conduit when *P*_{op}>*P*_{g} and enhances pressurization, but it also means easy backflow when diffusion-controlled pressurization is efficient and there is a smaller chance to move the plug. At low *K*_{mf}, flow into and out of the conduit is limited and degassing remains the major contributor to pressure build-up in the conduit. Following Melnik & Sparks (2005), we assume that *P*_{op} does not change significantly over the period of a cycle (hours). Integrating for the total mass of the magma extruding the conduit per cycle, and approximating *P*_{op}−*P*_{g}∼*P*_{mo}−*P*_{ar} of the order of 1 MPa, enables estimation of the value of *K*_{mf}. For SHV, we take 300 m ascent during a cycle of about 10 hrs, which gives *K*_{mf}∼300 m/ (3.6×10^{4} s)/10^{6} Pa∼10^{−8} m/Pa/s. For *K*_{mf} value estimated using equation 10, the characteristic time of the flow-controlled pressurization (equation 11) is the ratio between the gas volume in the conduit *H*α and overpressure during magma flow, τ=*H*α/*K*_{mf}/*P*_{g}.

Figure 4b demonstrates the flow-controlled pressurization for different values of magmatic overpressure, *P*_{op}, neglecting the diffusive-controlled pressurization. For these simulations, *P*_{op}=*P*_{ar}+*a*(*P*_{mo}−*P*_{ar}) is chosen to fall between *P*_{mo} and *P*_{ar} with *a*=0.8, or equal to *P*_{mo} (*a*=1), or above *P*_{mo} (*a*=1.2). Higher *P*_{op} means higher final pressure and higher chance for onset of plug motion. In addition, the figure demonstrates the role of the vesicularity, α, on the rate of pressurization due to magma flow. The higher vesicularity is, the slower the pressurization. When α approaches zero, the pressurization time reduces significantly. This may explain the very short pressurization time that Iverson *et al.* (2006) calculated in their model, which neglects vesicularity of magma.

### Deformable conduit (diffusion, magma flow, and crystallization)

When the plug sticks and pressure increases, the conduit walls deform. During the slip phase, the volcanic edifice relaxes due to magma extrusion. In this end member case, we assume that the reversible volume change associated with the tilt of the edifice during the stick stage is equivalent to the volume of magma extruded during the slip stage. If the expansion of the edifice (Δ*V*_{e}) is proportional to pressure, then:
12
where Δ*V*_{em} is the maximum volume change of the volcanic edifice. The ratio between Δ*V*_{e} and the volume of the bubbly region in the conduit (*V*_{co}), at *P*_{ar} prior to its expansion, is the volumetric strain of the magma:
13
and the expansion strain rate is:
14
Combining equations 14 and 9, we obtain:
15
Figure 4c presents the pressurization of magma due to volatile diffusion into the bubbles (no magma inflow) in a deformable conduit. The solution is very sensitive to the amount of tilt. Low tilt (i.e. low deformation) means that the gas flux into bubbles results in pressurization rather than expansion. In this case, a tilt of 1 microradian (μrad) is low enough so that diffusion-driven pressurization alone can initiate plug motion. However, more deformable walls (tilt of 10 µrad) leads to exhausting of volatiles before *P*_{mo} is reached (unless supersaturation is maintained by magma crystallization).

Figure 4d presents the pressure build-up due to magma inflow from the source region to a deformable conduit, neglecting the effect of diffusive influx of volatiles. The solution is very sensitive to the amount of deformation (tilt). Higher tilt leads to slower pressurization, since in addition to compressing the bubbles, the inflowing magma now fills the extra volume created by the deformation of the conduit walls. Since we assume that the overpressure at the source is not affected by the tilt at the top of the conduit, then all solutions will potentially reach that overpressure. If *P*_{mo} is lower than that overpressure, then the plug will move after a few hours. In this solution, the timescale is longer than the previous timescale (compare Fig. 4d with 4a–c).

Petrographic studies of the dome samples and ejecta reveal extensive groundmass crystallization in the upper conduit, particularly at pressures of 30 MPa or less (Clarke *et al.* 2007; Couch *et al.* 2003*a*; Sparks 1997; Sparks *et al.* 2000). Crystallization increases water concentration in the residual melt and enhances supersaturation and degassing. In order to quantify the role of crystallization, we estimate the amount of extra water that can be delivered during crystallization and compare the timescales of crystallization and diffusion: if the former is much faster, then the diffusion timescale remains as is; if it is slower, then the crystallization timescale will govern the pressurization.

As the crystallizing assemblage carries negligible water, the mass of water released due to crystallization (per volume of magma) is Δ*M*_{g}=*C*_{0}ρ_{m}Δβ, where Δβ is the fraction of magma that crystallized. The water concentration of the melt grows by Δ*C*=Δ*M*_{g}/*M*_{melt}=*C*_{0}ρ_{m}Δβ/(ρ_{m}(1−β)), or in other terms Δ*C*/*C*_{0}∼Δβ/(1−β). Assuming β∼0.5, further crystallization of 10% leads to a considerable increase of volatile concentration by 20%. The potential contribution due to crystallization is demonstrated in Figure 4a which shows the effect of the initial water content on degassing. Crystallization shifts the melt to higher curves at it proceeds.

Couch *et al.* (2003*b*) have determined the timescales of crystal nucleation and growth in melt with composition similar to that of SHV. They carried out experiments where samples were decompressed in steps to simulate the episodic magma ascent that would be expected in the stick-slip extrusion cycles. The experiments suggest that during decompression and degassing, the actual crystal growth lags somewhat behind equilibrium crystallization. We assume a similar behaviour and model the growth of the crystal mass during the slip phase as exponential. The equilibrium crystal fraction, β_{e}, follows β_{e}(*P*_{g})=0.115exp((1.2×10^{8}−*P*_{g})/6.67×10^{7}), while the dynamic crystallization lags behind: β_{d}(*P*_{g})=0.115exp((1.2×10^{8}−*P*_{g})/8×10^{7}). During the stick stage, the crystal content caches up with the equilibrium value according to: β_{stick}=β_{d}+(β_{e}(*P*_{slip})−β_{d}) · (1−e^{−t/τ}). Finally, at each step of the numerical solution, the water content, *C*_{S} in equations 9 and 15 is adjusted by factor 1/(1−β) according to the volume decrease of the residual melt during crystallization (ignoring the small water content of crystallizing assemblages).

### Depressurization

#### Plug motion

The pressurization stage comes to its end with the onset of plug motion, which is when the stress conditions meet the criterion for the onset of slip. During the slip stage, the motion of the plug is modelled using the rate- and state-dependent friction approach (e.g. Dieterich 1979; Ruina 1983) that is commonly used to describe stick-slip motion along pre-defined sliding surfaces. The model accounts for the evolution of frictional strength as a function of slip, slip-velocity and state variables that characterize the properties of the sliding surfaces. This simplified approach is widely used in seismology to simulate important aspects of earthquake cycles (e.g. Ben-Zion & Rice 1993; Burridge & Knopoff 1967).

Various laboratory experiments have been conducted in order to study variations of the dynamic friction during block sliding (for review see Marone 1998). An abrupt change in the slip velocity from *v*_{0} to *v* results in an almost instantaneous change of the dynamic friction coefficient. The friction changes from its steady-state value, *f*_{0}, to its dynamic value with an amplitude proportional to the logarithm of the velocity ratio, *a* · ln(*v*/*v*_{0}) (Marone 1998). This change is followed by an evolutionary effect with an opposite sign and amplitude *b* · ln(*v*/*v*_{0}) over a characteristic slip distance *L*, which is less than a millimetre (Marone 1998). The overall change in the quasi-static friction, *f* (also termed steady-state friction), is controlled by the value of *a*−*b*:
16

Stesky (1978), Blanpied *et al.* (1991, 1995) Kilgore *et al.* (1993), Dieterich & Kilgor (1996) and others measured the values of *a* and *b* for granite samples at various pressures and temperatures. If (*a*−*b*) is negative, there is an overall change of velocity-weakening and the frictional response favours unstable sliding (e.g. Scholz 1998). This regime is used to model earthquake dynamics (e.g. Ben-Zion & Rice 1993). On the other hand, if (*a*−*b*) is positive, the overall change is referred to as velocity-strengthening and the frictional response favours stable aseismic sliding. Blanpied *et al.* (1991) reported changes of the (*a*−*b*) sign from negative to positive for wet crystalline granite at about 300 °C and its increase up to 0.04 at 600 °C and 100 MPa. Positive (*a*−*b*) values for semi-brittle behaviour of granites at high temperatures are attributed to the onset of crystal plasticity of quartz (Scholz 1998). Similar semi-brittle behaviour of sufficiently crystalline magma in SHV noted by Sparks *et al.* (2000) and other dacite domes (Smith *et al.* 2001), suggests that positive (*a*−*b*) value should be adopted for the SHV conditions. Extrapolating the experimental results of Blanpied *et al.* (1991) to the typical temperatures in an active dome, high positive values of about (*a*−*b*)=0.1 are expected.

The ratio between the driving stress (equation 3) and the normal stress (equation 4) is now related to the plug velocity using the steady-state friction (equation 16): 17

Rearranging, we obtain the rate of motion: 18

The criterion for plug arrest is the decrease of τ_{dr}/σ_{N} below the dynamic friction, which in a seismological modelling of the tectonic earthquake source is taken as 0.2–0.3 (Ben-Zion & Rice 1993).
19

Motion continues as long as *P*_{g}>*P*_{ar}. Lower dynamic friction leads to a lower arrest pressure, larger pressure amplitude and longer relaxation period.

Plug motion with velocity *v* leads to magma expansion which reduces gas pressure and introduces an additional term to equation 15. The decompression due to expansion is equal to *dP*_{g}/*dt*=−*vP*_{g}/α*H*; this relation reflects the change in gas volume due to plug motion, similar to equations 11 and 14. Accounting for this motion and for the diffusive flux, magma flow and conduit deformation (equation 15), the rate of pressure change is:
20

#### Gas escape

Magma degassing in a conduit is generally attributed to permeable gas flow laterally into the conduit wall (Eichelberger 1995; Jaupart & Allegre 1991), or vertically through the magma column (Boudon *et al.* 1998). Gonnermann & Manga (2003) estimated the timescale for permeable gas flow through the vesicular magma within the conduit of the order of weeks to years. They also speculated that shear-induced fragmentation may create an interconnected fractured network of high permeability. This assumption is consistent with the observed correlation between gas emission and tilt cycles (Watson *et al.* 2000) and the suggested mechanism of crack opening related to plug movement in lava domes (Tuffen & Dingwell 2005; Tuffen *et al.* 2003). Hence, we assume that gas escape is ineffective during the stick stage and is very efficient during slip stage. Two different regions along the conduit should be considered: permeable gas flow horizontally through the magma below the plug and vertical gas flow along plug walls. There are also two permeabilities related to each region, namely the permeability of magma and the plug margin. It is unclear which mechanism might be rate limiting for overall gas escape. However, in both regions the gas flux is described by the Darcy's law:
21
where *q* is flux of the gas, η is the gas viscosity, and *k* is the effective permeability representing the whole system. If the rate of the gas escape is limited by the horizontal flow through the magma below the slipping plug, the permeability is vesicularity dependent: *k*(α)=5×10^{−12}α^{3.5} m^{2} (Costa 2006). Permeability of the plug walls strongly depends on the aspect ratio of the open cracks and could hardly be estimated. Using the conservation of gas mass, the rate of gas volume change is: *e*=−*div*(*q*). Finally, we use the Biot relation between effective pressure and fluid volume change:
22
where *M*_{b} is the Biot modulus of magma (approximately the bulk modulus of melt). We approximate the second derivative of *P*_{g} by the ratio between pressure drop and cross-sectional area, characterized by a squared length scale. The pressure drop is the difference between the gas and atmospheric pressure; the latter is negligible (relative to *P*_{g}). Introducing some representative length scale *L* (e.g. Edmonds *et al.* 2003) for the whole system leads to the approximation ∇^{2}*P*_{g}≅*P*_{g}/*L*^{2}. This length scale is related to the radius of the conduit below the plug, or to the height of the plug, or both, depending of which mechanism is rate-limiting. With these, equation 22 is simplified to:
22a
where *B* (1/s) is the pressure diffusion coefficient *B*=*M*_{b}*k*/η*L*^{2} controlling the timescale of the gas loss due to permeability. Equation 22a is the generalization of the complex problem of the gas escape from the vesiculating magma, which detailed description is out of the scope of this study. Combining equations 22a and 20 yields:
23

Before we discuss the coupled solution of equations 18 and 23, we search for the steady‐ state solution of equation 23, for which *dP*_{g}/*dt*=0:
24
where *P*_{sr} is the steady-state pressure and *C*_{R} is related to *P*_{sr} through the solubility law (equation 5). The value of the steady-state pressure, obtained from equation 24, distinguishes between two modes of the system: steady-state extrusion and cyclic motion. If depressurization is not efficient, *P*_{sr}>*P*_{ar}, then a steady motion of the plug is expected. However, under conditions leading to efficient decompression, the gas pressure meets the arrest value before the steady-state conditions are achieved, *P*_{g}=*P*_{ar}>*P*_{sr}. In this case, the depressurization stage is arrested and a new cycle of pressurization begins. The efficiency of depressurization is mostly controlled by frictional parameter (*a*−*b*) and the initial slip velocity. Figure 5 allows distinguishing between steady-state extrusion and cyclic motion for a given values of (*a*−*b*) and *v*_{0}. For reasonable range of the initial slip velocity *v*_{0}=0.1–10 m/s, corresponding to the plug acceleration 10^{−2}−10^{2} g over a few mm slip, the cyclic motion is realized for the (*a*−*b*) values of about 0.1 or higher (middle curve). In the case of effective diffusion, higher (*a*−*b*) values are required for the onset of cyclic motion (upper curve). On the other hand, efficient gas loss due to permeability enhances the depressurization and lower (*a*−*b*) values are sufficient for the onset of cyclic motion (lower curve). As discussed above, high positive values of about (*a*−*b*)=0.1 are expected. With this value, cyclic activity is expected at initial slip velocity *v*_{0}>0.15 m/s.

Under conditions of the transition between steady-state motion and cyclic activity, the duration of the depressurization cycle is extremely sensitive to the friction parameter (*a*−*b*). Figure 6 presents the depressurization during the plug motion for (*a*−*b*)=0.09–0.11. Higher (*a*−*b*) values shorten the duration of depressurization from of *c*. 4 hrs for (*a*−*b*)=0.1 to *c*. 1 hr (*a*−*b*)= 0.11. Permeable gas escape contributes to depressurization and further reduces the duration of the depressurization. With high permeability (*B*=5–7×10^{−5} s^{−1}) the relaxation time is *c*. 0.5 hr.

Figure 7 combines together the pressurization and depressurization stages into a single cycle for reasonable values of the model parameters (*a*−*b*=0.1, *v*_{0}=0.1 m/s, tilt=10^{−5} rad). Accounting for permeability provides efficient mechanism for the depressurization in addition to the plug motion. We should note here that the physical parameters (permeability, length scale and others forming B-value in equation 22a) controlling the gas lost by permeable flow along the conduit walls in the slip phase are not well constrained. Low gas flux associated with low *B*-values (low permeability and long pathways for gas escape) leads to negligible permeability-controlled pressure release. On the other hand, enhanced permeability-controlled pressure release may significantly shorten the depressurization stage, as shown above.

## Discussion and conclusions

We presented a model of cyclic volcanic activity that includes pressure accumulation and relaxation during stick and slip stages. Below, we discuss the role of the various processes controlling the rate of pressure change during these stages.

### Pressurization mechanisms

The present model takes into account two pressurization mechanisms: diffusion-controlled exsolution of water from a supersaturated melt into bubbles and magma inflow from a deeper source into the conduit. The supersaturation is the result of non-equilibrium decompression and crystallization of the magma. As long as the pressure in the conduit does not exceed the static friction, it builds up towards the equilibrium pressure which is controlled by the saturation pressure (under diffusion control, Fig. 4a) or the magmatic overpressure (under magma flow control, Fig. 4b).

If pressurization is driven only by diffusion from an initially supersaturated melt, then the pressure asymptotically approaches equilibrium pressure, *P*_{0}, depending on the number density of bubbles, the initial supersaturation and conduit deformation. The timescale of this process follows from equations 9 and 14 and is:

The duration of the diffusive controlled pressurization from *P*_{g}=*P*_{ar} to *P*_{g}=*P*_{mo} may be estimated as .

Using the reasonable values for the various parameters (Table 1), *t*_{mo} is of the order of an hour for an undeformable conduit and is about twice as long, once wall deformation is included. The minimum water content in the melt needed for the onset of the plug motion in the deformable conduit by a diffusive mechanism as presented in Figure 8. The calculations are done based on the water mass balance considerations. For the typical range of water contents estimated for the SHV magma (1.5–1.6%, Rutherford & Devine 2003; Sparks *et al.* 2000) the critical pressure for the onset of motion could be achieved by diffusive mechanism only for small values of tilt (<10^{−6} µrad). For the observed tilt values (above 5×10^{−6} µrad), the required water content is 1.75%. Higher water contents may in fact be plausible, since estimates of bulk water content based on melt inclusions and phase equilibrium are typically minimum values as they do not include any excess exsolved water in the magma chamber. Alternatively, an additional source of pressure is needed for onset of plug motion.

Crystallization of the magma is a source for extra degassing. It acts to maintain supersaturation, to increase the release time and to enhance pressure build-up until *P*_{mo} is reached. Couch *et al.* (2003*a*, *b*) reported timescales of a few hours to days for crystal growth under continuous decompression compatible with the SHV decompression rate. When pressure was kept constant, the crystals approached equilibrium within a few hours. Hence, the crystallization process ensures that enough water is released and extends the time for pressurization and surpassing *P*_{mo} to a few hours.

The other pressurization mechanism is related to the magma flow from the feeding source region into the conduit, which is blocked at its top by the plug. Magma flux is assumed to be proportional (with flow rate coefficient) to the pressure difference between the magmatic overpressure in the source and the ambient pressure. In the case of efficient diffusion-controlled pressurization, when high ambient pressure builds up in the conduit, the model accounts for backflow into the source, reducing the ambient pressure to the magmatic overpressure. The characteristic time of the flow-controlled pressurization is the ratio between the gas volume in the conduit *H*α, and the flow rate of the overpressurized magma, τ=*H*α/*K*_{mf}*P*_{g}, where *P*_{g} is a typical value of the gas pressure. For typical SHV conditions (*H*=10^{3} m, α=0.2, *K*_{mf}=10^{−8} m/Pa/s, *P*_{g}=10^{7} Pa), the pressurization timescale for a rigid conduit is *c*. 1 hr, compatible with the numerical solution in Figure 4b. This timescale is extended by an order of magnitude when the deformation of the SHV conduit is taken into account (Fig. 4d). In agreement with the derived timescale of the flow-driven pressurization, the numerical solutions show that the increase in vesicularity reduces the rate of pressurization. We note that when α approaches zero, the pressurization time reduces significantly. This may explain the very short pressurization time that Iverson *et al.* (2006) achieved in their model, which ignores vesicularity of magma.

Figure 9 presents the pressure evolution by diffusion (upper curve), magma flow (lower curve) and their coupled solution (middle curve with maximum). When the two processes are taken into account, the equilibrium pressure is controlled by the magmatic overpressure. At the initial stage, the two processes act to increase pressurization rate, with the timescale for both diffusion and flow of *c*. 1 hr (assuming tilt of only 10^{−6} μrad). After one hour, the pressure achieves its maximum and from then it decreases and asymptotically approaches the magma overpressure within a few hours. In other words, the long-term pressure is dictated by the magma overpressure and in the short term (hour) the diffusion can contribute to pressure build-up. If decompression-induced degassing is sufficient for reaching *P*_{mo}, the time needed to set the plug in motion may be approximated by the diffusion timescale discussed above. In the SHV, the observed tilt is *c*. 10^{−5} μrad, timescales should be longer than an hour and the relative contribution of decompression-induced degassing becomes less important. Pressurization is controlled by the rate of crystallization or by the flow rate from the feeding source region, which act over extended timescale.

### Depressurization mechanisms

The presented model utilizes the frictional framework assuming sliding on pre-existing surfaces. When the difference between magma pressure and ambient pressure exceeds the static friction between the plug and the host rock, the plug starts to extrude, and consequently bubbles expand and overpressure is relaxed. Figure 7 extends the pressurization stage up to the onset of plug motion with abrupt transition to the depressurization stage. This simplified approach leaves the detailed description of the micro-crack nucleation and very initial stage of slip out of the scope of the presented study. These processes are expected to smooth the transition from stick- to slip-mode, shown schematically in Fig. 7.

The detailed description of the very initial stage of slip is out of the scope of the presented model. During this stage, we assume a short-term (<1 second) inertia controlled acceleration over a distance of a few millimetres (these conditions were analysed by Iverson *et al.* 2006). The present model accounts for the decay of the accelerated motion followed by a steady-state slip during which the velocity slowly decreases. In reality, the onset of the motion may be more gradual, and involve competition between the formation of slip surfaces, their healing and the failure and generation of asperities. These processes lead to smoothing of the transition between the stick and the slip stages. For example, plug motion is associated with low-frequency earthquakes (Hydayat *et al.* 2000; Neuberg *et al.* 2006; Tuffen & Dingwell 2005; Tuffen *et al.* 2003), which are probably generated at the stage of the formation of the slip surfaces and asperities failure during the plug motion. In the case of lava domes, temperatures are high and fracture healing is fast (Hydayat *et al.* 2000; Neuberg *et al.* 2006; Tuffen & Dingwell 2005; Tuffen *et al.* 2003), so it is clear that a model of motion on smoothed slip surfaces is a simplification. While macroscopic motion is continuous, it may consist of many microscopic slip events (e.g. the repetitive drumbeats described by Iverson, 2006). Such geometrical complexities and detailed earthquake mechanism are not considered in this model. However, it does account for the dependence of frictional slip on pressure and is used here to express the pressure release during dome extrusion and plug motion in order to get a first order model approximation of the process.

As the plug extrudes, the bubbly magma expands and the development of permeability in the magma allows gas escape from the system, probably along connected shear planes along the plug walls (Tuffen & Dingwell 2005; Tuffen *et al.* 2003). Both plug motion and gas escape lead to gas pressure relaxation with a timescale that is highly sensitive to the frictional parameter (*a*−*b*) and permeability (see Fig. 6). When magma overpressure drops to the dynamic strength of the slip surfaces, the plug sticks to the conduit walls and blocks the vent.

### Different regimes

Three different dynamic regimes of volcanic eruptions are described by the model: stick solution, constant slip, and the intermediate stick-slip cyclic activity. Additional cases not discussed here are failure of the wall-rock during pressurization that may lead to either gas escape or to dyke and sill injection, and the case where *P*_{mo} exceeds the internal strength of the magma. In this latter case, the initiation of motion and gas escape along the plug-wallrock surfaces may lead to sudden decompression, a rarefaction wave, magma disruption and vulcanian explosions.

The style of eruption is controlled by the relative values of four pressures. Two of them are the mechanical limits for plug motion and arrest, *P*_{mo} and *P*_{ar}; the other two are the steady states for gas pressurization during the stick stage and for relaxation during the slip stage, *P*_{sp}, *P*_{sr}.

In the stick solution, if *P*_{sp}<*P*_{mo} then the conduit remains sealed, no magma extrudes and the dome does not grow.

For constant slip: *P*_{sp}>*P*_{mo}, *P*_{sr}>*P*_{ar} the gas pressure increases along the pressurization curve until *P*_{mo} is achieved and the plug is set in motion (Fig. 10a). Cycle starts with *P*_{g}=*P*_{ar} and high pressurization rate (point a). When *P*_{g}=*P*_{mo}, the plug is set in motion and the fast transition from the pressurization to relaxation regime is presented by the dotted line. Plug motion releases the pressure, but contributions from continued crystallization and water transfer to the bubbles, magma flow into the conduit and relaxation of the chamber walls may compensate and lead to a steady-state extrusion with *P*_{sr}>*P*_{ar}. The extrusion rate for the steady-state plug motion is expressed by equation 18 with *P*_{g}=*P*_{sr}. Periods of steady extrusion with no discernible pulsations at SHV (Sparks *et al.* 1998; Wadge *et al.* 2006) can be attributed to this flow regime.

For cyclic activity when *P*_{sr}<*P*_{ar}, stick-slip behaviour is expected (Fig. 10c, d). When *P*_{g}=*P*_{ar}, the plug arrests and a new pressurization stage begins. The dotted line at *P*_{ar} in Figure 10c represents the transition from pressure decrease along the relaxation curve to its increase along the pressurization curve. Thus, cyclic pattern is now expected (Fig. 10d) with gas pressure oscillating between *P*_{mo} and *P*_{ar}; *dP*_{g}/*dt* is never zero and the steady-state solutions are never reached, as long as *P*_{sp}>*P*_{mo}. This condition demands recharge of the conduit with fresh magma.

### Limitations of the model

The presented model is an heuristic one and considers certain processes. It splits the whole system into two cells: the first stands for the magma evolution in the conduit, while the second represents the plug and the brittle failure associated with its motion. The model ignores any kind of thermal effects or heat exchange between these virtual cells and surroundings. In this sense, it is isothermal and the only interaction with the host rock is mechanical response of edifice to the magma pressure change. The adopted relation between conduit deformation and magma pressure mimics linear elasticity of the edifice and leaves out any irreversible deformation components (e.g. Cayol 2003; Widiwijayanti *et al.* 2005) or effect of shear stress along the conduit wall (discussed by Green *et al.* 2006). The amount of deformation, or coefficient of proportionality between pressure and volumetric deformation, is not well constrained.

The sharp transition between the ductile magma and brittle plug disregards the semi-brittle behaviour observed in the laboratory experiments, where dome materials at low strain rate act as non-Newtonian fluids while acting as brittle solids at high strain rates. (e.g. Dingwell 1998; Lavalle *et al.* 2007) and discussed by Neuberg *et al.* (2006) and Tuffen & Dingwell (2005) in context of SHV seismic activity. In this context, a low rate of pressure build-up might result in aseismic flow of the magma without reaching the brittle-ductile transition. It remains to be established whether such an extrusion regime is capable of causing cycles, and the occurrence of seismic swarms during inflation cycles indicates a brittle mechanism of the kind we have modelled. The model assumes frictional plug sliding along pre-defined surfaces and uses rate- and state-dependent friction to produce a first order approximation for average rate of plug motion. It does not account for local details such as asperity failures and gradual nucleation and formation of slip surfaces and associated generation of low-frequency seismicity. After the arrest of plug motion, the model considers fast healing and full recovery of strength during the stick. This is similar to Tuffen *et al.* (2003) and reasonable under the high temperature conditions typical for SHV.

Detailed kinetics of the fracture nucleation, development and healing is expected to introduce smoothed transition from stick to slip stages. In order to simulate the details of the volcanic activity, including stochastic temporal cyclic behaviour, evolving seismicity pattern, as well as details of magma flow from the source region into the conduit, one has to consider a 2- or 3-D model incorporating the physical processes discussed in this paper.

The formation of the plug itself is beyond the scope of this model. A plausible mechanism is that gas escape from the foamy magma below the overlying plug results in densification of the magma by bubble collapse as the internal pressure in the bubbles reduces below the magma pressure itself. This collapse involves volume reduction and will act against prolonged movement of the plug. Rapid loss of gas from the melt to the underpressured bubbles in this stage could also induce further microlite formation, large increases in magma viscosity and strength. This densified and rheologically stiffened magma will then add a zone of impermeable magma to the base of the plug, which will assist in the start of a new cycle of pressure build-up. A more advanced model of the cycles would seek to include the compaction and gas pressure evolution in bubbles in the relaxation stage of the stick-slip cycles.

### Concluding remarks

The present model demonstrates that the timescale of water diffusion is too fast to account for the pressurization of the dome. Edifice deformation increases the pressurization time, but it is still of the order of an hour or two and shorter than the 6–10 hours of the Montserrat Dome. When edifice deformation is included, the amount of water released due to decompression and degassing is not sufficient to set the plug in motion. Crystallization and magma flow into and out of the conduit may lengthen the time and produce higher pressures. Gas escape during the slip stage ensures that pressure drops and the plug sticks. Together, they produce a physical model that accounts for the major processes that must take place and reproduces the right timescale for pressurization and relaxation.

## Acknowledgments

We acknowledge the reviewers H. Tuffen and A. Longo for their constructive reviews which greatly improved the manuscript. We thank O. Melnik for fruitful discussions. R. S. J. Sparks would like to acknowledge a Royal Society Wolfson Merit award and NERC research grants.

- © The Geological Society of London 2008