Long-term rheology , heat budget and dynamic permeability of deforming and reacting rocks: from laboratory to geological scales

The goal of this session is to reconcile short-time/small-scale and long-time/large-scale observations, including geodynamic processes such as subduction, collision, rifting or mantle lithosphere interactions. Despite the remarkable advances in experimental rock mechanics, the implications of rock-mechanics data for large temporal and spatial scale tectonic processes are still not straightforward, since the latter are strongly controlled by local lithological stratification of the lithosphere, its thermal structure, fluid content, tectonic heritage, metamorphic reactions and deformation rates.

Mineral reactions have mechanical effects that may result in the development of pressure variations and thus are critical for interpreting microstructural and mineral composition observations. Such effects may fundamentally influence element transport properties and rheological behavior.
Here, we encourage presentations focused on the interplay between metamorphic processes and deformation on all scales, on the rheological behavior of crustal and mantle rocks and time scales of metamorphic reactions in order to discuss
(1) how and when up to GPa-level differential stress and pressure variations can be built and maintained at geological timescales and modelling of such systems,
(2) deviations from lithostatic pressure during metamorphism: fact or fiction?,
(3) the impact of deviations from lithostatic pressure on geodynamic reconstructions.
(4) the effect of porous fluid and partial melting on the long-term strength.
We therefore invite the researchers from different domains (rock mechanics, petrographic observations, geodynamic and thermo-mechanical modelling) to share their views on the way forward for improving our knowledge of the long-term rheology and chemo-thermo-mechanical behavior of the lithosphere and mantle.

Co-organized by EMRP1/GMPV7/TS2
Convener: Yury Podladchikov | Co-conveners: Shun-ichiro Karato, Evangelos Moulas, Magdalena Scheck-Wenderoth, Lucie Tajcmanova
vPICO presentations
| Fri, 30 Apr, 09:00–10:30 (CEST)

vPICO presentations: Fri, 30 Apr

Chairpersons: Shun-ichiro Karato, Lucie Tajcmanova, Magdalena Scheck-Wenderoth
Haobo Wang, Shuyun Cao, Franz Neubauer, Junyu Li, Xuemei Cheng, and Meixia Liu

Studies of crustal anatexis have given valuable insights into the evolution of metamorphism–deformation and the tectonic processes at convergent plate margins during orogeny. The transition of metatexite to diatexite migmatite records crucial information about the tectono–thermal evolution and rheology of the deep crust. Along the Ailao Shan–Red River shear zone, metatexite migmatites, diatexite migmatites and leucogranites are widely distributed within the upper amphibolite and granulite facies zones of the Diancang Shan metamorphic complex. The high–pressure granulite–facies metamorphism with mineral assemblage comprising garnet + kyanite + K–feldspar + plagioclase + biotite + quartz + melt is first recognized from the patch metatexite migmatites in the complex. Detailed petrographic evidence and phase diagram reveal that the migmatite underwent nearly isothermal decompression metamorphism, presenting a clockwise P–T path. The peak metamorphic P–T conditions are constrained by phase diagram at ca. 11 kbar and 810 °C, and the amount of melt generated during heating is up to 18 mol%. The extraction and segregation of melts are evidenced by the presence of leucosomes within migmatites and leucogranite dikes, which record the melt flow network through the crust. Zircons and monazites from migmatites record the ages of the melting episode that began at ca. 36 Ma and lasted to ca. 20 Ma. All these results are in accord with orogenic crust thickening accompanied by pervasive anatexis during the Later Eocene to the early Oligocene in the Ailao Shan–Red River shear zone. Combined with available data related to the other continental–exhumed shear zone, we propose that the crustal anatexis has an important effect on the thermal–state of deep–seated shear zones, is thus controlling the rheological behavior of the lithosphere and plays the essential role in the initial localizing of shearing in the lower crust.

How to cite: Wang, H., Cao, S., Neubauer, F., Li, J., Cheng, X., and Liu, M.: High–pressure granulite–facies metamorphism and anatexis of deep continental crust: new insights from the Cenozoic Ailaoshan–Red River shear zone, SE Asia, EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-313,, 2021.

Jonas Vanardois, Pierre Trap, and Didier Marquer

