Effect of Thermoelastic Properties of the Pyrope-Almandine Solid Solutions on the Entrapment Pressure of Garnet-Related Elastic Geobarometer

The pyrope (Prp)–almandine (Alm) solid solutions are the most fundamental garnet components on the Earth, and both the quartz inclusions in garnet (QuiG) barometry and the garnet inclusions in diamond barometry need to be constrained by the thermoelastic parameters of Prp-Alm solid solution garnets. Here, we report the thermoelastic properties of a series of synthetic Prp-Alm solid solutions based on the high-pressure and high-temperature (HP–HT) in situ synchrotron single-crystal x-ray diffraction (SCXRD) experiments up to ∼20 GPa and 700 K, using diamond anvil cell (DAC). Fitting the SCXRD data by the Birch-Murnaghan equation of state (BM-EoS) and the thermal-pressure EoS, we obtain the thermoelastic parameters of Prp-Alm solid solution garnets, including bulk modulus (K 0), its pressure derivative (K′0), and the thermal expansion coefficient ( α 0 ). The K 0 along the Prp-Alm solid solution changes linearly with Prp content within their uncertainties and can be expressed by K 0 (GPa) = 181.0(8) – 0.11(1) X prp (R 2 = 0.91, X prp is the Prp mole fraction and K′0 fixed at 4). Our result indicates that the compressibility of the Prp-Alm solid solution increases with the increasing Prp content. However, the thermal expansion coefficient of Prp-Alm solid solution at ambient pressure shows a non-linear trend with Prp content and can be expressed by α 0 (10−5 K−1) = 2.7 (1) + 3.0 (5) X Prp − 3.2 (4) X 2 Prp (R 2 = 0.985). It shows that the Prp-Alm solid solution with intermediate composition has a larger thermal expansion coefficient than those close to the endmembers at ambient conditions. Furthermore, we also evaluated the influence of thermoelastic properties of the Prp-Alm solid solution on the entrapment pressure (P e ) estimation for two types of elastic geobarometers. Our results indicate that the garnet component may significantly influence entrapment pressure, and among the thermoelastic parameters of garnet, the thermal expansion coefficient has the main effect on the estimation of P e .


