Self-Consistent Thermodynamic Parameters of Diopside at High Temperatures and High Pressures: Implications for the Adiabatic Geotherm of an Eclogitic Upper Mantle

: Using an iterative numerical approach, we have obtained the self-consistent thermal expansion, heat capacity, and Grüneisen parameters of diopside (MgCaSi 2 O 6 ) over wide pressure and temperature ranges based on experimental data from the literature. Our results agree well with the published experimental and theoretical data. The determined thermodynamic parameters exhibit nonlinear dependences with increasing pressure. Compared with other minerals in the upper mantle, we found that the adiabatic temperature gradient obtained using the thermodynamic data of diopside is larger than that of garnet while lower than that of olivine, when ignoring the Fe incorporation. Combining our results with thermodynamic parameters of garnet obtained in previous studies, we have estimated the adiabatic temperature gradient and geotherm of an eclogitic upper mantle in a depth range of 200–450 km. The results show that the estimated adiabatic temperature gradient of the eclogite model is ~16% and ~3% lower than that of the pyrolite model at a depth of 200 km and 410 km, respectively. However, the high mantle potential temperature of the eclogite model leads to a similar temperature as the pyrolite model in a depth range of 200–410 km.


Introduction
Temperature is one of the fundamental quantities for determining the physical and chemical properties of the Earth's interior. At present, however, it is impossible to measure the temperature of the deep Earth directly. Since the seismic discontinuities in the mantle are usually related to mineral phase transitions [1], the temperatures at the seismic discontinuities (i.e., the 410 km, 520 km, 660 km discontinuities and the mantle-core boundary) can be determined by comparing the depth of the seismic velocity discontinuities with the phase transition pressures of mantle minerals [2]. In addition, since the radiative transfer of heat is regarded to be negligible in the mantle, the temperature gradient of the mantle is expected to be nearly adiabatic [3]. As a consequence, the adiabatic temperature gradient ((∂T/∂z) S ) between two seismic discontinuities can be estimated using thermodynamic parameters at high temperatures and high pressures (P-Ts) via Equation (1) [4,5]: modynamic data on diopside to estimate the adiabatic temperature gradient and geotherm of an eclogitic upper mantle in a depth range of 200-450 km, based on the assemblage of eclogitic mineral compositions, and compare our predictions with previous studies.

Theory
The approach used in this study has been described in published papers [19][20][21]; the fundamental equations are shown below. Note that the uncertainties of the derived parameters are estimated from the experimental measuring error propagations using a Taylor series expansion [20].
The thermal expansion at constant pressure is related to volume (V) as α = (1/V) × (∂V/∂T) P ; the integration of this expression yields [36] V(T) = V T0 exp T T 0 α(T)dT (2) where V T0 is the volume at a reference temperature T 0 . The isothermal derivative of V with respect to pressure can be written as where v B refers to the bulk velocity, which is related to the P-wave velocity (v P ) and Swave velocity (v S ) as Moreover, the isothermal derivative of the heat capacity with respect to pressure can be evaluated by To start the calculation, first, we use the experimental volume data at high temperatures to obtain α(P 0 ,T) at ambient pressure by Equation (2). With α(P 0 ,T) and the relationship between C P and temperature at ambient pressure C P (P 0 ,T), the approximate V at an arbitrary reference pressure can be estimated using Equation (3). Then, the resulting V is used to update the value of α and C P at the same pressure as in Equations (2) and (5), respectively. Hence, iterations of this loop lead to converged values of V, α, and C P as a function of temperature and pressure based on experimental elastic wave velocity. Furthermore, the self-consistent elastic parameters, including the adiabatic bulk modulus (K S ), isothermal bulk modulus (K T ), and shear modulus (G), as well as the Grüneisen parameter (γ), can be consequently derived (Equations (6)-(9)). The elastic wave velocities of diopside have been obtained using Brillouin spectroscopy [37], ultrasonic interferometry [23], and first-principle calculations [24] under high-pressure conditions. The experimental data presented by Li and Neuville [23] indicated that v P and v S of diopside increase linearly with pressure and temperature within their measured P-T range. Meanwhile, the results by Sang and Bass [37] and Zou et al. [24] suggested that both v P and v S show noticeable nonlinear dependence on pressure over~8 GPa. Hence, we use the data from Li and Neuville [23] and fit it with a two-dimensional equation: The fitting coefficients are listed in Table 1, with P in GPa, T in K, and v P and v S in m/s. The adjusted R 2 values are 0.979 and 0.953 for v P and v S , respectively. The calculated elastic wave velocities for Equation (10) agree well with the available data [24,37]. The largest differences are less than 0.5% and 1.2% for v P and v S , respectively, within our calculated P-T range. The ambient pressure unit-cell volumes of diopside were measured using X-ray diffraction at high temperatures [28,29,38,39]. In this study, we fitted the experimental volume data from Pandolfo et al. [28] to Equation (2) with the program EosFit7 GUI [40], thereby obtaining the thermal expansion as a function of temperature using the formula suggested by Fei [41]: Here, we take T 0 = 303 K and V T0 = 428.70 Å 3 , which are the values suggested by Pandolfo et al. [28]. The volumes calculated using Equations (2) and (11) show good agreement with the previous experimental results [28,29,38,39], with the largest difference being less than 0.3%.

