BAND Fragment Analysis Recipe¶
Use the PLAMS BAND fragment recipe to run two fragment calculations plus a full-system calculation and extract energy decomposition terms.
Initialization¶
from scm.base import ChemicalSystem
from scm.inputs import AMS, BAND
from scm.plams import view
from scm.plams.recipes.bandfragment import BANDFragmentJob
# build the surface
from ase.build import fcc111, add_adsorbate
Build Surface & Fragments¶
We first build a gold surface and add the hydrogen adsorbate.
mol = fcc111("Au", size=(2, 2, 3))
add_adsorbate(mol, "H", 1.5, "ontop")
mol.center(vacuum=10.0, axis=2)
view(ChemicalSystem.from_ase_atoms(mol), padding=-4, width=400)
The system is then split into two fragments based on the species.
surface = mol.copy()
symbols = surface.get_chemical_symbols()
del surface[[i for i in range(len(symbols)) if symbols[i] != "Au"]]
adsorbate = mol.copy()
del adsorbate[[i for i in range(len(symbols)) if symbols[i] == "Au"]]
If available, optimized fragments can also be loaded.
# from ase import io
# surface_opt = io.read("surface_opt.xyz")
# adsorbate_opt = io.read("adsorbate_opt.xyz")
# assert len(surface_opt) == len(surface)
# assert len(adsorbate_opt) == len(adsorbate)
Set Up & Run Job¶
For efficiency in this example, we use a minimal basis and reduced computational details to run the job.
base_settings = AMS(Task="SinglePoint", Engine=BAND())
base_settings.Engine.Basis.Type = "SZ"
base_settings.Engine.Basis.Core = "Large"
base_settings.Engine.DOS.CalcDOS = False
base_settings.Engine.KSpace.Regular.NumberOfPoints = [3, 3]
base_settings.Engine.BeckeGrid.Quality = "Basic"
base_settings.Engine.ZlmFit.Quality = "Basic"
base_settings.Engine.UseSymmetry = False
base_settings.Engine.XC.gga = "PBE"
base_settings.Engine.XC.dispersion = "Grimme4"
eda_settings = AMS(Engine=BAND(PEDA=True))
eda_job = BANDFragmentJob(
fragment1=ChemicalSystem.from_ase_atoms(surface),
fragment2=ChemicalSystem.from_ase_atoms(adsorbate),
settings=base_settings,
full_settings=eda_settings,
# fragment1_opt=ChemicalSystem.from_ase_atoms(surface_opt),
# fragment2_opt=ChemicalSystem.from_ase_atoms(adsorbate_opt),
)
eda_job.run();
[05.11|17:14:26] JOB plamsjob STARTED
[05.11|17:14:26] JOB plamsjob RUNNING
[05.11|17:14:26] JOB plamsjob/frag1 STARTED
[05.11|17:14:26] JOB plamsjob/frag1 RUNNING
[05.11|17:15:04] JOB plamsjob/frag1 FINISHED
[05.11|17:15:04] JOB plamsjob/frag1 SUCCESSFUL
[05.11|17:15:04] JOB plamsjob/frag2 STARTED
[05.11|17:15:04] JOB plamsjob/frag2 RUNNING
[05.11|17:15:07] JOB plamsjob/frag2 FINISHED
[05.11|17:15:07] JOB plamsjob/frag2 SUCCESSFUL
[05.11|17:15:07] JOB plamsjob/full STARTED
[05.11|17:15:07] JOB plamsjob/full RUNNING
[05.11|17:15:42] JOB plamsjob/full FINISHED
[05.11|17:15:42] JOB plamsjob/full SUCCESSFUL
[05.11|17:15:42] JOB plamsjob FINISHED
[05.11|17:15:42] JOB plamsjob SUCCESSFUL
Energy Decomposition Results¶
Finally, we extract the results of the energy decomposition:
import pandas as pd
results = eda_job.results
eda_res = eda_job.results.get_energy_decomposition(unit="eV")
df = pd.DataFrame([{"Term": t, "Energy [eV]": e} for t, e in eda_res.items()])
df
Term | Energy [eV] | |
|---|---|---|
0 | E_int | -1.886021 |
1 | E_int_disp | -0.094696 |
2 | E_Pauli | 11.635861 |
3 | E_elstat | -5.387854 |
4 | E_orb | -8.039060 |
5 | E_1 | -19.288663 |
6 | E_2 | -0.018088 |
See also¶
Python Script¶
#!/usr/bin/env python
# coding: utf-8
# ## Initialization
from scm.base import ChemicalSystem
from scm.inputs import AMS, BAND
from scm.plams import view
from scm.plams.recipes.bandfragment import BANDFragmentJob
# build the surface
from ase.build import fcc111, add_adsorbate
# ## Build Surface & Fragments
# We first build a gold surface and add the hydrogen adsorbate.
mol = fcc111("Au", size=(2, 2, 3))
add_adsorbate(mol, "H", 1.5, "ontop")
mol.center(vacuum=10.0, axis=2)
view(ChemicalSystem.from_ase_atoms(mol), padding=-4, width=400, picture_path="picture1.png")
# The system is then split into two fragments based on the species.
surface = mol.copy()
symbols = surface.get_chemical_symbols()
del surface[[i for i in range(len(symbols)) if symbols[i] != "Au"]]
adsorbate = mol.copy()
del adsorbate[[i for i in range(len(symbols)) if symbols[i] == "Au"]]
# If available, optimized fragments can also be loaded.
# from ase import io
# surface_opt = io.read("surface_opt.xyz")
# adsorbate_opt = io.read("adsorbate_opt.xyz")
# assert len(surface_opt) == len(surface)
# assert len(adsorbate_opt) == len(adsorbate)
# ## Set Up & Run Job
# For efficiency in this example, we use a minimal basis and reduced computational details to run the job.
base_settings = AMS(Task="SinglePoint", Engine=BAND())
base_settings.Engine.Basis.Type = "SZ"
base_settings.Engine.Basis.Core = "Large"
base_settings.Engine.DOS.CalcDOS = False
base_settings.Engine.KSpace.Regular.NumberOfPoints = [3, 3]
base_settings.Engine.BeckeGrid.Quality = "Basic"
base_settings.Engine.ZlmFit.Quality = "Basic"
base_settings.Engine.UseSymmetry = False
base_settings.Engine.XC.gga = "PBE"
base_settings.Engine.XC.dispersion = "Grimme4"
eda_settings = AMS(Engine=BAND(PEDA=True))
eda_job = BANDFragmentJob(
fragment1=ChemicalSystem.from_ase_atoms(surface),
fragment2=ChemicalSystem.from_ase_atoms(adsorbate),
settings=base_settings,
full_settings=eda_settings,
# fragment1_opt=ChemicalSystem.from_ase_atoms(surface_opt),
# fragment2_opt=ChemicalSystem.from_ase_atoms(adsorbate_opt),
)
eda_job.run();
# ## Energy Decomposition Results
# Finally, we extract the results of the energy decomposition:
import pandas as pd
results = eda_job.results
eda_res = eda_job.results.get_energy_decomposition(unit="eV")
df = pd.DataFrame([{"Term": t, "Energy [eV]": e} for t, e in eda_res.items()])
print(df)