这是一篇面向计算生物与分子模拟入门者的分子动力学（MD）方法导读。文章承接前文蒙特卡洛（MC）方法，指出MC因无真实时间维度无法计算动力学性质，由此引出MD的核心——基于牛顿第二定律对时间做离散化积分来复现分子运动轨迹，并阐明多体系统混沌特性导致长时轨迹必然偏离，但统计力学的系综平均保证计算结果的有效性。

文章系统梳理了MD全流程要点：基础层介绍Verlet家族积分器演化，点明Velocity-Verlet因可显式监控温度成为主流软件默认算法；工程层讲解Nose-Hoover链恒温器和Parrinello-Rahman恒压器以复现NVT/NPT系综，对比指出早期Andersen与Berendsen恒温器的缺陷及1992年Nose-Hoover链问世是MD动力学速率计算可信的分水岭，同时归纳截断半径、邻居列表、能量最小化等算力优化手段。应用层说明绝对自由能不可算而只算自由能差，解析FEP、伞形采样、元动力学的适用场景，强调集体变量（CV）选择的关键性，并对比弦方法、Path-MetaD、TAPS及过渡路径采样（TPS）在反应路径搜索中的优劣。

文末回归分子模拟本质——用有限算力在高维相空间逼近真实物理规律，适合作为MD入门认知框架搭建参考。

从赌博到真正的动力学

很多人第一次接触蒙特卡洛（MC）都会被它的采样思路惊艳，但一定要先记牢它的先天缺陷：MC的世界里没有真实时间。你看到的步数只是数学迭代，对应不到真实物理里的飞秒、皮秒，所以它永远算不了蛋白质折叠、扩散系数、反应速率这些和动力学相关的性质。

那为什么还要学MC？它在相空间探索上的效率是MD比不了的：不需要算复杂受力，只靠势能函数就能通过大步跳跃跨过能垒，在构象搜索、全局最优解寻找、特定系综转换场景里，至今仍是首选工具。

从能量拆分的角度看，系统总能量分为和温度绑定的动能、和原子位置绑定的势能，二者相互独立。MC只处理势能项，动能部分可以用理想气体模型轻松求解，这也是它轻量高效的根本原因。如果想看分子真实的演化过程，就必须回到牛顿力学框架。

确定性的舞蹈

MD的核心逻辑非常简单：给定初始位置和速度，通过牛顿第二定律做时间积分，预测体系的运动轨迹。因为多体体系不存在解析解，我们只能把时间切成极小的步长，靠积分器带着体系一步步走。

一个合格的积分器必须满足两个硬指标：严格守恒能量和动量，且具备时间可逆性——微观物理定律本身是可逆的，你把速度反向，积分器就该带你回到过去。至于速度快、步长大、易实现都是附加要求。

你永远不可能得到和真实物理完全一致的长时轨迹，多体系统的混沌特性会带来李雅普诺夫不稳定性，初始条件哪怕只有10^-6量级的速度偏差，也会在模拟中被指数级放大。但统计力学已经帮我们兜了底：个别轨迹的偏离不影响群体平均值的正确性，这正是MD成立的前提。

Verlet家族

MD界事实上的标准积分器全部来自Verlet家族。最早的原始Verlet算法靠泰勒展开消去了速度项，实现极其简单，天生满足时间可逆，但速度是隐含的，想算动能、温度都得额外做差分，精度不够高。

现在几乎所有主流MD软件（GROMACS、LAMMPS等）默认用的都是Velocity-Verret，它显式地分两步更新位置和速度：先按当前速度走半步位置，根据新位置算加速度，再补上剩下的半步速度。这种“左右横跳”的策略既保住了精度，还能让你实时监控体系温度，完全适配生产级模拟的需求。

在计算机上的实验室

第一篇里提过的系综，到这里才真正落地到计算里。纯牛顿力学对应的是孤立系统NVE系综，但真实的生物实验大多在恒温恒压环境下做，我们需要给模拟加上“控温控压”的算法，在计算机里搭出和试管一样的实验环境。

