Very long period oscillations of Mount Erebus Volcano

1. Introduction

[2] Mount Erebus is a large stratovolcano (elevation 3794 m, volume ∼1670 km³) overlying very thin continental crust (∼20 km [Bannister et al., 2000]) in the Terror Rift at the western boundary of the West Antarctic rift system [Behrendt, 1999]. Erebus, along with its subsidiary volcanoes, forms Ross Island, a key research and infrastructure site for Antarctic science hosting McMurdo (United States) and Scott (New Zealand) bases. Above ∼3500 m, Mount Erebus has a summit plateau composed of interbedded phonolitic pyroclastic, bomb deposits and lava flows, surrounding a summit cone. This summit configuration is believed to have been assembled during the past 37,000 years (C. Harpel, personal communication, 2001). The summit cone hosts an approximately 250-m-radius, 120-m-deep, Main Crater, which, in turn hosts the near vertically walled, roughly circular, 100-m-deep Inner Crater with a radius of ∼80 m. The Inner Crater contains a convecting phonolitic lava lake containing abundant anorthoclase phenocrysts [Kyle, 1977; Dunbar et al., 1994], which over the last 10 years has varied from 5 to 15 m in radius. Except for short-lived shallow burial by eruptive debris in September through December 1984 [Keinle et al., 1985; Caldwell and Kyle, 1994; Kyle, 1986; Kaminuma, 1994], the lava lake has been a persistent feature for three decades. Most recently (through October 2003), the lava lake has resided within the northern half of the Inner Crater (Figure 1). The observation of what were probably large Erebus Strombolian eruptions by J. C. Ross in January 1841 [Kyle et al., 1982; Kyle, 1994] is consistent with a similar summit magmatic system existing at that time. Open volcanic systems [e.g., Albarède, 1985] featuring long-lived conduit lava lakes are uncommon [Tazieff, 1994a, 1994b; Rothery and Oppenheimer, 1994], with only a handful of other worldwide cases, e.g., Erta 'Ale, Ethiopia [Le Guern et al., 1979], Nyiragongo, Zaire [Tazieff, 1984], Kilauea, Hawaii [e.g., Tilling, 1987; Wolfe et al., 1987], and Masaya, Nicaragua [e.g., Stoiber et al., 1986].

[3] Building on seasonal volcanological studies in the early 1970s [e.g., Giggenbach et al., 1973], short-period seismic monitoring at Erebus has been performed since December 1980, initially as a component of the International Mount Erebus Seismic Study [Dibble et al., 1984; Kaminuma, 1994]. Since 1992, seismic recording has been a key component of the interdisciplinary Mount Erebus Volcano Observatory (MEVO [Rowe et al., 2000]). Short-period seismic and infrasonic observations have also been supplemented with video recording [Dibble et al., 1988; Dibble, 1999; R. Aster et al., New instrumentation delivers multidisciplinary real-time data from Mount Erebus, Antarctica, submitted to Eos Transactions of the American Geophysical Union, 2003, hereinafter referred to as Aster et al., submitted manuscript, 2003] and with GPS and tilt measurements [Kyle and Johns, 2000; Kyle et al., 2002; Aster et al., submitted manuscript, 2003]. Erebus eruptive activity is predominantly Strombolian, caused by gas slugs up to 5 m in radius (directly observable in video recordings during clear weather) decompressing explosively at the lava lake surface. Slugs are usually single, but a few percent of the time are noted to be multiple detonations with separation times of less than 1 s, as observed both on video and in infrasound records [Rowe et al., 2000].

[4] Infrasound observations indicate bubble overpressures at the lava lake surface of around 10⁵ Pa [Johnson, 2000], placing an upper bound on slug sequestration depth of ∼6 m, assuming negligible decompression on ascent and the Erebus conduit magmatic density estimates of Dibble [1994]. Eruptive activity during the past three decades has also included infrequent phreatic events, notably on 19 October 1993, when a new crater was explosively excavated southwest of the Inner Crater [Dibble et al., 1995]. Strombolian eruptive activity has varied through unusually elevated (September through December 1984 [Kaminuma, 1994]) and depressed periods (December 1990 to December 1995). Small ash eruptions with durations of up to several minutes, and lava flows of several tens of meters (arising from two distinct Inner Crater vents within 100 m of the lava lake) were also noted during December 2000 to January 2001. Harmonic tremor, previously rare on the mountain, has recently become common (December 2000 to January 2001; March 2002 to July 2002; January 2003 to present) with some large amplitude events visible at the Scott Base Global Seismic Network station SBA, 35 km from the summit [Ruiz et al., 2002]. This tremor does not appear on infrasound records and appears to be several kilometers deep, originate below the several kilometers summit, and to be decoupled from eruptive activity. The most recent tremor activity (frequent episodes from January through July 2003) has coincided with a period of very low eruptive activity (M. Ruiz, personal communication, 2003). Interpretation of these tremor episodes is at present ambiguous, as some (and perhaps all) of them are generated near the coast of Ross Island by large icebergs [Talandier et al., 2002]. Erebus produces only a very small number of internal volcanic long-period or volcanotectonic earthquakes [Rowe et al., 2000], although isolated swarms have been associated with a possible episodic dike injection north of the summit [Rowe, 1988].

