Eruption plumes extended more than 30 km in altitude in both phases of the Millennium eruption of Paektu (Changbaishan) volcano

Large explosive eruptions disperse tephra widely and inject volatiles into the atmosphere, potentially perturbing climate. However, climate impacts vary with sulfur load, latitude, and timing. There are few well-constrained large eruptions to empirically link eruption parameters to climate effects. The 946 CE Millennium eruption (ME) of Paektu (Changbaishan) produced widespread tephra, including in Greenland (>7000 km away), implying a large magnitude (~M6.7), yet climate proxies show limited impact. This discrepancy raises questions about whether eruption parameters—volume, plume height, duration, and volatile release—have been accurately constrained. This study integrates tephra thickness and composition data with a dual-step inversion approach using the FALL3D model to reconstruct the dynamics and volumes of the two compositionally distinct phases (comendite and trachyte) of the ME, and to reassess their volatile emissions and potential climate impact.

- Distal tephra (B-Tm) from Paektu was first identified in northern Japan and linked to Paektu by geochemical fingerprinting; occurrences are widespread across Japan, NE China, Russia, and as cryptotephra in Greenland ice cores. - Volcanic history of Paektu includes multiple large explosive eruptions over the last 80 ka; the 946 CE ME was Plinian. - Dating: Wiggle-matched radiocarbon and cosmogenic markers (AD 775 event) in trees buried by PDCs give 946 ± 3 CE; Greenland ice core chronology (after correction) is consistent; historical records (Nov 946 and Feb 947 CE) suggest two phases or extended activity. - Proximal stratigraphy shows two phases: comenditic (light grey) and trachytic (dark), each producing fallout and PDCs; distal records often contain both phases indicated by glass chemistry. - Previous volume estimates: Fallout 13.4–37.4 km³ DRE for entire deposit; PDCs ~4.1–5.7 km³ DRE extra-caldera plus ~2.1 km³ intra-caldera. Earlier works suggested comendite dominance, but proximal thicknesses and distal glass proportions indicate comparable fallout volumes for both phases. - Climate impact proxies (tree rings, Greenland nssS ~30 ppb) indicate limited cooling and modest sulfur injection compared to eruptions like Tambora. This suggests either low syn-eruptive volatile release or tropospheric plumes, warranting re-evaluation of plume heights and volatile budgets.

Data and assumptions: - Compiled tephra thicknesses and compositions across >700,000 km² east of Paektu (83 total thickness points; 23 points where comendite vs. trachyte phase contributions can be separated by glass chemistry). For cryptotephra, 25,000 shards/g dry sediment was equated to ~1 mm thickness based on Lake Suigetsu calibrations. - Proximal sites (< ~40 km) excluded in inversions due to model limitations; deposit bulk density assumed 1000 kg/m³. Dual-step inversion strategy: 1) First-step inversion (semi-analytical; PARFIT-2.3.1 with embedded HAZMAP): - Searched >18,000 daily-averaged ERA5 wind profiles (1961–2010) to identify compatible meteorology and ranges for Eruption Source Parameters (ESP). - Weighted least squares with proportional and statistical weights; constrained TEM, column height, and total grain size distribution (TGSD) using limited data. TGSD parameterised following Costa et al. based on rhyolitic viscosity (~10^7 Pa s) and plume height ~30 km; densities and shapes per Bonadonna & Phillips; ash aggregation parameterised (Cornell et al.). Explored ranges: plume height 15–45 km; Suzuki coefficient 1–9; diffusion 1000–30,000 m²/s; aggregate density 100–1000 kg/m³. - Outcomes (examples of optimal cases): comendite ~16–20 km³ DRE with 22–26 km columns; trachyte ~17–36 km³ DRE with 26–35 km columns; combined solution ~15–19 km³ DRE with 38–45 km columns (noting large uncertainty and model simplifications). 2) Second-step ensemble-based inversion (FALL3D-8.2 with GNC method): - For each optimal meteorological configuration from step 1, ran six ensemble simulations per phase (ensemble size 128) using FALL3D-8.2 over a 3-D domain (grid 150×150×60; resolution ~0.23°×0.16°). Emission profile: Suzuki distribution; 10 ash bins + 1 aggregate bin; TGSD per Costa et al.; aggregate density central value 400 kg/m³ with perturbations. - Perturbations included: column height (±25%), source start time window, source duration (±3 h), Suzuki A (±50%), wind components (Gaussian ±25%). MER parameterised per Degruyter & Bonadonna. - Data assimilation: used 23 phase-specific observations to reconstruct deposit mass loading fields; validated total deposit against 83 observations. Metrics included RMSE (weighted by observation errors), Bias, MAE, SMAPE, and composite indices. Optimal meteorology for comendite: start 23 Sep 1986; for trachyte: 30 Oct 1992 (alternatives 21 Apr 1994; 29 Nov 1978 depending on metric). - Inversion provided weight factors for ensemble members to reconstruct source terms and time-varying emission rate profiles; durations and peak mass eruption rates derived from integrated source terms. Model evaluation and uncertainties: - RMSE close to 1 indicates deviations comparable to observation errors; small positive bias suggests slight underestimation within error bounds; >68% observations within a factor ~2.4 of analyses. Mean volume uncertainty ~factor 2.4 for each phase.

