diff --git a/burnman/calibrants/Anderson_1989.py b/burnman/calibrants/Anderson_1989.py new file mode 100644 index 00000000..29868441 --- /dev/null +++ b/burnman/calibrants/Anderson_1989.py @@ -0,0 +1,49 @@ +# This file is part of BurnMan - a thermoelastic and thermodynamic toolkit for +# the Earth and Planetary Sciences +# Copyright (C) 2012 - 2024 by the BurnMan team, released under the GNU +# GPL v2 or later. + +import numpy as np +from burnman.eos.birch_murnaghan import BirchMurnaghanBase as BM3 +from burnman.classes.calibrant import Calibrant +from burnman.utils.unitcell import molar_volume_from_unit_cell_volume + +""" +Anderson_1989 +^^^^^^^^^^^^^ +""" + + +class Au(Calibrant): + """ + The Au pressure standard reported by + Anderson et al. (1989; https://doi.org/10.1063/1.342969). + """ + + def __init__(self): + def _pressure(volume, temperature, params): + + # Isothermal pressure (GPa) + pressure_model = BM3() + P0 = pressure_model.pressure(params["T_0"], volume, params) + + # Thermal pressure + Pth_rel = (7.14e6 + params["dKTdT"] * np.log(params["V_0"] / volume)) * ( + temperature - 300.0 + ) + + return P0 + Pth_rel + + Z = 4.0 + _params = { + "V_0": molar_volume_from_unit_cell_volume(67.850, Z), + "K_0": 166.65e9, + "Kprime_0": 5.4823, + "dKTdT": -11.5e6, + "n": 1.0, + "T_0": 300.0, + "P_0": 0.0, + "Z": Z, + } + + Calibrant.__init__(self, _pressure, "pressure", _params) diff --git a/burnman/calibrants/Dorogokupets_2007.py b/burnman/calibrants/Dorogokupets_2007.py new file mode 100644 index 00000000..21149066 --- /dev/null +++ b/burnman/calibrants/Dorogokupets_2007.py @@ -0,0 +1,478 @@ +# This file is part of BurnMan - a thermoelastic and thermodynamic toolkit for +# the Earth and Planetary Sciences +# Copyright (C) 2012 - 2024 by the BurnMan team, released under the GNU +# GPL v2 or later. + +import numpy as np +from burnman.classes.calibrant import Calibrant +from burnman.constants import gas_constant +from burnman.eos.vinet import Vinet + +""" +Dorogokupets_2007 +^^^^^^^^^^^^^^^^^ +""" + +_materials_data = { + "Ag": { + "E_0": 0.0, + "P_0": 0.0, + "T_0": 298.15, + "n": 1.0, + "V_0": 10.272e-6, + "K_0": 99.65e9, + "Kprime_0": 6.11, + "theta_B1": 130.6, + "d_B1": 8.572, + "m_B1": 0.121, + "theta_B2": 103.6, + "d_B2": 5.326, + "m_B2": 0.449, + "theta_E1": 111.9, + "m_E1": 0.766, + "theta_E2": 189.12, + "m_E2": 1.664, + "gamma_0": 2.376, + "gamma_inf": 1.481, + "beta": 2.507, + "a": 6.70e-6, + "m": 3.44, + "e": 25.9e-6, + "f": -1.0, + "g": 0.666, + "h": -2.0, + "H": 15239, + "S": 0.732, + }, + "Al": { + "E_0": 0.0, + "P_0": 0.0, + "T_0": 298.15, + "n": 1.0, + "V_0": 9.999e-6, + "K_0": 72.67e9, + "Kprime_0": 4.62, + "theta_B1": 245.8, + "d_B1": 5.575, + "m_B1": 0.987, + "theta_B2": 300.0, # Not used as m = 0 + "d_B2": 1.0, # Not used as m = 0 + "m_B2": 0.0, + "theta_E1": 240.2, + "m_E1": 1.000, + "theta_E2": 356.2, + "m_E2": 1.013, + "gamma_0": 2.144, + "gamma_inf": 1.017, + "beta": 3.942, + "a": 5.14e-6, + "m": 3.44, + "e": 54.1e-6, + "f": -1.0, + "g": 1.8, + "h": -2.0, + "H": 8679, + "S": 0.998, + }, + "Au": { + "E_0": 0.0, + "P_0": 0.0, + "T_0": 298.15, + "n": 1.0, + "V_0": 10.215e-6, + "K_0": 166.70e9, + "Kprime_0": 6.00, + "theta_B1": 95.7, + "d_B1": 8.290, + "m_B1": 0.681, + "theta_B2": 106.4, + "d_B2": 3.239, + "m_B2": 0.417, + "theta_E1": 170.6, + "m_E1": 1.063, + "theta_E2": 105.2, + "m_E2": 0.839, + "gamma_0": 2.965, + "gamma_inf": 1.142, + "beta": 3.030, + "a": 25.33e-6, + "m": 3.79, + "e": 18.92e-6, + "f": -1.0, + "g": 0.66, + "h": -2.0, + "H": 11.69e3, # note typo + "S": 1.067, + }, + "Cu": { + "E_0": 0.0, + "P_0": 0.0, + "T_0": 298.15, + "n": 1.0, + "V_0": 7.113e-6, + "K_0": 133.41e9, + "Kprime_0": 5.37, + "theta_B1": 123.7, + "d_B1": 3.776, + "m_B1": 0.115, + "theta_B2": 175.4, + "d_B2": 10.372, + "m_B2": 0.711, + "theta_E1": 187.4, + "m_E1": 0.756, + "theta_E2": 286.9, + "m_E2": 1.418, + "gamma_0": 1.974, + "gamma_inf": 1.554, + "beta": 4.647, + "a": 3.50e-6, + "m": 3.46, + "e": 27.698e-6, + "f": -1.0, + "g": 0.666, + "h": -2.0, + "H": 11687, + "S": 1.407, + }, + "Pt": { + "E_0": 0.0, + "P_0": 0.0, + "T_0": 298.15, + "n": 1.0, + "V_0": 9.091e-6, + "K_0": 276.07e9, + "Kprime_0": 5.30, + "theta_B1": 95.2, + "d_B1": 8.199, + "m_B1": 0.329, + "theta_B2": 148.4, + "d_B2": 4.005, + "m_B2": 0.383, + "theta_E1": 214.6, + "m_E1": 1.211, + "theta_E2": 140.8, + "m_E2": 1.077, + "gamma_0": 2.802, + "gamma_inf": 1.538, + "beta": 5.550, + "a": 160.9e-6, + "m": 4.06, + "e": 260.0e-6, + "f": -1.0, + "g": 2.4, + "h": -2.0, + "H": 32572, + "S": 0.631, + }, + "Ta": { + "E_0": 0.0, + "P_0": 0.0, + "T_0": 298.15, + "n": 1.0, + "V_0": 10.851e-6, + "K_0": 191.39e9, + "Kprime_0": 3.81, + "theta_B1": 72.6, + "d_B1": 5.536, + "m_B1": 0.117, + "theta_B2": 101.8, + "d_B2": 24.513, + "m_B2": 0.396, + "theta_E1": 144.0, + "m_E1": 1.118, + "theta_E2": 214.9, + "m_E2": 1.369, + "gamma_0": 1.714, + "gamma_inf": 1.241, + "beta": 6.825, + "a": 61.9e-6, + "m": 4.00, + "e": 167.0e-6, + "f": -1.0, + "g": 1.3, + "h": -2.0, + "H": 36278, + "S": 4.910, + }, + "W": { + "E_0": 0.0, + "P_0": 0.0, + "T_0": 298.15, + "n": 1.0, + "V_0": 9.545e-6, + "K_0": 306.00e9, + "Kprime_0": 4.17, + "theta_B1": 182.8, + "d_B1": 13.270, + "m_B1": 0.513, + "theta_B2": 172.5, + "d_B2": 3.305, + "m_B2": 0.174, + "theta_E1": 287.6, + "m_E1": 1.166, + "theta_E2": 213.8, + "m_E2": 1.145, + "gamma_0": 1.553, + "gamma_inf": 0.694, + "beta": 3.698, + "a": -39.3e-6, + "m": 2.67, + "e": 40.4e-6, + "f": -1.0, + "g": 0.2, + "h": -2.0, + "H": 14714, + "S": 0.672, + }, + "MgO": { + "E_0": 0.0, + "P_0": 0.0, + "T_0": 298.15, + "n": 2.0, + "V_0": 11.248e-6, + "K_0": 160.31e9, + "Kprime_0": 4.18, + "theta_B1": 447.3, + "d_B1": 11.248, + "m_B1": 1.429, + "theta_B2": 384.0, + "d_B2": 3.593, + "m_B2": 0.276, + "theta_E1": 703.8, + "m_E1": 2.570, + "theta_E2": 466.0, + "m_E2": 1.725, + "gamma_0": 1.522, + "gamma_inf": 1.111, + "beta": 4.509, + "a": 13.56e-6, + "m": 5.23, + "e": 0, + "f": 0.0, + "g": 0, + "h": 0.0, + "H": 0, + "S": 0, + }, + "Diamond": { + "E_0": 0.0, + "P_0": 0.0, + "T_0": 298.15, + "n": 1.0, + "V_0": 3.417e-6, + "K_0": 443.16e9, + "Kprime_0": 3.777, + "theta_B1": 1202.1, + "d_B1": 9.604, + "m_B1": 1.163, + "theta_B2": 1135.1, + "d_B2": 3.380, + "m_B2": 0.218, + "theta_E1": 1687.2, + "m_E1": 1.396, + "theta_E2": 1033.7, + "m_E2": 0.223, + "gamma_0": 0.820, + "gamma_inf": 0.615, + "beta": 10.121, + "a": -23.85e-6, + "m": 1.22, + "e": 0, + "f": 0.0, + "g": 0, + "h": 0.0, + "H": 0, + "S": 0, + }, +} + + +do_cold_eos = Vinet() + + +def _F_qhB(m, d, theta, temperature): + # Equation 8a and equations just after Equation 7 + g = d * np.log(1.0 + theta / (temperature * d)) + b = 1.0 / (np.exp(g) - 1.0) + return m * ((d - 1.0) / (2.0 * d) * theta - temperature * np.log(1.0 + b)) + + +def _F_qhE(m, theta, temperature): + # Equation 8b + return m * (theta / 2.0 + temperature * np.log(1.0 - np.exp(-theta / temperature))) + + +def _F_anh_factor(m, theta, temperature): + # Equation 12 + zeta = np.exp(theta / temperature) + f = np.power(0.5 * theta + theta / (zeta - 1.0), 2.0) + 2.0 * theta * theta * ( + zeta / np.power(zeta - 1.0, 2.0) + ) + return m * f + + +def _F_th(temperature, x, params): + # Equation 10 + F_theta = np.power(x, -params["gamma_inf"]) * np.exp( + (params["gamma_0"] - params["gamma_inf"]) + / params["beta"] + * (1.0 - np.power(x, params["beta"])) + ) + theta_B1 = params["theta_B1"] * F_theta + theta_B2 = params["theta_B2"] * F_theta + theta_E1 = params["theta_E1"] * F_theta + theta_E2 = params["theta_E2"] * F_theta + + # Equation 8 + F_qh = ( + _F_qhB(params["m_B1"], params["d_B1"], theta_B1, temperature) + + _F_qhB(params["m_B2"], params["d_B2"], theta_B2, temperature) + + _F_qhE(params["m_E1"], theta_E1, temperature) + + _F_qhE(params["m_E2"], theta_E2, temperature) + ) + + # Equation 12 + F_anh = ( + params["a"] + * np.power(x, params["m"]) + / 6.0 + * ( + _F_anh_factor(params["m_B1"], theta_B1, temperature) + + _F_anh_factor(params["m_B2"], theta_B2, temperature) + + _F_anh_factor(params["m_E1"], theta_E1, temperature) + + _F_anh_factor(params["m_E2"], theta_E2, temperature) + ) + ) + + # Equation 13 + F_el = ( + -1.5 + * params["n"] + * params["e"] + * np.power(x, params["g"]) + * temperature + * temperature + ) + + # Equation 14 + F_def = ( + -1.5 + * params["n"] + * temperature + * np.exp( + params["S"] * np.power(x, params["f"]) + - params["H"] * np.power(x, params["h"]) / temperature + ) + ) + return F_qh + F_anh + F_el + F_def + + +def _helmholtz_energy(volume, temperature, params): + # Equation 6 + E = do_cold_eos.molar_internal_energy(0.0, 0.0, volume, params) + + # Just after Equation 6b + Vrel = volume / params["V_0"] + Fth = _F_th(temperature, Vrel, params) + Fth0 = _F_th(params["T_0"], Vrel, params) + + # Equation 5 + return E + gas_constant * (Fth - Fth0) + + +def _pressure(volume, temperature, params): + dV = params["V_0"] * 1.0e-7 + + F1 = _helmholtz_energy(volume + dV / 2.0, temperature, params) + F0 = _helmholtz_energy(volume - dV / 2.0, temperature, params) + return -(F1 - F0) / dV + + +class Ag(Calibrant): + """ + The Ag pressure standard reported by + Dorogokupets and Oganov (2007; https://doi.org/10.1103/PhysRevB.75.024115). + """ + + def __init__(self): + Calibrant.__init__(self, _pressure, "pressure", _materials_data["Ag"]) + + +class Al(Calibrant): + """ + The Al pressure standard reported by + Dorogokupets and Oganov (2007; https://doi.org/10.1103/PhysRevB.75.024115). + """ + + def __init__(self): + Calibrant.__init__(self, _pressure, "pressure", _materials_data["Al"]) + + +class Au(Calibrant): + """ + The Au pressure standard reported by + Dorogokupets and Oganov (2007; https://doi.org/10.1103/PhysRevB.75.024115). + """ + + def __init__(self): + Calibrant.__init__(self, _pressure, "pressure", _materials_data["Au"]) + + +class Cu(Calibrant): + """ + The Cu pressure standard reported by + Dorogokupets and Oganov (2007; https://doi.org/10.1103/PhysRevB.75.024115). + """ + + def __init__(self): + Calibrant.__init__(self, _pressure, "pressure", _materials_data["Cu"]) + + +class diamond(Calibrant): + """ + The diamond pressure standard reported by + Dorogokupets and Oganov (2007; https://doi.org/10.1103/PhysRevB.75.024115). + """ + + def __init__(self): + Calibrant.__init__(self, _pressure, "pressure", _materials_data["Diamond"]) + + +class MgO(Calibrant): + """ + The MgO pressure standard reported by + Dorogokupets and Oganov (2007; https://doi.org/10.1103/PhysRevB.75.024115). + """ + + def __init__(self): + Calibrant.__init__(self, _pressure, "pressure", _materials_data["MgO"]) + + +class Pt(Calibrant): + """ + The Pt pressure standard reported by + Dorogokupets and Oganov (2007; https://doi.org/10.1103/PhysRevB.75.024115). + """ + + def __init__(self): + Calibrant.__init__(self, _pressure, "pressure", _materials_data["Pt"]) + + +class Ta(Calibrant): + """ + The Ta pressure standard reported by + Dorogokupets and Oganov (2007; https://doi.org/10.1103/PhysRevB.75.024115). + """ + + def __init__(self): + Calibrant.__init__(self, _pressure, "pressure", _materials_data["Ta"]) + + +class W(Calibrant): + """ + The W pressure standard reported by + Dorogokupets and Oganov (2007; https://doi.org/10.1103/PhysRevB.75.024115). + """ + + def __init__(self): + Calibrant.__init__(self, _pressure, "pressure", _materials_data["W"]) diff --git a/burnman/calibrants/__init__.py b/burnman/calibrants/__init__.py index 3ceca19f..a9564d5f 100644 --- a/burnman/calibrants/__init__.py +++ b/burnman/calibrants/__init__.py @@ -7,7 +7,7 @@ """ Calibrant database - - :mod:`~burnman.calibrants.Decker_1971` + - :mod:`~burnman.calibrants.Anderson_1989` - :mod:`~burnman.calibrants.Armentrout_2015` - :mod:`~burnman.calibrants.Campbell_2009` - :mod:`~burnman.calibrants.Chidester_2021` @@ -16,6 +16,7 @@ - :mod:`~burnman.calibrants.Dewaele_2012` - :mod:`~burnman.calibrants.Dewaele_2013` - :mod:`~burnman.calibrants.Dewaele_2020` + - :mod:`~burnman.calibrants.Dorogokupets_2007` - :mod:`~burnman.calibrants.Dorogokupets_2017` - :mod:`~burnman.calibrants.Dubrovinsky_1998` - :mod:`~burnman.calibrants.Fei_2007` @@ -46,6 +47,7 @@ - :mod:`~burnman.calibrants.Zhao_2000` """ +from . import Anderson_1989 from . import Armentrout_2015 from . import Campbell_2009 from . import Chidester_2021 @@ -54,6 +56,7 @@ from . import Dewaele_2012 from . import Dewaele_2013 from . import Dewaele_2020 +from . import Dorogokupets_2007 from . import Dorogokupets_2017 from . import Dubrovinsky_1998 from . import Fei_2007 diff --git a/misc/benchmarks/calibrant_benchmarks.py b/misc/benchmarks/calibrant_benchmarks.py index b22bfaf5..81390118 100644 --- a/misc/benchmarks/calibrant_benchmarks.py +++ b/misc/benchmarks/calibrant_benchmarks.py @@ -63,6 +63,23 @@ def check_figures(): ) +def check_Anderson_1989(): + Au = calibrants.Anderson_1989.Au() + name = str(Au).split(" ")[0][1:] + compression = 0.34 + print(f"Checking {name} at compression {compression}...") + V = (1 - compression) * Au.params["V_0"] + T = 300.0 + print( + f"Temperature: {T} K, Pressure: {Au.pressure(V, T)/1.e9:.2f} (reported to be 216.06 GPa)" + ) + T = 3000.0 + print( + f"Temperature: {T} K, Pressure: {Au.pressure(V, T)/1.e9:.2f} (reported to be 222.44 GPa)" + ) + print() + + def check_Decker_1971(): NaCl = calibrants.Decker_1971.NaCl_B1() V_0 = NaCl.params["V_0"] @@ -74,7 +91,7 @@ def check_Decker_1971(): "Pressures from Decker 1971 calibrant " "vs. tabulated data in original paper (given in GPa)" ) - print(f"\nV={V:.4e} m^3/mol (standard state volume):") + print(f"V={V:.4e} m^3/mol (standard state volume):") T_C = [25.0, 100.0, 200.0, 300.0, 500.0, 800.0] P_kbar = [0.00, 2.13, 5.00, 7.89, 13.72, 22.48] for i, T in enumerate(T_C): @@ -97,6 +114,38 @@ def check_Decker_1971(): print("") +def check_Dorogokupets_Oganov_2007(): + tests = [ + [calibrants.Dorogokupets_2007.Ag(), 0.7, 100.371, 104.735, 111.146], + [calibrants.Dorogokupets_2007.Al(), 0.7, 56.986, 60.046, 64.905], + [calibrants.Dorogokupets_2007.Au(), 0.7, 164.83, 169.00, 175.03], + [calibrants.Dorogokupets_2007.Cu(), 0.7, 118.576, 124.032, 132.235], + [calibrants.Dorogokupets_2007.Pt(), 0.7, 242.676, 247.403, 254.730], + [calibrants.Dorogokupets_2007.Ta(), 0.7, 130.887, 133.891, 138.468], + [calibrants.Dorogokupets_2007.W(), 0.7, 222.444, 224.854, 228.550], + [calibrants.Dorogokupets_2007.MgO(), 0.7, 116.670, 120.904, 128.081], + [calibrants.Dorogokupets_2007.diamond(), 0.7, 301.514, 303.890, 309.757], + ] + + temperatures = [298.15, 1000, 2000] + + for test in tests: + calibrant = test[0] + Vrel = test[1] + pressures_GPa = test[2:] + + name = str(calibrant).split(" ")[0][1:] + print(f"Checking {name}...") + for i in range(len(temperatures)): + P = calibrant.pressure(calibrant.params["V_0"] * Vrel, temperatures[i]) + print( + f"V={Vrel}*V0, P={P/1.e9:.2f} GPa, T={temperatures[i]} K: {P/1.e6 - pressures_GPa[i]*1.e3:.3f} MPa difference" + ) + print() + + if __name__ == "__main__": + check_Anderson_1989() check_Decker_1971() + check_Dorogokupets_Oganov_2007() check_figures() diff --git a/misc/ref/calibrant_benchmarks.py.out b/misc/ref/calibrant_benchmarks.py.out index c6797f7f..7a5c8516 100644 --- a/misc/ref/calibrant_benchmarks.py.out +++ b/misc/ref/calibrant_benchmarks.py.out @@ -1,5 +1,8 @@ -Pressures from Decker 1971 calibrant vs. tabulated data in original paper (given in GPa) +Checking burnman.calibrants.Anderson_1989.Au at compression 0.34... +Temperature: 300.0 K, Pressure: 216.06 (reported to be 216.06 GPa) +Temperature: 3000.0 K, Pressure: 222.44 (reported to be 222.44 GPa) +Pressures from Decker 1971 calibrant vs. tabulated data in original paper (given in GPa) V=2.7013e-05 m^3/mol (standard state volume): 25.0 C: 0.00, 0.00 100.0 C: 0.21, 0.21 @@ -26,6 +29,43 @@ V=1.9044e-05 m^3/mol (last row): 500.0 C: 20.77, 20.78 800.0 C: 21.66, 21.67 +Checking burnman.calibrants.Dorogokupets_2007.Ag... +V=0.7*V0, P=100.37 GPa, T=298.15 K: 0.250 MPa difference +V=0.7*V0, P=104.73 GPa, T=1000 K: -0.662 MPa difference +V=0.7*V0, P=111.15 GPa, T=2000 K: -0.997 MPa difference +Checking burnman.calibrants.Dorogokupets_2007.Al... +V=0.7*V0, P=56.97 GPa, T=298.15 K: -11.383 MPa difference +V=0.7*V0, P=60.03 GPa, T=1000 K: -12.426 MPa difference +V=0.7*V0, P=64.86 GPa, T=2000 K: -43.175 MPa difference +Checking burnman.calibrants.Dorogokupets_2007.Au... +V=0.7*V0, P=164.83 GPa, T=298.15 K: -0.482 MPa difference +V=0.7*V0, P=169.00 GPa, T=1000 K: -4.064 MPa difference +V=0.7*V0, P=175.03 GPa, T=2000 K: 4.056 MPa difference +Checking burnman.calibrants.Dorogokupets_2007.Cu... +V=0.7*V0, P=118.65 GPa, T=298.15 K: 78.072 MPa difference +V=0.7*V0, P=124.11 GPa, T=1000 K: 77.859 MPa difference +V=0.7*V0, P=132.31 GPa, T=2000 K: 77.234 MPa difference +Checking burnman.calibrants.Dorogokupets_2007.Pt... +V=0.7*V0, P=242.66 GPa, T=298.15 K: -13.987 MPa difference +V=0.7*V0, P=247.39 GPa, T=1000 K: -15.448 MPa difference +V=0.7*V0, P=254.71 GPa, T=2000 K: -16.100 MPa difference +Checking burnman.calibrants.Dorogokupets_2007.Ta... +V=0.7*V0, P=130.95 GPa, T=298.15 K: 60.149 MPa difference +V=0.7*V0, P=133.95 GPa, T=1000 K: 58.981 MPa difference +V=0.7*V0, P=138.53 GPa, T=2000 K: 58.331 MPa difference +Checking burnman.calibrants.Dorogokupets_2007.W... +V=0.7*V0, P=222.43 GPa, T=298.15 K: -17.322 MPa difference +V=0.7*V0, P=224.83 GPa, T=1000 K: -19.063 MPa difference +V=0.7*V0, P=228.53 GPa, T=2000 K: -20.886 MPa difference +Checking burnman.calibrants.Dorogokupets_2007.MgO... +V=0.7*V0, P=116.72 GPa, T=298.15 K: 52.963 MPa difference +V=0.7*V0, P=120.96 GPa, T=1000 K: 51.888 MPa difference +V=0.7*V0, P=128.13 GPa, T=2000 K: 51.510 MPa difference +Checking burnman.calibrants.Dorogokupets_2007.diamond... +V=0.7*V0, P=301.53 GPa, T=298.15 K: 13.945 MPa difference +V=0.7*V0, P=303.90 GPa, T=1000 K: 13.360 MPa difference +V=0.7*V0, P=309.77 GPa, T=2000 K: 10.797 MPa difference + Checking burnman.calibrants.Fei_2007.Au... Checking burnman.calibrants.Fei_2007.Pt... Checking burnman.calibrants.Holmes_1989.Pt...