Diamond-inclusion system recording old deep lithosphere conditions at Udachnaya (Siberia)

Diamonds and their inclusions are unique fragments of deep Earth, which provide rare samples from inaccessible portions of our planet. Inclusion-free diamonds cannot provide information on depth of formation, which could be crucial to understand how the carbon cycle operated in the past. Inclusions in diamonds, which remain uncorrupted over geological times, may instead provide direct records of deep Earth’s evolution. Here, we applied elastic geothermobarometry to a diamond-magnesiochromite (mchr) host-inclusion pair from the Udachnaya kimberlite (Siberia, Russia), one of the most important sources of natural diamonds. By combining X-ray diffraction and Fourier-transform infrared spectroscopy data with a new elastic model, we obtained entrapment conditions, Ptrap = 6.5(2) GPa and Ttrap = 1125(32)–1140(33) °C, for the mchr inclusion. These conditions fall on a ca. 35 mW/m2 geotherm and are colder than the great majority of mantle xenoliths from similar depth in the same kimberlite. Our results indicate that cold cratonic conditions persisted for billions of years to at least 200 km in the local lithosphere. The composition of the mchr also indicates that at this depth the lithosphere was, at least locally, ultra-depleted at the time of diamond formation, as opposed to the melt-metasomatized, enriched composition of most xenoliths.

www.nature.com/scientificreports www.nature.com/scientificreports/ The results from the elastic geobarometry applied to diamond-mchr host-inclusion system in this work provide the first direct determination of provenance depth for a mchr-bearing diamond and new constraints on the thermal state and composition of the deep Siberian cratonic lithosphere at the time of diamond formation.
pressure-Volume-temperature equation of State of Reference mchr. The best fit to the P-V data  Table 2) for Ac139 is obtained for a second-order Birch-Murnaghan (BM2) EoS 27 , which gives the following parameters at room conditions: unit-cell volume V 0 = 576. 19(3) Å 3 , Reuss isothermal bulk modulus K RT0 = 182.8(6) GPa (first pressure derivative K RT0 ′ fixed to 4).
The full thermoelastic behavior of a mineral can be expressed through a P-V-T EoS that incorporates the concept of thermal pressure 28 . We therefore estimated the P-V-T EoS of Ac139 by combining the BM2 EoS at room conditions with the T-V data available in the literature for similar compositions. To constrain the high-T thermal behavior of our Ac139, we selected data 29 Figure 2b shows that the normalized unit-cell volumes as a function of T reported in these two studies are in very good agreement, suggesting that relatively large substitutions of Mg for Fe and of Cr for Al have negligible effects on the thermal behavior. Moreover, the literature data 30 clearly show no cation order-disorder process, implying no order-disorder also for the mixed composition 29 . Since the chemical composition of sample Ac139 lies between the two chemical compositions investigated in the literature 29,30 , these two datasets can be considered representative also of its high-T behavior. To further constrain the thermal behavior, we also used the low-T data 31 obtained on a synthetic polycrystalline MgAl 2 O 4 spinel from −196 °C to room-T. The incorporation of the low-T data 31 , even referring to a chemical composition that is different from Ac139, does not affect the final fitting, but provides a further constraint (i.e. lower uncertainties) for the thermal EoS parameters. On a cubic mineral, like mchr, in which the thermal expansion and the axial moduli are equal in all directions, a thermal pressure model based on the quasi-harmonic approximation 28 can be applied. Therefore, in order to determine the P-V-T EoS for mchr, we fitted P-V data for Ac139 and the selected T-V data from the literature by combining our BM2 EoS with the thermal pressure model 28 29,30 , we are confident that our P-V-T EoS parameters can be applied to a wide range of mchr compositions, including the sample here investigated MgCr2.
Relatively to the effect of inversion of the spinel component on the inclusion volume at elevated temperature, it is negligible. In fact, for MgAl 2 O 4 spinel from inversion (x) = 0, corresponding to an ideally ordered phase, to x = 0.29, corresponding to disorder at 1100 °C, volume change is only 0.3%, that is from 529.475 to 528.082 A 3 (ref. 32 ).
The effect of pressure is also negligible, as up to 7.5 GPa pressure has little or no influence on the degree of order in the MgAl 2 O 4 spinel 33 .
Moreover, the influence of inversion on the elastic properties of MgAl 2 O 4 spinel is also very limited. In fact, from x = 0 to x = 0.30, the decrease of bulk modulus K RT0 is about 1%, as shown by comparing results of density functional theory and Brillouin scattering technique 34 Table 4). In Fig. 3 we show one of the FTIR spectra we have collected on MgCr2 diamond host (in Supplementary Data 1 to 6 all FTIR spectra were available). These values are typical for Udachnaya diamonds 35 and allow to classify the diamond as  [29][30][31] . The black line represents the fitting to the P-V and T-V data obtained from a BM2 EoS combined with the thermal-P model 28 Table 4). These two ages cover the range of reported ages for Udachnaya diamonds [37][38][39] . The retrieved T res values differ by less than their respective uncertainties and lie within the typical range of T res for peridotitic diamonds from Udachnaya (ca. 1125-1175 °C, ref. 37 ). Since a short residence at higher T during the early history of the diamond would have not produced significant effects on the final N-aggregation state, the estimated T res of 1125-1140 °C can be considered as a minimum T estimate for diamond formation.
The SC-XRD measurements performed on inclusion MgCr2 before and after its release from the diamond provided unit-cell volumes (V 0 ) of 573.3(2) and 576.63(2) Å 3 , respectively. We used our P-V EoS (Supplementary Table 2) to convert the residual volume strain recorded by the inclusion still trapped in the host diamond into the hydrostatic component of the residual stress, i.e. the residual pressure (P inc ), and obtained a value of 1.073(72) GPa at room conditions. The esd on P inc was calculated by propagation of uncertainties on the EoS parameters (see Supplementary Table 3) and on V 0 measurements. Using this P inc , the P-V-T EoS determined for the compositionally similar mchr Ac139, and the P-V-T EoS for diamond 7 , we can now calculate the entrapment conditions for the MgCr2 inclusion (shape of the inclusion and its distance from the diamond surface did not require further corrections 16 ). Since both the host and the inclusion have cubic symmetry, we can apply the geobarometric model developed for isotropic systems and calculate the entrapment isomeke for the specific P inc using the EosFit-Pinc software 10 . The entrapment isomeke is the line that defines all the possible entrapment conditions for a given P inc (refs 6,8 ). A specific entrapment pressure (P trap ) can be obtained if the T trap is independently constrained. Since the isomeke has a small ∂P/∂T, the T-dependency of the elastic geobarometer is low. Assuming for the MgCr2 inclusion a minimum T trap equal to the estimated T res of 1125-1140 °C, we obtain a minimum P trap of 6.5(2) GPa (the uncertainty derives from the error propagation of all parameters used to calculate it), corresponding to a minimum depth of about 205(6) km ( Fig. 4; the data used to plot the entrapment isomeke in figure are reported in Supplementary Table 5; the input files to run EosFit-Pinc software for our diamond -MgCr2 pair are available in Supplementary Data 7 and 8).

