《Entropy》：Computational Methods for Molecular Dynamics of Supercooled Water Between 200 and 273 K
编辑推荐：
水在接近并低于其熔点时表现出异常的热力学和结构行为。对这一行为进行计算建模仍具有挑战性：随着温度降低，核量子效应（NQEs）变得愈加显著，采样效率下降，计算方法的选择对预测结构和动力学性质的准确性产生关键影响。本综述旨在为从业者提供在低温下进行可靠水分子动力学
水在接近并低于其熔点时表现出异常的热力学和结构行为。对这一行为进行计算建模仍具有挑战性：随着温度降低，核量子效应（NQEs）变得愈加显著，采样效率下降，计算方法的选择对预测结构和动力学性质的准确性产生关键影响。本综述旨在为从业者提供在低温下进行可靠水分子动力学模拟的实用指南，重点放在过冷区（200–273 K）。研究人员考察了常用的从头算分子动力学（AIMD）方法，并批判性评估了通过路径积分分子动力学（PIMD）、环聚合物分子动力学（RPMD）和质心分子动力学（CMD）纳入核量子效应的策略。在过冷区冷却过程中，核量子效应的相对重要性稳步增加，因此将原子核经典处理所引入的误差增大，并影响所关注的结构和动力学性质。密度泛函近似对水结构的已知局限性是在其与NQE处理相互作用的背景下评估，而非孤立地评估。通过讨论实际考虑因素，包括采样效率和路径积分珠子数随温度的标度，以及方法本身，本综述旨在为过冷水的可靠AIMD研究提供参考。
1. 引言
水在接近并低于熔点时表现出异常的热力学与结构行为。过冷水在云化学、低温保存、食品加工和行星科学中具有重要作用。经典分子动力学（MD）无法描述涉及电子结构的化学反应，刚性非极化力场又忽略电子极化与核量子效应（NQEs）。从头算分子动力学（AIMD）在低温下面临采样效率下降、核量子效应增强以及密度泛函精度不足三重挑战。本综述聚焦200–273 K过冷区，为可靠AIMD模拟提供方法选择、采样收敛与实验验证指南。
2. 水在低温下的结构
2.1 过冷液态水
过冷水可在微米尺度液滴中过冷至约235 K的均质冻结极限，快速冷却纳米液滴和超快X射线探测已将亚稳液态水延伸至227 K。冷却过程中第一氧–氧径向分布函数（RDF）峰升高，约3.3 ?处存在温度无关的等吸收点，第二峰增强，表明四面体配位增强。这些结构变化常以低密度液体（LDL）和高密度液体（HDL）双态框架解释，并涉及液–液转变以及“无主之地”温区等未决问题。
2.2 结晶冰相
冰Ih是200–273 K、常压下的热力学稳定相。冰Ic在实际样品中多为含六方层错堆垛的堆垛无序冰Isd，其密度和结构参数因此可变。小液滴冷冻时，较小液滴倾向于形成含堆垛层错的立方冰，较大液滴则更倾向六方冰。
3. 水的从头算分子动力学方法
3.1 经典原子核方法
3.1.1 Born–Oppenheimer分子动力学（BOMD）
BOMD在每个核构型下自洽收敛电子结构，核在基态势能面上按牛顿方程演化，可精确守恒总能量。由于原子核被经典处理，该方法无法描述零点能与隧穿效应，温度越低、热de Broglie波长越大，误差越严重。BOMD已广泛实现于CP2K、VASP和Quantum ESPRESSO等软件。
3.1.2 Car–Parrinello分子动力学（CPMD）
CPMD将电子波函数作为虚质量动力学变量，避免逐电子步自洽场收敛，但依赖电子与核运动的绝热解耦。水中O–H伸缩频率很高，迫使采用小时间步；低温下核运动变慢，维持绝热解耦更加困难，容易偏离Born–Oppenheimer势能面。对于需要考虑核量子效应的低温水模拟，标准CPMD通常不是首选。
3.2 路径积分分子动力学（PIMD）
PIMD利用量子统计力学与环聚合物之间的同构关系，将每个量子粒子表示为P个由谐波弹簧连接的珠子。收敛要求珠子数P随β增大而增大，因此低温下所需珠子数上升。标准PIMD对每个珠子独立计算力，成本随P线性增加；从300 K冷却至200 K时，固定精度下的最低成本约增加1.5倍。基于q-TIP4P/F模型的测试表明，密度、等温压缩率和静态介电常数对240–300 K范围内的珠子数不敏感，但等压热容对珠子数敏感，且随P增大反而偏离参考值，说明不同可观测量的珠子数收敛行为不能简单通用。
3.2.2 加速路径积分方法
环聚合物收缩（RPC）只在收缩后的环上计算昂贵势能力；PI-GLE和PIGLET通过有色噪声热浴加速路径积分收敛，可在300 K仅用约6个珠子、230 K用约8个珠子获得收敛的核量子效应。增强采样方法如伞形采样、元动力学、seeding和forward-flux sampling可用于加速构型采样和计算成核速率。
3.3 实时近似量子动力学
环聚合物分子动力学（RPMD）近似Kubo变换的实时相关函数，但会产生虚假共振峰；热浴化RPMD（TRPMD）可抑制该伪影。质心分子动力学（CMD）基于Feynman路径质心框架，但对O–H伸缩带产生人工红移和展宽，即曲率问题；准质心分子动力学（QCMD）通过曲线坐标消除该问题。环聚合物瞬子理论用于低温隧穿速率常数计算。
3.4 机器学习原子间势
机器学习原子间势（MLPs）以从头算数据训练，可在极低计算成本下替代势能面，并与BOMD、PIMD、RPMD或CMD结合，使纳秒级路径积分轨迹成为可能。但多数MLPs在环境温度附近训练，迁移至过冷区并不保证可靠，必须在目标温度和密度下用实验或从头算参考进行明确验证。
3.5 AIMD方法比较与方法选择
方法选择取决于目标可观测量、所需精度和温度区间。经典核AIMD适用于核量子效应修正小于其他建模误差的情形；PIMD或PIGLET用于需要量子统计时；RPMD、CMD/QCMD和瞬子方法用于需要实时相关函数、振动谱或速率常数时。各方法边界并不截然，应针对具体性质进行基准验证。
3.6 低温下的实际采样与收敛
冷却使结构弛豫时间迅速增长，轨迹长度是否足以代表平衡构型必须用分块平均等方法检验。过冷区还存在结晶与平衡的竞争，但结构弛豫时间通常远短于结晶时间，因此可在监测是否出现冰核的前提下实现亚稳液态平衡采样。
3.7 结构与热力学验证目标
仅用氧–氧RDF不足以验证结构，应同时检验氧–氢和氢–氢偏分布函数，并直接在倒空间比较结构因子S(Q)与实验值。环统计可用于诊断过冷液体中的氢键网络。该温区内关键实验参考点是等压热容最大值，超快量热法将其定位于约229 K。
4. 水的密度泛函选择
密度泛函选择是DFT-AIMD中影响最大的决定之一。目前多数泛函基准在环境条件附近获得，过冷区内的直接DFT验证仍稀疏，因而泛函排名主要依赖环境温度数据。
4.1 GGA泛函
BLYP无色散校正时过度结构化液体并低估扩散和密度，加入D3色散校正后显著改善。PBE也倾向过度结构化，D3校正的改进有限。revPBE-D3在环境温度经典AIMD中表现良好，部分来自泛函固有误差与忽略核量子效应之间的抵消，其在200–273 K的适用性仍属推断而非直接验证。
4.2 meta-GGA泛函
SCAN满足多个精确约束，能较好描述冰的多形体相对稳定性，但将液态水密度高估约5%，且存在数值网格敏感性。r2SCAN数值稳定性改善，但仍需在约350 K的高温下模拟以匹配300 K实验结构。密度校正的DC-SCAN可接近耦合簇精度，但通常需借助多体势而非直接AIMD。
4.3 杂化泛函
PBE0含25%精确交换，能改善反应势垒和自相互作用误差，但经典核下易过度结构化液态水。revPBE0-D3结合路径积分模拟可定量再现液态水的O–O RDF、自扩散系数和红外/拉曼光谱，是少数与实验定量一致的方案之一。B97X-V和B97M-V在水团簇基准中误差最小，相关meta-GGA泛函B97M-rV以meta-GGA成本达到接近杂化泛函的精度。
4.4 色散校正
常用方案包括D3(BJ)、D4和VV10非局域关联泛函。D3(BJ)因验证充分、开销低且代码普及，目前是水模拟的默认选择。
5. 核量子效应
5.1 核量子效应的物理起源
Born–Oppenheimer近似将电子与核运动分离，标准MD进一步将核当作经典粒子，因而忽略零点能、量子隧穿和交换效应。O–H伸缩模式的零点能与氢键强度相当，导致质子显著离域。沿共价键轴的量子涨落增强氢键，垂直方向的涨落则削弱氢键，二者竞争使净效应因体系和性质而异。水中质子局域在深共价势阱中，交换效应可忽略。
5.2 同位素取代作为NQE的实验探针
H/D同位素效应是验证量子核模拟的主要实验途径。液态水和冰Ih中，H/D取代对氧–氧距离的次级几何同位素效应仅约0.001 ?，表明竞争性量子效应几乎完全抵消。同位素分辨的几何、质子动能和分馏比可同时严格检验电子结构方法和核处理方案。
5.3 NQE的结构后果
经典核模拟会过度结构化O–H?O氢键；路径积分模拟中质子离域展宽分子内O–H峰并降低第一O–O峰，使其更接近衍射实验。中子Compton散射测得的质子动量分布也支持路径积分描述。NQE还使O–H伸缩振动相对经典模拟产生约200 cm
?1的红移。
5.4 经典与量子原子核：与实验的匹配比较
在q-TIP4P/F势能面上，PIMD密度在平衡态和过冷区均与实验接近；经典与量子核密度差最多约0.5%。自扩散系数方面，量子核在260、280和300 K更接近实验，250 K时经典核反而更接近；量子核在整个温区偏差不超过约8%，经典核在260–280 K低估约15%。在revPBE0-D3势能面上，密度对核处理约偏离1%，说明核量子效应的抵消程度具有模型依赖性。因此，若目标仅为密度且势能面已验证，经典核在200–273 K可用；但针对RDF、质子动量分布、同位素效应、自扩散系数和等压热容，经典核不是最优选择。
6. 结论与展望
可靠低温模拟需将核处理精度与电子结构精度相匹配。BOMD是核量子效应较小时的合理基线；CPMD在低温精度上较弱；PIMD及其加速方案用于量子统计；RPMD、CMD/QCMD用于振动谱和速率但各有伪影。泛函选择不能与核处理脱钩，revPBE-D3在环境温度的良好表现部分依赖误差抵消，meta-GGA模拟常需升高温度补偿，revPBE0-D3结合路径积分是目前较可靠的组合。过冷区的特殊性在于零点能与热能之比随降温增大，错误抵消在环境温度校准后未必成立。机器学习势使纳秒级路径积分模拟成为可能，但其过冷区迁移性必须验证。这一方法路径有望为云中冰核形成、低温保存和复杂有机分子生成等应用提供定量预测能力。