INTRODUCTION
Diamond is the only direct sample obtained from inaccessible portions of Earth  and the subduction zone metamorphism controls many global geochemical cycles and the lithosphere (Ashley et al., 2016). Estimating the pressure and temperature of these minerals or rocks provides fundamental data for studying many such geological processes, but their determination remains extremely controversial (Bebout, 2007;Tajčmanová et al., 2021). To date, various methods have been developed to address the geological history of rocks. Previous studies have used mineral-mineral equilibrium based on Gibbs free energy minimization or partitioning major or trace elements between phases to recover these geological processes (De Capitani and Petrakakis, 2010;Holland and Powell, 2011). Despite advances in analytical techniques, geothermobarometers still suffer from problems due to alteration processes, such as the erasure of peak metamorphic mineral assemblages (e.g., Korsakov et al., 2009;Jamtveit et al., 2016), chemical re-equilibration, diffusion, and kinetic limitations (e.g., Anzolini et al., 2019;Gonzalez et al., 2019). Additionally, part of the difficulty in determining the pressure of rocks is the lack of typical pressure-dependent mineral equilibria, especially for the diamond-inclusion system (Ashley et al., 2016;Nestola et al., 2017).
Ideally, geobarometers should yield accurate pressure estimates using only commonly found minerals and not be susceptible to significant resetting during exhumation. Recently, an alternative method based on the mechanical equilibrium between entrapped mineral inclusions (e.g., quartz and garnet) and host mineral (e.g., garnet and diamond) has received significant attention and could vastly expand the range of barometers available to petrologists (e.g., Enami et al., 2007;Kohn, 2014;Milani et al., 2015;Cisneros et al., 2020). Inclusions encapsulated in host minerals such as garnet or diamond have different elastic properties; they produce a residual pressure (P inc ) following exhumation, thereby allowing us to determine the entrapment pressure (Rosenfeld and Chase, 1961;Zhang, 1998). The P inc and the equations of state (EoS) of the hostinclusion system can be used to calculate unique P-T curves (called isomekes) along which the pressure of the host and the inclusion are equal despite the changes in P-T conditions. Detailed discussions of this application are given by Angel et al. (2015) and Gonzalez et al. (2019).
Much of the previous work has focused on the factors affecting the calculation accuracy of the elastic geobarometer, such as the viscosity/plasticity relaxation properties of the inclusions (Zhong et al., 2020), the elastic anisotropy of minerals (Murri et al., 2018;Mazzucchelli et al., 2019), and the non-ideal geometry of inclusions . Moreover, precisely determining the P inc in inclusions has attracted many studies, as measured by x-ray diffractometry (Angel et al., 2014b) or Raman spectroscopy (Enami et al., 2007;Korsakov et al., 2009;Ashley et al., 2016;Murri et al., 2018;Bonazzi et al., 2019). However, researchers rarely explored the effect of thermoelastic parameters of minerals in entrapment pressure (P e ) calculations, which is crucial to improving elastic geobarometry precision (Angel et al., 2015;Moulas et al., 2020). Only Angel et al. (2015) considered the influence of thermoelastic parameters of diamond in determining the pressure of formation of diamondinclusion systems. Moreover, as garnet rarely has a homogeneous composition, Moulas et al. (2020) considered the possible propagation of errors caused by the uncertainty of garnet composition in entrapment pressure (P e ) calculations. Therefore, in this study, we evaluated the influence of the thermoelastic parameters of the concerned minerals on the calculated results of the two types of elastic geobarometers, garnet (inclusion)-diamond (host) system and quartz (inclusion)-garnet (host) system (QuiG). The garnet-diamond system is very appropriate for estimating the formation pressure of diamond because the garnet is the most abundant occurring inclusion mineral in the diamond (Milani et al., 2015). Similarly, QuiG is also suitable for geobarometers because quartz is highly compressible relative to garnet, and the garnet host can maintain a large stress difference caused by changes in P-T conditions (Ashley et al., 2016). In these two types of elastic geobarometers, diamond and quartz are generally virtually pure endmember, and their thermoelastic parameters have been accurately determined (Angel et al., 2015;Angel et al., 2017a). However, the effect of composition of garnet on its thermoelastic parameters is still not very clear Fan et al., 2018;Xu et al., 2019).
Understanding the thermoelastic properties of garnet, especially the Prp-Alm solid solution, is not only important to model the structure of the subducted slab but also a primary step in quantifying the P-T conditions of subduction and exhumation (Bass et al., 2008;Ashley et al., 2014). There are already extensive studies about the thermoelastic properties of Prp and Alm endmember component (Zhang L. et al., 1998;Gwanmesia et al., 2006;Fan et al., 2009;Li et al., 2011;Zou et al., 2012;Arimoto et al., 2015;Chantel et al., 2016;Hu et al., 2016;Fan et al., 2019b;Hartwig and Galkin, 2021). However, there is still a lack of tight constraints on the variation behavior of thermoelastic parameters developed when Prp-Alm solid solution is varying with composition (Takahashi and Liu, 1970;Huang and Chen, 2014;Milani et al., 2015). It is known that, for different mineral compositions, a variety of different cations incorporated in the garnets will affect their thermoelastic properties (Zhang et al., 1999;Kuang et al., 2019). In the 1970s, Takahashi and Liu (1970) conducted the pressure-volume (P-V) EoS study on four garnet samples with different compositions (Prp 100 , Prp 60 Alm 31 , Prp 22 Alm 72 , and Alm 100 ) up to 32.8 GPa using single-crystal x-ray diffraction (SCXRD) and diamond anvil cell (DAC). They inferred that the bulk modulus (K 0 ) of the Prp-Alm solid solution decreases with the increasing Alm content. However, the conclusion from Takahashi and Liu (1970) needs to be further verified because they did not use any gasket in their high-pressure SCXRD experiments, which may significantly influence their experimental results. Similarly, Huang and Chen (2014) also performed the P-V EoS study on three synthetic Prp-Alm solid solutions (Prp 83 Alm 17 , Prp 54 Alm 46 , and Prp 30 Alm 70 ) up to~21 GPa using in situ XRD combined with DAC, but obtained the opposite conclusion from Takahashi and Liu (1970), that is, the K 0 increases with the increasing Alm content. Furthermore, Huang and Chen (2014) also obtained the relationship between the K 0 of the Prp-Alm solid solutions and Alm content: K 0 = 170.5 (2.6) + 0.12 (4) X Alm (X Alm is the mole fraction of Alm). Subsequently, Milani et al. (2015) studied the P-V EoS on three synthetic Prp-Alm solid solutions (Prp 100 , Prp 60 Alm 40 , and Alm 100 ) and obtained the same conclusion as that of Huang and Chen (2014). Nevertheless, to date, the thermoelastic properties of the Prp-Alm solid solutions at simultaneous high-pressure and high-temperature (HP-HT) conditions are still not entirely well constrained.
In this study, we selected eight synthetic Prp-Alm solid solution samples ( ) for the synchrotron SCXRD measurements up to~20 GPa and 700 K. Based on the obtained thermoelastic parameters of the systematic Prp-Alm solid solutions, we discussed the effects of composition (Mg 2+ -Fe 2+ substitution) on the thermoelastic parameters of Prp-Alm solid solutions. Finally, we used these results to calculate the P e in two host-inclusion systems and discussed the effects of the thermoelastic parameters of Prp-Alm-rich garnet on the results of P e calculations.

