Deep magma mobilization years before the 2021 CE Fagradalsfjall eruption, Iceland

A major challenge for volcanology is the recognition of when a volcanic system is evolving toward eruption. The deeper parts of magmatic systems play a key role in supply of magma to the surface and, therefore, transitions from unrest to eruption (Sparks and Cashman, 2017). Geophysical techniques do not have sufficient spatial resolution to image deep crustal melt distributions on a scale relevant to physical models of melt migration (Maclennan, 2019). Uncovering the behavior of the deeper part of the magmatic system can be achieved through temporal records of magmatic processes derived from crystals carried by the magma, which then can be interpreted in combination with geophysical data (Kahl et al., 2011, 2013; Saunders et al., 2012; Pankhurst et al., 2018; Rasmussen et al., 2018). Linking magmatic processes to geophysical monitoring records requires robust and internally consistent time-scale estimates to be determined from erupted products. Modeling the diffusive relaxation of compositional zoning in mineral phases at magmatic temperatures can provide such estimates. This method, diffusion chronometry, is a key tool for recovering time scales of magma residence, mixing, and mobilization (Costa et al., 2003; Cooper and Kent, 2014; Kahl et al., 2015, 2017; Hartley et al., 2016; Lynn et al., 2017; Mutch et al., 2019a, 2019b; Caracciolo et al., 2021). Insights from such studies have supported the concept that volcanic systems are dominated by magmatic mushes: solid frameworks of crystals with interconnected melt. Mush disaggregation plays a key role in assembly of eruptible magmas that exist in a high-meltfraction state (Maclennan, 2019). Such disruption and upheaval provide the shifts in crystal composition that initiate diffusion chronometers within carried crystals. Few studies have considered diffusion in multiple mineral phases in basaltic rocks to test the robustness of mineral geospeedometers (Mutch et al., 2021). We report the results of studying olivine (Fe-Mg) and plagioclase (Mg) diffusion chronometers in parallel, thus examining deep magma mobilization and mush disaggregation prior to the 2021 CE Fagradalsfjall eruption in southwestern Iceland. The well-documented geophysical precursors of the eruption present an ideal opportunity for comparing geophysical and petrological data sets and exploring the links between magma mobilization and eruption priming.

Over the past 4000 yr, eruptive activity along the Reykjanes Peninsula has been episodic, with 200–500-yr-long eruption periods separated by 600–1200 yr of quiescence (Sӕmundsson et al., 2020). During each eruption period, volcanic activity has taken place within all of the volcanic systems of the Reykjanes Peninsula but at different times (Fig. S1 in the Supplemental Material¹). After a repose period of 781 yr, new eruptive activity commenced on 19 March 2021, with a basaltic eruption at Geldingadalir in the Fagradalsfjall hyaloclastite complex (Fig. 1). Geothermobarometric, melt inclusion, volcanic gas, and geodetic data provide strong evidence that the Fagradalsfjall magma was sourced from a compositionally diverse set of primary mantle melts that resided and mixed close to the Moho (depth of ~15–20 km; Weir et al., 2001; Halldórsson et al., 2022; Sigmundsson et al., 2022). The latest Icelandic eruption known to have been directly sourced from a near-Moho storage zone was the Borgarhraun eruption in the Þeistareykir volcanic system, which occurred 8000–10,000 yr ago (Maclennan et al., 2003). The Fagradalsfjall eruption is thus not only the first known on the Reykjanes Peninsula since the 13^th century but also the first deep-sourced (15–20 km; Halldórsson et al., 2022) eruption on a spreading-ridge system to be monitored with modern instrumentation. Crystal records preserve a rare history of near-Moho magmatic processes for comparison with contemporaneous geophysical observations.

