约束的构建

在物理引擎中，约束 (Constraint)是一个用来限制物体运动的规则。无论是防止一个物体穿透另一个物体（接触约束），还是用一个铰链将两个物体连接在一起（关节约束），其本质都是在描述一个系统不应进入的“禁止状态”。现代物理引擎，特别是基于冲量的引擎，已经发展出一种极其强大和统一的框架来处理所有这些看似不同的问题：将它们都表述为约束方程，然后通过一个通用的约束求解器 (Constraint Solver) 来求解。

1. 约束方程 C(x) = 0

所有约束的出发点都是一个描述几何关系的方程。我们希望这个方程的值在任何时候都等于（或大于等于）零。

• 等式约束 (Equality Constraint): $C(\mathbf{x})=0$。例如，一个球关节约束要求两个物体上的两个锚点${\mathbf{p}}_A$和${\mathbf{p}}_B$必须重合。其约束方程就是： $$C(\mathbf{x})={\mathbf{p}}_A-{\mathbf{p}}_B=0$$

• 不等式约束 (Inequality Constraint): $C(\mathbf{x})>=0$。例如，一个非穿透约束要求两个物体之间的距离（或穿透深度的相反数）必须大于等于零。 $$C(\mathbf{x})=distance(A,B)>=0$$

这里的$\mathbf{x}$代表了系统中所有物体的广义坐标（位置和朝向）。

2. 速度级约束：C_dot = 0

直接在位置层面$C(x)=0$，求解约束非常困难，因为它是一个非线性方程组。一种更有效的方法是求解速度级 (velocity-level)约束。我们对约束方程$C(x)$求时间导数，得到$\dot{C}$，并要求$\dot{C}=0$。这意味着我们要求物体的速度必须满足一个条件，使得约束在下一瞬间仍然被满足。

$$\dot{C}=\frac{\partial C}{\partial x}\frac{dx}{dt}=\mathbf{Jv}=0$$

其中，$\mathbf{v}$：系统的广义速度向量，包含了所有物体的线速度 $v_i$ 和角速度 ${\omega }_i$。$\mathbf{J}$: 雅可比矩阵 (Jacobian Matrix)。这是整个约束动力学的核心。它描述了当物体的速度发生变化时，约束方程 $C$ 的值会如何变化。$J$的每一行对应一个约束，每一列对应一个物体的某个速度分量。

3. 雅可比矩阵 J 的构建

雅可比矩阵 J 是约束 C 对广义坐标 x 的偏导数。对于一个只涉及两个物体 A 和 B 的约束，其雅可比矩阵的一行通常有12列（物体A的3个线速度分量+3个角速度分量，物体B的3个线速度分量+3个角速度分量）：

$$J=[J_{v}A,J_{\omega }A,J_{v}B,J_{\omega }B]$$

$$v=[v_A,{\omega }_A,v_B,{\omega }_B{]}^T$$

$$\dot{C}=J_{v}A\cdot v_A+J_{\omega }A\cdot {\omega }_A+J_{v}B\cdot v_B+J_{\omega }B\cdot {\omega }_B$$

示例：接触约束 (Contact Constraint)

考虑一个在${\mathbf{p}}_A$和${\mathbf{p}}_B$点接触的非穿透约束，法线为$\mathbf{n}$。约束方程为

$$C=(p_A-p_B)\cdot n>=0。$$

我们要求$\dot{C}>=0$。

$$\dot{C}=(v_A+{\omega }_A\times r_A-(v_B+{\omega }_B\times r_B))\cdot n$$

其中 $r_A$ 和 $r_B$ 是从质心到接触点的向量。

$$\dot{C}=v_A\cdot n+{\omega }_A\cdot (r_A\times n)-v_B\cdot n-{\omega }_B\cdot (r_B\times n)$$

通过将这个表达式与$J\cdot v$的形式进行比较，我们可以直接读出雅可比矩阵的各个部分： $$J_{v}A=n$$

$$J_{\omega }A=r_A\times n$$

$$J_{v}B=-n$$

$$J_{\omega }B=-(r_B\times n)$$