The causes for heterogeneous deformation with strain partitioning into kilometre-scale shear zone within the partially molten crust and the spatiotemporal feedback relationships between strain localization and melt organization still remain unclear. In order to tackle these questions and unravel the strain localization in a partially molten crustal scale shear zone, we used field observations and thermodynamic modelling in the Eastern Variscan Shear Zone (EVSZ) located in the Aiguille-Rouge massif (Western Alps). The EVSZ is an orogen scale, 10 km wide and 600 km long, transpressional high strain corridor recognized in the French External Crystalline Massifs. The EVSZ affected the partially-molten late Variscan crust during late Carboniferous times (340-300 Ma). In this contribution we present a detailed field-map survey of the mid- and lower crusts focussed on the partitioning and strain pattern in the Aiguille-Rouge EVSZ. Detailed mapping revealed that high-strain deformation domains and orthogneiss occurrences are spatially related. New petrological, thermobarometrical and LA-ICP-MS dating also better constrain the P-T-t-D evolution of the partially molten crust along the EVSZ. Field observations and P-T pseudosection calculations show that among the three dominant lithologies forming the mid- and lower-crusts, i.e. metapelite, metagreywacke and orthogneiss, the latter is the most fertile if considering H2O-fluid-saturated melting. During prograde evolution at pressure between-12-15 kbar, orthogneisses reached the solidus at lower temperature and produced higher melt fraction than the metasedimentary rocks.  The water-present melting in the orthogneisses may have initiate strain localization at the end of the prograde evolution. Thus, the favoured localization of the shear zone within the metagranites is explained by a higher melt fraction than in the metapelites and metagreywackes. PTDt path and thermobarometrical modelling suggest that these transpressional deformation conditions occurred under suprasolidus conditions from at least 12 kbar to 4 kbar during a near isothermal decompression. During this cooling path, while crystallization of anatectic melts might have provoked strain hardening in the orthogneisses, a strength decrease might be controlled by a higher proportion of micas in metapelites and metagreywackes as suggested by forward modelling of modal proportion of mica. This change in the nature of the weakest phase, starting with melt in metagranites and followed by micas in metasedimentary rocks, seems to control the progressive localization and broadening of the crustal scale shear zone during clockwise P-T-t path. Our results suggest that H2O-fluid-saturated melting of metagranites has a first order rheological impact on the birth and growth of the orogen scale shear zone in the lower continental crust.

How to cite: Vanardois, J., Trap, P., and Marquer, D.: Orogen scale shear zone development in partially Variscan molten crust. The critical role of water-present melting in metagranite inducing strain localization., EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-1360,, 2021.

Konstantin Huber, Johannes C. Vrijmoed, and Timm John

Serpentinite dehydration in subduction zones plays an important role in Earth’s deep water cycle. In order to keep this water cycle in balance, an efficient rock dehydration mechanism at depth is needed to keep pace with loss of ocean water due to subduction of hydrated oceanic lithosphere. Field observations in non-deformed meta-serpentinites in Erro Tobbio, Ligurian Alps, show that serpentinite dehydration at depth occurs by a channelized vein network rather than pervasive flow. The mineral assemblage in the veins is characterized by a high abundance of metamorphic olivine. Plümper et al. (2017) showed that on small scales (μm-mm) the formation of these veins is controlled by intrinsic chemical heterogeneities in the rock. Field observations suggest that on larger scales the fluid escape is governed by mechanical processes such as hydraulic fracturing. On small scales, where dehydration is chemically controlled, reactive fluid flow is an important process because changes in the fluid chemistry may trigger or hinder further dehydration reactions in the rock. Because of its high solubility and high abundance as a rock forming component, Si might be a key metasomatic agent for first-order effects on the dehydration process.

Following the approach of Beinlich et al. (2020) we extended the model of Plümper et al. (2017) to a reactive fluid flow model for serpentinite dehydration that accounts for the Si content of the fluid. As input for our model we use mineral chemical data of non-dehydrated serpentinites from the Mirdita ophiolite in Albania that are representative for serpentinized oceanic lithosphere that enters a subduction zone, hence has not experienced any subduction-related metamorphic processes. The results of our model suggest that the high abundance of metamorphic olivine observed in the Erro Tobbio meta-serpentinites hence the purification towards a olivine-dominated assemblage is the result of interaction with an external fluid in the veins after they have been formed from the intrinsic chemical heterogeneities.


  • Beinlich, A. et al. (2020). “Instantaneous rock transformations in the deep crust driven by
    reactive fluid flow”. In: Nature Geoscience 13.4, pp. 307–311. doi: 10.1038/s41561-
  • Plümper, O. et al. (2017). “Fluid escape from subduction zones controlled by channel-
    forming reactive porosity”. In: Nature Geoscience 10.2, pp. 150–156. doi: 10.1038/

How to cite: Huber, K., Vrijmoed, J. C., and John, T.: Olivine enrichment in dehydration veins in serpentinites by reactive fluid flow, EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-9813,, 2021.

junyu Li, shunyun Cao, Xuemei Cheng, Haobo Wang, and Wenxuan Li

Adakite‐like potassic rocks are widespread in post-collisional settings and provide potential insights into deep crustal or crust-mantle interaction processes including asthenosphere upwelling, partial melting, lower crustal flow, thickening and collapse of the overthickened orogen. However, petrogenesis and compositional variation of these adakite‐like potassic rocks and their implications are still controversial. Potassic magmatic rocks are abundant developed in the Jinshajiang–Ailaoshan tectono-magmatic belt that stretches from eastern Tibet over western Yunnan to Vietnam. Integrated studies of structure, geochronology, mineral compositions and geochemistry indicate adakite-like potassic rocks with different deformation are exposed along the Ailaoshan-Red River shear zone. The potassic felsic rocks formed by mixing and partial melting between enriched mantle-derived ultrapotassic and thickened ancient crust-derived magmas. The mixing of the mafic and felsic melts and their extended fractional crystallization of plagioclase, K-feldspar, hornblende and biotite gave rise to the potassic magmatic rocks. Zircon geochronology provide chronological markers for emplacement at 35–37 Ma of these adakite-like potassic rocks along the shear zone. Temperature and pressure calculated by amphibole-plagioclase thermobarometry range from 3.5 to 5.9 kbar and 650 to 750 ℃, respectively, and average emplacement depths of ca. 18 km for granodiorite within this suite. In combination with the results of the Cenozoic potassic magmatism in the Jinshajiang–Ailaoshan tectono-magmatic belt, we suggest that in addition to partial melting of the thickened ancient continental crust, magma underplating and subsequent crust-mantle mixing beneath the ancient continental crust have also played an important role in crustal reworking and strongly affected the rheological properties and density of rocks. The exhumation underlines the role of lateral motion of the Ailaoshan-Red River shear zone initiation by potassic magma-assisted rheological weakening and exhumation at high ambient temperatures within the shear zone.

