MolCryst-MLIPs：分子晶体机器学习原子间势数据库

《Journal of Chemical Theory and Computation》：MolCryst-MLIPs: A Machine-Learned Interatomic Potentials Database for Molecular Crystals

　　

研究人员提出了一个名为MolCryst-MLIPs的开放分子晶体（MC）机器学习原子间势（MLIPs）数据库。首个版本包含针对九个分子晶体体系（苯甲酰胺、苯甲酸、香豆素、均四甲苯、异烟酰胺、烟酸、烟酰胺、吡嗪酰胺和间苯二酚）的微调MACE模型，这些模型是利用自动化机器学习流水线（AMLP）开发的，该流水线将整个MLIP开发工作流程（从参考数据生成到模型训练和验证）简化为可复现且用户友好的流程。模型从MACE-MH-1基础模型（omol头）进行微调，在所有体系中实现了平均能量MAE为0.141 kJ·mol^–1·atom^–1、平均力MAE为0.648 kJ·mol^–1·?^–1。在DFT标记的多晶型集合上对三个最先进的基础模型进行基准测试，只有微调后的模型能够可靠地识别最稳定的多晶型，平均Kendall τ = 0.397表明与完整DFT能景存在中等秩相关性。通过能量守恒、P2取向序参数和径向分布函数评估的动力学稳定性和结构完整性，通过分子动力学模拟进行了评价。所发布的模型和数据集构成一个不断增长的经过验证的MLIPs开放数据库，可用于在不同热力学条件下对目标化合物的多晶型能景进行生产级MD模拟。

**MolCryst-MLIPs：面向分子晶体多晶型研究的机器学习原子间势数据库解读**

**一、研究背景与意义**

多晶型现象是指同一化合物结晶形成多种不同结构的现象，在分子晶体（MCs）研究领域具有核心重要性，尤其对药物开发影响深远——同一化合物的不同多晶型在溶解度、熔点和生物利用度等方面可能存在显著差异。分子晶体依靠非共价相互作用维系，包括π–π堆积、C–H···π接触和氢键，这些相互作用共同决定了晶体所采用的堆积排列方式。准确模拟多晶型现象具有内在挑战性，因为晶体稳定性主要由这些分子间非共价力主导，必须以足够精度加以描述，才能解析决定竞争晶型相对稳定性的分子内与分子间相互作用的微妙关联。多晶型之间的晶格能差异通常仅为几kJ·mol^–1，且互变体系可能随温度变化出现稳定性交叉，这使得多晶型排序成为检验势函数精度的特别敏感测试。经典力场虽计算高效，但往往缺乏区分此类能量相近晶型或在变化热力学条件下可靠捕捉动力学行为所需的分辨率。量子力学（QM）方法如密度泛函理论（DFT）原则上可获得高保真结构和能量数据，但DFT高昂的计算成本——尤其是多晶型可靠排序所需的严格收敛标准——限制了其在从头算分子动力学（AIMD）中仅能处理小体系和短时间尺度。

机器学习原子间势（MLIPs）作为弥合这一鸿沟的有力替代方案应运而生：通过从QM参考数据中学习，其精度可与DFT媲美，同时以极低计算成本实现分子动力学（MD）模拟。MLIP开发面临三个关键选择：模型架构的选择、训练数据集的组成与覆盖范围、以及从头训练还是微调现有基础模型。这些选择深度关联：架构决定数据需求和精度上限，而所需训练数据的数量和多样性又取决于从头训练还是利用现有基础模型。从头训练可完全控制势能面（PES），但需要大型精心策划的数据集和大量计算投入。微调则利用已在广泛化学空间上训练的基础模型所习得的表征，大幅降低数据需求和训练成本，同时保持目标体系的高精度。在近期基础模型中，UMA和Orb旨在广泛覆盖分子与材料体系，而MACE多头基础模型（MACE-MH-1）及其omol头——在涵盖分子、有机和有机金属体系的OMOL数据集上训练——是分子晶体应用的自然选择。尽管MLIPs在原子和分子体系中取得显著成功，分子晶体仍构成独特挑战：精确建模需要同时捕捉强分子内相互作用和支配多晶型稳定性的微妙非共价力。大多数现有基础模型主要基于气相分子或无机固体训练，周期性晶体环境代表性不足。OMC25数据集提供了超过2700万个DFT标记的分子晶体结构，UMA纳入了源自该数据集的专用分子晶体训练子集，MACE-MH-1通过跨域微调在X23分子晶体基准上取得现有基础模型中最低MAE。

