于红英，赵彦微，许栋铭

(哈尔滨工业大学 机电工程学院，150001哈尔滨)

轨迹综合是平面连杆机构综合的3类基本问题之一，其主要目的是找到一个连杆机构，使其连杆平面上某一点能依次通过某一预先给定的轨迹.平面铰链四杆机构的轨迹综合一般有3种方法:图解法、解析法和图谱法.

传统的图解法存在作图误差大，综合精度低等缺点，目前已经很少使用.解析法通常与各种优化方法结合运用，通过寻求最小目标函数来获得最优解［1-4］，求解精度较高，但通常要受到初值选择以及寻优方法的影响而难于收敛.随着计算机技术的发展，计算机的海量存储能力和快速检索能力使得图谱法在轨迹综合中应用得越来越广泛.文献［5］通过建立已知机构的连杆转角曲线数据库，将要实现的封闭轨迹曲线也转化成连杆转角曲线，再与数据库中已有转角曲线比较，进而获得满足要求的机构.文献［6-7］中提出利用神经网络对平面连杆机构进行轨迹综合，通过运动学仿真建立了一个大样本库，然后对神经网络进行训练，可以通过训练好的神经网络找到综合问题的近似解.文献［8-9］中采用小波分析方法提取连杆轨迹的特征参数，构建连杆曲线图谱库进行轨迹综合.文献［10］中采用均匀B样条曲线来拟合连杆轨迹，提取B样条曲线控制顶点的坐标作为轨迹特征参数进行轨迹综合.文献［11］中将连杆曲线进行快速傅里叶变换，提取傅里叶级数作为描述连杆曲线的谐波特征参数，然后生成连杆曲线的电子图谱库，再从库中检索出满足要求的机构.上述各种方法中有些方法图谱库中存储数据较多，有些方法需要进行复杂后处理才能得到机构的实际尺寸和安装尺寸，导致不易用编程的方法快速找到合适的机构.

本文以曲柄摇杆机构为例，采用三次非均匀B样条曲线来拟合平面连杆曲线，通过程序自动获取连杆曲线的型值点进行B样条曲线拟合，将B样条曲线的控制多边形各相邻两边间的夹角作为表征连杆曲线形状特征的参数，连同对应的机构尺寸参数建立连杆曲线的电子图谱库，再利用人工神经网络进行特征参数匹配，快速获得合适的平面铰链四杆机构，最终达到轨迹综合的目的.

1 原始连杆曲线的生成

如图1为曲柄摇杆机构ABCD，其中假设AB为曲柄，各杆实际杆长为l1、l2、l3、l4和l5，BP与连杆BC沿逆时针方向所成的角为δ(0°≤δ＜360°)．

图1 曲柄摇杆机构的参数

为研究方便，本文采用量纲一的相对杆长，即取曲柄r1=l1/l1=1，r2=l2/l1，r3=l3/l1，r4=l4/l1，r5=l5/l1．众所周知，只有r1、r2、r3及r4满足一定的长度关系时，才能构成以r1为曲柄的曲柄摇杆机构．现以r2=a(a≥1)为例说明如何确定r3、r4的取值范围．如图2所示，建立以相对杆长r2、r3、r4为坐标轴的三维直角坐标系Or2r3r4．过r2轴上的一点(a，0，0)作一平行于Or3r4的平面M，在平面M内过点(a，1，a)和点(a，a，1)作两条与r3轴正向成45°的射线EE'和FF'，连接E、F两点．可以证明，线段EF与两射线EE'、FF'所围成的条形区域为r3、r4的取值范围．如图2中G点的坐标为(a，b，c)，则代表r1=1，r2=a，r3=b，r4=c的曲柄摇杆机构尺寸型．当r2变化时，r3、r4在条形区域内的值应满足如下关系:

图2 曲柄摇杆机构杆长尺寸分布

