平面G²连续的三次PH过渡曲线的构造

Construction of Cubic PH Transition Curve of Plane G² Continuity

摘要: 本文用三次PH曲线来构造互相包含的两圆之间G2连续的C型过渡曲线。同时证明了互相包含的两圆之间的过渡曲线的唯一性条件。最后，给出过渡曲线的生成步骤，并且通过数值例子证明了该方法的有效性。

Abstract: In this paper, the cubic PH curves are used to construct the G2 continuous C type transition curves between two circles that contain each other. At the same time, the uniqueness condition of the transition curve between two circles containing each other is proved. Finally, the generating steps of the transition curve are given, and the effectiveness of the method is proved by a numerical example.

1. 引言

在CAGD中，通常通过构造G²连续的过渡曲线来实现几何拼接的需要，这要求过渡曲线在端点处保持与原曲线切向量平行，曲率相等，并且过渡曲线最好内部曲率单调或者有尽可能少的极值点，这类曲线被称为Spiral线，在道路轨迹设计、卫星运行轨迹设计以及机器人的行走路线等方面有着广泛用途。

在Spiral线的构造上，很多学者已经有了许多研究成果。传统的几何中采用的是clothoid [1] 曲线(也称cornu回旋线)，这是一类曲率单调的曲线，不仅如此，其曲率值与弧长具有线性函数关系，常用于道路轨迹设计。但是，clothoid曲线的表达式所用的Frensel积分与CAD中的NURBS系统不兼容。鉴于此缺点，Walton等人提出了用三次Bezier曲线 [2] 和五次PH曲线 [3] [4] [5] [6] 来构造满足曲率单调的曲线段。此后，三次Beizer曲线成为CAD及CAGD中最常用的曲线，其次数低且形状易控，并且在计算机中便于实施。用三次Bezier曲线构造的过渡曲线曲率单调，克服了表达形式与CAD系统的不兼容性这一缺点。而五次PH曲线相对于三次Bezier曲线，在构造过渡曲线时曲线的弧长以及等距线的表达式是有理多项式的形式，这一优点不仅为几何设计带来了便利，也解决了以往构造过渡曲线时存在的问题。

由Kneser定理 [7] 可知，由于Spiral曲线的曲率单调性，其上的任何一点的曲率圆一定包含一个较小的曲率圆，同时也被一个较大的曲率圆所包含。因此，Walton等人结合Kneser定理总结了五种最常见的过渡曲线类型：1) 直线与圆之间的J型过渡曲线；2) 圆与圆之间的C型过渡曲线；3) 圆与圆之间的S型过渡曲线；4) 直线与直线之间的过渡曲线；5) 圆与圆之间的C型过渡曲线，其中一圆包含在另一圆内。前四类过渡曲线可以通过构造五次PH过渡曲线，解简单的线性方程得到；而第五类过渡曲线则要通过解复杂非线性方程组的数值解来实现，而又无法预知解的存在性，同时，过渡曲线的次数较高，不便于实际应用。文献 [8] [9] 分别提出用三次PH曲线构造两曲率圆内含和两曲率圆不互相包含情况下的C型G²连续过渡曲线。

本文在五次PH曲线的基础上，进一步研究，提出用较低次的三次PH曲线来构造过渡曲线。讨论了三次PH曲线的曲率的单调性，给出保证曲率单调的充分条件。该方法通过解四次方程，比较便于计算，最终实现用三次PH曲线构造第五类过渡曲线。本文最后给出过渡曲线的生成方法，并且通过数值例子说明了该方法的有效性。

2. 问题背景

Farouki和Sakkali在1990年提出了Pythagorean-hodograph平面曲线，简称PH曲线。PH曲线的主要特征是具有有理性质的等距线及多项式弧长，且其弧长可用含参数的多项式表示。

定义1：对于多项式参数曲线 $P\left(t\right)=\left(x\left(t\right),y\left(t\right)\right)$，若存在多项式 $\sigma \left(t\right)$，使得 ${x^{\prime }}^2\left(t\right)+{y^{\prime }}^2\left(t\right)={\sigma }^2\left(t\right)$，则称 $P\left(t\right)$ 为PH曲线。

对于参数曲线方程 $P\left(t\right)=\left(x\left(t\right),y\left(t\right)\right)$，其切向量为 $P^{\prime }\left(t\right)=\left(x^{\prime }\left(t\right),y^{\prime }\left(t\right)\right)$，则曲线 $P\left(t\right)$ 的曲率为：

$k\left(t\right)=\frac{\left|P^{\prime }\left(t\right)\times P^{″}\left(t\right)\right|}{{\left|P^{\prime }\left(t\right)\right|}^3}$ (1)

其中

$P^{″}\left(t\right)=\left(x^{″}\left(t\right),y^{″}\left(t\right)\right);$