尽管有这些进展，仅凭基础模型可能无法达到多晶型排序或自由能计算等苛刻任务所需的精度。例如，在吖啶多晶型家族的近期研究中，MACE-MPA基础模型未能重现实验确立的稳定性排序，而体系特异性微调模型则恢复了正确排序。在体系特异性高质量DFT数据上进行微调提供了弥合这一差距的实用途径，可提高目标化合物多晶型能景上的精度和迁移性。然而，组装必要的训练数据、选择适当的超参数以及验证跨多个多晶型的动力学稳定性仍是密集流程，通常需要大量专业知识——这对更广泛社区构成实际障碍。为此，研究人员利用自动化机器学习流水线（AMLP）系统地对MACE-MH-1进行九个高多晶型分子晶体体系的微调，证明从数据集生成到模型验证的完整工作流程可以可复现且易操作的方式执行，进而建立了MolCryst-MLIPs开放数据库。

**二、关键技术方法**

研究人员采用AMLP流水线统筹所有工作流程，涵盖DFT输入生成、主动学习、数据集策展、MACE训练配置和训练后验证。实验晶体结构（.cif文件）来自剑桥结构数据库（CSD）。DFT参考数据使用维也纳从头算模拟包（VASP）生成，采用Perdew–Burke–Ernzerhof（PBE）交换关联泛函结合Grimme的D4色散校正，平面波截断能650 eV，电子收敛标准EDIFF = 10^–7 eV。AIMD模拟在25、300、400和500 K四个温度下进行，辅以主动学习循环识别模拟失败结构并补充至训练集。模型微调采用MACE-MH-1基础模型（omol头），遵循包含随机权重平均（SWA）的两阶段训练协议，使用Adam优化器配合AMSGrad和早停机制。验证方案涵盖静态与动态两个范畴：结构保真度通过DFT与MACE优化几何对比和RMSD15堆积相似性评估；动力学鲁棒性通过微正则（NVE）能量守恒测试、正则（NVT）系综模拟、以及恒压（NPT）全柔性晶胞模拟评估。NPT模拟选取苯甲酸、烟酸和间苯二酚三个代表性体系，分别代表羧酸二聚体、氢键两性离子网络和多羟基三维网络等不同氢键 motifs。

**三、研究结果**

**3.1 参考数据生成**

基于CSD检索的65个晶胞松弛任务均通过AMLP自动生成一致DFT参数输入。物理有效性验证显示：九个体系的多晶型相对晶格能范围为0.09 kJ mol^–1（均四甲苯）至4.64 kJ mol^–1（间苯二酚），平均值2.19 kJ mol^–1，在预期的小于10 kJ mol^–1范围内；DFT优化产生的体积收缩为?3.0%至?6.4%（均值?4.7%），应变张量范数为0.025至0.081，符合0 K计算相对室温实验结构的物理预期。AIMD和主动学习共生成113,953个结构，能量与力间呈现强线性相关（Pearson r = 0.907–0.949）。数据集经DBSCAN聚类分析、力幅值10.0 eV/?截断和能量阈值0.2 eV/atom三重过滤后，以85/15比例随机分割为训练集和验证集。

**3.2 MACE模型微调与评估**

