随机太阳系模拟碰撞检测问题

Randomised solar system simulation collision detection issue

因此,我使用基本的运动方程和牛顿力学并使用 OpenGL 和 C 开发了一个相当稳定的模拟太阳系的随机行星。

当行星相互碰撞或太阳碰撞时,我也有一些碰撞检测。

我的问题是,我当前的代码会检测到碰撞并相应地移除行星,但是,在查看模拟行星时,行星可以 overlap/collide 减少大约一半的半径(只是根据我的猜测)看)在碰撞之前是 'detected'.

我有理由相信我的碰撞检测方法是足够的,但如果我错了请纠正我。

我唯一能想到的另一个问题是 'units' 在世界坐标、眼睛坐标等方面存在一些差异。

非常感谢您提供的所有帮助 and/or 建议。

代码:

float* vecSub( float *pV0, float *pV1, float *pVRes )
{
    if (pV0 && pV1 && pVRes)
    {
        pVRes[0] = pV0[0] - pV1[0];
        pVRes[1] = pV0[1] - pV1[1];
        pVRes[2] = pV0[2] - pV1[2];
        return pVRes;
    }
    return 0;
}

float vecLength( float *pV )
{
    if(pV) return sqrtf(pV[0]*pV[0]+pV[1]*pV[1]+pV[2]*pV[2]);
    return 0.0f;
}

float vecDistance( float *pV1, float *pV2 )
{
    float fLen=0.0f;

    if(pV1 && pV2)
    {
        float av[4];

        vecSub(pV2, pV1, av);
        fLen=vecLength(av);
    }

    return fLen;
}

if (vecDistance(pPlanet->afPosition, pPlanet1->afPosition) <= pPlanet->fRadius + pPlanet1->fRadius)
{
    //Collision resolution code
}

我假设你有一个离散时间模拟,具有固定的时间步长。

您的碰撞检测似乎还不够。如果行星的速度足够高,并且它们的半径足够小,您可能会让它们彼此擦肩而过,而不会检测到碰撞。如果碰撞不是正面碰撞,而是切向碰撞,情况会变得更糟。

更谨慎的做法是不观察行星,而是观察半径更大的球(动态取决于行星的速度),当两个 碰撞时,切换到更精细的时间步长。

如果你走这条路,请小心地将模拟与可视化分开。

假设您使用一个小结构来描述坐标,另一个用于轴对齐的边界框:

typedef struct {
    double  x;
    double  y;
    double  z;
} vec3d;

typedef struct {
    vec3d   min;
    vec3d   max;
} box3d; /* Axis-aligned bounding box */

假设您模拟 N 个球形物体,半径在 double radius[] 中描述,当前坐标在 vec3d curr[](或 vec3d *curr)中,之前的坐标在 vec3d prev[](或vec3d *prev)。

我们想检查是否有任何对象发生碰撞或相交,但这样做效率很高。 为避免做不必要的工作,请为轴对齐的边界框使用一个数组,以便该框包含先前坐标和当前坐标处的球形对象。如果两个对象的轴对齐边界框相交,则两个对象只能碰撞或相交:

static inline double dmin(const double d1, const double d2)
{
    return (d1 <= d2) ? d1 : d2;
}

static inline double dmax(const double d1, const double d2)
{
    return (d1 >= d2) ? d1 : d2;
}

void update_boxes(box3d *box, const size_t count,
                  const vec3d *curr, const vec3d *prev,
                  const double *radius)
{
    size_t  i;
    for (i = 0; i < count; i++) {
        box[i].min.x = dmin(curr[i].x, prev[i].x) - radius[i];
        box[i].max.x = dmax(curr[i].x, prev[i].x) + radius[i];
        box[i].min.y = dmin(curr[i].y, prev[i].y) - radius[i];
        box[i].max.y = dmax(curr[i].y, prev[i].y) + radius[i];
        box[i].min.z = dmin(curr[i].z, prev[i].z) - radius[i];
        box[i].max.z = dmax(curr[i].z, prev[i].z) + radius[i];
    }
}

为了计算方框,我们采用中心坐标,并按半径扩展方框。 (轴对齐的边界框不需要精确,它们只需要在两个位置覆盖球体。如果它们更大,则意味着我们将做一些不必要的工作。我们使用轴对齐的边界框,因为检查其中两个是否相交非常快。)

如果我们假设两个球形物体在每个时间步长中都具有恒定的速度,那么在时间步长中时间t0 <= t && t <= 1)的位置是

xi(t) = (1-t)*prev[i].x + t*curr[i].x;
yi(t) = (1-t)*prev[i].y + t*curr[i].y;
zi(t) = (1-t)*prev[i].z + t*curr[i].z;

这在两个位置之间的 3D 中很简单 linear interpolationt = 0 在前一个时间步,t = 1 在当前时间步。

如果我们在索引 ik 处构建两个此类对象之间的平方距离方程,我们得到

L(t) = SQUARE( ((1-t)*prev[i].x + t*curr[i].x) - ((1-t)*prev[k].x + t*prev[k].x) )
     + SQUARE( ((1-t)*prev[i].y + t*curr[i].y) - ((1-t)*prev[k].y + t*prev[k].y) )
     + SQUARE( ((1-t)*prev[i].z + t*curr[i].z) - ((1-t)*prev[k].z + t*prev[k].z) )