先说恒温器，这是最容易踩坑的部分：如果你要算扩散系数这类动力学性质，选错恒温器等于直接废掉数据。Andersen恒温器的逻辑是随机给原子赋予玻尔兹曼分布的新速度，相当于原子经历了瞬时剧烈碰撞，直接抹除了之前的运动记忆，彻底破坏了速度自相关函数，算出来的动力学参数完全没有物理意义。Berendsen恒温器靠缩放速度平滑调温，本质是一种全局约束，只保证总动能符合目标，但体系自身的涨落不符合正则系综的分布特征。

• Berendsen Thermostat：虽然它看起来温和一点，通过缩放速度（Velocity rescaling）来平滑调整，但它本质上是一种全局约束。它像是一个强权政府，只在乎“平均财富”（总动能），却无视个体的波动。这种人为的干预会导致系统的涨落不符合 Canonical Ensemble

现在学界公认靠谱的是Nose-Hoover链：它不直接改原子速度，而是给体系加一串嵌套的虚拟动力学变量，第一个变量控制体系温度，第二个控制第一个变量的变化，一般套5层就能让体系的动能涨落完全符合正态分布，既不破坏微观动力学相关性，又能精准复现NVT系综。业内甚至有共识：1992年Nose-Hoover链问世之前的MD动力学速率相关文章，基本都可以跳过。

Nose-Hoover Chain. 可以看到虚线几乎看不出了，和正态重合。上面的甜甜圈也变成了俯瞰正态分布的形状。

恒压器的逻辑更直观：让模拟盒子变成“会呼吸”的，体系压力大就放大体积，压力小就缩小体积，经典的Parrinello-Rahman算法还支持盒子变形，用来研究固态相变非常合适，本质是改坐标系尺度平衡内外压力，复现NPT系综。

算得快才是硬道理

MD的生产级模拟动辄上百纳秒，哪怕用GPU加速也要跑几天，效率永远是核心痛点。最耗时的部分是原子间相互作用的计算：如果暴力算两两之间的力，10万原子体系单步就要做100亿次运算，完全不可行。

第一步优化是截断半径：范德华这类短程力衰减极快，超过一定距离的力直接忽略，先把计算量降下来。但计算机每步还是要扫所有原子对的距离，依然慢。

Neighbor list.用一个 linked list 来表征。

第二步是Verlet列表：给截断半径加一层缓冲带，不用每步都更新邻居列表，只有当原子走出缓冲带才重新算，中间十几步只需要算列表里已经存好的原子对。

体系再大就用元胞链表：把整个模拟盒子切成和截断半径差不多大的小格子，每个原子找邻居只需要看自己所在格子和周围26个相邻格子，直接把计算复杂度从O(N²)降到O(N)，是大规模MD的必备优化。

跑正式模拟前还必须做能量最小化：你手动搭的体系难免有原子重叠，在伦纳德-琼斯势能曲线里，距离过近意味着巨大的排斥力，直接跑MD会直接让体系“爆炸”。一般先用最速下降法快速把重叠的原子挤开，再用共轭梯度法精细调整到稳定的亚稳态，相当于让体系先站稳再起跑。

圣杯——自由能计算

MD的终极目标几乎都是算自由能。统计力学里所有物理量都来自配分函数，但绝对自由能在分子模拟里是算不出来的：重要性采样的设计让我们只采样了相空间里能量低的极小一部分区域，根本不知道没采样的“荒原”有多大，配分函数对整个相空间的积分永远得不到。

所以我们只算自由能差ΔG。所有可算的量本质都是两个积分的比值，类似量子场论里的重整化：我们不关心那个发散的绝对自由能，只关心两个状态的相对差异，把不可测的配分函数归一化消掉就行。不管是药物结合自由能，还是沿某个维度的平均力势（PMF），全是这个逻辑。

算自由能差有三个主流方法：自由能微扰（FEP）适合两个状态非常接近的场景，比如把氨基酸的甲基突变成羟基，步子大了公式直接失效；伞形采样是人为加偏置势把体系“摁”在能垒顶采样，采完再用WHAM或者MBAR把偏置去掉，专门解决“翻不过能垒”的问题；现在最流行的是元动力学（Metadynamics），每采样一个位置就加一个高斯势“填坑”，逼着体系探索没去过的地方，最后把加的所有势倒过来就是自由能地貌。