How to cite: Li, J., Cao, S., Cheng, X., Wang, H., and Li, W.: Structure and spatial-temporal relationship of potassic magmatism linked to a continental-scale strike-slip shear zone and post-collisional Cenozoic extension, EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-3884,, 2021.

Benjamin Hess and Jay Ague

Thermodynamic modeling in active tectonic settings typically makes the assumption that stress is equal in all directions. This allows for the application of classical equilibrium thermodynamics. In contrast, geodynamic modeling indicates that differential, or non-hydrostatic, stresses are widespread. Non-hydrostatic equilibrium thermodynamics have been developed by past workers [1], but their application to geological systems has generated controversy in recent years [2-5]. Therefore, we seek to clarify how stress influences the chemical potential of non-hydrostatically stressed elastic solids. To quantify this, we consider the effects of stress variation on the equilibrium between the single-component polymorph pairs of kyanite/sillimanite, quartz/coesite, calcite/aragonite, and diamond/graphite.

The stress on each interface of a solid can be decomposed into components normal to the interface and parallel to the interface. In our work, we determine the shift in the temperature of equilibrium on fixed interfaces between polymorph pairs as a function of both interface-normal and interface-parallel stress variation. We find that the influence of normal stress variation on the equilibrium temperature of polymorphs is approximately two orders of magnitude greater than interface-parallel stress variation. Thus, at a fixed temperature, normal stress determines the chemical potential of a given interface to first order. Consequently, high-pressure polymorphs will preferentially form normal to the maximum stress, while low-pressure polymorphs, normal to the minimum stress.

Nonetheless, interface-parallel stress variations can meaningfully affect the stability of phases that are at or near equilibrium. We demonstrate the surprising result that for a given polymorph pair, a decrease in interface-parallel stresses can make a high-pressure polymorph more stable relative to a low-pressure polymorph on the given interface.

The effects of non-hydrostatic stress on mineral assemblages are most likely to be seen in dry systems. Many reactions in metamorphic systems are fluid-mediated, and fluids cannot sustain non-hydrostatic stress. Consequently, in systems with interconnected, fluid-filled porosity, mineral assemblages will tend to form at a pressure approximately equal to the fluid pressure. In contrast, in dry systems all reactions occur directly between solids which can sustain non-hydrostatic stress. This facilitates the application of non-hydrostatic thermodynamics. Consequently, dry rocks containing polymorphs such as such as quartzites, marbles, and peridotites represent ideal lithologies for the testing and application of these concepts. By influencing the chemical potential of solid interfaces, non-hydrostatic stress alters the thermodynamic driving force and subsequent kinetics of polymorphic reactions. This likely results in preferential orientations of polymorphs which could influence seismic anisotropy and potentially generate seismicity.

[1] Larché, F., & Cahn, J. W. (1985). Acta Metallurgica, 33(3), 331-357.

[2] Hobbs, B. E., & Ord, A. (2016). Earth-Science Reviews, 163, 190-233.

[3] Powell, R., Evans, K. A., Green, E. C. R., & White, R. W. (2018). Journal of Metamorphic Petrology, 36(4), 419-438.

[4] Tajčmanová, L., Podladchikov, Y., Powell, R., Moulas, E., Vrijmoed, J. C., & Connolly, J. A. D. (2014). Journal of Metamorphic Petrology, 32(2), 195-207.

[5] Wheeler, J. (2018). Journal of Metamorphic Petrology, 36(4), 439-461.

How to cite: Hess, B. and Ague, J.: Quantifying the influence of non-hydrostatic stress on polymorph equilibria, EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-8029,, 2021.

Reinier van Noort, Lawrence Hongliang Wang, and Viktoriya Yarushina

Understanding fluid flow patterns in the shallow and deep earth is one of the major challenges of modern earth sciences. Fluid flow may be slow and pervasive, or fast and focused. In the deep earth, focused fluid flow may result in, for example, dikes, veins, volcanic diatremes and gas venting systems. In the shallow Earth, focused fluid flow can be found in the form of fluid escape pipes and gas conducting chimneys, mud volcanoes, sand injectites, pockmarks, hydrothermal vent complexes, etc.

