A little more than 760 ka ago, a supervolcano on the eastern edge of California (United States) underwent one of North America’s largest Quaternary explosive eruptions. Over this ∼6-day-long eruption, pyroclastic flows blanketed the surrounding ∼50 km with more than 1400 km3 of the now-iconic Bishop Tuff, with ashfall reaching as far east as Nebraska. Collapse of the volcano’s magma reservoir created the restless Long Valley Caldera. Although no rhyolitic eruptions have occurred in 100 k.y., beginning in 1978, ongoing uplift suggests new magma may have intruded into the reservoir. Alternatively, the reservoir could be approaching final crystallization, with present-day uplift related to the expulsion of fluid from the last vestiges of melt. Despite 40 years of diverse investigations, the presence of large volumes of melt in Long Valley’s magma reservoir remain unresolved. Here we show, through full waveform seismic tomography, a mid-crustal zone of low shear-wave velocity. We estimate the reservoir contains considerable quantities of melt, >1000 km3, at melt fractions as high as ∼27%. While supervolcanoes like Long Valley are rare, understanding the volume and concentration of melt in their magma reservoirs is critical for determining their potential hazard.

Volcanoes capable of explosive caldera-forming supereruptions are exceedingly rare, yet are arguably the most globally catastrophic natural process on the planet. In a single eruption, these volcanoes can erupt >1000× the volume that erupted from Mount St. Helens (Washington State, USA) in 1980 (1 km3) (Crosweller et al. 2012). Of the 13 Quaternary-active supervolcanoes in the world, three are in the continental United States: Long Valley (California), Valles (New Mexico), and Yellowstone (Wyoming) (Crosweller et al. 2012).

At approximately 767 ka, Long Valley volcano in eastern California (Fig. 1) erupted >1400 km3 of rhyolitic ash and pyroclastics in an ∼6-day-long Plinian eruption (tephra equivalent; Hildreth and Wilson, 2007). With the exception of the 639 ka eruption of Yellowstone, this is North America’s most recent supereruption (Crosweller et al. 2012). Unlike Yellowstone, Long Valley has no deep-mantle hotspot source, and regional mafic volcanism is relatively young (<4 Ma) (Bailey et al., 1976). Tectonically, the volcano sits in a right-stepping offset in the dextral Walker Lane shear zone (Fig. 1). The region is bounded by the Sierra Nevada to the west and the Basin and Range to the east (Fig. 1), with volcanism focusing at Long Valley at ca. 2.5 Ma at a left-stepping offset in the Sierra range-front fault system (Bailey et al., 1976).

Long Valley’s caldera-forming supereruption was followed by 120 k.y. of subPlinian activity and 400 m of resurgent dome uplift (Hildreth, 2017). While minor eruptions continued to 100 ka, there have been no eruptions on this resurgent dome in 500 k.y. (Hildreth, 2017). This quiescence, along with a perceived migration of volcanism westward toward compositionally distinct systems (Hildreth, 2017), a post–ca. 300 ka decrease in high-temperature hydrothermal activity, and the current absence of magmatic gases (CO2, 3He/4He), have led to the perspective that the caldera-forming magma reservoir is now near wholly crystallized (Hildreth, 2017). This assessment holds crucial implications for the interpretation of ∼30 yr of ongoing uplift centered on the post-caldera resurgent dome, and the long-term hazard of the volcano (Hill and Montgomery-Brown, 2015).

Is uplift at Long Valley driven by the exsolution of magmatic fluids from the last few vestiges of melt from the caldera-forming reservoir; i.e., “second boiling” (Hildreth, 2017)? Or could the uplift be related to the intrusion of new magma (Battaglia et al., 1999)? While North America’s other two Quaternary supervolcanoes, Valles and Yellowstone, likely contain significant quantities of melt (Huang et al., 2015; Steck et al., 1998), Long Valley remains enigmatic.

Geophysical Overview

The search for a modern magma reservoir beneath Long Valley led to more than 20 geophysical studies over the past 40 yr, most commonly based on local-earthquake tomography (LET) of P-wave and/or S-wave traveltime (VP, VS) or S-wave attenuation. While these studies agree there is likely no melt within the upper ∼6 km (Romero et al., 1993; Foulger et al., 2003; Seccia et al., 2011; Lin, 2015) their ability to resolve a deeper magma reservoir (>8 km) is limited by the focal depths and raypaths of local earthquakes. Beneath the caldera, LET studies can, at most, resolve to depths of ∼6–10 km. Despite this limitation, these studies have found evidence for a low VP and VS or high S-wave attenuating zone near the base of their models (Romero et al., 1993; Foulger et al., 2003; Seccia et al., 2011; Lin, 2015). However, given the proximity to the relatively low-temperature base (100 °C) of the Long Valley Exploration Well (Sorey et al., 2000), low VP/VS, and low resistivity, many suggest these anomalies reflect magmatically derived fluids and/or hydrothermal alteration (Romero et al., 1993; Lin, 2015; Peacock et al., 2016).