[5] An increasing diversity of very long period (VLP; T > 5 s period pulses or tremor) seismic signals associated with active volcanoes have been discovered during the past decade, spurred by the ready availability of portable broadband seismic instrumentation for volcano research through IRIS PASSCAL and other programs [Aster et al., 2000; Chouet, 2003]. Coupled with increasingly interdisciplinary instrumentation [e.g., Ripepe et al., 2001; Aster et al., submitted manuscript, 2003], this rapidly growing data archive has produced new observations of volcanic processes. The volcanological setting for VLP signal generation is diverse, ranging from low-to-intermediate viscosity systems such as Erebus [Rowe et al., 1998, 2000], Miyakejima [Kumagai et al., 2001; Kobayashi et al., 2003], and Stromboli [Neuberg et al., 1994; Ripepe et al., 2001; Chouet, 2003; Chouet et al., 2003], to higher viscosity systems such as Popocatepetl, Mexico [Arciniega-Ceballos et al., 1999]; Merapi, Indonesia [Hidayat et al., 2000, 2002]; Aso, Japan [Kawakatsu et al., 1994; Kaneshima et al., 1996; Kawakatsu et al., 2000]; and Usu, Japan [Kawakatsu and Yamamoto, 2001]. Conversely, other active volcanoes such as Karymsky, Sangay [Johnson, 2000; Johnson et al., 2003], and Arenal [Hagerty et al., 2000] have not been observed to generate such signals. VLP events can be associated with eruptions (e.g., Erebus, Stromboli, Merapi, Popocatepetl), or with noneruptive internal processes (e.g., Aso). Some VLP observations have been associated with magma transport at depths of several kilometers, for example, at Long Valley Caldera and Kilauea [Ohminato et al., 1998; Hill and Swanson, 2001; Hill et al., 2002] and Bandai, Japan [Nishimura et al., 2001]. In this paper we summarize observations and consider source processes for the VLP signals observed at Erebus.

2. Broadband Seismic Observations

[6] Monitoring of Erebus activity prior to 1996 was accomplished using a network of short period (SP; f > 1 Hz) seismometers augmented by a near-summit infrasonic sensor. During November 1996 to January 1997, three IRIS PASSCAL 30-s corner seismometers were deployed in a pilot experiment on the summit plateau at sites NKB, LVA, and HUT, located 0.7, 1.4, and 1.9 km from the lava lake, respectively [Rowe et al., 1998] (Figure 1). The 1996–1997 deployment showed that VLP signals with modal periods as grave as 20 s are unfailingly associated with Strombolian eruptions from the lava lake [Rowe et al., 1998, 2000]. In December 1999 to January 2000, MEVO and PASSCAL personnel temporarily deployed broadband (30-s Guralp CMG-3ESP and 120-s CMG-3T) and infrasound [Johnson, 2000, 2003; Johnson et al., 2003] stations at eight sites on the summit plateau across an azimuth range of ∼180°. In 2000–2001, a Guralp 40-T (30 s) sensor was permanently installed at site E1S in conjunction with a rim-sited, telemetered, Inner Crater surveillance camera. Most recently, in 2002–2003, five stations featuring seismometers, infrasound, tilt, GPS, weather, and other sensors were installed on the volcano (Aster et al., submitted manuscript, 2003).

[7] Broadband displacement seismograms from Erebus Strombolian explosions are dominated by VLP signals with highly similar waveforms across the network (Figures 2 and 3). VLP signals are associated solely with Strombolian events, and have never been observed except in association with eruptions from the lava lake. Ash eruptions from other vents within the Inner Crater produce no signal at frequencies below the SP band. The repeatability of VLP observations indicates a nondestructive or self-reconstructing source mechanism stimulated by the eruptive disruption of the lava lake from least late 1996 through at least early 2003.

[8] A characteristic Strombolian eruption sequence observed on the video camera is shown in Figure 4 and is interpreted in Figure 5. The corresponding seismogram from station E1S (∼700 m from the lava lake; Figure 1) is shown in Figure 6. Because VLP signals are similar from event to event (except for amplitude, which ranges over approximately a factor of 25 in displacement), and second-order variability in the preeruptive onset (see below), seismogram stacking can be applied to enhance signal to noise. Such stacking is especially valuable for reducing the strong Ross Island oceanic microseism near 7 s.

[9] Seismograms show that the VLP onset typically leads the SP signal by several seconds and shows significant variability from event to event. On some high signal-to-noise events (Figure 6) the VLP onset can be observed 5 s or more before the explosion. Individual event seismograms show that the onset varies in polarity, amplitude, and duration (Figure 7). The onset shown in stacks is thus not necessarily representative. Particle motions are observed to steepen systematically as the signal evolves (Figure 6) [Rowe et al., 1998]. Figure 8 shows VLP stacked displacement seismogram particle motions, where each stack has been properly normalized to preserve mean relative VLP amplitudes. Particle motions project radially inward toward rough convergence ∼220 m west-northwest of the lava lake (Figure 8a). Projected particle motion dips intersect a vertical axis underneath the lava lake at inconsistent depths ranging from 0.6 km (NKB) to 3.3 km (E1S) (Figure 8b). Figure 9 shows corresponding stacked VLP seismograms at each site for the particle motions of Figure 8.

[10] Figure 10 shows vertical and radial component displacement power spectral densities estimated from the stacked seismograms of Figure 9. VLP spectra exhibit distinct spectral peaks at about 20.7, 11.3, and 7.8 s, which we denote T₀, T₁, and T₂, respectively. The VLP displacement power band at f ≤ 0.2 Hz is spectrally distinct from, and of much larger amplitude than, the SP band at f ≥ 1 Hz, and the two frequency ranges are separated by a band of very low spectral amplitude. Power spectra featuring such “dead bands” are unusual in transient seismic events (appearing occasionally, for example, in tsunamigenic earthquakes [e.g., Schindele et al., 1995]). Recordings made in 1996–1997, compared to 1999–2000, suggest that the amplitude (but not period) ratio of the VLP spectral peaks is a variable feature of the system, with a higher proportion of 20.7 s energy apparent in spectra recorded in 1996 relative to later years (Figure 10). Figure 11 shows a vertical component spectrogram of stacked VLP signals recorded at E1S. Mode amplitudes show no systematic dependence on event size (across an observed amplitude range of about a factor of 25), or during individual events codas. This is surprising behavior, given that the VLP source is associated with the extremely nonlinear surface explosion and evisceration of the lava lake. Signal partitioning between the SP and VLP bands is approximately conserved from event to event (Figure 12). Figure 11 also shows very low-Q shorter-period VLP modes not visible in individual events.