SAMPLES AND EXPERIMENTAL METHODS
The single crystals of Prp-Alm solid solutions (Prp 9 Alm 91 , Prp 14 Alm 86 , Prp 23 Alm 76 , Prp 31 Alm 68 , Prp 48 Alm 52 , Prp 57 Alm 43 , Prp 67 Alm 32 , Prp 78 Alm 22 , and Prp 87 Alm 12 ) used in this study were synthesized at HP-HT conditions using a multi-anvil pressure apparatus (YJ-3000t) at the Institute of Geochemistry, Chinese Academy of Sciences, Guiyang, China. The detailed information about the garnet sample synthesis and subsequent compositional and structural analysis are described by Kuang et al. (2019). In general, with the increasing Alm content, the typical size of granular Prp-Alm crystals increases from~60 μm to more than 100 μm. The Prp-Alm solid solutions with good euhedral crystals, high quality, and without inclusions were selected for the SCXRD experiments. We successively double-side polished our sample pieces to~20-30 μm thickness, using the 3M diamond lapping films with grain size from coarse to a final 1 μm. The polished platelets were then cleaved into several square pieces of the desired size (~40-60 μm) for SCXRD experiments.
High pressure and room temperature (HP-RT) SCXRD experiments were conducted on eight Prp-Alm solid solution samples using three short symmetrical DACs (cell-1: Prp 9 Alm 91 , Prp 31 Alm 68 , Prp 48 Alm 52 ; cell-2: Prp 57 Alm 43 , Prp 67 Alm 32 , Prp 78 Alm 22 , and Prp 87 Alm 12 , and cell-3: Prp 23 Alm 76 ). Rhenium (Re) gaskets of~250 μm thickness were pre-indented to~65 μm thickness using a pair of 500-μm culet size diamond anvils. Subsequently, a~300-μm-diameter cylindrical hole was drilled in the pre-indented area as the sample chamber. The selected single-crystal platelets with diameters of~40-60 μm were then placed into the sample chamber together with gold powder, which served as the pressure standard (Fei et al., 2007). The diffraction patterns of gold were collected before and after sample data collection for each pressure, and the average pressure values were used for the EoS fitting. A small ruby sphere of~10 μm was also loaded into each DAC sample chamber and used as the pressure indicator during the neon gas-loading. The neon pressure medium was loaded into the sample chamber using the gas-loading system at GeoSoilEnviroCARS (GSECARS) (Rivers et al., 2008). After the neon gas loading, the diameter of the sample chamber was~220 μm.
HP-HT SCXRD experiments were carried out using an externally heated DAC (EHDAC) equipped with an alumina ceramic heater coiled with two Pt wires of 200 μm in diameter and 48 cm length (Kantor et al., 2012;Fan et al., 2019a). A rhenium foil was used as the gasket material and pre-indented tõ 60 μm thickness using a pair of diamond anvils with 500-μm culet size and then a 300-μm-diameter sample chamber was drilled at the center of pre-indentation. The selected singlecrystal platelets (Prp 14 Alm 86 , Prp 31 Alm 68 , Prp 48 Alm 52 , and Prp 67 Alm 32 ) with a diameter of~40 μm were loaded into the sample chamber together with gold powder, which served as the pressure standard at HP-HT conditions. For each P-T point, the diffraction patterns of gold were collected before and after sample data collection, and the average pressure values were used. Likewise, a small ruby sphere of~10 μm was loaded as the pressure indicator for the gas-loading with neon as the pressure transmitting medium using the GSECARS gasloading system. Heating was carried out by a resistanceheating system (Fan et al., 2010). A Pt 90 Rh 10 -Pt 100 thermocouple was attached to one of the diamond surfaces, approximately 500 μm away from its culet, and clad with a ceramic adhesive (Resbond 920) for temperature measurements (Kantor et al., 2012). To minimize temperature instability for each heating run, we first heated the sample chamber to an expected temperature and then kept it at this temperature for 15 min. The temperatures of the sample chamber were actively stabilized within ±1 K using the temperature-power feedback program with a remotely controlled Tektronix Keithley DC power supply during the experiments (Sinogeikin et al., 2006).
Both of the HP-RT and HP-HT synchrotron SCXRD experiments were conducted at the 13-BMD beamline of the GeoSoilEnviroConsortium for Advanced Radiation Sources (GSECARS) of Advanced Photon Source (APS), Argonne National Laboratory (ANL). An incident x-ray beam of 0.3344 Å wavelength focused on a 3 × 7 μm 2 area was used to determine the unit-cell volume of Prp-Alm solid solution samples in the DACs. Diffraction images were acquired on a stationary Perkin-Elmer area detector. Tilting and rotation of the detector and the sample-to-detector distance were calibrated using ambient LaB 6 as the diffraction standard. Wide and stepped φ-rotation exposure were collected for the single-crystal samples at each loading run, with an exposure time of 2 s/deg. The φ-rotation (the opening angle of the DAC is ±15°) axis was horizontal and perpendicular to the incident x-ray direction. The diffraction images collected at each P-T point were analyzed using the GSE_ADA/RSV software (Dera et al., 2013).

