加速球体之间的碰撞检测

Question

我正在编写一个物理引擎/模拟器，其中包含 3D 太空飞行、行星/恒星引力、船舶推力和相对论效应。到目前为止，一切进展顺利，但是，我需要帮助的一件事是碰撞检测算法的数学。

我正在使用的运动迭代模拟基本上如下：（

注意：3D向量全部大写。）

For each obj

    obj.ACC = Sum(all acceleration influences)

    obj.POS = obj.POS + (obj.VEL * dT) + (obj.ACC * dT^2)/2     (*EQ.2*)

    obj.VEL = obj.VEL + (obj.ACC * dT)

Next

其中：

obj.ACC is the acceleration vector of the object
obj.POS is the position or location vector of the object
obj.VEL is the velocity vector of the object

obj.Radius is the radius (scalar) of the object

dT is the time delta or increment

我基本上需要做的是找到一些有效的公式，该公式源自（EQ.2) 上面的两个对象 (obj1, obj2) 并告诉它们是否发生碰撞，如果发生，发生在什么时间。我需要准确的时间，以便我可以确定它是否在这个特定的时间增量中（因为加速度在不同的时间增量下会有所不同），并且还可以定位准确的位置（我知道如何做，因为给定时间）

对于这个引擎，我将所有对象建模为球体，所有这个公式/算法需要做的就是找出在哪些点：

(obj1.POS - obj2.POS).Distance = (obj1.Radius + obj2.Radius)

其中 .Distance 是正标量值。 （如果更容易的话，您也可以对两边进行平方，以避免 .Distance 计算中隐含的平方根函数）。

（是的，我知道很多很多其他碰撞检测问题，但是，他们的解决方案似乎都非常适合其引擎和假设，并且似乎没有一个符合我的条件：3D、球体和在模拟增量中应用的加速度如果我错了，请告诉我。）

---

一些澄清：

1）对于我来说，检查时间增量之前和之后两个球体的交点是不够的。在许多情况下，它们的速度和位置变化将远远超过它们的半径。

2）回复：效率，我不需要帮助（无论如何在这一点上）来确定可能的碰撞候选者，我认为我已经涵盖了。

---

出现了很多：

3）我的增量运动方程（EQ.2）是一个同时应用速度和加速度的二次方程：

obj.POS = obj.POS + (obj.VEL * dT) + (obj.ACC * dT^2)/2

另一个澄清，似乎 我见过的物理引擎（当然还有我听说过的每个游戏引擎）仅线性增量运动方程仅应用速度：

obj.POS = obj.POS + (obj.VEL * dT)

这就是为什么我不能使用 StackOverflow、维基百科和整个网络上常见的碰撞检测解决方案，例如查找两条线段的交点/最近路径。我的模拟处理对结果至关重要的可变加速度，因此我需要的是两个抛物线段的交集/最接近方法。

Answer 1

在 AShelley 提到的网页上，最近点接近方法是针对两个物体匀速运动的情况开发的。然而，我相信相同的矢量微积分方法可用于在两个物体都以恒定的非零加速度（二次时间依赖性）移动的情况下得出结果。

在这种情况下，距离平方函数的时间导数是三阶（三次）而不是一阶。因此，“最近接近时间”将有 3 个解决方案，这并不奇怪，因为两个对象的路径都是弯曲的，因此可能存在多个交叉点。对于此应用程序，您可能希望使用当前模拟步骤定义的间隔内的 t 的最早值（如果存在这样的时间）。

我计算出了导数方程，它应该给出最接近的时间：

0 = |D_ACC|^2 * t^3 + 3 * dot(D_ACC, D_VEL) * t^2 + 2 * [ |D_VEL|^ 2 + 点(D_POS, D_ACC) ] * t + 2 * 点(D_POS, D_VEL)

其中：

D_ACC = ob1.ACC-obj2.ACC

D_VEL = ob1.VEL -obj2.VEL（更新前）

D_POS = ob1.POS-obj2.POS（也在更新前）

和dot(A, B) = Ax*Bx + Ay*通过 + Az*Bz

（注意，幅值的平方|A|^2可以使用dot(A, A)计算）

来解决这个问题对于t，您可能需要使用找到的公式维基百科。

当然，这只会给您最接近的时刻。此时您需要测试距离（使用类似于方程 2 的方法）。如果它大于(obj1.Radius + obj2.Radius)，则可以忽略它（即不发生碰撞）。然而，如果距离较小，则意味着球体在此时刻之前发生碰撞。然后，您可以使用迭代搜索来测试早期的距离。也可能提出另一种（甚至更复杂的）考虑尺寸的推导，或者可能找到其他一些解析解，而无需诉诸迭代求解。

编辑：由于阶数较高，方程的某些解实际上是最远距离的矩。我相信在所有情况下，3 个解决方案中的 1 个或 3 个解决方案中的 2 个将是分离最远的时间。您可以通过评估相对于时间的二阶导数（通过将一阶导数设置为零而找到的 t 值）来分析测试您是否处于最小值或最大值：

D''(t) = 3 * |D_ACC|^2 * t^2 + 6 * dot(D_ACC, D_VEL) * t + 2 * [ |D_VEL|^2 + dot(D_POS, D_ACC) ]

如果二阶导数计算为一个正数，那么您就知道在给定时间 t 内距离是最小值，而不是最大值。

Answer 2