$\left|P^{\prime }\left(t\right)\right|=\sqrt{{\left(x^{\prime }\left(t\right)\right)}^2+{\left(y^{\prime }\left(t\right)\right)}^2}.$

如果参数曲线，它的切向量连续变化，曲率也连续变化，则称该曲线是二阶几何连续，记为G²连续。若该过渡曲线内部曲率值变化是由正到负或者由负到正，则该过渡曲线称为两圆弧间的S型G²连续过渡曲线，如果该过渡曲线内部曲率值不变号，则称为C型G²连续过渡曲线。

定义2 对于三次多项式Bezier曲线

$P\left(t\right)=\left(x\left(t\right),y\left(t\right)\right)={\displaystyle \underset{i=0}{\overset{3}{\sum }}{P}_i{B}_i^3\left(t\right)}$ (2)

其中 $P_i=\left(x_i,y_i\right)$ 为三次PH曲线的控制顶点； $B_i^3\left(t\right)=\left(\begin{array}{c}3\\ i\end{array}\right){\left(1-t\right)}^{3-i}{t}^i\left(0\le t\le 1\right)$ 为Bernstein基函数， $i=0,1,2,3$ 。根据定义1，我们设

$P^{\prime }\left(t\right)=\left(x^{\prime }\left(t\right),y^{\prime }\left(t\right)\right)=(u^2\left(t\right)-v^2\left(t\right),2u\left(t\right)v(t))$

其中 $u\left(t\right),v\left(t\right)$ 满足

