University of Bologna

Pressurized cavities are commonly used to compute ground deformation in volcanic areas: the set of available solutions is limited and in some cases the moment tensors inferred from inversion of geodetic data cannot be associated with any of the available models. Two different source models (pure tensile source, TS and mixed tensile/shear source, MS) are studied using a boundary element approach for rectangular dislocations buried in a homogeneous elastic medium employing a new C/C++ code which provides a new implementation of the dc3d Okada fortran code. Pressurized triaxial cavities are obtained assigning the overpressure in the middle of each boundary element distributed over the cavity surface. The MS model shows a moment domain very similar to triaxial ellipsoidal cavities. The TS and MS models are also compared in terms of the total volume increment limiting the analysis to cubic sources: the observed discrepancy (~10%) is interpreted in terms of the different deformation of the source interior which provides significantly different internal contributions (~30%). Comparing the MS model with a Mogi source with the some volume, the overpressure of the latter must be ~37% greater than the former, in order to obtain the same surface deformation; however the outward expansion and the inner contraction separately differ by ~±10% and the total volume increments differ only by ~2%. Thus, the density estimations for the intrusion extracted from the MS model and the Mogi model are nearly identical.

In the present paper we introduce a numerical model for the representation of displacement, strain and stress due to single forces embedded in a layered elastic half-space. The code EFGRN/EFCMP (Elastic Forces GReeN functions/Elastic Forces CoMPutation) is able to represent the mechanical effects due to pre-assigned distributions of single forces. Even if internal deformation sources can be described by distributions of equivalent body forces with vanishing resultant and moment, single forces are employed in geophysics to represent hydraulic and/or lithostatic loads, effects of internal density anomalies, and even some kind of seismic events. A distribution of single forces is also used to describe the effects of an inelastic inclusion located inside an elastic medium. In fact, the recent literature shows that poro-elastic and thermo-elastic inclusions can be represented using single forces distributed on their boundaries. EFGRN/EFCMP shares the benefits of rapid and semi-analytical calculation offered by the parent code, EFGRN/EFCMP , which is instead suitable for the representation of extended dislocation sources, as seismic faults. The present code also provides an option for computing the effects of a distribution of single forces embedded in a homogeneous half-space, by using the analytical solutions of Mindlin. Accordingly, EFGRN/EFCMP can be a valid support both for the representation of forward models of deformation sources and for the procedures of inversion of geodetic data in a layered medium. We show some applications of the code and we provide several scripts in MATLAB language which help the user to quickly start using EFGRN/EFCMP

The influence of the hydrothermal circulation on seismicity and uplift observed at the Campi Flegrei caldera (Italy) is a topic of great interest to the scientific community. Recently, Thermo-Poro-Elastic (TPE) inclusions were proposed as likely deformation sources. They are suitable to explain the mechanical effects induced by hot and pressurized hydrothermal fluids, possibly exsolved from underlying magma, and pervading an overlying brittle layer. Recent works show that a TPE inclusion located at approximately 2 km depth below the Campi Flegrei caldera significantly contributed to the large and rapid soil uplift observed during the ‘82-’84 unrest phase. In the present work we demonstrate that such a source of deformation is likely playing a role even in the current unrest phase, which is characterized by a much lower uplift-rate with respect to the one occurred in the previous unrest phase. We will show that the time-series of soil uplift observed in the last 18 years can be reproduced by assuming the reactivation of the same deformation source responsible of the ‘82-’84 unrest located within a shallow brittle layer at about 2 km depth. The presence of a brittle layer has been evidenced in the past by tomographic studies and is confirmed by a sharp variation of the b-value at the corresponding depth.We believe that our results provide very important insights and evidences, supporting the existence and the importance of an active thermo-poro-elastic deformation source, which can be useful for understanding the unrest of the Campi Flegrei caldera, from both a scientific and geohazard perspective.

