一.三阶 Bézier 曲线定义
三阶 Bézier 曲线由四个控制点 P₀, P₁, P₂, P₃ 定义，其参数方程为：
B(t) = (1 - t)³P₀ + 3(1 - t)²tP₁ + 3(1 - t)t²P₂ + t³P₃, t ∈ [0,1]
二.输入条件
设：
- 起点：P₀ = (x₀, y₀)
- 终点：P₃ = (x₃, y₃)
- 起点切线方向（单位向量）：T₀ = (cos θ₀, sin θ₀)
- 终点切线方向（单位向量）：T₃ = (cos θ₃, sin θ₃)
- 起点曲率：κ₀
- 终点曲率：κ₃
三.Bezier 曲线的导数与几何条件
Bézier 曲线的一阶导数为：
B′(t) = 3(1 - t)²(P₁ - P₀) + 6(1 - t)t(P₂ - P₁) + 3t²(P₃ - P₂)
在端点处：
- B′(0) = 3(P₁ - P₀) ⇒ P₁ = P₀ + (B′(0))/3
- B′(1) = 3(P₃ - P₂) ⇒ P₂ = P₃ - (B′(1))/3
设起点和终点的切向速度大小分别为 s₀ 和 s₃，则：
- B′(0) = s₀ T₀
- B′(1) = s₃ T₃
因此： P₁ = P₀ + (s₀/3) T₀
P₂ = P₃ - (s₃/3) T₃
四.利用曲率确定参数 s₀ 和 s₃
平面曲线在参数 t 处的曲率公式为：
κ(t) = |B′(t) × B″(t)| / |B′(t)|³
其中二维叉积：a × b = aₓbᵧ - aᵧbₓ
Bezier 曲线的二阶导数为：
B″(t) = 6(1 - t)(P₂ - 2P₁ + P₀) + 6t(P₃ - 2P₂ + P₁)
在 t = 0 处：
- B′(0) = s₀ T₀
- B″(0) = 6(P₂ - 2P₁ + P₀) = 6[(P₃ - (s₃/3)T₃) - 2(P₀ + (s₀/3)T₀) + P₀]
= 6[P₃ - P₀ - (2s₀/3)T₀ - (s₃/3)T₃]
但为简化，我们使用曲率在端点的近似表达式。
对于三阶 Bezier 曲线，起点曲率 κ₀ 的近似公式为：
κ₀ ≈ |T₀ × (P₂ - 2P₁ + P₀)| / (|P₁ - P₀|²)
但更实用的方法是使用以下曲率与参数速度的关系：
在 t = 0 处，有：
κ₀ = |B′(0) × B″(0)| / |B′(0)|³
计算：
- B′(0) = s₀ T₀
- B″(0) = 6(P₁ - 2P₀ + P₂)，但 P₁ 和 P₂ 依赖于 s₀, s₃
代入 P₁ 和 P₂ 的表达式：
令 ΔP = P₃ - P₀
令 V₀ = T₀, V₃ = T₃
则： P₁ = P₀ + (s₀/3) V₀
P₂ = P₃ - (s₃/3) V₃
代入 B″(0)： B″(0) = 6[(P₁ - 2P₀ + P₂)]
= 6[(P₀ + (s₀/3)V₀) - 2P₀ + (P₃ - (s₃/3)V₃)]
= 6[-P₀ + P₃ + (s₀/3)V₀ - (s₃/3)V₃]
= 6[ΔP + (s₀/3)V₀ - (s₃/3)V₃]
则： B′(0) × B″(0) = (s₀ V₀) × {6[ΔP + (s₀/3)V₀ - (s₃/3)V₃]}
= 6s₀ [V₀ × ΔP + (s₀/3)(V₀ × V₀) - (s₃/3)(V₀ × V₃)]
= 6s₀ [V₀ × ΔP - (s₃/3)(V₀ × V₃)] （因 V₀ × V₀ = 0）
|B′(0)|³ = (s₀)³
因此：
κ₀ = |6s₀ [V₀ × ΔP - (s₃/3)(V₀ × V₃)]| / s₀³
= 6 |V₀ × ΔP - (s₃/3)(V₀ × V₃)| / s₀²
同理，在终点 t = 1：
κ₃ = 6 |V₃ × (ΔP) - (s₀/3)(V₃ × V₀)| / s₃²
注意：V₃ × V₀ = - (V₀ × V₃)
五.求解 s₀ 和 s₃
为简化，常假设 s₀ = s₃ = s（对称情况），或使用启发式方法。
一种常用近似方法是忽略高阶耦合项，采用：
s₀ ≈ √( |2 (T₀ × (P₃ - P₀))| / |κ₀| ) ，若 κ₀ ≠ 0
s₃ ≈ √( |2 (T₃ × (P₀ - P₃))| / |κ₃| ) ，若 κ₃ ≠ 0
更稳健的方法是通过数值求解以下两个方程：
κ₀ = 6 |V₀ × ΔP - (s₃/3)(V₀ × V₃)| / s₀²
κ₃ = 6 |V₃ × (-ΔP) - (s₀/3)(V₃ × V₀)| / s₃²
注意：二维叉积为标量，可用行列式计算： V₀ × V₃ = cosθ₀ sinθ₃ - sinθ₀ cosθ₃ = sin(θ₃ - θ₀)
六.最终控制点计算步骤
- 计算 ΔP = P₃ - P₀
- 计算叉积：
- C₀₃ = T₀ × T₃ = sin(θ₃ - θ₀)
- D₀ = T₀ × ΔP
- D₃ = T₃ × (-ΔP)
- 求解以下方程组（可用数值方法如牛顿法）： s₀² = 6 |D₀ - (s₃/3) C₀₃| / |κ₀|
s₃² = 6 |D₃ + (s₀/3) C₀₃| / |κ₃| - 计算控制点： P₁ = P₀ + (s₀/3) T₀
P₂ = P₃ - (s₃/3) T₃
七.特殊情况处理
- 若 κ₀ = 0，则起点附近近似直线，可设 s₀ 较大或基于切线对齐设定。
- 若切线方向与弦方向接近平行，叉积接近零，需正则化处理。
此方法可用于路径规划、字体设计、机器人轨迹生成等需要 G² 连续（位置、切线、曲率）拼接的场景。