Focused fluid flow has been reproduced in visco-plastic models of flow through porous materials. However, the mechanisms that cause fluid flow to focus along such relatively narrow channels, with transiently elevated permeability, have not been investigated thoroughly in experiments. We have carried out experiments in a transparent Hele-Shaw cell. In our experiments, a hydrous fluid is injected into an aggregate of viscous grains, and the mechanisms by which this injected fluid flows are recorded using a digital camera. Our experiments demonstrate a dependence of fluid flow mechanisms on the injection rate. At low injection rate, we observe the formation of a slowly-rising diapir. As this diapir slowly rises through the porous medium, it is fed by transient, focused fluid flow following the path of the rising diapir. Once the diapir escapes through the surface of our aggregate, continued fluid flow through the porous aggregate is focused and transient. At high injection rate, instead of a diapir fed by focused fluid flow, an open channel forms as a result of local fluidization of the granular material.

Our experimental observations are interpreted through visco-plastic models simulating the experimental conditions. These numerical models can reproduce the diapirs observed in our experiments at low flow rate by assuming flow through a layered porous aggregate, with a layer with relatively high bulk viscosity overlying a layer with relatively low bulk viscosity. For low injection rates, such a model reproduces focused fluid flow in the low-viscosity layer, that feeds into a slowly rising diapir in the high-viscosity layer. This model observation thus suggests that the passage of the rising diapir in our experiments leaves a trail, where the aggregate bulk viscosity is lowered and along which ongoing fluid flow can focus transiently.

How to cite: van Noort, R., Wang, L. H., and Yarushina, V.: Visualisation of fluid flow mechanisms through a viscous-porous rock-analogue medium – experiment and model results, EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-8341,, 2021.

Lucie Tajcmanova, Yury Podladchikov, and Evangelos Moulas

Quantifying natural processes that shape our planet is a key to understanding the geological observations. Many phenomena in the Earth are not in thermodynamic equilibrium. Cooling of the Earth, mantle convection, mountain building are examples of dynamic processes that evolve in time and space and are driven by gradients. During those irreversible processes, entropy is produced. In petrology, several thermodynamic approaches have been suggested to quantify systems under chemical and mechanical gradients. Yet, their thermodynamic admissibility has not been investigated in detail. Here, we focus on a fundamental, though not yet unequivocally answered, question: which thermodynamic formulation for petrological systems under gradients is appropriate – mass or molar?  We provide a comparison of both thermodynamic formulations for chemical diffusion flux, applying the positive entropy production principle as a necessary admissibility condition. Furthermore, we show that the inappropriate solution has dramatic consequences for understanding the key processes in petrology, such as chemical diffusion in the presence of stress gradients.

How to cite: Tajcmanova, L., Podladchikov, Y., and Moulas, E.: The choice of a thermodynamic formulation dramatically affects modelled chemical zoning in minerals, EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-15826,, 2021.

Stefan Markus Schmalholz, Evangelos Moulas, and Yuri Podladchikov

Melting is a major process of plate tectonics, affecting divergent and convergent plate boundaries. Melting of rock is also a typical example of a coupled geological process, in which the associated transformation affects the heat transfer via the latent heat of fusion and the rock deformation via the volume change. However, petrological studies on melting usually focus on chemical aspects, such as differentiation of involved components, thermal studies usually focus only on the impact of latent heat on heat transfer, such as done in the classical Stefan problem of solidification. Similarly, studies focusing on lithosphere and mantle deformation usually only consider the impact on the effective viscosity, such as weakening due to partial melting, or the impact on buoyancy due to density changes. Many studies do, therefore, not consider coupling of melting, heat transfer and rock deformation. Indeed, a common assumption is that rock pressure, or mean stress, remains lithostatic during melting. While this assumption is attractive due to its simplicity, it is against the common knowledge derived from physical experiments and the well-established mechanical theories. Furthermore, theoretical models of melt migration would not work if pressure is everywhere lithostatic, or hydrostatic, because melt migration is driven by local deviations from the static stress state.

Here, we present simple mathematical models based on the fundamental laws of physics and thermodynamics (e.g. conservation of mass, momentum and energy) to study the fundamental coupling of melting, heat transfer and rock deformation, and to quantify dynamic pressure variations due to melting. We show both analytical and numerical solutions for these models. We discuss applications of these solutions to experiments and geological observations and estimate magnitudes of dynamic pressure resulting from melting under natural conditions.

How to cite: Schmalholz, S. M., Moulas, E., and Podladchikov, Y.: Melting and dynamic pressure: coupling of reactions, heat transfer and deformation, EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-10210,, 2021.

Denis Anuprienko, Viktoriya Yarushina, and Yury Podladchikov

Understanding interactions between rock and fluids is important for many applications including CO2 storage in the subsurface. Today significant effort is aimed at research on CO2 flow through low-permeable shale formations. In some experiments, CO2 is injected in a shale sample at a constant rate, and the upstream pressure exhibits rise until a certain moment followed by a decline, representing the so called breakthrough phenomenon. After the breakthrough, downstream flux significantly rises. This behavior was thought to be the result of fracture occurence or mechanical effects. 

Here, we present a 3D numerical model of flow through experiments in shale. Our model accounts for poroelastic compaction/decompaction of shale, its time-dependent permeability, and two-phase flow, the fluid phases being CO2 and air. The model also accounts for a capillary entry pressure threshold observed in experiments. The key feature of the model are saturation-based relative permeabilities which result in sharp overall permeability increases as the CO2 moves through the shale sample. The model is implemented for 3D calculations with the finite volume method. Our results show that CO2 breakthrough is a natural outcome of two-phase fluid flow dynamics and does not need a fracture to exhibit pressure behavior observed in experiments.