本文按照一定的步长变化r2、r3、r4的值，编写程序生成了一定数量的曲柄摇杆机构的尺寸型．对于每一个尺寸型的曲柄摇杆机构，随连杆点P位置的不同(即r5和δ发生变化)产生的连杆曲线也不同．文献［12］在连杆平面上取了5行10列共50个点作为连杆点P的50个位置(见图3)．

图3 连杆点位置分布

本文为丰富图谱库中连杆曲线，选取了11行21列共231个点作为P的231个位置，其坐标值可表示为

根据曲柄摇杆机构的尺寸参数变化，利用计算机可批量完成连杆曲线的数据采集，即每一个尺寸型的曲柄摇杆机构的每一个连杆点在曲柄旋转一周时所形成的连杆曲线上N个点Qi(i=1，2，3，…，N)的坐标．通过试验，每条连杆曲线上取N=720时即可代表该曲线．

2 连杆曲线的归一化

为便于连杆曲线特征匹配，需要对连杆曲线进行归一化处理，使所有曲线的长、短和放置方式统一，无论待识别的曲线是正放、竖放、斜放，或大或小都能被识别出来.归一化具体方法如下.

1)曲柄摇杆机构A0B0C0D0P0的起始位置和曲柄沿逆时针旋转时生成的连杆曲线如图4(a)所示，找出连杆曲线上相距最远的两个点G0、H0．

2)将机构A0B0C0D0P0和连杆曲线平移，使偏左的点G0与坐标原点重合，平移变换公式如式(1)所示:

式中T1为平移变换矩阵，i代表归一化过程中需要进行坐标变换的各点的标号.变换后机构和连杆曲线如图4(b)所示.

图4 连杆曲线归一化过程

3)如图4(b)所示，计算G01H01与x轴正向逆时针所成的角α0，将机构A01B01C01D01P01和连杆曲线顺时针旋转α0使G01H01与x轴重合，旋转变换公式如式(2)所示，变换后机构和连杆曲线如图4(c)所示.

式中R1为旋转变换矩阵．

4)如图4(c)所示，为使曲线凸起部分靠近y轴，找到当前曲线上距x轴最远的点Ⅰ02，Ⅰ02到x轴的垂足为点J02．计算G02J02和J02H02的长度，若，则将机构A02B02C02D02P02连同连杆曲线绕y轴做反射变换，再平移，使H02与坐标原点重合，变换公式如式(3)所示，变换结果如图4(d)所示;若，则保持原图，不进行变换，但交换G02、H02两点代号，即使H02点位于坐标原点．

式中:S1为反射变换矩阵，T2为平移变换矩阵，且当xⅠ02＞xH02-xⅠ02时，

当xⅠ02≤xH02-xⅠ02时，

5)如图4(d)所示，为保证曲线凸起部分在x轴的上方较大，需判断Ⅰ03点纵坐标的正负，然后进行式(4)的变换．若yⅠ03＜0，则绕x轴做反射变换;若yⅠ03≥0，则不变化．变换结果如图4(e)所示.

式中S2为反射变换矩阵，且当yⅠ03＜0时，

当yⅠ03≥0时，

6)如图4(e)所示，将机构A04B04C04D04P04连同连杆曲线缩小倍，比例变换公式如式(5)所示，变换结果如图4(f)所示.

式中S为比例变换矩阵．

机构A0B0C0D0P0连同连杆曲线经过归一化处理，得到归一化后的机构A'0B'0C'0D'0P'0和归一化后的连杆曲线以及在归一化过程中相关的6个归一化参数:xG0、yG0、α0、xⅠ02、yⅠ03、k0．归一化变换矩阵是以上6个变换矩阵的乘积，即

3 连杆曲线型值点的选择

在使用B样条曲线对连杆曲线进行拟合时，连杆曲线上型值点的选择对曲线的特征参数提取及其最终的拟合精度具有重大影响.由于本文利用神经网络进行曲线特征参数的匹配，而神经网络的输入特征参数的数目必须是固定的，连杆曲线又形态各异，初始型值点数目可能并不相等，因此必须设计一套型值点自动选择算法，从连杆曲线的720个点中找到m个能够描述曲线特征的点qi(i=0，1，…，m-1)．具体过程如下:

