ASEによる表面モデルの作成方法

fcc111など、典型的な表面の作成

from ase.build import fcc111
surf = fcc111("Cu", size=[3,3,5], a=3.2, vacuum=10.0)

fcc100, fcc110, bcc100, bcc110, bcc111, fcc111, hcp0001が使える

一般的な表面の作成

from ase.build import bulk, surface
bulk = bulk("Cu", "fcc", a=3.2)
surf = surface(bulk, (1,1,1), 5, vacuum=10.0)  # lattice, Miller indices, layers
surf = surf*[3, 3, 1]  # if you want to replicate
surf = surf.translate([0, 0, -10]) # if you want to shift

CIF/POSCARをロードして表面を作る例

from ase import Atoms, Atom
from ase.build import surface
from ase.visualize import view
from ase.io import read
from ase.io import write
bulk = read("1000053.cif") # or read("oldPOSCAR")
surf = surface(lattice=bulk, indices=(1,1,1), layers=5, vacuum=10.0)
surf.translate([0,0,-10.0+0.1]) # shift slightly for easier to view
view(surf*[5,5,1])
write("POSCAR", surf)

原子・分子を表面に吸着させる

add_adsorbateを使う

from ase import Atoms
from ase.build import fcc111, add_adsorbate

surf = fcc111("Cu", size=[3,3,6], vacuum=10.0)
mol  = Atoms("CO", [(0,0,0),(0,0,1.2)])
add_adsorbate(surf, mol, 1.6, position="ontop")
# slab, adsorbate, height, position=(0,0), offset=None, mol_index=0)

offsetはfloatでよい(e.g. offset=(1, 1), offset=(0.5, 0))。unit cellのfractional coordinateになっている
メジャーなsurfaceに対してはstringでの吸着位置がある

surface	available adsorption sites
fcc100	ontop, hollow, bridge
fcc110	ontop, hollow, longbridge, shortbridge
bcc100	ontop, hollow, bridge
fcc111	ontop, bridge, fcc, hcp
hcp0001	ontop, bridge, fcc, hcp
bcc110	ontop, shortbridge, longbridge, hollow
bcc111	ontop, hollow

表面の一部を固定した計算

Atomsのconstraintを設定することでPOSCARに原子固定情報が追加される
ASEで表面を構築した場合、tagが自動的に設定されるのでこれとase.constraint.FixAtomsを用いてconstraint条件を設定する(下記スクリプトのconstraints)
- adsorbateはtag=0, 表面層はtag=1など。深いほどtagの値が増える
スクリプトの例

from ase.calculators.vasp import Vasp
from ase.optimize import QuasiNewton
from ase.constraints import FixAtoms
from ase.build import fcc111

surf = fcc111(symbol='Pd', size=[2, 2, 4], a=3.9, vacuum=10.0)
constraints = FixAtoms(indices=[atom.index for atom in surf if atom.tag >= 3])
surf.constraints = constraints

# using ASE optimizer
calc = Vasp(prec='Normal', xc='PBE', encut=400.0,
            ibrion=-1, nsw=0, kpts=[2, 2, 1], potim=0.5)
surf.set_calculator(calc)
dyn = QuasiNewton(surf, trajectory='pd-relax.traj')
dyn.run(fmax=0.05)

energy = surf.get_potential_energy()
printf("Energy after relaxation: {energy} eV")