the p-t conditions of the Diamond-mchr System: implications for the Mantle Under
Udachnaya. Thermobarometric data obtained for mantle xenoliths from the Udachnaya kimberlite show relatively large scatter compared to many other cratonic mantle sections, making it difficult to outline a unique geotherm for this locality 40 (Fig. 4). Existing data suggest that the mantle under Udachnaya experienced very cold conditions during at least part of its history, as testified by several peridotites and peridotitic garnets falling close to a ca. 35 mW/m 2 conductive geotherm 40,41 . The lithosphere under Udachnaya extended to ca. 230 km depth at the time of kimberlite eruption, but the depth extension of the cold geotherm is poorly constrained due to scatter of P-T data. In fact, the deepest portion of the lithosphere records significant thermal perturbation, with the majority of high-P xenoliths (>5 GPa) being shifted to ca. 200-300 °C higher T (ref. 40 ) (Fig. 4).
The estimated P trap − T res conditions for diamond -MgCr2 pair (6.5 GPa at ca. 1130 °C; see the P -T data relative to the isomeke in Fig. 4 provided between 700 and 1400 °C in Supplementary Data 9) fall right on the 35 mW/m 2 geotherm (Fig. 4). This indicates that cold conditions have persisted for a long geological time in the Udachnaya lithosphere to at least 200 km depth. Heating processes shortly before or during kimberlite eruption might lead to overestimation of the real T res in the mantle and, consequently, of P trap . However, this effect cannot be significant, otherwise P-T estimates along the entrapment isomeke would become incompatible with even the coldest cratonic geotherm (Fig. 4). We cannot exclude that our diamond has formed during a short thermal pulse at conditions somewhat hotter than the 35 mW/m 2 geotherm, as this would not be recorded in its N-aggregation state 1 . In this case, however, the real P trap would be slightly higher (e.g. ca. 7 GPa for a T trap of 1300 °C) and the  www.nature.com/scientificreports www.nature.com/scientificreports/ depth extension of the long-lived cold geotherm would only be greater (i.e., to ca. 1130 °C and 7 GPa, or 220 km). Much higher P trap , however, would imply long-term residence at conditions, again, too cold for a cratonic lithosphere (Fig. 4). For the same reasons, significant P trap underestimation due to plastic deformation and P inc relaxation can be excluded. Disregarding the inevitable small discrepancies between P-T estimates obtained with different thermobarometric methods, our results provide independent confirmation of literature hypothesis 41 , based on scanty xenolith (megacrystalline peridotites) and garnet xenocryst data, that a cold 35 mW/m 2 geotherm originally extended to ca. 1200 °C.
Abundant experimental data and thermodynamic modeling 2,42-46 showed that the stability of mchr in peridotite is strongly dependent on the peridotite 'fertility' and that in the garnet-spinel stability field Cr# spinel increases rapidly with increasing P and slightly with decreasing T. At given P-T, Cr# spinel depends on the spinel modal abundance and is thus lower in more spinel-rich, more depleted peridotites 2 . Thermodynamic modeling in natural peridotitic systems 2 predicts that, along a 35 mW/m 2 geotherm in very depleted harzburgite with bulk Cr# = 0.32, spinel is only stable to 5.5 GPa and 1000 °C, where it reaches a Cr# spinel of ca. 0.9. These conditions are not compatible with our P-T estimates (Fig. 4). Along hotter geotherms or in less depleted peridotites, the stability of spinel would be restricted to even lower P. This suggests that our mchr belonged to a more depleted, possibly dunitic and spinel-rich assemblage. The persistence of a Cr# ∼ 0.9 (namely, 0.86) in our mchr at 6.5 GPa further supports this hypothesis.
Compared with classical xenolith geothermobarometry, diamond elastic geothermobarometry has a considerable potential to yield ancient lithospheric thermal states, because diamonds can be very old, and their inclusions remain exceptionally preserved from chemical exchanges with surrounding minerals, fluids and melts even under dramatically changing physico-chemical conditions. With our unprecedented application of elastic geothermobarometry to a mchr inclusion in a gem-quality octahedral diamond from the Udachnaya kimberlite, we provided evidence that cold conditions corresponding to a ca. 35 mW/m 2 geotherm persisted for billions of years to a depth of at least 200 km in the Siberian craton lithosphere, and that at this depth the lithosphere was, at least locally, ultra-depleted.

