C++求两条贝塞尔曲线的交点

求两条贝塞尔曲线的交点，可以通过求解方程组来实现。假设两条贝塞尔曲线分别为 $B_1(t)$ 和 $B_2(t)$，则它们的交点满足以下方程组：

$$
\begin{cases}
B_1(t_1) - B_2(t_2) = 0\\
B_1'(t_1) \cdot B_2''(t_2) - B_2'(t_2) \cdot B_1''(t_1) = 0
\end{cases}
$$

其中 $B_i(t)$ 表示第 $i$ 条贝塞尔曲线在参数 $t$ 处的点，$B_i'(t)$ 表示第 $i$ 条贝塞尔曲线在参数 $t$ 处的切向量，$B_i''(t)$ 表示第 $i$ 条贝塞尔曲线在参数 $t$ 处的二阶导向量。

具体实现可以采用牛顿迭代法来求解方程组的根。

以下是一个简单的 C++ 实现：

#include <iostream>
#include <cmath>

using namespace std;

const double EPSILON = 1e-6; // 停止迭代的精度要求

struct Point {
    double x;
    double y;
};

struct BezierCurve {
    Point p0;
    Point p1;
    Point p2;
    Point p3;
};

Point bezier(BezierCurve curve, double t) {
    double u = 1 - t;
    double uu = u * u;
    double uuu = uu * u;
    double tt = t * t;
    double ttt = tt * t;
    Point p;
    p.x = uuu * curve.p0.x + 3 * uu * t * curve.p1.x + 3 * u * tt * curve.p2.x + ttt * curve.p3.x;
    p.y = uuu * curve.p0.y + 3 * uu * t * curve.p1.y + 3 * u * tt * curve.p2.y + ttt * curve.p3.y;
    return p;
}

Point bezierDerivative(BezierCurve curve, double t) {
    double u = 1 - t;
    double uu = u * u;
    double tt = t * t;
    Point p;
    p.x = 3 * uu * (curve.p1.x - curve.p0.x) + 6 * u * t * (curve.p2.x - curve.p1.x) + 3 * tt * (curve.p3.x - curve.p2.x);
    p.y = 3 * uu * (curve.p1.y - curve.p0.y) + 6 * u * t * (curve.p2.y - curve.p1.y) + 3 * tt * (curve.p3.y - curve.p2.y);
    return p;
}

Point bezierSecondDerivative(BezierCurve curve, double t) {
    double u = 1 - t;
    Point p;
    p.x = 6 * u * (curve.p2.x - 2 * curve.p1.x + curve.p0.x) + 6 * t * (curve.p3.x - 2 * curve.p2.x + curve.p1.x);
    p.y = 6 * u * (curve.p2.y - 2 * curve.p1.y + curve.p0.y) + 6 * t * (curve.p3.y - 2 * curve.p2.y + curve.p1.y);
    return p;
}

double distance(Point a, Point b) {
    double dx = a.x - b.x;
    double dy = a.y - b.y;
    return sqrt(dx * dx + dy * dy);
}

Point findIntersection(BezierCurve curve1, BezierCurve curve2) {
    double t1 = 0.5; // 初始化迭代参数
    double t2 = 0.5;
    double dt = 1.0;
    while (dt > EPSILON) { // 迭代求解方程组
        Point p1 = bezier(curve1, t1);
        Point p2 = bezier(curve2, t2);
        double d = distance(p1, p2);
        if (d < EPSILON) { // 当两条曲线已经相交时，直接返回交点
            return p1;
        }
        Point v1 = bezierDerivative(curve1, t1);
        Point v2 = bezierDerivative(curve2, t2);
        Point a1 = bezierSecondDerivative(curve1, t1);
        Point a2 = bezierSecondDerivative(curve2, t2);
        double f1 = p1.x - p2.x;
        double f2 = p1.y - p2.y;
        double f3 = v1.x * a2.y - v1.y * a2.x - v2.x * a1.y + v2.y * a1.x;
        double f4 = a1.x * a2.y - a1.y * a2.x;
        double d1 = v1.x * a2.y - v1.y * a2.x;
        double d2 = v2.x * a1.y - v2.y * a1.x;
        double d3 = a1.x * a2.y - a1.y * a2.x;
        double h1 = f1 - d1 / d3 * f3;
        double h2 = f2 - d2 / d3 * f3;
        double h3 = f4;
        double det = h1 * h2 - h3 * h3;
        double t1_new = t1 - (h2 * f3 - h3 * f2) / det;
        double t2_new = t2 - (h1 * f3 - h3 * f1) / det;
        if (t1_new < 0 || t1_new > 1 || t2_new < 0 || t2_new > 1) { // 当迭代参数越界时，重新初始化
            t1_new = 0.5;
            t2_new = 0.5;
        }
        dt = max(abs(t1_new - t1), abs(t2_new - t2));
        t1 = t1_new;
        t2 = t2_new;
    }
    return bezier(curve1, t1);
}

int main() {
    BezierCurve curve1 = {{0, 0}, {1, 2}, {3, -1}, {4, 0}};
    BezierCurve curve2 = {{2, 0}, {3, 1}, {1, 3}, {0, 2}};
    Point p = findIntersection(curve1, curve2);
    cout << "Intersection point: (" << p.x << ", " << p.y << ")" << endl;
    return 0;
}

在上述代码中，`BezierCurve` 结构体表示一条贝塞尔曲线，包含四个控制点。`bezier` 函数用来计算曲线在参数 $t$ 处的点，`bezierDerivative` 和 `bezierSecondDerivative` 函数分别用来计算曲线在参数 $t$ 处的一阶导向量和二阶导向量。`distance` 函数用来计算两个点之间的距离。`findIntersection` 函数则用来求解两条贝塞尔曲线的交点，采用牛顿迭代法求解方程组的根。