多体系统动力学方程数值求解技巧

由多体力学方程的结构形式可知，多体力学仿真数值求解的核心通常是对常微分方程初值问题的处理。其公式如下：

求解式（1-3）的基本途径有以下3种：

（1）化导数为差商的方法，即用差商来近似代替导数，从而得到数值解序列。具有代表性的是各种欧拉方法。

（2）数值积分法，将方程化成积分形式，利用梯形、龙贝格、高斯等数值积分方法得到解序列。(www.xing528.com)

（3）利用泰勒公式的近似求解。典型的方法是各阶龙格-库塔公式。另外，为了充分利用有用信息，进一步提高计算结果的精度，还提出了线形多步法来代替单步法的思想，典型的如亚当姆斯（Adams）法和哈明（Hamming）法。

由于在方程求解时经常会遇到系统的特征值在数值上相差若干个数量级的情况，描述这种系统的微分方程，称为刚性（Stiff）方程。对这种方程的处理必须采用特殊的方法，现在常用的方法有隐式或半隐式Kunge-Kutta法、自动变阶变步长的Gear法、隐式或显式Adams法等，而且对于线性病态系统，还可以用增广矩阵法和蛙跳算法等。

多柔体系统在数值计算时，慢变大幅变量与快变微幅变量的耦合会导致方程严重的病态，这个问题已成为多柔体系统动力学发展的一个“瓶颈”，引起了学者们的普遍关注。Gear方法被认为是求解刚性微分方程的很有效的办法，但用到多柔体系统动力学方程上有很大的不便，即Gear方法需计算方程右端项的Jacobian矩阵，这对复杂多柔体系统动力学方程而言，几乎是难以做到的。