The best evidence to date for a deep magma reservoir comes from limited lower-resolution teleseismic studies, which show an 8–30% reduction in VP between 7 km and 20 km depth (Dawson et al., 1990; Steck and Prothero, 1994; Weiland et al., 1995). The most recent of these studies imaged two isolated systems: a mid-to-upper crustal reservoir at ∼10–18 km depth, and a second, deeper, low-velocity zone, centered at ∼25 km depth (Weiland et al., 1995). This type of multilevel storage is similarly observed at the Valles (Steck et al., 1998) and Yellowstone calderas (Huang et al., 2015), and is supported by the latest perspectives on the transcrustal development of large silicic systems (Cashman et al., 2017).

To address the longstanding uncertainties of Long Valley’s magma reservoir, we solved for the caldera’s crustal shear-wave velocity structure using three-dimensional (3-D) full-waveform tomography. We invert for traveltime differences between source and forward-modeled seismograms using a multi-step iterative process that includes 3-D finite-difference wave propagation simulations and the calculation of 3-D finite-frequency sensitivity kernels (Flinders and Shen, 2017). Our method allows us to account for the complex 3-D spatial sensitivity of wave propagation, scattering of short-period waves by topography, and P/S-wave velocity cross-dependence, phenomena typically ignored in traditional tomography but potentially significant at volcanic settings at short periods (2–30 s).

Source seismograms are Rayleigh-wave, ambient noise cross correlations, derived from all seismic stations within 150 km of the caldera over the past 26 yr (Fig. 1; Figs. DR1 and DR2 in the GSA Data Repository1). These data are independent of earthquake locations and provide sensitivity to seismic structure deeper than previous LET studies at spatial resolutions higher than teleseismic body-wave studies. In later iterations, we supplement these data with Rayleigh-wave records from 11 large regional earthquakes (Fig. 1; Tables DR1 and DR2).

We image a spheroidal low shear-wave velocity (VS) zone (30 km in diameter) that underlies the entire caldera (Fig. 2). Near-surface (<3 km depth) low velocity is indicative of volcanoclastic caldera infill, with the primary reservoir low-VS zone extending from ∼5 km depth, near the base of down-dropped Paleozoic metasedimentary roof (Bailey et al., 1976), to more than 20 km (Fig. 2D). On average, this zone is ∼20% slower (2780 m/s) than expected for a mid-crustal granite (3460 m/s) near its solidus temperature (350 MPa, ∼680 °C; Evans et al., 2016; Ji et al., 2002). The magnitude of the reduced velocity and the vertical extents agree with what has been observed teleseismically (Dawson et al., 1990; Steck and Prothero, 1994; Weiland et al., 1995) and with the depths of previously imaged P-wave reflectors from active-source seismic refraction experiments (Hill et al., 1985). The zone also correlates with an ∼10% slow VP zone imaged using regional-scale earthquake tomography (Thurber et al., 2009). Our low-VS zone encompasses a previously imaged two-station receiver function reflector at 7–11 km depth (Fig. 2E) attributed to 30–60% melt (Seccia et al., 2011). While hydrothermal alteration contributes to reduced velocity at shallow depths (<7 km) (Peacock et al., 2016), it does not likely have significant effects at greater depths (Barnes, 1997). Similarly, the velocity is too low to be explained by heat remaining from a crystallized magma reservoir, and implies the presence of residual melt.

Melt Estimates

To constrain the melt fraction in the low-VS zone, we correct for variations in VS in granite from temperature and pressure (Fig. DR10), assuming the entire zone is held at lithostatic pressure and at the mean Bishop Tuff eruptive temperature (∼750 °C; Evans et al., 2016; Ji et al., 2002). While this temperature correction will overestimate the average reservoir temperature, it will provide a conservative/minimum estimate of the magnitude of the VS perturbation attributed to melt fraction. We limit these estimates to regions where the minimum contrast between the reservoir velocity and regional Walker Lane crust (WLC) are >5% (Fig. 2E). Using two independent partial derivatives of VS with respect to melt fraction, one derived experimentally (Caricchi et al., 2008) (δVS/δM% = −26 m/s) and one we derived specific to an average mineral composition of the Bishop Tuff (δVS/δM% = −23 m/s) at 750 °C and 350 MPa (Fig. DR11), we estimate an average melt content of 23% ± 4%.

