hydrogen – water

$$\mathrm{2H_2 + O_2 = 2H_2O}$$

$$f\mathrm{O}_2=\left(\frac{X\mathrm{H_2O}}{X\mathrm{H_2}}\right)^2\times\exp\left(\frac{\Delta_rG}{RT}\right)$$ $$\Delta_rG = 2\left(\Delta_fH_\mathrm{H_2O} + H^{\circ}_\mathrm{H_2O} - H^{\circ}_\mathrm{298.15,H_2O} - TS^{\circ}_\mathrm{H_2O}\right) - 2\left(\Delta_fH_\mathrm{H_2} + H^{\circ}_\mathrm{H_2} - H^{\circ}_\mathrm{298.15,H_2} - TS^{\circ}_\mathrm{H_2}\right) - \left(H^{\circ}_\mathrm{O_2} - H^{\circ}_\mathrm{298.15,O_2} - TS^{\circ}_\mathrm{O_2}\right)$$ Values for enthalpy and entropy are evaluated using a Shomate equation: $$H^{\circ} - H^{\circ}_\mathrm{298.15} = AT + \frac{BT^2}{2}+ \frac{CT^3}{3} + \frac{DT^4}{4} - \frac{E}{T} +F - H$$ $$S^{\circ} = A\ln{T} + BT + \frac{CT^2}{2} + \frac{DT^3}{3} - \frac{E}{2T^2} + G$$ All parameters are given in the NIST Chemistry WebBook “Gas phase thermochemistry data” pages for water, hydrogen, and oxygen.

This buffer is calculated at atmospheric pressures (1 bar), regardless of the pressure setting. It is only directly comparable to other buffers when they are likewise calculated at 1 bar. It is valid from 500 K to 6000 K, but does not take into account gas dissociation into monoatomic molecules.

sulfur dioxide – hydrogen sulfide

$$\mathrm{\frac{2}{3}H_2S + O_2 = \frac{2}{3}SO_2 + \frac{2}{3}H_2O}$$

$$f\mathrm{O}_2=\left(\frac{X\mathrm{SO_2}}{X\mathrm{H_2S}}\times f\mathrm{H_2O}\right)^{2/3}\times\exp\left(\frac{\Delta_rG}{RT}\right)$$ $$\Delta_rG = \frac{2}{3}\left(\Delta_fH_\mathrm{SO_2} + H^{\circ}_\mathrm{SO_2} - H^{\circ}_\mathrm{298.15,SO_2} - TS^{\circ}_\mathrm{SO_2}\right) + \frac{2}{3}\left(\Delta_fH_\mathrm{H_2O} + H^{\circ}_\mathrm{H_2O} - H^{\circ}_\mathrm{298.15,H_2O} - TS^{\circ}_\mathrm{H_2O}\right)$$ $$ - \frac{2}{3}\left(\Delta_fH_\mathrm{H_2S} + H^{\circ}_\mathrm{H_2S} - H^{\circ}_\mathrm{298.15,H_2S} - TS^{\circ}_\mathrm{H_2S}\right) - \left(H^{\circ}_\mathrm{O_2} - H^{\circ}_\mathrm{298.15,O_2} - TS^{\circ}_\mathrm{O_2}\right)$$ Values for enthalpy and entropy are evaluated using a Shomate equation: $$H^{\circ} - H^{\circ}_\mathrm{298.15} = AT + \frac{BT^2}{2}+ \frac{CT^3}{3} + \frac{DT^4}{4} - \frac{E}{T} +F - H$$ $$S^{\circ} = A\ln{T} + BT + \frac{CT^2}{2} + \frac{DT^3}{3} - \frac{E}{2T^2} + G$$ All parameters are given in the NIST Chemistry WebBook “Gas phase thermochemistry data” pages for water, hydrogen sulfide, and sulfur dioxide, and oxygen.

This buffer is calculated at atmospheric pressures (1 bar), regardless of the pressure setting. It is only directly comparable to other buffers when they are likewise calculated at 1 bar. It is valid from 500 K to 6000 K, but does not take into account gas dissociation into monoatomic molecules. All gases are considered as ideal gases.

nickel – nickel oxide (bunsenite)

$$\mathrm{2Ni + O_2 = 2NiO}$$

$$\log_{10}f\mathrm{O}_2=\frac{\mu\mathrm{O}_2 + 2\left(\int_1^P V^\mathrm{NiO}\,dP - \int_1^P V^\mathrm{Ni}\,dP\right)}{RT\ln10}$$

Parameterisation at 1 bar from O'Neill and Pownceby (1993): $$\mu\mathrm{O}_2=-478967+248.514T-9.7961T\ln{T}$$ Calibrated from 700 K to 1700 K.

