在视频游戏中模拟物理现象是非常普遍的，因为大多数游戏的灵感都来源于现实世界中的事物。刚体动力学(即固体、不可弯曲物体的运动与相互作用)是迄今为止游戏中模拟得最多的一种物理效果。
在本系列教程中，我们将探索刚体模拟，从本文中的简单刚体运动开始，然后在后续部分中通过碰撞和约束来探讨物体间的相互作用。
这是我们关于游戏物理的三部分系列的第一部分。本系列的其余部分包括：
- Part 2：固体物体的碰撞检测
- Part 3：受约束的刚体模拟
物理模拟是计算机科学中的一个领域，旨在利用计算机重现物理现象。通常，这些模拟通过将数值方法应用于现有理论，以获得尽可能接近我们在现实世界中观察到的结果。作为高级游戏开发者，这使我们能够预测并仔细分析某物在实际构建之前的行为方式，而这几乎总是更简单且成本更低的做法。
物理模拟的应用范围非常广泛。最早的计算机已经被用来进行物理模拟——例如，在军事中预测弹道运动。它也是土木和汽车工程中的重要工具，能够揭示某些结构在地震或车祸等事件中的行为。不仅如此，我们还可以模拟天体物理学、相对论以及许多其他在自然界奇迹中能够观察到的复杂现象。
在视频游戏中模拟物理现象非常普遍，因为大多数游戏的灵感都来源于现实世界中的事物。许多游戏完全依赖于物理模拟来增加趣味性。因此，物理模拟游戏需要一个稳定的模拟系统，不会崩溃或变慢，而这通常并不容易实现。
在任何游戏中，只有某些特定的物理效果是关注的焦点。刚体动力学（即固体、不可弯曲物体的运动与相互作用）是迄今为止游戏中模拟得最多的一种物理效果。这是因为我们在现实生活中与之交互的大多数物体都是相对刚性的，而且模拟刚体相对简单（尽管我们会看到，这也并不意味着它很轻而易举）。不过，一些其他游戏需要模拟更复杂的实体，如可变形体、流体、磁性物体等。
在本系列游戏物理教程中，我们将探索刚体物理模拟，从本文中的简单刚体运动开始，然后在后续部分中通过碰撞和约束来探讨物体间的相互作用。我们将介绍并解释现代游戏物理引擎（如Box2D、Bullet Physics和Chipmunk Physics）中最常用的方程。
刚体动力学
在游戏物理中，我们希望让屏幕上的物体动起来，并赋予它们逼真的物理行为。这是通过基于物理的程序化动画实现的，即通过数值计算并应用物理理论定律生成的动画。
动画是通过连续显示一系列图像来生成的，物体在每张图像之间略微移动。当图像快速连续显示时，效果就是物体看起来平滑且连续的运动。因此，为了在物理模拟中让物体动起来，我们需要根据物理定律每秒多次更新物体的物理状态（例如位置和方向），并在每次更新后重新绘制屏幕。
物理引擎是执行物理模拟的软件组件。它接收要模拟的物体的规格说明以及一些配置参数，然后可以逐步推进模拟。每一步都将模拟向前推进几分之一秒，之后可以将结果显示在屏幕上。需要注意的是，物理引擎只执行数值模拟。如何处理结果可能取决于游戏的需求。并非每一步的结果都需要绘制到屏幕上。
刚体的运动可以使用牛顿力学来建模，牛顿力学基于艾萨克·牛顿著名的三大运动定律：
- 惯性定律：如果物体上没有施加力，其速度（运动的速度和方向）不会改变。
- 力、质量和加速度：作用在物体上的力等于物体的质量乘以其加速度（速度的变化率）。公式为 F = ma。
- 作用与反作用定律：“对于每一个作用力，都有一个大小相等、方向相反的反作用力。”换句话说，当一个物体对另一个物体施加力时，第二个物体会对第一个物体施加一个大小相等、方向相反的力。
基于这三条定律，我们可以制作一个物理引擎，能够重现我们熟悉的动态行为，从而为玩家创造沉浸式体验。
下图展示了物理引擎的一般流程的高级概述：
向量
为了理解物理模拟的工作原理，掌握向量及其运算的基础知识至关重要。如果你已经熟悉向量数学，可以继续阅读。但如果你不熟悉，或者想复习一下，请花点时间阅读游戏中的物理学：向量。
粒子模拟
理解刚体模拟的一个很好的起点是从粒子开始。模拟粒子比模拟刚体更简单，我们可以使用相同的原理来模拟后者，只需为粒子添加体积和形状即可。此外，你可以在网上轻松找到关于粒子物理引擎的优质教程，涵盖各种场景和真实感/复杂度的层次。
粒子只是空间中的一个点，具有位置向量、速度向量和质量。根据牛顿第一定律，只有当力作用在粒子上时，其速度才会改变。当其速度向量的长度不为零时，其位置会随时间变化。
要模拟粒子系统，我们首先需要创建一个具有初始状态的粒子数组。每个粒子必须具有固定的质量、空间中的初始位置和初始速度。然后，我们需要启动模拟的主循环，对于每个粒子，计算当前作用在其上的力，根据力产生的加速度更新其速度，然后根据刚刚计算的速度更新其位置。
力的来源可能因模拟类型而异。它可以是重力、风力、磁力等，或者是这些力的组合。它可以是全局力，例如恒定的重力，也可以是粒子之间的力，例如吸引力或排斥力。
为了使模拟以逼真的速度运行，我们“模拟”的时间步长应与自上次模拟步骤以来经过的实际时间相同。然而，可以缩放此时间步长以使模拟运行得更快，或缩小它以使其以慢动作运行。
假设我们有一个粒子，质量为 m，在时间 ti 的位置为 p(ti)，速度为 v(ti)。此时作用在该粒子上的力为 f(ti)。该粒子在未来时间 ti + 1 的位置 p(ti + 1) 和速度 v(ti + 1) 可以通过以下公式计算：
从技术上讲，我们在这里所做的是使用半隐式欧拉法对粒子的运动常微分方程进行数值积分。这是大多数游戏物理引擎使用的方法，因为它简单且在小时间步长（dt）下具有可接受的精度。如果你想更深入地了解这种方法的数值计算原理，可以参考由 Andy Witkin 博士和 David Baraff 博士编写的《基于物理的建模：原理与实践》课程中的《微分方程基础》讲义，这是一篇非常不错的文章。
以下是一个用 C 语言编写的粒子模拟示例：
#define NUM_PARTICLES 1
// Two dimensional vector.
typedef struct {
float x;
float y;
} Vector2;
// Two dimensional particle.
typedef struct {
Vector2 position;
Vector2 velocity;
float mass;
} Particle;
// Global array of particles.
Particle particles[NUM_PARTICLES];
// Prints all particles' position to the output. We could instead draw them on screen
// in a more interesting application.
void PrintParticles() {
for (int i = 0; i < NUM_PARTICLES; ++i) {
Particle *particle = &particles[i];
printf("particle[%i] (%.2f, %.2f)\n", i, particle->position.x, particle->position.y);
}
}
// Initializes all particles with random positions, zero velocities and 1kg mass.
void InitializeParticles() {
for (int i = 0; i < NUM_PARTICLES; ++i) {
particles[i].position = (Vector2){arc4random_uniform(50), arc4random_uniform(50)};
particles[i].velocity = (Vector2){0, 0};
particles[i].mass = 1;
}
}
// Just applies Earth's gravity force (mass times gravity acceleration 9.81 m/s^2) to each particle.
Vector2 ComputeForce(Particle *particle) {
return (Vector2){0, particle->mass * -9.81};
}
void RunSimulation() {
float totalSimulationTime = 10; // The simulation will run for 10 seconds.
float currentTime = 0; // This accumulates the time that has passed.
float dt = 1; // Each step will take one second.
InitializeParticles();
PrintParticles();
while (currentTime < totalSimulationTime) {
// We're sleeping here to keep things simple. In real applications you'd use some
// timing API to get the current time in milliseconds and compute dt in the beginning
// of every iteration like this:
// currentTime = GetTime()
// dt = currentTime - previousTime
// previousTime = currentTime
sleep(dt);
for (int i = 0; i < NUM_PARTICLES; ++i) {
Particle *particle = &particles[i];
Vector2 force = ComputeForce(particle);
Vector2 acceleration = (Vector2){force.x / particle->mass, force.y / particle->mass};
particle->velocity.x += acceleration.x * dt;
particle->velocity.y += acceleration.y * dt;
particle->position.x += particle->velocity.x * dt;
particle->position.y += particle->velocity.y * dt;
}
PrintParticles();
currentTime += dt;
}
}
如果你调用 RunSimulation 函数（针对单个粒子），它将输出类似以下内容：
particle[0] (-8.00, 57.00)
particle[0] (-8.00, 47.19)
particle[0] (-8.00, 27.57)
particle[0] (-8.00, -1.86)
particle[0] (-8.00, -41.10)
particle[0] (-8.00, -90.15)
particle[0] (-8.00, -149.01)
particle[0] (-8.00, -217.68)
particle[0] (-8.00, -296.16)
particle[0] (-8.00, -384.45)
particle[0] (-8.00, -482.55)
正如你所看到的，粒子从 (-8, 57) 的位置开始，然后它的 y 坐标开始越来越快地下降，因为它在重力作用下向下加速。
以下动画直观地展示了单个粒子模拟的三个步骤序列：
最初，在 t = 0 时，粒子位于 p0。经过一步后，它沿着其速度向量 v0 指向的方向移动。在下一步中，力 f1 作用在粒子上，速度向量开始变化，仿佛被拉向力的方向。在接下来的两步中，力的方向发生变化，但继续将粒子向上拉动。
刚体物理模拟
刚体是一种不能变形的固体。这种固体在现实世界中并不存在（即使是最硬的材料在受到力作用时也会发生至少非常微小的形变）但刚体是游戏开发者常用的物理模型，它简化了对固体动力学的研究，尤其是在我们可以忽略形变的情况下。
刚体类似于粒子的扩展，因为它也具有质量、位置和速度。此外，它还具有体积和形状，因此可以旋转。这比听起来要复杂得多，尤其是在三维空间中。
刚体自然围绕其质心旋转，而刚体的位置被认为是其质心的位置。我们定义刚体的初始状态为质心位于原点且旋转角度为零。在任何时刻 t，刚体的位置和旋转都是初始状态的偏移量。
质心是物体质量分布的中点。如果你想象一个质量为 M 的刚体由 N 个微小粒子组成，每个粒子的质量为 mi，位置为 ri（微小粒子位于刚体内部），则质心可以通过以下公式计算：
这个公式表明，质心是粒子位置的加权平均值，权重为它们的质量。如果物体的密度在整个体积内是均匀的，那么质心与物体形状的几何中心（也称为形心）相同。游戏物理引擎通常只支持均匀密度，因此可以使用几何中心作为质心。
然而，刚体并不是由有限数量的离散粒子组成的，它们是连续的。因此，我们应该使用积分而不是有限求和来计算质心，如下所示：
其中，r 是每个点的位置向量，𝜌 (rho) 是一个函数，表示物体内每个点的密度。本质上，这个积分与有限求和的作用相同，但它是在连续体积 V 中进行的。
由于刚体可以旋转，我们必须引入其角属性，这些属性类似于粒子的线性属性。在二维空间中，刚体只能绕垂直于屏幕的轴旋转，因此我们只需要一个标量来表示其方向。我们通常使用弧度（一个完整的圆周为 0 到 2π）作为单位，而不是角度（一个完整的圆周为 0 到 360），因为这可以简化计算。
为了旋转，刚体需要一定的角速度，这是一个标量，单位为弧度/秒，通常用希腊字母 ω (omega) 表示。然而，刚体要获得角速度，需要受到某种旋转力的作用，我们称之为扭矩，用希腊字母 τ (tau) 表示。因此，应用于旋转的牛顿第二定律为：
其中，α (alpha) 是角加速度，I 是转动惯量。
对于旋转，转动惯量类似于线性运动中的质量。它定义了改变刚体角速度的难度。在二维空间中，它是一个标量，定义如下：
其中 V 表示该积分应在物体体积内的所有点上进行，r 是每个点相对于旋转轴的位置向量，r² 实际上是 r 与自身的点积，𝜌 是一个函数，表示物体内每个点的密度。
例如，一个质量为 m、宽度为 w、高度为 h 的二维矩形绕其质心的转动惯量为：
在这里，你可以找到一系列公式，用于计算不同形状的物体绕不同轴的转动惯量。
当力作用在刚体的某一点时，它可能会产生扭矩。在二维空间中，扭矩是一个标量，作用在刚体某点上的力 f 产生的扭矩 τ 可以通过以下公式计算，其中 r 是从质心到作用点的偏移向量：
其中 θ (theta) 是 f 和 r 之间的最小夹角。
前面的公式正是 r 和 f 的叉积长度的公式。因此，在三维空间中，我们可以这样做(上面的τ是标量，下面的是向量)：
二维模拟可以看作是一种三维模拟，其中所有刚体都是薄而平的，并且所有运动都发生在 xy 平面上，这意味着在 z 轴上没有运动。这意味着 f 和 r 始终位于 xy 平面内，因此 τ 的 x 和 y 分量始终为零，因为叉积始终垂直于 xy 平面。这反过来意味着它始终平行于 z 轴。因此，只有叉积的 z 分量是重要的。这导致二维空间中的扭矩计算可以简化为：
令人难以置信的是，仅仅在模拟中增加一个维度就会使事情变得复杂得多。在三维空间中，刚体的方向必须用四元数表示，这是一种四元素向量。转动惯量由一个 3x3 矩阵表示，称为惯性张量，它不是恒定的，因为它依赖于刚体的方向，因此随着刚体的旋转而随时间变化。要了解 3D 刚体模拟的所有细节，你可以查看优秀的《刚体模拟I：无约束刚体动力学》，这也是 Witkin 和 Baraff 的《基于物理的建模：原理与实践》课程的一部分。
模拟算法与粒子模拟非常相似。我们只需要添加刚体的形状和旋转属性：
#define NUM_RIGID_BODIES 1
// 2D 盒子形状。物理引擎通常有一些不同类别的形状，
// 例如圆形、球体（3D）、圆柱体、胶囊体、多边形、多面体（3D）等...
typedef struct {
float width;
float height;
float mass;
float momentOfInertia;
} BoxShape;
// 计算盒子形状的转动惯量并存储在 momentOfInertia 变量中。
void CalculateBoxInertia(BoxShape *boxShape) {
float m = boxShape->mass;
float w = boxShape->width;
float h = boxShape->height;
boxShape->momentOfInertia = m * (w * w + h * h) / 12;
}
// 二维刚体
typedef struct {
Vector2 position; // 位置
Vector2 linearVelocity; // 线速度
float angle; // 角度
float angularVelocity; // 角速度
Vector2 force; // 力
float torque; // 扭矩
BoxShape shape; // 形状
} RigidBody;
// 全局刚体数组。
RigidBody rigidBodies[NUM_RIGID_BODIES];
// 在输出中打印每个刚体的位置和角度。
// 我们也可以在屏幕上绘制它们。
void PrintRigidBodies() {
for (int i = 0; i < NUM_RIGID_BODIES; ++i) {
RigidBody *rigidBody = &rigidBodies[i];
printf("body[%i] p = (%.2f, %.2f), a = %.2f\n", i, rigidBody->position.x, rigidBody->position.y, rigidBody->angle);
}
}
// 初始化刚体，随机位置和角度，线速度和角速度为零。
// 它们都被初始化为具有随机尺寸的盒子形状。
void InitializeRigidBodies() {
for (int i = 0; i < NUM_RIGID_BODIES; ++i) {
RigidBody *rigidBody = &rigidBodies[i];
rigidBody->position = (Vector2){arc4random_uniform(50), arc4random_uniform(50)};
rigidBody->angle = arc4random_uniform(360)/360.f * M_PI * 2;
rigidBody->linearVelocity = (Vector2){0, 0};
rigidBody->angularVelocity = 0;
BoxShape shape;
shape.mass = 10;
shape.width = 1 + arc4random_uniform(2);
shape.height = 1 + arc4random_uniform(2);
CalculateBoxInertia(&shape);
rigidBody->shape = shape;
}
}
// 在刚体的某个点施加力，产生扭矩。
void ComputeForceAndTorque(RigidBody *rigidBody) {
Vector2 f = (Vector2){0, 100}; // 施加的力
rigidBody->force = f;
// r 是从质心到力作用点的“力臂向量”
Vector2 r = (Vector2){rigidBody->shape.width / 2, rigidBody->shape.height / 2};
rigidBody->torque = r.x * f.y - r.y * f.x; // 计算扭矩
}
void RunRigidBodySimulation() {
float totalSimulationTime = 10; // 模拟将运行 10 秒。
float currentTime = 0; // 累积已过去的时间。
float dt = 1; // 每个时间步长为 1 秒。
InitializeRigidBodies();
PrintRigidBodies();
while (currentTime < totalSimulationTime) {
sleep(dt);
for (int i = 0; i < NUM_RIGID_BODIES; ++i) {
RigidBody *rigidBody = &rigidBodies[i];
ComputeForceAndTorque(rigidBody);
Vector2 linearAcceleration = (Vector2){rigidBody->force.x / rigidBody->shape.mass, rigidBody->force.y / rigidBody->shape.mass};
rigidBody->linearVelocity.x += linearAcceleration.x * dt;
rigidBody->linearVelocity.y += linearAcceleration.y * dt;
rigidBody->position.x += rigidBody->linearVelocity.x * dt;
rigidBody->position.y += rigidBody->linearVelocity.y * dt;
float angularAcceleration = rigidBody->torque / rigidBody->shape.momentOfInertia;
rigidBody->angularVelocity += angularAcceleration * dt;
rigidBody->angle += rigidBody->angularVelocity * dt;
}
PrintRigidBodies();
currentTime += dt;
}
}
调用 RunRigidBodySimulation 函数（针对单个刚体）将输出其位置和角度，如下所示：
body[0] p = (36.00, 12.00), a = 0.28
body[0] p = (36.00, 22.00), a = 15.28
body[0] p = (36.00, 42.00), a = 45.28
body[0] p = (36.00, 72.00), a = 90.28
body[0] p = (36.00, 112.00), a = 150.28
body[0] p = (36.00, 162.00), a = 225.28
body[0] p = (36.00, 222.00), a = 315.28
body[0] p = (36.00, 292.00), a = 420.28
body[0] p = (36.00, 372.00), a = 540.28
body[0] p = (36.00, 462.00), a = 675.28
body[0] p = (36.00, 562.00), a = 825.28
由于施加了一个向上的 100 牛顿的力，只有 y 坐标发生变化。这个力并没有直接作用在质心上，因此它产生了扭矩，导致旋转角度增加（逆时针旋转）。
图形化结果与粒子动画类似，但现在我们有形状在移动和旋转。
这是一个无约束的模拟，因为刚体可以自由移动，并且它们之间没有相互作用(它们不会发生碰撞)。没有碰撞的模拟相当无聊，也不太实用。在本系列的下一部分中，我们将讨论碰撞检测，以及一旦检测到碰撞后如何解决碰撞问题。