3. VLP Source Considerations

3.1. Moment

[11] Assuming a shallow source, the shortest period principal VLP mode (T₂ = 7.8 s) observed at the furthest broadband station (∼2 km) is approximately 1/12 of a wavelength from the eruption source, and its character, apart from polarization, is similar across all stations (Figure 2). The VLP signal is thus an approximately quasi-static, source time function dominated, and near-field deformation. An order-of-magnitude estimate for the seismic moment necessary to produce observed VLP displacements can be obtained by considering a point source embedded in a half-space. Mogi [1958] gives a useful analytic formula for the static displacement field from a step time function point source isotropic moment tensor in a half-space [e.g., Johnson, 1974]. Expressed in terms of an equivalent pressure change Δp in a spherical cavity of radius a,

where (u_r(r), u_z(r)) is the (radial, vertical) displacement field at the surface (z = 0) as a function of source-receiver radial distance r, z₀ is the source depth, and a Poisson solid is assumed. The equivalent isotropic seismic moment of a Mogi source is [e.g., Aki and Richards, 2002]

[12] From equation (1), the term Δpa³ required for a spherical cavity to produce the observed near-field displacements associated with the large explosions (∼20 μm of vertical displacement at a range of 700 m) can be estimated under the quasi-static assumption and assuming no large effects from topography or heterogeneity. For a nominal shallow source depth z₀ = 100 m, and using a Poisson solid model for the Erebus summit region (μ = 3.9 × 10⁹ Pa [Dibble, 1994]), we obtain Δpa³ ≈ 5.3 × 10¹⁰ N m, corresponding (equation (3)) to a time-integrated moment rate of M₀ ≈ 5.0 × 10¹¹ N m, comparable to the scalar moments observed at Merapi volcano [Hidayat et al., 2002]. As the summit region elastic properties of Dibble [1999] were determined using seismic first arrival times, this moment estimate is probably an upper bound, as the bulk elastic compliance of the summit region, including its complex stratigraphy and magmatic system, is presumably somewhat less. The analogous microearthquake moment magnitude [e.g., Abercrombie, 1996] is

Erebus Strombolian explosions, on the other hand, typically have SP seismic energies corresponding to local magnitudes of order 1 or less, with the very largest eruptions of the last three decades (1984) having magnitudes near 2.5 [Dibble et al., 1988]. The VLP moment is thus anomalously large (by more than an order of magnitude) relative to the SP moment, as can also be noted in the displacement power spectral densities (Figures 10 and 12). In contrast, simple explosion spectra produce flat or decreasing displacement spectral amplitudes at low frequencies [e.g., Brune, 1970; Boatwright, 1980; Walter et al., 1988]. The spectral isolation of the short-period and VLP bands, and the SP/VLP moment discrepancy, suggests that the eruption process excites two distinct source mechanisms.

3.2. Effect of Tilt on Seismograms

[13] Long-period observations on horizontal seismograms can be susceptible to surface tilt [Wielandt and Forbriger, 1999; Wilson et al., 2002; Hidayat et al., 2002]. This arises because a seismometer in a gravity field cannot distinguish between tilt and horizontal acceleration. The Mogi source tilt in radians is

and the apparent tilt-induced horizontal displacement signal at a period of T is the twice-integrated apparent horizontal acceleration

where g is gravitational acceleration. The magnitude of equation (6) becomes quadratically larger with increasing period and produces outward tilt and outward deformation. An outward tilt will reduce apparent radial amplitudes and correspondingly increase particle motion dip. For a nominal shallow source depth of z₀ = 100 m, the effect for a Mogi source in a half-space is negligible at periods ≤ T₀ = 20.7 s (Figure 10) at distances greater than several times z₀, using the Erebus summit region elastic parameters of Dibble [1994].

[14] A Mogi source is only useful here as a first-order model to assess a best fit isotropic source size. Aside from the dip discrepancy described above, the Mogi source also fails to adequately model observed VLP displacements in other ways. First, the observed particle motions become increasingly elliptical with distance and have progressively greater inclinations with advancing time (Figure 8a) than the source-directed inclinations predicted by the Mogi model. Second, displacement amplitude fall-off with distance predicted by a shallow Mogi source is more rapid and regular than that observed. Finally, although observed VLP particle motions are roughly radial with respect to the lava lake, clear amplitude departures exist (e.g., the amplitude at station NOE; Figure 9). These discrepancies suggest that the source deviates significantly from isotropy, principally in that single force and/or force couple source components significantly exceed horizontal ones, particularly later in the VLP signal (for example, as in an oscillating finite subhorizontal crack [Fialko et al., 2001]), and that the source may furthermore be significantly spatially extended. Finally, significant discrepancies between Mogi model predictions and observations may also be expected because the elastic half-space is likely to be a poor model for characterizing the topography and internal structure of a complex volcano. For this reason the conclusion that tilt is unimportant is provisional.

3.3. Source Q

[15] A characteristic feature of the Erebus VLP source is its oscillatory nature (e.g., Figure 11). The decay of the modal components can be quantified using a decaying linear oscillator model. For a damped linear oscillator with a time domain amplitude

where Q is the source quality factor, the spectral amplitude is

Fitting normalized Erebus VLP spectral peaks to equation (8) using the −10 dB spectral widths (approximately 0.014, 0.015, and 0.068 Hz, for the spectral peaks near 0.048, 0.085, and 0.12 Hz, respectively) give Q values of approximately 11, 18, and 4.

4. VLP Source Models

[16] We next investigate repeatable source processes for generating Erebus VLP signals consistent with observed amplitudes, duration, period, Q, and particle motion characteristics. In the summit region of an active volcano, such source processes may generally include elastic/acoustic, gravity, and/or flow/inertial processes within an elastic summit edifice.

