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.

Downloads: Notebook | Script ?

Requires: AMS2027 or later

Related documentation
../_images/band_frag_4_0_031c89db.png

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)
image generated from notebook

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)