所有这些方法都有一个前提：你得选对集体变量（CV）——也就是你认为能表征反应过程的那个坐标，比如两个蛋白的距离、键的扭转角。如果选错了CV，比如蛋白折叠的关键是拓扑结构你却只看了距离，算出来的结果完全无效。

何以知道CV？

选对CV之后还要找反应路径。物理学家给了过渡态的严格定义：分选系数（Committor）等于0.5的等能面，就是真正的分水岭——从这个点出发，一半概率到产物，一半概率到反应物。但这个方法太费算力，没人会对每个点都跑几百条轨迹。

路径搜索的核心是找最小自由能路径（MFEP），要在“下山快、路径连续、适配高维”三个互斥的要求里找平衡。最早的弦方法把路径当绳子在势能面演化，但必须预先指定CV，选错了就找不到正确路径。Path-MetaD不用预先降维，但高维采样的时候路径点很容易乱序，“绳子”直接断开。

竺老师团队的TAPS算法解决了这个问题：把路径点重排转化为旅行商问题，用找最短路径的逻辑保证点序永远不乱，既保留了Path-MetaD的高效，又不会断路径。如果不想预设CV，还可以用过渡路径采样（TPS），直接在相空间做扰动抽轨迹，不需要任何先验假设，但遇到多能垒体系很容易卡在一个能垒附近反复横跳，没法探索完整路径。

Recap

整个分子模拟的脉络其实非常清晰：从热力学到统计力学，再到系综理论，发现配分函数算不出来，于是转去做计算机采样，分出MC和MD两条路线；MD为了适配真实实验环境，不断优化积分器、恒温恒压算法、计算效率，又发展出各种增强采样方法找自由能和反应路径。所有设计的出发点都是同一个：用有限的算力，在不可穷尽的高维相空间里，还原最接近真实的物理规律。

如果你刚入门MD，这篇内容可以帮你搭起完整的认知框架；如果是已经在一线做研究的同行，也可以回头补一下各个方法的设计初衷——很多时候我们踩的坑，本质上都是忘了方法本身的假设前提。

原文链接：https://mp.weixin.qq.com/s/yPTyg7EMQ1DFR3xspgdXsA

这是一篇面向计算生物与分子模拟入门者的分子动力学（MD）方法导读。文章承接前文蒙特卡洛（MC）方法，指出MC因无真实时间维度无法计算动力学性质，由此引出MD的核心——基于牛顿第二定律对时间做离散化积分来复现分子运动轨迹，并阐明多体系统混沌特性导致长时轨迹必然偏离，但统计力学的系综平均保证计算结果的有效性。

文章系统梳理了MD全流程要点：基础层介绍Verlet家族积分器演化，点明Velocity-Verlet因可显式监控温度成为主流软件默认算法；工程层讲解Nose-Hoover链恒温器和Parrinello-Rahman恒压器以复现NVT/NPT系综，对比指出早期Andersen与Berendsen恒温器的缺陷及1992年Nose-Hoover链问世是MD动力学速率计算可信的分水岭，同时归纳截断半径、邻居列表、能量最小化等算力优化手段。应用层说明绝对自由能不可算而只算自由能差，解析FEP、伞形采样、元动力学的适用场景，强调集体变量（CV）选择的关键性，并对比弦方法、Path-MetaD、TAPS及过渡路径采样（TPS）在反应路径搜索中的优劣。

文末回归分子模拟本质——用有限算力在高维相空间逼近真实物理规律，适合作为MD入门认知框架搭建参考。

从赌博到真正的动力学

很多人第一次接触蒙特卡洛（MC）都会被它的采样思路惊艳，但一定要先记牢它的先天缺陷：MC的世界里没有真实时间。你看到的步数只是数学迭代，对应不到真实物理里的飞秒、皮秒，所以它永远算不了蛋白质折叠、扩散系数、反应速率这些和动力学相关的性质。

那为什么还要学MC？它在相空间探索上的效率是MD比不了的：不需要算复杂受力，只靠势能函数就能通过大步跳跃跨过能垒，在构象搜索、全局最优解寻找、特定系综转换场景里，至今仍是首选工具。

从能量拆分的角度看，系统总能量分为和温度绑定的动能、和原子位置绑定的势能，二者相互独立。MC只处理势能项，动能部分可以用理想气体模型轻松求解，这也是它轻量高效的根本原因。如果想看分子真实的演化过程，就必须回到牛顿力学框架。