How to cite: Anuprienko, D., Yarushina, V., and Podladchikov, Y.: Effect of capillary pressure and geomechanics on multiphase fluid flow in rocks, EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-10244,, 2021.

Yury Podladchikov, Viktoriya Yarushina, and Benjamin Malvoisin

Deformation, chemical reactions, and fluid flow in the geological materials are coupled processes. While some reactions are thought to be a consequence of fluid assisted dissolution on the stressed mineral surfaces and precipitation on the free surface, other reactions are caused by mineral replacement wherein a less stable mineral phase is replaced by a more stable phase, involving a change in solid volume and build-up of stresses on grain contacts, also known as a force of crystallization. Most of the existing models of chemical reactions coupled with fluid transport either assume dissolution-precipitation process or mineral growth in rocks. However, dissolution-precipitation models used together with fluid flow modelling predict a very limited extent of reaction hampered by pore clogging and blocking of reactive surfaces, which will stop reaction progress due to the limited supply of fluid to reactive surfaces. Yet, field observations report that natural rocks can undergo 100% hydration/carbonation. Mineral growth models, on the other hand, preserve solid volume but do not consider its feedback on porosity evolution. In addition, they predict the unrealistically high force of crystallization on the order of several GPa that must be developed in minerals during the reaction. Here, using a combination of effective media theory and irreversible thermodynamics approaches, we propose a new model for reaction-driven mineral expansion, which preserves porosity and limits unrealistically high build-up of the force of crystallization by allowing inelastic failure processes at the pore scale. To fully account for the coupling between reaction, deformation, and fluid flow we derive macroscopic poroviscoelastic stress-strain constitute laws, that account for chemical alteration and viscoleastic deformation of porous rocks. These constitutive equations are then used to simulate the reactive transport in porous rocks.

How to cite: Podladchikov, Y., Yarushina, V., and Malvoisin, B.: On the constitutive equations for coupled chemical reaction and deformation of porous rocks, EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-15602,, 2021.

Annelore Bessat, Sébastien Pilet, Stefan M. Schmalholz, and Yuri Podladchikov

The formation of alkaline magmas observed worldwide requires that low degree-melts, potentially formed in the asthenosphere, were able to cross the overlying lithosphere. Fracturing in the upper, brittle part of the lithosphere may help to extract this melt to the surface. However, the mechanism of extraction in the lower, ductile part of the lithosphere is still contentious. Metasomatic enrichment of the lithospheric mantle demonstrates that such low-degree melts interact with the lithosphere, but the physical aspect of this process remains unclear.

Here, we aim to better understand, first, the percolation of magma in a porous viscous medium at pressure (P) and temperature (T) conditions relevant for the base of the lithosphere, and second, the impact of chemical differentiation on melt migration. We investigate theoretically the process of melt migration employing the fundamental laws of physics and thermodynamics. We simulate melt percolation numerically with a one-dimensional (1-D) Thermo-Hydro-Mechanical-Chemical (THMC) model of porosity waves coupled with thermodynamic results obtained from numerical Gibbs energy minimisation calculations. We perform THMC modelling and Gibbs energy minimisations with self-developed numerical algorithms using MATLAB and linear programming routines. We employ a simple ternary system of Forsterite/Fayalite/Enstatite for the solid and melt. Model variables, such as solid and melt densities or mass concentrations of MgO and SiO in solid and melt, are a function of pressure (P), temperature (T) and total silica concentration of the system (X). These variables are pre-computed with Gibbs energy minimisation and implemented in the THMC porosity wave transport code via parameterized equations, determining the T-P-X dependence of the model variables.

First results show that the total silica concentration and the temperature gradient are important parameters to consider in melt migration by reactive porosity waves. We discuss results of a systematic series of 1-D simulations and we present preliminary results form a 2-D reactive porosity wave model.

How to cite: Bessat, A., Pilet, S., Schmalholz, S. M., and Podladchikov, Y.: Melt migration by reactive porosity waves, EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-2786,, 2021.

Dániel Kiss, Evangelos Moulas, Lisa Rummel, and Boris Kaus

A recent focus of studies in geodynamic modeling and magmatic petrology is to understand the coupled behavior between deformation and magmatic processes. Here, we present a 2D numerical model of an upper crustal magma (or mush) chamber in a visco-elastic host rock, with coupled thermal, mechanical and chemical processes, accounting for thermodynamically consistent material parameters. The magma chamber is isolated from deeper sources of magma (at least periodically) and it is cooling, and thus shrinking. We quantify the changes of pressure and stress around a cooling magma chamber and a warming host rock, using a compressible visco-elastic formulation, considering both simplified idealized and more complex and realistic geometries of the magma chamber.

We present solutions based on a self-consistent system of the conservation equations for coupled thermo-mechanical-chemical processes, under the assumptions of slow (negligible inertial forces), visco-elastic deformation and constant chemical bulk composition. The thermodynamic melting/crystallization model is based on a pelitic melting model calculated with Perple_X, assuming a granitic composition and is incorporated as a look-up table. We will discuss the numerical implementation, show the results of systematic numerical simulations, and illustrate the effect of volume changes due to temperature changes (including the possibility melting and crystallization) on stress and pressure evolution in magmatic systems.