4.1. Trapped Body Waves

[17] A starting point for considering possible VLP source processes is trapped acoustic or seismic body waves bounded by a strong impedance contrast. Dibble [1994] investigated elastic and density properties of the highly vesiculated shallow magma at Erebus and concluded that the physical properties of the upper few hundred meters are primarily affected by the growth of H₂O-dominated foam [Crespo, 1969; Murase and McBirney, 1973; Miksis and Ting, 1987]. Dibble [1994] and others [e.g., Garces, 2000] have noted that P wave velocities in a lava foam can be as low as ∼20 m/s, increasing rapidly and approaching nonvesiculated magma velocities at depths of a few hundred meters (Figure 13).

4.1.1. Longitudinal Resonance

[18] Travel time contributions for vertically traveling waves in the highly vesiculated magma at z < 300 m from slowness integration of the models of Figure 13 range from approximately 0.68 to 1.64 s. Velocities increase rapidly below this highly vesiculated region. A vertical conduit resonating in a fundamental longitudinal acoustic conduit mode (featuring an approximate displacement node at depth and an approximate displacement antinode at the surface) would have to have a length comparable to a substantial fraction of the thickness of the crust (∼11.5 km) to generate a T₀ = 20.7 s fundamental longitudinal mode as a λ/4 wave. A longitudinal resonator with open-closed (antinode-node) boundary conditions would display period ratios for its three lowest harmonics of approximately 1:0.33:0.20, contrary to the closer-spaced observed values T₀:T₁:T₂ ≈ 1:0.56:0.39. On the other hand, an extensive conduit resonator having open-open or closed-closed boundary conditions would have T₀:T₁:T₂ ≈ 1:0.50:0.33, which is much closer to the observed modal ratios. However, pure acoustic excitation of VLP modes from the lava lake explosion seems unviable given the VLP/SP moment discrepancy discussed in 3, and the observation that the VLP signal is precursory to the excitation of the seismoacoustic SP (explosion) signal. Thus, while vertically traveling conduit P wave modes are likely important contributors to the SP spectrum [e.g., Garces, 2000], explosively excited longitudinal acoustic modes can be dismissed as a cause of the VLP signal. For similar reasons, we also dismiss vertically propagating wave types such as tube waves [e.g., Marzetta and Schoenberg, 1986].

4.1.2. Radial Resonance

[19] An alternative VLP source possibility is radial acoustic or seismic standing waves with a quasi-cylindrical or other conduit-associated geometry. Two elastic possibilities for such modes with potentially suitably slow phase velocities are acoustic modes in the low-velocity magma foam, or slow crack waves trapped in a large fluid-filled crack system [Chouet, 1986; Ferrazzini and Aki, 1987]. The eigenmodes of a cylindrical resonator satisfying Laplace's equation for wave propagation at phase velocity c(ω) are [e.g., Blackstock, 2000]

where J_m, (m ≥ 0) is the mth-order Bessel function of the first kind, a is the resonator radius, z_n,m^{(1, 2)} (n ≥ 0) is the value of the nth zero (1) or zero of the derivative (2) of J_m (depending on whether the perimeter displacement boundary condition is nodal or antinodal), and the eigenfrequencies are

[20] Because of the observed approximate radial symmetry of the displacement field (most notably its lack of phase reversal with azimuth; Figures 2, 8a, 8b, and 9), we restrict consideration to the purely radial (m = 0) modes. A prediction of the radial acoustic mode model for a cylindrical resonator with boundaries defined by an abrupt large increase in seismic impedance is that modal period ratios should nearly conform to those of a nodal boundary condition. Predicted and observed modal period ratios for nondispersive waves (n ≤ 8) are shown in Table 1.

Table 1. Normalized Periods of Linear Radial Modes in a Cylindrical Conduit, m = 0

n	zn,0(1)	Tn, 0/T0 (Nodal Boundary)	zn,0(2)	Tn,0/T1 (Antinodal Boundary)
0	2.4048	1	0	∞
1	5.5201	0.4357	3.8317	1
2	8.6537	0.2779	7.0156	0.5462
3	11.7915	0.2039	10.1735	0.3776
4	14.9309	0.1611	13.3237	0.2876
5	18.0711	0.1331	16.4706	0.2326
6	21.2116	0.1134	19.6159	0.1953
7	24.3525	0.0988	22.7601	0.1684

[21] Observed period ratios for the three prominent spectral peaks of the VLP oscillator are approximately 1:0.56:0.39. The antinodal boundary condition thus produces a much better fit to the period data than does the nodal boundary condition appropriate for acoustic modes. Additionally, for acoustic radial modes in the shallow summit magma chamber at Erebus, using the minimum near-surface acoustic velocity (20.1 m/s) of Dibble [1994], and assuming that the gravest observed VLP resonance, T₀ = 20.7 s, represents u_0,0⁽¹⁾ (9, Table 1), the required resonator radius is

(Figure 14), which substantially exceeds the lava lake (5–15 m) and Inner Crater radii (∼80 m; Figure 1). As with the longitudinal acoustic case, radial acoustic mode models fail to account for the VLP/SP moment discrepancy discussed in 3. Additionally, it is difficult to see how a stable foam acoustic velocity could be sustained in the face of lava lake-eviscerating explosions. Finally, the anticipation of the SP signal by the VLP signal and the modal period ratios thus suggest that a purely acoustic explanation for the VLP signal is inadequate and that SP and VLP excitation mechanisms and energy sources must be distinct, although clearly genetically related.

