ase

ASE - Atomic Simulation Environment

ASE は Atoms オブジェクトを中心に構築されており、位置、原子番号、単位格子を持つ原子の集合を表します。エネルギーと力を計算するさまざまな「Calculator」（外部コード）とやり取りするための共通インターフェースを提供します。

使用する場面

• 分子、結晶、表面、ナノ粒子などの複雑な原子構造を構築する場合。
• 構造最適化（エネルギーが最小となる構造の探索）を実行する場合。
• さまざまなアンサンブルで分子動力学（MD）シミュレーションを行う場合。
• ポテンシャルエネルギー面（PES）や遷移状態（NEB 法）を計算する場合。
• 原子構造ファイル形式（CIF、XYZ、POSCAR、PDB など）を相互変換する場合。
• 状態密度（DOS）やバンド構造などの電子特性を計算する場合。
• 複数のソフトウェアパッケージを使用するシミュレーションのワークフローを自動化する場合。

リファレンスドキュメント

公式ドキュメント：https://wiki.fysik.dtu.dk/ase/
Calculator の一覧：https://wiki.fysik.dtu.dk/ase/ase/calculators/calculators.html
検索パターン：ase.Atoms、ase.build、ase.optimize、ase.io.read、ase.calculators

基本原則

Atoms オブジェクト

ASE の中核です。次の情報を格納します。

• positions：オングストローム単位のデカルト座標（Nx3 配列）。
• symbols：化学元素（例：'H'、'Fe'）。
• cell：単位格子ベクトル（3x3 行列、または 3 個／6 個の長さ・角度）。
• pbc：周期境界条件（各軸に対するブール値）。

Calculator

ASE 自体はエネルギーを計算しません。Calculator（例：GPAW、VASP、EMT、Lennard-Jones）を Atoms オブジェクトに割り当てます。atoms.get_potential_energy() を呼び出すと、ASE は Calculator に計算を依頼します。

単位

重要なルール：ASE はエネルギーに eV（電子ボルト）、距離にオングストロームを使用します。時間には ASE の単位を使用します（多くの場合、fs に変換されます）。

クイックリファレンス

インストール

pip install ase

標準的なインポート

from ase import Atoms
from ase.build import molecule, bulk, surface
from ase.io import read, write
from ase.optimize import BFGS
from ase.visualize import view

基本パターン - 分子の作成と最適化

from ase.build import molecule
from ase.optimize import BFGS
from ase.calculators.emt import EMT # Simple effective medium theory calculator

# 1. Build structure
atoms = molecule('H2O')

# 2. Attach a calculator
atoms.calc = EMT()

# 3. Optimize geometry
opt = BFGS(atoms, trajectory='opt.traj')
opt.run(fmax=0.05) # Converge until forces are < 0.05 eV/Ang

print(f"Final Energy: {atoms.get_potential_energy():.3f} eV")

重要なルール

✅ すること

• 周期境界条件を指定する - 結晶には pbc=[True, True, True] を設定し、孤立分子には False を設定します。
• NumPy のベクトル化アクセスを使用する - NumPy を使って atoms.positions または atoms.get_positions() を直接操作します。
• 収束を確認する - オプティマイザーが目標の fmax に到達したことを必ず確認します。
• 単位格子を正しく設定する - 周期系では、自己相互作用を防げるだけの大きさの単位格子が必要です。
• ase.build を使用する - スラブ（surface）やナノチューブの作成など、複雑なタスクには組み込み関数を活用します。
• 軌跡を監視する - .traj ファイルを使って、ASE-GUI で最適化や MD の進行状況を確認します。

❌ しないこと

• 原子を手動でループ処理する - 原子を移動する際に Python のループを使うのは避け、atoms.positions += displacement を使用します。
• ASE の単位を忘れる - from ase.units import Hartree, kcal, mol を使って明示的に変換せずに、kcal/mol や Hartree を混在させないでください。
• 真空を無限と見なす - スラブの非周期方向には、十分な真空層（例：10-15 Å）を確保します。
• 原子番号をハードコードする - 可読性を高めるため、79 の代わりに化学記号（'Au'）を使用します。

アンチパターン（絶対に避ける）

# ❌ BAD: Moving atoms in a Python loop
for atom in atoms:
    atom.position += [0.1, 0, 0]

# ✅ GOOD: Vectorized movement
atoms.positions += [0.1, 0, 0]

# ❌ BAD: Not defining a cell for a "bulk" system
iron = Atoms('Fe', positions=[[0, 0, 0]]) # No cell = not a crystal!

# ✅ GOOD: Use the bulk constructor
from ase.build import bulk
iron = bulk('Fe', 'bcc', a=2.87)

# ❌ BAD: Ignoring forces after optimization
energy = atoms.get_potential_energy() # Optimization might have failed!

# ✅ GOOD: Check forces
forces = atoms.get_forces()
max_force = (forces**2).sum(axis=1).max()**0.5
if max_force > 0.05:
    print("Warning: Structure not converged!")

構造の構築（ase.build）

分子、結晶、表面

from ase.build import molecule, bulk, surface, add_adsorbate

