UNE source calibration

A UNE source-spectral model39,40,41 is implemented to calculate the seismic waveforms for a given UNE detonation size. The validity of the UNE source-spectral model is checked, and the apparent seismic moment in the model is calibrated for the UNE magnitude20. The spectral-amplitude ratios of seismic records of North Korean UNEs at common stations yield source-spectral amplitude ratios of the UNEs because raypath effects such as geometrical spreading and seismic attenuation are naturally corrected.

The spectral-amplitude ratios between the 2009 and 2013 UNEs are calculated for six azimuthal ranges, 70–95°, 95–120°, 120–145°, and 145–178° for stations on Japanese islands, and 176–200° and 200–227° for stations on the Korean Peninsula (Fig. 1(a–c)). The spectral-amplitude ratios are found to be similar for all azimuths (Fig. 1(c)). The comparison between the observed and theoretical spectral-amplitude ratios verifies the theoretical UNE source-spectral model. The overshoot parameter in the source-spectral model is ξ = 1.0520.

Quasi-observed seismic waveforms

The seismic waveforms for a certain UNE size are calculated by using the waveforms of previous UNEs after correction of the source spectra. Seismic waveforms for an m b 5.1 UNE are calculated based on those of the 2009 m b 4.7 North Korea UNE. The synthesized waveforms are compared with the observed waveforms of the 2013 m b 5.1 North Korean UNE to verify the method (Fig. 1(d)). The synthesized waveforms match the observed waveforms well. The change in frequency content with the UNE size is modeled well. Seismic waveforms for m b 7.0 UNE are calculated from both the 2009 m b 4.7 and 2013 m b 5.1 UNEs. The calculated seismic waveforms based on the 2009 UNE waveforms match those based on the 2013 UNE waveforms well to confirm the validity of the method (Fig. 1(e)). Seismic waveforms for UNEs with magnitudes of 5.0 to 7.6 are synthesized (see, supplementary materials).

Peak ground motions

The strong ground motions induced by the hypothetical UNEs around Baekdusan are inferred from strong motion attenuation curves that are calibrated to quasi-observed seismic waveforms at known stations. The strong motions attenuate with distance from the sources. The strong motion at a certain distance can be inferred from a reference strong motion attenuation equation. The strong motion attenuation equation is determined by using the seismic records of 55 earthquakes that occurred around the Korean Peninsula from 2001 to 2013. The magnitudes are M L 3.0–5.2, and the focal depths are 0.26–19.8 km. Seismic waveforms with signal-to-noise ratios of greater than 2 are selected for analysis. We analyze 7973 horizontal and 3623 vertical record sections from 163 stations. The hypocentral distances are 4–630 km.

The seismic waveforms are corrected for instrumental responses and bandpass-filtered between 0.01 and 30 Hz. The zero-to-peak ground motions (accelerations, velocities) are measured. A reference strong-motion attenuation curve is determined as a function of distance from the observed peak ground motions. The event-strength-calibration constants adjusting the levels of peak ground motions for hypothetical UNEs are determined from the quasi-observed waveforms (Fig. 2(c)) (see supplementary materials). The strong ground motions on the surface of Baekdusan volcano are determined from the strong motion attenuation curves (Fig. 3(a)).

Figure 2: Strong motion attenuation curves. (a) An example of the observed horizontal PGVs of an M L 4.8 earthquake and their fitted attenuation curve. (b) Difference between the observed and predicted horizontal PGVs, suggesting good representation of observed PGVs with the PGV attenuation curves. (c) Determined event-strength-calibration constants for PGA and PGV attenuation curves for North Korean UNEs as a function of magnitude, presenting strong linear relationships. Full size image