在每个球体的起始位置和结束位置之间绘制一条线。如果生成的线段相交，则球体肯定在某个点发生碰撞，并且一些聪明的数学可以找到碰撞发生的时间。还要确保检查线段之间的最小距离（如果它们不相交）是否小于 2*半径。这也表明发生了碰撞。

从那里您可以倒退您的增量时间以准确地发生碰撞，以便您可以正确计算力。

您是否考虑过使用已经可以实现此功能的物理库？许多图书馆使用更先进、更稳定（更好的积分器）的系统来求解您正在使用的方程组。我想到了子弹物理。

Answer 3

op询问碰撞时间。稍微不同的方法将精确计算它......

请记住，位置投影方程为：

NEW_POS=POS+VEL*t+(ACC*t^2)/2

如果我们替换 POS< /code> 与 D_POS=POS_A-POS_B、VEL 与 D_VEL=VEL_A-VEL_B 和 ACC=ACC_A-ACC_B对于对象 A 和 B，我们得到：

$D_NEW_POS=D_POS+D_VEL*t+(D_ACC*t^2)/2

这是对象之间的矢量距离公式。为了得到它们之间的标量距离的平方，我们可以取这个方程的平方，展开后如下：

distsq(t) = D_POS^2+2*dot(D_POS,D_VEL)*t + ( dot(D_POS, D_ACC)+D_VEL^2)*t^2 + dot(D_VEL,D_ACC)*t^3 + D_ACC^2*t^4/4

为了找到碰撞发生的时间，我们可以将方程设置为半径总和的平方并求解t：

0 = D_POS^2-(r_A+r_B)^2 + 2*dot(D_POS,D_VEL )*t + (dot(D_POS, D_ACC)+D_VEL^2)*t^2 + dot(D_VEL,D_ACC)*t^3 + D_ACC^2*t^4/4

现在，我们可以求解使用四次公式计算方程。

四次公式将产生 4 个根，但我们只对实数根感兴趣。如果存在双实根，则两个对象恰好在一个时间点接触边缘。如果有两个真实的根，则对象在根1和根2之间连续重叠（即根1是碰撞开始的时间，根2是碰撞停止的时间）。四个实根意味着物体在根对 1,2 和 3,4 之间连续碰撞两次。

在R中，我使用polyroot()解决如下：

# initial positions
POS_A=matrix(c(0,0),2,1)
POS_B=matrix(c(2,0),2,1)
# initial velocities
VEL_A=matrix(c(sqrt(2)/2,sqrt(2)/2),2,1)
VEL_B=matrix(c(-sqrt(2)/2,sqrt(2)/2),2,1)
# acceleration
ACC_A=matrix(c(sqrt(2)/2,sqrt(2)/2),2,1)
ACC_B=matrix(c(0,0),2,1)
# radii
r_A=.25
r_B=.25
# deltas
D_POS=POS_B-POS_A
D_VEL=VEL_B-VEL_A
D_ACC=ACC_B-ACC_A


# quartic coefficients
z=c(t(D_POS)%*%D_POS-r*r, 2*t(D_POS)%*%D_VEL, t(D_VEL)%*%D_VEL+t(D_POS)%*%D_ACC, t(D_ACC)%*%D_VEL, .25*(t(D_ACC)%*%D_ACC))
# get roots
roots=polyroot(z)
# In this case there are only two real roots...
root1=as.numeric(roots[1])
root2=as.numeric(roots[2])

# trajectory over time
pos=function(p,v,a,t){
  T=t(matrix(t,length(t),2))
  return(t(matrix(p,2,length(t))+matrix(v,2,length(t))*T+.5*matrix(a,2,length(t))*T*T))
}

# plot A in red and B in blue
t=seq(0,2,by=.1) # from 0 to 2 seconds.
a1=pos(POS_A,VEL_A,ACC_A,t)
a2=pos(POS_B,VEL_B,ACC_B,t)
plot(a1,type='o',col='red')
lines(a2,type='o',col='blue')

# points of a circle with center 'p' and radius 'r'
circle=function(p,r,s=36){
  e=matrix(0,s+1,2)
  for(i in 1:s){
    e[i,1]=cos(2*pi*(1/s)*i)*r+p[1]
    e[i,2]=sin(2*pi*(1/s)*i)*r+p[2]
  }
  e[s+1,]=e[1,]
  return(e)
}

# plot circles with radius r_A and r_B at time of collision start in black
lines(circle(pos(POS_A,VEL_A,ACC_A,root1),r_A))
lines(circle(pos(POS_B,VEL_B,ACC_B,root1),r_B))
# plot circles with radius r_A and r_B at time of collision stop in gray
lines(circle(pos(POS_A,VEL_A,ACC_A,root2),r_A),col='gray')
lines(circle(pos(POS_B,VEL_B,ACC_B,root2),r_B),col='gray')

对象 A 遵循红色轨迹从左下到右上。对象B遵循从右下角到左上角的蓝色轨迹。两个物体在时间 0.9194381 和时间 1.167549 之间连续碰撞。两个黑色圆圈刚刚接触，表明重叠的开始 - 并且重叠会及时继续，直到物体到达灰色圆圈的位置。

Answer 4

看来您想要 最近接近点 (CPA)。如果它小于半径之和，就会发生碰撞。链接中有示例代码。您可以使用当前速度计算每一帧，并检查 CPA 时间是否小于您的刻度大小。您甚至可以缓存 CPA 时间，并且仅在对任一项目应用加速时进行更新。