1)通过数值构造的方式计算出连杆曲线上各点的一阶导数、二阶导数、曲率，选取轨迹上的极值点、拐点以及曲率极值点作为初始型值点，如图5(a)中“+”所示.

2)删去步骤1)中较为密集的型值点，遵循删去曲率较小点的原则，剩余的型值点如图5(b)中“o”所示.

3)步骤2)中剩余型值点数为s，设最终有m+t个特征参数作为神经网络的输入，其中t为连杆曲线精确拟合时再次插入的型值点数．为了保证神经网络输入特征参数数目的一致性，应使初始型值点数由s变到m．若s＞m，则将另外的s-m个型值点依次插入到两个距离最远的型值点中;若s＜m，则在型值点密集处删去m-s个点．经过反复试验，m=12时在型值点数目较小的情况下拟合效果较好，如图5(c)中“◁”所示.

4)对初始型值点进行调整，步骤3)中获得的初始型值点并不能完全满足拟合精度的要求.因此必须要对初始型值点进行调整.已知连杆曲线的N个点Qi(i=1，2，3，…，N)，结合下文中的B样条曲线拟合方法，由m个初始型值点得到初始拟合曲线上的N个点Wi(i=1，2，3，…，N)．为了便于度量拟合质量，定义最大拟合误差esmax=max(‖Qi-Wi‖)(i=1，2，3，…，N)、平均拟合误差，其中‖Qi-Wi‖表示点Qi到点Wi的距离．最大拟合误差反应了拟合曲线和连杆曲线最大分离程度，而平均拟合误差则反应了整体的拟合效果．规定esmax＜0.010和es＜0.008时满足初始拟合精度要求．分别检查esmax和es，若满足要求则可以不再调整型值点，若达不到要求，则找到连杆曲线上出现最大拟合误差点的曲线段的两个型值点，按照一定步长调整左右型值点位置，重新拟合曲线．重复以上步骤，保证型值点数目不变情况下调整型值点的位置，使拟合误差满足预期的精度要求，得到m个型值点qi(i=0，1，…m-1)．最终得到的型值点如图5(d)中“*”所示．

图5 初始型值点选择的过程

4 曲线特征提取与图谱库建立

给定连杆轨迹曲线(通常都是封闭的曲线)上的一组首尾重合的数据点，可构造出一条通过该组数据点且满足误差要求的Ck-1连续(具有k-1阶参数连续性)的k次封闭的B样条插值曲线［13-15］．本文采用具有二阶参数连续性的三次非均匀B样条闭合的插值曲线对连杆曲线上的数据点进行拟合，将提取的型值点qi(i=0，1，…，m-1)作为已知的数据点，拟合的B样条曲线由n+1个控制顶 点di(i=0，1，…，n)与 节 点 矢 量U=［u0，u1，…，un+k+1］来定义．

首先，要确定节点矢量U．本文对于封闭的三次非 均 匀B样 条曲线定义域［uk，un+1］=［u3，un+1］=［0，1］内的节点采用规范积累弦长参数化法来确定，该方法能够如实反映数据点按弦长的分布情况，定义域外节点确定为u0=un-2-1，u1=un-1-1，u2=un-1，un+2=1+u4，un+3=1+u5，un+4=1+u6．

对于闭合的连杆曲线，B样条曲线的首末端点连接处同样要求具有C2连续性，令首末3个控制顶点依次相重合d0=dn-2，d1=dn-1，d2=dn，故未知控制顶点共有n-2个．对于插值曲线，反算过程中应使曲线的分段连接点依次与B样条曲线定义域内节点一一对应起来，则

用于拟合m+1个数据点qi(i=0，1，…，m)的三次B样条插值曲线方程为