Equations of State of Prp-Alm Solid Solutions at High Pressure and Room Temperature
The unit-cell volumes change with pressure at room temperature for eight Prp-Alm solid solution garnets and are presented in Supplementary Table S1 and Figure 1, showing that the unit-cell volumes of Prp-Alm solid solution with increasing pressure garnets decrease systematically. Using the console program EosFit7c (Angel et al., 2014a;Gonzalez-Platas et al., 2016), the EoS parameters, including the room P-T volume, the isothermal bulk modulus, and its pressure derivative of these Prp-Alm solid solutions, were refined by fitting their P-V data using two different EoSs (Birch-Murnaghan EoS and Vinet EoS).
Birch-Murnaghan equation of state. The P-V relations of these Prp-Alm solid solutions have been determined by fitting their HP-RT data to the third-order Birch-Murnaghan EoS, which is given in the following form: where V 0 , K 0 , and K′ 0 are the unit-cell volume, isothermal bulk modulus, and its pressure derivative at ambient conditions, respectively. We obtained the V 0 , K 0 , and K′ 0 of Prp-Alm solid solutions, as shown in Table 1. The refined value of V 0 was within 1σ or so of the measured V 0 by XRD at ambient conditions, indicating the accuracy of the refined results (Angel, 2000). It can be found that the bulk moduli of Prp-Alm solid solutions gradually increase with increasing Alm content, especially when K′ 0 was fixed at 4 (Table 1). Furthermore, this trend is also in agreement with the conclusion of Huang and Chen (2014) and Milani et al. (2015). Frontiers in Earth Science | www.frontiersin.org January 2022 | Volume 9 | Article 833405 Vinet equation of state. We also analyzed the P-V data using the Vinet EoS (Vinet et al., 1986), which is derived from a "universal equation" for solids and is represented as follows: where y = x 1/3 , x = V/V 0 , and ƞ 0 = (3/2) (K′ 0-1). Analyses of Equation 2 yielded V 0 , K 0, and K′ 0 for the Prp-Alm solid solutions (  Milani et al. (2015) reported the least value (163.7 GPa) due to their largest K′ 0 value (6.4) than those of other studies (4.3-4.6). There is a trade-off between the bulk modulus and its pressure derivative, which negatively correlates (Gatta et al., 2011). So, we obtained the confidence ellipses of K 0 and K′ 0 for pyrope from previous XRD studies, as shown in Figure 2. From Figure 2, we can see that the previous XRD results of pyrope all fall within the 68.3% (1σ) confidence interval, with mean values of 171 (8) GPa and 4 (2) for K 0 and K′ 0 , respectively. Meanwhile, the obtained K 0 of pyrope [170.0 (2)-175 (2) GPa] by using Brillouin light spectroscopy (BLS) and ultrasonic interferometry (UI) methods are also consistent with the results from the XRD method (Supplementary Table S2 and Figure 2). For the Alm endmember, the K 0 value [172.6 (15) GPa] obtained by Milani et al. (2015) was also smaller than other XRD results [179 (3)-185 (3) GPa]. The main reason remains that the K′ 0 value (5.6) obtained by Milani et al. (2015) is higher than those of other studies (4-4.2). Moreover, Arimoto et al. (2015) reported that the K 0 value [174.2 (12) GPa] of Alm endmember using the UI method is slightly smaller than the result [179 (3) GPa] using the XRD method.
As a comparison, until now, the previous EoS studies on the Prp-Alm solid solutions are relatively limited. Huang and Chen (2014) measured three synthetic Prp-Alm solid solutions (Prp 83 Alm 17 , Prp 54 Alm 46 , and Prp 30 Alm 70 ) by the XRD method up to 7, 21, and 19 GPa, respectively. With fixed K′ 0 at 4.3, the K 0 obtained for these three samples were 172 (4), 174 (2), and 183 (2) GPa, respectively. Subsequently, Milani et al. (2015) reported a significantly smaller value [167.2 (17) GPa], also using the XRD method for synthetic Prp 60 Alm 40 . As a consequence, the previous experimental studies indicated that the K 0 range along the Prp-Alm solid solutions from Prp to Alm endmembers is expected to be 163.7 (1)-185 (3) Table S2). Our results show that the bulk moduli of our eight Prp-Alm solid solutions ( Table 1) are within this range. Because the bulk modulus and its pressure derivative have a negative correlation of bulk modulus and its pressure derivative, the K′ 0 is permanently fixed at 4 for the following discussion of the relationship between the bulk moduli of Prp-Alm solid solutions and their compositions. Figure 3 shows the bulk moduli of Prp-Alm solid solutions as a function of composition, and we can observe that the compressibilities of Prp-Alm solid solutions increase with increasing Prp concentrations. Furthermore, by linear fitting of the results in this study, the bulk moduli of the Prp-Alm solid solutions as a function of Prp content can be expressed by K 0 (GPa) = 181.0 (8)-0.11(1) X prp (R 2 = 0.91, K′ 0 fixed at 4). It is worth noting that the K 0 of most previous studies are in harmony with our fitting line, except for Milani et al. (2015) (Figure 3). It may be attributed to the    Table S2).

