CREST构象集合工具：RMSD、拆分、比较与原子顺序

搜索结束得到多帧XYZ后，通常还要比较结构、挑选构象、拆成独立目录，并保持原子编号与能量的对应关系。本篇用三份真实GFN2-xTB优化的乙醇结构，加一份完全重复的记录，演示这些操作。重复只用作文件处理对照。

下载输入、原生工具调用与绘图程序

1. 先跑通完整小工作流#

将资料包保存为crest-tools.zip：

unzip crest-tools.zip -d crest-tools
cd crest-tools
sinfo
sbatch job.slurm

模板使用CPU分区，提交前确认资源可用。程序最多并发8个单线程工具。环境为官方GNU12预编译CREST3.0.2、xTB6.7.1；绘图用sisl0.16.4环境中的NumPy和Matplotlib。输入与结果分别在：

ls cases
cat analysis/rmsd.csv
ls figures

下面按用途说明每个命令。需要单独重跑时使用新的工作目录，并在已申请的计算资源内执行，避免覆盖整包原结果。

2. RMSD：整体、重原子与刚体运动#

module load crest/3.0.2
cd cases/rmsd-conformer
crest --rmsd reference.xyz variant.xyz
crest --rmsdheavy reference.xyz variant.xyz
cd ../rmsd-rotated
crest --rmsd reference.xyz rotated.xyz
cd ../..

两份乙醇构象的整体RMSD为0.41092745Å，重原子RMSD为0.03390029Å；对同一结构只做旋转和平移，结果为0Å（输出精度内）。重原子RMSD小，说明这里更明显的变化来自氢原子及O–H方向。

RMSD在最佳刚体叠合后计算，不能把文件中坐标的逐项差直接当作RMSD。本包用独立Kabsch算法核对。3.0.2的这两个命令要求原子顺序一致；重新排列原子后应先恢复编号对应关系。

3. 拆分XYZ时保留能量对应表#

mkdir split-example
cp cases/split/ensemble.xyz split-example/
cd split-example
crest --splitfile ensemble.xyz
ls SPLIT
python ../split_energies.py
cat split-energies.csv
cd ..

本例得到SPLIT/STRUC1/struc.xyz至STRUC4/struc.xyz，坐标与原始四帧逐一对应。该版原生拆分文件的注释行为空，不会保留原XYZ注释里的能量；能量仍在输入ensemble.xyz中。配套split_energies.py按原帧序号写CSV，将结构文件、原始能量/Eh和原注释关联起来，没有重新计算或填造能量。

如果只想拆指定范围，可以在文件名后加两个序号。对自己的数据，先检查输入的能量注释格式；本包采用xTB的energy:格式，不把任意注释中的数字都当作能量。

4. 比较两个集合，而不只比较两份XYZ#

cd cases/compare
crest reference.xyz --compare ensemble.xyz subset.xyz
cd ../..

本例A集合有4条结构记录，B集合有2条；CREGEN按本例阈值分别归为2个构象家族，随后得到两组对应关系。输出包含家族间RMSD矩阵，对角元素为0，非对角约0.41093Å。记录数、旋转异构体数和构象家族数含义不同，不能只看文件里有多少帧。

这里需要把参考坐标放在命令开头，两个集合文件放在--compare后。RMSD阈值、本例的能量阈值与旋转常数判据会影响分组；“没有被归为同一家族”并不自动代表发生化学反应。

5. ZSORT与resortensemble究竟做什么#

cd cases/zsort-resort
crest reference.xyz --zsort -T 1
ls zcoord zmatrix1.zmat
cd ../..

ZSORT按分子内部连接重排原子，便于建立一致的Z-matrix。3.0.2写出Turbomole格式的zcoord，没有自动保存.atorder。本包的build_atom_map.py根据元素和原子间距离，从原XYZ与实际zcoord恢复对应关系，再明确生成.atorder。映射采用该版原生换算1 bohr=0.52917726Å，坐标单位不能混用。

--resortensemble的作用是按已有映射恢复原子顺序，不是按能量排序。整包已准备好一个按ZSORT次序排列的四帧集合和映射，可单独体验恢复：

mkdir restore-example
cp cases/zsort-resort/.atorder restore-example/
cp cases/zsort-resort/ordered-ensemble.xyz restore-example/
cd restore-example
crest --resortensemble ordered-ensemble.xyz
head ensemble_tmp.xyz
cd ..

本例恢复后四帧的元素顺序、坐标与原能量注释逐一对应，坐标差小于1e-9 Å。换成自己的数据时，映射必须来自同一套原子编号；不能拿另一体系的.atorder套用。

6. 对称性与拓扑#

cd cases/symmetries
crest --symmetries ensemble.xyz
cat symmetries
cd ../topology
crest reference.xyz --testtopo reference.xyz
cd ../..

本例第三帧反式乙醇识别为Cs，其余记录为C1。对称性对几何微小畸变和算法容差敏感，处理自己的体系时应保留实际坐标。

拓扑输出列出9个原子的连续配位数和邻居编号，本例共有8条共价键、1个分子片段。连续配位数可能是4.0010、3.9969之类，不是整数形式价态，也不是原子电荷。

7. 重画图并换成自己的数据#

module load sisl/0.16.4
python analyze.py

先保持原子编号和元素组成一致，再依次做比较、拆分和结果整理。需要重新计算单点、优化结构或汇总电子能量权重，继续看CREST单点、优化与权重。命令功能参考CREST官方关键词说明；其中标为3.1的新功能不能直接当作3.0.2已有功能。