Heat Capacity under High-Temperature and Ambient Pressure Conditions
The heat capacity of diopside at high temperatures below 1000 K has been determined by Stebbins et al. [32] and Krupka et al. [30] by using calorimetry. Then, Richet and Fiquet [31] extended the measuring temperature to~1600 K, and the experimental data were fitted to different forms. Here, we use Equation (12), which was suggested by Richet and Fiquet [31] to be the best for describing their data: C P = −186.908 + 57.729 ln(T) + 5.5483 × 10 4 T −1 − 2.2557 × 10 7 T −2 + 2.4436 × 10 9 T −3 (12) where C P is in J/mol K and T is in K. The quality of the fit is excellent, with an adjusted R 2 value of 0.993 and the largest error being 0.12%. Figure 1 presents comparisons between our calculated unit-cell volume of diopside (Supplementary Materials) and previously published results. At 300 K, our data shows good agreement with the results measured using X-ray diffraction [42,43]. The excellent consistency with available experimental results also holds at high temperatures up to 1073 K [23,25], with the largest differences being less than~0.3%. Though our unit-cell volume is slightly lower than the first-principle calculation result under high-temperature conditions [24], the uncertainty is within 1.2%.  [24]. The squares and triangles represent the measurements at 300 K by Tribaudino et al. (2000) [43] and Thompson and Downs (2008) [42], respectively. The diamonds and circles correspond to the measurements at various temperatures by Zhao et al. (1998) [25] and Li and Neuville (2010) [23], respectively.