九个体系的微调模型在验证集上平均能量MAE为0.141 kJ mol^–1 atom^–1，力MAE为0.648 kJ mol^–1 ?^–1。在相对晶格能（ΔE_latt）排序测试中，MACE-MH-1基础模型（omol头）对香豆素和吡嗪酰胺呈现近乎平坦的能景，对烟酰胺则给出定性错误的稳定性排序；Orb模型在多数体系中无法重现DFT排序；UMA模型虽恢复多数体系的正确排序，但多晶型能量差距落于其自身MAE范围内而无法真正区分。相比之下，微调模型在全部四个代表性体系中恢复正确的DFT稳定性排序，包括纠正异烟酰胺中omol模型完全反转的排序。跨九个体系的综合指标显示，微调模型平均ΔE_latt MAE为0.98 kJ mol^–1，Kendall τ = 0.397，最稳定多晶型罚分为0.51 kJ mol^–1，均优于三个基础模型（UMA τ = 0.145，omol τ = ?0.044，orb τ = ?0.250）。微调模型在9个体系中的5个实现零罚分，即精确识别DFT最稳定多晶型。对105个CSD多晶型（含40个因晶胞过大而从训练集排除的晶型）的迁移性测试表明，平均RMSD15为0.16 ?，94%结构低于0.3 ?。

**3.3 分子动力学动力学稳定性**

NVE模拟中，全部九个模型展现优异能量守恒，累计漂移保持在10^–7量级。NVT模拟中，通过P2取向序参数监测：吡嗪酰胺的15个晶型（含训练外样本PYRZIN05）全部保持取向有序，P2均值?0.5对应垂直排列（PYRZIN14、PYRZIN18、PYRZIN23），近零对应随机取向无序（PYRZIN17），0.1–0.2对应人字形堆积（PYRZIN05），0.4以上对应平行堆积（其余晶型）；苯甲酸在600 K时BENZAC20和BENZAC22失去结构完整性，符合其热稳定性极限；香豆素的COUMAR14在600 K发生向不同堆积 motif 的重新取向。径向分布函数（RDFs）分析确认C–N和C–O对的第一峰保持尖锐（分子内键合完整性），O–O对峰随温度升高而展宽但整体峰位和结构 motif 保持（分子间堆积保持）。

**3.4 恒压稳定性**

NPT模拟（300 K，环境压力，全柔性晶胞）中，苯甲酸、烟酸和间苯二酚的全部多晶型在100 ps生产动力学中保持稳定。晶胞体积围绕稳定均值波动且无系统性漂移，生产运行前半段与后半段的平均V/V0差异不超过0.5%。三个晶胞角保持在初始值几度以内，无剪切失稳或堆积转变。时间平均体积接近0 K参考值：苯甲酸膨胀2.5–3.3%，烟酸和间苯二酚保持在V0的0.5%以内，该膨胀小于DFT晶胞相对室温实验结构的平均收缩（?4.7%）。

**四、讨论与结论**

讨论部分指出，微调基础模型对于分子晶体多晶型的体系特异性应用至关重要：基础模型缺乏区分同一化合物多晶型的分辨率，而针对DFT数据的定向微调可大幅提升判别能力，可靠识别最稳定多晶型并改善全能景的稳定性排序。微调模型可靠地迁移至训练中未包含的实验晶体结构——包括DFT水平上因成本过高而无法松弛的大晶胞多晶型——跨105个多晶型以平均RMSD15 = 0.16 ?重现实验堆积。需强调的是，这是化学体系内而非跨体系的迁移性：对化学上不同的化合物部署需要微调新模型，AMLP为此提供了协议。NVE、NVT和NPT三层动力学验证共同确认了模型在固定体积和自由变形条件下的结构稳定性。残余局限在于势函数精度受PBE + D4参考水平约束：对近简并晶型，交换关联泛函的不确定性与待解析能量差相当，因此模型重现的是DFT层级而非精确实验层级。配套DFT数据集与势函数一同发布，为未来任何新型基础模型提供了可直接用于分子晶体域微调的资源，无需重新生成昂贵的DFT参考数据。MolCryst-MLIPs数据库作为首个经过验证的、可直接用于有机晶体多晶型大规模MD模拟的势函数版本，为晶体结构探索和有限温度自由能计算提供了有力工具，其中相对AIMD降低的计算成本带来决定性实用优势。