运动点与运动线段的碰撞时刻及位置求解（连续碰撞检测）

Your initial modeling of the problem is spot-on—linear trajectories for all points is the right approach here. The issue with SymPy's massive solution is that you're solving the full system of two equations directly, but we can simplify this drastically by leveraging 2D geometric properties (collinearity via cross product) to reduce the problem to a single quadratic equation, which is trivial to solve in real-time engine code.

Core Insight: Collision as Collinearity + Projection #

For the dynamic point ( P(t) ) to collide with the dynamic segment ( UV(t) ) at time ( t \in [0,1] ), two conditions must hold:

1. ( P(t) ) is collinear with ( U(t) ) and ( V(t) ) (cross product of ( P(t)-U(t) ) and ( V(t)-U(t) ) equals 0)
2. ( P(t) ) lies between ( U(t) ) and ( V(t) ) (projection parameter ( a \in [0,1] ))

Step 1: Simplify Trajectories with Displacement Vectors #

First, rewrite all motion using displacement increments instead of start/end points to clean up the math:

• ( \Delta P = eP - sP ), so ( P(t) = sP + \Delta P \cdot t )
• ( \Delta U = eU - sU ), so ( U(t) = sU + \Delta U \cdot t )
• ( \Delta V = eV - sV ), so ( V(t) = sV + \Delta V \cdot t )

Define helper vectors:

• Initial position differences: ( PU_0 = sP - sU ), ( VU_0 = sV - sU )
• Displacement differences: ( \Delta PU = \Delta P - \Delta U ), ( \Delta VU = \Delta V - \Delta U )

Step 2: Build the Collinearity Quadratic Equation #

Collinearity means the cross product of ( P(t)-U(t) ) and ( V(t)-U(t) ) is zero. Expanding this gives a quadratic equation in ( t ):
[ A t^2 + B t + C = 0 ]
Where:

• ( A = \text{cross}(\Delta PU, \Delta VU) ) (2D cross product: ( x_1 y_2 - x_2 y_1 ))
• ( B = \text{cross}(PU_0, \Delta VU) + \text{cross}(\Delta PU, VU_0) )
• ( C = \text{cross}(PU_0, VU_0) )

Step 3: Solve the Quadratic Equation #

Solve for ( t ) using the quadratic formula:

• Discriminant: ( D = B^2 - 4AC )
  • If ( D < 0 ): No real roots → no collision
  • If ( D = 0 ): One real root ( t = -B/(2A) )
  • If ( D > 0 ): Two real roots ( t_1 = (-B - \sqrt{D})/(2A) ), ( t_2 = (-B + \sqrt{D})/(2A) )

Step 4: Validate Collision Candidates #

For each candidate ( t \in [0,1] ):

1. Compute ( U(t) ), ( V(t) ), ( P(t) )
2. Check if ( P(t) ) lies on segment ( UV(t) ):
   • Calculate vector ( VU = V(t) - U(t) ) and ( PU = P(t) - U(t) )
   • If ( |VU|^2 \approx 0 ): Segment is a point—check if ( P(t) ) equals this point
   • Else, compute projection parameter ( a = \frac{\text{dot}(PU, VU)}{|VU|^2} )
   • If ( a \in [0,1] ) (with small epsilon for floating-point error), this is a valid collision

Step 5: Handle Edge Cases #

• Initial collision: Check if ( P(0) ) is on ( UV(0) ) first (return ( t=0 ))
• Segment collapses to a point: Solve ( U(t) = V(t) ) (linear equation) and check if ( P(t) ) matches that point
• Floating-point precision: Use a small epsilon (e.g., ( 1e-6 )) when comparing values for equality or range checks

Real-Time Implementation (Pseudocode) #

struct Vec2 { float x, y; };
Vec2 operator-(Vec2 a, Vec2 b) { return {a.x - b.x, a.y - b.y}; }
Vec2 operator+(Vec2 a, Vec2 b) { return {a.x + b.x, a.y + b.y}; }
Vec2 operator*(Vec2 v, float s) { return {v.x*s, v.y*s}; }
float cross(Vec2 a, Vec2 b) { return a.x*b.y - a.y*b.x; }
float dot(Vec2 a, Vec2 b) { return a.x*b.x + a.y*b.y; }
float clamp(float val, float min, float max) { return std::max(min, std::min(max, val)); }

const float EPS = 1e-6;