Elastic Properties at High Temperatures and High Pressures
With the combination of the determined unit-cell volume and previously published v P and v S values at P-Ts [23], the K T , K S , and G of diopside can be calculated via Equations (6)- (8); the results are shown in Table 2    As shown in Table 2, the K T , K S , and G values under ambient conditions have been determined as 114.8 GPa, 116 GPa, and 72.8 GPa, respectively, which generally agree with the experimental [23,25,37,[42][43][44][45] and theoretical [24,46] results. All three elastic moduli show a nonlinear relationship with pressure and a linear relationship with temperature. At room temperature, the first pressure derivatives of the elastic moduli have been determined as ∂K T /∂P = 5.17, ∂K S /∂P = 5.18, and ∂G/∂P = 1.67, and the second pressure derivatives of elastic moduli have been determined as ∂ 2 K T /∂P 2 = −0.066/GPa, ∂ 2 K S /∂P 2 = −0.067/GPa, and ∂ 2 G/∂P 2 = −0.0106/GPa. Also, the temperature derivatives of elastic moduli have been determined as ∂K T /∂T = −0.0142 GPa/K, ∂K S /∂T = −0.011 GPa/K, and ∂G/∂T = −0.0114 GPa/K under ambient pressure conditions; these values are consistent with previously presented results [23][24][25]. Table 2. Elastic moduli of diopside and their first and second derivatives with respect to pressure and temperature. The pressure and temperature dependences of these elastic moduli are shown in Figures 2 and 3, respectively. In Figure 2b,c, the calculated K S and G values are in good agreement with the available experimental measurements under the relevant conditions [23,37]. Also, in Figure 2b, our ambient and high temperature K S values generally agree with the first-principle calculations performed by Zou et al. [24] for under 10 GPa. Then, the pressure derivative becomes smaller than that from Zou et al. [24], leading to a significant separation of~10% at 20 GPa. Then, in Figure 2c, our calculated G value is in good agreement with the result from Zou et al. [24] at 1100 K. Compared to Zou et al. [24], our result is larger than that from Zou et al. [24] at 300 K, and smaller than that from Zou et al. [24] at 2000 K.
In Figure 3, at ambient pressure, our elastic moduli are generally consistent with the experimental study by Isaak et al. [33], which presented the ambient pressure elastic properties of chrome-diopside. The largest differences between our and Isaak's results are 2.1%, 0.6%, and 2.8% for K T , K S , and G, respectively. In addition, our K S and G values are generally consistent with the results from Zou et al. [24] below 10 GPa. Then, the separations get more significant along with the pressure due to the differences in the pressure derivatives of K S and G.

Thermodynamic Properties at High Temperatures and High Pressures
Our calculated thermal expansion, heat capacity, and Grüneisen parameters of diopside are illustrated in Figures 4-6, respectively. Figure 4a shows the temperature dependence of the thermal expansion of diopside at several pressure conditions, as well as previous results for comparisons. At ambient pressure, the temperature derivative of our calculated thermal expansion ∂α/∂T is~5.0 × 10 -4 /K 2 over 500 K. This value is lower than that presented by Isaak et al. [33], which was obtained from the combination of previous studies [38,39,47]. Moreover, the difference in ∂α/∂T causes a 17% deviation at 1300 K. On the other hand, our result is in good agreement with the recent data obtained via an X-ray diffraction study by Hovis et al. [29], with the largest difference being 4.7% in the temperature range of 500-1500 K.   Currently, the published heat capacities of diopside are limited under ambient conditions. In Figure 5a, since the heat capacity of diopside for our calculation was taken from previous experimental studies [31], our obtained heat capacity is in good agreement with the experimental data presented by Krupka et al. [30], with a difference of 1.9%, and generally agrees with the result presented by Stebbins et al. [32] below~850 K, with a difference of 1.5%.
The Grüneisen parameter is a valuable thermodynamic parameter in geophysics which can set limitations on the pressure and temperature dependence of the thermal properties of the mantle and core and constrain the adiabatic temperature gradient [48]. In Figure 6a, our calculated Grüneisen parameter shows a rapid decrease to~500 K and then tends to nearly a constant, with a very small temperature derivative of -5.0 × 10 -5 /K at 2000 K. The result is 15% lower than that from Isaak et al. [33] at 1300 K, which is probably caused by the difference in the chemical composition of minerals. Also, the Grüneisen parameter proposed by Isaak et al. [33] demonstrated a rapid increase below~500 K; in contrast, our determined data shows an opposite trend, exhibiting a similar pattern to garnet and olivine under ambient pressure conditions [20,49].
Additionally, to investigate the effects of pressure on the thermodynamic properties of diopside, the calculated thermal expansion, heat capacity, and Grüneisen parameters as a function of pressure at different temperatures are illustrated in Figures 4b, 5b and 6b, respectively. All α, C P , and γ values for diopside show nonlinear pressure dependences. Hence, we fitted the thermodynamic parameters into the equation N = N 0 + ∂N/∂P × P + ∂ 2 N/∂P 2 × P 2 at different temperatures; where N refers to the thermodynamic properties; N 0 refers to N ambient pressure; and ∂N/∂P and ∂ 2 N/∂P 2 refer to the first and second pressure derivatives of N, respectively. The fitting coefficients are listed in Table 3. In Table 3, all α, C P , and γ values show negative pressure dependences, and the effects of pressure on α and C P increase along with the temperature. At room temperature, the first pressure derivative of α is ∂α/∂P = -0.92 × 10 -6 /K GPa. This value is larger than that suggested by Zhao et al. [25] (-1.52(2) × 10 -6 /K GPa), which was determined in the P-T range of 8.2 GPa and 1280 K with a linear fitting equation. Then, increasing the temperature to 2000 K, ∂α/∂P decreases to ∂α/∂P = -1.37 × 10 -6 /K GPa. Similarly, the pressure derivative of C P (∂C P /∂P) decreases from -0.417 J/mol K GPa to -0.688 J/mol K GPa in the temperature range of 300-2000 K. In particular, the first pressure derivative of the Grüneisen parameter is ∂γ/∂P = −3.4 × 10 -3 /GPa at room temperature. Then, the value increases with temperature to −2.6 × 10 -3 /GPa at~500 K because of the variations at the low-temperature range, and decreases with the temperature again to -4.8 × 10 -3 /GPa at 2000 K.

Effect of Mineral Composition on the Adiabatic Temperature Gradient
Based on the obtained thermodynamic data on diopside, we have estimated the adiabatic temperature gradient ((∂T/∂z) S ) by using Equation (1); the result is shown in Figure 7. To investigate the effect of mineral composition on the adiabatic temperature gradient, the (∂T/∂z) S values obtained using thermodynamic data on Fe-bearing olivine  Figure 7 for comparison [19][20][21]. Note that the mantle potential temperature (T P ), which refers to the temperature that it would have upon ascending adiabatically to the Earth's surface without undergoing partial melting [50], is employed as T P = 1610 K [18]. Additionally, the pressure and gravitational acceleration variations with depth are from the study performed by Stacey and Davis [4]. The (∂T/∂z) S values obtained using thermodynamic data on diopside are 0.46 K/km, 0.39 K/km, and 0.32 K/km at depths of 0 km, 200 km, and 410 km, respectively; the value is~14% lower than that obtained by Stacey and Davis [4] in the depth range of 220-410 km. For the Fe-free minerals, at a depth of 200 km, the (∂T/∂z) S value obtained using thermodynamic data on diopside is 14% and 5% higher than that of pyrope and grossular, respectively, whilst it is 17% lower than that of forsterite. Also, this value is 33% and 23% lower than that of the Fe-bearing olivine and Fe-bearing garnet (almandine), respectively. Moreover, a previous study suggested that the depth derivatives of (∂T/∂z) S of those minerals from the same group are similar [20]. In this study, the depth derivative of (∂T/∂z) S obtained from thermodynamic data on diopside is -3.0 × 10 -4 K/km 2 in the depth range of 200-410 km, which is smoother than that of olivine (-5.0 × 10 -4 K/km 2 ) and steeper than that of garnet (-2.0 × 10 -4 K/km 2 ) [20].

Adiabatic Geotherm of the Eclogite Model
It is widely accepted that the bulk mantle mineral assemblages can be represented by the pyrolite model [11,51]. The seismic discontinuities at depths of 410 km, 520 km, and 660 km are associated with the phase transitions from (Mg,Fe) 2 SiO 4 olivine (α-phase) to wadsleyite (β-phase) and from wadsleyite to ringwoodite (γ-phase) and the dissociation of ringwoodite into bridgmanite and periclase, respectively, for the pyrolite model [52,53]. Hence, the temperature at these seismic discontinuities can be estimated by comparing the depth of these seismic discontinuities with the transition pressures in the system Mg 2 SiO 4 -Fe 2 SiO 4 [2,18]. With the combination of the thermal equations of states of (Mg, Fe) 2 SiO 4 olivine and its high-pressure phase minerals, as well as Equation (1), the adiabatic temperature gradient (∂T/∂z) S and adiabatic geotherm of a pyrolitic mantle was presented for~3000 km by Katsura et al. [18], which has been widely used in later investigations.
Meanwhile, the upper mantle under the subducted oceanic crust is generally considered eclogitic in composition [12,54]. The subducting slab has a layered structure: overlying crust of mid-ocean ridge basalt (MORB), followed by residual harzburgite and depleted pyrolite from top to bottom [55,56]. When the subducting slab is subducted into the mantle, the MORB transforms directly into eclogite, which mainly contains an assemblage of eclogitic garnet and clinopyroxene, and becomes denser than the harzburgite underlying it [57,58]. Then, the clinopyroxene starts to dissolve progressively into stishovite and garnet, which then completely forms majoritic garnet by a depth of~450 km [59][60][61]. Hence, the mineral composition of an eclogitic upper mantle is quite different from a pyrolitic upper mantle, which raises the question: Will the geotherm of the eclogite model also be different from the pyrolite model?
In order to estimate the adiabatic geotherm of an eclogitic upper mantle with Equation (1), it is necessary to determine the T P of the oceanic mantle first. Previous research on the generation of the oceanic upper mantle proposed the T P of mid-ocean ridges as 1533 K [62]; then, the temperature was suggested to be hotter [63,64]. Recently, with the consideration of water content of the peridotite in MORB, the T P of mid-ocean ridges was constrained to 1623-1683 K [65]. Also, the composition at the top of an eclogitic upper mantle (200 km) is 70% clinopyroxene, and~30% eclogitic garnet [66]. Therefore, with the high temperature and high pressure thermodynamic data on diopside and eclogitic garnet, the temperature at a depth of 200 km can be estimated as 1699-1762 K via Equation (13) [50]. It is worth noting that the eclogitic garnet is composed of approximately 40% pyrope, 30% almandine, and 30% grossular [67,68].
Based on the variations of mineral volume fraction with depth presented by Irifune et al. [66], our calculated (∂T/∂z) S value of the eclogite model in the depth range of 200-450 km is shown in Figure 8, together with the result obtained by using the seismic velocity [4] and the result for a pyrolitic upper mantle [18]. The (∂T/∂z) S value of the upper mantle was presented as a constant of 0.3 K/km in an early study [5]; later, a higher value was suggested at the top part of the upper mantle (0.4-0.5 K/km) [4,18]. In Figure 8, the (∂T/∂z) S value of the eclogitic model decreases from 0.39 K/km to 0.34 K/km at a depth range of 200-410 km when T P = 1623 K, which is~12% lower than that presented by Stacey and Davis [4]. Also, the (∂T/∂z) S value of the eclogite model is~17% lower than that of the pyrolite model suggested by Katsura et al. [18] at 200 km, but increasing the depth to 410 km, the values become closer, with a difference of 4%. Under 410 km, the (∂T/∂z) S value of the eclogite model maintains the same trend as the depth and decreases to 0.33 K/km at 450 km, while the (∂T/∂z) S value of the pyrolite model shows an increase because of the olivine-wadsleyite transition. Moreover, the value of (∂T/∂z) S for the eclogite model shows a~5% increase with T P increasing from 1623 K to 1683 K; the values are 0.41 K/km, 0.35 K/km, and 0.34 K/km at depths of 200 km, 410 km, and 450 km, respectively.
Our obtained adiabatic geotherm of an eclogitic mantle from 200 km to 450 km is shown in Figure 9. At depths of 200 km, 410 km, and 450 km, the temperatures of the eclogite model are 1699 K, 1775 K, and 1789 K when T P = 1623 K, respectively, and 1762 K, 1841 K, and 1855 K when T P = 1683 K, respectively. Under the condition of T P = 1683 K, the estimated temperature is generally consistent with the results proposed by Stacey and Davis [4] and Trubitsyn and Trubitsyn [69] in the depth range of~200-410 km, with the differences being less than 0.8% and 1%, respectively, but~80 K higher than the results obtained by Anderson [8] and Brown and Shankland [9] at a depth of 410 km. By comparing the depth of the 410 km discontinuity with the forsterite-wadsleyite phase transition pressure, the temperature at a depth of 409 km was determined to be 1830 ± 48 K for the pyrolite model, and consequently, the T P was presented as 1610 ± 35 K [18]. Because of the existence of olivine in the pyrolitic mantle, though the T P of the pyrolite model is lower than that of the eclogite model, the higher (∂T/∂z) S value of the pyrolite model leads to similar temperatures to the eclogite model in the depth range of 200-410 km. However, the temperature of the eclogite model does not show the increase demonstrated in the pyrolite model at the 410 km discontinuity, but maintains the same trend at a depth of 450 km.  It is worth noting that, except for olivine, the peridotitic garnet, clinopyroxene, and orthopyroxene are essential compositions in a pyrolite mantle [67,70]. Additionally, the majority of the clinopyroxene in an eclogite mantle is considered to be omphacite, which is mainly composed of diopside and jadeite (NaAlSi 2 O 6 ) [71]. Also, the majoritic garnet becomes dominant over a depth of~450 km [59][60][61]. Furthermore, Fe is a non-negligible content in both olivine and clinopyroxene [72,73]. Previous studies have suggested that Fe incorporation could increase the estimated adiabatic temperature gradient [20]. Hence, the differences between the geotherm of the pyrolite and eclogite models are probably not the same as what we described if we consider more complicated mineral compositions for the pyrolite model. More thermodynamic parameters of mantle minerals are still required to improve the geotherm model. Overall, this study shows that the adiabatic temperature gradient and adiabatic geotherm of different mantle compositions exhibit various trends with depth. As a consequence, suitable geotherms need to be taken into account when investigating different tectonic settings.

Conclusions
In summary, following a numerical iterative procedure, we obtained self-consistent thermodynamic parameters: thermal expansion, heat capacity, and Grüneisen parameters of diopside in a wide P-T range (20 GPa, 2000 K). All these thermodynamic parameters are nonlinearly and negatively correlated with pressure. Also, the pressure effects on thermal expansion and heat capacity increase with temperature, while the effect of pressure on the Grüneisen parameter decreases to a minimum at~500 K, then increases with temperature again. Based on thermodynamic data on major mantle minerals from this and previous studies, we have estimated the adiabatic temperature gradient and geotherm of an eclogitic upper mantle. The results suggest that the adiabatic temperature gradient of the eclogite model is lower than that of the pyrolite model at the same depth, which causes a slower rise of temperature than that of a pyrolitic upper mantle and leads to a similar adiabatic geotherm for both pyrolitic and eclogitic upper mantle in the depth range of 200-410 km.