- Reconstructed fallout volumes (DRE): comendite phase ~7.2 km³ (range 3–17); trachyte phase ~9.3 km³ (range 4–22). - Total magma volume including fallout and previously estimated PDCs: ~23 km³ DRE (range ~13–47), corresponding to eruption magnitude M6.5–M7.0; comparable to estimated caldera subsidence volume (~22 km³). - Eruption dynamics: two major Plinian phases each <1 day (comendite ~15 h; trachyte ~20 h); peak mass flow rates up to ~4×10^8 kg/s. - Plume heights: both phases reached ~30–40 km, well into the stratosphere (tropopause ~12 km at this latitude), implying efficient stratospheric injection of eruptive mixture during the climatic phases. - Phase contributions: fallout volumes of comendite and trachyte phases were comparable, aligning with distal glass-shard proportions and proximal thickness observations. - Volatile release (scaled using comendite data): syn-eruptive sulfur ~0.5–2.8 Tg S (≈1.0–5.6 Tg SO2), consistent with Greenland ice-core nssS and muted climate signals; pre-eruptive losses estimated at 5–30 Tg S, 6–32 Tg F, 2–15 Tg Cl for comendite, likely not transported effectively to the stratosphere. - Climate impact: despite very high plumes, limited syn-eruptive halogen and sulfur injection, high latitude and seasonal timing likely contributed to negligible hemispheric cooling; impacts likely regional (ash loading, lahars, acid rain, fluoride hazards).

The study resolves key parameters for the ME by explicitly separating and modelling the comenditic and trachytic fallout phases. Comparable fallout volumes for both phases, brief high-intensity durations, and extreme plume heights indicate classic Plinian behaviour with umbrella clouds and efficient stratospheric injection. However, the reconstructed syn-eruptive sulfur budgets are modest, explaining the limited climate response in proxy records (tree rings, ice-core sulfur). At high latitudes, especially in winter, stratospheric sulfate residence and radiative forcing are influenced by insolation and circulation patterns, further dampening detectable cooling signals. The findings reconcile widespread tephra dispersal and caldera formation with low climatic impact by showing that total erupted volume was large but syn-eruptive sulfur injection was relatively low and halogens likely scavenged in lower plumes. The explicit quantification of each phase supports refining volatile budgets—particularly once trachyte-phase volatile data become available—and improves validation of dispersal models and climate impact assessments for extratropical large eruptions.

A dual-step inversion combining PARFIT (semi-analytical) and ensemble-based FALL3D modelling reconstructed the two main fallout phases of the 946 CE Paektu ME. Best estimates indicate ~7 km³ DRE (comendite) and ~9 km³ DRE (trachyte) fallout, with columns reaching 30–40 km and peak mass flow rates ~4×10^8 kg/s. Incorporating prior PDC estimates yields a total erupted magma volume of ~23 km³ DRE (M6.7), consistent with caldera collapse volumes. Scaled petrological constraints imply low syn-eruptive sulfur release, aligning with minimal climatic perturbation in proxies. The work demonstrates that ensemble-based inversions can recover time-evolving source terms and uncertainties for past events. Future research should acquire volatile data specific to the trachyte phase to refine sulfur and halogen budgets, improve grain-size constraints, and explore seasonal and latitudinal dependencies of aerosol–climate interactions for extratropical eruptions.

- Data sparsity and uneven spatial coverage, especially for phase-resolved thicknesses, limit inversion constraints. - Grain-size data are scarce; TGSD was parameterised, introducing model dependency and potential bias. - First-step PARFIT uses simplified physics and cannot reliably constrain column height and wind intensity from thickness alone; proximal data (<40 km) excluded due to model limitations. - Conversion of cryptotephra shard concentrations to thickness (25,000 shards/g ≈ 1 mm) introduces uncertainty. - Deposit density assumed uniform (1000 kg/m³); aggregation and particle properties parameterised. - Meteorological drivers derived from modern ERA5 reanalysis and daily averages; selected historical analog days may not perfectly represent 946–947 CE conditions. - Volatile budgets for the trachyte phase are unconstrained due to lack of melt inclusion/matrix glass data. - Resulting volume uncertainties are large (≈factor 2.4 per phase).