[22] The acoustic resonators considered above utilize the low acoustic velocities of compressional waves in foamy magma. However, when contacts between solid and fluid media are considered, other wave types with slower phase velocities may exist. An intriguing possibility is a very low velocity crack or crack system wave [Chouet, 1986; Hidayat et al., 2002; Ferrazzini and Aki, 1987]. In its simplest geometry, a fluid crack between two half-spaces, Ferrazzini and Aki [1987] showed that such waves can become excited at long wavelengths with the phase velocity decreasing as the inverse square root of wavelength. Using such a model, the 7-s-period VLP signals of Aso Volcano can be modeled by a magma-filled crack with a dimension of ∼0.5 km [Kubotera, 1974]. A slow crack wave resonator model is alluring in that it offers the possibility of a conduit-enveloping resonator that could be excited by forces generated during eruption and recovery while being insensitive to the nonlinearities of magma and gas transport.

4.2. Gravity Waves

[23] In Table 1 the first three antinodal boundary condition modal periods fit the observed principal VLP resonances to a few percent (1:0.56:0.39 observed versus 1:0.55:0.38 theoretical). As the antinodal boundary condition is appropriate to the reflection condition for gravity waves, this leads to the consideration of possible magma gravity waves. Gravity-linked phenomena, being driven by gravitational energy release associated with mass transport prior, during, and after eruptions, have the appeal of offering a possible potential mechanism for VLP excitation that is genetically linked to the eruption, while being distinct from the elastic energy released in the explosive decompression of the SP source.

4.2.1. Surface Gravity Waves

[24] The theory of surface seiches in lakes or harbors has long been of interest in ocean and limnological systems [e.g., Jeffreys, 1925]. Nondispersive surface gravity waves in a half-space of depth h [Lighthill, 1978] propagate at phase velocities

for wavelengths λ ≫ h. A radial surface standing gravity wave with a fundamental (n = 1) modal period of 20 s in an a = 5 m radius lava lake requires a phase velocity of (equation (10); Table 1)

which requires an extremely shallow lava lake for shallow gravity waves (h on the order of centimeters using equation (12)).

[25] For λ ≪ h, dispersive surface gravity waves in a fluid half-space propagate at phase velocities

For a 20 s period wave, equation (14) gives c_d ≈ 32 m/s, which requires (10, Table 1) a lava lake of radius

which greatly exceeds the Erebus lava lake (5–15 m) and Inner Crater radii (∼80 m; Figure 1). We can also rule out surface gravity waves as candidates for the VLP resonance based on the consistent destruction of the lava lake surface in eruptions, the changing character of the lava lake surface during refill, and the stability of the observed modal periods across many field seasons despite significant changes in the dimensions of the lava lake.

4.2.2. Internal Gravity Waves

[26] A more general gravity wave phenomenon to be considered is internal gravity (buoyancy) waves. Such waves are capable of having very low phase velocities. If nondispersive gravity waves of suitable phase velocity and amplitude were to occur in the near Erebus summit magmatic system, the modal period of the three largest VLP signals would be well fit by the successive modes of a cylindrical gravity wave resonator (Table 1). The general behavior of small-displacement internal gravity waves in fluids for a one-dimensionally stratified density column is controlled by the density-versus-depth function (pycnocline), ρ(z), and by the free surface [Phillips, 1969; Lighthill, 1978]. Internal wave behavior in a 1-D pycnocline is characterized by the buoyancy (Väisälä-Brunt) frequency. If the bulk modulus is large enough so that the acoustic velocity

gravity and acoustic waves decouple, in which case, the buoyancy frequency is

where ρ′(z) = dρ/dz. N(z) is real for the common case of a stratified fluid with ρ′(z) decreasing and ρ(z) increasing with depth, the situation encountered in an increasingly vesiculated magma approaching the surface.

[27] The Sturm-Liouville equation describing the vertical fluid velocity field,

for linear horizontally propagating internal waves in a low-viscosity system bounded by displacement nodes at top and bottom at z = 0 and z = −D (as in a magma chamber with subhorizontal upper and lower geometry), is

where W = 0 at z = −D and z = 0, k is the wave number, and c is the phase velocity.

[28] There are no traveling wave solutions, W(z), to equation (19) for ω > N_max. Solutions at lower frequencies have motions that are approximately confined in a waveguide defined by a depth range where ω < N(z). For a given ω, eigensolutions W_n(z) exist for an infinite number of eigenwave numbers, k_n [Phillips, 1969], where each eigensolution has n + 1 extrema and n zero crossings. The lowest-order (n = 0) internal mode is purely sinuous, in the sense that the vertical velocity field throughout a vertical profile of the internal waveguide (ω < N(z)) is in phase. The n = 1 mode is purely varicose, in the sense that vertical velocities will be out of phase near the waveguide boundary, and feature an internal node. Higher modes incorporate a combination of sinuous and varicose features.

[29] For a specified N(z) and c, successively lower eigenvalues, k², of equation (19) correspond to longer wavelengths and lower frequencies associated with successively higher modes. The phase velocity dispersion, c(ω), for solutions to equation (19) coupled with equation (10) allows for the prediction of eigenfrequencies for horizontally traveling trapped internal gravity wave modes in a cylindrical chamber of radius a for an arbitrary buoyancy frequency depth function, N(z).

[30] To investigate internal gravity mode behavior expected for a stratified summit magma chamber, we spline interpolated the four pycnocline models of Dibble [1994] (Figure 15) and, to facilitate the calculation of smooth derivatives, fit each of these models to a scaled analytic formula

where ρ_min and ρ_max are the density extrema for each model and

where ξ_min and ξ_max are the extremal values of

Corresponding buoyancy frequencies (equation (17)) are shown in Figure 15b.

[31] Solutions for the eigenvalue-eigenfunction equation (19) for N(z) corresponding to each model were investigated for a bounded (0 ≥ z ≥ −D) fluid layer, such as could exist beneath the floor of the Inner Crater, and for a range of phase velocities, c, using a Fourier Grid Hamiltonian method [Marston and Balint-Kurti, 1989]. Figure 16 shows example eigenfunctions, and Figure 17 shows dispersion curves for internal gravity waves in each pycnocline model.

