Metamorphic rocks are the records of plate tectonic processes whose reconstruction relies on correct estimates of the pressures and temperatures (P-T) experienced by these rocks through time. Unlike chemical geothermobarometry, elastic geobarometry does not rely on chemical equilibrium between minerals, so it has the potential to provide information on overstepping of reaction boundaries and to identify other examples of non-equilibrium behavior in rocks. Here we introduce a method that exploits the anisotropy in elastic properties of minerals to determine the unique P and T of entrapment from a single inclusion in a mineral host. We apply it to preserved quartz inclusions in garnet from eclogite xenoliths hosted in Yakutian kimberlites (Russia). Our results demonstrate that quartz trapped in garnet can be preserved when the rock reaches the stability field of coesite (the high-pressure and high-temperature polymorph of quartz) at 3 GPa and 850 °C. This supports a metamorphic origin for these xenoliths and sheds light on the mechanisms of craton accretion from a subducted crustal protolith. Furthermore, we show that interpreting P and T conditions reached by a rock from the simple phase identification of key inclusion minerals can be misleading.

The mechanisms attending the downward transport of crustal material into the mantle and its return back to Earth’s surface (exhumation) are still a matter of vigorous debate. Chemical information only allows the interpretation of the measurements on mineral and rock composition in terms of pressure (P) for perfectly lithostatic systems under ideal chemical equilibrium. However, significant overstepping of reaction boundaries (Spear and Pattison, 2017) as well as the presence of non-lithostatic stresses might prevent the correct interpretation of P, and in turn depths, reached by crustal rocks during subduction and metamorphism. Host-inclusion geobarometry provides an alternative and complementary method to determine pressures and temperatures (P-T) attained during the history of rocks (Zhang, 1998; Angel et al., 2014b, 2015). A mineral trapped as an inclusion within another host mineral is not free to expand or contract as would a free crystal but is constrained by the host mineral. This results in the development of stress in the inclusion that differs from the external stress or pressure applied to the host mineral, both while it is in the earth and afterwards when we examine the rock at room pressure in the laboratory. The stress state of the inclusion arises from the change in P and T from the time of its entrapment, so measurement of the stress state of the still-entrapped inclusion while in the laboratory enables the conditions of entrapment to be calculated, provided no plastic or brittle deformation occurred upon exhumation after entrapment. However, the current state of the art is based upon the assumption that both the host and the inclusion are elastically isotropic. But, no mineral is isotropic in elastic properties and this may cause errors in calculated P and T. This also means that a measurement of a single inclusion pressure while the host is at room pressure provides only one constraint on the entrapment conditions. As a consequence, one can only calculate a line in P-T space (the entrapment isomeke) which represents possible entrapment conditions of the inclusion (Rosenfeld and Chase, 1961; Angel et al., 2014b). Here we describe how the anisotropy of mineral inclusions can be exploited to determine unique P-T conditions last recorded by the rock. The basic idea behind this approach is that an anisotropic inclusion will exhibit different stresses and strain along different crystallographic directions. By measuring these strains from a single inclusion, we obtain two or three independent data which, in combination with the known P-T variation of the unit-cell parameters of the inclusion mineral and the host, enable both the P and T of entrapment or elastic equilibration of the inclusion to be determined.

The Mir pipe (Yakutian kimberlites, Russia) is a relatively young kimberlite (360 Ma) for which eruption temperatures of ∼1000 °C and very fast ascent rates or short residence time (<0.1 m.y.) have been estimated (Korsakov et al., 2009; Zhukov and Korsakov, 2015). The Mir pipe kimberlite carried to the surface eclogite xenoliths made of omphacite, garnet, and rutile (Taylor et al., 2003; Tomilenko et al., 2005; Zhukov and Korsakov, 2015). The major and trace element bulk composition of these eclogite xenoliths suggests that they may be derived from subducted oceanic crust (Shimizu and Sobolev, 1995; Taylor et al., 2003). The studied xenolith contains coarse homogeneous garnet hosting relatively large primary quartz inclusions (Fig. 1). Previous models show that the short eruption time prevents any significant resetting of inclusion pressures by plastic flow of the garnet during upward transport from the mantle (Zhong et al., 2018). This is confirmed by the estimated residual pressures of 1.0–1.2 GPa for the quartz inclusions, which are significantly higher than 0.1–0.6 GPa reported on other quartz inclusions in garnets from coesite-grade and diamond-grade ultrahigh-pressure (UHP) rocks (see Korsakov et al., 2009). The relatively large sizes and high pressures of the inclusions and the homogeneous composition of the garnet host, coupled with the fast kimberlite ascent which may have been insufficient to reset the rock-forming minerals, make these eclogites an ideal case to test and prove the potential of anisotropic elastic geobarometry.