确定性的舞蹈

MD的核心逻辑非常简单：给定初始位置和速度，通过牛顿第二定律做时间积分，预测体系的运动轨迹。因为多体体系不存在解析解，我们只能把时间切成极小的步长，靠积分器带着体系一步步走。

一个合格的积分器必须满足两个硬指标：严格守恒能量和动量，且具备时间可逆性——微观物理定律本身是可逆的，你把速度反向，积分器就该带你回到过去。至于速度快、步长大、易实现都是附加要求。

你永远不可能得到和真实物理完全一致的长时轨迹，多体系统的混沌特性会带来李雅普诺夫不稳定性，初始条件哪怕只有10^-6量级的速度偏差，也会在模拟中被指数级放大。但统计力学已经帮我们兜了底：个别轨迹的偏离不影响群体平均值的正确性，这正是MD成立的前提。

Verlet家族

MD界事实上的标准积分器全部来自Verlet家族。最早的原始Verlet算法靠泰勒展开消去了速度项，实现极其简单，天生满足时间可逆，但速度是隐含的，想算动能、温度都得额外做差分，精度不够高。

现在几乎所有主流MD软件（GROMACS、LAMMPS等）默认用的都是Velocity-Verret，它显式地分两步更新位置和速度：先按当前速度走半步位置，根据新位置算加速度，再补上剩下的半步速度。这种“左右横跳”的策略既保住了精度，还能让你实时监控体系温度，完全适配生产级模拟的需求。

在计算机上的实验室

第一篇里提过的系综，到这里才真正落地到计算里。纯牛顿力学对应的是孤立系统NVE系综，但真实的生物实验大多在恒温恒压环境下做，我们需要给模拟加上“控温控压”的算法，在计算机里搭出和试管一样的实验环境。

先说恒温器，这是最容易踩坑的部分：如果你要算扩散系数这类动力学性质，选错恒温器等于直接废掉数据。Andersen恒温器的逻辑是随机给原子赋予玻尔兹曼分布的新速度，相当于原子经历了瞬时剧烈碰撞，直接抹除了之前的运动记忆，彻底破坏了速度自相关函数，算出来的动力学参数完全没有物理意义。Berendsen恒温器靠缩放速度平滑调温，本质是一种全局约束，只保证总动能符合目标，但体系自身的涨落不符合正则系综的分布特征。

• Berendsen Thermostat：虽然它看起来温和一点，通过缩放速度（Velocity rescaling）来平滑调整，但它本质上是一种全局约束。它像是一个强权政府，只在乎“平均财富”（总动能），却无视个体的波动。这种人为的干预会导致系统的涨落不符合 Canonical Ensemble

现在学界公认靠谱的是Nose-Hoover链：它不直接改原子速度，而是给体系加一串嵌套的虚拟动力学变量，第一个变量控制体系温度，第二个控制第一个变量的变化，一般套5层就能让体系的动能涨落完全符合正态分布，既不破坏微观动力学相关性，又能精准复现NVT系综。业内甚至有共识：1992年Nose-Hoover链问世之前的MD动力学速率相关文章，基本都可以跳过。

Nose-Hoover Chain. 可以看到虚线几乎看不出了，和正态重合。上面的甜甜圈也变成了俯瞰正态分布的形状。

恒压器的逻辑更直观：让模拟盒子变成“会呼吸”的，体系压力大就放大体积，压力小就缩小体积，经典的Parrinello-Rahman算法还支持盒子变形，用来研究固态相变非常合适，本质是改坐标系尺度平衡内外压力，复现NPT系综。

算得快才是硬道理

MD的生产级模拟动辄上百纳秒，哪怕用GPU加速也要跑几天，效率永远是核心痛点。最耗时的部分是原子间相互作用的计算：如果暴力算两两之间的力，10万原子体系单步就要做100亿次运算，完全不可行。

第一步优化是截断半径：范德华这类短程力衰减极快，超过一定距离的力直接忽略，先把计算量降下来。但计算机每步还是要扫所有原子对的距离，依然慢。

Neighbor list.用一个 linked list 来表征。

