AI力场二次开发教程(03):第一个AI力场分子——SMILES 到 OpenMM System 全流程
AI力场二次开发教程3第一个AI力场分子——SMILES 到 OpenMM System 的全流程复现版本声明本文运行环境为第 02 篇所建的espaloma0.3.2openff-toolkit0.19.0OpenMMPython 3.10本文完整复现 espaloma 官方部署示例锚点 A所有代码以 espaloma 仓库 README 的公开脚本为真源不涉及任何未公开接口。目标是让你不靠“看教程”而靠“跑通代码”获得第一个由 AI 力场参数化的 OpenMM System。一句话结论将Molecule.from_smiles(CN1CNC2C1C(O)N(C(O)N2C)C)依次交给esp.Graph(molecule)、esp.get_model(...)、model(...)前向、esp.graphs.deploy.openmm_system_from_graph(molecule_graph)即可在数行代码内拿到 OpenMMSystem其原子由 GNN图神经网络推理出的键长、键角、二面角与非键参数定义。〇、本篇要解决的认知问题esp.Graph(molecule)到底把分子变成了什么样的 “异构图”heterograph里面装了什么字段esp.get_model(latest)加载的模型执行model(heterograph)时GNN 完成了怎样的“图消息传递”openmm_system_from_graph是如何把网络输出“落回”经典势函数、组装成 OpenMMSystem的生成的System里通常有哪些能量项Bond/Constrain/Angle/Proper/RB/Nonbonded/GBSA…怎么数它们从 SMILES 到 System哪一个环节可能因为分子化学如单价异常、环标注中断如何处置一、机制解析1.1esp.Graph从化学分子到 DGL 异构图Espaloma 依赖 DGLDeep Graph Library把分子建模为异构图heterograph。之所以“异构”是因为同一张图里面有多种边类型原子-原子键、原子-原子-原子角、三中心二面角等分别对应不同的遁址类型。esp.Graph(molecule)会为每个原子生成n1节点、为每条键生成n2、为每个角生成n3、为每条二面角生成n4级别的初始特征边则包括n1-n2、n2-n3等映射。molecule_graph.heterograph即该 DGL 图对象可直接传给神经网络前向。具体的字段名如n1、n2、f1、f2是 espaloma 内部约定的图 schema你应通过print(molecule_graph.heterograph.ntypes)与print(molecule_graph.heterograph.etypes)现场验证不要硬背。1.2 GNN 前向消息传递 读出新参数model(graph)是EspalomaModel对象的前向。它进行若干轮消息传递每个节点根据邻居隐向量更新自身表示最终用readout层为每条键/角/二面角/原子输出一套数值。关键点它不枚举原子类型而是把化学环境编码进连续向量——这正是第 01 篇反复强调的“无缺参范式”的直接体现。1.3openmm_system_from_graph把推理结果变成 Systemesp.graphs.deploy.openmm_system_from_graph(molecule_graph)是部署层函数。它读取molecule_graph.nodes上已推理出的参数用 OpenMM 的ForceField/SystemBuilder语义组装出一个System。你会在 System 上看到 standard 力学项键合、角度、二面角、非键 vdW/静电其力常数与平衡值均来自 GNN 输出而非预定义查表。内置的OpenMM项与我们期望的势函数一一对应HarmonicBondForce键、HarmonicAngleForce/PeriodicTorsionForce角/二面角、NonbondedForcevdW点电荷。因为 Espaloma 0.3.x 主要输出全原子 Lennard-Jones 与键合项多数生成的 System 里能量项数量一般在 4–8 个之间含可选的 PositionalRestraint/GBSA。真实数量以运行环境与所用内置模板为准配合下方计数代码即可获得精确数字。二、完整代码与逐行剖析2.1 端到端咖啡因全流程锚点 A 原样复现# espaloma 官方部署示例锚点A原文首次运行会下载权重fromopenff.toolkit.topologyimportMoleculeimportespalomaasesp# step1 : 从 SMILES 构造分子moleculeMolecule.from_smiles(CN1CNC2C1C(O)N(C(O)N2C)C)# 咖啡因print(分子原子数,molecule.n_atoms)# step2 : 构造异构图molecule_graphesp.Graph(molecule)print(图节点类型,list(molecule_graph.heterograph.ntypes))print(图边类型 ,list(molecule_graph.heterograph.etypes))# step3 : 加载模型并前向espaloma_modelesp.get_model(latest)# 网络拉取官方权重espaloma_model(molecule_graph.heterograph)# 触发GNN推理print(推理完成。二面体张量字段读取示例,list(molecule_graph.nodes.keys()))# step4 : 部署成 OpenMM Systemopenmm_systemesp.graphs.deploy.openmm_system_from_graph(molecule_graph)print(System 已生成。)逐行剖析Molecule.from_smiles来自openff.toolkit.topology它把 SMILES 解析为带价态/立体/环信息的Molecule详见第 04 篇。molecule.n_atoms让你直观核对图规模。esp.Graph(molecule)内部用 RDKit/OpenFF 信息构建 DGL 异构图。打印ntypes/etypes能就地确认图 schema如是否含n1/n2/n3与f1/f2。esp.get_model(latest)拉取默认权重若内网受限可用本地espaloma-0.3.2.pt路径代替见第 02 篇。openmm_system_from_graph签名按官方deploy模块公开接口若签名有差异请inspect.signature(esp.graphs.deploy.openmm_system_from_graph)现场核实。2.2 数一数 System 里的能量项# 在 2.1 之后续跑枚举 System 中的 Forcen_forcesopenmm_system.getNumForces()foriinrange(n_forces):fopenmm_system.getForce(i)print(i,force_to_name(f),f 类型{type(f).__name__})defforce_to_name(f):# 用 OpenMM 类的属性尽量给出人类可读标签returngetattr(f,getName,lambda:type(f).__name__)()逐行剖析openmm_system.getNumForces()/getForce(i)是 OpenMMSystem的标准遍历方式。每个Force子类如HarmonicBondForce、NonbondedForce常通过getName()自带标签若调用失败退化为type(f).__name__。你应看到类似HarmonicBondForce / HarmonicAngleForce / PeriodicTorsionForce / NonbondedForce的一串把它们与论文/文档中的“键、角、二面角、非键”一一对应。说明上文force_to_name里的getName与getNumForces等均为 OpenMM 官方公开 API不同 OpenMM 小版本对getName的支持以官方文档为准。2.3 换分子快速跑复制即跑打印能量项个数fromopenff.toolkit.topologyimportMoleculeimportespalomaasespdefparam_smiles(smiles,weightlatest):molMolecule.from_smiles(smiles)gesp.Graph(mol)modelesp.get_model(weight)model.eval()model(g.heterograph)# 推理sysesp.graphs.deploy.openmm_system_from_graph(g)returng,sysforsmiin[c1ccccc1,CC(O)O,CCO]:# 苯/乙酸/乙醇_,sysparam_smiles(smi)print(f{smi}: 能量项个数 {sys.getNumForces()})逐行剖析param_smiles封装了“图→推理→System”三步示范如何复用。model.eval()见第 02 篇保证推理口径一致。循环三种常见分子验证 Espaloma 对不同官能团都给出相同路径的结果直观体现“无缺参”。三、常见报错与排查报错现象可能原因处理Molecule.from_smiles抛InvalidValenceTemplateExceptionSMILES 价态异常改用 RDKit 生成合理 SMILES或用带环符的标准表示KeyError: n1或TypeError于 heterograph图 schema 与所用权重版本不匹配list(g.heterograph.ntypes)核对字段确认权重与依赖匹配openmm_system_from_graph报ValueError: Bonds span topologyvalence 或显式 H 数与图不一致用Molecule.from_pdb/from_rdkit清洗后再to_openff推理后graph.nodes无参数张量前向顺序被跳过或模型未 eval确认执行了model(g.heterograph)且model.eval()OpenMMValueError: Non-zero preconditioner非键项 NaN/负 LJ检查分子几何合理性必要时先generate_conformers再验以上错误处理是接口层面的通用排查Espaloma 具体异常类名以官方仓库 issue 与文档为准。四、动手练习三分子对比依次用 §2.3 处理c1ccccc1、CC(O)O、CCO打印各自的能量项列表与个数写一张 3 行对比表。图 schema 画图选一个分子用print(molecule_graph.heterograph.ntypes)与etypes列出所有节点/边类型再手工画出它对应的 DGL 异构图ASCII 即可标出n1→n2→n3的路径。本地权重版把espaloma-0.3.2.pt下载到本地改写 §2.1 用本地路径调用esp.get_model验证与latest输出在能量项个数上一致。合成一个分子触发失败尝试Molecule.from_smiles([Li])或含过渡金属的表示观察在哪一步出错把异常文本写进你的笔记并说明为何这类体系需要专门处理。五、小结与下一篇预告你现在已经亲手跑通了“SMILES → 异构图 → GNN 推理 → OpenMM System”的完整链路理解了esp.Graph的异构字段与openmm_system_from_graph的部署语义。更重要的是你验证了不同官能团分子都走同一路径——这正是 AI 力场“无缺参”的关键证据。下一篇我们要把视野从 Espaloma 拉回到更通用的“分子准备”层第 04 篇《SMIRNOFF 与 OpenFF Toolkit 入门——分子准备》深入 SMIRKS 子结构语法与Molecule从 SMILES/SDF/PDB 的读入、原子键二面角属性与find_matches的使用为自建力场与覆盖率分析打基础。本篇认知问题回显FAQQesp.Graph(molecule) 是如何把分子构造成异构图其中有哪些字段A它把分子抽象为 DGL 异构图用 n1/n2/n3 等节点类型对应原子、键、角等化学语境f1/f2 等张量承载特征边类型映射相邻节点具体的 ntypes/etypes 可用g.heterograph.ntypes和etypes打印查验。Qesp.get_model 加载的模型执行 model(heterograph) 时 GNN 完成了什么图消息传递A各节点沿连边反复接收并聚合邻居隐向量更新自身表示传 N 轮后再由 readout 头为键、角、二面角、原子分别输出数值参数最终用一次前向完成对整张图的参数化。Qopenmm_system_from_graph 是如何把网络输出落回经典势函数生成 OpenMM System 的A它读取异构图节点上已推理出的键长、键角、二面角与非键参数按 OpenMM 力场语义组装出 System其中各物理项分别对应 HarmonicBondForce、HarmonicAngleForce、PeriodicTorsionForce 与 NonbondedForce。Q生成的 OpenMM System 中通常有哪些能量项怎么统计个数A多为键合项、角度、二面角与非键 vdW/静电项共约 4–8 个用system.getNumForces()再遍历getForce(i)逐一打印type(f).__name__即可统计种类与数量。Q从 SMILES 到 System 全流程中哪一步会因分子化学特性中断如何处置AMolecule.from_smiles 或因价态异常抛 InvalidValenceTemplateException异构图构造或因显式氢与环标注不一致而报错可用 RDKit 清洗 SMILES 并保证价态/环符规范后再生成 System。