# Molecule from database
h2o = molecule('H2O')

# Bulk crystal (Copper fcc)
cu = bulk('Cu', 'fcc', a=3.6)

# Surface slab (Al 111 surface, 3 layers)
slab = surface('Al', (1, 1, 1), layers=3)
slab.center(vacuum=10, axis=2) # Add 10A vacuum on Z axis

# Adding an adsorbate (CO on a surface)
co = molecule('CO')
add_adsorbate(slab, co, height=2.0, position='ontop')

最適化とダイナミクス（ase.optimize、ase.md）

幾何構造の最小化

from ase.optimize import QuasiNewton

# BFGS, LBFGS, GPMin, QuasiNewton are common choices
opt = QuasiNewton(atoms, trajectory='relax.traj', logfile='relax.log')
opt.run(fmax=0.01)

分子動力学

from ase.md.langevin import Langevin
from ase import units

# MD at 300K with Langevin thermostat
dyn = Langevin(atoms, 
               timestep=1.0 * units.fs, 
               temperature_K=300, 
               friction=0.01)

# Run for 1000 steps
dyn.run(1000)

ファイル入出力（ase.io）

読み込みと書き込み

from ase.io import read, write

# Read from XYZ or CIF
atoms = read('structure.cif')

# Read multiple frames from a trajectory
frames = read('simulation.traj', index=':') # All frames
last_5 = read('simulation.traj', index='-5:') # Last 5 frames

# Write to different formats
write('output.poscar', atoms, format='vasp')
write('movie.gif', frames, rotation='10x,10y,10z') # Create animated gif

応用：遷移状態（NEB）

ナッジド・エラスティック・バンド

from ase.mep import NEB
from ase.optimize import BFGS

# Initial and Final states (must have same atoms/order)
initial = read('initial.traj')
final = read('final.traj')

# Create 5 images between initial and final
images = [initial]
for i in range(5):
    images.append(initial.copy())
images.append(final)

# Interpolate positions
neb = NEB(images)
neb.interpolate()

# Run optimization on all images
opt = BFGS(neb, trajectory='neb.traj')
opt.run(fmax=0.05)

実践的なワークフロー

1. 格子定数のスキャン（状態方程式）

import numpy as np
from ase.build import bulk
from ase.calculators.emt import EMT

def find_opt_lattice():
    volumes = []
    energies = []
    for a in np.linspace(3.5, 3.7, 5):
        atoms = bulk('Cu', 'fcc', a=a)
        atoms.calc = EMT()
        volumes.append(atoms.get_volume())
        energies.append(atoms.get_potential_energy())
    
    # You can then fit these to Birch-Murnaghan EOS
    return volumes, energies

2. 振動解析（熱力学）

from ase.vibrations import Vibrations

# Calculate vibrations (requires a calculator)
vib = Vibrations(atoms)
vib.run()
vib.summary()

# Get Helmholtz free energy at 300K
from ase.thermochemistry import IdealGasThermo
thermo = IdealGasThermo(vib_energies=vib.get_frequencies(), 
                        geometry='nonlinear', 
                        potentialenergy=atoms.get_potential_energy())
F = thermo.get_helmholtz_free_energy(temperature=300)

3. PySCF との連携（量子化学）

from ase.calculators.pyscf_calc import PySCF
from ase.build import molecule

atoms = molecule('H2')
# Define PySCF calculator
atoms.calc = PySCF(mol_options={'basis': '6-31g', 'spin': 0},
                   method='RKS.xc = "b3lyp"')

energy = atoms.get_potential_energy()

パフォーマンスの最適化

対称性の利用

大きな結晶を構築する場合、ASE 内の spglib 連携を使うと対称性の特定に役立ちます。ただし、特定の拘束条件を追加しない限り、ASE のコア構造は最適化中に対称性を自動的に維持しません。

制約

原子を固定すると（例：スラブの下層）、緩和を大幅に高速化できます。

from ase.constraints import FixAtoms

# Fix all atoms with Z-coordinate < 5.0
c = FixAtoms(mask=[atom.position[2] < 5.0 for atom in atoms])
atoms.set_constraint(c)

よくある落とし穴と解決策

「計算機が見つからない」エラー

# ❌ Problem: atoms.get_potential_energy() fails
# ✅ Solution: Ensure .calc is set and external code is in PATH
import os
os.environ['VASP_COMMAND'] = 'mpirun vasp_std' # Example for VASP

セルベクトルの符号規約

ASE では通常、右手系の座標系が使われます。異なる規約を使っている可能性がある古いコードからインポートする際は注意してください。

Atoms オブジェクトの純粋性

リストから原子を削除すると、インデックスが変わります。

# ❌ Problem: Deleting atoms in a loop by index
# ✅ Solution: Use a mask or the .pop() method carefully
del atoms[[0, 5, 10]] # Bulk delete by index list

ASE は、原子シミュレーションの世界における標準的な「接着剤」です。コードを 1 行（atoms.calc）変更するだけで、単純な経験的ポテンシャルから高コストな DFT 計算に切り替えられるため、ハイスループット計算材料科学に欠かせません。