第二步是Verlet列表：给截断半径加一层缓冲带，不用每步都更新邻居列表，只有当原子走出缓冲带才重新算，中间十几步只需要算列表里已经存好的原子对。

体系再大就用元胞链表：把整个模拟盒子切成和截断半径差不多大的小格子，每个原子找邻居只需要看自己所在格子和周围26个相邻格子，直接把计算复杂度从O(N²)降到O(N)，是大规模MD的必备优化。

跑正式模拟前还必须做能量最小化：你手动搭的体系难免有原子重叠，在伦纳德-琼斯势能曲线里，距离过近意味着巨大的排斥力，直接跑MD会直接让体系“爆炸”。一般先用最速下降法快速把重叠的原子挤开，再用共轭梯度法精细调整到稳定的亚稳态，相当于让体系先站稳再起跑。

圣杯——自由能计算

MD的终极目标几乎都是算自由能。统计力学里所有物理量都来自配分函数，但绝对自由能在分子模拟里是算不出来的：重要性采样的设计让我们只采样了相空间里能量低的极小一部分区域，根本不知道没采样的“荒原”有多大，配分函数对整个相空间的积分永远得不到。

所以我们只算自由能差ΔG。所有可算的量本质都是两个积分的比值，类似量子场论里的重整化：我们不关心那个发散的绝对自由能，只关心两个状态的相对差异，把不可测的配分函数归一化消掉就行。不管是药物结合自由能，还是沿某个维度的平均力势（PMF），全是这个逻辑。

算自由能差有三个主流方法：自由能微扰（FEP）适合两个状态非常接近的场景，比如把氨基酸的甲基突变成羟基，步子大了公式直接失效；伞形采样是人为加偏置势把体系“摁”在能垒顶采样，采完再用WHAM或者MBAR把偏置去掉，专门解决“翻不过能垒”的问题；现在最流行的是元动力学（Metadynamics），每采样一个位置就加一个高斯势“填坑”，逼着体系探索没去过的地方，最后把加的所有势倒过来就是自由能地貌。

所有这些方法都有一个前提：你得选对集体变量（CV）——也就是你认为能表征反应过程的那个坐标，比如两个蛋白的距离、键的扭转角。如果选错了CV，比如蛋白折叠的关键是拓扑结构你却只看了距离，算出来的结果完全无效。

何以知道CV？

选对CV之后还要找反应路径。物理学家给了过渡态的严格定义：分选系数（Committor）等于0.5的等能面，就是真正的分水岭——从这个点出发，一半概率到产物，一半概率到反应物。但这个方法太费算力，没人会对每个点都跑几百条轨迹。

路径搜索的核心是找最小自由能路径（MFEP），要在“下山快、路径连续、适配高维”三个互斥的要求里找平衡。最早的弦方法把路径当绳子在势能面演化，但必须预先指定CV，选错了就找不到正确路径。Path-MetaD不用预先降维，但高维采样的时候路径点很容易乱序，“绳子”直接断开。

竺老师团队的TAPS算法解决了这个问题：把路径点重排转化为旅行商问题，用找最短路径的逻辑保证点序永远不乱，既保留了Path-MetaD的高效，又不会断路径。如果不想预设CV，还可以用过渡路径采样（TPS），直接在相空间做扰动抽轨迹，不需要任何先验假设，但遇到多能垒体系很容易卡在一个能垒附近反复横跳，没法探索完整路径。

Recap

整个分子模拟的脉络其实非常清晰：从热力学到统计力学，再到系综理论，发现配分函数算不出来，于是转去做计算机采样，分出MC和MD两条路线；MD为了适配真实实验环境，不断优化积分器、恒温恒压算法、计算效率，又发展出各种增强采样方法找自由能和反应路径。所有设计的出发点都是同一个：用有限的算力，在不可穷尽的高维相空间里，还原最接近真实的物理规律。

如果你刚入门MD，这篇内容可以帮你搭起完整的认知框架；如果是已经在一线做研究的同行，也可以回头补一下各个方法的设计初衷——很多时候我们踩的坑，本质上都是忘了方法本身的假设前提。

原文链接：https://mp.weixin.qq.com/s/yPTyg7EMQ1DFR3xspgdXsA