Atomic Simulation Environmentの使い方: エネルギー計算と構造最適化

Last updated at 2025-03-10Posted at 2024-01-29

Atomic Simulation Environment (ASE)とは

第一原理計算/量子化学計算を行うPythonのライブラリ
- https://wiki.fysik.dtu.dk/ase/
第一原理計算の主要部分はcalculatorとして条件を設定し、実際の計算のインプットを編集することなく計算を実行できる
さまざまな計算パッケージをcalculatorに用いることができる

ASEを使う利点

第一原理計算や分子動力学計算を行う計算プログラムは非常に数が多く、それぞれが独自のインプット・アプトプットを持つ
これらのプログラムを用いて計算を行うためには、インプットの作成方法やアウトプットの解析方法に習熟していることが望ましく、学習コストが高い
ASEを用いることで、calculatorの設定を通して計算プログラムのインプット作成処理ができる --> 計算プログラムのインプットの詳細を知っていなくても計算ができる(もちろん知っていた方が利点は多い)
Pythonをベースとして用いるため、for文によるループ処理やif文による分岐が簡単に導入できる

環境設定

ASEのインストール: pip install ase
ASEに組み込まれているポテンシャル(EMT)を用いる場合はこれ以外の設定は不要

VASPを用いる場合

第一原理計算としてメジャーなソフトウェアであるVASPを用いる場合は設定が必要
VASPがインストールされていることが前提(実行ファイルは~/vasp.X.X.X/bin/vasp_stdとする)
ASEから参照するVASP実行用のPythonスクリプト(run_vasp.pyとする)を用意する

import os
exitcode = os.system('~/vasp.X.X.X/bin/vasp_std')
# exitcode = os.system('mpirun -np 2 ~/vasp.X.X.X/bin/vasp_std') # MPI並列の場合はこちらを使う

~/.bashrcに上記のVASP_SCRIPTの場所と擬ポテンシャルの場所(VASP_PP_PATH)を設定する

$ export VASP_SCRIPT=$HOME/somewhere/run_vasp.py
$ export VASP_PP_PATH=$HOME/vasp/potentials/

run_vasp.pyや擬ポテンシャルの場所は適宜設定すること

エネルギー計算

ASEを用いた計算の最も基本的な使い方としてエネルギー計算を行う
Pt結晶のself-consistent field (SCF)計算を行う
ここでは以下のcalculatorを用いる。それぞれ異なる計算手法となっており、計算コスト・計算精度が違ってくる

EMT (多体ポテンシャル)
DFT (VASP)
M3GNet (ニューラルネットワーク・ポテンシャル)
CHGNet (ニューラルネットワーク・ポテンシャル)

1. EMT

Effective medium theory (EMT): 多体ポテンシャルの一種
対応している元素: Al, Cu, Ag, Au, Ni, Pd, Pt, H, C, N, O

from ase.build import bulk
from ase.calculators.emt import EMT

solid = bulk(name="Pt", crystalstructure="fcc", a=3.9)

solid.calc = EMT()

energy = solid.get_potential_energy()
print(f"energy = {energy:5.3} eV")

2. DFT

Vienna ab initio simulation package (VASP)によるDFT計算からエネルギーを求める

from ase.build import bulk
from ase.calculators.vasp import Vasp

solid = bulk(name="Pt", crystalstructure="fcc", a=3.9)

solid.calc = Vasp(xc="pbe", pp="potpaw_PBE.64")

energy = solid.get_potential_energy()
print(f"energy = {energy:5.3} eV")

Vaspのxcは汎函数, ppは擬ポテンシャルの種類を指定する
RuntimeError: Looking for PP for potpaw_PBE/Pt/POTCARが出力される場合はVASP_PP_PATHもしくはppの指定がうまくいっていない
上記だと$VASP_PP_PATH/potpaw_PBE.64以下が参照される

3. M3GNet

ニューラルネットワーク・ポテンシャルの1種
DFTでフィットされたポテンシャルで、DFTと同等の精度であるが計算が早い
M3GNetのインストール: pip install matgl
- matglはdgl(deep graph library)を使うが、PyTorchのバージョンが新しいと対応していないことがある。torch==2.3.0, dgl==2.2.1なら動いた(2025/2/22)