这组向量$[n,r_A\times n,-n,-(r_B\times n)]$就是这个接触约束的雅可比“行”。它精确地描述了两个物体的线速度和角速度如何影响它们在接触点沿法线方向的分离速度。

4. 冲量与雅可比的关系

约束是通过施加力或冲量来满足的。一个约束力/冲量 λ 会导致物体的速度发生变化 Δv。这个变化与施加的冲量之间，通过雅可比矩阵的转置 J^T 联系起来：

$$\mathbf{M}\Delta \mathbf{v}={\mathbf{J}}^T\lambda$$

• λ: 拉格朗日乘子 (Lagrange Multiplier)，在基于冲量的求解器中，它代表了满足约束所需的冲量大小。
• M: 广义质量矩阵，一个包含了所有物体质量和惯性张量的对角矩阵。
• ${\mathbf{J}}^T$雅可比矩阵的转置。它将一个标量冲量 λ “分配”到各个物体的线动量和角动量上。

从上式可以解出速度变化 Δv：

$$\Delta v=M^{-1}{J}^T\lambda$$

5. 最终的约束方程

我们的目标是找到一个冲量 λ，使得施加该冲量后的新速度 $v_{n}ew=v_{o}ld+\Delta v$能够满足速度级约束 $J\cdot v_{n}ew=0$（或更通用的$J\cdot v_{n}ew=b$，其中$b$是一个目标速度，例如用于模拟反弹）。

$$J(v_{old}+\Delta v)=b$$

$$J\cdot v_{old}+J\cdot \Delta v=b$$

$$J(M^{-1}{J}^T\lambda )=b-J\cdot v_{old}$$

令 $K=J{M}^{-1}{J}^T$，我们得到一个关于未知冲量 λ 的线性方程：

$$\mathbf{K}\lambda =\mathbf{b}-\mathbf{J}{\mathbf{v}}_{old}$$

• K: 这个标量（或小矩阵）被称为有效质量 (Effective Mass)。它代表了在约束方向上，系统抵抗速度变化的惯性。K 的计算 J M⁻¹ J^T。
• J·v_old: 碰撞前的相对速度。
• b: 目标相对速度。对于简单的接触，b=0。为了处理弹性碰撞（恢复系数 e），我们可以设置 b = -e * (J·v_old)。为了修正位置误差（穿透），我们可以引入一个 Baumgarte 稳定化项 b = -β/Δt * C(x)。

伪代码示例 (构建接触约束的雅可比和有效质量):

void buildContactConstraint(Constraint& c, RigidBody& A, RigidBody& B, ContactPoint& p) {
    Vector3 rA = p.position - A.position;
    Vector3 rB = p.position - B.position;
    Vector3 n = p.normal;

    // 构建雅可比矩阵
    c.J_vA = n;
    c.J_wA = rA.cross(n);
    c.J_vB = -n;
    c.J_wB = -rB.cross(n);

    // 计算有效质量的倒数 (K = J * M_inv * J^T)
    float invEffectiveMass = A.inverseMass + B.inverseMass + 
                             c.J_wA.dot(A.inverseInertiaTensorWorld * c.J_wA) + 
                             c.J_wB.dot(B.inverseInertiaTensorWorld * c.J_wB);
    c.effectiveMass = 1.0f / invEffectiveMass;

    // 计算 Baumgarte 稳定化项 (用于修正穿透)
    float beta = 0.2f;
    float slop = 0.01f;
    c.baumgarte_bias = (beta / dt) * max(0.0f, p.penetration - slop);
}

总结

约束的构建是将物理问题转化为数学问题的关键一步。通过定义一个几何约束方程$C(x)$，并对其求导，我们得到了一个线性的速度级约束 $J\cdot v=0$。雅可比矩阵 J 在此过程中扮演了核心角色，它像一座桥梁，连接了物体的运动速度和约束的变化。最终，我们得到了一个关于未知约束冲量$\lambda$的线性方程$K\lambda =b-Jv$。这个方程是统一的、普适的，无论是接触、摩擦还是关节，都可以被转化为这个形式。如何高效地求解这个由成百上千个此类方程构成的巨大系统，就是下一章——迭代求解器——将要解决的问题。