Monitoring of quiescent volcanoes, such as Campi Flegrei (Italy), involves the measurement of geochemical and geophysical parameters that are expected to change as eruptive conditions approach. Some of these changes are associated with the hydrothermal activity that is driven by the release of heat and magmatic fluids. This work focuses on the properties of the porous medium and on their effects on the signals generated by the circulating fluids. The TOUGH2 porous media flow model is applied to simulate a shallow hydrothermal system fed by a source of magmatic fluids. The simulated activity of the source, with periods of increased fluid discharge, generates changes in gas composition, gravity, and ground deformation. The same boundary conditions and source activity were applied to simulate the evolution of homogeneous and heterogeneous systems, characterized by different rock properties. Phase distribution, fluid composition, and the related signals depend on the nature and properties of the rock sequence through which the fluids propagate. Results show that the distribution of porosity and permeability affects all the observable parameters, controlling the timing and the amplitude of their changes through space and time. Preferential pathways for fluid ascent favor a faster evolution, with larger changes near permeable channels. Slower changes over wider areas characterize less permeable systems. These results imply that monitoring signals do not simply reflect the evolution of the magmatic system: intervening rocks leave a marked signature that should be taken into account when monitoring data are used to infer system conditions at depth.

Dikes and sills are the moving building blocks of the plumbing system of volcanoes and play a fundamental role in the accretionary processes of the crust. They nucleate, propagate, halt, resume propagation, and sometimes change trajectory with drastic implications for the outcome of eruptions (Sigmundsson et al., 2010). Their dynamics is still poorly understood, in particular when different external influencing factors are interacting. Here we apply a boundary element model to study dike and sill formation, propagation and arrest in different scenarios. We model dikes as finite batches of compressible fluid magma, propagating quasi-statically in an elastic medium, and calculate their trajectories by maximising the energy release of the magma-rock system. We consider dike propagation in presence of density layering, of density plus rigidity layering, of a weakly welded interface between layers, under the action of an external stress field (of tectonic or topographic origin). Our simulations predict sill formation in several situations: i) when a horizontal weak interface is met by a propagating dike; ii) when a sufficiently high compressive tectonic environment is experienced by the ascending dike and iii) in case a dike, starting below a volcanic edifice, propagates away from the topographic load with a low dip angle. We find that dikes halt and stack when they become negatively buoyant and when they propagate with low overpressure at their upper tip toward a topographic load. Neutral buoyancy by itself cannot induce dikes to turn into sills, as previously suggested.

Campi Flegrei caldera (Italy) was affected by a new unrest phase during 2011-13. We exploit two COSMO-SkyMed datasets to map the deformation field, obtaining displacement rates reaching 9 cm/yr in 2012 in the caldera center. The resulting dataset is fitted in a geophysical inversion framework using finite element forward models to account for the 3D heterogeneous medium. The best-fit model is a North dipping mixed-mode dislocation source lying at ~5 km depth. The driving mechanism is ascribable to magma input into the source of the large 1982-84 unrest (since similar source characteristics were inferred) that generates initial inflation followed by additional shear slip accompanying the extension of crack tips. The history and the current state of the system indicate that Campi Flegrei is able to erupt again, and the advanced techniques adopted provide useful information for short-term forecasting.

We develop a mathematical model describing dyke propagation in proximity of an elastic discontinuity of the embedding medium. The dyke is modelled as a fluid-filled crack in plane strain configuration employing the boundary element method. The pressure gradient along the crack is assumed proportional to the difference between the densities of the host rock and the fluid. Mass conservation is imposed during propagation and fluid compressibility is taken into account. The path followed by the crack is found by maximizing the total energy release, given by the sum of the elastic and gravitational contributions. The mathematical simulations provide a sort of ‘refraction phenomenon’, that is a sudden change in the direction of propagation when the crack crosses the boundary separating different rigidities: if the dyke enters a softer medium, its path deviates towards the vertical, if the dyke enters a harder medium its path deviates away from the vertical and may even become arrested as a horizontal sill along the interface, if the rigidity contrast is large. Gravitational energy plays a major role during propagation; in particular, in proximity of layer boundaries, this role is enhanced by the shift of the centre of mass due to changes of dyke shape. Mathematical results were validated by laboratory experiments performed injecting tilted air-filled cracks through gelatin layers with different rigidities.