from ase.build import bulk
import matgl
from matgl.ext.ase import PESCalculator

solid = bulk(name="Pt", crystalstructure="fcc", a=3.9)

potential = matgl.load_model("M3GNet-MP-2021.2.8-PES")
solid.calc = PESCalculator(potential=potential)

energy = solid.get_potential_energy()
print(f"energy = {energy:5.3} eV")

4. CHGNet

M3GNet同様でニューラルネットワーク・ポテンシャルの1種
CHGNetのインストール: pip install chgnet
M3GNetより精度が高い気がする

from ase.build import bulk
from chgnet.model.dynamics import CHGNetCalculator
from chgnet.model.model import CHGNet
    
solid = bulk(name="Pt", crystalstructure="fcc", a=3.9)

chgnet = CHGNet.load()
potential = CHGNetCalculator(potential=chgnet, properties="energy")
solid.calc = potential
        
energy = solid.get_potential_energy()
print(f"energy = {energy:5.3} eV")

構造最適化

Ptのfcc(111)面を構築して構造最適化計算を行う
原子層は4層とし、上2層のみを最適化(下2層は固定)

1. EMT, M3GNet, CHGNet

from ase.calculators.emt import EMT
from ase.build import fcc111
from ase.optimize import FIRE
from ase.constraints import FixAtoms

# M3Gnet
import matgl
from matgl.ext.ase import PESCalculator

# CHGNet
from chgnet.model.dynamics import CHGNetCalculator
from chgnet.model.model import CHGNet

surf = fcc111(symbol="Pt", size=[1, 1, 4], a=3.9, vacuum=12.0)
surf.pbc = True
c = FixAtoms(indices=[atom.index for atom in surf if atom.tag >= 3])
surf.constraints = c

# --- EMT
calc = EMT()

# --- M3GNet
# potential = matgl.load_model("M3GNet-MP-2021.2.8-PES")
# calc = PESCalculator(potential=potential)

# --- CHGNet
# potential = CHGNet.load()
# calc = CHGNetCalculator(potential=potential, properties="energy")

surf.calc = calc
        
opt = FIRE(surf, trajectory="pt-relax.traj")
opt.run(fmax=0.05, steps=100)

fcc111はlayer numberに応じてAtomのtagを自動的に決める
- 第一層: tag=1, 第二層: tag=2, 吸着子: tag=0, などとなる。このtagを用いて原子層を選択する。tagを確認するにはprint(atoms.get_tags())
ASEでよく使われるoptimizerはBFGSだが、FIREやFIRE2のほうが素性がよいかもしれない
.trajファイルの可視化: ase gui pt-relax.traj

2. VASP

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

surf = fcc111(symbol="Pt", size=[1, 1, 4], a=3.9, vacuum=12.0)
surf.pbc = True
c = FixAtoms(indices=[atom.index for atom in surf if atom.tag >= 3])
surf.constraints = c

surf.calc = Vasp(xc="pbe", pp="potpaw_PBE.64", ibrion=2, nsw=10,potim=0.1, directory="pt_relax")
surf.get_potential_energy()

EMTのようにASEのoptimizerを使うのでも良いが、VASPの最適化を用いた方が収束が早い(と思う)
Vaspのibrion, nsw, potimはVASPのINCARで指定するものと同じ
directoryは必須ではないが、設定した場合そのディレクトリで(ない場合は作成される)VASPの計算が行われるので設定しておくとよい。以下のどちらかで設定できる
- Atoms.calc.directory = "some_directory"
- Vasp(direcoty="some_directory")

おわりに

ここではASEの基本的な使い方としてエネルギー計算と構造最適化の例を示した
CalculatorにはEMTとVASPを用いた。EMTはASEに含まれておりそのまま利用でき計算も早いが、計算精度が相当悪いのと、利用できる元素の種類が非常に少ない(H, C, N, O, Al, Cu, Ag, Au, Ni, Pd, Ptのみ)が欠点であり、コードの確認だけに用いた方がよい
NNPを使うとより多くの元素を扱うことができ
他の機能は別ページで説明する予定

You get articles that match your needs
You can efficiently read back useful information
You can use dark theme

What you can do with signing up