The existence of large scale partially molten regions in the crust where mantle derived melts can accumulate, differentiate, mix, and assimilate to eventually ascend towards the surface or crystallize to build plutonic bodies is a widely accepted concept in geology. With the impact that volcanic eruptions can have on economic infrastructure (e.g., Ajtai et al., 2010) and the security of the population (e.g., Naranjo et al., 1986) as well as the importance of magmatic systems in the creation of ore deposits (e.g., Hedenquist & Lowenstern, 1994), the comprehension of such systems is highly relevant. Yet, due to our inability to directly observe the dynamics within magmatic systems or to recreate lab scale models with realistic pressure and temperature conditions, our understanding remains incomplete. Numerical models that obey the laws of physics can however go beyond the limits of observations or laboratory experiments and can shed light on previously unexplored processes. Instead they are limited by the availability and quality of data as well as the complexity of a system. Here, we want to demonstrate the advantages of using joint interpretations of different data sets and physically consistent thermomechanical models.

Geodynamic modeling has become a powerful tool to investigate how different mechanical and thermodynamical parameters influence and control geological systems such as orogens, subduction zones, magmatic systems, basins, and other terrestrial bodies (e.g.,

<p>For the past decades, several numerical studies have successfully reproduced the concentric uplift pattern observed above the Altiplano-Puna Magma Body (APMB) in the central Andes. However, the temperature- and strain rate-dependent viscoelastoplastic rheology of rocks, the buoyancy of magma, the effects of modelling in 3D as well as the shape of the magma body have often been simplified or neglected.</p> <p>Here, we use a joint interpretation of seismic imaging and gravity anomalies to constrain location, 3D shape and density of the magma body. With the help of the thermo-mechanical finite difference code LaMEM, we then model the surface deformation and test our results against observations made by Interferometric Synthetic-Aperture Radar (InSAR) missions. This way, we gain insights into the dynamics and rheology of the present-day magmatic system and can test how a change to the current conditions (e.g., magma influx) could impact it.</p> <p>We find that only an APMB with a maximum thickness of 14 to 18 km and a corresponding density contrast to the surrounding host rock of 100 to 175 kg/m<sup>3</sup> satisfies both tomography and Bouguer data. Based on that and the chemistry of eruption products, we estimate the melt content of the APMB to be on the order of 20 - 25%. We also find that the observed uplift can be reproduced by magma-induced buoyancy forces without the need for an additional pressure source or magma injection within the mush, and that the geometry of the top of the magma body exerts a major control on the deformation pattern at the surface.</p>

<p>Advanced numerical methods and increasing computational power have enabled us to incorporate numerical forward models into geodynamic inverse frameworks. We now have several strategies to constrain the rheological properties of the crust and lithosphere. Yet, the initial geometry of geological formations (e.g., salt bodies, magma bodies, subducting slabs) and associated uncertainties are, in most cases, excluded from the inverse problem and assumed to be part of the a priori knowledge. Usually, geometrical properties remain constant, or we employ simplified bodies like planes, spheres or ellipsoids for their parameterization.</p><p>Here, we present a simple method to parameterize complex three-dimensional bodies and incorporate them into geodynamic inverse problems. Our approach enables us to automatically create an entire ensemble of initial geometries within the uncertainty limits of geophysical imaging data. This not only allows us to account for geometric uncertainties, but also enables us to investigate the sensitivity of geophysical data to the geometrical properties of the geological structures.</p><p>We present 3 areas of application for our method, covering salt diapirs, magmatic systems and subduction zones, using both synthetic and real data.</p>

Volcano deformation is most frequently interpreted in terms of models of surface deformation due to processes in magma bodies of various geometries. The most widely applied model is that of a point source of pressure embedded within a uniform elastic half space (Mogi, 1958), but a range of more advanced models and approaches exist (e.g., Fialko, Khazan, et al., 2001;Hickey et al., 2016). As liquid magma flows in/out of the these "deformation sources," they expand/contract. Most often, such magma flow is considered to cause uniform pressure change on the boundary of the magma body, and the density difference between magma and host rock is not considered specifically. It has, however, been demonstrated in a number of studies that magma buoyancy can cause significant stresses in volcano roots and contribute to failure of magma bodies (e.g., Sigmundsson et al., 2020). A particular phenomenon not considered by traditional volcano deformation models is the effect of accumulated exsolved volatiles in volcano roots and their release during eruptions.During major explosive eruptions an excess of gas may be observed, beyond that which can be explained by a petrological calculation of the original dissolved volatile amounts and the volume of erupted lavas. Excess gas was observed in the 1991 eruption of Pinatubo, Philippines and an analysis from Wallace and Gerlach (1994) showed that this could be explained by a pre-existing gas/volatile phase representing 0.7-1.3 wt% of the erupted magma. Volatile accumulation was proposed to occur in the roof zone of the system.