Adsorption Energies#
Pre-trained ODAC models are versatile across various MOF-related tasks. To begin, we’ll start with a fundamental application: calculating the adsorption energy for a single CO2 molecule. This serves as an excellent and simple demonstration of what you can achieve with these datasets and models.
For predicting the adsorption energy of a single CO2 molecule within a MOF structure, the adsorption energy (\(E_{\mathrm{ads}}\)) is defined as:
Each term on the right-hand side represents the energy of the relaxed state of the indicated chemical system. For a comprehensive understanding of our methodology for computing these adsorption energies, please refer to our paper.
Loading Pre-trained Models#
Need to install fairchem-core or get UMA access or getting permissions/401 errors?
Install the necessary packages using pip, uv etc
Get access to any necessary huggingface gated models
Get and login to your Huggingface account
Request access to https://huggingface.co/facebook/UMA
Create a Huggingface token at https://huggingface.co/settings/tokens/ with the permission “Permissions: Read access to contents of all public gated repos you can access”
Add the token as an environment variable using
huggingface-cli loginor by setting the HF_TOKEN environment variable.
A pre-trained model can be loaded using FAIRChemCalculator. In this example, we’ll employ UMA to determine the CO2 adsorption energies.
from fairchem.core import FAIRChemCalculator, pretrained_mlip
predictor = pretrained_mlip.get_predict_unit("uma-s-1")
calc = FAIRChemCalculator(predictor, task_name="odac")
/home/runner/work/_tool/Python/3.12.11/x64/lib/python3.12/site-packages/torchtnt/utils/version.py:12: UserWarning: pkg_resources is deprecated as an API. See https://setuptools.pypa.io/en/latest/pkg_resources.html. The pkg_resources package is slated for removal as early as 2025-11-30. Refrain from using this package or pin to Setuptools<81.
import pkg_resources
WARNING:root:device was not explicitly set, using device='cuda'.
Adsorption in rigid MOFs: CO2 Adsorption Energy in Mg-MOF-74#
Let’s apply our knowledge to Mg-MOF-74, a widely studied MOF known for its excellent CO2 adsorption properties. Its structure comprises magnesium atomic complexes connected by a carboxylated and oxidized benzene ring, serving as an organic linker. Previous studies consistently report the CO2 adsorption energy for Mg-MOF-74 to be around -0.40 eV [1] [2] [3].
Our goal is to verify if we can achieve a similar value by performing a simple single-point calculation using UMA. In the ODAC23 dataset, all MOF structures are identified by their CSD (Cambridge Structural Database) code. For Mg-MOF-74, this code is OPAGIX. We’ve extracted a specific OPAGIX+CO2 configuration from the dataset, which exhibits the lowest adsorption energy among its counterparts.
import matplotlib.pyplot as plt
from ase.io import read
from ase.visualize.plot import plot_atoms
mof_co2 = read("structures/OPAGIX_w_CO2.cif")
mof = read("structures/OPAGIX.cif")
co2 = read("structures/co2.xyz")
fig, ax = plt.subplots(figsize=(5, 4.5), dpi=250)
plot_atoms(mof_co2, ax)
ax.set_axis_off()
The final step in calculating the adsorption energy involves connecting the FAIRChemCalculator to each relaxed structure: OPAGIX+CO2, OPAGIX, and CO2. The structures used here are already relaxed from ODAC23. For simplicity, we assume here that further relaxations can be neglected. We will show how to go beyond this assumption in the next section.
mof_co2.calc = calc
mof.calc = calc
co2.calc = calc
E_ads = (
mof_co2.get_potential_energy()
- mof.get_potential_energy()
- co2.get_potential_energy()
)
print(f"Adsorption energy of CO2 in Mg-MOF-74: {E_ads:.3f} eV")
Adsorption energy of CO2 in Mg-MOF-74: -0.481 eV
Adsorption in flexible MOFs#
The adsorption energy calculation method outlined above is typically performed with rigid MOFs for simplicity. Both experimental and modeling literature have shown, however, that MOF flexibility can be important in accurately capturing the underlying chemistry of adsorption [1] [2] [3]. In particular, uptake can be improved by treating MOFs as flexible. Two types of MOF flexibility can be considered: intrinsic flexibility and deformation induced by guest molecules. In the Open DAC Project, we consider the latter MOF deformation by allowing the atomic positions of the MOF to relax during geometry optimization [4]. The addition of additional degrees of freedoms can complicate the computation of the adsorption energy and necessitates an extra step in the calculation procedure.
The figure below shows water adsorption in the MOF with CSD code WOBHEB with added defects (WOBHEB_0.11_0) from a DFT simulation. A typical adsorption energy calculation would only seek to capture the effects shaded in purple, which include both chemisorption and non-bonded interactions between the host and guest molecule. When allowing the MOF to relax, however, the adsorption energy also includes the energetic effect of the MOF deformation highlighted in green.
To account for this deformation, it is vital to use the most energetically favorable MOF geometry for the empty MOF term in Eqn. 1. Including MOF atomic coordinates as degrees of freedom can result in three possible outcomes:
The MOF does not deform, so the energies of the relaxed empty MOF and the MOF in the adsorbed state are the same
The MOF deforms to a less energetically favorable geometry than its ground state
The MOF locates a new energetically favorable geoemtry relative to the empty MOF relaxation
The first outcome requires no additional computation because the MOF rigidity assumption is valid. The second outcome represents physical and reversible deformation where the MOF returns to its empty ground state upon removal of the guest molecule. The third outcome is often the result of the guest molecule breaking local symmetry. We also found cases in ODAC in which both outcomes 2 and 3 occur within the same MOF.
To ensure the most energetically favorable empty MOF geometry is found, an addition empty MOF relaxation should be performed after MOF + adsorbate relaxation. The guest molecule should be removed, and the MOF should be relaxed starting from its geometry in the adsorbed state. If all deformation is reversible, the MOF will return to its original empty geometry. Otherwise, the lowest energy (most favorable) MOF geometry should be taken as the reference energy, \(E_{\mathrm{MOF}}\), in Eqn. 1.
H2O Adsorption Energy in Flexible WOBHEB with UMA#
The first part of this tutorial demonstrates how to perform a single point adsorption energy calculation using UMA. To treat MOFs as flexible, we perform all calculations on geometries determined by geometry optimization. The following example corresponds to the figure shown above (H2O adsorption in WOBHEB_0.11_0).
In this tutorial, \(E_{x}(r_{y})\) corresponds to the energy of \(x\) determined from geometry optimization of \(y\).
First, we obtain the energy of the empty MOF from relaxation of only the MOF: \(E_{\mathrm{MOF}}(r_{\mathrm{MOF}})\)
import ase.io
from ase.optimize import BFGS
mof = ase.io.read("structures/WOBHEB_0.11.cif")
mof.calc = calc
relax = BFGS(mof)
relax.run(fmax=0.05)
E_mof_empty = mof.get_potential_energy()
print(f"Energy of empty MOF: {E_mof_empty:.3f} eV")
Step Time Energy fmax
BFGS: 0 20:20:47 -1077.292857 0.165906
BFGS: 1 20:20:47 -1077.295713 0.138830
BFGS: 2 20:20:47 -1077.302078 0.137076
BFGS: 3 20:20:48 -1077.304994 0.122721
BFGS: 4 20:20:48 -1077.309243 0.101695
BFGS: 5 20:20:48 -1077.311390 0.083808
BFGS: 6 20:20:48 -1077.313546 0.078572
BFGS: 7 20:20:48 -1077.315553 0.090393
BFGS: 8 20:20:49 -1077.317810 0.099523
BFGS: 9 20:20:49 -1077.319836 0.084360
BFGS: 10 20:20:49 -1077.321704 0.076825
BFGS: 11 20:20:49 -1077.323775 0.092502
BFGS: 12 20:20:50 -1077.326227 0.083269
BFGS: 13 20:20:50 -1077.328768 0.080299
BFGS: 14 20:20:50 -1077.331000 0.083239
BFGS: 15 20:20:50 -1077.332912 0.077169
BFGS: 16 20:20:51 -1077.334902 0.082378
BFGS: 17 20:20:51 -1077.337360 0.094597
BFGS: 18 20:20:51 -1077.340262 0.083951
BFGS: 19 20:20:51 -1077.343059 0.092127
BFGS: 20 20:20:52 -1077.345413 0.083144
BFGS: 21 20:20:52 -1077.347398 0.075780
BFGS: 22 20:20:52 -1077.349435 0.079734
BFGS: 23 20:20:52 -1077.351919 0.095029
BFGS: 24 20:20:53 -1077.354742 0.127323
BFGS: 25 20:20:53 -1077.357638 0.110035
BFGS: 26 20:20:53 -1077.360337 0.076491
BFGS: 27 20:20:53 -1077.362893 0.072993
BFGS: 28 20:20:54 -1077.365629 0.118983
BFGS: 29 20:20:54 -1077.368700 0.136169
BFGS: 30 20:20:54 -1077.371404 0.086550
BFGS: 31 20:20:54 -1077.373486 0.060053
BFGS: 32 20:20:54 -1077.375340 0.076341
BFGS: 33 20:20:55 -1077.377237 0.106647
BFGS: 34 20:20:55 -1077.379239 0.098119
BFGS: 35 20:20:55 -1077.381107 0.078262
BFGS: 36 20:20:55 -1077.382691 0.068215
BFGS: 37 20:20:56 -1077.384179 0.095465
BFGS: 38 20:20:56 -1077.385719 0.105318
BFGS: 39 20:20:56 -1077.387182 0.065271
BFGS: 40 20:20:56 -1077.388445 0.061075
BFGS: 41 20:20:57 -1077.389693 0.063910
BFGS: 42 20:20:57 -1077.390921 0.081666
BFGS: 43 20:20:57 -1077.392084 0.072794
BFGS: 44 20:20:57 -1077.393100 0.037364
Energy of empty MOF: -1077.393 eV
Next, we add the H2O guest molecule and relax the MOF + adsorbate to obtain \(E_{\mathrm{MOF+H2O}}(r_{\mathrm{MOF+H2O}})\).
mof_h2o = ase.io.read("structures/WOBHEB_H2O.cif")
mof_h2o.calc = calc
relax = BFGS(mof_h2o)
relax.run(fmax=0.05)
E_combo = mof_h2o.get_potential_energy()
print(f"Energy of MOF + H2O: {E_combo:.3f} eV")
Step Time Energy fmax
BFGS: 0 20:20:58 -1091.588456 1.130484
BFGS: 1 20:20:58 -1091.607674 0.292276
BFGS: 2 20:20:58 -1091.612768 0.244457
BFGS: 3 20:20:58 -1091.631192 0.232970
BFGS: 4 20:20:59 -1091.637789 0.244934
BFGS: 5 20:20:59 -1091.648793 0.185920
BFGS: 6 20:20:59 -1091.655536 0.183822
BFGS: 7 20:20:59 -1091.663495 0.164533
BFGS: 8 20:21:00 -1091.671228 0.176775
BFGS: 9 20:21:00 -1091.679115 0.167616
BFGS: 10 20:21:00 -1091.686678 0.182172
BFGS: 11 20:21:00 -1091.694618 0.176638
BFGS: 12 20:21:01 -1091.704442 0.192324
BFGS: 13 20:21:01 -1091.715585 0.175239
BFGS: 14 20:21:01 -1091.726366 0.136995
BFGS: 15 20:21:01 -1091.735699 0.162123
BFGS: 16 20:21:02 -1091.744603 0.161332
BFGS: 17 20:21:02 -1091.754438 0.176043
BFGS: 18 20:21:02 -1091.764819 0.152776
BFGS: 19 20:21:02 -1091.774011 0.117210
BFGS: 20 20:21:03 -1091.781095 0.153273
BFGS: 21 20:21:03 -1091.787279 0.168392
BFGS: 22 20:21:03 -1091.793470 0.147293
BFGS: 23 20:21:03 -1091.800675 0.158254
BFGS: 24 20:21:04 -1091.808516 0.146093
BFGS: 25 20:21:04 -1091.815110 0.208642
BFGS: 26 20:21:04 -1091.822059 0.115493
BFGS: 27 20:21:04 -1091.827561 0.229688
BFGS: 28 20:21:05 -1091.833614 0.134056
BFGS: 29 20:21:05 -1091.838681 0.142427
BFGS: 30 20:21:05 -1091.844496 0.090061
BFGS: 31 20:21:05 -1091.847854 0.143851
BFGS: 32 20:21:05 -1091.852878 0.123032
BFGS: 33 20:21:06 -1091.856760 0.300291
BFGS: 34 20:21:06 -1091.861755 0.128652
BFGS: 35 20:21:06 -1091.866178 0.161280
BFGS: 36 20:21:06 -1091.872053 0.157821
BFGS: 37 20:21:07 -1091.877129 0.207711
BFGS: 38 20:21:07 -1091.885248 0.168115
BFGS: 39 20:21:07 -1091.889387 0.582541
BFGS: 40 20:21:07 -1091.896758 0.149674
BFGS: 41 20:21:08 -1091.901791 0.121329
BFGS: 42 20:21:08 -1091.914629 0.500383
BFGS: 43 20:21:08 -1091.920202 0.166124
BFGS: 44 20:21:08 -1091.930444 0.239970
BFGS: 45 20:21:09 -1091.942839 0.420396
BFGS: 46 20:21:09 -1091.951166 0.162663
BFGS: 47 20:21:09 -1091.964936 0.540581
BFGS: 48 20:21:09 -1091.985874 0.344845
BFGS: 49 20:21:10 -1092.006594 0.383737
BFGS: 50 20:21:10 -1092.025834 0.203954
BFGS: 51 20:21:10 -1092.027268 1.287413
BFGS: 52 20:21:10 -1092.063286 0.285507
BFGS: 53 20:21:11 -1092.084171 0.255845
BFGS: 54 20:21:11 -1092.101609 0.548366
BFGS: 55 20:21:11 -1092.125817 0.305899
BFGS: 56 20:21:11 -1092.139135 0.236076
BFGS: 57 20:21:12 -1092.154027 0.223446
BFGS: 58 20:21:12 -1092.164332 0.264231
BFGS: 59 20:21:12 -1092.178087 0.307265
BFGS: 60 20:21:12 -1092.187986 0.371080
BFGS: 61 20:21:13 -1092.201334 0.336026
BFGS: 62 20:21:13 -1092.211310 0.301372
BFGS: 63 20:21:13 -1092.221330 0.222567
BFGS: 64 20:21:13 -1092.230946 0.149221
BFGS: 65 20:21:14 -1092.239634 0.116998
BFGS: 66 20:21:14 -1092.246706 0.101949
BFGS: 67 20:21:14 -1092.252686 0.116666
BFGS: 68 20:21:14 -1092.258605 0.125347
BFGS: 69 20:21:15 -1092.264686 0.142368
BFGS: 70 20:21:15 -1092.269889 0.118685
BFGS: 71 20:21:15 -1092.274174 0.080362
BFGS: 72 20:21:15 -1092.278093 0.125676
BFGS: 73 20:21:16 -1092.281937 0.169223
BFGS: 74 20:21:16 -1092.285520 0.180595
BFGS: 75 20:21:16 -1092.288644 0.132724
BFGS: 76 20:21:16 -1092.291663 0.091648
BFGS: 77 20:21:17 -1092.294983 0.086504
BFGS: 78 20:21:17 -1092.298289 0.085589
BFGS: 79 20:21:17 -1092.301025 0.077517
BFGS: 80 20:21:17 -1092.303246 0.064986
BFGS: 81 20:21:17 -1092.305127 0.084747
BFGS: 82 20:21:18 -1092.306707 0.082859
BFGS: 83 20:21:18 -1092.308032 0.075575
BFGS: 84 20:21:18 -1092.309340 0.063983
BFGS: 85 20:21:18 -1092.310827 0.090306
BFGS: 86 20:21:19 -1092.312496 0.097845
BFGS: 87 20:21:19 -1092.314104 0.095230
BFGS: 88 20:21:19 -1092.315532 0.058497
BFGS: 89 20:21:19 -1092.316911 0.042171
Energy of MOF + H2O: -1092.317 eV
We can now isolate the MOF atoms from the relaxed MOF + H2O geometry and see that the MOF has adopted a geometry that is less energetically favorable than the empty MOF by ~0.2 eV. The energy of the MOF in the adsorbed state corresponds to \(E_{\mathrm{MOF}}(r_{\mathrm{MOF+H2O}})\).
mof_adsorbed_state = mof_h2o[:-3]
mof_adsorbed_state.calc = calc
E_mof_adsorbed_state = mof_adsorbed_state.get_potential_energy()
print(f"Energy of MOF in the adsorbed state: {E_mof_adsorbed_state:.3f} eV")
Energy of MOF in the adsorbed state: -1077.116 eV
H2O adsorption in this MOF appears to correspond to Case #2 as outlined above. We can now perform re-relaxation of the empty MOF starting from the \(r_{\mathrm{MOF+H2O}}\) geometry.
relax = BFGS(mof_adsorbed_state)
relax.run(fmax=0.05)
E_mof_rerelax = mof_adsorbed_state.get_potential_energy()
print(f"Energy of re-relaxed empty MOF: {E_mof_rerelax:.3f} eV")
Step Time Energy fmax
BFGS: 0 20:21:20 -1077.115541 0.933867
BFGS: 1 20:21:20 -1077.146944 0.841690
BFGS: 2 20:21:20 -1077.193260 0.831750
BFGS: 3 20:21:20 -1077.233594 0.573527
BFGS: 4 20:21:21 -1077.253517 0.357214
BFGS: 5 20:21:21 -1077.268459 0.261913
BFGS: 6 20:21:21 -1077.280001 0.254889
BFGS: 7 20:21:21 -1077.290050 0.281306
BFGS: 8 20:21:22 -1077.299127 0.249011
BFGS: 9 20:21:22 -1077.307287 0.180263
BFGS: 10 20:21:22 -1077.312033 0.149457
BFGS: 11 20:21:22 -1077.316245 0.166602
BFGS: 12 20:21:23 -1077.319909 0.178935
BFGS: 13 20:21:23 -1077.325187 0.197543
BFGS: 14 20:21:23 -1077.329561 0.187120
BFGS: 15 20:21:23 -1077.334411 0.155985
BFGS: 16 20:21:23 -1077.338673 0.186094
BFGS: 17 20:21:24 -1077.343280 0.172124
BFGS: 18 20:21:24 -1077.347679 0.153904
BFGS: 19 20:21:24 -1077.351017 0.108817
BFGS: 20 20:21:24 -1077.353406 0.108391
BFGS: 21 20:21:25 -1077.355602 0.105608
BFGS: 22 20:21:25 -1077.357813 0.135340
BFGS: 23 20:21:25 -1077.360345 0.132757
BFGS: 24 20:21:25 -1077.362868 0.104287
BFGS: 25 20:21:26 -1077.365160 0.095923
BFGS: 26 20:21:26 -1077.367029 0.075546
BFGS: 27 20:21:26 -1077.368624 0.074004
BFGS: 28 20:21:26 -1077.370274 0.091634
BFGS: 29 20:21:27 -1077.371997 0.100105
BFGS: 30 20:21:27 -1077.373641 0.066669
BFGS: 31 20:21:27 -1077.374806 0.047987
Energy of re-relaxed empty MOF: -1077.375 eV
The MOF returns to its original empty reference energy upon re-relaxation, confirming that this deformation is physically relevant and is induced by the adsorbate molecule. In Case #3, this re-relaxed energy will be more negative (more favorable) than the original empty MOF relaxation. Thus, we take the reference empty MOF energy (\(E_{\mathrm{MOF}}\) in Eqn. 1) to be the minimum of the original empty MOF energy and the re-relaxed MOf energy:
E_mof = min(E_mof_empty, E_mof_rerelax)
# get adsorbate reference energy
h2o = mof_h2o[-3:]
h2o.calc = calc
E_h2o = h2o.get_potential_energy()
# compute adsorption energy
E_ads = E_combo - E_mof - E_h2o
print(f"Adsorption energy of H2O in WOBHEB_0.11_0: {E_ads:.3f} eV")
Adsorption energy of H2O in WOBHEB_0.11_0: -0.701 eV
This adsorption energy closely matches that from DFT (–0.699 eV) [1]. The strong adsorption energy is a consequence of both H2O chemisorption and MOF deformation. We can decompose the adsorption energy into contributions from these two factors. Assuming rigid H2O molecules, we define \(E_{\mathrm{int}}\) and \(E_{\mathrm{MOF,deform}}\), respectively, as
\(E_{\mathrm{int}}\) describes host host–guest interactions for the MOF in the adsorbed state only. \(E_{\mathrm{MOF,deform}}\) quantifies the magnitude of deformation between the MOF in the adsorbed state and the most energetically favorable empty MOF geometry determined from the workflow presented here. It can be shown that
For H2O adsorption in WOBHEB_0.11, we have
E_int = E_combo - E_mof_adsorbed_state - E_h2o
print(f"E_int: {E_int}")
E_int: -0.9784625172610912
E_mof_deform = E_mof_adsorbed_state - E_mof_empty
print(f"E_mof_deform: {E_mof_deform}")
E_mof_deform: 0.2775588035583496
E_ads = E_int + E_mof_deform
print(f"E_ads: {E_ads}")
E_ads: -0.7009037137027416
\(E_{\mathrm{int}}\) is equivalent to \(E_{\mathrm{ads}}\) when the MOF is assumed to be rigid. In this case, failure to consider adsorbate-induced deformation would result in an overestimation of the adsorption energy magnitude.