pymatgen结构、表面与界面建模

下载输入、实测结果和完整脚本。在解压目录运行 bash run.sh，结果写入 work/；reference/保留原始结构与CSV，脚本不会自动提交新的电子结构计算。

1 加载环境，建立Si#

module load pymatgen/2026.9.24
python - <<'BUILD'
from pymatgen.core import Structure, Lattice
s = Structure.from_spacegroup('Fd-3m', Lattice.cubic(5.43),
                              ['Si'], [[0, 0, 0]])
s.to(filename='POSCAR-Si')
s.to(filename='Si.cif')
s.to(filename='Si.json')
print(len(s), s.volume)
BUILD

结果为8原子、160.103007 ų。资料包逐一读取CIF/POSCAR/JSON，并用 StructureMatcher(primitive_cell=False) 核对同一结构；比较坐标时允许晶格等价的原子重排，不要求文本行顺序相同。

pmg structure -f POSCAR-Si -s 0.001

标准Si空间群227。对其中一个原子施加(0.001,0.002,0.0015) Å位移，收紧位置容差得到低对称群，放宽至0.01 Å时恢复227。容差越大不代表识别结果越准确，要与结构实际精度对应。

Cu Kα条件下，本例第一个衍射峰约28.47°；这是结构计算的理想粉末峰位，不含实验仪器展宽、织构和样品缺陷。

2 Cu(111)切面与CO吸附初始模型#

将下面内容保存为 surface.py：

from pymatgen.core import Structure, Lattice, Molecule
from pymatgen.core.surface import SlabGenerator
from pymatgen.analysis.adsorption import AdsorbateSiteFinder
cu = Structure.from_spacegroup('Fm-3m', Lattice.cubic(3.61),
                               ['Cu'], [[0, 0, 0]])
slab = SlabGenerator(cu, (1,1,1), 8, 15,
                     center_slab=True).get_slabs()[0]
finder = AdsorbateSiteFinder(slab)
sites = finder.find_adsorption_sites()
co = Molecule(['C','O'], [[0,0,0], [0,0,1.15]])
model = finder.add_adsorbate(co, sites['ontop'][0])
model.to(filename='POSCAR-Cu111-CO')
print(len(slab), len(model), {k:len(v) for k,v in sites.items()})

python surface.py

本例切面4个Cu原子，加CO后6个原子；有1类顶位、1类桥位、2类空位吸附点，总计4类。这里的空位指表面几何hollow位点，不是删除晶格原子的缺陷空位。

下图为便于观察，将输出晶胞在面内重复3×3显示；原始计算输入仍保留最小晶胞。

最小切片厚度8 Å和最小真空15 Å会按整层取整；本次实际真空沿表面法向约18.758 Å。斜晶胞应按表面法向与体积/面内面积求真空，不能直接用c轴长度减z坐标跨度。

3 Cu/Ni共格界面#

from pymatgen.analysis.interfaces.coherent_interfaces import CoherentInterfaceBuilder
from pymatgen.analysis.interfaces.zsl import ZSLGenerator
ni = Structure.from_spacegroup('Fm-3m', Lattice.cubic(3.52),
                               ['Ni'], [[0,0,0]])
b = CoherentInterfaceBuilder(
    substrate_structure=cu, film_structure=ni,
    substrate_miller=(1,1,1), film_miller=(1,1,1),
    zslgen=ZSLGenerator(max_area=60,
                        max_length_tol=0.04, max_angle_tol=0.02))
print(b.terminations)
x = next(b.get_interfaces(b.terminations[0], gap=2,
    vacuum_over_film=15, film_thickness=3,
    substrate_thickness=3, in_layers=True))
x.to(filename='POSCAR-Cu-Ni-interface')

在前述Cu定义后接上这段代码运行。本例得到Cu3Ni3，选择的终止面和匹配变换记录在RESULTS中。两者原始晶格常数约差2.56%，共格化会引入应变；要结合实际晶面、面积、层数和应变阈值选择匹配。

这些图对应初始几何，不能直接给出最稳定吸附位或界面结合能。继续做DFT前，需要松弛、增加层数和面内尺寸，并用一致参照计算能量。

4 生成缺陷结构#

from pymatgen.core import Structure
s = Structure.from_file('POSCAR-Si')
s.make_supercell([2,2,2])
s.remove_sites([0])
s.to(filename='POSCAR-vacancy')

得到Si63。另将8原子常规胞的第一个Si替为Ge，得到Si7Ge1。缺陷结构构造没有自动指定合理的电荷态、化学势或稀释极限；形成能计算还需原始完整超胞和对应修正。