Samples.
We have studied two mchr samples: the first sample was a released crystal used to determine the EoS of mchr, while the second sample was still included within a diamond and was used to determine the depth of formation of the diamond-mchr host-inclusion pair using the thermoelastic properties determined on the first one. The first mchr, labeled Ac139, is a mchr single crystal with a size of 42 × 25 × 15 μm 3 and high crystal quality, i.e. sharp diffraction profiles and absence of twinning and inclusions. It was extracted from a lithospheric diamond from the Aikhal kimberlite (Siberia, Russia) by crushing the diamond and then studied by in-situ high-pressure single-crystal X-ray diffraction. The second mchr, labeled MgCr2, was included in a gemstone octahedral diamond (longest dimension = 2.15 mm) from the Udachnaya kimberlite, Siberia, Russia (Fig. 1), and belongs to the same set of diamonds previously studied 22,47 . The MgCr2 inclusion was alone within its host and located in intermediate position between  61 with the two-pyroxene thermometer 62 . For the three highest-P megacrystalline peridotites (clinopyroxene-free) the orthopyroxene-garnet thermometer 63 was used. Note that this thermometer is less reliable as it is sensitive to redox conditions 64 . Xenolith data are from the compilation of well-equilibrated lherzolites 63 , with additional data on megacrystalline peridotites 40,65 . Conductive geotherms for 35 and 40 mW/ m 2 are after ref. 66 ; the diamond-graphite boundary is after ref. 67 . the gravimeter center and external surface of the diamond. Its morphology is nearly equant cube-octahedral, with an average crystal radius of 0.34 mm, and no optically visible fractures were observed. The mchr inclusion was extracted from its host by crushing the diamond only after a first set of in-situ X-ray analyses. electron microprobe analyses. Chemical analyses were carried out on sample Ac139 using a CAMECA SX50 electron microprobe (WDS mode, 20 kV, 20 nA, 2 μm beam diameter) at CNR-IGG of Padova. Standards used were Kakanui pyrope (New Zealand) from the Smithsonian Museum (TAP analyser crystal, MgKα element emission line); Amelia albite (Virginia) (TAP, NaKα); diopside (TAP, SiKα); Al 2 O 3 (TAP, AlKα); MnTiO 3 (PET, TiKα; LIF, MnKα); Cr 2 O 3 (LIF, CrKα); Fe 2 O 3 (LIF, FeKα); sphalerite (LIF, ZnKα); NiO (LIF, NiKα). Raw data were reduced with the PAP-type correction software provided by CAMECA. The resulting oxides wt% are reported in Supplementary Table 1. Chemical analyses on sample MgCr2 were carried out by standardless EDS analysis on a FEG scanning-electron microscope, using an Oxford Instruments SDD detector and the AZtec TruQ software installed at the Oxford Instruments in Wiessbaden (Germany).
fourier transform infrared spectroscopy. Fourier Transform Infrared (FTIR) were acquired using a Thermo Scientific Nicolet Centaurμs FT-IR Microscope (located at the Department of Industrial Engineering, University of Padova) operated via OMNIC software. The spectra were collected in the 4000-650 cm −1 wavelength range, with 64 scans of sample exposure and 0.5 cm −1 of spectral resolution. The background spectra were run for 120 s before the analyses and the sample spectra were subtracted by them. The measurements were performed on a diamond plate collected from the diamond crushed material after the MgCr2 crystal release. The plate was mounted in wax to allow the analyses in transmission mode. We acquired different FTIR spectra over six different locations on the plate with an FTIR spot of about 100 μm, covering most of its area. The absorbance spectra were processed using DiaMap 48,49 to obtain the N content and aggregation state 50 . This program fits the individual N absorption band and deconvolutes the spectrum in the 1001 to 1350 cm −1 range. Then, the concentrations of N impurities aggregated in pairs (A defects) and of N impurities aggregated in clusters plus vacancy (B defects) were calculated using literature absorption coefficients 51,52 . Errors in calculation of N content and aggregation state are considered of ±10% 48 and depend on the quality of the FTIR spectra. The original FTIR spectra collected in these analyses are reported in the supplementary data (1 to 6), while the data on N content and aggregation state are reported in Supplementary Table 4. A representative FTIR spectrum for the diamond host of the sample MgCr2 (data from such spectrum are relative to Ud_4 in Supplementary Table 4) is shown in Fig. 2.
We used the constants from the literature 53,54 (i.e. E/R = 81160 K; ln(A) = 12.59; where: E = activation energy; R = gas constant; A = Arrhenius constant) to calculate the residence temperature, T res , as a function of time for the MgCr2 diamond. The age of our diamond is unknown, but previous geochronological work on Udachnaya diamonds suggested two major peaks of diamond formation, one at 3.1-3.5 Ga (Re-Os dating on sulphide inclusions 35,36 ) and one at ca. 2 Ga (Sm-Nd dating on garnet inclusions 37 ).
The calculations indicate a mantle residence temperature, T res , of 1125(32) °C for the oldest age (3.5 Ga) and of 1140(33) °C for the youngest age (2 Ga). Uncertainties reported in parentheses (estimated standard deviations, esd) were calculated from the propagation of the uncertainties on the FTIR measurement and on constants taken from literature 53,54 . X-ray diffraction measurements. In-situ high-pressure single-crystal X-ray diffraction (SC-XRD) experiments were performed on mchr Ac139, in order to determine its thermoelastic properties and Pressure − Volume Equation of State. The crystal was loaded in an ETH-type diamond-anvil cell (DAC) equipped with Be-backing plates 55 . We used a steel gasket pre-indented to a thickness of about 90 μm with a spark-eroded hole of 250 μm in diameter and a mixture of methanol:ethanol with ratio 4:1 to ensure hydrostatic pressures up to the maximum P reached in this work (7.852 GPa) (it is well established that for such pressure values the mixture used behaves hydrostatically 56 ). The diameter of the diamond culet was 600 μm. The sample was loaded together with a single crystal of quartz used as an internal pressure standard 57 . The unit-cell edge was measured on a four-circle STOE STADI IV diffractometer equipped with a point detector and controlled by the SINGLE software 58 and installed at the Department of Geosciences, University of Padova. The SINGLE software is based on the 8-position centering method 59 and vector least-squares refinement of the unit-cell parameters 60 . The data were acquired using MoKα radiation at 50 kV and 40 mA. Unit-cell parameters were determined from the centering of an average of 19 reflections for each high-pressure measurement, in the 2θ range between 8° and 34°. The unit-cell edges and unit-cell volumes at different pressures are reported in Supplementary Table 2.
The determination of the residual pressure (P inc ) on inclusion MgCr2 was obtained by comparing the unit-cell parameters measured on the crystal when still trapped in its host and, then, in air after extraction from the diamond. The SC-XRD measurements were performed by using a Rigaku-Oxford Diffraction Supernova diffractometer equipped with a micro-source MoKα radiation, working at 50 kV and 0.8 mA and with a sample-to-detector distance of 68 mm, and a Pilatus 200 K Dectris detector. The diffractometer was installed at the Department of Geosciences, University of Padova. The unit-cell edges before and after the release from the diamond host were measured by collecting data up to 2θ = 80°, using a thin slicing in ω of 0.1° (see values at the section Diamond residence temperature and depth of formation of the diamond-mchr pair).