其中 SQUARE(expr) = (expr)*(expr)。当其导数为零时,它达到最小值。如果我们求解 t,我们发现恰好有一个实根,这意味着两个物体以恒定速度沿着这些线性路径,在时间 t:[=37= 彼此最接近]

t = ( (prev[i].x - prev[k].x) * ( (prev[i].x - prev[j].x) - (curr[i].x - curr[k].x) )
    + (prev[i].y - prev[k].y) * ( (prev[i].y - prev[j].y) - (curr[i].y - curr[k].y) )
    + (prev[i].z - prev[k].z) * ( (prev[i].z - prev[j].z) - (curr[i].z - curr[k].z) )
    ) / ( SQUARE( (prev[i].x - prev[k].x) - (curr[i].x - curr[k].x) )
        + SQUARE( (prev[i].y - prev[k].y) - (curr[i].y - curr[k].y) )
        + SQUARE( (prev[i].z - prev[k].z) - (curr[i].z - curr[k].z) ) )

这仅在除数非零时有效(它永远不会为负,因为它是平方和)。

我们只对 t >= 0t <= 1 的情况感兴趣;也就是说,两个对象在前一个时间步和当前时间步之间彼此最接近的情况。

如果确实发生了这种情况,我们需要做的就是将 t 代入方程式,并与 SQUARE(radius[i] + radius[k]) 进行比较以确定两个物体是否发生碰撞。

让我们看一个示例函数,它也计算轴对齐的边界框,使用它们进行快速选择,并且应该正确检测碰撞:

void handle_collisions(const size_t        count,
                       box3d       *const  box,
                       vec3d       *const  curr,
                       const vec3d *const  prev,
                       const double *const radius)
{
    size_t  i, k;

    for (k = 0; k < count; k++) {
        box[k].min.x = dmin(prev[k].x, curr[k].x) - radius[k];
        box[k].max.x = dmax(prev[k].x, curr[k].x) + radius[k];
        box[k].min.y = dmin(prev[k].y, curr[k].y) - radius[k];
        box[k].max.y = dmax(prev[k].y, curr[k].y) + radius[k];
        box[k].min.z = dmin(prev[k].z, curr[k].z) - radius[k];
        box[k].max.z = dmax(prev[k].z, curr[k].z) + radius[k];

        for (i = 0; i < k; i++) {
            if (box[k].min.x <= box[i].max.x &&
                box[k].min.y <= box[i].max.y &&
                box[k].min.z <= box[i].max.z &&
                box[k].max.x >= box[i].min.x &&
                box[k].max.y >= box[i].max.y &&
                box[k].max.z >= box[i].max.z) {

                /* A collision is possible, since the axis-aligned
                   bounding boxes intersect. Check. */

                const vec3d p = { prev[i].x - prev[k].x,
                                  prev[i].y - prev[k].y,
                                  prev[i].z - prev[k].z };
                const vec3d d = { p.x - curr[i].x + curr[k].x,
                                  p.y - curr[i].y + curr[k].y,
                                  p.z - curr[i].z + curr[k].z };
                const double tn = p.x * d.x + p.y * d.y + p.z * d.z;
                const double td = d.x * d.x + d.y * d.y + d.z * d.z;
                if (tn >= 0.0 && tn <= td) {
                    const double  t1 = tn / td;
                    const double  t0 = 1.0 - t1;
                    const vec3d   loc_k = { t0*prev[k].x + t1*curr[k].x,
                                            t0*prev[k].y + t1*curr[k].y,
                                            t0*prev[k].z + t1*curr[k].z };
                    const vec3d   loc_i = { t0*prev[i].x + t1*curr[i].x,
                                            t0*prev[i].y + t1*curr[i].y,
                                            t0*prev[i].z + t1*curr[i].z };
                    const vec3d   delta = { loc_i.x - loc_k.x,
                                            loc_i.y - loc_k.y,
                                            loc_i.z - loc_k.z };
                    if (delta.x*delta.x + delta.y*delta.y + delta.z*delta.z <=
                        (radius[i] + radius[k])*(radius[i] + radius[k])) {

                        /* Collision occurs at time t  (0 <= t && t <= 1),
                           between object k (at loc_k) and object i (at loc_i).
                        */
                    }
                }
            }
        }
    }
}

从技术上讲,两个以上的物体可以在同一时间步内发生碰撞,尽管这种情况极为罕见。问题是,如果您 "remove" 在上述循环中考虑的任何碰撞对象,您可能会错过它们。 (我怀疑不担心的情况很少见,但我是一个偏执的人,不喜欢不担心的事情。:)

为了避免这种情况,我只是将上述循环中的冲突保存到某种数组或 disjoint-set data structure 中。然后,在上面的代码之后,合并任何碰撞的对象。 (请注意,不相交集数据结构使这更容易,因为如果你有冲突 A-B、B-C 和 C-D,不相交集实际上将其解析为 A-B、A-C 和 A-D。否则,如果你将 B 连接到 A,并销毁B,将 C 连接到 B 是不可能的。同样,这是一个极端情况,但不担心这种极端情况是拥有可靠的模拟器与每千次模拟或十亿分之一崩溃一次的模拟器之间的区别,出于一个不为人知的不明原因。)