Figure 3: Strong motions induced by UNEs. (a) Estimation of PGV on the surface of Baekdusan volcano induced by an m b 7.0 UNE. (b) Variation of PGVs and equivalent peak dynamic stress changes as a function of magnitude. Spatial distribution of PGVs induced by (c) m b 6.0 and (d) m b 7.0 UNEs. The maps were created by using the software Generic Mapping Tools (http://gmt.soest.hawaii.edu/). Full size image

Hypothetical UNEs with magnitudes (m b ) of 5.0–7.6 are expected to produce peak ground velocities (PGVs) on the surface of Baekdusan volcano of 0.00040–0.01610 m/s with a logarithmic 95% confidence range of 0.765 in the horizontal direction, and 0.00025–0.00922 m/s with a logarithmic 95% confidence range of 0.752 in the vertical direction (Fig. 3(b); see, supplementary materials). Also, the peak ground accelerations (PGAs) are 0.02689–0.28030 m/s2 with a logarithmic 95% confidence range of 0.868 in the horizontal direction, and 0.01616–0.15308 m/s2 with a logarithmic 95% confidence range of 0.860 in the vertical direction. The horizontal and vertical PGVs reach 0.0017 and 0.0010 m/s for an m b 6.0 UNE, and 0.0069 and 0.0040 m/s for an m b 7.0 UNE (Fig. 3(c,d)). The horizontal and vertical PGAs reach 0.0683 and 0.0398 m/s2 for an m b 6.0 UNE, and 0.1684 and 0.0917 m/s2 for an m b 7.0 UNE.

Peak dynamic stress change in magma chamber

We compute the dynamic stress changes induced in the magma chamber by an m b 7.0 nuclear explosion with PyLith, which is a finite-element code for dynamic and quasi-static simulations of crustal deformation42. We consider an impulsive plane wave incident to the magma chamber in order to estimate the upper bound of the peak dynamic stress change. The plane wave is generated in the left margin of the domain. The strength of the input signal is set to produce transient waves of which peak ground motion is equivalent to that expected from quasi-observed seismic waveforms on the volcano surface (Fig. 4).

Figure 4: Dynamic stress changes in magma chamber. (a) A spherical and three spheroidal magma chamber models in a 3-D domain with annotated boundary conditions and dimensions. A mesh is presented for the magma chamber interface. (b) Snapshots of the ground-motion velocities in the crust and horizontal dynamic stress changes (σ xx ) in the spherical magma chamber with a V P /V S ratio of 1.85 at lapse times of 0.8, 1.3 and 2.0 s. (c) Peak horizontal dynamic stress changes (σ xx,max ) and peak dynamic pressure changes (p max ) in the spherical magma chamber with a V P /V S ratio of 1.85 and density of 2500 kg/m3 as a function of time. (d) Peak horizontal dynamic stress changes and pressure changes as a function of the V P /V S ratio in spherical magma chambers with three different densities. The peak horizontal dynamic stress changes are found to be nearly constant, while the peak pressure changes increased with the V P /V S ratio. Full size image

A spherical model with a radius of 3 km and three spheroidal models with a flattening of 0.5 are considered for the magma chamber (Fig. 4(a)). We assess the variation in the induced stress depending on the geometry and density of the magma chamber. The densities of the magma chamber are set to be 2500, 2580 and 2660 kg/m3 for three different compositions of the magma chamber, representing 50, 30 and 10 vol% of silicic magma compared to the typical continental crust (2700 kg/m3)43.

To represent various possible states of the magma chamber during an explosive volcanic eruption44, we consider four V P /V S ratios of 1.65, 1.75, 1.85, and 1.95 (Table 1), representing various magma compositions and pore-fluid pressures. The V P /V S ratio is low in felsic rock and high in mafic rock45. Also, the V P /V S ratio increases with pore-fluid pressure. These V P /V S ratios are equivalent to Poisson’s ratios of 0.209, 0.258, 0.293, and 0.321.

Table 1 Medium properties of four magma chamber models: V P /V S ratio, Poisson’s ratio (ν), seismic-wave velocities (V P , V S ), and bulk and shear moduli (K, μ). Full size table

The transient wavefield comprised compressional pulses satisfying the observed peak ground motions on the volcano surface to induce dynamic stress changes in the magma chamber (Fig. 4(b–d)). The radial peak dynamic stress change in the magma chamber (Δσ xx,max ) are found to be constant for V P /V S ratios, while the peak dynamic pressure change in the magma chamber (p max ) increase with the V P /V S ratio. On the other hand, Δσ xx,max and p max generally increase with the density. The spherical magma chamber model produces larger Δσ xx,max and p max than the three spheroidal models.