bool dynamicPointDynamicSegmentCCD(
    Vec2 sP, Vec2 eP,
    Vec2 sU, Vec2 eU,
    Vec2 sV, Vec2 eV,
    float& out_t, Vec2& out_pos
) {
    Vec2 deltaP = eP - sP;
    Vec2 deltaU = eU - sU;
    Vec2 deltaV = eV - sV;

    Vec2 PU0 = sP - sU;
    Vec2 VU0 = sV - sU;
    Vec2 deltaPU = deltaP - deltaU;
    Vec2 deltaVU = deltaV - deltaU;

    float A = cross(deltaPU, deltaVU);
    float B = cross(PU0, deltaVU) + cross(deltaPU, VU0);
    float C = cross(PU0, VU0);

    float earliest_t = 2.0f; // Initialize to value >1 (no collision)

    // Handle linear case (A ≈ 0)
    if (fabs(A) < EPS) {
        if (fabs(B) < EPS) {
            // Check initial collinearity
            if (fabs(C) < EPS) {
                Vec2 VU = sV - sU;
                float len_sq = dot(VU, VU);
                if (len_sq < EPS) {
                    if (dot(sP - sU, sP - sU) < EPS) {
                        out_t = 0.0f;
                        out_pos = sP;
                        return true;
                    }
                } else {
                    float a = dot(sP - sU, VU) / len_sq;
                    if (a >= -EPS && a <= 1.0f + EPS) {
                        out_t = 0.0f;
                        out_pos = sP;
                        return true;
                    }
                }
            }
            return false;
        } else {
            float t = -C / B;
            if (t >= -EPS && t <= 1.0f + EPS) {
                t = clamp(t, 0.0f, 1.0f);
                Vec2 U = sU + deltaU * t;
                Vec2 V = sV + deltaV * t;
                Vec2 P = sP + deltaP * t;
                Vec2 VU = V - U;
                float len_sq = dot(VU, VU);
                if (len_sq < EPS) {
                    if (dot(P - U, P - U) < EPS && t < earliest_t) {
                        earliest_t = t;
                        out_pos = P;
                    }
                } else {
                    float a = dot(P - U, VU) / len_sq;
                    if (a >= -EPS && a <= 1.0f + EPS && t < earliest_t) {
                        earliest_t = t;
                        out_pos = P;
                    }
                }
            }
        }
    } else {
        // Solve quadratic equation
        float D = B*B - 4*A*C;
        if (D < -EPS) return false;
        D = std::max(D, 0.0f);
        float sqrtD = sqrt(D);
        float t1 = (-B - sqrtD) / (2*A);
        float t2 = (-B + sqrtD) / (2*A);

        auto check_t = [&](float t) {
            if (t < -EPS || t > 1.0f + EPS) return;
            t = clamp(t, 0.0f, 1.0f);
            Vec2 U = sU + deltaU * t;
            Vec2 V = sV + deltaV * t;
            Vec2 P = sP + deltaP * t;
            Vec2 VU = V - U;
            float len_sq = dot(VU, VU);
            if (len_sq < EPS) {
                if (dot(P - U, P - U) < EPS && t < earliest_t) {
                    earliest_t = t;
                    out_pos = P;
                }
            } else {
                float a = dot(P - U, VU) / len_sq;
                if (a >= -EPS && a <= 1.0f + EPS && t < earliest_t) {
                    earliest_t = t;
                    out_pos = P;
                }
            }
        };

        check_t(t1);
        check_t(t2);
    }

    // Check if segment collapses to a point at any t
    Vec2 deltaUV = deltaU - deltaV;
    Vec2 SV0 = sV - sU;
    if (fabs(cross(deltaUV, SV0)) < EPS) {
        float dot_uv = dot(deltaUV, deltaUV);
        if (fabs(dot_uv) < EPS) {
            // U and V are always the same point
            Vec2 uv_point = sU;
            if (fabs(cross(deltaP, uv_point - sP)) < EPS) {
                float dot_p = dot(deltaP, deltaP);
                if (fabs(dot_p) < EPS) {
                    if (dot(sP - uv_point, sP - uv_point) < EPS && 0.0f < earliest_t) {
                        earliest_t = 0.0f;
                        out_pos = sP;
                    }
                } else {
                    float t = dot(uv_point - sP, deltaP) / dot_p;
                    if (t >= -EPS && t <= 1.0f + EPS) {
                        t = clamp(t, 0.0f, 1.0f);
                        if (t < earliest_t) {
                            earliest_t = t;
                            out_pos = uv_point;
                        }
                    }
                }
            }
        } else {
            float t0 = dot(SV0, deltaUV) / dot_uv;
            if (t0 >= -EPS && t0 <= 1.0f + EPS) {
                t0 = clamp(t0, 0.0f, 1.0f);
                Vec2 uv_point = sU + deltaU * t0;
                Vec2 p_point = sP + deltaP * t0;
                if (dot(p_point - uv_point, p_point - uv_point) < EPS && t0 < earliest_t) {
                    earliest_t = t0;
                    out_pos = uv_point;
                }
            }
        }
    }

    if (earliest_t <= 1.0f + EPS) {
        out_t = earliest_t;
        return true;
    }
    return false;
}

Why This Works for Your Requirements #

• Dynamic motion: Handles all linear motion of the point and segment (including rotation, stretching, and translation)
• Exact solution: Uses algebraic quadratic equation solving—no iterative approximations or sweep volumes
• Collision details: Returns the earliest collision time ( t ) and exact collision position
• Edge case coverage: Handles segments collapsing to points, initial collisions, and floating-point precision issues

内容的提问来源于stack exchange，提问作者PawelBoe