How to cite: Kiss, D., Moulas, E., Rummel, L., and Kaus, B.: 2D thermo-mechanical-chemical coupled numerical models of interactions between a cooling magma chamber and a visco-elastic host rock, EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-11661,, 2021.

Salavat Ishbulatov, Viktoriya Yarushina, and Yury Podladchikov

The reliability of geomechanical and petrophysical laboratory experiments depends on coring operation. One of the steps where the core material undergoes critical loads is decompression during the core retrieval operation. Currently, a few numerical and analytical models simulate that process only with critical simplifications. The analytical solution considers only homogeneous media that neglects micro defects. FEM methods calculate slower than FDM up to several orders, simulating lifting processes with dynamic boundary conditions.

We present an axisymmetric cylindrical model of fully coupled fluid flow and elastic deformation solution by pseudo-transient numerical method. Calculation in the physical domain allows for high efficiency of parallelization on GPU, making it possible to simulate with high resolution of loading a core sample.

How to cite: Ishbulatov, S., Yarushina, V., and Podladchikov, Y.: A Numerical Simulation of Poroelastic Cylinder Decompression Problem on CUDA in an Axisymmetric Domain, EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-11955,, 2021.

Evangelos Moulas, Konstantin Zingerman, Anatoly Vershinin, Vladimir Levin, and Yury Podladchikov

The recent improvement in spectroscopic methods has allowed the detailed characterization of minerals with very high spatial resolution. Such methods allow the accurate estimation of residual pressures in mineral inclusions from exhumed metamorphic rocks. The residual inclusion pressures can be used to recast the pressure conditions during metamorphic recrystallization (e.g. Moulas et al., 2020). The most common assumptions in the aforementioned models are that 1) the rheology of the host-inclusion system is elastic, 2) the pressure was the same in the host and the inclusion phase at the time of recrystallization and, 3) the host and the inclusion can be treated as elastically isotropic phases.

In this work we focus on isotropic host-inclusion systems. Such solutions appear to be sufficient even for anisotropic minerals such as the Quartz-in-Garnet system (Bonazzi et al., 2019; Moulas et al., 2020; Thomas and Spear, 2018). In addition, numerous experimental studies show that mineral Equations-of-State (EoS) are non-linear and, therefore, mechanical solutions which consider linear-elastic host-inclusion systems may be inadequate. We present two analytical solutions for the host-inclusion problem which can be applied in systems under large strain. In the first approach we consider that the volumetric deformation of minerals is non-linear and the deviatoric stresses can be approximated by linear elasticity. The resulting solution is:

(Eq. 1)

where ΔP is the residual pressure difference, G is the shear modulus, V are the mineral volumes, superscripts h,i indicate host/inclusion and, “ini”/”fin” indicate initial and final P-T conditions respectively. This result is similar but not identical with a previously published solution (Guiraud and Powell, 2006). However, we demonstrate that Guiraud and Powell’s (2006) solution is a linearization of this formulation and its accuracy decreases with increasing pressure range.  Finally, we discuss our results in the framework of a newly-derived, fully-non-linear elastic solution that considers the effects of large finite strain in Neo-Hookean materials (Levin et al., 2020). We conclude that, for common mineral barometry applications, the effects of geometrical non linearity are minor and the application of Eq. 1 is sufficient.



Bonazzi, M., Tumiati, S., Thomas, J.B., Angel, R.J., Alvaro, M., 2019. Assessment of the reliability of elastic geobarometry with quartz inclusions. Lithos 350–351, 105201.

Guiraud, M., Powell, R., 2006. P–V–T relationships and mineral equilibria in inclusions in minerals. Earth and Planetary Science Letters 244, 683–694.

Levin, V.A., Podladchikov, Y.Y., Zingerman, K.M., 2020. An exact solution to the Lame problem for a hollow sphere for new types of nonlinear elastic materials in the case of large deformations. European Journal of Mechanics A / Solids (under revision).

Moulas, E., Kostopoulos, D., Podladchikov, Y., Chatzitheodoridis, E., Schenker, F.L., Zingerman, K.M., Pomonis, P., Tajčmanová, L., 2020. Calculating pressure with elastic geobarometry: A comparison of different elastic solutions with application to a calc-silicate gneiss from the Rhodope Metamorphic Province. Lithos 378–379, 105803.

Thomas, J.B., Spear, F.S., 2018. Experimental study of quartz inclusions in garnet at pressures up to 3.0 GPa: evaluating validity of the quartz-in-garnet inclusion elastic thermobarometer. Contributions to Mineralogy and Petrology 173, 42.

How to cite: Moulas, E., Zingerman, K., Vershinin, A., Levin, V., and Podladchikov, Y.: Large strain formulations for host-inclusion systems and their applications to mineral elastic geobarometry, EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-14498,, 2021.

Yury Alkhimenkov, Beatriz Quintal, and Yury Podladchikov