Equations of State of Prp-Alm Solid Solutions at High Pressure and High Temperature
The unit-cell volumes of Prp 14 Alm 86 , Prp 31 Alm 68 , Prp 48 Alm 52 , and Prp 67 Alm 32 at HP-HT conditions up to~19.67 GPa and 700 K are shown in Supplementary Table S3. Figure 4 also shows the pressure-volume-temperature (P-V-T) data of these Prp-Alm solid solution samples. In this study, we applied the thermal-pressure model EoS to describe the P-V-T behavior of Prp-Alm solid solutions. There are two ways to calculate the P-V-T properties of a mineral (as shown in Figure 1 of Milani et al., 2017). One way to characterize the P-V-T behavior of minerals is to describe how the V T0 , K T0 , and K' T0 vary with temperature (at room pressure) and then use these parameters to calculate the isothermal compressibility at the high temperature (as shown by the rough path in Figure 1 of Milani et al., 2017). A potential weakness of this method is that the variation of K' T0 with temperature has rarely been measured, where K' T0 should increase slightly with increasing temperature, rather than not change with temperature, as the common assumption. Additionally, this assumption, coupled with the approximation that dK T0 /dT is a constant, has often given rise to a non-physical negative coefficient of thermal expansion for many materials at reasonably moderate pressures (Helffrich and Connolly, 2009).
Furthermore, the maximum experimental temperature of only 700 K in this study and the relatively limited number of hightemperature data points make the fitting of the thermal expansion coefficient at high temperatures more inaccurate. Therefore, we employ the concept of thermal pressure (Anderson, 1995) to describe the P-V-T behavior of Prp-Alm solid solutions. Then, as shown by the thin path in Figure 1 of Milani et al., 2017, the P at a given V and T can be expressed as the sum of two terms: The P(V, T ref ) is the isothermal EoS for minerals at the reference temperature, up to the volume at P and T ref . The thermal-pressure Pth(T) is the pressure generated by raising the temperature from T ref to T at a constant volume, which is the isochor of the mineral passing through the final P-T point. Details of the employment of a thermal-pressure EoS and any necessity assumptions have been reviewed by Milani et al. (2017). In this study, we used the Birch-Murnaghan isothermal EoS in combination with the Holland and Powell (2011) model of thermal-pressure EoS to determine the thermoelastic parameters of four Prp-Alm solid solutions, as shown in Table 2. All described in this study have been performed with the EosFit7c program (Angel et al., 2014a;Gonzalez-Platas et al., 2016). The Einstein temperature, θ E , of Prp-Alm solid solutions in this study are calculated from the values of the pyrope (320 K) and almandine (400 K) and according to the ideal solution model (Milani et al., 2015). The V 0 and K T0 derived from HP-HT data were roughly consistent with that from P-V data at room temperature, especially as the K' T0 was fixed at 4 ( Table 1 and Table 2). Figure 5 shows the thermal expansion coefficient at ambient conditions (α 0 ) of four Prp-Alm solid solutions in this study, combined with the results from previous studies. As shown in

Mg 2+ -Fe 2+ Substitution Effect on Bulk Moduli of Prp-Alm Solid Solutions
Among previous studies about the compositional effect on bulk modulus of different mineral groups, minerals with Mg 2+ -Fe 2+ solid solutions have been extensively investigated because of their importance in understanding the composition and dynamic property of the Earth's interior (Takahashi and Liu, 1970;Zhang J. et al., 1998;Zhang et al., 1999;Speziale et al., 2004;Huang and Chen, 2014). The Fe and Mg component dependence of bulk modulus in mantle minerals is interesting for modeling the mineralogical composition of the mantle and mineral crystal chemistry. Anderson and Anderson (1970) proposed a simple bulk modulus (K 0 )-volume (V 0 ) systematics, where K 0 V 0 = constant. However, for most of the ferromagnesian silicates and carbonates that have been studied so far, the bulk modulus appears to increase with increasing Fe content, such as garnet (Huang and Chen, 2014), olivine (Speziale et al., 2004), and siderite (Zhang J. et al., 1998), which is inconsistent with bulk modulus-volume systematics. To sum up, it is necessary to understand the Mg 2+ -Fe 2+ substitution effect on the bulk modulus of garnet. Figure 3 shows the bulk moduli of Prp-Alm solid solutions as functions of compositions in this study and compared with those of previous studies. Mg 2+ -Fe 2+ substitution results in a nearly linear decrease of the bulk modulus. The reduction of bulk moduli from Prp 9 Alm 91 to Prp 87 Alm 12 is~3%, mainly consistent with the~6% decrease from Prp 30 Alm 70 to Prp 83 Alm 17 reported by Huang and Chen (2014).
The bulk moduli of minerals are mainly controlled by their constituent polyhedral compressibilities, primarily determined by the mean cation-oxygen bond length in the polyhedral (Hazen et al., 1994). Many previous high-pressure studies have shown that the SiO 4 tetrahedrons of silicate minerals are essentially rigid units (Robinson et al., 1971;Zhang et al., 1997;Xu et al., 2017). The selected bond distances of Prp-Alm solid solutions at the ambient condition (data from Kuang et al., 2019) as functions of composition are shown in Figure 6. Inspection of Figure 6 demonstrates that the d <M-O> (M is the divalent cations in X-site) decreases significantly with the increasing Prp content, mainly due to the smaller effective ionic radius of Mg 2+ (0.89 Å) in eightfold coordination than that of Fe 2+ (0.92 Å) (Shannon, 1976). In contrast, the average d <Si-O> of the SiO 4 tetrahedron and d <Al-O> of the AlO 6 octahedron change negligibly with increasing Prp mol%. These results imply that the compressibility of the MO 8 dodecahedron mainly controls the bulk modulus of the Prp-Alm solid solution. However, other factors, such as the compressibilities of the SiO 4 tetrahedron and the AlO 6 octahedron, may be secondary important to the bulk modulus of Prp-Alm solid solution.
Besides, this study shows that FeO 8 dodecahedron is less compressible than MgO 8 dodecahedron. Similar conclusions were also obtained about the structures with sixfold coordinated Mg 2+ and Fe 2+ cations (Zhang et al., 1997), implying that the bonding character plays an essential role in such anomalous compression, whereas structure type and polyhedral coordination may be less critical. Furthermore, the electronegativity of the element, which describes the ability of an atom to attract electrons, may be an essential factor for the bonding character. The electronegativity of Mg (1.31) is distinctly smaller than that of Fe (1.83) (Allred, 1961), resulting in the less compressible FeO 8 dodecahedron than MgO 8 dodecahedron.
Moreover, Hazen et al. (1994) suggested that the SiO 4 octahedra and AlO 6 tetrahedra framework may also dictate the bulk modulus of garnet. Subsequently, Zhang L. et al. (1998) showed that the kinking degree of the Si-O-Al angle decreased continuously with increasing pressure, and the Si-O-Al angle in Alm is less kinked than that in Prp over the same pressure range (Zhang et al., 1999). Furthermore, the Raman spectroscopy results of Kuang et al. (2019) show that both the stretching motions of Si-O and the rotation motions of SiO 4 decreased linearly with increasing Alm content. Thus, the degree of distortion and rotation of SiO 4 may also affect the bulk moduli of the Prp-Alm solid solutions.
To sum up, the bulk moduli of Prp-Alm solid solutions are primarily governed by the compressibility of MO 8 dodecahedra. The difference in the compressibilities of MO 8 can be attributed to the different effective ionic radius and electronegativities of the Mg and Fe atoms. Also, the kinking degree of the Si-O-Al angle and the distortion and rotation of SiO 4 may affect the bulk moduli of Prp-Alm solid solutions.

