网格弹簧质点系统模拟（Spring-Mass System by Euler Integration）

　　弹簧质点模型是利用牛顿运动定律来模拟物体变形的方法。如下图所示，该模型是一个由m×n个虚拟质点组成的网格，质点之间用无质量的、自然长度不为零的弹簧连接。其连接关系有以下三种：

　　1.连接质点[i, j]与[i+1,j]，[i, j]与[i,j+1]的弹簧，称为“结构弹簧”；

　　2.连接质点[i, j]与[i+1,j+1]，[i+1,j]与[i,j+1]的弹簧，称为“剪切弹簧”；

　　3.连接质点[i, j]与[i+2,j]，[i, j]与[i, j+2]的弹簧，称为“弯曲弹簧”。

　　这三种弹簧分别用于与结构力（拉力或压力）、剪力和弯矩相关的计算。

　　弹簧质点运动时受到内力和外力和影响，内力包括弹簧的弹性力和阻尼力，外力包括重力以及空气阻力等。

　　弹簧的弹性力遵从Hooke定律，质点i和质点j之间的弹性力可以表示为：

其中：，，k_s表示弹性系数，x_i表示质点i的位置，l₀表示弹簧的自然长度。

　　弹簧的阻尼力与质点相对速度在弹簧方向上的分量成正比，可以表示为：

其中：v_ij = v_j – v_i，k_d表示阻尼系数，v_i表示质点i的速度。

　　牛顿运动定律是弹簧质点模型的关键，将该二阶微分方程改写成两个一阶微分方程：

（1）

　　利用显式欧拉积分求解上述微分方程可以得到：

（2）

　　整个过程的伪代码如下：

　　显式欧拉积分表示简单，求解快速，但其对积分步长有限制，一旦超过某个临界值T₀时，积分结果就会发散。而临界值T₀主要与弹簧的弹性系数有关，弹簧的弹性系数越大，临界值T₀就越小，对应的步长也越小，这就意味着两个可视帧之间需要更多次的迭代才能使结果稳定。

　　一种提高积分精度的方法是采用四阶Runge-Kutta积分，其精度能达到O(h⁴)：

　　整个过程的伪代码如下：

　　隐式欧拉积分比显式欧拉积分具有更高的稳定性，弹簧质点模型的隐式欧拉积分表达式如下：

（3）

　　将f^t+h一阶泰勒展开，代入式（3）第一个微分方程得到：

　　整理得：

　　两边乘以M：

　　上式中M是以3×3对角矩阵为子矩阵的对角矩阵，两个Jacobian矩阵和都是子矩阵为3×3的对称稀疏矩阵，只有当质点i和质点j之间有弹簧连接时对应的子矩阵才非零（对角线也非零）。隐式欧拉积分需要求解线性方程组来得到下一时刻状态，但是稳定性好，能够取较大的步长进行计算。

参考文献：

[1] David Baraff and Andrew Witkin. 1998. Large steps in cloth simulation. In Proceedings of the 25th annual conference on Computer graphics and interactive techniques (SIGGRAPH ‘98). ACM, New York, NY, USA, 43-54.

[2] Matthias Müller, Jos Stam, Doug James, and Nils Thürey. 2008. Real time physics: class notes. In ACM SIGGRAPH 2008 classes (SIGGRAPH ‘08). ACM, New York, NY, USA, , Article 88 , 90 pages.

附录

标量函数对向量的导数：

，其中：x∈Rⁿ。