Deformation sources in volcanic areas are generally modeled in terms of pressurized tri-axial ellipsoids or pressurized cracks with simple geometrical shapes, embedded in a homogeneous half-space. However, the assumption of a particular source mechanism and the neglect of medium heterogeneities bias significantly the estimate of source parameters. A more general approach describes the deformation source in terms of a suitable moment tensor. Ratios between moment tensor eigenvalues are shown to provide a strong diagnostic tool for the physical interpretation of the deformation source and medium heterogeneities may be accounted for through 3D finite element computations. Leveling and EDM data, collected during the 1982–84 unrest episode at Campi Flegrei (Italy), are employed to retrieve the complete moment tensor according to a Bayesian inversion procedure, considering the heterogeneous elastic structure of the volcanic area. Best fitting moment tensors are found to be incompatible with any pressurized ellipsoid or crack. Taking into account the deflation of a deeper magma reservoir, which accompanies the inflation of a shallower source, data fit improves considerably but the retrieved moment tensor of the shallow source is found to be incompatible with pressurized ellipsoids, still. Looking for alternative physical models of the dislocation source, we find that the best fit moment tensor can be best interpreted in terms of a mixed mode (shear and tensile) dislocation at 5.5 km depth, striking EW and dipping by ~25°–30° to the North. Gravity changes are found to be compatible with the intrusion of ~60–70·10^6 m^3 of volatile rich magma with density ~2400 kg/m^3.

The 2012 Emilia Romagna (Italy) seismic sequence has been extensively studied given the occurrence of two mainshocks, both temporally and spatially close to each other. The recent literature accounts for several fault models, obtained with different inversion methods and different datasets. Several authors investigated the possibility that the second event was triggered by the first mainshock with elusive results. In this work, we consider all the available InSAR and GPS datasets and two planar fault geometries, which are based on both seismological and geological constraints. We account for a layered, elastic half-space hosting the dislocation and compare the slip distribution resulting from the inversion and the related changes in Coulomb Failure Function (CFF) obtained with both a homogeneous and layered half-space. Finally, we focus on the interaction between the two main events, discriminating the contributions of coseismic and early postseismic slip of the mainshock on the generation of the second event and discuss the spatio-temporal distribution of the seismic sequence. When accounting for both InSAR and GPS geodetic data we are able to reproduce a detailed coseismic slip distribution for the two mainshocks that is in accordance with the overall aftershock seismicity distribution. Furthermore, we see that an elastic medium with depth dependent rigidity better accounts for the lack of the shallow seismicity, amplifying, with respect to the homogeneous case, the mechanical interaction of the two mainshocks

The ground deformation produced by a spherical overpressure source in a heterogeneous elastic and/or viscoelastic medium is investigated by numerical models based on the finite element method. Sources are assumed to be located at different depths beneath Mount Etna, Sicily, Italy, the structure of which is approximated as axially symmetric. Finite element modelling allows to incorporate in the analysis realistic features such as topographic relief and the laterally heterogeneous multi-layered structure inferred from seismic tomography. In order to avoid introducing artifacts in the solution, great care was taken to calibrate the computational domain necessary to reproduce analytical results accurately. An elastic analysis, performed initially, shows significant changes of the deformation field with respect to homogeneous half-space solutions: topography induces slight but detectable changes in the deformation field; in particular the maximum value of the vertical component is shifted away from the symmetry axis. When introducing the elastic heterogeneities, the ground deformation is found to be more confined to the proximity of the axis and its amplitude is mostly sensitive to the presence of low rigidity layers above the source. The ratio of maximum radial to vertical deformation is significantly larger for deeper sources. A further development of the model includes the study of inelastic properties assuming a Maxwell viscoelastic rheology for different layers. If the viscoelastic rheology is applied only to layers deeper than the source, the solutions are affected in different ways according to the distance of the source from the viscoelastic layer. If a viscoelastic layer is present above the source, a very large amplification (by more than 100%) of the surface deformation is predicted by the model; moreover, uplift transients are found to be followed by subsidence, without invoking any decrease in source overpressure. The most striking effects are observed when the source is embedded within a viscoelastic layer: in this case a static equilibrium configuration is not attained and, in the long term, both components of deformation reverse their signs in proximity to the axis. Furthermore, the surface deformation becomes nearly independent of source depth, in the long term. Simple physical explanations are proposed for the different cases.