Equation of State¶
Set up environment (optional)¶
These steps are required for Google Colab, but may work on other systems too:
[1]:
# import locale
# locale.getpreferredencoding = lambda: "UTF-8"
# ! pip uninstall torch torchaudio torchvision numpy -y
# ! uv pip install janus-core[all] data-tutorials torch==2.5.1 --system
# get_ipython().kernel.do_shutdown(restart=True)
[2]:
from ase import Atoms
from ase.io import read
from weas_widget import WeasWidget
from data_tutorials.data import get_data
from janus_core.calculations.eos import EoS
Use data_tutorials to get the data required for this tutorial:
[3]:
get_data(
url="https://raw.githubusercontent.com/stfc/janus-tutorials/main/data/",
filename=["beta_quartz.cif"],
folder="data",
)
try to download beta_quartz.cif from https://raw.githubusercontent.com/stfc/janus-tutorials/main/data/ and save it in data/beta_quartz.cif
saved in data/beta_quartz.cif
Equation of state for α-quartz¶
Build the structure:
[4]:
α_quartz = Atoms(
symbols=(*["Si"] * 3, *["O"] * 6),
scaled_positions=[
[0.469700, 0.000000, 0.000000],
[0.000000, 0.469700, 0.666667],
[0.530300, 0.530300, 0.333333],
[0.413500, 0.266900, 0.119100],
[0.266900, 0.413500, 0.547567],
[0.733100, 0.146600, 0.785767],
[0.586500, 0.853400, 0.214233],
[0.853400, 0.586500, 0.452433],
[0.146600, 0.733100, 0.880900],
],
cell=[
[4.916000, 0.000000, 0.000000],
[-2.45800, 4.257381, 0.000000],
[0.000000, 0.000000, 5.405400],
],
pbc=True,
)
v=WeasWidget()
v.from_ase(α_quartz)
v.avr.model_style = 1
v.avr.show_hydrogen_bonds = True
v
[4]:
Calculate the equation of state using the MACE-MP potential:
[5]:
mace_eos = EoS(
struct=α_quartz.copy(),
arch="mace_mp",
device="cpu",
model_path="small",
calc_kwargs={"default_dtype": "float64"},
minimize_kwargs={"filter_func": None},
min_volume=0.75,
max_volume=1.25,
n_volumes=20,
).run()
/home/runner/work/janus-core/janus-core/.venv/lib/python3.12/site-packages/e3nn/o3/_wigner.py:10: UserWarning: Environment variable TORCH_FORCE_NO_WEIGHTS_ONLY_LOAD detected, since the`weights_only` argument was not explicitly passed to `torch.load`, forcing weights_only=False.
_Jd, _W3j_flat, _W3j_indices = torch.load(os.path.join(os.path.dirname(__file__), 'constants.pt'))
cuequivariance or cuequivariance_torch is not available. Cuequivariance acceleration will be disabled.
Downloading MACE model from 'https://github.com/ACEsuit/mace-mp/releases/download/mace_mp_0/2023-12-10-mace-128-L0_energy_epoch-249.model'
Cached MACE model to /home/runner/.cache/mace/20231210mace128L0_energy_epoch249model
Using Materials Project MACE for MACECalculator with /home/runner/.cache/mace/20231210mace128L0_energy_epoch249model
Using float64 for MACECalculator, which is slower but more accurate. Recommended for geometry optimization.
/home/runner/work/janus-core/janus-core/.venv/lib/python3.12/site-packages/mace/calculators/mace.py:139: UserWarning: Environment variable TORCH_FORCE_NO_WEIGHTS_ONLY_LOAD detected, since the`weights_only` argument was not explicitly passed to `torch.load`, forcing weights_only=False.
torch.load(f=model_path, map_location=device)
Step Time Energy fmax
LBFGS: 0 17:11:09 -71.020250 0.406340
LBFGS: 1 17:11:10 -71.032943 0.324324
LBFGS: 2 17:11:10 -71.057296 0.251886
LBFGS: 3 17:11:10 -71.064274 0.188059
LBFGS: 4 17:11:10 -71.066728 0.112945
LBFGS: 5 17:11:10 -71.070293 0.166075
LBFGS: 6 17:11:10 -71.074647 0.187184
LBFGS: 7 17:11:10 -71.078118 0.141462
LBFGS: 8 17:11:10 -71.079346 0.074487
[6]:
mace_eos["eos"].plot(show=True)
[6]:
<Axes: title={'center': 'birchmurnaghan: E: -71.091 eV, V: 114.353 Å$^3$, B: 185.324 GPa'}, xlabel='volume [Å$^3$]', ylabel='energy [eV]'>
Equation of state for β-quartz¶
Perform the same calculation for β-quartz:
[7]:
β_quartz = read("data/beta_quartz.cif")
v=WeasWidget()
v.from_ase(β_quartz)
v.avr.model_style = 1
v.avr.show_hydrogen_bonds = True
v
[7]:
[8]:
mace_eos_beta = EoS(
struct=β_quartz.copy(),
arch="mace_mp",
device="cpu",
model_path="small",
calc_kwargs={"default_dtype": "float64"},
minimize_kwargs={"filter_func": None},
min_volume=0.75,
max_volume=1.25,
n_volumes=20,
).run()
Using Materials Project MACE for MACECalculator with /home/runner/.cache/mace/20231210mace128L0_energy_epoch249model
Using float64 for MACECalculator, which is slower but more accurate. Recommended for geometry optimization.
/home/runner/work/janus-core/janus-core/.venv/lib/python3.12/site-packages/mace/calculators/mace.py:139: UserWarning: Environment variable TORCH_FORCE_NO_WEIGHTS_ONLY_LOAD detected, since the`weights_only` argument was not explicitly passed to `torch.load`, forcing weights_only=False.
torch.load(f=model_path, map_location=device)
Step Time Energy fmax
LBFGS: 0 17:11:12 -70.779216 0.490905
LBFGS: 1 17:11:12 -70.798248 0.413527
LBFGS: 2 17:11:12 -70.843094 0.020448
[9]:
mace_eos_beta["eos"].plot(show=True)
[9]:
<Axes: title={'center': 'birchmurnaghan: E: -71.024 eV, V: 124.056 Å$^3$, B: 172.905 GPa'}, xlabel='volume [Å$^3$]', ylabel='energy [eV]'>
Combining plots for α-quartz and β-quartz:
[10]:
import matplotlib.pyplot as plt
ax = plt.gca()
data_alpha = mace_eos["eos"].getplotdata()
data_beta = mace_eos_beta["eos"].getplotdata()
ax.plot(data_alpha[4], data_alpha[5], ls="-", color="C3", label="α-quartz")
ax.plot(data_alpha[6], data_alpha[7], ls="", marker="x", color="C4", mfc="C4")
ax.plot(data_beta[4], data_beta[5], ls="-", color="C0", label="β-quartz")
ax.plot(data_beta[6], data_beta[7], ls="", marker="x", color="C2", mfc="C2")
ax.set_xlabel("volume [Å$^3$]")
ax.set_ylabel("energy [eV]")
ax.legend()
plt.show()
Comparing MACE to CHGNeT and SevenNet¶
[11]:
sevennet_eos = EoS(
struct=α_quartz.copy(),
arch="sevennet",
device="cpu",
minimize_kwargs={"filter_func": None},
).run()
sevennet_eos["eos"].plot(show=True)
Step Time Energy fmax
LBFGS: 0 17:11:18 -71.035088 0.385548
LBFGS: 1 17:11:18 -71.046516 0.309889
LBFGS: 2 17:11:18 -71.069389 0.239566
LBFGS: 3 17:11:18 -71.075714 0.174659
LBFGS: 4 17:11:18 -71.077728 0.100454
LBFGS: 5 17:11:19 -71.080536 0.149301
LBFGS: 6 17:11:19 -71.084152 0.173767
LBFGS: 7 17:11:19 -71.087196 0.135812
LBFGS: 8 17:11:19 -71.088341 0.074211
[11]:
<Axes: title={'center': 'birchmurnaghan: E: -71.098 eV, V: 114.510 Å$^3$, B: 183.486 GPa'}, xlabel='volume [Å$^3$]', ylabel='energy [eV]'>
[12]:
chgnet_eos = EoS(
struct=α_quartz.copy(),
arch="chgnet",
device="cpu",
minimize_kwargs={"filter_func": None},
).run()
chgnet_eos["eos"].plot(show=True)
/home/runner/work/janus-core/janus-core/.venv/lib/python3.12/site-packages/nvidia_smi.py:810: SyntaxWarning: invalid escape sequence '\A'
mem = 'N\A'
/home/runner/work/janus-core/janus-core/.venv/lib/python3.12/site-packages/nvidia_smi.py:831: SyntaxWarning: invalid escape sequence '\A'
maxMemoryUsage = 'N\A'
CHGNet v0.3.0 initialized with 412,525 parameters
CHGNet will run on cpu
Step Time Energy fmax
LBFGS: 0 17:11:24 -74.947889 0.359846
LBFGS: 1 17:11:24 -74.958481 0.298901
LBFGS: 2 17:11:24 -74.981046 0.196990
LBFGS: 3 17:11:24 -74.985586 0.135577
LBFGS: 4 17:11:24 -74.986985 0.078268
[12]:
<Axes: title={'center': 'birchmurnaghan: E: -75.009 eV, V: 115.249 Å$^3$, B: 176.634 GPa'}, xlabel='volume [Å$^3$]', ylabel='energy [eV]'>
[13]:
print(f"MACE energy [eV]: {mace_eos['e_0']}")
print(f"SevenNet energy [eV]: {sevennet_eos['e_0']}")
print(f"CHGNeT energy [eV]: {chgnet_eos['e_0']}")
print()
print(f"MACE volume [Å^3]: {mace_eos['v_0']}")
print(f"SevenNet volume [Å^3]: {sevennet_eos['v_0']}")
print(f"CHGNeT volume [Å^3]: {chgnet_eos['v_0']}")
print()
print(f"MACE bulk_modulus [GPa]: {mace_eos['bulk_modulus']}")
print(f"SevenNet bulk_modulus [GPa]: {sevennet_eos['bulk_modulus']}")
print(f"CHGNeT bulk_modulus [GPa]: {chgnet_eos['bulk_modulus']}")
MACE energy [eV]: -71.09060046602642
SevenNet energy [eV]: -71.09814252551035
CHGNeT energy [eV]: -75.00900973428655
MACE volume [Å^3]: 114.35325408752064
SevenNet volume [Å^3]: 114.51025209091137
CHGNeT volume [Å^3]: 115.2490925284205
MACE bulk_modulus [GPa]: 185.3243187627527
SevenNet bulk_modulus [GPa]: 183.48614457865406
CHGNeT bulk_modulus [GPa]: 176.63355666500502