Mg 2+ -Fe 2+ Substitution Effect on Thermal Expansion Coefficient of Prp-Alm Solid Solutions
Contrary to the effect of Mg 2+ -Fe 2+ substitution on the bulk modulus of Prp-Alm solid solutions, the α 0 of Prp-Alm solid solutions varies with the composition forming a parabola. As shown in Figure 5, at ambient temperatures, the composition of the Prp-Alm solid solution near the intermediate composition displays a larger thermal expansion coefficient than the composition near the endmember. While the causes of this result are yet to be revealed, we think that the most likely reason is the local arrangement of Mg 2+ and Fe 2+ . For the garnet, the different-sized divalent cations substitution on the dodecahedral site can make its thermodynamic properties (e.g., volume, enthalpy, entropy, free energy, etc.) deviate from the ideal mixing results (Bosenick et al., 2001;Vinograd et al., 2004;Dachs et al., 2014;Du et al., 2016). Besides, the magnitude of the non-ideal mixing is correlated with the size difference between Frontiers in Earth Science | www.frontiersin.org January 2022 | Volume 9 | Article 833405 the substituting cations on the dodecahedral site (Geiger and Rossman, 1994). The study on the local structure of the pyropegrossular garnet solid solution shows that the Mg/Ca substitution can produce significant amounts of short-range order rather than long-range order in garnets (Bosenick et al., 2000). Furthermore, the compositions near Prp 50 Grs 50 depart most strongly from random mixing, and all the compositions tend to random disorder at high temperatures. Although the size difference between Mg 2+ (0.89Å) and Fe 2+ (0.92Å) is smaller than that between Mg 2+ and Ca 2+ (1.12 Å) (Shannon, 1976), the influence on the local arrangement by the Mg/Fe substitution cannot be overlooked. It is speculated that local disorder caused by the Mg/ Fe substitution is the main reason for a non-linear relationship between thermal expansion coefficient and the composition of the Prp-Alm solid solution.