$\{\begin{array}{l}v\left(t\right)=u_0\left(1-t\right)+u_{1}t\\ v\left(t\right)=v_0\left(1-t\right)+v_{1}t\end{array}$

根据文献 [9] 中三次PH曲线控制多边形的几何特征，则该三次Bezier控制顶点可表示为：

$\{\begin{array}{l}{P}_1=P_0+\frac{u_0^2-v_0^2,2u_0{v}_0}{3}\\ P_2=P_1+\frac{u_0{u}_1-v_0{v}_1,u_0{v}_1+u_1{v}_0}{3}\\ P_3=P_2+\frac{u_1^2-v_1^2,2u_1{v}_1}{3}\end{array}$ (3)

则，由式(1)可得：

$\{\begin{array}{l}k\left(t\right)=\frac{2\left(u{v}^{\prime }-u^{\prime }v\right)}{{\left(u^2+v^2\right)}^2}\\ k\left(0\right)=\frac{2\left(u_0{v}_1-u_1{v}_0\right)}{{\left(u_0^2+v_0^2\right)}^2}\\ k\left(1\right)=\frac{2\left(u_0{v}_1-u_1{v}_0\right)}{{\left(u_1^2+v_1^2\right)}^2}\\ k^{\prime }\left(t\right)=\left\{2\left[\left(u{v}^{″}-u^{″}v\right)\left(u^2+v^2\right)-4\left(u{v}^{\prime }-u^{\prime }v\right)\left(u{u}^{\prime }+v{v}^{\prime }\right)\right]\right\}/{\left(u^2+v^2\right)}^3\end{array}$ (4)

为了方便计算，取三次PH曲线的起点向量 $P_0=\left(0,0\right)$，则该控制多边形的边向量 $P_1-P_0=\left(1,0\right)$，即 $P\left(t\right)$ 在初始端点的切向量 $P^{\prime }\left(0\right)\left|\right|\left(1,0\right)$，带入式(3)可得：

$\{\begin{array}{l}{P}_0=\left(0,0\right)\\ P_1=P_0+\frac{1}{3}\left(u_0^2,0\right)\\ P_2=P_1+\frac{1}{3}\left(u_0{u}_1,u_0{v}_1\right)\\ P_3=P_2+\frac{1}{3}\left(u_1^2-v_1^2,2u_1{v}_1\right)\end{array}$ (5)

由式(4)可知，曲线在端点 $t=0$ 和 $t=1$ 处的曲率值分别为

$\{\begin{array}{l}k\left(0\right)=\frac{2u_0{v}_1}{u_0^4}\\ k\left(1\right)=\frac{2u_0{v}_1}{{\left(u_1^2+v_1^2\right)}^2}\end{array}$

3. 三次PH曲线构造过渡曲线

首先，由曲线为G²连续可得：

$\{\begin{array}{l}{P}_0=\left(0,0\right)\\ P^{\prime }\left(0\right)\left|\right|\left(1,0\right)\\ \left(k\left(0\right),k\left(1\right)\right)=\left(\frac{1}{r_0},\frac{1}{r_1}\right)\end{array}$ (6)

设三次PH曲线的两端点分别为 $P_0$ 和 $P_{\text{3}}$， $P_0$ 处的曲率圆为 ${\Omega }_0$，其的圆心及半径分别为 $C_0$ 和 $r_0$ ；同样， $P_{\text{3}}$ 处的曲率圆为 ${\Omega }_{\text{1}}$，其圆心及半径分别为 $C_{\text{1}}$ 和 $r_{\text{1}}$ 。三次PH曲线及两曲率圆在坐标系中的设定如图1所示。若设圆心距为r，则有 $\left|C_1-C_0\right|=r$ 。

令

$\lambda =\sqrt[4]{\frac{r_0}{r_1}}$ (7)

不妨设 $r_1>r_0$，则有 $\text{0}<\lambda <1$ 。同时，设三次PH曲线的控制多边形的边向量 $P_1-P_0$ 与 $P_2-P_1$ 的夹角为 $\theta$，且 $1<\theta <\frac{\text{π}}{2}$，则

$\frac{v_1}{u_1}=\tan \theta$ (8)

可得：

$\{\begin{array}{c}{u}_0=\sqrt{2}{r}_0^{\frac{3}{8}}{r}_1^{\frac{1}{8}}\sqrt{\sin \theta }\\ u_1=\sqrt{2}{r}_0^{\frac{1}{8}}{r}_1^{\frac{3}{8}}\cos \theta \sqrt{\sin \theta }\\ v_1=\sqrt{2}{r}_0^{\frac{1}{8}}{r}_1^{\frac{3}{8}}{\left(\sin \theta \right)}^{\frac{3}{2}}\end{array}$ (9)

定理1 [8] ：由图1所示的控制多边形所决定的三次PH曲线能够保证曲率严格单调递增的充分条件为

Figure 1. Coordinate system setting and Endpoint curvature Circle of cubic PH Curve

图1. 三次PH曲线的坐标系设定及端点曲率圆

$0<\theta <\arccos \lambda$ (10)

则式(2)为C型过渡曲线。

定理2：两端点处的曲率圆为 ${\Omega }_{\text{0}}$ 、 ${\Omega }_{\text{1}}$，对应曲率分别为 $k_0$ 、 $k_1$，曲率半径为 $r_0$ 、 $r_1$，令 $\lambda =\sqrt[\text{4}]{\frac{r_0}{r_1}}$ 。两曲率圆的圆心距r满足

$\sigma \left(\lambda \right)\left(r_{\text{1}}-r_{\text{0}}\right)<r<\left(r_{\text{1}}-r_{\text{0}}\right)$ (11)

其中

$\sigma \left(\lambda \right)=\frac{\sqrt{17{\lambda }^4+10{\lambda }^2+9}}{3\left(1+{\lambda }^2\right)}$ (12)

则三次PH过渡曲线是唯一的。

证明：设三次PH曲线的两端点分别为 $P_0$ 和 $P_3$，起始点 $P_0$ 的曲率圆的圆心及曲率半径分别为 ${\Omega }_0$ 和 $r_0$，末端点 $P_3$ 的曲率圆的圆心及曲率半径分别为 ${\Omega }_{\text{1}}$ 和 $r_{\text{1}}$，三次PH曲线的坐标系设定及端点曲率圆如图1所示。

由图1可知， $C_0=\left(0,r_0\right)$， $C_{\text{1}}=\left(x_3-r_1\sin 2\theta ,y_3+r_{\text{1}}\cos 2\theta \right)$，则圆心距为 $C_1-C_0=\left(x_3-r_1\sin 2\theta ,y_3+r_1\cos 2\theta -r_0\right)$ 。

构造函数 $\psi \left(\theta \right)={\left|C_1-C_0\right|}^2-r^2$ 。这里令 $h=\cos \theta$，则有 $\text{0}\le h<1$ 。将 $r_0={\lambda }^4{r}_1$ 带入 $\psi \left(\theta \right)$ 可得关于h的四次方程：

$\begin{array}{l}\psi \left(h\right)=-\frac{16}{9}{\lambda }^4{r}_0^2{h}^4+\frac{16}{9}{\lambda }^3\left(1+{\lambda }^2\right)r_0^2{h}^3+\frac{8}{9}{\lambda }^2\left(-1+{\lambda }^2\right)r_0^2{h}^2\\ -\frac{16}{9}{\lambda }^3\left(1+{\lambda }^2\right)r_0^{2}h+\frac{1}{9}\left(-9r^2+\left(9-8{\lambda }^2+14{\lambda }^4-8{\lambda }^6+9{\lambda }^8\right)r_0^2\right)=0\end{array}$ (13)

由定理1可知， $\lambda <h=\cos \theta <\text{1}$ 。根据零点定理，四次方程 $\psi \left(\theta \right)=\text{0}$ 在 $\left(\lambda ,\text{1}\right)$ 内有根的充分条件为： $\psi \left(1\right)\psi \left(\lambda \right)<0$ 。由于

$\psi \left(\lambda \right)=-r^2+\frac{{\left(-1+{\lambda }^2\right)}^2\left(9+10{\lambda }^2+17{\lambda }^4\right)}{\text{9}}{r}_1^2<0$

则

$\psi \left(1\right)=-r^2+{\left(r_0-r_1\right)}^2>0$

则有 $\psi \left(1\right)>0$，即有 $r<r_0-r_1$ 。因而，只需圆心距r的范围满足 $\sigma \left(\lambda \right)\left(r_{\text{1}}-r_{\text{0}}\right)<r<\left(r_{\text{1}}-r_{\text{0}}\right)$，其中 $\sigma \left(\lambda \right)=\frac{\sqrt{17{\lambda }^4+10{\lambda }^2+9}}{3\left(1+{\lambda }^2\right)}$，且 $0<\sigma \left(\lambda \right)<1$ 。

下面，我们讨论式(13)的根的唯一性。首先对式(13)求导，可得：

${\psi }^{\prime }\left(h\right)=-\frac{16}{9}{\lambda }^3\left(1+{\lambda }^2\right)r_1^2+\frac{16}{9}{\lambda }^2{\left(-1+{\lambda }^2\right)}^2{r}_1^{2}h+\frac{16}{3}{\lambda }^3\left(1+{\lambda }^2\right)r_1^2{h}^2-\frac{64}{9}{\lambda }^4{r}_1^2{h}^3,$

且

${\psi }^{\prime }\left(0\right)=-\frac{16}{9}{\lambda }^3\left(1+{\lambda }^2{r}_1^2\right)<0,{\psi }^{\prime }\left(\lambda \right)=0.$

再对式(13)二次求导，可得：

${\psi }^{″}\left(h\right)=\frac{16}{9}{\lambda }^2{\left(-1+{\lambda }^2\right)}^2{r}_1^2+\frac{32}{3}{\lambda }^3\left(1+{\lambda }^2\right)r_1^{2}h-\frac{64}{3}{\lambda }^4{r}_1^2{h}^2,$

由于式(13)是一个关于h的开口向下的二次函数， ${\psi }^{″}\left(0\right)=\frac{16}{9}{\lambda }^2{\left(-1+{\lambda }^2\right)}^2{r}_1^2>0$，且 ${\psi }^{″}\left(h\right)=\frac{16}{9}{\lambda }^2\left(1-{\lambda }^4+4{\lambda }^2{\left(1-{\lambda }^2\right)}^2\right)>0$，所以 ${\psi }^{\prime }\left(h\right)$ 是一个单调递增的二次函数。又由于 ${\psi }^{\prime }\left(h\right)=0$，当 $h\in \left(0,\lambda \right)$ 时， ${\psi }^{\prime }\left(h\right)<0$，因此 $\psi \left(h\right)$ 在 $\left(0,\lambda \right)$ 上单调递减。因此，式(13)在 $\left(\lambda ,\text{1}\right)$ 内的根是唯一的。

4. 三次PH过渡曲线的生成算法

1) 输入所要生成的过渡曲线的初始端点 $P_0$，根据式(7)的要求，输入半径 $r_0$ 、 $r_{\text{1}}$ 及 $\lambda$ 的值。

2) 根据上一步骤中半径 $r_0$ 、 $r_{\text{1}}$ 的值，得到圆心距r，由式(13)解出 $\theta$ 和h的值。

3) 将 $\theta$ 带入式(9)，得出 $u_0$ 、 $u_{\text{1}}$ 、 $u_{\text{2}}$ 的值。

4) 将上一步得出的 $u_0$ 、 $u_{\text{1}}$ 、 $u_{\text{2}}$ 的值代入式(5)和式(2)，因而输出C型三次PH过渡曲线。

5. 数值例子

例：设起点 $P_0=\left(0,0\right)$ 处的曲率圆圆心为 $C_0=\left(0,1\right)$， $r_{\text{0}}=1$ 。按照三次PH过渡曲线的生成算法，取 $r_1=\text{2}$，则 $\lambda =\text{0}\text{.8409}$， $\theta =\text{0}\text{.3728}$，由此确定该三次PH过渡曲线如图2。

Figure 2. Cubic PH transition curve of G2 continuity between two circles containing each other

图2. 互相包含的两圆之间C型G2连续的三次PH过渡曲线

6. 结束语

本文使用较低次数的三次PH曲线，构造两圆互相包含情况下的过渡曲线。由于过渡曲线次数较低，因此便于计算，且克服了以往方法与CAD系统的不兼容性。同时，本文给出了求两圆互相包含情况下生成PH过渡曲线的方法。