[32] Phase velocities for gravity waves in a pycnocline with a characteristic length scale are dispersive and also possess a shortest allowable period that makes it impossible to obtain a good fit to the observed VLP periods using a simple model of radial resonance for the lowest-order (sinuous) mode. Figure 18 shows radial mode period solutions for a sub-Inner Crater reservoir of radius a = 50 m, where successive values of a/z_n,0⁽²⁾ (Table 1) are plotted as horizontal lines, and the ω-normalized phase velocities, c(T)/ω = Tc (T)/2π, are overlain as a function of period, T. Intersections between the ω-normalized phase velocities and the a/z_n,0⁽²⁾ give period solutions satisfying equation (10) and can be compared to observed VLP periods. The fit between the observed T₀, T₁, and T₂ and the periods predicted by simple internal gravity wave radial modes in a quasi-cylindrical reservoir is seen to be poor for the Dibble [1994] model-derived pycnoclines (Figure 15). Other solutions for this system show that this is true for higher modes as well. Furthermore, reducing dispersion by invoking a shallow reservoir (even a very shallow one of less than 10 m thickness) is also unsuccessful in producing a suitably nondispersive phase velocity curve (e.g., Figure 14).

[33] An additional complication for gravity waves is the damping effect of the kinematic viscosity, ν = η/ρ (where η is the viscosity). For Erebus near-surface phonolite, a reasonable viscosity range of 10³ to 10⁴ Pa s [Dibble et al., 1984] gives a minimum kinematic viscosity of order 0.5–5 Pa s m³/kg, assuming a nominal density of 2 × 10³ kg/m³ for the magma. Viscous energy attenuation for a linear internal wave occurs at a rate [e.g., Lighthill, 1978] of k²v per unit time, where k is the wave number. A wavelength of λ₀ = 19.7 m (Figure 16), with a corresponding wave number of 0.32 m⁻¹, predicts energy dissipation per second of greater than 5%. After one 10 s period as a radial standing wave in a cylindrical chamber, the energy in such a mode is dissipated by an amount greater than of order 0.95¹⁰ ≈ 0.60, or is equivalently decreased in amplitude by ∼0.77. This lower attenuation value is in the right range for the observed modal decay (the observed Q values of ∼10 corresponds to an energy decrease of ∼0.62 per period).

[34] The above linearized decay calculation, however, is sensitive to the value of ν, and assumes perfect conduit wall reflection and Newtonian viscosity. If the viscosity of Erebus phonolite is at the upper end of the estimated range (∼10⁴ Pa s), viscous forces will be strong relative to inertial forces, and highly oscillatory gravity wave modes would be prohibited. Unambiguous and long-standing visual observations of lava lake convection, with velocities of a few m/s, further suggest that convection may significantly disrupt the pycnocline and internal waves phenomena. A final difficulty with the internal gravity wave model is the amplitudes necessary for producing the observed displacement amplitudes. The simple Mogi model considered earlier requires Δpa³ ≈ 4 × 10¹⁰ N m. For a spherical source comparable in radius to the Inner Crater (a = 80 m), the corresponding excess pressure is Δp = 0.8 × 10⁵ Pa, or ∼0.8 times atmospheric pressure. This order of pressure variation would require substantial gravity wave amplitudes, even for the minimal amplitude case of a surface wave. The necessary surface wave height required to produce this level of pressure variation would be on the order of 3 m for a density contrast of ρ = 2600 kg/m³. The reduced density contrast for internal waves would require proportionately larger amplitudes to exert similar forces. Because of these difficulties, a pure gravity wave process appears untenable as the principal cause of the VLP signal.

4.3. Oscillatory Recharge

[35] Another candidate source for the VLP signal is oscillatory gravity-driven inertial effects in the lava lake and upper conduit system following its explosive evisceration, including conduit oscillations and oscillatory recharge. Two possibilities are oscillatory conduit wall forces during recharge through a constriction and inertially induced recharge oscillations from a summit magma reservoir due to eruptive disruption of the system equipotential.

4.3.1. Constricted Flow

[36] Julian [1994, 2000] investigated nonlinear vibrations induced by constricted flow in a fluid transport system. In the Julian model, oscillations ensue from feedback between pressure and elastic constriction variations in response to the Bernoulli effect. Such oscillations can mimic linear resonances, but show additional nonlinear effects such as period doubling. The Julian constriction oscillator, however, still produces harmonically related modal periods, which are not observed in the Erebus VLP signal. A feature of the Julian model which adds to its attractiveness for explaining features observed in volcanic tremor is that small changes in source constriction conditions can drive significant changes (such as spectral gliding [e.g., Hagerty et al., 2000]) in spectral line frequencies, as well as transitions between harmonic and chaotic behavior. However, this sensitivity argues against it being important here because long-term stability of the Erebus VLP signal indicates very stable oscillator conditions.

4.3.2. Oscillatory Conduit-Reservoir Flow

[37] Oscillatory recharge of the lava lake system might follow eruptive evisceration from an underlying, and perhaps surrounding, summit magma reservoir. In this case the removal of the upper few tens of meters of the lava lake in eruptions initiates what is essentially a slug test of the magmatic conduit, and the VLP signal arises from pressure variations and fluid flow reaction forces applied to the summit region during the subsequent (∼5 min) refill of the system. The response of such a system is of interest in aquifer characterization and has been modeled in this context by numerous authors [e.g., Rehbinder and Apazidis, 1994; McElwee and Zenner, 1998; Zlotnik and McGuire, 1998]. Highly oscillatory well level responses result if the dimensionless damping coefficient, C_D, (which is inversely proportional to the radial hydraulic conductivity) is sufficiently low [Van der Kamp, 1976]. For example, the oscillatory linearized level response, w_d, of a fluid column in a transmissive medium to a slug test, assuming radial Darcian flow, consists of a decaying harmonic oscillator response in which the oscillating fluid column and aquifer constitute the inertial and spring components of the system, respectively [Zlotnik and McGuire, 1998]