Studies on the generation of ignimbrite eruptions and silicic batholiths provide constraints on the plutonic/volcanic ratios necessary to generate mid-to-upper crustal rhyolitic magmas (Lipman and Bachmann, 2015). By increasing the minimum VS contrast between the reservoir and WLC, we can calculate a range of reservoir and melt volumes across these geologically reasonable plutonic/volcanic ratios. Ratios are estimated using our calculated reservoir volume compared to the 767 ka Bishop Tuff dense-rock-equivalent volume (650 km3) (Hildreth and Wilson, 2007). For a minimum velocity contrast of 12%, the reservoir volume is 3300 km3, equivalent to an ∼5:1 plutonic/volcanic ratio, and contains ∼900 km3 of rhyolitic melt (27%) (Figs. 2E and 3). Using a smaller velocity contrast (8.5%), equivalent to an ∼10:1 ratio, the reservoir volume increases to 6600 km3 and contains ∼1400 km3 of melt (22%) (Figs. 2E and 3). For comparison, based on P-wave non-full waveform methods, Yellowstone’s upper-crustal reservoir has been estimated to contain a similar 900–1400 km3 of rhyolitic melt (Huang et al., 2015; Chu et al., 2010). While the volumes of melt at these two supervolcanoes appear equivalent, Yellowstone’s total reservoir volume is considerably larger (4300–10000 km3) (Huang et al., 2015; Chu et al., 2010). However, at Long Valley, melt is likely significantly more concentrated.

Irrespective of the thermal history of the reservoir over the past 500 k.y., our estimate of ∼23% melt suggests that a substantial volume of fluid has exsolved from the reservoir. Crystallization of a 4 wt% H2O pure melt (Hildreth and Wilson, 2007) to a 0.77 crystal fraction is approximately equivalent to exsolution of half of the original water content (Botcharnikov et al., 2005). Although hydrothermal pathways are smaller than our tomographic resolution, following the waning of the paleo-hydrothermal system (300 ka; Hildreth, 2017), fluid could be accumulating near the roof of the modern reservoir (Fig. 4). However, sensitivity testing indicates that it is not possible to reproduce our imaged low-velocity zone from vertical smearing of a shallow low-velocity or fluid-rich zone (Figs. DR8 and DR9).

Thermal Implications

Our melt estimates imply significant intrusions occurred following the last resurgent-dome eruptions at 500 ka. Without continued heat input, the caldera-forming reservoir, at an assumed initial 70% melt fraction, would have crystallized to <5% melt in 225 ± 100 k.y. through conductive cooling alone (Fig. DR12A). To account for ∼23% melt solely from the remains of crystallization would require the reservoir to be at a 70% melt fraction as recently as 160 ka ± 50 ka. This seems unlikely given the significant decrease in high-temperature hydrothermal activity in the central caldera after ca. 300 ka (Hildreth, 2017), and implies even more-recent intrusion. Assuming intrusions reactivated the reservoir, the thermal front could take a minimum of ∼90 k.y. to emplace an appreciable 25 °C perturbation above the background geotherm at the depth of the Long Valley Exploration Well (Fig. DR12B). Any reactivation more recent than this would likely not yet be observable.

Recently, the inability of geophysical methods to consistently identify melt-rich magma reservoirs (>10%) in the upper crust has been used to argue that large volumes of melt are ephemeral (Cashman et al., 2017). Although controversial, lithium diffusion experiments argue that these reservoirs only approach a melt fraction required for eruptability (>35%, i.e., above rheological lock-up) (Cashman et al., 2017) for, at most, centuries after rapid heating from new intrusions (i.e., “cold storage”; Rubin et al., 2017). So, while the reactivation necessary to produce a contemporary 23% melt fraction at Long Valley likely occurred no later than 90 ka, these “cold storage” hypotheses could suggest significantly more-recent changes to the reservoir.

Although we cannot discriminate between magmatic intrusion and mobilization of exsolved fluids as the driver of recent uplift at Long Valley, we can conclude the mid-crustal reservoir is still melt-rich. We estimate the reservoir currently contains enough melt to support another supereruption comparable in size to the caldera-forming eruption at 767 ka. However, this volume and a relatively high melt fraction in no way ensures that the magma is eruptible. Equally important is how that melt is distributed within the reservoir, a characteristic that remains beyond our tomographic resolution. As tomography provides average solutions smoothed over spatial scales larger than individual melt-filled dikes, sills, and fissures, it can overestimate volumes and underestimate melt fractions. Melt at Long Valley could be concentrated in smaller zones, at melt fractions above the rheological lock-up window (35%). Future research focusing on mid-crustal seismic anisotropy, dense-station receiver functions, expanded magnetotelluric studies, or continued scientific drilling would help address the question of melt distribution, and is crucial to progressing our understanding of this volcano.

We thank those who contributed data to the Incorporated Research Institutions for Seismology Data Management Center (U.S. National Science Foundation) and Northern California Earthquake Data Center. W. Hildreth and J. Fierstein introduced Flinders to the geology of the Long Valley area. C.F. Boomer provided critical last-minute support. Resources supporting this work were provided by the NASA High-End Computing Program through the NASA Advanced Supercomputing Division at Ames Research Center, California. This work was funded by the U.S. Geological Survey Mendenhall Research Fellowship Program. Any use of trade, firm, or product names is for descriptive purposes only and does not imply endorsement by the U.S. Government.

1GSA Data Repository item 2018290, methods and data, Figures DR1–DR12, and Tables DR1–DR3, is available online at http://www.geosociety.org/datarepository/2018/ or on request from [email protected].