IMPLICATIONS FOR ELASTIC GEOBAROMETER
It is ordinary to use host-inclusion systems to infer the geological history of a mineral or a rock (Guiraud and Powell, 2006). The history of study on elastic geobarometers can be traced back to the 1980s. Roedder and Bodnar (1980) made the first estimation for the formation pressure of minerals from the fluid inclusions in them. The basic principle is that when the host mineral captures the fluid, it experiences a different P-T path than the host mineral because its volume is limited by the space inside the mineral and is subjected to the pressure exerted on it by the host mineral during subsequent changes in temperature and pressure conditions. In turn, this difference and its associated physicochemical properties can be used to estimate the P-T conditions experienced by the host mineral. This principle is not only applicable to fluid/melt inclusions but also to other solid mineral inclusions within the mineral that differ from its thermoelastic properties. This method underwent significant developments in the past decade and is gradually gaining the interest of geologists (e.g., Murri et al., 2018;Alvaro et al., 2019;Anzolini et al., 2019;Bonazzi et al., 2019;Gonzalez et al., 2019;Mazzucchelli et al., 2019;Moulas et al., 2020). The pressure estimation by an elastic geobarometer has some prerequisite assumptions (Bonazzi et al., 2019): (1) The inclusions are in the same pressure environment as the host mineral at the time of capture, and the inclusions fill the vacancies entirely within the host mineral. (2) The host mineral and the inclusions are always reversibly and elastically deformed after capture. In this way, the captured pressure can be estimated by measuring the residual pressure of the inclusions due to the difference in elastic properties between the two minerals under normal temperature and pressure conditions and by combining the equation of the state of the two minerals. In this study, the method for estimating the pressure is based on entrapment isomekes proposed by Angel et al. (2017b). Entrapment isomekes is a single P-T path along which the fractional volume change in host and inclusion are equal (Gonzalez et al., 2019). The calculation of P e was performed with EosFit-Pinc, and the details of this method were described in Angel et al. (2017b).
As mentioned in the introduction section of this study, much of the previous work has focused on the factors affecting the calculation accuracy of the elastic geobarometer except the influence of the thermoelastic parameters of minerals. So, we selected two garnet-related host-inclusion pairs to evaluate the influence of the thermoelastic parameters of the concerned minerals on the calculated results of elastic geobarometers. Frontiers in Earth Science | www.frontiersin.org January 2022 | Volume 9 | Article 833405 9 Quartz is generally a virtually pure endmember, and it is stable over a wide range of metamorphic conditions (Korsakov et al., 2009). Garnet has also served as a central focus of metamorphic studies to discern P-T paths because it can store the entire history of metamorphism by its inclusion suites such as quartz, graphite, and rutile (Ague and Carlson, 2013;Baxter et al., 2013;Caddick and Kohn, 2013). Hence, QuiG is considered a particularly useful geobarometer and has been widely used to infer the geological history of the rock from the collisional orogenic belts (Wolfe and Spear, 2018), the subduction complex (Bayet et al., 2018;Alvaro et al., 2019), and ultrahigh-pressure metamorphic complex (Korsakov et al., 2009;Gonzalez et al., 2019;Zhong et al., 2019;Moulas et al., 2020). Similarly, as the primeval samples from inaccessible regions of Earth, diamonds and their inclusions provide direct information of deep Earth's evolution. This information is precious if it can be combined with depth estimates (Anzolini et al., 2019). An estimate of their formative pressure (depth of provenance) can be used to constrain and understand the environment in which they formed. Therefore, the thermoelastic parameters of the minerals play a crucial role in determining the accuracy of the pressure calculation. The garnet, a widespread mineral solid solution in nature, can have varying elastic properties depending on its composition (Du et al., 2015;Fan et al., 2017;Fan et al., 2018;Xu et al., 2019). In most cases, however, the thermoelastic properties used for calculation are the results of the endmember garnet (Milani et al., 2015;Gonzalez et al., 2019;Zhong et al., 2019), or assumed that the solid solution is isotropic and ideal (Johnson et al., 2021). Based on the results of this study, the thermoelastic properties of the Prp-Alm solid solutions may not follow the ideal mode (e.g., Figure 5). Therefore, the thermoelastic properties of garnet estimated from the endmember components may cause divergence. It is necessary to discuss the influence of garnet composition on the accuracy of these two elastic geobarometers.
These calculations utilized the thermoelastic coefficients of quartz with a curved α-β transition (Angel et al., 2017a) and the thermoelastic coefficients of the diamond from Angel et al. (2015). The shear modulus and its P and T derivative of Prp-Alm garnet were calculated by combining the results of the endmembers pyrope and almandine in proportion to their measured abundances, assumed as an ideal solid solution (Arimoto et al., 2015;Fan et al., 2019b). The final full set of thermoelastic coefficients for the P e estimation are shown in Table 3.
Firstly, the uncertainties of the thermoelastic properties of garnets in the QuiG application should be considered. Here, we take the Prp 14 Alm 86 garnet to demonstrate the effect of the uncertainties in thermoelastic properties. We chose 0.65 GPa as the P inc obtained by laser Raman spectroscopy for quartz inclusions in garnet from eclogite for our P e estimation (Zhong et al., 2019). Supplementary Table S4 lists the calculated P e for Prp 14 Alm 86 considering the uncertainties of the thermoelastic properties. Compared to the bulk modulus, the thermal expansion coefficient in this study has the most considerable effect in calculated P e ( Figure 7A). A 1.7% uncertainty in the value of α 0 of garnet  Arimoto et al. (2015) and Fan et al. (2019b), assuming it is an ideal solid solution. b Values were calculated by α 300 (10 −5 K −1 ) = 2.7 (1) + 3.0 (5) X Prp − 3.2 (4) X 2 Prp . c Values were calculated by the endmember Prp and Alm result of Milani et al. (2015), assuming it is an ideal solid solution.
changes the P e by~0.007 GPa at 500 K, and with the temperature rising, it makes the P e change~0.016 GPa at 1300 K. The accuracy of bulk modulus obtained in this study has almost no effect on the P e estimation. A 0.22% uncertainty in the value of K 0 of garnet changes the P e by only 0.001 GPa. In conclusion, the maximal error caused by the thermoelastic properties of garnet on P e estimation is less than 0.02 GPa ( Figure 7A). Besides, as shown in Figure 7A, we have calculated the P e of the other two garnets with different components (Prp 48 Alm 52 and Prp 87 Alm 12 ), and the results obtained are consistent with the Prp 14 Alm 86 .
For the P e calculations of the diamond-garnet pair, we assume a residual pressure P inc on the garnet inclusion of 0.2 (1) GPa, and we use the thermoelastic properties given in Table 3 Figure 7B are that a 0.22% uncertainty of K 0 may contribute an error of~0.03 GPa to the calculated P e . Same as in QuiG, the thermal expansion coefficient in this study has the most considerable effect in the calculated P e of the diamond-garnet system ( Figure 7B). A 1.7% uncertainty in the value of α 0 of garnet changes the P e by~0.092 GPa at 800 K, and with the temperature rising, it makes the P e change~0.254 GPa at 1600 K. The maximal error caused by the thermoelastic properties of garnet in the diamond-garnet system is~0.3 GPa at 1,600 K ( Figure 7B). It can be found that for the same garnet composition, the P e change caused by the uncertainties of its thermoelastic properties is unequal in the two types of elastic geobarometers. The P e change in the diamond-garnet pair is almost ten times as large as it is in QuiG.
Here, the Pyp-Alm garnets of different compositions (Table 3) are used to calculate the P e , and the results of two elastic geobarometers are shown in Figures 8A,B. The effect of garnet composition on the predicted entrapment pressure is relatively large, especially for the diamond-garnet system ( Figure 8B). At 1,400 K, the different compositions of garnet can cause a variation of P e up to 4.5 GPa in the diamond-garnet system. In the QuiG, the difference of inferred entrapment pressure may reach ca. 0.35 GPa at 1,060 K. As mentioned  Frontiers in Earth Science | www.frontiersin.org January 2022 | Volume 9 | Article 833405 above, this is mainly due to the distinct thermal expansion coefficient among them. To sum up, even though the procedure of applying an elastic barometer to host-inclusion mineral systems requires several steps that involve various assumptions (Angel et al., 2015), and there are lots of possible factors that could influence the precision and accuracy of calculated P e value as mentioned above, its advantages could not be ignored. Much previous work suggested that QuiG, and diamond-garnet elastic barometer, may preserve different aspects of garnet history in complex terranes (Wolfe and Spear, 2018;Alvaro et al., 2019;Spear and Wolfe, 2020;Johnson et al., 2021). The advantages of QuiG or other inclusion-based methods are relatively easy to apply and require far less data acquisition than conventional methods. Our results suggest that the entrapment pressure may be significantly influenced by the thermoelastic parameters of garnet, which is selected to have a different component from the research sample, while the thermal expansion coefficient is the main influencing factor for the elastic barometer to P e estimation, compared with other thermoelastic parameters. Therefore, the effect of the composition of garnet on the P e estimation could not be ignored. Furthermore, one factor that must not be overlooked is the effect of viscous relaxation, which occurs very rapidly at temperatures >110 K  and the preserved P e may then be a reflection of re-equilibration (Zhong et al., 2020). In the future application of elastic barometer, as the most common garnet component in the Earth, the thermoelastic parameters (especially the thermal expansion coefficient) of the pyrope-almandine-grossular ternary solid solutions must also be considered.