characterized by a single oscillation frequency

where

and L_e is the effective length of the fluid column.

[38] As the linearized equation (23) produces only a single spectral component, it does not provide an explanation for the observed relative periods and amplitudes of the individual VLP modal components (e.g., Figure 11). One possibility is multiple recharge paths producing a superposition of responses, although this is an ad hoc explanation for what would be more satisfactorily and simply explicable as a multimodal phenomenon. A perhaps better model is thus the contribution of nonlinear effects to the recharge response [e.g., McElwee and Zenner, 1998] such as non-Darcian flow, nonlinear frictional losses, or other nonlinearities to the system. The observed linearity of the VLP response with respect to event amplitude indicates that any such nonlinearity would have to be weak to keep the VLP modal frequencies constant throughout each event and to keep the VLP amplitude proportional to eruption size, as indicated by the SP energy. An inertial and/or viscous recharge delay would conveniently explain the time between the eruptive evisceration of the lava lake and the maximum VLP amplitude, which is typically around 5 s later. A further appealing aspect of the recharge model is that, although it requires the summit magmatic system to be suitably well connected for rapid recharge, it does not require the presence of an excitable long-lived resonator with reflecting boundary conditions necessary to trap strong acoustic or gravity wave modes. Observational support for underdamped, posteruptive recharge at Erebus comes from the video records of Dibble et al. [1988] during 1987–1988 and from small bubble bursts observed in 2001, in which the lava lake surface was observed to show surging behavior with a period of 8.8 ± 1.6 s, close to T₂, the most energetic of the VLP modal periods (Figure 11). As noted above, such grave periods are inconsistent with surface gravity wave modes, indicating that these oscillations in the lava lake surface could be forced by surging recharge. A schematic showing the lava lake eruption and inferred refill processes is shown in Figure 19.

5. Discussion

[39] VLP signals at Erebus are a repeating response of the volcano that presumably require stable conditions to be maintained in the near-summit conduit system and perhaps adjoining structures, throughout hundreds of Strombolian explosions. The extent of the Erebus lava lake has varied seasonally (between approximately 5 and 15 m radius) during several decades of observation, and it has apparently maintained its open conduit Strombolian activity, while exhibiting little inflationary/deflationary activity associated with magma storage or depletion at shallow depths, even during exceptional Strombolian eruptive activity in 1984 [Otway et al., 1994]. This is consistent with the lava lake being the exposed top of a more extensive summit magma reservoir that provides sufficient thermal mass and stable conditions to maintain the lava lake, sustain the convective circulation and crystallization of anorthoclase phenocrysts [Dunbar et al., 1994], and sustain characteristic Strombolian activity on at least a decadal timescale.

[40] Variations in the lava lake radius appear to represent a change in the extent to which the underlying magma system is capped by the Inner Crater floor. The size of the Inner Crater may be somewhat coupled with the extent of the near-summit magma system, and has also been observed to increase in size from mass wasting events over the past several decades (P. Kyle, personal communication, 2001). A working model of the summit region (Figure 20) incorporates a summit magma reservoir with a radius comparable to that of the Inner Crater (∼80 m). Variability in timing and polarity of preeruptive VLP signals associated with gas slug rise (Figure 7) suggests a well-connected underlying plumbing system featuring multiple sites for slug coalescence and differing corresponding buoyant ascent paths, risetimes, and perhaps depths. Such variability is indicated by the modeling of general Strombolian systems [Vergniolle and Brandeis, 1996; Vergniolle et al., 1996], where even subtle conduit ruggedness is shown to be capable of sequestering gas bubbles at multiple sites [e.g., Jaupart and Vergniolle, 1988]. Alternatively, differences in the slug ascent phase between events could conceivably also be affected by variable magmatic viscoelasticity.

