Khalili Mey, Liwo Adam, Rakowski Franciszek, Grochowski Paweł, Scheraga Harold A
Baker Laboratory of Chemistry and Chemical Biology, Cornell University, Ithaca, New York 14853-1301, USA.
J Phys Chem B. 2005 Jul 21;109(28):13785-97. doi: 10.1021/jp058008o.
The Lagrange formalism was implemented to derive the equations of motion for the physics-based united-residue (UNRES) force field developed in our laboratory. The C(alpha)...C(alpha) and C(alpha)...SC (SC denoting a side-chain center) virtual-bond vectors were chosen as variables. The velocity Verlet algorithm was adopted to integrate the equations of motion. Tests on the unblocked Ala(10) polypeptide showed that the algorithm is stable in short periods of time up to the time step of 1.467 fs; however, even with the shorter time step of 0.489 fs, some drift of the total energy occurs because of momentary jumps of the acceleration. These jumps are caused by numerical instability of the forces arising from the U(rot) component of UNRES that describes the energetics of side-chain-rotameric states. Test runs on the Gly(10) sequence (in which U(rot) is not present) and on the Ala(10) sequence with U(rot) replaced by a simple numerically stable harmonic potential confirmed this observation; oscillations of the total energy were observed only up to the time step of 7.335 fs, and some drift in the total energy or instability of the trajectories started to appear in long-time (2 ns and longer) trajectories only for the time step of 9.78 fs. These results demonstrate that the present U(rot) components (which are statistical potentials derived from the Protein Data Bank) must be replaced with more numerically stable functions; this work is under way in our laboratory. For the purpose of our present work, a nonsymplectic variable-time-step algorithm was introduced to reduce the energy drift for regular polypeptide sequences. The algorithm scales down the time step at a given point of a trajectory if the maximum change of acceleration exceeds a selected cutoff value. With this algorithm, the total energy is reasonably conserved up to a time step of 2.445 fs, as tested on the unblocked Ala(10) polypeptide. We also tried a symplectic multiple-time-step reversible RESPA algorithm and achieved satisfactory energy conservation for time steps up to 7.335 fs. However, at present, it appears that the reversible RESPA algorithm is several times more expensive than the variable-time-step algorithm because of the necessity to perform additional matrix multiplications. We also observed that, because Ala(10) folds and unfolds within picoseconds in the microcanonical mode, this suggests that the effective (event-based) time unit in UNRES dynamics is much larger than that of all-atom dynamics because of averaging over the fast-moving degrees of freedom in deriving the UNRES potential.
采用拉格朗日形式主义来推导我们实验室开发的基于物理的联合残基（UNRES）力场的运动方程。选择C(α)…C(α)和C(α)…SC（SC表示侧链中心）虚拟键向量作为变量。采用速度Verlet算法对方程进行积分。对无阻碍的Ala(10)多肽的测试表明，该算法在短时间内直至时间步长为1.467 fs时是稳定的；然而，即使时间步长缩短至0.489 fs，由于加速度的瞬间跳跃，总能量仍会出现一些漂移。这些跳跃是由UNRES的U(rot)分量产生的力的数值不稳定性引起的，U(rot)分量描述了侧链旋转异构体状态的能量学。对Gly(10)序列（其中不存在U(rot)）以及用简单数值稳定的谐振势代替U(rot)的Ala(10)序列进行的测试运行证实了这一观察结果；仅在时间步长为7.335 fs时观察到总能量的振荡，并且仅在时间步长为9.78 fs时，在长时间（2 ns及更长）轨迹中才开始出现总能量的一些漂移或轨迹的不稳定性。这些结果表明，当前的U(rot)分量（它们是从蛋白质数据库导出的统计势）必须用数值更稳定的函数来代替；我们实验室正在进行这项工作。为了我们目前的工作目的，引入了一种非辛变时间步长算法来减少规则多肽序列的能量漂移。如果加速度的最大变化超过选定的截止值，该算法会在轨迹的给定时刻缩小时间步长。通过该算法，对无阻碍的Ala(10)多肽进行测试时，总能量在时间步长为2.445 fs时得到合理守恒。我们还尝试了辛多时间步长可逆RESPA算法，并在时间步长达到7.335 fs时实现了令人满意的能量守恒。然而，目前看来，由于需要进行额外的矩阵乘法，可逆RESPA算法比变时间步长算法贵几倍。我们还观察到，由于Ala(10)在微正则模式下在皮秒内折叠和展开，这表明由于在推导UNRES势时对快速移动的自由度进行了平均，UNRES动力学中的有效（基于事件）时间单位比全原子动力学中的大得多。
