ase

ASE - 原子模拟环境

ASE 围绕 Atoms 对象构建，该对象表示由原子组成的集合，其中包含原子的位置、原子序数和晶胞。它提供了一个通用接口，用于与各种“计算器”（计算能量和力的外部代码）交互。

适用场景

• 构建复杂的原子结构：分子、晶体、表面和纳米颗粒。
• 运行几何优化（寻找能量最低的结构）。
• 在各种系综中执行分子动力学（MD）模拟。
• 计算势能面（PES）和过渡态（NEB 方法）。
• 在原子文件格式之间转换（CIF、XYZ、POSCAR、PDB 等）。
• 计算电子性质，例如态密度（DOS）和能带结构。
• 自动化涉及多个软件包的模拟工作流。

参考文档

官方文档：https://wiki.fysik.dtu.dk/ase/
计算器列表：https://wiki.fysik.dtu.dk/ase/ase/calculators/calculators.html
搜索模式：ase.Atoms、ase.build、ase.optimize、ase.io.read、ase.calculators

核心原则

Atoms 对象

ASE 的核心。它存储：

• positions：以埃为单位的笛卡尔坐标（N×3 数组）。
• symbols：化学元素（例如，'H'、'Fe'）。
• cell：晶胞矢量（3×3 矩阵，或 3/6 个长度/角度）。
• pbc：周期性边界条件（每个轴对应一个布尔值）。

计算器

ASE 本身不计算能量。你需要将计算器（例如，GPAW、VASP、EMT、Lennard-Jones）附加到 Atoms 对象上。当你调用 atoms.get_potential_energy() 时，ASE 会请求计算器执行计算。

单位

重要规则：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 单位 - 不要混用 kcal/mol 或 Hartree，除非明确进行转换（from ase.units import Hartree, kcal, mol）。
• 假设真空无限大 - 对于薄片中的非周期性方向，确保留有足够的真空填充（例如，10-15 Å）。
• 硬编码原子序数 - 为了提高可读性，应使用化学符号（'Au'），而不是 79。

反模式（绝对不要）

# ❌ 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)

文件 I/O（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 的核心结构在优化期间不会自动强制满足对称性。

约束

固定原子（例如 slab 的底层）可以显著加快弛豫速度。

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 是原子模拟领域的标准“粘合剂”。只需更改一行代码（atoms.calc），就能让你从简单的经验势切换到昂贵的 DFT 计算，因此它对于高通量计算材料科学不可或缺。