Campbell et al (2009) use a Mie-Grüneisen–Birch-Murnaghan equation of state: $$P=\frac{3}{2}K_0\left(\left(\frac{V_0}{V}\right)^{7/3}-\left(\frac{V_0}{V}\right)^{5/3}\right)\left(1-\frac{3}{4}\left(4-K'_0\right)\left(\left(\frac{V_0}{V}\right)^{2/3}-1\right)\right)+\gamma_0\frac{V^{q-1}}{V_0^q}\left(E(\theta_D,T) - E(\theta_D,295)\right)$$ \(E\left(\theta_D,T\right)=3nRTD_3\left(\frac{\theta_D}{T}\right)\) where \(D_3()\) is the Debye function.

Direct calculation of \(\int_1^P V\,dP\) requires \(V=f(P)\) which is not available, so instead it is calculated using integration by parts: $$\int_{P_0}^{P_1} V\,dP = \Delta PV - \int_{V_0}^{V_1} P\,dV=P(V_1)V_1 - P(V_0)V_0 - \int_{V_0}^{V_1} P\,dV$$ $$\int P\,dV=\left(E(\theta_D,T) - E(\theta_D,295)\right)\frac{\gamma_0}{q}\left(\frac{V}{V_0}\right)^q-\frac{9K_0\left(\left(K_0'-4\right)V_0^3+\left(14-3K_0'\right)V^{2/3}V_0^{7/3}+\left(3K_0'-16\right)V^{4/3}V_0^{5/3}\right)}{16V^2}$$ Values for \(V_n\) are found numerically by solving for \(V\) which satisfy \(P(V)-P=0\).

All required parameters are given in Table 1 of Campbell et al (2009) and are calibrated to about 100 GPa and 2500 K.

iron – iron(II) oxide (wüstite)

$$\mathrm{2Fe + O_2 = 2FeO}$$

$$\log_{10}f\mathrm{O}_2=\frac{\mu\mathrm{O}_2 + 2\left(\int_1^P V^\mathrm{FeO}\,dP - \int_1^P V^\mathrm{Fe}\,dP\right)}{RT\ln10}$$

Parameterisation at 1 bar from O'Neill and Pownceby (1993): $$\mu\mathrm{O}_2=\begin{cases} -605568+1366.42T-182.7955T\ln{T}+0.10359T^2, & 833 < T < 1042 \\[2ex] -519113+59.129T+8.9276T\ln{T}, & 1042 \leq T \leq 1184 \\[2ex] -550915+269.106T-16.9484T\ln{T}, & 1184 < T < 1644\end{cases}$$ Calibrated at 1 bar from 833 K to 1644 K.

For FeO, Campbell et al (2009) use a Mie-Grüneisen–Birch-Murnaghan equation of state: $$P=\frac{3}{2}K_0\left(\left(\frac{V_0}{V}\right)^{7/3}-\left(\frac{V_0}{V}\right)^{5/3}\right)\left(1-\frac{3}{4}\left(4-K'_0\right)\left(\left(\frac{V_0}{V}\right)^{2/3}-1\right)\right)+\gamma_0\frac{V^{q-1}}{V_0^q}\left(E(\theta_D,T) - E(\theta_D,295)\right)$$ \(E\left(\theta_D,T\right)=3nRTD_3\left(\frac{\theta_D}{T}\right)\) where \(D_3()\) is the Debye function and \(n=2\).

Direct calculation of \(\int_1^P V\,dP\) requires \(V=f(P)\) which is not available, so instead it is calculated using integration by parts: $$\int_{P_0}^{P_1} V\,dP = \Delta PV - \int_{V_0}^{V_1} P\,dV=P(V_1)V_1 - P(V_0)V_0 - \int_{V_0}^{V_1} P\,dV$$ $$\int P\,dV=\left(E(\theta_D,T) - E(\theta_D,295)\right)\frac{\gamma_0}{q}\left(\frac{V}{V_0}\right)^q-\frac{9K_0\left(\left(K_0'-4\right)V_0^3+\left(14-3K_0'\right)V^{2/3}V_0^{7/3}+\left(3K_0'-16\right)V^{4/3}V_0^{5/3}\right)}{16V^2}$$ Values for \(V_n\) are found numerically by solving for \(V\) which satisfy \(P(V)-P=0\).

All required parameters are given in Table 1 of Campbell et al (2009) and are calibrated to about 100 GPa and 2500 K.

For Fe, the equation of state described by Dorogokupets et al (2017) is used. The equation is given in terms of the Helmholz free energy \(F\). Expressions for \(F(V,T)\) and \(P(V)\) are given in the Supplementary Files of Dorogokupets et al (2017). \(V(P,T)\) is calculated by numerically solving \(P(V)-P=0\) for the pressure of interest. \(G(P,T)\) is calculated using: $$G = F+PV$$ This Gibbs energy is used to find phase transitions, and the resulting pressures are used as the boundaries for integration. The \(\int_1^P V^\mathrm{Fe}\,dP\) integral is evaluated numerically.

This equation of state is calibrated to >300 GPa and >3000 K. The BCC α-Fe and γ-Fe phases are stable at low pressures. At high pressures, the FCC γ-Fe exists at high temperatures whereas the HCP ε-Fe is stable at lower temperatures. Evaluation of \(G\) leads to prediction of FCC γ-Fe at extreme pressures (>100–150 GPa). This phase transition is spurious and results from extrapolation beyond the calibrated range, and as such is not considered here.

Note relationship between three reactions in the Fe-O system. IW is only stable at temperature of about 830 K and higher. The invariant intersection point with the other buffers (IM and WM) shifts to lower temperature and higher oxygen fugacity with increasing pressure. See schematic diagram:

Wüstite is non-stoichiometric, and its non-stoichiometry is a function of oxygen fugacity, and is determined by the coexistence of iron as a buffering assemblage. The two studies from which wüstite data is used here contained iron, ensuring that wüstite had the appropriate stoichiometry for the IW buffer.

iron(II,III) oxide – iron(III) oxide
(magnetite – hematite)

$$\mathrm{4Fe_3O_4 + O_2 = 6Fe_2O_3}$$

$$\log_{10}f\mathrm{O}_2=\log_{10}f_{\mathrm{1\,bar}}\mathrm{O}_2+\log_{10}f_P\mathrm{O}_2$$

Parameterisation at 1 bar from Fegley (2013) Table 10-16, from 950 K to 1870 K: $$\log_{10}f_{1 bar}\mathrm{O}_2=-25.7139 - \frac{19375}{T} + 11.4191\log_{10}T$$

$$\log_{10}f_P\mathrm{O}_2 = \frac{6\int_1^P V^\mathrm{hm}\,dP - 4\int_1^P V^\mathrm{mt}\,dP}{RT\ln{10}}$$ Fugacity at pressure is calculated using a modified Tait equation of state including a thermal pressure term as described by Holland and Powell (2011): $$V=V_0\left(1-a\left(1-\left(1+b\left(P-P\mathrm{th}\right)\right)^{-c}\right)\right)$$ $$\int_{P_0}^P V\,dP=V_0\left(\left(P-P_0\right)\left(1-a\right) + \frac{a\left(\left(1+b\left(P-P\mathrm{th}\right)\right)^{1-c} - \left(1+b\left(P_0-P\mathrm{th}\right)\right)^{1-c}\right)}{b(1-c)}\right)$$ $$a=\frac{1+K'_0}{1+K'_0+K_TK''_0}$$ $$b=\frac{K'_0}{K_0}-\frac{K''_0}{1+K'_0}$$ $$c=\frac{1+K'_0+K_0 K''_0}{{K'_0}^2 + K'_0 - K_0 K''_0}$$ $$P\mathrm{th}=\alpha_0K_0\frac{\theta}{\xi_0}\left(\frac{1}{\mathrm{e}^u-1} - \frac{1}{\mathrm{e}^{u_0}-1}\right)$$ $$\theta = \frac{10636}{S/n} + 6.44 \quad\xi=\frac{u^2\mathrm{e}^u}{\left(\mathrm{e}^u-1\right)^2} \quad u=\frac{\theta}{T}$$ $$n=\begin{cases} 7, & \mathrm{magnetite} \\ 5, & \mathrm{hematite} \end{cases}$$ $$\xi_o = \xi(T_0) \quad u_0=u(T_0) \quad T_0 = 298.15\,\mathrm{K}$$

All required values can be obtained from the Holland and Powell (2011) dataset version 6.33 (tc-633.txt).

The valid pressure range of this buffer is limited by the breakdown of magnetite at high pressures: “…magnetite breaks down via a disproportionation reaction to Fe₄O₅ and hematite (Fe₂O₃)…this reaction occurs at ~9.5–11 GPa and 973–1673 K.” (Woodland et al, 2012)

fayalite – magnetite + quartz

$$3\mathrm{Fe_2SiO_4 + O_2 = 3SiO_2 + 2Fe_3O_4}$$

$$\log_{10}f\mathrm{O}_2=\frac{\mu\mathrm{O}_2 + \left(3\left(\int_1^P V^\mathrm{qz}\,dP + G^\mathrm{ex}_{P,T}\right) + 2\int_1^P V^\mathrm{mt}\,dP - 3\int_1^P V^\mathrm{fa}\,dP\right)}{RT\ln10}$$

Parameterisation at 1 bar from O'Neill (1987):$$\mu\mathrm{O}_2=-587474+1584.427T-203.3164T\ln{T}+0.09271T^2$$ Calibrated at 1 bar from 900 K to 1420 K.

Fugacity at pressure is calculated using a modified Tait equation of state including a thermal pressure term as described by Holland and Powell (2011): $$V=V_0\left(1-a\left(1-\left(1+b\left(P-P\mathrm{th}\right)\right)^{-c}\right)\right)$$ $$\int_{P_0}^P V\,dP=V_0\left(\left(P-P_0\right)\left(1-a\right) + \frac{a\left(\left(1+b\left(P-P\mathrm{th}\right)\right)^{1-c} - \left(1+b\left(P_0-P\mathrm{th}\right)\right)^{1-c}\right)}{b(1-c)}\right)$$ $$a=\frac{1+K'_0}{1+K'_0+K_TK''_0}$$ $$b=\frac{K'_0}{K_0}-\frac{K''_0}{1+K'_0}$$ $$c=\frac{1+K'_0+K_0 K''_0}{{K'_0}^2 + K'_0 - K_0 K''_0}$$ $$P\mathrm{th}=\alpha_0K_0\frac{\theta}{\xi_0}\left(\frac{1}{\mathrm{e}^u-1} - \frac{1}{\mathrm{e}^{u_0}-1}\right)$$ $$\theta = \frac{10636}{S/n} + 6.44 \quad\xi=\frac{u^2\mathrm{e}^u}{\left(\mathrm{e}^u-1\right)^2} \quad u=\frac{\theta}{T}$$ $$n=\begin{cases} 7, & \mathrm{fayalite\,and\,magnetite} \\ 3, & \mathrm{quartz} \end{cases}$$ $$\xi_o = \xi(T_0) \quad u_0=u(T_0) \quad T_0 = 298.15\,\mathrm{K}$$

α-quartz to β-quartz transition considered using Landau theory. The transition at 1 bar is included in the expression for \(\mu\mathrm{O}_2\), hence only the pressure dependent Landau term is sought by subtracting the pressure independent Landau contribution from the total using: $$G^\mathrm{ex}_{P,T} = G^\mathrm{ex}_{P,T}(V_\mathrm{max}) - G^\mathrm{ex}_{P,T}(V_\mathrm{max}=0)$$ $$G^\mathrm{ex}_{P,T} = S_\mathrm{max}\left(T_{c0}\left(Q^2+\frac{1}{3}\left(Q^6 - Q^6_0\right)\right) - T_{cP} Q^2 - T \left(Q^2_0 - Q^2 \right) \right)+P V_\mathrm{max} Q_0^2$$ $$T_{cP} = T_{c0}+P\frac{V_\mathrm{max}}{S_\mathrm{max}} \quad Q^2_0 = \sqrt{1-\frac{T_0}{T_{c0}}}$$ $$Q^2=\begin{cases} 0, & \mathrm{if}\,T > T_{cP} \\[2ex] \sqrt{\frac{T_{cP}-T}{T_{c0}}} & \mathrm{otherwise} \end{cases}$$

All required values can be obtained from the Holland and Powell (2011) dataset version 6.33 (tc-633.txt).

Note limited stability range for the actual phase assemblage. Although the fO₂ curve can be extrapolated beyond the stability range and be used as a reference curve (as ΔFMQ), it is best to avoid doing so.

fayalite – ferrosilite + magnetite

Including high pressure phase transitions for Fe-ringwoodite, ferrite-type fayalite, spinel-type ferrosilite, akimotoite, and perovskite-type ferrosilite.

$$\mathrm{6Fe_2SiO_4+O_2=3Fe_2Si_2O_6+2Fe_3O_4}$$

$$\log_{10}f\mathrm{O}_2=\frac{3\Delta G^\mathrm{Fe_2Si_2O_6} + 2\Delta G^\mathrm{Fe_3O_4} - 6\Delta G^\mathrm{Fe_2SiO_4} - \Delta G^\mathrm{O_2}}{RT\ln{10}}$$

$$\Delta G^\mathrm{solids} = \Delta_fH^\circ + \int_{298.15}^T Cp\,dT - T\left(S^\circ + \int_{298.15}^T \frac{Cp}{T}\,dT\right) + G^\mathrm{excess} + \int_1^P V\,dP$$ $$\Delta G^\mathrm{O_2} = H^{\circ}_\mathrm{O_2} - H^{\circ}_\mathrm{298.15,O_2} - TS^{\circ}_\mathrm{O_2}$$

\(\Delta G\) for solids is evaluated using: $$\int Cp\,dT = aT+0.5bT^2-cT^{-1}+2dT^{0.5}$$ $$\int\frac{Cp}{T}\,dT = a\ln{T}+bT-0.5cT^{-2}-2dT^{-0.5}$$ \(\Delta G^\mathrm{O_2}\) is evaluated using a Shomate equation and data from the NIST Chemistry WebBook “Gas phase thermochemistry data” pages. $$H^{\circ} - H^{\circ}_\mathrm{298.15} = AT + \frac{BT^2}{2}+ \frac{CT^3}{3} + \frac{DT^4}{4} - \frac{E}{T} +F - H$$ $$S^{\circ} = A\ln{T} + BT + \frac{CT^2}{2} + \frac{DT^3}{3} - \frac{E}{2T^2} + G$$

Fugacity at pressure is calculated using a modified Tait equation of state including a thermal pressure term as described by Holland and Powell (2011): $$V=V_0\left(1-a\left(1-\left(1+b\left(P-P\mathrm{th}\right)\right)^{-c}\right)\right)$$ $$\int_{P_0}^P V\,dP=V_0\left(\left(P-P_0\right)\left(1-a\right) + \frac{a\left(\left(1+b\left(P-P\mathrm{th}\right)\right)^{1-c} - \left(1+b\left(P_0-P\mathrm{th}\right)\right)^{1-c}\right)}{b(1-c)}\right)$$ $$a=\frac{1+K'_0}{1+K'_0+K_TK''_0}$$ $$b=\frac{K'_0}{K_0}-\frac{K''_0}{1+K'_0}$$ $$c=\frac{1+K'_0+K_0 K''_0}{{K'_0}^2 + K'_0 - K_0 K''_0}$$ $$P\mathrm{th}=\alpha_0K_0\frac{\theta}{\xi_0}\left(\frac{1}{\mathrm{e}^u-1} - \frac{1}{\mathrm{e}^{u_0}-1}\right)$$ $$\theta = \frac{10636}{S/n} + 6.44 \quad\xi=\frac{u^2\mathrm{e}^u}{\left(\mathrm{e}^u-1\right)^2} \quad u=\frac{\theta}{T}$$ $$n=\begin{cases} 7, & \mathrm{Fe_2SiO_4} \\ 10, & \mathrm{Fe_2Si_2O_6} \\ 7, & \mathrm{Fe_3O_4}\end{cases}$$ $$\xi_o = \xi(T_0) \quad u_0=u(T_0) \quad T_0 = 298.15\,\mathrm{K}$$

Landau theory is used to calculate excess Gibbs energy for magnetite phase transitions: $$G^\mathrm{excess} = S_\mathrm{max}\left(T_{c0}\left(Q^2+\frac{1}{3}\left(Q^6 - Q^6_0\right)\right) - T_{cP} Q^2 - T \left(Q^2_0 - Q^2 \right) \right)+P V_\mathrm{max} Q_0^2$$ $$T_{cP} = T_{c0}+P\frac{V_\mathrm{max}}{S_\mathrm{max}} \quad Q^2_0 = \sqrt{1-\frac{T_0}{T_{c0}}}$$ $$Q^2=\begin{cases} 0, & \mathrm{if}\,T > T_{cP} \\[2ex] \sqrt{\frac{T_{cP}-T}{T_{c0}}} & \mathrm{otherwise} \end{cases}$$

All required values can be obtained from the Holland and Powell (2011) dataset version 6.33 (tc-633.txt).

The valid pressure range of this buffer is limited by the breakdown of magnetite at high pressures: “…magnetite breaks down via a disproportionation reaction to Fe₄O₅ and hematite (Fe₂O₃)…this reaction occurs at ~9.5–11 GPa and 973–1673 K.” (Woodland et al, 2012)

Also note complex phase relations in real high pressure assemblages containing FeO and SiO₂, making FFM more of a theoretical exercise in metastable buffer curves.

rhenium – rhenium oxide

$$\mathrm{Re + O_2 = ReO_2}$$

$$\log_{10}f\mathrm{O}_2=\frac{\mu\mathrm{O}_2 + \left(\int_1^P V^\mathrm{ReO_2}\,dP - \int_1^P V^\mathrm{Re}\,dP\right)}{RT\ln10}$$

Parameterisation at 1 bar from Pownceby and O'Neill (1994): $$\mu\mathrm{O}_2=-451020+297.595T-14.6585T\ln{T}$$ Calibrated between 850 K and 1250 K.

Pressure effect is estimated using ambient conditions molar volumes as \(\int_1^P \Delta V \,dP = 0.9917P\), and is unlikely to be reliable at pressures higher than several kbar.

ruthenium – ruthenium oxide

$$\mathrm{Ru + O_2 = RuO_2}$$

Parameterisation at 1 bar from O'Neill and Nell (1997): $$\mu\mathrm{O}_2=-324563+344.151T-22.1155T\ln{T}$$ Calibrated 700 K to 1800 K.

Armstrong et al (2020) measured the properties of ruthernium metal and oxides at high temperature and pressure. For Ru metal, they use a modified Tait equation of state: $$V=V_{0,T}\left(1-a\left(1-\left(1+bP\right)^{-c}\right)\right)$$ $$\int_{P_0}^{P_1} V\,dP=V_{0,T}\left(1-a\right)\left(P_1 - P_0\right)+V_{0,T}\frac{a}{b(1-c)}\left(\left(1+bP_1\right)^{1-c} - \left(1+bP_0\right)^{1-c}\right)$$ $$V_{0,T}=V_0\exp\left(\int_{298.15}^T \alpha_0 + \alpha_1 dT\right) = V_0\exp\left(\alpha_0\left(T-298.15\right) + 0.5\alpha_1\left(T^2 - 298.15^2\right)\right)$$ $$a=\frac{1+K'_T}{1+K'_T+K_TK''_T}$$ $$b=\frac{K'_T}{K_T}-\frac{K''_T}{1+K'_T}$$ $$c=\frac{1+K'_T+K_T K''_T}{{K'_T}^2 + K'_T - K_T K''_T}$$ $$K_T=K_0+\left(\frac{\partial K}{\partial T}\right)_P\left(T-T_0\right) \quad K'_T = 4 \quad K'_T = - \frac{K'_T}{K_T}$$

Ruthenium oxide is described using a modified Tait equation of state including a thermal pressure term: $$V=V_0\left(1-a\left(1-\left(1+b\left(P-P\mathrm{th}\right)\right)^{-c}\right)\right)$$ Assuming \(P_0=0\approx 1\,\mathrm{bar}\), then: $$\int_0^P V\,dP=V_0P\left(1-a + \frac{a}{b(1-c)P}\left(\left(1+b\left(P-P\mathrm{th}\right)\right)^{1-c} - \left(1-bP\mathrm{th}\right)^{1-c}\right)\right)$$ $$a=\frac{1+K'_0}{1+K'_0+K_TK''_0}$$ $$b=\frac{K'_0}{K_0}-\frac{K''_0}{1+K'_0}$$ $$c=\frac{1+K'_0+K_0 K''_0}{{K'_0}^2 + K'_0 - K_0 K''_0}$$ $$K'_0 = 4 \quad K''_0 = - \frac{K'_0}{K_0}$$ $$P\mathrm{th}=\alpha_0K_0\frac{\theta}{\xi_0}\left(\frac{1}{\mathrm{e}^u-1} - \frac{1}{\mathrm{e}^{u_0}-1}\right)$$ $$\theta = \frac{10636}{S/n} + 6.44 \quad\xi=\frac{u^2\mathrm{e}^u}{\left(\mathrm{e}^u-1\right)^2} \quad u=\frac{\theta}{T}\quad n=3$$ $$\xi_o = \xi(T_0) \quad u_0=u(T_0) \quad T_0 = 298.15\,\mathrm{K}$$

Ruthenium oxide undergoes two phase transitions. It is tetragonal at low pressures, orthorhombic at intermediate pressures, and cubic at high pressures. The tetragonal to orthorhombic phase change is modelled using Landau Theory: $$G^\mathrm{ex}_{P,T} = S_\mathrm{max}\left(T_{c0}\left(Q^2+\frac{1}{2}\left(Q^4 - Q^4_0\right)\right) - T_{cP} Q^2 - T \left(Q^2_0 - Q^2 \right) \right)+P V_\mathrm{max} Q_0^2$$ $$T_{cP} = T_{c0}+P\frac{V_\mathrm{max}}{S_\mathrm{max}} \quad Q^2_0 = 1-\frac{T_0}{T_{c0}}$$ $$Q^2=\begin{cases} 0, & \mathrm{if}\,T > T_{cP} \\[2ex] \frac{T_{cP}-T}{T_{cP}} & \mathrm{otherwise} \end{cases}$$ The \(G^\mathrm{ex}_{P,T}\) term is added to the Landau-independent Gibbs energy of the tetragonal/orthogonal phase: \(G^\mathrm{total}_{P,T} = G_{P,T} + G^\mathrm{ex}_{P,T}\), where: $$G_{P,T} = \Delta_fH + \int_{T_0}^{T} Cp dT - T\left(S + \int_{T_0}^{T} \frac{Cp}{T} dT\right) $$ $$Cp = a+bT+cT^{-2}+dT^{-0.5}$$ Once the phase transition at a given \(T\) is found by solving for \(P\) which satisfies \(G_{P,T}^\mathrm{tet-orth} = G_{P,T}^\mathrm{cubic}\), oxygen fugacity is calculated by:. $$\log_{10}f\mathrm{O}_2=\frac{\mu\mathrm{O}_2 + \left(\int_1^P V^\mathrm{RuO_2}\,dP - \int_1^P V^\mathrm{Ru}\,dP\right)}{RT\ln10}$$ If the maximum pressure is below the transition pressure, then: $$\int_1^P V^\mathrm{RuO_2}\,dP = \int_1^P V^\mathrm{RuO_2(tet-orth)}\,dP + G^\mathrm{ex}_{P,T}$$ If the maximum pressure is above the transition pressure, then: $$\int_1^P V^\mathrm{RuO_2}\,dP = \int_1^{P_\mathrm{trans}} V^\mathrm{RuO_2(tet-orth)}\,dP + \int_{P_\mathrm{trans}}^P V^\mathrm{RuO_2(cubic)}\,dP + G^\mathrm{ex}_{P_{trans},T} = $$ $$\int_1^{P_\mathrm{trans}} V^\mathrm{RuO_2(tet-orth)}\,dP + \int_1^P V^\mathrm{RuO_2(cubic)}\,dP - \int_1^{P_\mathrm{trans}} V^\mathrm{RuO_2(cubic)}\,dP + G^\mathrm{ex}_{P_{trans},T}$$

All required parameters are given in Tables 5 and 6 of Armstrong (2020). Calibrated to 25 GPa between 773 and 2500 K.

graphite/diamond – CO fluid

$$\mathrm{C + O_2 = CO_2 (\pm CO)}$$

$$\log_{10}f\mathrm{O}_2=\frac{\Delta G^\circ_\mathrm{CO_2} - \left(\Delta G^\circ_\mathrm{C} + \int_1^P V_\mathrm{C}\,dP\right) - \Delta G^\circ_\mathrm{O_2} }{RT\ln{10}} + \frac{\ln f\mathrm{CO}_2}{\ln{10}}$$

Oxygen fugacity values for the CCO buffer are strongly dependent on the fugacity of carbon dioxide (\(f\mathrm{CO}_2\)). There is lack of consensus between experimental studies and molecular dynamics simulations. Simulations results in consistently higher \(f\mathrm{CO}_2\) than is observed in experimental studies, which is then propagated to higher \(f\mathrm{O}_2\).

Experimental work from Frost and Wood (1997) which also summarises previous work shows that the CCO buffer follows the same slope of the FMQ buffer in (\(f\mathrm{O}_2-P\)) space, with (\(f\mathrm{O}_2\)) of the CCO buffer being roughly equal to that of FMQ at low temperatures (~800 °C), and becoming more reduced relative to FMQ at higher temperatures (ΔFMQ=–2 at ~1400 °C and 1 GPa).

In contrast, molecular dynamics simulations show a steeper slope, resulting in higher \(f\mathrm{O}_2\) values, especially at lower temperatures. For example, the simulated data of Li et al (2021) is 1 \(\log_{10}f\mathrm{O}_2\) higher than that obtained using the experimental data of Frost and Wood (1997) at 800 °C and 3 GPa. Other studies (e.g., Fu et al (2017) and Duan and Zhang (2006)) have a similar slope, although the discrepancy is of a slightly lower magnitude. The different models converge at higher temperatures, with all experimental and simulation studies agreeing to within 0.5 \(\log_{10}f\mathrm{O}_2\) unit at at 1400 °C and 3 GPa, but then deviating again at higher pressures.

As the reasons for the discrepancies are speculative without any additional studies, there is little justification to pick one model over the other. Here, two models are considered: the experimental study of Frost and Wood (1997) and the moluecular dynamics simulation of Li et al (2021). The true value is most likely lies somewhere between the two curves.

Gas species in the standard state are evaluated using a Shomate equation: $$\Delta G^\circ_\mathrm{CO_2,\,O_2} = H^{\circ} - H^{\circ}_\mathrm{298.15} - TS^{\circ}$$ $$H^{\circ} - H^{\circ}_\mathrm{298.15} = AT + \frac{BT^2}{2}+ \frac{CT^3}{3} + \frac{DT^4}{4} - \frac{E}{T} +F - H$$ $$S^{\circ} = A\ln{T} + BT + \frac{CT^2}{2} + \frac{DT^3}{3} - \frac{E}{2T^2} + G$$ All parameters are given in the NIST Chemistry WebBook “Gas phase thermochemistry data” pages for carbon dioxide, and oxygen.

Elemental carbon requires two terms: standard state \(\Delta G^\circ_\mathrm{C}\) and high pressure \(\int_1^P V_\mathrm{C}\,dP\). Equations are described in Holland and Powell (2011) and parameters are available in their dataset version 6.33 (tc-633.txt).

$$\Delta G^\circ_\mathrm{C} = \Delta_fH^\circ + \int_{298.15}^T Cp\,dT - T\left(S^\circ + \int_{298.15}^T \frac{Cp}{T}\,dT\right) + \int_1^P V\,dP$$ $$\int Cp\,dT = aT+0.5bT^2-cT^{-1}+2dT^{0.5}$$ $$\int\frac{Cp}{T}\,dT = a\ln{T}+bT-0.5cT^{-2}-2dT^{-0.5}$$

\(\int_1^P V_\mathrm{C}\,dP\) is calculated using a modified Tait equation of state including a thermal pressure term as described by : $$V=V_0\left(1-a\left(1-\left(1+b\left(P-P\mathrm{th}\right)\right)^{-c}\right)\right)$$ $$\int_{P_0}^P V\,dP=V_0\left(\left(P-P_0\right)\left(1-a\right) + \frac{a\left(\left(1+b\left(P-P\mathrm{th}\right)\right)^{1-c} - \left(1+b\left(P_0-P\mathrm{th}\right)\right)^{1-c}\right)}{b(1-c)}\right)$$ $$a=\frac{1+K'_0}{1+K'_0+K_TK''_0}$$ $$b=\frac{K'_0}{K_0}-\frac{K''_0}{1+K'_0}$$ $$c=\frac{1+K'_0+K_0 K''_0}{{K'_0}^2 + K'_0 - K_0 K''_0}$$ $$P\mathrm{th}=\alpha_0K_0\frac{\theta}{\xi_0}\left(\frac{1}{\mathrm{e}^u-1} - \frac{1}{\mathrm{e}^{u_0}-1}\right)$$ $$\theta = \frac{10636}{S/n} + 6.44 \quad\xi=\frac{u^2\mathrm{e}^u}{\left(\mathrm{e}^u-1\right)^2} \quad u=\frac{\theta}{T}$$ $$\xi_o = \xi(T_0) \quad u_0=u(T_0) \quad T_0 = 298.15\,\mathrm{K}\quad n=1$$

\(\Delta G^\circ_\mathrm{C} + \int_1^P V_\mathrm{C}\,dP\) is evaluated twice for the two carbon phases (graphite and diamond). The lower value for a given \(P\) and \(T\) is used in the expression for \(\log_{10}f\mathrm{O}_2\) given above.

Experimental carbon dioxide fugacity

Frost and Wood (1997) parameterise oxygen fugacity, with the equation provided in an erratum: $$\log_{10}f\mathrm{O}_2=4.3927-\frac{21234-117.9P+0.891T^2}{T}+7.679\times10^{-3}P+3.627\times10^{-4}P^2$$ This parameterisation is strongly overestimating \(f\mathrm{O}_2\) at low pressure (< 0.5 GPa), and underestimating \(f\mathrm{O}_2\) at high pressures (> 10 GPa). Instead, we use the parameters fitted by Frost and Wood (1997) an expanded Redlich-Kwong equation of state: $$V=\frac{RT}{P}+b-\frac{aR\sqrt{T}}{(RT+bP)(RT+2bP)}+c\sqrt{P}+dP$$ This equation is integrated to obtain: $$\ln f\mathrm{CO}_2 = \ln{P} + \frac{bP}{RT} + \frac{a\left(\ln(RT+bP)-\ln(RT+2bP)\right)}{RTb\sqrt{T}} + \frac{2cP\sqrt{P}}{3RT} + \frac{dP^2}{2RT} $$ This equations results in the correct \(f\mathrm{O}_2\) values for low pressures, and a slope consisent with simulated data at high pressures, up to about 20 GPa.

Simulated carbon dioxide fugacity

Calculated \(f\mathrm{CO}_2\) using equations of state and parameters derived from molecular dynamics simulations tend to be higher than those obtained experimentally. Here, the study of Li et al (2021) is used, as it results in the highest \(f\mathrm{O}_2\) values and can therefore be used as an upper bound for the range of simulated \(f\mathrm{CO}_2\). They used a modified Lee and Kesler equation of state: $$P = \frac{RT}{V}\left( 1+\frac{B}{V_r}+\frac{C}{V_r^2}+\frac{D}{V_r^4} + \frac{E}{V_r^5} + \frac{F}{V_r^2}\left(\beta+\frac{\gamma}{V_r^2}\right)\exp\left(-\frac{\gamma}{V_r^2}\right)\right)$$

Direct calculation of \(\int_1^P V\,dP\) requires \(V=f(P)\) which is not available, so instead it is calculated using integration by parts: $$\ln f\mathrm{CO}_2 = \int_{P_0}^{P_1} V\,dP = \Delta PV - \int_{V_0}^{V_1} P\,dV=P(V_1)V_1 - P(V_0)V_0 - \int_{V_0}^{V_1} P\,dV$$ \(\int_{V_0}^{V_1} P\,dV\) is integrated numerically, and values for \(V_n\) are found numerically by solving for \(V\) which satisfy \(P(V)-P=0\).

forsterite + graphite/diamond – magnesite + enstatite

$$\mathrm{Mg_2SiO_4+C+O_2=MgCO_3+MgSiO_3}$$

forsterite + diopside + graphite/diamond – dolomite + enstatite

$$\mathrm{Mg_2SiO_4+\frac{1}{2}CaMgSi_2O_6+C+O_2=\frac{1}{2}CaMg(CO_3)_2+2MgSiO_3}$$

diopside + graphite/diamond – coesite + dolomite

$$\mathrm{\frac{1}{2}CaMgSi_2O_6 +C+O_2=\frac{1}{2}CaMg(CO_3)_2 + SiO_2}$$

tungsten – tungsten(IV) oxide

$$\mathrm{W + O_2 = WO_2}$$

$$\log_{10}f\mathrm{O}_2=\frac{\mu\mathrm{O}_2 + \left(\int_1^P V^\mathrm{WO_2}\,dP - \int_1^P V^\mathrm{W}\,dP\right)}{RT\ln10}$$

Parameterisation at 1 bar from O'Neill and Pownceby (1993): $$\mu\mathrm{O}_2=-569087+300.479T-15.9697T\ln{T}$$ Calibrated from 700 K to 1700 K.

Shofner et al (2016) use a Mie-Grüneisen–Birch-Murnaghan equation of state: $$P=\frac{3}{2}K_0\left(\left(\frac{V_0}{V}\right)^{7/3}-\left(\frac{V_0}{V}\right)^{5/3}\right)\left(1-\frac{3}{4}\left(4-K'_0\right)\left(\left(\frac{V_0}{V}\right)^{2/3}-1\right)\right)+\gamma_0\frac{V^{q-1}}{V_0^q}\left(E(\theta_D,T) - E(\theta_D,295)\right)$$ \(E\left(\theta_D,T\right)=3nRTD_3\left(\frac{\theta_D}{T}\right)\) where \(D_3()\) is the Debye function, and \(n=\begin{cases} 1, & \mathrm{W} \\ 3, & \mathrm{WO_2} \end{cases}\).

Direct calculation of \(\int_1^P V\,dP\) requires \(V=f(P)\) which is not available, so instead it is calculated using integration by parts: $$\int_{P_0}^{P_1} V\,dP = \Delta PV - \int_{V_0}^{V_1} P\,dV=P(V_1)V_1 - P(V_0)V_0 - \int_{V_0}^{V_1} P\,dV$$ $$\int P\,dV=\left(E(\theta_D,T) - E(\theta_D,295)\right)\frac{\gamma_0}{q}\left(\frac{V}{V_0}\right)^q-\frac{9K_0\left(\left(K_0'-4\right)V_0^3+\left(14-3K_0'\right)V^{2/3}V_0^{7/3}+\left(3K_0'-16\right)V^{4/3}V_0^{5/3}\right)}{16V^2}$$ Values for \(V_n\) are found numerically by solving for \(V\) which satisfy \(P(V)-P=0\).

\(\mathrm{WO_2}\) undergoes two phase transitions at 4 GPa and 32 GPa, and are independent of temperature in the Shofner et al (2016) model. Thus in cases where the required pressure is greater than 4 GPa or 32 GPa, the \(\int V\,dP\) integral is broken into two or three integrals, respectively, e.g.: $$\int_{P_0}^{P_1} V^{\mathrm{WO_2}}\,dP = \int_{P_0}^{4\,\mathrm{GPa}} V^{\mathrm{ap}}\,dP + \int_{4\,\mathrm{GPa}}^{32\,\mathrm{GPa}} V^{\mathrm{hp}}\,dP + \int_{32\,\mathrm{GPa}}^{P_1} V^{\mathrm{hpm}}\,dP$$

All required parameters are given in Table 1 of Shofner et al (2016) and are calibrated to about 70 GPa and 2300 K.

copper – copper(I) oxide (cuprite)

$$\mathrm{4Cu + O_2 = 2Cu_2O}$$

$$\log_{10}f\mathrm{O}_2=\frac{\mu\mathrm{O}_2 + 2\int_1^P V^\mathrm{cuprite}\,dP - 4\int_1^P V^\mathrm{copper}\,dP}{RT\ln10}$$

Parameterisation at 1 bar from O'Neill (1988): $$\mu\mathrm{O}_2= -347705+246.096T - 12.9053T\ln{T}$$ Calibrated at 1 bar from 750 K to 1330 K.

Fugacity at pressure is calculated using a modified Tait equation of state including a thermal pressure term as described by Holland and Powell (2011): $$V=V_0\left(1-a\left(1-\left(1+b\left(P-P\mathrm{th}\right)\right)^{-c}\right)\right)$$ $$\int_{P_0}^P V\,dP=V_0\left(\left(P-P_0\right)\left(1-a\right) + \frac{a\left(\left(1+b\left(P-P\mathrm{th}\right)\right)^{1-c} - \left(1+b\left(P_0-P\mathrm{th}\right)\right)^{1-c}\right)}{b(1-c)}\right)$$ $$a=\frac{1+K'_0}{1+K'_0+K_TK''_0}$$ $$b=\frac{K'_0}{K_0}-\frac{K''_0}{1+K'_0}$$ $$c=\frac{1+K'_0+K_0 K''_0}{{K'_0}^2 + K'_0 - K_0 K''_0}$$ $$P\mathrm{th}=\alpha_0K_0\frac{\theta}{\xi_0}\left(\frac{1}{\mathrm{e}^u-1} - \frac{1}{\mathrm{e}^{u_0}-1}\right)$$ $$\theta = \frac{10636}{S/n} + 6.44 \quad\xi=\frac{u^2\mathrm{e}^u}{\left(\mathrm{e}^u-1\right)^2} \quad u=\frac{\theta}{T}$$ $$n=\begin{cases} 1, & \mathrm{copper} \\ 3, & \mathrm{cuprite} \end{cases}$$ $$\xi_o = \xi(T_0) \quad u_0=u(T_0) \quad T_0 = 298.15\,\mathrm{K}$$

All required values can be obtained from the Holland and Powell (2011) dataset version 6.33 (tc-633.txt).

copper(I) oxide – copper(II) oxide
(cuprite – tenorite)

$$\mathrm{2Cu_2O + O_2 = 4CuO}$$

$$\log_{10}f\mathrm{O}_2=\frac{\mu\mathrm{O}_2 + 4\int_1^P V^\mathrm{tenorite}\,dP - 2\int_1^P V^\mathrm{cuprite}\,dP}{RT\ln10}$$

Parameterisation at 1 bar from O'Neill (1988): $$\mu\mathrm{O}_2= -292245+377.012T - 23.1976T\ln{T}$$ Calibrated at 1 bar from 800 K to 1300 K.

Fugacity at pressure is calculated using a modified Tait equation of state including a thermal pressure term as described by Holland and Powell (2011): $$V=V_0\left(1-a\left(1-\left(1+b\left(P-P\mathrm{th}\right)\right)^{-c}\right)\right)$$ $$\int_{P_0}^P V\,dP=V_0\left(\left(P-P_0\right)\left(1-a\right) + \frac{a\left(\left(1+b\left(P-P\mathrm{th}\right)\right)^{1-c} - \left(1+b\left(P_0-P\mathrm{th}\right)\right)^{1-c}\right)}{b(1-c)}\right)$$ $$a=\frac{1+K'_0}{1+K'_0+K_TK''_0}$$ $$b=\frac{K'_0}{K_0}-\frac{K''_0}{1+K'_0}$$ $$c=\frac{1+K'_0+K_0 K''_0}{{K'_0}^2 + K'_0 - K_0 K''_0}$$ $$P\mathrm{th}=\alpha_0K_0\frac{\theta}{\xi_0}\left(\frac{1}{\mathrm{e}^u-1} - \frac{1}{\mathrm{e}^{u_0}-1}\right)$$ $$\theta = \frac{10636}{S/n} + 6.44 \quad\xi=\frac{u^2\mathrm{e}^u}{\left(\mathrm{e}^u-1\right)^2} \quad u=\frac{\theta}{T}$$ $$n=\begin{cases} 3, & \mathrm{cuprite} \\ 2, & \mathrm{tenorite} \end{cases}$$ $$\xi_o = \xi(T_0) \quad u_0=u(T_0) \quad T_0 = 298.15\,\mathrm{K}$$

All required values can be obtained from the Holland and Powell (2011) dataset version 6.33 (tc-633.txt).

iron(II) – iron(II,III) oxide (wüstite – magnetite)

$$\mathrm{6FeO + O_2 = 2Fe_3O_4}$$

$$\log_{10}f\mathrm{O}_2=\frac{\mu\mathrm{O}_2 + 2\int_1^P V^\mathrm{Fe_3O_4}\,dP - 6\int_1^P V^\mathrm{FeO}\,dP}{RT\ln10}$$

Parameterisation at 1 bar from O'Neill (1988): $$\mu\mathrm{O}_2= -581927+65.618T + 38.741T\ln{T}$$ Calibrated at 1 bar from 833 K to 1270 K.

Fugacity at pressure is calculated using a modified Tait equation of state including a thermal pressure term as described by Holland and Powell (2011): $$V=V_0\left(1-a\left(1-\left(1+b\left(P-P\mathrm{th}\right)\right)^{-c}\right)\right)$$ $$\int_{P_0}^P V\,dP=V_0\left(\left(P-P_0\right)\left(1-a\right) + \frac{a\left(\left(1+b\left(P-P\mathrm{th}\right)\right)^{1-c} - \left(1+b\left(P_0-P\mathrm{th}\right)\right)^{1-c}\right)}{b(1-c)}\right)$$ $$a=\frac{1+K'_0}{1+K'_0+K_TK''_0}$$ $$b=\frac{K'_0}{K_0}-\frac{K''_0}{1+K'_0}$$ $$c=\frac{1+K'_0+K_0 K''_0}{{K'_0}^2 + K'_0 - K_0 K''_0}$$ $$P\mathrm{th}=\alpha_0K_0\frac{\theta}{\xi_0}\left(\frac{1}{\mathrm{e}^u-1} - \frac{1}{\mathrm{e}^{u_0}-1}\right)$$ $$\theta = \frac{10636}{S/n} + 6.44 \quad\xi=\frac{u^2\mathrm{e}^u}{\left(\mathrm{e}^u-1\right)^2} \quad u=\frac{\theta}{T}$$ $$n=\begin{cases} 7, & \mathrm{magnetite} \\ 2, & \mathrm{wustite} \end{cases}$$ $$\xi_o = \xi(T_0) \quad u_0=u(T_0) \quad T_0 = 298.15\,\mathrm{K}$$

All required values can be obtained from the Holland and Powell (2011) dataset version 6.33 (tc-633.txt), using the ferropericlase (“fper”) end-member for wüstite.

The valid pressure range of this buffer is limited by the breakdown of magnetite at high pressures: “…magnetite breaks down via a disproportionation reaction to Fe₄O₅ and hematite (Fe₂O₃)…this reaction occurs at ~9.5–11 GPa and 973–1673 K.” (Woodland et al, 2012)

Note relationship between three reactions in the Fe-O system. WM is only stable at temperature of about 830 K and lower. The invariant intersection point with the other buffers (IM and WM) shifts to lower temperature and higher oxygen fugacity with increasing pressure. See schematic diagram:

iron – iron(II,III) oxide (magnetite)

$$\mathrm{\frac{3}{2}Fe + O_2 = \frac{1}{2}Fe_3O_4}$$

$$\log_{10}f\mathrm{O}_2=\frac{\mu\mathrm{O}_2 + 0.5\int_1^P V^\mathrm{Fe_3O_4}\,dP - 1.5\int_1^P V^\mathrm{Fe}\,dP}{RT\ln10}$$

Parameterisation at 1 bar from O'Neill (1988): $$\mu\mathrm{O}_2= -607673+1060.994T - 132.3909T\ln{T}+0.06657T^2$$ Calibrated at 1 bar from 750 K to 833 K.

For Fe, the equation of state described by Dorogokupets et al (2017) is used. The equation is given in terms of the Helmholz free energy \(F\). Expressions for \(F(V,T)\) and \(P(V)\) are given in the Supplementary Files of Dorogokupets et al (2017). \(V(P,T)\) is calculated by numerically solving \(P(V)-P=0\) for the pressure of interest. \(G(P,T)\) is calculated using: $$G = F+PV$$ This Gibbs energy is used to find phase transitions, and the resulting pressures are used as the boundaries for integration. The \(\int_1^P V^\mathrm{Fe}\,dP\) integral is evaluated numerically.

This equation of state is calibrated to >300 GPa and >3000 K. The BCC α-Fe and γ-Fe phases are stable at low pressures. At high pressures, the FCC γ-Fe exists at high temperatures whereas the HCP ε-Fe is stable at lower temperatures. Evaluation of \(G\) leads to prediction of FCC γ-Fe at extreme pressures (>100–150 GPa). This phase transition is spurious and results from extrapolation beyond the calibrated range, and as such is not considered here.

For magnetite, fugacity at pressure is calculated using a modified Tait equation of state including a thermal pressure term as described by Holland and Powell (2011): $$V=V_0\left(1-a\left(1-\left(1+b\left(P-P\mathrm{th}\right)\right)^{-c}\right)\right)$$ $$\int_{P_0}^P V\,dP=V_0\left(\left(P-P_0\right)\left(1-a\right) + \frac{a\left(\left(1+b\left(P-P\mathrm{th}\right)\right)^{1-c} - \left(1+b\left(P_0-P\mathrm{th}\right)\right)^{1-c}\right)}{b(1-c)}\right)$$ $$a=\frac{1+K'_0}{1+K'_0+K_TK''_0}$$ $$b=\frac{K'_0}{K_0}-\frac{K''_0}{1+K'_0}$$ $$c=\frac{1+K'_0+K_0 K''_0}{{K'_0}^2 + K'_0 - K_0 K''_0}$$ $$P\mathrm{th}=\alpha_0K_0\frac{\theta}{\xi_0}\left(\frac{1}{\mathrm{e}^u-1} - \frac{1}{\mathrm{e}^{u_0}-1}\right)$$ $$\theta = \frac{10636}{S/n} + 6.44 \quad\xi=\frac{u^2\mathrm{e}^u}{\left(\mathrm{e}^u-1\right)^2} \quad u=\frac{\theta}{T}$$ $$\xi_o = \xi(T_0) \quad u_0=u(T_0) \quad T_0 = 298.15\,\mathrm{K}\quad n=7$$

All required values can be obtained from the Holland and Powell (2011) dataset version 6.33 (tc-633.txt).

The valid pressure range of this buffer is limited by the breakdown of magnetite at high pressures: “…magnetite breaks down via a disproportionation reaction to Fe₄O₅ and hematite (Fe₂O₃)…this reaction occurs at ~9.5–11 GPa and 973–1673 K.” (Woodland et al, 2012)

Note relationship between three reactions in the Fe-O system. IM is only stable at temperature of about 830 K and lower. The invariant intersection point with the other buffers (IM and WM) shifts to lower temperature and higher oxygen fugacity with increasing pressure. See schematic diagram:

tungsten carbide – tungsten oxide

$$\mathrm{WC + O_2 = WO_2 + C}$$

Calibration from Taylor and Foley (1989): $$\log_{10}f\mathrm{O}_2=11.82-\frac{33780}{T}+0.066\frac{P}{T}$$ Temperature in kelvin and pressure in bar.

Valid between 1000 °C and 1250 °C, and 9 kbar 35 kbar.

iron + quartz – fayalite

Including high pressure phase transitions for coesite, stishovite, and Fe-ringwoodite.

$$\mathrm{2Fe + SiO_2 + O_2 = Fe_2SiO_4}$$

$$\log_{10}f\mathrm{O}_2=\frac{\mu\mathrm{O}_2 + \int_1^P V^\mathrm{fayalite}\,dP - \left(\int_1^P V^\mathrm{quartz}\,dP + G^\mathrm{ex}_{P,T}\right) - 2\int_1^P V^\mathrm{Fe}\,dP}{RT\ln10}$$

Parameterisation at 1 bar from O'Neill (1987):$$\mu\mathrm{O}_2=\begin{cases} -542941-33.182T+22.446T\ln{T}, & 900 < T < 1042 \\[2ex] -562377+103.384T+5.4771T\ln{T}, & 1042 \leq T \leq 1184 \\[2ex] -602739+369.704T-27.3443T\ln{T}, & 1184 < T < 1420\end{cases}$$ Calibrated at 1 bar from 900 K to 1420 K.

\(\int_1^P V\,dP\) for quartz and fayalite is calculated using a modified Tait equation of state including a thermal pressure term as described by Holland and Powell (2011): $$V=V_0\left(1-a\left(1-\left(1+b\left(P-P\mathrm{th}\right)\right)^{-c}\right)\right)$$ $$\int_{P_0}^P V\,dP=V_0\left(\left(P-P_0\right)\left(1-a\right) + \frac{a\left(\left(1+b\left(P-P\mathrm{th}\right)\right)^{1-c} - \left(1+b\left(P_0-P\mathrm{th}\right)\right)^{1-c}\right)}{b(1-c)}\right)$$ $$a=\frac{1+K'_0}{1+K'_0+K_TK''_0}$$ $$b=\frac{K'_0}{K_0}-\frac{K''_0}{1+K'_0}$$ $$c=\frac{1+K'_0+K_0 K''_0}{{K'_0}^2 + K'_0 - K_0 K''_0}$$ $$P\mathrm{th}=\alpha_0K_0\frac{\theta}{\xi_0}\left(\frac{1}{\mathrm{e}^u-1} - \frac{1}{\mathrm{e}^{u_0}-1}\right)$$ $$\theta = \frac{10636}{S/n} + 6.44 \quad\xi=\frac{u^2\mathrm{e}^u}{\left(\mathrm{e}^u-1\right)^2} \quad u=\frac{\theta}{T}$$ $$n=\begin{cases} 7, & \mathrm{fayalite} \\ 3, & \mathrm{quartz} \end{cases}$$ $$\xi_o = \xi(T_0) \quad u_0=u(T_0) \quad T_0 = 298.15\,\mathrm{K}$$

α-quartz to β-quartz transition considered using Landau theory. The transition at 1 bar is included in the expression for \(\mu\mathrm{O}_2\), hence only the pressure dependent Landau term is sought by subtracting the pressure independent Landau contribution from the total using: $$G^\mathrm{ex}_{P,T} = G^\mathrm{ex}_{P,T}(V_\mathrm{max}) - G^\mathrm{ex}_{P,T}(V_\mathrm{max}=0)$$ $$G^\mathrm{ex}_{P,T} = S_\mathrm{max}\left(T_{c0}\left(Q^2+\frac{1}{3}\left(Q^6 - Q^6_0\right)\right) - T_{cP} Q^2 - T \left(Q^2_0 - Q^2 \right) \right)+P V_\mathrm{max} Q_0^2$$ $$T_{cP} = T_{c0}+P\frac{V_\mathrm{max}}{S_\mathrm{max}} \quad Q^2_0 = \sqrt{1-\frac{T_0}{T_{c0}}}$$ $$Q^2=\begin{cases} 0, & \mathrm{if}\,T > T_{cP} \\[2ex] \sqrt{\frac{T_{cP}-T}{T_{c0}}} & \mathrm{otherwise} \end{cases}$$

All required values for quartz and fayalite can be obtained from the Holland and Powell (2011) dataset version 6.33 (tc-633.txt).

\(\int_1^P V\,dP\) for Fe is calculated using the equation of state described by Dorogokupets et al (2017) is used. The equation is given in terms of the Helmholz free energy \(F\). Expressions for \(F(V,T)\) and \(P(V)\) are given in the Supplementary Files of Dorogokupets et al (2017). \(V(P,T)\) is calculated by numerically solving \(P(V)-P=0\) for the pressure of interest. \(G(P,T)\) is calculated using: $$G = F+PV$$ This Gibbs energy is used to find phase transitions, and the resulting pressures are used as the boundaries for integration. The \(\int_1^P V^\mathrm{Fe}\,dP\) integral is evaluated numerically.

This equation of state is calibrated to >300 GPa and >3000 K. The BCC α-Fe and γ-Fe phases are stable at low pressures. At high pressures, the FCC γ-Fe exists at high temperatures whereas the HCP ε-Fe is stable at lower temperatures. Evaluation of \(G\) leads to prediction of FCC γ-Fe at extreme pressures (>100–150 GPa). This phase transition is spurious and results from extrapolation beyond the calibrated range, and as such is not considered here.

iron + rutile – ilmenite

$$\mathrm{2Fe + 2TiO_2 + O_2 = 2FeTiO_3}$$

$$\log_{10}f\mathrm{O}_2=\frac{\mu\mathrm{O}_2 + 2\left(\int_1^P V^\mathrm{ilmenite}\,dP + G^\mathrm{ex}_{P,T}\right) - 2\int_1^P V^\mathrm{rutile}\,dP - 2\int_1^P V^\mathrm{Fe}\,dP}{RT\ln10}$$

Parameterisation at 1 bar from O'Neill et al (1988):$$\mu\mathrm{O}_2=\begin{cases} -570745+97.278T+5.0278T\ln{T}, & 800 < T < 1042 \\[2ex] -578438+139.612T, & 1042 \leq T \leq 1184 \\[2ex] -574852+136.606T, & 1184 < T < 1420\end{cases}$$ Calibrated at 1 bar from 800 K to 1340 K. Note typographical error in the original paper–missing decimal point in the first expression.

\(\int_1^P V\,dP\) for rutile and ilmenite is calculated using a modified Tait equation of state including a thermal pressure term as described by Holland and Powell (2011): $$V=V_0\left(1-a\left(1-\left(1+b\left(P-P\mathrm{th}\right)\right)^{-c}\right)\right)$$ $$\int_{P_0}^P V\,dP=V_0\left(\left(P-P_0\right)\left(1-a\right) + \frac{a\left(\left(1+b\left(P-P\mathrm{th}\right)\right)^{1-c} - \left(1+b\left(P_0-P\mathrm{th}\right)\right)^{1-c}\right)}{b(1-c)}\right)$$ $$a=\frac{1+K'_0}{1+K'_0+K_TK''_0}$$ $$b=\frac{K'_0}{K_0}-\frac{K''_0}{1+K'_0}$$ $$c=\frac{1+K'_0+K_0 K''_0}{{K'_0}^2 + K'_0 - K_0 K''_0}$$ $$P\mathrm{th}=\alpha_0K_0\frac{\theta}{\xi_0}\left(\frac{1}{\mathrm{e}^u-1} - \frac{1}{\mathrm{e}^{u_0}-1}\right)$$ $$\theta = \frac{10636}{S/n} + 6.44 \quad\xi=\frac{u^2\mathrm{e}^u}{\left(\mathrm{e}^u-1\right)^2} \quad u=\frac{\theta}{T}$$ $$n=\begin{cases} 3, & \mathrm{rutile} \\ 5, & \mathrm{ilmenite} \end{cases}$$ $$\xi_o = \xi(T_0) \quad u_0=u(T_0) \quad T_0 = 298.15\,\mathrm{K}$$

Ilmenite phase transitions considered using Landau theory. The transition at 1 bar is included in the expression for \(\mu\mathrm{O}_2\), hence only the pressure dependent Landau term is sought by subtracting the pressure independent Landau contribution from the total using: $$G^\mathrm{ex}_{P,T} = G^\mathrm{ex}_{P,T}(V_\mathrm{max}) - G^\mathrm{ex}_{P,T}(V_\mathrm{max}=0)$$ $$G^\mathrm{ex}_{P,T} = S_\mathrm{max}\left(T_{c0}\left(Q^2+\frac{1}{3}\left(Q^6 - Q^6_0\right)\right) - T_{cP} Q^2 - T \left(Q^2_0 - Q^2 \right) \right)+P V_\mathrm{max} Q_0^2$$ $$T_{cP} = T_{c0}+P\frac{V_\mathrm{max}}{S_\mathrm{max}} \quad Q^2_0 = \sqrt{1-\frac{T_0}{T_{c0}}}$$ $$Q^2=\begin{cases} 0, & \mathrm{if}\,T > T_{cP} \\[2ex] \sqrt{\frac{T_{cP}-T}{T_{c0}}} & \mathrm{otherwise} \end{cases}$$

All required values for rutile and ilmenite can be obtained from the Holland and Powell (2011) dataset version 6.33 (tc-633.txt).

\(\int_1^P V\,dP\) for Fe is calculated using the equation of state described by Dorogokupets et al (2017) is used. The equation is given in terms of the Helmholz free energy \(F\). Expressions for \(F(V,T)\) and \(P(V)\) are given in the Supplementary Files of Dorogokupets et al (2017). \(V(P,T)\) is calculated by numerically solving \(P(V)-P=0\) for the pressure of interest. \(G(P,T)\) is calculated using: $$G = F+PV$$ This Gibbs energy is used to find phase transitions, and the resulting pressures are used as the boundaries for integration. The \(\int_1^P V^\mathrm{Fe}\,dP\) integral is evaluated numerically.

This equation of state is calibrated to >300 GPa and >3000 K. The BCC α-Fe and γ-Fe phases are stable at low pressures. At high pressures, the FCC γ-Fe exists at high temperatures whereas the HCP ε-Fe is stable at lower temperatures. Evaluation of \(G\) leads to prediction of FCC γ-Fe at extreme pressures (>100–150 GPa). This phase transition is spurious and results from extrapolation beyond the calibrated range, and as such is not considered here.

iron + ilmenite – ulvöspinel

$$\mathrm{2Fe + 2FeTiO_3 + O_2 = 2Fe_2TiO_4}$$

$$\log_{10}f\mathrm{O}_2=\frac{\mu\mathrm{O}_2 + 2\int_1^P V^\mathrm{usp}\,dP - 2\left(\int_1^P V^\mathrm{ilm}\,dP + G^\mathrm{ex}_{P,T}\right) - 2\int_1^P V^\mathrm{Fe}\,dP}{RT\ln10}$$

Parameterisation at 1 bar from O'Neill et al (1988):$$\mu\mathrm{O}_2=\begin{cases} -505563+4.99T+13.3151T\ln{T}, & 900 < T < 1042 \\[2ex] -521274+112.599T, & 1042 \leq T \leq 1184 \\[2ex] -587764+553.582T-54.37T\ln{T}, & 1184 < T < 1500\end{cases}$$ Calibrated at 1 bar from 800 K to 1340 K. Note typographical error in the original paper–missing decimal point in the first expression.

\(\int_1^P V\,dP\) for ulv&ouml;spinel and ilmenite is calculated using a modified Tait equation of state including a thermal pressure term as described by Holland and Powell (2011): $$V=V_0\left(1-a\left(1-\left(1+b\left(P-P\mathrm{th}\right)\right)^{-c}\right)\right)$$ $$\int_{P_0}^P V\,dP=V_0\left(\left(P-P_0\right)\left(1-a\right) + \frac{a\left(\left(1+b\left(P-P\mathrm{th}\right)\right)^{1-c} - \left(1+b\left(P_0-P\mathrm{th}\right)\right)^{1-c}\right)}{b(1-c)}\right)$$ $$a=\frac{1+K'_0}{1+K'_0+K_TK''_0}$$ $$b=\frac{K'_0}{K_0}-\frac{K''_0}{1+K'_0}$$ $$c=\frac{1+K'_0+K_0 K''_0}{{K'_0}^2 + K'_0 - K_0 K''_0}$$ $$P\mathrm{th}=\alpha_0K_0\frac{\theta}{\xi_0}\left(\frac{1}{\mathrm{e}^u-1} - \frac{1}{\mathrm{e}^{u_0}-1}\right)$$ $$\theta = \frac{10636}{S/n} + 6.44 \quad\xi=\frac{u^2\mathrm{e}^u}{\left(\mathrm{e}^u-1\right)^2} \quad u=\frac{\theta}{T}$$ $$n=\begin{cases} 7, & \mathrm{ulvospinel} \\ 5, & \mathrm{ilmenite} \end{cases}$$ $$\xi_o = \xi(T_0) \quad u_0=u(T_0) \quad T_0 = 298.15\,\mathrm{K}$$

Ilmenite phase transitions considered using Landau theory. The transition at 1 bar is included in the expression for \(\mu\mathrm{O}_2\), hence only the pressure dependent Landau term is sought by subtracting the pressure independent Landau contribution from the total using: $$G^\mathrm{ex}_{P,T} = G^\mathrm{ex}_{P,T}(V_\mathrm{max}) - G^\mathrm{ex}_{P,T}(V_\mathrm{max}=0)$$ $$G^\mathrm{ex}_{P,T} = S_\mathrm{max}\left(T_{c0}\left(Q^2+\frac{1}{3}\left(Q^6 - Q^6_0\right)\right) - T_{cP} Q^2 - T \left(Q^2_0 - Q^2 \right) \right)+P V_\mathrm{max} Q_0^2$$ $$T_{cP} = T_{c0}+P\frac{V_\mathrm{max}}{S_\mathrm{max}} \quad Q^2_0 = \sqrt{1-\frac{T_0}{T_{c0}}}$$ $$Q^2=\begin{cases} 0, & \mathrm{if}\,T > T_{cP} \\[2ex] \sqrt{\frac{T_{cP}-T}{T_{c0}}} & \mathrm{otherwise} \end{cases}$$

All required values for rutile and ilmenite can be obtained from the Holland and Powell (2011) dataset version 6.33 (tc-633.txt).

\(\int_1^P V\,dP\) for Fe is calculated using the equation of state described by Dorogokupets et al (2017) is used. The equation is given in terms of the Helmholz free energy \(F\). Expressions for \(F(V,T)\) and \(P(V)\) are given in the Supplementary Files of Dorogokupets et al (2017). \(V(P,T)\) is calculated by numerically solving \(P(V)-P=0\) for the pressure of interest. \(G(P,T)\) is calculated using: $$G = F+PV$$ This Gibbs energy is used to find phase transitions, and the resulting pressures are used as the boundaries for integration. The \(\int_1^P V^\mathrm{Fe}\,dP\) integral is evaluated numerically.

This equation of state is calibrated to >300 GPa and >3000 K. The BCC α-Fe and γ-Fe phases are stable at low pressures. At high pressures, the FCC γ-Fe exists at high temperatures whereas the HCP ε-Fe is stable at lower temperatures. Evaluation of \(G\) leads to prediction of FCC γ-Fe at extreme pressures (>100–150 GPa). This phase transition is spurious and results from extrapolation beyond the calibrated range, and as such is not considered here.

annite + almandine – muscovite + magnetite + quartz

$$\mathrm{KFe_3AlSi_3O_{10}(OH)_2 + Fe_3Al_2Si_3O_{12} + O_2 = KAl_3Si_3O_{10}(OH)_2 + 2Fe_3O_4 + 3SiO_2}$$

The BAMM buffer (biotite-almandine-muscovite-magnetite) was suggested by Zen (1985) as a representative redox indicator for peraluminous granites and metamorphic rocks. With only pure phases considered, it lies roughly between NNO and FMQ in upper crustal conditions. Zen's original formulation consists of polynomial fits to \(\Delta G\) calculated at fixed intervals using thermodynamic data available at the time, and using constant volume changes. Here, the buffer is calculated using accurate equations of state for the solids employing the more recently available thermodynamic data.

$$\log_{10}f\mathrm{O}_2=\frac{3\Delta G^\mathrm{quartz} + 2\Delta G^\mathrm{magnetite} + \Delta G^\mathrm{muscovite} - \Delta G^\mathrm{almandine} - \Delta G^\mathrm{annite} - \Delta G^\mathrm{O_2}}{RT\ln{10}}$$

$$\Delta G^\mathrm{quartz} = \Delta_fH^\circ + \int_{298.15}^T Cp\,dT - T\left(S^\circ + \int_{298.15}^T \frac{Cp}{T}\,dT\right) + G^\mathrm{excess} + \int_1^P V\,dP$$ $$\Delta G^\mathrm{magnetite} = \Delta_fH^\circ + \int_{298.15}^T Cp\,dT - T\left(S^\circ + \int_{298.15}^T \frac{Cp}{T}\,dT\right) + G^\mathrm{excess} + \int_1^P V\,dP$$ $$\Delta G^\mathrm{muscovite} = \Delta_fH^\circ + \int_{298.15}^T Cp\,dT - T\left(S^\circ + \int_{298.15}^T \frac{Cp}{T}\,dT\right) + \int_1^P V\,dP$$ $$\Delta G^\mathrm{almandine} = \Delta_fH^\circ + \int_{298.15}^T Cp\,dT - T\left(S^\circ + \int_{298.15}^T \frac{Cp}{T}\,dT\right) + \int_1^P V\,dP$$ $$\Delta G^\mathrm{annite} = \Delta_fH^\circ + \int_{298.15}^T Cp\,dT - T\left(S^\circ + \int_{298.15}^T \frac{Cp}{T}\,dT\right) + \int_1^P V\,dP$$ $$\Delta G^\mathrm{O_2} = H^{\circ}_\mathrm{O_2} - H^{\circ}_\mathrm{298.15,O_2} - TS^{\circ}_\mathrm{O_2}$$

\(\Delta G\) for solids is evaluated using: $$\int Cp\,dT = aT+0.5bT^2-cT^{-1}+2dT^{0.5}$$ $$\int\frac{Cp}{T}\,dT = a\ln{T}+bT-0.5cT^{-2}-2dT^{-0.5}$$ \(\Delta G^\mathrm{O_2}\) is evaluated using a Shomate equation and data from the NIST Chemistry WebBook “Gas phase thermochemistry data” pages. $$H^{\circ} - H^{\circ}_\mathrm{298.15} = AT + \frac{BT^2}{2}+ \frac{CT^3}{3} + \frac{DT^4}{4} - \frac{E}{T} +F - H$$ $$S^{\circ} = A\ln{T} + BT + \frac{CT^2}{2} + \frac{DT^3}{3} - \frac{E}{2T^2} + G$$

Fugacity at pressure is calculated using a modified Tait equation of state including a thermal pressure term as described by Holland and Powell (2011): $$V=V_0\left(1-a\left(1-\left(1+b\left(P-P\mathrm{th}\right)\right)^{-c}\right)\right)$$ $$\int_{P_0}^P V\,dP=V_0\left(\left(P-P_0\right)\left(1-a\right) + \frac{a\left(\left(1+b\left(P-P\mathrm{th}\right)\right)^{1-c} - \left(1+b\left(P_0-P\mathrm{th}\right)\right)^{1-c}\right)}{b(1-c)}\right)$$ $$a=\frac{1+K'_0}{1+K'_0+K_TK''_0}$$ $$b=\frac{K'_0}{K_0}-\frac{K''_0}{1+K'_0}$$ $$c=\frac{1+K'_0+K_0 K''_0}{{K'_0}^2 + K'_0 - K_0 K''_0}$$ $$P\mathrm{th}=\alpha_0K_0\frac{\theta}{\xi_0}\left(\frac{1}{\mathrm{e}^u-1} - \frac{1}{\mathrm{e}^{u_0}-1}\right)$$ $$\theta = \frac{10636}{S/n} + 6.44 \quad\xi=\frac{u^2\mathrm{e}^u}{\left(\mathrm{e}^u-1\right)^2} \quad u=\frac{\theta}{T}$$ $$n=\begin{cases} 22, & \mathrm{annite} \\ 20, & \mathrm{almandine} \\ 3, & \mathrm{quartz} \\ 7, & \mathrm{magnetite} \\ 21, & \mathrm{muscovite}\end{cases}$$ $$\xi_o = \xi(T_0) \quad u_0=u(T_0) \quad T_0 = 298.15\,\mathrm{K}$$

Landau theory is used to calculate excess Gibbs energy for quartz and magnetite phase transitions: $$G^\mathrm{excess} = S_\mathrm{max}\left(T_{c0}\left(Q^2+\frac{1}{3}\left(Q^6 - Q^6_0\right)\right) - T_{cP} Q^2 - T \left(Q^2_0 - Q^2 \right) \right)+P V_\mathrm{max} Q_0^2$$ $$T_{cP} = T_{c0}+P\frac{V_\mathrm{max}}{S_\mathrm{max}} \quad Q^2_0 = \sqrt{1-\frac{T_0}{T_{c0}}}$$ $$Q^2=\begin{cases} 0, & \mathrm{if}\,T > T_{cP} \\[2ex] \sqrt{\frac{T_{cP}-T}{T_{c0}}} & \mathrm{otherwise} \end{cases}$$

All required values can be obtained from the Holland and Powell (2011) dataset version 6.33 (tc-633.txt).

iron(II) oxide – iron(III) oxide
(wüstite – hematite)

$$\mathrm{4FeO + O_2 = 3Fe_2O_3}$$

At pressures below 10 GPa, the w&uuml;stite - hematite buffer is a hypothetical buffer, because the two minerals are not stable together and will react to form magnetite. It is, however, potentially useful as an indicator of the oxygen fugacity constrained by the metastable coexistence of simple end-member divalent and trivalent iron oxides (e.g. Fudali 1965). At higher pressures, magnetite is no longer stable and coexisting wüstite - hematite are stable. However, wüstite is non-stoichiomertric and data for the form of wüstite that exists in equilibrium with hematite does not exist. Instead, data for wüstite is taken from Campbell et al (2009) and is only an approximation. Real oxygen fugacity values are likely to be slightly higher than reported here.

$$\log_{10}f\mathrm{O}_2=\frac{2\Delta G^\mathrm{hematite} - 4\Delta G^\mathrm{wustite} - \Delta G^\mathrm{O_2}}{RT\ln{10}}$$

$$\Delta G^\mathrm{hematite} = \Delta_fH^\circ + \int_{298.15}^T Cp\,dT - T\left(S^\circ + \int_{298.15}^T \frac{Cp}{T}\,dT\right) + G^\mathrm{excess} + \int_1^P V\,dP$$ $$\Delta G^\mathrm{wustite} = \Delta_fH^\circ + \int_{298.15}^T Cp\,dT - T\left(S^\circ + \int_{298.15}^T \frac{Cp}{T}\,dT\right) + \int_1^P V\,dP$$ $$\Delta G^\mathrm{O_2} = H^{\circ}_\mathrm{O_2} - H^{\circ}_\mathrm{298.15,O_2} - TS^{\circ}_\mathrm{O_2}$$

\(\Delta G\) for solids is evaluated using: $$\int Cp\,dT = aT+0.5bT^2-cT^{-1}+2dT^{0.5}$$ $$\int\frac{Cp}{T}\,dT = a\ln{T}+bT-0.5cT^{-2}-2dT^{-0.5}$$ \(\Delta G^\mathrm{O_2}\) is evaluated using a Shomate equation and data from the NIST Chemistry WebBook “Gas phase thermochemistry data” pages. $$H^{\circ} - H^{\circ}_\mathrm{298.15} = AT + \frac{BT^2}{2}+ \frac{CT^3}{3} + \frac{DT^4}{4} - \frac{E}{T} +F - H$$ $$S^{\circ} = A\ln{T} + BT + \frac{CT^2}{2} + \frac{DT^3}{3} - \frac{E}{2T^2} + G$$

For hematite, fugacity at pressure is calculated using a modified Tait equation of state including a thermal pressure term as described by Holland and Powell (2011): $$V=V_0\left(1-a\left(1-\left(1+b\left(P-P\mathrm{th}\right)\right)^{-c}\right)\right)$$ $$\int_{P_0}^P V\,dP=V_0\left(\left(P-P_0\right)\left(1-a\right) + \frac{a\left(\left(1+b\left(P-P\mathrm{th}\right)\right)^{1-c} - \left(1+b\left(P_0-P\mathrm{th}\right)\right)^{1-c}\right)}{b(1-c)}\right)$$ $$a=\frac{1+K'_0}{1+K'_0+K_TK''_0}$$ $$b=\frac{K'_0}{K_0}-\frac{K''_0}{1+K'_0}$$ $$c=\frac{1+K'_0+K_0 K''_0}{{K'_0}^2 + K'_0 - K_0 K''_0}$$ $$P\mathrm{th}=\alpha_0K_0\frac{\theta}{\xi_0}\left(\frac{1}{\mathrm{e}^u-1} - \frac{1}{\mathrm{e}^{u_0}-1}\right)$$ $$\theta = \frac{10636}{S/n} + 6.44 \quad\xi=\frac{u^2\mathrm{e}^u}{\left(\mathrm{e}^u-1\right)^2} \quad u=\frac{\theta}{T} \quad n=5$$ $$\xi_o = \xi(T_0) \quad u_0=u(T_0) \quad T_0 = 298.15\,\mathrm{K}$$

Landau theory is used to calculate excess Gibbs energy for hematite phase transitions: $$G^\mathrm{excess} = S_\mathrm{max}\left(T_{c0}\left(Q^2+\frac{1}{3}\left(Q^6 - Q^6_0\right)\right) - T_{cP} Q^2 - T \left(Q^2_0 - Q^2 \right) \right)+P V_\mathrm{max} Q_0^2$$ $$T_{cP} = T_{c0}+P\frac{V_\mathrm{max}}{S_\mathrm{max}} \quad Q^2_0 = \sqrt{1-\frac{T_0}{T_{c0}}}$$ $$Q^2=\begin{cases} 0, & \mathrm{if}\,T > T_{cP} \\[2ex] \sqrt{\frac{T_{cP}-T}{T_{c0}}} & \mathrm{otherwise} \end{cases}$$

All required values can be obtained from the Holland and Powell (2011) dataset version 6.33 (tc-633.txt).

For wüstite, Campbell et al (2009) use a Mie-Grüneisen–Birch-Murnaghan equation of state: $$P=\frac{3}{2}K_0\left(\left(\frac{V_0}{V}\right)^{7/3}-\left(\frac{V_0}{V}\right)^{5/3}\right)\left(1-\frac{3}{4}\left(4-K'_0\right)\left(\left(\frac{V_0}{V}\right)^{2/3}-1\right)\right)+\gamma_0\frac{V^{q-1}}{V_0^q}\left(E(\theta_D,T) - E(\theta_D,295)\right)$$ \(E\left(\theta_D,T\right)=3nRTD_3\left(\frac{\theta_D}{T}\right)\) where \(D_3()\) is the Debye function and \(n=2\).

Direct calculation of \(\int_1^P V\,dP\) requires \(V=f(P)\) which is not available, so instead it is calculated using integration by parts: $$\int_{P_0}^{P_1} V\,dP = \Delta PV - \int_{V_0}^{V_1} P\,dV=P(V_1)V_1 - P(V_0)V_0 - \int_{V_0}^{V_1} P\,dV$$ $$\int P\,dV=\left(E(\theta_D,T) - E(\theta_D,295)\right)\frac{\gamma_0}{q}\left(\frac{V}{V_0}\right)^q-\frac{9K_0\left(\left(K_0'-4\right)V_0^3+\left(14-3K_0'\right)V^{2/3}V_0^{7/3}+\left(3K_0'-16\right)V^{4/3}V_0^{5/3}\right)}{16V^2}$$ Values for \(V_n\) are found numerically by solving for \(V\) which satisfy \(P(V)-P=0\).

All required parameters are given in Table 1 of Campbell et al (2009) and are calibrated to about 100 GPa and 2500 K.

fayalite – hematite + quartz

$$\mathrm{2Fe_2SiO_4 + O_2 = 2Fe_2O_3 + 2SiO_2}$$

The hypothetical FHQ buffer was suggested by Giggenbach (1987) as a “non-specific” redox buffer containing fayalite as a proxy for Fe²⁺ and hematite as a proxy for Fe³⁺. It is intended to approximate rock-buffering conditions in volcanic systems, and occasionally called the “rock” buffer. See Moretti and Stefánsson (2020) for more details.

$$\log_{10}f\mathrm{O}_2=\frac{2\Delta G^\mathrm{quartz} + 2\Delta G^\mathrm{hematite} - 2\Delta G^\mathrm{fayalite} - \Delta G^\mathrm{O_2}}{RT\ln{10}}$$

$$\Delta G^\mathrm{quartz} = \Delta_fH^\circ + \int_{298.15}^T Cp\,dT - T\left(S^\circ + \int_{298.15}^T \frac{Cp}{T}\,dT\right) + G^\mathrm{excess} + \int_1^P V\,dP$$ $$\Delta G^\mathrm{hematite} = \Delta_fH^\circ + \int_{298.15}^T Cp\,dT - T\left(S^\circ + \int_{298.15}^T \frac{Cp}{T}\,dT\right) + G^\mathrm{excess} + \int_1^P V\,dP$$ $$\Delta G^\mathrm{fayalite} = \Delta_fH^\circ + \int_{298.15}^T Cp\,dT - T\left(S^\circ + \int_{298.15}^T \frac{Cp}{T}\,dT\right) + \int_1^P V\,dP$$ $$\Delta G^\mathrm{O_2} = H^{\circ}_\mathrm{O_2} - H^{\circ}_\mathrm{298.15,O_2} - TS^{\circ}_\mathrm{O_2}$$

\(\Delta G\) for solids is evaluated using: $$\int Cp\,dT = aT+0.5bT^2-cT^{-1}+2dT^{0.5}$$ $$\int\frac{Cp}{T}\,dT = a\ln{T}+bT-0.5cT^{-2}-2dT^{-0.5}$$ \(\Delta G^\mathrm{O_2}\) is evaluated using a Shomate equation and data from the NIST Chemistry WebBook “Gas phase thermochemistry data” pages. $$H^{\circ} - H^{\circ}_\mathrm{298.15} = AT + \frac{BT^2}{2}+ \frac{CT^3}{3} + \frac{DT^4}{4} - \frac{E}{T} +F - H$$ $$S^{\circ} = A\ln{T} + BT + \frac{CT^2}{2} + \frac{DT^3}{3} - \frac{E}{2T^2} + G$$

Fugacity at pressure is calculated using a modified Tait equation of state including a thermal pressure term as described by Holland and Powell (2011): $$V=V_0\left(1-a\left(1-\left(1+b\left(P-P\mathrm{th}\right)\right)^{-c}\right)\right)$$ $$\int_{P_0}^P V\,dP=V_0\left(\left(P-P_0\right)\left(1-a\right) + \frac{a\left(\left(1+b\left(P-P\mathrm{th}\right)\right)^{1-c} - \left(1+b\left(P_0-P\mathrm{th}\right)\right)^{1-c}\right)}{b(1-c)}\right)$$ $$a=\frac{1+K'_0}{1+K'_0+K_TK''_0}$$ $$b=\frac{K'_0}{K_0}-\frac{K''_0}{1+K'_0}$$ $$c=\frac{1+K'_0+K_0 K''_0}{{K'_0}^2 + K'_0 - K_0 K''_0}$$ $$P\mathrm{th}=\alpha_0K_0\frac{\theta}{\xi_0}\left(\frac{1}{\mathrm{e}^u-1} - \frac{1}{\mathrm{e}^{u_0}-1}\right)$$ $$\theta = \frac{10636}{S/n} + 6.44 \quad\xi=\frac{u^2\mathrm{e}^u}{\left(\mathrm{e}^u-1\right)^2} \quad u=\frac{\theta}{T}$$ $$n=\begin{cases} 5, & \mathrm{hematite} \\ 7, & \mathrm{fayalite} \\ 3, & \mathrm{quartz}\end{cases}$$ $$\xi_o = \xi(T_0) \quad u_0=u(T_0) \quad T_0 = 298.15\,\mathrm{K}$$

Landau theory is used to calculate excess Gibbs energy for quartz and hematite phase transitions: $$G^\mathrm{excess} = S_\mathrm{max}\left(T_{c0}\left(Q^2+\frac{1}{3}\left(Q^6 - Q^6_0\right)\right) - T_{cP} Q^2 - T \left(Q^2_0 - Q^2 \right) \right)+P V_\mathrm{max} Q_0^2$$ $$T_{cP} = T_{c0}+P\frac{V_\mathrm{max}}{S_\mathrm{max}} \quad Q^2_0 = \sqrt{1-\frac{T_0}{T_{c0}}}$$ $$Q^2=\begin{cases} 0, & \mathrm{if}\,T > T_{cP} \\[2ex] \sqrt{\frac{T_{cP}-T}{T_{c0}}} & \mathrm{otherwise} \end{cases}$$

All required values can be obtained from the Holland and Powell (2011) dataset version 6.33 (tc-633.txt).

supercritical oxygen fluid

$$\mathrm{O_2}$$

This model uses a combination of two equations of state. The first equation (BS), by Belonoshko and Saxena (1991) is only calibrated for pressure greater than 5 kbar. The second equation is a Redlich-Kwong (RK) equation with parameters from Connolly (2016), but it overestimates oxyxgen fugacity at pressures greater than about 50 kbar.

Here, the RK equation is used at pressures below 5 kbar and the BS equation is used at pressures above 40 kbar. In between, \(\log_{10}f\mathrm{O}_2\) is varied linearly between 5 and 40 kbar: $$\log_{10}f\mathrm{O}_2 = \log_{10}f\mathrm{O}_2^\mathrm{BS}\times x + \log_{10}f\mathrm{O}_2^\mathrm{RK}\times (1-x)$$ $$x = \frac{P - 5000}{40000 - 5000}$$

BS

The relevant equations from Belonoshko and Saxena (1991) are $$f\mathrm{O}_2=f_\mathrm{5\,kbar}\times\exp\left(\frac{\int_5^P V\,dP}{RT}\right)$$ $$f_\mathrm{5\,kbar} = \exp\left(\frac{1000\left(f_1 + f_2T + f_3 \ln{T}\right)}{RT}\right)$$ $$P = \frac{a}{V} + \frac{b}{V^2} + \frac{c}{V^n}$$ $$a=\left(a_1 + a_2T/1000\right)\times10^4$$ $$b=\left(b_1 + b_2T/1000\right)\times10^6$$ $$c=\left(c_1 + c_2T/1000\right)\times10^9$$

Direct calculation of \(\int_1^P V\,dP\) requires \(V=f(P)\) which is not available, so instead it is calculated using integration by parts: $$\int_{P_0=\mathrm{5\,kbar}}^{P_1} V\,dP = \Delta PV - \int_{V_0}^{V_1} P\,dV=P(V_1)V_1 - P(V_0)V_0 - \int_{V_0}^{V_1} P\,dV$$ $$\int V\,dP = a\left(1-\ln V\right) + \frac{2b}{V} + \frac{cnV^{1-n}}{n-1}$$ The limits of integration are determined numerically by solving for \(V\) which satisfy \(P(V)-P=0\).

All parameters are given in Tables 4 and 5 of Belonoshko and Saxena (1991). Calibrated from 5 kbar to 1 Mbar, and from 400 K to 4000 K:

RK

Oxygen fugacity is calculated using a Redlich-Kwong equation of state: $$f\mathrm{O}_2=\exp\left(\frac{\int_1^P V\,dP}{RT}\right)$$ $$P=\frac{RT}{V-b} - \frac{a}{\sqrt{T}V(V-b)}$$ $$a=\frac{1}{9\left(\sqrt[3]{2}-1\right)}\times\frac{R^2T_\mathrm{c}^{2.5}}{P_\mathrm{c}}\quad b=\frac{\sqrt[3]{2}-1}{3}\times\frac{RT_\mathrm{c}}{P_\mathrm{c}}$$

Direct calculation of \(\int_1^P V\,dP\) requires \(V=f(P)\) which is not available, so instead it is calculated using integration by parts: $$\int_{P_0=\mathrm{1\,bar}}^{P_1} V\,dP = \Delta PV - \int_{V_0}^{V_1} P\,dV=P(V_1)V_1 - P(V_0)V_0 - \int_{V_0}^{V_1} P\,dV$$ $$\int V\,dP = RT\left(\frac{V}{V-b} - \ln(V-b)\right) + \frac{a}{\sqrt{T}}\left(\frac{\ln{\frac{V}{V+b}}}{b} - \frac{1}{V+b}\right)$$ The limits of integration are determined numerically by solving for \(V\) which satisfy \(P(V)-P=0\). All parameters are given in Table 2 (row: “vWRK”) of Connolly (2016). Also works at pressures below 5 kbar.