Phasing and climate forcing potential of the Millennium Eruption of Mt. Baekdu

The study addresses whether the Millennium Eruption (ME) of Mt. Baekdu (late 946 CE) comprised a continuous single event or two temporally separated phases (rhyolite then trachyte), and evaluates the eruption’s climate forcing potential. Prior proximal/medial tephra investigations and the presence of an intervening lahar, grain-size analyses, and historical chronicles suggest a potential hiatus of months between phases, but this remained unresolved. The authors leverage high-resolution Greenland ice-core records to reconstruct the timing and phasing of cryptotephra and sulfate deposition, infer aerosol formation altitudes from sulfur isotopes, and assess the residence time and seasonal timing of sulfate aerosols to understand why the ME left little evidence of Northern Hemisphere summer cooling despite its magnitude.

- Magnitude and timing: The ME (VEI 6–7) occurred in late 946 CE, with widespread tephra across East Asia and signals in Greenland ice cores. Historical Japanese chronicles report ashfall on 3 Nov 946 and thunderous sounds on 7 Feb 947. - Two-phase eruption hypothesis: Stratigraphy indicates an early rhyolite (Phase 1) followed by trachyte (Phase 2), with lahar deposits and distal B-Tm grain-size analyses hinting at a hiatus of months. - Sulfur emissions: Petrologic estimates range widely (2–11 Tg S syn-eruptive; up to ~45 Tg including hypothesized pre-eruptive S), though more recent work suggests total 2–7 Tg with substantial pre-eruptive sulfur loss. - Plume and injection heights: Proximal deposit inversions suggest 25–35 km plume heights; distal constraints indicated 15–20 km; recent modeling from thicknesses suggests 30–40 km. However, sulfur injection heights can be much lower than isopleth-derived plume tops. - Climate impacts: Extra-tropical NH winter eruptions can cool climate, yet tree-ring reconstructions and historical records show little to no NH summer cooling after the ME. Ice-core sulfur suggests modest stratospheric injection (~1.7 ± 0.6 Tg). These discrepancies motivate a detailed phasing and aerosol-altitude analysis using high-resolution ice-core chemistry and isotopes.

- Study material: A section of the NGRIP1 Greenland ice core spanning the ME deposition was analyzed using both discrete and continuous methods. Two longitudinal samples (~10 cm² cross-section) were taken; discrete subsampling at 2.5–5 cm across the volcanic interval and 10–15 cm for background. - Continuous-flow aerosol analyses (DRI): Employed a continuous flow melter with three streams (inner 10% for elemental chemistry including sulfur, Na, Ca by two Element2 ICP-MS; middle 20% for high-resolution liquid conductivity and size-resolved insoluble particle counts via Abakus; outer 70% discarded). Effective depth resolutions: ~1.8 mm (conductivity), ~20 mm (elemental chemistry), ~3 mm (insoluble particles). Sea-salt corrections computed ssNa and ssS; non-sea-salt sulfur (nssS) derived accordingly. - Cryptotephra extraction and geochemistry: Discrete melted samples were centrifuged; residues concentrated, dried on Kapton tape, mounted in epoxy, lightly polished, carbon-coated. Shards (4–38 µm) identified and analyzed by SEM-EDS and EPMA (JEOL iSP100; 15 kV; beam currents 5 nA/1 nA; beam sizes 5 µm/3 µm depending on shard size; counting times: 30 s for Si, Fe, Mg, Ca, K, P, Mn, Ti; 20 s Al; 10 s Na at 5 µm setup, 8 s at 3 µm). International glass standards (ML3B, StHs6/80-G, Lipari obsidian) used; accuracy generally within 2σ of references. - Annual net accumulation model: Based on digitized 1991–1995 north-central Greenland monthly net accumulation (four sites; Shuman et al. 2001). Generated 1000 realizations assuming uniform distributions within observed monthly ranges (January range inferred from neighboring winter months) to obtain monthly fractions summing to an annual mean of 191 ± 22 mm w.e. Fractions applied to convert depth separations into elapsed days using NGRIP average annual ice-layer thickness for 940–949 CE (0.183 ± 0.023 m for timing; 0.197 m ice yr⁻¹ for flux calculations). - Time-interval estimation between phases: Identified double spikes in large (4.5–9.5 µm) particle counts and liquid conductivity. Depth separations (e.g., ~11–18 mm depending on proxy) converted to time using the accumulation model and an assumed eruption start of 1 Nov 946. A forward-air-trajectory-based transport delay of 14 ± 7 days from source to Greenland was applied. Monte Carlo (n=1000) used to propagate uncertainties and derive probability distributions for spike timings and inter-spike interval. - Forward air trajectory analysis: NOAA HYSPLIT (READY) 240-h forward trajectories using GDAS 1° data (years 2005–2023) initialized at Mt. Baekdu from 8 and 14 km AMSL indicated likely Greenland arrival within ~10 days for 14 km starts and potentially longer/less certain for 8 km. Adopted a 14 ± 7 day deposition onset window (2σ). - Sulfur isotope analysis: Sulfate isolated via automated anion exchange (prepFAST MC; AG1X8 resin) and measured for δ34S, δ33S, δ32S on a Neptune Plus MC-ICP-MS. Results reported vs. V-CDT; Δ33S computed from measured ratios. Background correction applied when volcanic fraction fvol > 0.60 using background sulfate concentration and isotope ratios to isolate volcanic δ34S and Δ33S; uncertainties assessed via Monte Carlo. - Sulfate deposition estimation: Converted high-resolution nss-sulfate concentrations to flux using NGRIP accumulation (0.197 m ice yr⁻¹) and integrated over the event interval. Due to overlap with a shoulder peak, summed the first half of the main sulfur peak and doubled it, subtracting background flux similarly scaled. Noted that this simplification underrepresents the true exponential tail of aerosol decay. - Chronology and depth references: Microparticle onset at 218.5225 m defines initial deposition depth; depths for spikes and sulfur end (218.431 m) tied to accumulation fractions to estimate elapsed days from 1 Nov 946. - Quality control: Consistency between discrete and continuous sulfur concentrations verified; cross-checked seasonal markers (Na, Ca) and ECM H+ for event seasonality.

