Inverting ground displacement to resolve traction changes on fractures in a 3D heterogeneous crust

Summary

Over the past few decades, InSAR and GNSS have significantly improved the monitoring of volcanic and tectonic processes by providing frequent, high-resolution measurements over large areas. To track magma transport and better understand the seismic cycle, these data are typically analyzed using kinematic solutions that model fractures as dislocations within an elastic, homogeneous half-space. Such simplifying assumptions can lead to inaccuracies when estimating the location and strength of subsurface sources. To address these limitations and determine more physically meaningful quantities - such as stress vectors rather than dislocations - we propose a new method for inferring variable pressure or traction changes on fractures, through the inversion of noisy InSAR or GNSS ground displacement data. These fractures may represent faults or magma intrusions embedded in a three-dimensional, heterogeneous crust with prominent topography. Our approach is based on a fictitious domain method inspired by the extended finite element method (XFEM), which is particularly well-suited for modeling fractures that need not conform to volume meshes. The inversion framework minimizes a cost function comprising the misfit between modeled and observed displacements, along with two regularization terms applied to both the tractions and their gradients. We present an inversion algorithm based on the forward and adjoint problem. The method is reformulated for InSAR data, which is challenging due to missing data, measurement noise, and the projection of displacements into Earth-satellite line-of-sight (LOS) directions, resulting in an incomplete displacement field. A corresponding adaptation for GNSS data is also introduced. Synthetic tests demonstrate the effectiveness and robustness of the method even with a single LOS observation, or in the presence of spatial gaps and noise, as well as with few GNSS data. They show that, neglecting heterogeneities in the inversion may bias the distribution of pressure and their amplitude. Regularization parameters are determined using either L-curve analysis or a compromise function, both providing reliable results. We confirm that log-log L-curves, and not lin-lin L-curves, are required for appropriate regularization selection. As source depth or noise levels increase, the method favors stronger regularization, yielding smoother solutions. We apply the method to the May 2016 intrusion and eruption at Piton de la Fournaise, where InSAR and GNSS data are available for the beginning and end of the intrusion and eruption. Using a homogeneous elastic medium, the results agree well with those from a previous boundary element study, that assumed a homogeneous medium and pressurization over a circular portion of a dike. Inversion results for a heterogeneous elastic medium show greater complexity in the pressure distribution within the dike, but remain consistent with our previous findings. Furthermore, by inverting traction changes rather than pressure changes, we identify a shear displacement at depth–interpreted as a flank slip occurring simultaneously with the magma intrusion.