将曲线定义域内的节点值代入式(6)应满足插值条件，即

式(7)中含有m+1=n-1个方程，对于C2连续的三次非均匀B样条闭合曲线，由于首末数据点相重合q0=qm，则方程数目减少一个，剩m=n-2个．至此可由m=n-2个方程构成的方程组解出n-2个未知控制顶点．根据以上条件整理可得到以下矩阵方程:

其中

解矩阵方程即可求出未知的控制顶点，至此可得到全部控制顶点.由反算出的控制顶点通过德布尔算法［13］即可求得B样条曲线上的点.

已知连杆曲线的m个初始型值点，首先利用上述方法得到初始的拟合曲线上的N个点Wi(i=1，2，3，…，N)，然后对连杆曲线进行精确拟合．规定esmax＜0.008和es＜0.005时满足精确拟合精度要求．具体方法为:在连杆曲线上找到与拟合曲线的误差最大点处插入一个型值点，更新原有的型值点并重新进行曲线拟合，如此反复进行t次，经过多次的试验，取t=4已能满足曲线拟合的误差要求，不同形状原始曲线的拟合情况见图6，最大拟合误差、平均拟合误差如表1所示，可以看出拟合结果比较理想.

图6 不同形状的曲线拟合

表1 原始曲线的拟合误差

为使神经网络输入特征参数的数目尽可能少，本文选择将控制顶点顺序连接成的控制多边形各相邻两边间的夹角θi(i=1，2，…，m+t)作为连杆曲线的特征参数．

至此，连杆轨迹曲线的特征参数提取完成．将大量的连杆曲线的特征参数θi(i=1，2，…，m+t)和对应机构的6个尺寸参数r1、r2、r3、r4、r5、δ 存储，形成连杆曲线的电子图谱库．

5 寻找满足给定轨迹要求的机构

本文选用广义回归神经网络，利用已建立的电子图谱库中数据作为神经网络训练样本，即m+t=16个控制多边形各相邻两边间夹角θi(i=1，2，…，m+t)作为输入，对应6个机构尺寸参数r1、r2、r3、r4、r5、δ作为输出，对神经网络进行训练并保存．

给定f个有序离散点(离散点数目f不受限制，可在4≤f≤N范围内任意给定)，应用本文所述方法进行轨迹综合，寻找满足要求机构的步骤如下:

1)将给定的f个离散点进行B样条曲线拟合，得到一条封闭的期望曲线．

2)利用与图谱库中连杆曲线相同的处理方法，对期望曲线进行归一化处理，记录此过程中期望曲线的6个归一化参数:xG1、yG1、α1、xⅠ12、yⅠ13、k1，并对归一化后期望曲线进行型值点的选择以及特征参数的提取．

3)将期望曲线特征参数θ1，i(i=1，2，…，m+t)输入到已训练好并保存的神经网络中，输出与期望曲线对应的量纲一的机构尺寸参数(r2，1、r2，2、r2，3、r2，4、r2，5、δ2)，对应的机构为A2B2C2D2P2．

4)由于神经网络具有自学习能力，输出的满足给定轨迹要求的机构尺寸参数可能并不对应图谱库中任何一组数据．因此，需要利用与图谱库中连杆曲线相同的处理方法，将输出机构A2B2C2D2P2及其描绘的连杆曲线也进行一次归一化处理，得到归一化后机构A'2B'2C'2D'2P'2．

5)将归一化后的期望曲线连同步骤4)中的机构A'2B'2C'2D'2P'2，按照步骤2)中记录的6个归一化参数进行反归一化，得到满足要求的实际机构A1B1C1D1P1及其实际参数(l1，1、l1，2、l1、3、l1，4、

l1，5、δ1)．

6)计算期望曲线与实际机构描绘的连杆曲线的综合误差，包括误差平均值es和误差最大值esmax．并计算误差与期望曲线上最远两点间距离k1的比值，即Ks=es/k1和Ksmax=esmax/k1，将Ks、Ksmax作为度量在有量纲情况下轨迹综合的精度参数．