CONCLUSION
In summary, the thermoelastic parameters of systematic Prp-Alm solid solutions were obtained by fitting the P-V-T or P-V data to the EoS and compared with previous studies. The relationship between the bulk moduli of the Prp-Alm solid solutions and Prp content can be expressed as K 0 (GPa) = 181.0 (8) − 0.11 (1) X prp (R 2 = 0.91, K′ 0 fixed at 4). Our results demonstrate that the Prp content will decrease the incompressibilities of Prp-Alm solid solutions, which is inconsistent with bulk modulus-volume systematics: K 0 V 0 = constant. The possible reasons for this difference may be multiple and complicated. In our view, the most likely reason is the larger compressibility of MgO 8 than that of FeO 8 . Also, the kinking degree of the Si-O-Al angle and the degree of distortion and rotation of SiO 4 may affect the compressibilities of the Prp-Alm solid solutions. The thermal expansion coefficient with the component of Prp-Alm solid solution also has a relational expression, α 0 (10 −5 K −1 ) = 2.7 (1) + 3.0 (5) X Prp − 3.2 (4) X 2 Prp (R 2 = 0.985), where the Prp-Alm solid solutions with intermediate composition show a larger coefficient than the endmembers at ambient condition. Furthermore, for two types of elastic geobarometers, we evaluate the correlation of the uncertainties of thermoelastic properties with P e estimation accuracy. Compared with bulk modulus, the thermal expansion coefficient has the main effect on the estimation of P e . The P e change caused by the uncertainties of thermoelastic properties of garnets in the diamond-garnet pair is almost 10 times larger than it is in QuiG. Besides, by calculating P e with different garnet components, the result suggests that the garnet component may significantly influence entrapment pressure, especially in the diamond-garnet system and the pyrope-almandine-grossular ternary solid solutions EoS must also be considered in the future study.

DATA AVAILABILITY STATEMENT
The original contributions presented in the study are included in the article/Supplementary Material, further inquiries can be directed to the corresponding author.