[41] Eruptive differences between the various Erebus active Inner Crater vents and pronounced variability in the frequency of Strombolian eruptions further suggests that the propensity of the system to sequester eruptive size slugs is sensitive to the shallow conduit geometry. Although the lava lake delivers intact slugs with explosion times measured in tenths of seconds, the ash vent conduit is, in contrast, suitably constricted to prolong eruptions to several minutes. Because of this gradual energy release and mass transport, ash eruptions do not produce the rapid disturbance to the system apparently necessary to excite the VLP source. The lava flow vent (Werner's Fumerole) is a smaller lava lake, but has never been observed to host Strombolian eruptions, presumably because its conduit does not intersect a gas slug sequestration zone.

[42] After surveying a host of fundamental models in this paper, we feel that source details of VLP oscillations at Erebus still remain somewhat enigmatic. Nevertheless, we can point to some potentially key clues. First, the VLP excitation unambiguously precedes the surface explosion, and this initial VLP signal is much more variable than the extended coda. VLP oscillations are thus initially excited by the passage of the ascending gas slug and not by the surface explosion. This variable initial VLP signature likely reflects variable geometric and ascent characteristics of the gas slugs. This suggests conduit and/or gas transport complexity beyond that observed in laboratory Strombolian systems, where a consistent downward single force caused by magma flow around a rising slug dominates [Ripepe et al., 2001; Uhira and Takeo, 1994].

[43] The highly repeatable extended VLP coda begins during the rapid removal of upper conduit material and application of the eruptive reaction force, and persists through the (several minute) refill of the lava lake. This suggests that the VLP source is a relatively narrowband oscillator that is somewhat insensitive to the details of excitation during this period. The distinct spectral bands and amplitudes for the VLP and SP signals, together with a prolonged excitation period that matches the refill period of the lava lake system, suggests that the VLP source derives its energy primarily from gravitational and inertial forces. This gravitational energy is of two forms, with initial excitation being attributable to the buoyant ascent of the gas slug (and the complementary downward flow of surrounding conduit material), and the VLP coda being a response to the rapid removal of eruptive products, subsequent gravitational disequilibrium, and reestablishment of equilibrium via shallow magmatic transport leading to the reconstitution of the lava lake.

[44] Because of the eruptive destruction of the upper conduit system and stable VLP signals observed during posteruptive recovery, we suggest that the VLP signal is not primarily due to trapped modes or processes in the dynamically changing shallow magmatic system, but is rather the nearly linear response of an enveloping and more stable structure, such as the summit magma reservoir and/or a surrounding slow crack wave resonator to oscillatory excitation forces within the conduit. The two most viable processes for the extended VLP coda in this regard are resurgent recharge of the vesiculated magma, and a corresponding excitation of slow crack modes in a fractured near-summit zone, with the minutes-long coda sustained during the recharge of the lava lake. A slow crack wave resonator in a fractured hydrothermal or other fluid-saturated region, excited by gas or multiphase transport, has previously been invoked to explain 7-s-period oscillatory signals at Aso [Kawakatsu et al., 1994; Kaneshima et al., 1996; Kawakatsu et al., 2000] and for 20-s-period pulse-like signals at Merapi [Hidayat et al., 2002]. Hidayat et al. at Merapi and Chouet et al. [2003] at Stromboli further suggest that such an underlying crack/conduit system could play the dual role of facilitating both gas slug accumulation and VLP response, and is conceivable that such a geometry may be essential to generating a VLP signal. Comparisons between eruptive and noneruptive VLP signals and between VLP and non-VLP generating eruptive systems at the increasing number of relatively well instrumented volcanoes worldwide thus merits further investigation.

[45] Moment rate tensor modeling of VLP systems is an essential tactic in pursuing further insight into source properties and model testing and is a component of ongoing Erebus investigations. Despite improvements in data quality, source modeling of VLP systems remains a difficult problem which requires more detailed consideration of forces generated by multiphase fluid flow during and after eruptions, accurate near-field modeling of Green's functions in complex media, including realistic estimates of tilt, and further study of possibilities for the excitation and propagation of slow seismic phases in crack systems and in more generally fractured media. Nonradial aspects of VLP particle motions (Figure 8a) and steep summit topography suggest that accurate moment tensor estimation at Erebus will require proper incorporation of topography in Green's function calculations. Progressive steepening of VLP particle motion dip (Figure 6) suggests that the force system evolves temporally from the ascending slug phase through the eruption to the extended coda phase. Much of this steepening is attributable to the rapid decay of the more shallowly dipping, initially large amplitude, T₀ (7.8 s) mode and the relative persistence of the steeper-dipping, graver modes [Rowe et al., 1998] (Figure 11). From this we believe that moment tensor inversion methodologies allowing for a temporally varying equivalent force system will be necessary to adequately model the time evolution of the VLP signal.

6. Conclusions

[46] We present detailed observations spanning several field seasons and consider first-order source issues for Erebus VLP near-field deformations generated during Strombolian explosions and lava lake recovery. The Erebus VLP signal is characterized by linearity with respect to signal amplitude, repeatability, highly oscillatory character and near-radial particle motions. In this regard it appears similar to other eruptive (e.g., Stromboli) and noneruptive (e.g., Merapi, Aso) observations of VLP volcano signals. The facility with which complementary observations, such as video, can be made of the exposed conduit at close (several hundreds of meters) range makes Erebus a favorable locale for continued observation and investigation into these phenomena.

[47] A proposed model for Strombolian eruptions and repeatable VLP signal generation at Erebus consists of the following key elements: (1) a conduit system with appropriate geometry to sequester gas slugs; (2) a near-summit magmatic system of sufficient extent and thermal mass to maintain physical openness and stability and to rapidly recharge and reset the system; and (3), possibly, a surrounding narrowband and temporally stable VLP generating structure capable of being excited by gas and magma transport in the conduit system. A corresponding eruption/VLP scenario for Erebus consists of five phases:

[48] 1. Eruptions are preceded by a preparation period characterized by growth and sequestration of gas slugs (Figure 19a) with radii up to ∼5 m (volumes of ∼500 m³). If this phase is accompanied by a significant volumetric stress changes, the associated time constant is longer than is observable with 120 s broadband seismometers (at least several minutes).

[49] 2. The earliest indication of an impending eruption on the broadband seismograms is the preeruptive VLP signal generated by the buoyant ascent of substantially intact gas slugs (Figure 19b). Until a critical buoyancy threshold is reached, gas slugs are sequestered at shallow (between tens of meters and ∼6 m) depth. Seismic polarity and timing characteristics of the ascent phase show significant variability, both from event to event and from season to season, suggesting multiple gas slug source regions and associated variability in ascent characteristics and source forces. Source components during the buoyant ascent phase involve at least three processes acting in the upper conduit system walls: (1) expansion of the bubble with decompression, (2) leading expansive displacement, including precursory bulging of the lava lake, and (3) trailing void implosion.

[50] 3. The eruption (Figure 19c) is caused by the emergence of the overpressurized gas bubble at the surface of the lava lake and subsequent explosive evisceration of the upper conduit system (Figure 19d). The eruption also generates a reaction force that may further contribute to the excitation of the VLP signal.

[51] 4. The eruption is followed by an extended (several minutes) oscillatory VLP coda that persists as the lava lake reforms (Figure 19e). It is suggested that the VLP coda arises from of oscillatory magma recharge forces interacting with viscous gravity waves and/or driving a stable slow wave near-summit crack system resonator. The repeatability of the VLP coda indicates that the VLP source mechanism is temporally stable on at least a multi-year time scale (Figures 19c and 19d).

[52] 5. The cycle completes with the recovery of the lava lake and cessation of the VLP signal (Figure 19f). High signal-to-noise stacked VLP seismograms and video observations indicate that the recovery times for the lava lake and the duration of the VLP signal are similar (∼5 min), suggesting an association between the VLP source and the gravitational disequilibrium in the shallow conduit system.