Fluid injection is one of the main triggers of induced seismicity. Accurate numerical modeling of such processes is crucial for the safety of many affected regions. We propose a high-resolution numerical simulation of the strain localization in elasto-plastic and poro-visco-elasto-plastic media with a particular focus on the fluid pressure distribution. The resolution of our numerical model is 10000 by 10000 grid cells. The simulation is accelerated using graphical processing units (GPUs), thus, the total simulation time is in the order of a few minutes. We implement a pressure-dependent Mohr-Coulomb plastic law and study the influence of fluid pressure on the triggering of shear bands. Mean stress is partitioned between fluid pressure and total pressure. This study is particularly important since the effective stress law (the difference between fluid and total pressures) controls brittle failure. We vary viscosity and permeability as well as initial conditions for fluid pressure to explore the physics of shear bands nucleation. We show that fluid pressure in hydro-mechanically coupled media significantly affects the strain localization pattern compared to only elasto-plastic media. Permeability and viscosity are important parameters that control the fluid pressure distribution in the localized shear zones. This work is a preliminary study to model induced seismicity due to the fluid injection in fluid-saturated rocks described as fully coupled poro-visco-elasto-plastic media.

How to cite: Alkhimenkov, Y., Quintal, B., and Podladchikov, Y.: Fluid-total pressure partitioning in shear banding poro-visco-elasto-plastic media, EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-15459,, 2021.

Lyudmila Khakimova, Nikolai Belov, Artyom Myasnikov, Anatoly Vershinin, Kirill Krapivin, Anna Isaeva, Vladimir Dobrozhanskiy, and Yury Podladchikov

This work is devoted to developing the self-consistent thermo-hydro-chemo-mechanical reactive transport model to predict and describe natural and industrial petroleum processes at different scales.

We develop a version of the front tracking approach for multicomponent multiphase flow in order to treat spontaneous splitting of discontinuities. We revisit the solution for the Riemann problem and systematically classify all possible configurations as functions of initial concentrations on both sides of the discontinuity. We validate the algorithm against finite volume high-resolution technics and high-order spectral finite elements.

To calculate the parameters of phase equilibria, we utilize an approach based on the direct minimization of the Gibbs energy of a multicomponent mixture. This method ensures the consistency of the thermodynamic lookup tables. The core of the algorithm is the non-linear free-energy constrained minimization problem, formulated in the form of a linear programming problem by discretization in compositional space.

The impact of the complex rheological response of porous matrix on the morphology of fluid flow and shear deformation localization is considered. Channeling of porosity waves and shear bands morphology and their orientation is investigated for viscoelastoplastic both shear and bulk rheologies.

How to cite: Khakimova, L., Belov, N., Myasnikov, A., Vershinin, A., Krapivin, K., Isaeva, A., Dobrozhanskiy, V., and Podladchikov, Y.: Multicomponent multiphase reactive fluid flow in viscoelastoplastic porous media: localization patterns of fluid flow and strain, EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-15478,, 2021.

Anatoly Vershinin, Vladimir Levin, and Yury Podladchikov

The presentation describes an approach to solving problems of modeling the development of zones of localization of plastic deformations within the framework of a poroelastoplastic model generalizing Biot's model. A distinctive feature of this model is a two-way coupling between mechanical processes occurring in a porous elastoplastic matrix and a saturating viscous fluid.

  For the numerical solution of the problem, a variational formulation based on the Galerkin method and the isoparametric  spectral element method (SEM) is used to discretize the geometric model and PDEs on curvilinear unstructured SEM meshes. SEM orders up to the 15th were used for calculations.

  The software implementation of the developed algorithm based on SEM is performed using CUDA. A spectral element mesh is naturally mapped to a CUDA grid of SMs, and accordingly, each spectral element is mapped to a streaming block, within which individual nodes are processed by the corresponding threads within the block.

The research for this article is performed partially in Schmidt Institute of Physics of the Earth of the Russian Academy of Sciences and supported by the Russian Science Foundation under grant № 19-77-10062.

How to cite: Vershinin, A., Levin, V., and Podladchikov, Y.: Poroelastoplastic modeling of borehole shear bands on high order curvilinear meshes using CUDA technology, EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-15600,, 2021.

Michel Jaboyedoff, Emmanuel Wyser, and Yury Podladchikov

Strain localization plays an important role in the mechanical response of a slumping mass and defines the overall behaviour of such process. We study strain localization with the help of the Material Point Method (MPM), which is well-suited to simulate large deformation problem.

We implemented both mechanical and hydromechanical (i.e., we assume fully saturated conditions of the material) MPM-based solvers within a rate-dependent formulation framework under a GPU architecture. We selected an explicit MPM formulation enriched with the Generalized Interpolation Material Point (GIMP) variant, which fixes a major flaw of MPM, i.e., the cell-crossing error. To avoid spurious oscillation of the pressure field (due to the use of low-order elements) for both solid and liquid phase, we used an element-based averaging technique. This minimizes volumetric locking problems. This numerical framework allows to study high-resolution two-dimensional elasto-plastic (i.e., Mohr-Coulomb plasticity) problems in an affordable amount of time. The solvers were written in a CUDA C environment on a single Nvidia GPU. We report a speed-up factor of 500 compared to a similar MATLAB implementation.