- Two distinct depositional spikes: Large (4.5–9.5 µm) insoluble particle concentrations and high-resolution liquid conductivity show two well-separated spikes with depth separations of ~11–18 mm (exceeding instrumental resolution), mirrored in ECM H+ (20–30 mm separation). Tephra shard geochemistry indicates Spike 1 (sample T4) is dominantly rhyolitic and Spike 2 (T4–T5) is dominantly trachytic, matching proximal stratigraphy. - Hiatus duration: Using an annual net accumulation model and transport delay (14 ± 7 days after 1 Nov 946), Spike 1 occurred 46 ± 21 days after 1 Nov, Spike 2 at 77 ± 23 days, yielding an inter-spike interval of 31 ± 25 days (2σ). Probability the interval is within 1 week: 2.5%; within 4 weeks: 42.7%; within 8 weeks: 97.9%. - Seasonal context: Microparticle peak aligns with seasonal Na and Ca maxima, consistent with late 946 CE timing. Sulfate deposition persisted for 232 ± 42 days after 1 Nov (primarily winter–spring), implying little sulfate remained by summer 947 CE. - Sulfur isotopes: Δ33S values show negligible S-MIF with no positive-to-negative Δ33S evolution, indicating sulfate formation largely below the stratospheric ozone layer (UT/LS origin), contrasting with known stratospheric events (e.g., Tambora, Okmok). - Near-synchronous deposition: Insoluble particle and sulfur fallout in NGRIP1 are near-synchronous, consistent with UT/LS residence times (weeks–months) rather than mid–high stratospheric residence (1–3 years). - Sulfur burden and deposition: Estimated ME sulfate deposition at NGRIP ~11 kg km⁻² (NEEM ~9 kg km⁻²), far lower than Tambora (~40 kg km⁻²), reflecting lower sulfur yield and injection height for the ME. - Injection vs. plume height: Despite modeled maximum plume heights of 30–40 km, the sulfur injection likely remained in the UT/LS, potentially due to a very high mass eruption rate (~4 × 10⁶ kg s⁻¹) leading to >15 km differences between plume top and sulfur injection height, as observed in analogous datasets. - Historical constraint: The inferred 1–2 month phase separation makes it unlikely that thunderous sounds reported in Japan on 7 Feb 947 (~14 weeks after 1 Nov 946) were associated with a second ME phase.

The ice-core glaciochemical and cryptotephra records resolve the phasing of the Millennium Eruption, demonstrating two compositionally distinct depositional pulses separated by approximately one month. This supports a hiatus between rhyolitic and trachytic phases rather than a single continuous episode. The lack of sulfur mass-independent fractionation and the near-synchronous sulfate–particle deposition indicate sulfate aerosol formation largely within the upper troposphere/lower stratosphere, implying shorter atmospheric residence times and faster scavenging. Combined with the winter timing of the aerosol veil, these factors reduce radiative forcing during the critical NH summer season, explaining the absence of a robust summer cooling signal in tree-ring reconstructions and historical records. The findings reconcile the large eruptive magnitude with the muted climatic response by constraining the altitude and duration of sulfate aerosol presence and by quantifying the short hiatus between phases. They also suggest that historical February 947 acoustic reports are unlikely to record a second ME phase, given transport and deposition timescales inferred from the ice core.

High-resolution analyses of the NGRIP1 ice core reveal two discrete cryptotephra and glaciochemical spikes corresponding to rhyolitic and trachytic phases of the 946 CE Millennium Eruption, separated by about one month. Sulfur isotope evidence (negligible Δ33S) and near-synchronous microparticle–sulfate deposition indicate UT/LS aerosol formation with short atmospheric lifetimes. Most sulfate deposition occurred in boreal winter–spring, and total sulfur loading and deposition were modest relative to major stratospheric events, explaining the limited climatic impact. The refined eruption timeline and inferred sulfur injection heights improve understanding of caldera eruption hazards and their climatic implications. Future work should aim to: (i) better constrain monthly to sub-monthly accumulation histories in northern Greenland to reduce phasing uncertainties; (ii) integrate additional ice-core sites for spatial deposition patterns; and (iii) couple tephra and sulfur isotope constraints with aerosol–climate models to quantify radiative forcing and dynamical responses more precisely.

- Accumulation model constraints: Monthly fractions are derived from a short (1991–1995) observational record, assuming applicability to 946–947 CE; interannual variability and non-normality introduce uncertainty. - Assumed eruption start and transport time: Calculations adopt a 1 Nov 946 start date and a 14 ± 7 day transit to Greenland based on modern HYSPLIT trajectories; true 10th-century meteorology and transport may differ. - Depth-to-time conversion: Uses average annual layer thickness and model-derived monthly accumulation fractions; uncertainties propagate into spike timing and interval estimates. - Sulfate deposition integration: Shoulder overlap required halving-and-doubling of the main peak area and background subtraction, which does not fully capture the true exponential decay tail of volcanic aerosols. - Tephra identification: Limited shard sizes (4–38 µm) and small counts per sample may underrepresent full compositional variability; potential minor non-ME shards noted. - Injection height inference: UT/LS inference relies on S-isotope signatures and synchronicity; lack of direct observational constraints from 946 CE means some ambiguity remains regarding exact injection altitudes and hemispheric transport.