Fresh olivine and plagioclase macrocrysts were hand picked from crushed lava and tephra samples erupted on 20–23 March, 2 April, 28 April, and 6 May 2021, from Geldingadalir (Fig. 1B). See the Supplemental Material for details of the analytical methods and for supplemental data (Tables S3 and S4). The macrocryst cargo shows textural and chemical characteristics expected for storage in magmatic mush, including nodules (Fig. S2). Olivine crystals are mostly euhedral to subhedral, with crystal lengths between 0.3 and 1.9 mm, displaying normal, reverse, and complex zonation in forsterite (Fo) content (Fig. S3), with extensive core plateaus of variable composition (Fo_85–89; Fig. S4) and rims (Fo_85–87; Fig. S4) in contact with the groundmass glass (Figs. S5–S41). Plagioclase crystals occur both as macrocrysts and groundmass phase in lava and tephra samples. Crystals have an apparent aspect ratio of 1.2–2 in two-dimensional sections with long axes >4 mm. Overgrowth rims (anorthite [An] content of An78–80; Fig. S4) are present in all samples and preferentially occur on faces perpendicular to (010) as expected during rapid crystallization (Holness et al., 2007). Rims on (010) faces can be very narrow (<5 μm) or absent, meaning that crystal cores were exposed to the carrier liquid with little or no overgrowth. Plagioclase macrocryst cores are typically anorthitic (An₈₅-An₉₀) and uniform in composition. However, half of the crystals measured show internal variations of ~2 An units in zoneŝ100 μm wide. Most crystals show prominent gradients in Mg as far as 100 μm from the crystal edges (Figs. S42–S54).

We obtained a total of 231 time estimates from modeling the diffusive relaxation of Fe-Mg in olivine (n = 180) and Mg zoning in plagioclase (n = 51) using Autodiff (Couperthwaite et al., 2020), a one-dimensional finite-difference diffusion algorithm, and DFENS (Mutch et al., 2021), a diffusion chronometry method that combines finite-element models with Bayesian inversion (see the Supplemental Material and Tables S5 and S6 therein for model details). Examples of best-fit models are shown in Figure 2 (and in Figs. S5–S54). Olivine diffusion time scales range between 2 and 1800 days (~5 yr) with an average uncertainty on each time scale of ±0.28 log units (1σ). The median time scale of the olivine population is 115 days with 57% of the time scales being shorter than 6 mo (Fig. S55).

Plagioclase diffusion time scales range from 1 to 871 days (Fig. S55), with an average uncertainty on each time scale of ±0.29 log units (1σ). The median time scale of this plagioclase population is 91 days with 76% of time scales being shorter than 6 mo, which is slightly shorter than that observed in the corresponding olivine population but in accord within uncertainty.

Measured volcano-tectonic unrest along the Reykjanes Peninsula commenced in mid-December 2019 with an intense, week-long earthquake swarm at the southern margin of Fagradalsfjall (Sigmundsson et al., 2022; i in Fig. 3). Surface inflation at a rate of 3–4 mm/day associated with intense earthquake activity centered on the Svartsengi system (I in Fig. 3; locations in Fig. 1A) started in January 2020. Inflation continued until 6 February, marking the end of this first geodetic cycle (Cubuk-Sabuncu et al., 2021; Geirsson et al., 2021). Two more inflation episodes were detected at Svartsengi on 6 March–17 April and 15 May–22 July 2020 (II and III in Fig. 3). Further inflation episodes were observed at the Reykjanes (February 2020) and Krýsuvík (mid-July to November 2020; IV in Fig. 3) volcanic systems, culminating in a magnitude 5.6 earthquake on 20 October (Cubuk-Sabuncu et al., 2021; Geirsson et al., 2021). These inflation events have been attributed to intrusion of a sill at 3.2–4.3 km depth (Cubuk-Sabuncu et al., 2021; Geirsson et al., 2021) or the injection of magma-derived gas into the shallow crust (Flóvenz et al., 2022).

While the Svartsengi and Krýsuvík systems displayed uplift in 2020, the Fagradalsfjall system, located between Krýsuvík and Svartsengi (Fig. 1A), showed no detectable inflation or deflation (Geirsson et al., 2021). Unrest in the form of seismicity was, however, high at Fagradalsfjall during 2020, and seismic moment release increased notably after 19–20 July (ii in Fig. 3) with the rapid migration of seismic swarms to ~10 km northeast of the town of Grindavík and the occurrence of two magnitude 4.9 earthquakes (Cubuk-Sabuncu et al., 2021). These shallow earthquakes were accompanied by deep seismicity (10–12 km depth) at Fagradalsfjall (Greenfield et al., 2022). A new cycle of unrest initiated on 24 February 2021 with a magnitude 5.64 earthquake (iii in Fig. 3), followed by intense surface deformation related to dike formation and magma accumulation under the area between Keilir and Fagradalsfjall (Flóvenz et al., 2022; Sigmundsson et al., 2022) (Fig. 1A).

Combining complementary crystal diffusion age records with real-time geophysical observations (Fig. 3) and geothermobarometric calculations (see the Supplemental Material) enables a more complete overview of eruption run-up processes.