反归一化过程与归一化过程相反，对归一化后的机构和期望曲线进行反归一化的具体方法如下:

a)将归一化后的机构A'2B'2C'2D'2P'2和期望曲线放大k1倍，比例变换公式如式(8)所示，i代表反归一化过程中需要进行坐标变换的各点的标号．

式中S1为比例变换矩阵．

b)判断yⅠ13的正负性，若yⅠ13＜0，即期望曲线在归一化时绕x轴进行了反射变换，则要对当前图形做反射变换，反射变换公式如式(9)所示;若yⅠ13≥0，则不变化．

式中S'1为反射变换矩阵，且当yⅠ13＜0时，

当yⅠ13≥0时，

c)判断xⅠ12和k1-xⅠ12的大小，若xⅠ12＞k1-xⅠ12，即期望曲线归一化时先绕y轴进行了反射和平移变换，则将当前图形绕y轴先做反射变换再平移，变换公式如式(10)所示;若xⅠ12≤k1-xⅠ12，则不变化．

式中:S'2为反射变换矩阵，T'1为平移变换矩阵．且当xⅠ12＞k1-xⅠ12时，

d)将当前图形逆时针旋转α1角，旋转变换公式如式(11)所示:

式中R'1为旋转变换矩阵．

e)将当前图形平移，平移变换公式如式为

式中T'2为平移变换矩阵，得到实际机构A1B1C1D1P1．

反归一化变换矩阵是以上6个变换矩阵的乘积，即

6 算法实现与示例

下面给出计算实例，给定16个离散点的坐标值(见表2)进行轨迹综合，寻找近似通过这些离散点的曲柄摇杆机构.

表2 给定的离散点

应用本文所述的方法，首先，对给定的离散点通过B样条曲线拟合得到期望曲线，如图7(a).然后，对期望曲线经过归一化处理，得到归一化后的期望曲线如图7(b)，以及6个归一化参数(见表3).

图7 综合实例

表3 期望曲线的归一化参数

然后，进行型值点选择和特征参数提取，获取的型值点和控制多边形如图7(c)和图7(d)所示，得到期望曲线的特征参数 θ1，i(i=1，2，…，16)(见表4)．

表4 期望曲线的特征参数

将特征参数输入到已保存的神经网络中，输出与期望曲线对应的机构尺寸参数r2，1=1．000，r2，2=3.037，r2，3=3.461，r2，3=2.500，r2，5=1.108，δ2=57.009°，对应的机构为A2B2C2D2P2．将机构A2B2C2D2P2及其描绘的连杆曲线进行一次归一化，得到归一化后的机构A'2B'2C'2D'2P'2，其各点的坐标为A'2(0.278，-0.275)，B'2(0.374，-0.641)，C'2(1.470，-0.290)，D'2(0.517，-1.191)，P'2(0.484，-0.236)，如图7(e)所示.

最后，将机构A'2B'2C'2D'2P'2和归一化后的期望曲线，按照期望曲线的6个归一化参数xG1、yG1、α1、xⅠ12、yⅠ13和k1进行反归一化，得到的轨迹综合结果如图7(f)所示．实际机构A1B1C1D1P1尺寸参数l1，1=22.729 mm，l1，2=69.036 mm，l1，3=78.671 mm，l1，4=56.824mm，l1，5=25.188 mm，δ1=57.009°，其各点的坐标为A1(22.723，-5.826)，B1(38.742，-21.950)，C1(84.949，29.342)，D1(62.771，-46.138)，P1(32.230，2.382)．综合得到的机构连杆曲线与期望曲线的误差和综合精度参数分别为:误差平均值es=0.657 mm，误差最大值esmax=1.084 mm，Ks=1.10%，Ksmax=1.81%，能够满足轨迹综合的精度要求．为进一步提高轨迹综合的精度可以通过增加轨迹图谱的数量、提高特征参数的数目等方式得到，在此不再详述．