The birefringence haloes around the four selected inclusions (see Fig. 1) indicate the presence of significant residual stresses in the inclusions (Korsakov et al., 2007, 2009; Howell et al., 2010; Campomenosi et al., 2018;). The upshift of the Raman bands (Fig. 2A; for further details, see the GSA Data Repository1) at various positions across the inclusion confirms that the inclusion is under significant residual strain. Quantitative values of the strains have been determined from the wavenumber shifts of the Raman bands at 128, 206, 464, and 696 cm−1 by using the mode Grüneisen tensors of quartz (Murri et al., 2018, 2019; Angel et al., 2019; Bonazzi et al., 2019). The strains determined in this way from several inclusions in the same garnet are identical within estimated standard deviations, and agree with the measurements by single-crystal X-ray diffraction (Table 1).

The residual strains as determined by X-ray diffraction and by micro-Raman spectroscopy (MRS) cannot be directly used to back-calculate the residual pressure at entrapment conditions because they are the product of two processes: the contrast in the elastic properties of the host and inclusion over P and T that leads to the inclusions exhibiting an excess pressure, and the mutual elastic relaxation of the system driven by this excess pressure (Angel et al., 2014b). Before calculating entrapment conditions by using the equations of state of the host and inclusion minerals, the residual strains must be corrected for elastic relaxation (Angel et al., 2017b).

Correction for elastic relaxation of a spherical inclusion in an elastically isotropic system depends only on the elastic properties of the host and inclusion (Zhang, 1998; Angel et al., 2017b). For faceted or complex-shaped inclusions, and for all elastically anisotropic inclusions, the final stress state depends on both the geometry of the system and the elastic anisotropy of the inclusion-host pair (Eshelby, 1957; Zhukov and Korsakov, 2015; Campomenosi et al., 2018; Mazzucchelli et al., 2018, 2019). Therefore, for all of these cases, the elastic relaxation must be calculated numerically (Mazzucchelli et al., 2018, 2019) using the exact shape of the inclusion and its full elastic properties together with those for the host. We determined the three-dimensional (3-D) model of the entire sample (see Fig. 1E) from X-ray microtomography measurements at the TOMCAT (Tomographic Microscopy and Coherent Radiology Experiments) beamline of the Swiss Light Source (SLS; Paul Scherrer Institut) (Stampanoni et al., 2006). From these images, we created a meshed 3-D model in order to perform finite element (FE) analyses that allow us to calculate the amount of elastic relaxation that has occurred (Mazzucchelli et al., 2018). Because of the anisotropy of quartz, the amount of elastic relaxation is different in different directions. Only after correction for the elastic relaxation can we demonstrate that the measured quartz inclusions were subject to isotropic strain (Table 1; Fig. 2B) as a consequence of the change of the dimensions of the cubic host mineral with P and T. The corrected strains of all four inclusions are identical within the estimated uncertainties, consistent with them having been trapped under the same conditions and having experienced the same post-entrapment history.

Once the measured inclusion pressure or strain (as appropriate) has been corrected for the effects of mutual relaxation, the determination of possible entrapment conditions is a purely thermodynamic calculation (see the Data Repository). In the isotropic model (Angel et al., 2014b, 2015, 2017b), entrapment conditions are calculated as the conditions under which there are no stress gradients across the host and inclusion, which leads to an entrapment isomeke, which is a line in P-T space along which there would have been no strain gradients in the system. We exploit the anisotropy of quartz and the known variation with P and T of its unit-cell parameters (for details see the Data Repository) to calculate a line of possible entrapment conditions from each of the corrected strains calculated along a and c crystallographic axes (Table 1; see details in the Data Repository). The intersection of these two lines (e.g., crossing point of dashed and solid lines in Fig. 2C) provides a unique P and T of entrapment and /or equilibration for each inclusion, three of which are very similar and cluster around 3 GPa, close to the mantle geotherms corresponding to surface heat flow of 40–50 mW/m2 considered appropriate for these eclogites (Litasov et al., 2003).

