Slab + mol energy¶
触媒上での反応などでは、
触媒表面上に有機分子が吸着
表面上で反応が進行
反応終了後、生成物が触媒表面上から脱離
といった過程を扱います。 反応物・生成物が触媒上で吸着・脱離をするため、吸着構造(表面構造上に分子が存在する構造)を扱う必要があります。
本節では、吸着構造に関連するエネルギーとして、吸着エネルギーの計算を取り扱います。
吸着エネルギー
[1]:
from ase import Atoms
from ase.build import bulk
from ase.constraints import ExpCellFilter, StrainFilter
from ase.optimize import LBFGS, FIRE
import pfp_api_client
from pfp_api_client.pfp.calculators.ase_calculator import ASECalculator
from pfp_api_client.pfp.estimator import Estimator, EstimatorCalcMode
print(f"pfp_api_client: {pfp_api_client.__version__}")
estimator = Estimator(calc_mode=EstimatorCalcMode.CRYSTAL_U0, model_version="v3.0.0")
calculator = ASECalculator(estimator)
pfp_api_client: 1.19.0
吸着エネルギー¶
吸着エネルギー \(E_{\rm{adsorption}}\) は、ある分子が表面構造に吸着する前後のエネルギー差を表します。
\[E_{\rm{adsorption}} = E_{\rm{slab+mol}} - (E_{\rm{slab}} + E_{\rm{mol}})\]
\(E_{\rm{slab}}\)がSlab単体のエネルギー、\(E_{\rm{mol}}\)が分子単体のエネルギー、\(E_{\rm{slab+mol}}\)がSlab上に分子がついた吸着構造のエネルギーです。
CO on Pt¶
[2]:
from ase.optimize import LBFGS
from ase.build import fcc111, molecule, add_adsorbate
from ase.constraints import ExpCellFilter, StrainFilter
def get_opt_energy(atoms, fmax=0.001, opt_mode: str = "normal"):
atoms.set_calculator(calculator)
if opt_mode == "scale":
opt1 = LBFGS(StrainFilter(atoms, mask=[1, 1, 1, 0, 0, 0]))
elif opt_mode == "all":
opt1 = LBFGS(ExpCellFilter(atoms))
else:
opt1 = LBFGS(atoms)
opt1.run(fmax=fmax)
return atoms.get_total_energy()
まずは適切なPt のセルサイズを求めます。
Bulkを用意して最適化することでセルサイズを求めます。
[3]:
bulk_atoms = bulk("Pt", cubic=True)
bulk_atoms.cell
[3]:
Cell([3.92, 3.92, 3.92])
[4]:
bulk_atoms = bulk("Pt", cubic=True)
bulk_atoms.calc = calculator
E_bulk = get_opt_energy(bulk_atoms, fmax=1e-4, opt_mode="scale")
E_bulk
Step Time Energy fmax
LBFGS: 0 07:54:49 -21.749341 6.135577
LBFGS: 1 07:54:49 -21.296372 15.553638
/tmp/ipykernel_8561/1510506507.py:7: DeprecationWarning: Please use atoms.calc = calc
atoms.set_calculator(calculator)
/tmp/ipykernel_8561/1510506507.py:9: FutureWarning: Import StrainFilter from ase.filters
opt1 = LBFGS(StrainFilter(atoms, mask=[1, 1, 1, 0, 0, 0]))
LBFGS: 2 07:54:49 -21.809599 1.021796
LBFGS: 3 07:54:49 -21.811497 0.156920
LBFGS: 4 07:54:49 -21.811544 0.002066
LBFGS: 5 07:54:49 -21.811544 0.000010
[4]:
-21.811543944246235
[5]:
bulk_atoms.cell
[5]:
Cell([3.9674289205307174, 3.967428532149144, 3.967428615755973])
[6]:
import numpy as np
a = np.mean(np.diag(bulk_atoms.cell))
a
[6]:
3.967428689478611
ASE デフォルトの3.92A からは少し値が変わっていることがわかります。このセルサイズを用いて以下、吸着構造を作成していきます。
[7]:
slab = fcc111("Pt", a=a, size=(4, 4, 4), vacuum=40.0, periodic=True)
slab.calc = calculator
[8]:
E_slab = get_opt_energy(slab, fmax=1e-4, opt_mode="normal")
Step Time Energy fmax
LBFGS: 0 07:54:49 -329.391234 0.322358
LBFGS: 1 07:54:49 -329.444709 0.243355
/tmp/ipykernel_8561/1510506507.py:7: DeprecationWarning: Please use atoms.calc = calc
atoms.set_calculator(calculator)
LBFGS: 2 07:54:50 -329.511078 0.043360
LBFGS: 3 07:54:50 -329.512193 0.022307
LBFGS: 4 07:54:50 -329.512651 0.022152
LBFGS: 5 07:54:50 -329.514980 0.010735
LBFGS: 6 07:54:50 -329.515248 0.005379
LBFGS: 7 07:54:50 -329.515195 0.000678
LBFGS: 8 07:54:50 -329.515152 0.000035
[9]:
mol = molecule("CO")
mol.calc = calculator
E_mol = get_opt_energy(mol, fmax=1e-4)
Step Time Energy fmax
LBFGS: 0 07:54:50 -11.497234 0.820163
LBFGS: 1 07:54:50 -11.484524 1.958990
LBFGS: 2 07:54:50 -11.500194 0.029776
/tmp/ipykernel_8561/1510506507.py:7: DeprecationWarning: Please use atoms.calc = calc
atoms.set_calculator(calculator)
LBFGS: 3 07:54:50 -11.500197 0.001082
LBFGS: 4 07:54:50 -11.500199 0.000007
吸着構造の作成¶
吸着構造は add_adsorbate 関数を用いることで、Slab上に分子を載せる事ができます。
[10]:
%%time
slab_mol = fcc111("Pt", a=a, size=(4, 4, 4), vacuum=40.0)
slab_mol.pbc = True
mol2 = molecule("CO")
add_adsorbate(slab_mol, mol2, 3.0, "bridge")
E_slab_mol = get_opt_energy(slab_mol, fmax=1e-3)
Step Time Energy fmax
LBFGS: 0 07:54:50 -341.994046 3.707106
LBFGS: 1 07:54:51 -342.131358 4.426330
/tmp/ipykernel_8561/1510506507.py:7: DeprecationWarning: Please use atoms.calc = calc
atoms.set_calculator(calculator)
LBFGS: 2 07:54:51 -342.429021 3.370984
LBFGS: 3 07:54:51 -342.724384 0.943089
LBFGS: 4 07:54:51 -342.800038 0.793622
LBFGS: 5 07:54:51 -342.828912 0.437039
LBFGS: 6 07:54:51 -342.872906 0.237297
LBFGS: 7 07:54:51 -342.876484 0.111745
LBFGS: 8 07:54:51 -342.880215 0.075044
LBFGS: 9 07:54:51 -342.883296 0.095750
LBFGS: 10 07:54:51 -342.885186 0.100860
LBFGS: 11 07:54:51 -342.886260 0.065886
LBFGS: 12 07:54:51 -342.887206 0.042198
LBFGS: 13 07:54:51 -342.888046 0.035738
LBFGS: 14 07:54:51 -342.888617 0.043613
LBFGS: 15 07:54:52 -342.888869 0.031778
LBFGS: 16 07:54:52 -342.889159 0.022924
LBFGS: 17 07:54:52 -342.889444 0.022108
LBFGS: 18 07:54:52 -342.889710 0.026686
LBFGS: 19 07:54:52 -342.889830 0.021686
LBFGS: 20 07:54:52 -342.889972 0.017799
LBFGS: 21 07:54:52 -342.890118 0.016525
LBFGS: 22 07:54:52 -342.890202 0.016837
LBFGS: 23 07:54:52 -342.890292 0.011782
LBFGS: 24 07:54:52 -342.890313 0.010002
LBFGS: 25 07:54:52 -342.890382 0.010341
LBFGS: 26 07:54:52 -342.890426 0.012521
LBFGS: 27 07:54:53 -342.890447 0.009512
LBFGS: 28 07:54:53 -342.890484 0.007368
LBFGS: 29 07:54:53 -342.890478 0.007693
LBFGS: 30 07:54:53 -342.890493 0.008405
LBFGS: 31 07:54:53 -342.890537 0.007140
LBFGS: 32 07:54:53 -342.890515 0.006148
LBFGS: 33 07:54:53 -342.890559 0.005694
LBFGS: 34 07:54:53 -342.890563 0.005676
LBFGS: 35 07:54:53 -342.890609 0.005642
LBFGS: 36 07:54:53 -342.890605 0.005324
LBFGS: 37 07:54:53 -342.890591 0.005027
LBFGS: 38 07:54:53 -342.890633 0.004830
LBFGS: 39 07:54:53 -342.890661 0.004811
LBFGS: 40 07:54:54 -342.890648 0.004381
LBFGS: 41 07:54:54 -342.890633 0.004414
LBFGS: 42 07:54:54 -342.890678 0.004662
LBFGS: 43 07:54:54 -342.890645 0.005618
LBFGS: 44 07:54:54 -342.890670 0.004685
LBFGS: 45 07:54:54 -342.890705 0.003669
LBFGS: 46 07:54:54 -342.890699 0.003658
LBFGS: 47 07:54:54 -342.890704 0.005095
LBFGS: 48 07:54:54 -342.890697 0.004729
LBFGS: 49 07:54:54 -342.890707 0.003384
LBFGS: 50 07:54:54 -342.890703 0.003407
LBFGS: 51 07:54:54 -342.890725 0.003284
LBFGS: 52 07:54:54 -342.890736 0.003033
LBFGS: 53 07:54:54 -342.890737 0.002020
LBFGS: 54 07:54:55 -342.890725 0.001996
LBFGS: 55 07:54:55 -342.890712 0.001847
LBFGS: 56 07:54:55 -342.890735 0.001679
LBFGS: 57 07:54:55 -342.890737 0.000997
CPU times: user 377 ms, sys: 31.7 ms, total: 408 ms
Wall time: 4.37 s
[11]:
E_adsorp = E_slab_mol - (E_slab + E_mol)
print(f"E_adsorp {E_adsorp}, E_slab_mol {E_slab_mol}, E_slab {E_slab}, E_mol {E_mol}")
E_adsorp -1.8753863462617915, E_slab_mol -342.89073734826064, E_slab -329.5151524229501, E_mol -11.500198579048762
吸着エネルギーとして-1.87 eVが得られました。
[12]:
from pfcc_extras.visualize.view import view_ngl
view_ngl(slab_mol)
吸着サイトについて¶
add_adsorbate 関数の4番目に指定する position 引数として、座標値の代わりに吸着サイト名で指定することができます。
[13]:
%%time
slab_mol = fcc111("Pt", a=a, size=(4, 4, 4), vacuum=40.0)
slab_mol.pbc = True
mol2 = molecule("CO")
add_adsorbate(slab_mol, mol2, 3.0, "ontop")
E_slab_mol = get_opt_energy(slab_mol, fmax=1e-4)
Step Time Energy fmax
LBFGS: 0 07:55:06 -342.399910 0.902054
LBFGS: 1 07:55:06 -342.463817 1.697185
/tmp/ipykernel_8561/1510506507.py:7: DeprecationWarning: Please use atoms.calc = calc
atoms.set_calculator(calculator)
LBFGS: 2 07:55:06 -342.595038 0.461127
LBFGS: 3 07:55:06 -342.609559 0.374306
LBFGS: 4 07:55:06 -342.655779 0.619674
LBFGS: 5 07:55:06 -342.671951 0.635450
LBFGS: 6 07:55:06 -342.685604 0.410959
LBFGS: 7 07:55:06 -342.695938 0.185124
LBFGS: 8 07:55:06 -342.707631 0.303364
LBFGS: 9 07:55:06 -342.720453 0.393410
LBFGS: 10 07:55:06 -342.732852 0.282378
LBFGS: 11 07:55:07 -342.739609 0.113854
LBFGS: 12 07:55:07 -342.744490 0.145889
LBFGS: 13 07:55:07 -342.748745 0.254531
LBFGS: 14 07:55:07 -342.751913 0.220942
LBFGS: 15 07:55:07 -342.754120 0.085357
LBFGS: 16 07:55:07 -342.755419 0.053168
LBFGS: 17 07:55:07 -342.756486 0.118864
LBFGS: 18 07:55:07 -342.757801 0.145427
LBFGS: 19 07:55:07 -342.758954 0.095398
LBFGS: 20 07:55:07 -342.759748 0.036356
LBFGS: 21 07:55:07 -342.760303 0.043884
LBFGS: 22 07:55:08 -342.760789 0.069188
LBFGS: 23 07:55:08 -342.761185 0.053738
LBFGS: 24 07:55:08 -342.761434 0.028786
LBFGS: 25 07:55:08 -342.761635 0.022932
LBFGS: 26 07:55:08 -342.761780 0.031662
LBFGS: 27 07:55:08 -342.761886 0.023381
LBFGS: 28 07:55:08 -342.761952 0.010756
LBFGS: 29 07:55:08 -342.762002 0.011717
LBFGS: 30 07:55:08 -342.762022 0.017857
LBFGS: 31 07:55:08 -342.762073 0.012374
LBFGS: 32 07:55:08 -342.762078 0.005001
LBFGS: 33 07:55:08 -342.762066 0.004472
LBFGS: 34 07:55:08 -342.762079 0.006191
LBFGS: 35 07:55:09 -342.762057 0.003670
LBFGS: 36 07:55:09 -342.762067 0.001881
LBFGS: 37 07:55:09 -342.762076 0.002791
LBFGS: 38 07:55:09 -342.762088 0.003626
LBFGS: 39 07:55:09 -342.762083 0.002242
LBFGS: 40 07:55:09 -342.762098 0.001410
LBFGS: 41 07:55:09 -342.762075 0.001218
LBFGS: 42 07:55:09 -342.762077 0.001632
LBFGS: 43 07:55:09 -342.762120 0.001008
LBFGS: 44 07:55:09 -342.762066 0.000690
LBFGS: 45 07:55:09 -342.762096 0.000912
LBFGS: 46 07:55:10 -342.762103 0.001544
LBFGS: 47 07:55:10 -342.762073 0.001479
LBFGS: 48 07:55:10 -342.762093 0.000675
LBFGS: 49 07:55:10 -342.762060 0.000407
LBFGS: 50 07:55:10 -342.762065 0.000547
LBFGS: 51 07:55:10 -342.762073 0.000606
LBFGS: 52 07:55:10 -342.762052 0.000299
LBFGS: 53 07:55:10 -342.762062 0.000176
LBFGS: 54 07:55:10 -342.762070 0.000289
LBFGS: 55 07:55:10 -342.762041 0.000301
LBFGS: 56 07:55:10 -342.762044 0.000155
LBFGS: 57 07:55:10 -342.762075 0.000072
CPU times: user 381 ms, sys: 47.8 ms, total: 429 ms
Wall time: 4.87 s
この吸着サイトの名前は、 .info を表示することで “sites” 以下に定義されている値として確認ができます。
[14]:
slab_mol.info
[14]:
{'adsorbate_info': {'cell': array([[2.80539573, 0. ],
[1.40269787, 2.42954397]]),
'sites': {'ontop': (0, 0),
'bridge': (0.5, 0),
'fcc': (0.3333333333333333, 0.3333333333333333),
'hcp': (0.6666666666666666, 0.6666666666666666)},
'top layer atom index': 48}}
今回の例fcc111だと、 ontop, bridge, fcc, hcp があるようです。それぞれ実際に見てみましょう。
[15]:
orthogonal = False
slab_mol_ontop = fcc111("Pt", a=a, size=(4, 4, 4), vacuum=40.0, orthogonal=orthogonal, periodic=True)
mol = molecule("CO")
add_adsorbate(slab_mol_ontop, mol, 3.0, "ontop")
slab_mol_bridge = fcc111("Pt", a=a, size=(4, 4, 4), vacuum=40.0, orthogonal=orthogonal, periodic=True)
mol = molecule("CO")
add_adsorbate(slab_mol_bridge, mol, 3.0, "bridge")
slab_mol_fcc = fcc111("Pt", a=a, size=(4, 4, 4), vacuum=40.0, orthogonal=orthogonal, periodic=True)
mol = molecule("CO")
add_adsorbate(slab_mol_fcc, mol, 3.0, "fcc")
slab_mol_hcp = fcc111("Pt", a=a, size=(4, 4, 4), vacuum=40.0, orthogonal=orthogonal, periodic=True)
mol = molecule("CO")
add_adsorbate(slab_mol_hcp, mol, 3.0, "hcp")
[16]:
view_ngl([slab_mol_ontop, slab_mol_bridge, slab_mol_fcc, slab_mol_hcp])
[17]:
from ase.io import write
from IPython.display import Image
import matplotlib.pyplot as plt
import matplotlib.image as mpimg
write("output/pt_ontop.png", slab_mol_ontop, rotation="0x,0y,90z")
write("output/pt_bridge.png", slab_mol_bridge, rotation="0x,0y,90z")
write("output/pt_fcc.png", slab_mol_fcc, rotation="0x,0y,90z")
write("output/pt_hcp.png", slab_mol_hcp, rotation="0x,0y,90z")
fig, axes = plt.subplots(1, 4, figsize=(12, 3))
for i, adsorbate_name in enumerate(["ontop", "bridge", "fcc", "hcp"]):
ax = axes[i]
ax.imshow(mpimg.imread(f"output/pt_{adsorbate_name}.png"))
ax.set_axis_off()
ax.set_title(adsorbate_name)
fig.show()
上記のように、
ontop: ある１つの原子の上bridge: ある２つの原子の間fcc: ３つの原子の間で、2層目には原子が存在しない位置hcp: ３つの原子の間で、2層目の原子の直上の位置
という形となっています。
分子が表面構造のどの部分に吸着するのか(吸着サイト)は、表面・分子系によって異なり自明ではありません。 吸着エネルギーが最小(吸着構造のエネルギーが最小)となる場所を探す必要があります。
[18]:
E_slab_mol_ontop = get_opt_energy(slab_mol_ontop, fmax=1e-3)
E_slab_mol_bridge = get_opt_energy(slab_mol_bridge, fmax=1e-3)
E_slab_mol_fcc = get_opt_energy(slab_mol_fcc, fmax=1e-3)
E_slab_mol_hcp = get_opt_energy(slab_mol_hcp, fmax=1e-3)
Step Time Energy fmax
LBFGS: 0 07:55:12 -342.399913 0.902056
LBFGS: 1 07:55:12 -342.463811 1.697190
LBFGS: 2 07:55:12 -342.595038 0.461123
/tmp/ipykernel_8561/1510506507.py:7: DeprecationWarning: Please use atoms.calc = calc
atoms.set_calculator(calculator)
LBFGS: 3 07:55:13 -342.609555 0.374318
LBFGS: 4 07:55:13 -342.655791 0.619672
LBFGS: 5 07:55:13 -342.671952 0.635448
LBFGS: 6 07:55:13 -342.685588 0.410953
LBFGS: 7 07:55:13 -342.695934 0.185126
LBFGS: 8 07:55:13 -342.707676 0.303394
LBFGS: 9 07:55:13 -342.720449 0.393400
LBFGS: 10 07:55:13 -342.732849 0.282384
LBFGS: 11 07:55:13 -342.739609 0.113871
LBFGS: 12 07:55:13 -342.744498 0.145904
LBFGS: 13 07:55:13 -342.748734 0.254533
LBFGS: 14 07:55:13 -342.751916 0.220931
LBFGS: 15 07:55:13 -342.754138 0.085339
LBFGS: 16 07:55:14 -342.755417 0.053165
LBFGS: 17 07:55:14 -342.756478 0.118874
LBFGS: 18 07:55:14 -342.757792 0.145426
LBFGS: 19 07:55:14 -342.758954 0.095402
LBFGS: 20 07:55:14 -342.759746 0.036355
LBFGS: 21 07:55:14 -342.760329 0.043886
LBFGS: 22 07:55:14 -342.760766 0.069188
LBFGS: 23 07:55:14 -342.761198 0.053732
LBFGS: 24 07:55:14 -342.761433 0.028790
LBFGS: 25 07:55:14 -342.761631 0.022934
LBFGS: 26 07:55:14 -342.761783 0.031655
LBFGS: 27 07:55:14 -342.761892 0.023374
LBFGS: 28 07:55:14 -342.761949 0.010753
LBFGS: 29 07:55:15 -342.761996 0.011726
LBFGS: 30 07:55:15 -342.762020 0.017852
LBFGS: 31 07:55:15 -342.762072 0.012363
LBFGS: 32 07:55:15 -342.762076 0.004995
LBFGS: 33 07:55:15 -342.762060 0.004476
LBFGS: 34 07:55:15 -342.762074 0.006189
LBFGS: 35 07:55:15 -342.762062 0.003668
LBFGS: 36 07:55:15 -342.762069 0.001881
LBFGS: 37 07:55:15 -342.762082 0.002786
LBFGS: 38 07:55:15 -342.762088 0.003632
LBFGS: 39 07:55:15 -342.762079 0.002243
LBFGS: 40 07:55:15 -342.762090 0.001409
LBFGS: 41 07:55:16 -342.762078 0.001220
LBFGS: 42 07:55:16 -342.762077 0.001632
LBFGS: 43 07:55:16 -342.762085 0.001006
LBFGS: 44 07:55:16 -342.762064 0.000690
Step Time Energy fmax
LBFGS: 0 07:55:16 -341.994047 3.707106
LBFGS: 1 07:55:16 -342.131367 4.426340
LBFGS: 2 07:55:16 -342.429021 3.370986
LBFGS: 3 07:55:16 -342.724371 0.943083
LBFGS: 4 07:55:16 -342.800045 0.793617
LBFGS: 5 07:55:16 -342.828909 0.437020
LBFGS: 6 07:55:21 -342.872912 0.237315
LBFGS: 7 07:55:21 -342.876481 0.111741
LBFGS: 8 07:55:22 -342.880221 0.075046
LBFGS: 9 07:55:22 -342.883296 0.095753
LBFGS: 10 07:55:22 -342.885190 0.100857
LBFGS: 11 07:55:22 -342.886260 0.065888
LBFGS: 12 07:55:22 -342.887198 0.042198
LBFGS: 13 07:55:22 -342.888047 0.035742
LBFGS: 14 07:55:22 -342.888611 0.043609
LBFGS: 15 07:55:22 -342.888874 0.031778
LBFGS: 16 07:55:22 -342.889163 0.022925
LBFGS: 17 07:55:22 -342.889449 0.022104
LBFGS: 18 07:55:22 -342.889715 0.026693
LBFGS: 19 07:55:22 -342.889827 0.021689
LBFGS: 20 07:55:22 -342.889971 0.017798
LBFGS: 21 07:55:23 -342.890115 0.016524
LBFGS: 22 07:55:23 -342.890204 0.016821
LBFGS: 23 07:55:23 -342.890311 0.011784
LBFGS: 24 07:55:23 -342.890312 0.010002
LBFGS: 25 07:55:23 -342.890379 0.010338
LBFGS: 26 07:55:23 -342.890426 0.012522
LBFGS: 27 07:55:23 -342.890446 0.009515
LBFGS: 28 07:55:23 -342.890470 0.007352
LBFGS: 29 07:55:23 -342.890488 0.007702
LBFGS: 30 07:55:23 -342.890495 0.008391
LBFGS: 31 07:55:23 -342.890538 0.007104
LBFGS: 32 07:55:23 -342.890515 0.006136
LBFGS: 33 07:55:23 -342.890560 0.005694
LBFGS: 34 07:55:24 -342.890564 0.005681
LBFGS: 35 07:55:24 -342.890595 0.005638
LBFGS: 36 07:55:24 -342.890613 0.005325
LBFGS: 37 07:55:24 -342.890593 0.005027
LBFGS: 38 07:55:24 -342.890638 0.004830
LBFGS: 39 07:55:24 -342.890628 0.004811
LBFGS: 40 07:55:24 -342.890655 0.004379
LBFGS: 41 07:55:24 -342.890638 0.004417
LBFGS: 42 07:55:24 -342.890677 0.004661
LBFGS: 43 07:55:24 -342.890642 0.005589
LBFGS: 44 07:55:24 -342.890674 0.004677
LBFGS: 45 07:55:24 -342.890695 0.003663
LBFGS: 46 07:55:24 -342.890697 0.003655
LBFGS: 47 07:55:24 -342.890697 0.005094
LBFGS: 48 07:55:25 -342.890694 0.004730
LBFGS: 49 07:55:25 -342.890706 0.003387
LBFGS: 50 07:55:25 -342.890704 0.003405
LBFGS: 51 07:55:25 -342.890722 0.003285
LBFGS: 52 07:55:25 -342.890737 0.003030
LBFGS: 53 07:55:25 -342.890744 0.002017
LBFGS: 54 07:55:25 -342.890721 0.001996
LBFGS: 55 07:55:25 -342.890721 0.001847
LBFGS: 56 07:55:25 -342.890743 0.001673
LBFGS: 57 07:55:25 -342.890738 0.000996
Step Time Energy fmax
LBFGS: 0 07:55:25 -341.816867 3.637716
LBFGS: 1 07:55:25 -341.917219 4.362194
LBFGS: 2 07:55:26 -342.143875 3.080004
LBFGS: 3 07:55:26 -342.591196 2.065531
LBFGS: 4 07:55:26 -342.676578 1.254161
LBFGS: 5 07:55:26 -342.853274 0.396045
LBFGS: 6 07:55:26 -342.861501 0.281914
LBFGS: 7 07:55:26 -342.875726 0.297499
LBFGS: 8 07:55:26 -342.883727 0.291233
LBFGS: 9 07:55:26 -342.887155 0.137761
LBFGS: 10 07:55:26 -342.890183 0.139746
LBFGS: 11 07:55:26 -342.894227 0.239277
LBFGS: 12 07:55:26 -342.897257 0.203881
LBFGS: 13 07:55:26 -342.898728 0.084899
LBFGS: 14 07:55:26 -342.899574 0.060572
LBFGS: 15 07:55:27 -342.900480 0.124482
LBFGS: 16 07:55:27 -342.901278 0.125108
LBFGS: 17 07:55:27 -342.901797 0.063454
LBFGS: 18 07:55:27 -342.902055 0.022635
LBFGS: 19 07:55:27 -342.902177 0.051499
LBFGS: 20 07:55:27 -342.902406 0.063136
LBFGS: 21 07:55:27 -342.902540 0.040770
LBFGS: 22 07:55:27 -342.902595 0.010407
LBFGS: 23 07:55:27 -342.902655 0.014829
LBFGS: 24 07:55:27 -342.902664 0.021484
LBFGS: 25 07:55:27 -342.902694 0.015904
LBFGS: 26 07:55:27 -342.902721 0.005062
LBFGS: 27 07:55:27 -342.902700 0.008799
LBFGS: 28 07:55:27 -342.902728 0.013920
LBFGS: 29 07:55:28 -342.902768 0.011519
LBFGS: 30 07:55:28 -342.902771 0.004520
LBFGS: 31 07:55:28 -342.902798 0.005936
LBFGS: 32 07:55:28 -342.902790 0.011713
LBFGS: 33 07:55:28 -342.902779 0.011379
LBFGS: 34 07:55:28 -342.902786 0.005375
LBFGS: 35 07:55:28 -342.902779 0.002169
LBFGS: 36 07:55:28 -342.902765 0.004994
LBFGS: 37 07:55:28 -342.902789 0.005137
LBFGS: 38 07:55:28 -342.902790 0.002454
LBFGS: 39 07:55:28 -342.902776 0.000694
Step Time Energy fmax
LBFGS: 0 07:55:28 -341.805736 3.630377
LBFGS: 1 07:55:28 -341.906011 4.360864
LBFGS: 2 07:55:28 -342.132995 3.070995
LBFGS: 3 07:55:29 -342.575795 2.031545
LBFGS: 4 07:55:29 -342.659040 1.238490
LBFGS: 5 07:55:29 -342.830706 0.387560
LBFGS: 6 07:55:29 -342.838962 0.292390
LBFGS: 7 07:55:29 -342.850377 0.238325
LBFGS: 8 07:55:29 -342.861258 0.283475
LBFGS: 9 07:55:29 -342.864656 0.159977
LBFGS: 10 07:55:29 -342.868207 0.123694
LBFGS: 11 07:55:29 -342.872319 0.219206
LBFGS: 12 07:55:29 -342.875880 0.196974
LBFGS: 13 07:55:29 -342.877610 0.088268
LBFGS: 14 07:55:29 -342.878570 0.084757
LBFGS: 15 07:55:30 -342.879534 0.111292
LBFGS: 16 07:55:30 -342.880512 0.115116
LBFGS: 17 07:55:30 -342.881138 0.059364
LBFGS: 18 07:55:30 -342.881456 0.022330
LBFGS: 19 07:55:30 -342.881648 0.046044
LBFGS: 20 07:55:30 -342.881885 0.053383
LBFGS: 21 07:55:30 -342.882048 0.032017
LBFGS: 22 07:55:30 -342.882155 0.014687
LBFGS: 23 07:55:30 -342.882240 0.021256
LBFGS: 24 07:55:30 -342.882318 0.023533
LBFGS: 25 07:55:30 -342.882350 0.013021
LBFGS: 26 07:55:30 -342.882363 0.008467
LBFGS: 27 07:55:30 -342.882427 0.016556
LBFGS: 28 07:55:30 -342.882411 0.018341
LBFGS: 29 07:55:31 -342.882437 0.010138
LBFGS: 30 07:55:31 -342.882491 0.007715
LBFGS: 31 07:55:31 -342.882510 0.011601
LBFGS: 32 07:55:31 -342.882539 0.014853
LBFGS: 33 07:55:31 -342.882517 0.009812
LBFGS: 34 07:55:31 -342.882496 0.004526
LBFGS: 35 07:55:31 -342.882508 0.004830
LBFGS: 36 07:55:31 -342.882536 0.006192
LBFGS: 37 07:55:31 -342.882526 0.003496
LBFGS: 38 07:55:31 -342.882528 0.001782
LBFGS: 39 07:55:36 -342.882551 0.002238
LBFGS: 40 07:55:37 -342.882544 0.003163
LBFGS: 41 07:55:37 -342.882551 0.002128
LBFGS: 42 07:55:37 -342.882529 0.000808
[19]:
import pandas as pd
series = pd.Series({
"ontop": E_slab_mol_ontop,
"bridge": E_slab_mol_bridge,
"fcc": E_slab_mol_fcc,
"hcp": E_slab_mol_hcp,
})
E_adsorp_series = series - (E_slab + E_mol)
E_adsorp_series.plot(style="o-", xlabel="Adsorption site", ylabel="Energy [eV]", title="Adsorption energy for each adsorption site")
plt.show()
このように、吸着サイトによって吸着構造のエネルギーが異なることがわかります。 今回の計算では、fcc hollowサイトへの吸着構造が一番安定で、吸着エネルギーが一番小さいという結果が得られました。
この結果を以下の文献のTable 1.やTable S2. に記載のPBE 交換相関汎関数を用いて行われた計算結果と比較してみましょう。吸着エネルギーがTop、Bridge、HCP Hollowサイトでそれぞれ1.80, 1.88, 1.89 eV と報告されています。 (ちなみに、実験ではTopサイトが最安定となることが知られており、PBEを用いて行われた計算では実験結果と一致する結果が出ないと報告されています。)
[20]:
view_ngl([slab_mol_ontop, slab_mol_bridge, slab_mol_fcc, slab_mol_hcp])
FixAtoms constraints を使用する例¶
上記の計算例では、すべての原子座標を構造緩和しましたが、DFTを用いて吸着構造を求める場合、構造の安定性を高めるためにSlabの下層の原子座標を固定して構造緩和を行うことがあります。
この方法での計算も行ってみましょう。ASEのconstraintであるFixAtomsを用いることで、指定した index の原子座標を固定することができます。
[21]:
%%time
slab_mol_fixed = fcc111("Pt", a=a, size=(4, 4, 4), vacuum=40.0)
slab_mol_fixed.pbc = True
mol2 = molecule("CO")
add_adsorbate(slab_mol_fixed, mol2, 3.0, "bridge")
view_ngl(slab_mol_fixed)
CPU times: user 104 ms, sys: 3.72 ms, total: 108 ms
Wall time: 141 ms
今回は下から2層の原子を指定する際に、z座標でしきい値を設定します。
上記の可視化結果をマウスホバーしてもわかりますが、今回は以下のようにz軸の分布を可視化して確認してみます。
[22]:
import pandas as pd
atoms = slab_mol_fixed
df = pd.DataFrame({
"x": atoms.positions[:, 0],
"y": atoms.positions[:, 1],
"z": atoms.positions[:, 2],
"symbol": atoms.symbols,
})
df
[22]:
| x | y | z | symbol | |
|---|---|---|---|---|
| 0 | 0.000000 | 0.000000 | 40.000000 | Pt |
| 1 | 2.805396 | 0.000000 | 40.000000 | Pt |
| 2 | 5.610791 | 0.000000 | 40.000000 | Pt |
| 3 | 8.416187 | 0.000000 | 40.000000 | Pt |
| 4 | 1.402698 | 2.429544 | 40.000000 | Pt |
| ... | ... | ... | ... | ... |
| 61 | 7.013489 | 7.288632 | 46.871788 | Pt |
| 62 | 9.818885 | 7.288632 | 46.871788 | Pt |
| 63 | 12.624281 | 7.288632 | 46.871788 | Pt |
| 64 | 1.402698 | 0.000000 | 49.871788 | O |
| 65 | 1.402698 | 0.000000 | 48.721448 | C |
66 rows × 4 columns
[23]:
import plotly.express as px
coord = "z"
df_sorted = df.sort_values(coord).reset_index().rename({"index": "atom_index"}, axis=1)
fig = px.scatter(df_sorted, x=df_sorted.index, y=coord, color="symbol", hover_data=["x", "y", "z", "atom_index"])
fig.show()
2層目と3層目の境界が、42.29 - 44.58 の間であることがわかりました。
FixAtoms は indicesに固定したい原子のindexを指定するか、maskにそれぞれの原子を固定するかしないかの配列を指定します。 今回はmaskを使用してみます。
[24]:
from ase.constraints import FixAtoms
thresh = 43.0
constraint = FixAtoms(mask=slab_mol_fixed.positions[:, 2] < thresh)
constraint
[24]:
FixAtoms(indices=[0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31])
作成したconstraintは、 atoms.set_constraint methodで指定することで制約がかけられます。
[25]:
slab_mol_fixed.set_constraint(constraint)
可視化して、意図した場所を固定できているか確認してみます。
[26]:
v = view_ngl(slab_mol_fixed)
v.show_fix_atoms_constraint()
v
[27]:
E_slab_mol_fixed = get_opt_energy(slab_mol_fixed, fmax=1e-3)
Step Time Energy fmax
LBFGS: 0 07:55:39 -341.994038 3.707107
/tmp/ipykernel_8561/1510506507.py:7: DeprecationWarning:
Please use atoms.calc = calc
LBFGS: 1 07:55:39 -342.107696 4.426356
LBFGS: 2 07:55:39 -342.373290 3.317400
LBFGS: 3 07:55:40 -342.579909 2.000243
LBFGS: 4 07:55:40 -342.756639 0.870773
LBFGS: 5 07:55:40 -342.790093 0.374673
LBFGS: 6 07:55:40 -342.809422 0.291128
LBFGS: 7 07:55:40 -342.823264 0.266580
LBFGS: 8 07:55:40 -342.833950 0.124070
LBFGS: 9 07:55:40 -342.835572 0.088253
LBFGS: 10 07:55:40 -342.836857 0.075296
LBFGS: 11 07:55:40 -342.838347 0.068667
LBFGS: 12 07:55:40 -342.839204 0.064898
LBFGS: 13 07:55:41 -342.839661 0.064189
LBFGS: 14 07:55:41 -342.840044 0.051451
LBFGS: 15 07:55:41 -342.840505 0.050056
LBFGS: 16 07:55:41 -342.840820 0.032795
LBFGS: 17 07:55:41 -342.841054 0.030308
LBFGS: 18 07:55:41 -342.841195 0.027926
LBFGS: 19 07:55:41 -342.841300 0.032563
LBFGS: 20 07:55:41 -342.841448 0.026804
LBFGS: 21 07:55:41 -342.841483 0.013311
LBFGS: 22 07:55:41 -342.841538 0.014673
LBFGS: 23 07:55:41 -342.841598 0.019561
LBFGS: 24 07:55:41 -342.841640 0.017394
LBFGS: 25 07:55:41 -342.841670 0.009116
LBFGS: 26 07:55:41 -342.841695 0.009561
LBFGS: 27 07:55:42 -342.841708 0.011721
LBFGS: 28 07:55:42 -342.841772 0.014140
LBFGS: 29 07:55:42 -342.841791 0.009850
LBFGS: 30 07:55:42 -342.841842 0.006582
LBFGS: 31 07:55:42 -342.841835 0.008020
LBFGS: 32 07:55:42 -342.841845 0.010491
LBFGS: 33 07:55:42 -342.841836 0.009795
LBFGS: 34 07:55:42 -342.841877 0.007181
LBFGS: 35 07:55:42 -342.841859 0.006186
LBFGS: 36 07:55:42 -342.841893 0.009587
LBFGS: 37 07:55:42 -342.841898 0.010452
LBFGS: 38 07:55:42 -342.841930 0.007519
LBFGS: 39 07:55:43 -342.841952 0.006661
LBFGS: 40 07:55:43 -342.841949 0.008566
LBFGS: 41 07:55:43 -342.841968 0.010709
LBFGS: 42 07:55:43 -342.842006 0.008396
LBFGS: 43 07:55:43 -342.842020 0.005409
LBFGS: 44 07:55:43 -342.842028 0.004772
LBFGS: 45 07:55:43 -342.842030 0.006176
LBFGS: 46 07:55:43 -342.842019 0.005571
LBFGS: 47 07:55:43 -342.842023 0.002580
LBFGS: 48 07:55:43 -342.842029 0.002196
LBFGS: 49 07:55:43 -342.842009 0.002595
LBFGS: 50 07:55:43 -342.842016 0.002606
LBFGS: 51 07:55:43 -342.842025 0.001621
LBFGS: 52 07:55:44 -342.842046 0.001115
LBFGS: 53 07:55:44 -342.842026 0.000947
[28]:
E_adsorp = E_slab_mol_fixed - (E_slab + E_mol)
print(f"E_adsorp {E_adsorp}, E_slab_mol {E_slab_mol_fixed}, E_slab {E_slab}, E_mol {E_mol}")
E_adsorp -1.8266751227401983, E_slab_mol -342.84202612473905, E_slab -329.5151524229501, E_mol -11.500198579048762
下層を固定する方法を用いた場合、吸着エネルギーとして-1.83 eVが得られました。
参考文献
CO on Pd¶
同様に、Pd上でのCO吸着に関しても計算してみましょう。
[29]:
bulk_atoms = bulk("Pd", cubic=True)
np.mean(np.diag(bulk_atoms.cell))
[29]:
3.89
[30]:
bulk_atoms.calc = calculator
E_bulk = get_opt_energy(bulk_atoms, fmax=1e-4, opt_mode="scale")
a_pd = np.mean(np.diag(bulk_atoms.cell))
a_pd
Step Time Energy fmax
LBFGS: 0 07:55:44 -14.872189 4.690001
LBFGS: 1 07:55:44 -14.761179 7.138004
/tmp/ipykernel_8561/1510506507.py:7: DeprecationWarning:
Please use atoms.calc = calc
/tmp/ipykernel_8561/1510506507.py:9: FutureWarning:
Import StrainFilter from ase.filters
LBFGS: 2 07:55:44 -14.923395 0.707045
LBFGS: 3 07:55:44 -14.924615 0.152445
LBFGS: 4 07:55:44 -14.924664 0.001121
LBFGS: 5 07:55:44 -14.924666 0.000014
[30]:
3.9414397766391507
[31]:
from ase.build import fcc111, molecule, add_adsorbate
slab = fcc111("Pd", a=a_pd, size=(4, 4, 4), vacuum=40.0, periodic=True)
slab.calc = calculator
E_slab = get_opt_energy(slab, fmax=1e-4, opt_mode="normal")
Step Time Energy fmax
LBFGS: 0 07:55:44 -220.564856 0.115019
LBFGS: 1 07:55:44 -220.571327 0.092847
/tmp/ipykernel_8561/1510506507.py:7: DeprecationWarning:
Please use atoms.calc = calc
LBFGS: 2 07:55:44 -220.582634 0.010980
LBFGS: 3 07:55:44 -220.582743 0.001272
LBFGS: 4 07:55:44 -220.582707 0.000075
[32]:
mol = molecule("CO")
mol.calc = calculator
E_mol = get_opt_energy(mol, fmax=1e-4)
Step Time Energy fmax
LBFGS: 0 07:55:44 -11.497234 0.820158
LBFGS: 1 07:55:45 -11.484527 1.958979
LBFGS: 2 07:55:45 -11.500196 0.029784
/tmp/ipykernel_8561/1510506507.py:7: DeprecationWarning:
Please use atoms.calc = calc
LBFGS: 3 07:55:45 -11.500199 0.001079
LBFGS: 4 07:55:45 -11.500198 0.000010
[33]:
%%time
slab_mol = fcc111("Pd", a=a_pd, size=(4, 4, 4), vacuum=40.0)
slab_mol.pbc = True
mol2 = molecule("CO")
add_adsorbate(slab_mol, mol2, 3.0, "ontop")
E_slab_mol = get_opt_energy(slab_mol, fmax=1e-3)
Step Time Energy fmax
LBFGS: 0 07:55:45 -233.487857 1.079220
LBFGS: 1 07:55:45 -233.491375 1.894393
/tmp/ipykernel_8561/1510506507.py:7: DeprecationWarning:
Please use atoms.calc = calc
LBFGS: 2 07:55:45 -233.515930 0.671719
LBFGS: 3 07:55:45 -233.535165 0.525847
LBFGS: 4 07:55:45 -233.540258 0.339795
LBFGS: 5 07:55:45 -233.543925 0.176446
LBFGS: 6 07:55:45 -233.548123 0.344875
LBFGS: 7 07:55:45 -233.552475 0.420478
LBFGS: 8 07:55:45 -233.555398 0.288876
LBFGS: 9 07:55:46 -233.557775 0.141699
LBFGS: 10 07:55:46 -233.560570 0.271896
LBFGS: 11 07:55:46 -233.564281 0.313258
LBFGS: 12 07:55:46 -233.568309 0.284822
LBFGS: 13 07:55:46 -233.571348 0.197042
LBFGS: 14 07:55:46 -233.573282 0.113801
LBFGS: 15 07:55:46 -233.574606 0.148960
LBFGS: 16 07:55:46 -233.575601 0.152804
LBFGS: 17 07:55:46 -233.576371 0.098615
LBFGS: 18 07:55:46 -233.576925 0.083781
LBFGS: 19 07:55:46 -233.577258 0.069899
LBFGS: 20 07:55:46 -233.577481 0.047896
LBFGS: 21 07:55:46 -233.577654 0.040826
LBFGS: 22 07:55:46 -233.577835 0.067757
LBFGS: 23 07:55:47 -233.578010 0.069439
LBFGS: 24 07:55:47 -233.578084 0.034979
LBFGS: 25 07:55:47 -233.578156 0.019567
LBFGS: 26 07:55:47 -233.578243 0.037782
LBFGS: 27 07:55:47 -233.578312 0.055524
LBFGS: 28 07:55:47 -233.578400 0.041176
LBFGS: 29 07:55:47 -233.578460 0.014200
LBFGS: 30 07:55:47 -233.578464 0.009982
LBFGS: 31 07:55:47 -233.578483 0.020392
LBFGS: 32 07:55:47 -233.578480 0.019457
LBFGS: 33 07:55:47 -233.578483 0.008402
LBFGS: 34 07:55:47 -233.578493 0.002447
LBFGS: 35 07:55:52 -233.578524 0.008048
LBFGS: 36 07:55:53 -233.578501 0.008466
LBFGS: 37 07:55:53 -233.578520 0.004002
LBFGS: 38 07:55:53 -233.578512 0.000913
CPU times: user 284 ms, sys: 30.7 ms, total: 315 ms
Wall time: 7.95 s
[34]:
E_adsorp = E_slab_mol - (E_slab + E_mol)
print(f"E_adsorp {E_adsorp}, E_slab_mol {E_slab_mol}, E_slab {E_slab}, E_mol {E_mol}")
E_adsorp -1.4956079430033924, E_slab_mol -233.57851244191073, E_slab -220.58270687257885, E_mol -11.500197626328491
[35]:
from pfcc_extras.visualize.view import view_ngl
view_ngl(slab_mol)
このように、様々な表面や分子の組み合わせに対して吸着エネルギーが計算できます。 これは5章で説明する、効率の良い触媒探索に際して欠かせないプロセスとなります。
被覆率¶
本Tutorialで取り上げた例では、まっさらな表面構造上に、分子が１つだけ吸着するというものでした。 現実では、たくさんの分子が同時に表面と反応を起こし、同時に複数の分子が吸着しているような場合や、ほとんど吸着した分子で覆われているような場合も存在します。
その際には、分子がすでに複数個吸着している表面構造から、あらたに1つの分子が吸着する場合の吸着エネルギーなどにも興味が出てきます。
上記で上げたPd表面上でのCO分子の吸着エネルギーの被覆率依存性の解析に関しては以下の文献で報告されています。
複数箇所への吸着¶
また、今回扱えていない応用トピックとして、2箇所以上に吸着する例があります。 分子が大きく、複雑になってくると複数箇所が表面上に吸着することで安定となるようなものも存在します。