The longest diffusion time scales associated with deep magmatic unrest beneath Fagradalsfjall preceded geophysical eruption precursors on the Reykjanes Peninsula by at least 1 yr: a priming phase of deep magma accumulation and mush disaggregation occurred without detection through monitoring of shallow seismicity or geodetic inflation. In the year before eruption, however, the diffusion records show increased frequency of disequilibrium onset (increased gradient in Figure 3, in early 2020) synchronous with the onset of shallow geophysical unrest at neighboring volcanic systems along the Reykjanes Peninsula. After the extended priming phase, the deep near-Moho (15–20 km; Halldórsson et al., 2022) magmatic system may have supplied magma or magma-derived CO₂-rich fluids to the shallow system (4–7 km, depth of inflation sources and earthquakes; Flóvenz et al., 2022). This supply resulted in the widespread onset of geophysical unrest consisting of recurrent cycles of inflation and deflation and associated seismicity within the upper, brittle crust. Further localization and escalation of geophysical and petrological signals in the weeks and days prior to eruption mark the onset of the acute run-up phase (defined by the occurrence of intense seismicity and high rates of deformation beginning 24 February; iii in Fig. 3). Combined with the thermobarometric record presented by Halldórsson et al. (2022), this could be taken to represent full magma mobilization and transport through the crust.

During initial dike propagation, the magma may undergo a punctuated and protracted ascent as it tries to fracture the overlying crust or exploit preexisting weakness (Lister and Kerr, 1991). High rates of surface deformation and seismicity detected between 24 February and mid-March 2021 confirm the emplacement of a vertical, 9-km-long, two-segmented dike between the surface and 8 km depth in the area between Keilir and Fagradalsfjall (Fig. 1A). Three days before eruption, the seismicity was largely concentrated in two clusters, at depths of 1–2.5 km and ~0.5 km under the eruption sites, suggesting that magma had already started to accumulate near the surface (Sigmundsson et al., 2022). Petrological observations from early-erupted products retain a record of magma stalling or heat loss associated with ascent. In the 2021 magma, we observe relatively cool (1185–1210 °C) and shallow (0.5–2.5 kbar) rim zone formation in olivine from early-erupted samples (20 and 23 March), consistent (within uncertainty) with the observed depth of dike emplacement (0–8 km) and a consequent loss of heat during either an ascent hiatus or a decreased rate of magma influx into the dike (Halldórsson et al., 2022; Sigmundsson et al., 2022).

Petrological observations collected from samples erupted between 28 April and 6 May also reveal transitions in ascent dynamics throughout the eruption itself. These samples show evidence of early disruption at depth in their zoned olivine and plagioclase cargoes but retain higher eruption temperatures (1220–1230 °C) and more primitive glasses (MgO 8.4–8.6 wt%) and have little record of shallow (4.4–5.2 kbar) crystallization (Halldórsson et al., 2022). This indicates that they were flushed directly from the deeper (15–20 km) system and ascended unimpeded to the surface through a fully established conduit.

Our data emphasize how combining comprehensive petrological and geophysical records and a near-real-time sampling approach can provide powerful insights into the pre-eruptive state of magmatic systems in the critical time before eruption run-up begins. In particular, the high temporal resolution of our approach has allowed us to identify different pre-eruptive states in the underlying magma reservoir and pinpoint the transition to full-scale magma mobilization associated with dike propagation and later eruption. Our findings provide new insights into the dynamics of magmatic systems and may offer a template for interpreting pre-eruptive signals at basaltic volcanoes worldwide.

Kahl acknowledges funding by Deutsche Forschungsgemeinschaft (DFG) grant KA 3532/2-1. Mutch acknowledges funding from a Lamont-Doherty Earth Observatory (New York, USA) postdoctoral fellowship. Maclennan acknowledges funding from UK Research and Innovation (UKRI) grant NE/W004690/1. Bali acknowledges support from Icelandic Research Fund grant 228933-051. We thank A. Schmitt and G. Lloyd for their comments on an early draft of this manuscript, and I. Fin, O. Wienand, J. Wyn Williams, and N. Groschopf for sample preparation, polishing, and assistance with EMP analysis. We are grateful to J.S. Pálmadóttir, B.V. Óskarsson, and J.M.C. Belart for assistance with Figure 1. We gratefully acknowledge constructive comments provided by Ó. Flóvenz and C. Petrone.