Modelling the 2018 Kīlauea Caldera Collapse with a joint Finite Volume Method and Distinct Element Method approach

Kīlauea’s 2018 collapse represents one of the best-monitored caldera-forming events recorded. A dense network of geodetic and seismic instrument, complimented by still and satellite-based imagery, captured the full temporal evolution of summit deformation and clearly defined distinct collapse phases. An initial pre-collapse phase was characterised by lava-lake drainage and small, elastic surface displacements, followed by a three non-elastic collapse phases, captured on GNSS stations NPIT and CALS. This detailed and well-resolved, multiphase evolution makes Kīlauea an exceptional case for testing mechanical models of caldera collapse.

Analytical and continuum-based numerical models are commonly used to relate these surface displacements to deformation sources at depth. However, their elastic or viscoelastic material assumptions limit the representation of large-strain discontinuous deformation, such as fracturing and faulting, typical of caldera collapse events. To overcome this, we use 3D Discrete Element Method (DEM) modelling, in conjunction with Finite Volume Method (FVM), to capture a transition from elastic to non-elastic (fictional-plastic) behaviour similar to that during the 2018 Kīlauea collapse event.

Sentinel-1 acquisitions between the 5^th and 14^th of May 2018 were used to compute surface displacements during the elastic, pre-collapse subsidence phase. The resulting summit subsidence provided constraints on subsurface source characteristics and were used to test a range of chamber geometries, depths and pressure states using FVM models. This approach allowed for a rapid and systematic exploration of the trade-offs among these parameters and demonstrates the non-unique elastic surface displacement solutions, consistent with the observed elastic, pre-collapse deformation at Kīlauea.

The “best-fitting” parameter combinations were then used to inform forward modelling within the 3D DEM solutions. The initial source geometry, as constrained by FVM models, had a depth of 500m, vertical axis of 2000m and horizontal axes of 1500m. As underpressure was progressively increased to 6-8 MPa, deformation transitioned from elastic into non-elastic, as characterised by host-rock fracturing and accelerated summit subsidence. The DEM model subsidence curve mimics closely that measured by GNSS at Kīlauea. Furthermore, model fracture population characteristics through time show similarly with observed earthquake magnitude distributions. This study thus highlights the capacity of 3D DEM models for capturing structural, geodetic and seismic observations during large-strain discontinuous events at volcanoes.