在主频为2.54 GHz、奔腾处理器及内存为2 GB的计算机上利用本文所述方法综合10条期望轨迹，只需要1.35 s就能给出与期望轨迹相对应的实际的曲柄摇杆机构，可见本方法具有非常快的综合速度.

7 结 论

1)给出了基于三次非均匀B样条曲线理论的连杆曲线特征提取方法，利用B样条曲线控制多边形各相邻两边间的夹角作为连杆曲线的特征参数建立电子图谱库，减少了人工神经网络输入数据的数目.

2)通过连杆曲线归一化方法，使重点集中在轨迹曲线的形状上，忽略曲线的具体姿态，消除了由于坐标轴的选择对轨迹综合的影响，大大减少了电子图谱库中数据的冗余量.

3)通过程序自动获取连杆曲线的型值点，避免了对大量形状各异的连杆曲线的分类处理，大大简化了轨迹综合的过程.

4)本方法并不限制给定有序离散点的数目，将离散点拟合的期望曲线的特征参数输入神经网络进行特征匹配，试验表明了该方法的有效性.

5)本文所提方法也适用于其他类型的平面四杆机构的轨迹综合，如双曲柄机构、双摇杆机构和曲柄滑块机构.

［1］CHANEKAR P V，GHOSAL A.Optimal synthesis of adjustable planar four-bar crank-rocker type mechanisms for approximate multi-path generation［J］.Mechanism and Machine Theory，2013，69:263-277.

［2］TARI H，SU H.A complex solution framework for the kinetostatic synthesis of a compliant four-bar mechanism［J］.Mechanism and Machine Theory，2011，46(8):1137-1152.

［3］MATEKAR S B，GOGATE G R.Optimum synthesis of path generating four-bar mechanisms using differential evolution and a modified error function［J］.Mechanism and Machine Theory，2012，52:158-179.

［4］CABRERA J A，ORTIZ A，NADAL F，et al.An evolutionary algorithm for path synthesis of mechanisms［J］.Mechanism and Machine Theory，2011，46(2):127-141.

［5］HONGYING Y，DEWEI T，ZHIXING W.Study on a new computer path synthesis method of a four-bar linkage［J］.Mechanism and machine theory，2007，42(4):383-392.

［6］VASILIU A，YANNOU B.Dimensional synthesis of planar mechanisms using neural networks:application to path generator linkages［J］.Mechanism and Machine Theory，2001，36(2):299-310.

［7］ULLAH I，KOTA S.Optimal synthesis of mechanisms for path generation using Fourier descriptors and global search methods［J］.Journal of Mechanical Design，1997，119(4):504-510.

［8］王成志，纪跃波.小波分析在平面四杆机构轨迹综合中的应用研究［J］.机械工程学报，2004，40(8):34-39.

［9］GALÁN-MARÍN G，ALONSO F J，DEL-CASTILLO J M.Shape optimization for path synthesis of crank-rocker mechanisms using a wavelet-based neural network［J］.Mechanism and Machine Theory，2009，44(6):1132-1143.

［10］张新歌.基于神经网络的铰链四杆机构复演轨迹设计软件的开发［D］.苏州:苏州大学，2008.

［11］褚金奎，孙建伟.基于傅里叶级数理论的连杆机构轨迹综合方法［J］.机械工程学报，2010，46(13):31-41.

［12］郎内斯J A，纳耳桑G L.四连杆机构分析图谱［M］.李学英，译.北京:机械工业出版社，1966.

［13］PIEGL L A，TILLER W.The NURBS book［M］.New York:Springer，1997.

［14］PIEGL L A，TILLER W.Surface approximation to scanned data［J］.The visual computer，2000，16(7):386-395.

［15］PIEGL L A，TILLER W.Biarc approximation of NURBS curves［J］.Computer-Aided Design，2002，34(11):807-814.