Our results showcase a contribution of pore water pressures over shear banding. In particular, we report a significant influence of the liquid phase over the steady thickness of the shear bands and their location. Pore pressures add a viscous contribution to the elasto-plastic rheological model we choose, i.e., Mohr-Coulomb.

As a future perspective, even high resolution could be achieved considering the extension of the actual implementation toward a multi-GPU solver using MPI.

How to cite: Jaboyedoff, M., Wyser, E., and Podladchikov, Y.: Shear banding in elasto-plastic slumping process on a GPU: mechanical and hydromechanical MPM solver, EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-15861,, 2021.

Taras Gerya, Claudio Petrini, and Luca Dal Zilio

We present a newly developed marker in cell staggered finite difference poro-visco-elasto-plastic numerical model for spontaneous seismic cycle along fluid-bearing fault structures. The fully coupled hydro-mechanical multi-physics model includes poro-elastic compressibility of the solid matrix together with experimentally calibrated rate-dependent strength laws and strain-stress dependent dilation. Localised brittle/plastic deformation is treated accurately through global Picard iterations. To simulate deformation on both long- and short-time scale, an adaptive time stepping is used allowing the resolution of large seismic events with time steps on the order of milliseconds.

Our new numerical modelling tool allows to explore how the presence of pressurised fluids in the pore space of subduction interface and strike-slip zones triggers poro-elastic stress accumulation and release in form of various seismic cycles. The model is capable of simulating spontaneous quasi-periodic seismic events along self-consistently forming highly localized self-pressurised ruptures accommodating shear displacement between the plates. The generated elastic rebound events show slip velocities ranging from the order of Nm/s to m/s, covering the entire range of seismic and slow slip phenomena. The governing strength decrease along the propagating fracture is related mainly to the significant increase of fluid pressure generated by deformation induced plasto-elastic collapse of pores. The reduction of the effective pressure decreases the brittle/plastic strength of fluid-bearing rocks along the rupture, thus providing a dynamic feedback mechanism for the accumulated elastic stress release at the fault interface.  It is remarkable that the seismic behaviours for both slow slip and ordinary earthquakes can be generated within the same self-consistent poro-visco-elasto-plastic rheological framework without any involvement of rate- and state-dependent friction commonly used for seismicity modelling. We furthermore analyse how this process and the seismic cycle are affected by poro-elastic, rate weakening and dilation parameters.

How to cite: Gerya, T., Petrini, C., and Dal Zilio, L.: Roles of poro-elastic compressibility, rate-dependent strength and strain-stress dependent dilation for spontaneous generation of seismicity along fluid-bearing fault structures, EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-16148,, 2021.

Andrey Frendak, Yury Podladchikov, and Ludovic Räss

The ongoing warming of Earth's climate is about to trigger significant change in global sea-level and temperature distribution. Among first evidence is the thawing of hydrate-rich subsurface sediments leading to the potential release of large amounts of greenhouse gases into the oceans and the atmosphere. We hypothesise sea-level and water temperature variations to trigger abrupt changes in the stability of natural systems leading to the spontaneous localisation of flow and deformation. The sedimentary stacks we consider represent fluid saturated deformable porous environments described by interacting thermal, hydrological, chemical and mechanical processes. Resolving the multi-physics interactions or coupling is vital in order to accurately predict the rapid, non-linear and non-trivial evolution of natural systems in fragile equilibrium.

We here investigate how interactions among chemical, hydrological and mechanical processes lead to the spontaneous localisation of flow. We employ a novel numerical modelling framework based on the iterative implicit pseudo-transient method to understand how the relative role of reactions and pore fluid distribution impacts the local deformation in saturated porous media. Resolving strongly coupled multi-physical systems is challenging because accurate results require high resolution calculation in space and time.

Our study aims at better understanding how external forcing parameters such as e.g. pressure and temperature may lead to abrupt changes in the dynamic of complex systems. Ultimately, such investigations should permit further assessment of the longer term evolution and stability of natural systems.

How to cite: Frendak, A., Podladchikov, Y., and Räss, L.: Modelling reactive two-phase flow: challenges and implications, EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-16157,, 2021.

Ivan Utkin, Yury Podladchikov, and Oleg Melnik

One of the mechanisms of magma generation in the Earth's crust is the reaction of dehydration during subduction process. Water is released from subducting lithosphere which leads to the lowering of the melting temperature of mantle rock by hundreds of degrees.

In this work, we present a numerical study of the formation and rise of magma to the Earth's surface, considering partial melting and crystallization of rocks and chemical differentiation of magma. We develop a coupled model of the filtration flow of melt and magmatic fluid through deformable permeable rocks and a thermodynamic model of plagioclase melting based on Gibbs energy minimization approach. The formation of regions with a high melt concentration due to spontaneous focusing of filtration flow being the result of viscoplastic (de)compaction of the pore space is shown. The influence of mechanical properties of rocks and chemical composition of the system on the dynamics of the process is investigated.

How to cite: Utkin, I., Podladchikov, Y., and Melnik, O.: Mathematical modeling of the process of magma formation in the Earth’s crust, EGU General Assembly 2021, online, 19–30 Apr 2021, EGU21-16346,, 2021.