We show that the quartz inclusions in the garnet of the eclogite xenoliths from the Mir kimberlite record P-T conditions of elastic equilibration within the coesite stability field at pressures lower than the stability field of diamond (Fig. 2C). We can exclude the possibility that the inclusions were originally trapped as coesite and then inverted to quartz because the inclusions are perfect single crystals (see X-ray diffraction data in the Data Repository) and do not exhibit the palisade texture typical of quartz inclusions inverted from coesite (Schertl et al., 1991; Kotková et al., 2011). Neither is there evidence (Fig. 1) of cracks that would have allowed access of fluids to catalyze the otherwise very slow inversion from coesite. Fluid infiltration would also be required to remove silica from the inclusion to allow the pressure to drop below the pressure of the coesite-quartz equilibrium phase boundary, at which the pressure of an isolated inclusion consisting of a mixture of coesite and quartz should be buffered (Ferrero and Angel, 2018). The absence of coesite from this fragment of xenolith, which consists mostly of garnet and clinopyroxene instead, indicates that the garnet grew at relatively low P-T along the prograde subduction path of a crustal protolith, which produced minor quartz that became entrapped as inclusions. Subduction then proceeded toward UHP conditions to mantle depths of ∼100 km and the pressures in the quartz inclusions, much softer than the garnet host crystals, will lag behind the external pressures by ∼1.5 GPa (Angel et al., 2015). Therefore, at rock pressures of 3.5 GPa, the quartz experiences only 2.0 GPa, well within its stability field, and would not transform to coesite. A significant period of residence of the studied garnet at high temperatures would be consistent with the complete absence of compositional zoning and gradients in the eclogitic garnet host (Shimizu and Sobolev, 1995). We suggest that at the same time that the compositional gradients in the garnet were eliminated by diffusion, the pressure difference between the quartz inclusions and the garnet host drove plastic flow that relaxed the stress difference while maintaining the isolation of the inclusions so that no fluids were available to catalyze the quartz to coesite transformation. Similar resetting of inclusion pressures in garnets has also been documented in metamorphic settings (Ferrero and Angel, 2018). Subsequent eruption of the kimberlite would have been sufficiently rapid to prevent further significant resetting of the inclusion stress state (Zhong et al., 2018). Therefore, the remanent pressures that we measured on these inclusions reflect resetting at mantle pressures and temperatures following subduction (see Figs. 2 and 3). We note that the resetting temperature and pressure cluster around a mantle geotherm associated with surface heat flow of 45 mW/m2, appropriate for ambient conditions at the time of the kimberlite eruption, so they may indicate that the xenolith was incorporated at the time of eruption (Litasov et al., 2003). On the other hand, the measured remanent inclusion pressures definitely exclude the possibility that the studied mantle eclogite xenoliths are the product of direct high-pressure magma crystallization at mantle depths, as already supported by geochemical evidence for other xenoliths (Griffin and O’Reilly, 2007; Aulbach et al., 2017). Direct crystallization at depth would result in either coesite inclusions, or inclusions whose remanent pressures would indicate entrapment within the stability field of quartz.

This example illustrates that the anisotropic behavior of quartz is not just a complication that may invalidate an isotropic analysis of inclusion data, but actually provides additional information that can be used in certain circumstances to provide P and T estimates that are completely independent of chemical equilibrium in the rock. The same principles and methodology can certainly be applied to any other uniaxial mineral trapped in almost isotropic host minerals such as garnets.

This project received funding from the European Research Council under the European Union’s Horizon 2020 research and innovation program grant agreement 714936 (ERC-STG TRUE DEPTHS) to M. Alvaro. Alvaro and Scambelluri have been supported by the Ministero dell’Istruzione dell’Università e della Ricerca (MIUR)Progetti di Ricerca di Interesse Nazionale (PRIN)Bando PRIN 2017 - Prot. 2017ZE49E7_005. We acknowledge the Paul Scherrer Institut, Villigen, Switzerland, for provision of synchrotron radiation beamtime at the TOMCAT beamline of the Swiss Light Source. We thank Catherine Mottram, Jay Thomas, and one anonymous reviewer.

1GSA Data Repository item 2020006, theoretical backgound for the calculations, is available online at, or on request from
Gold Open Access: This paper is published under the terms of the CC-BY license.