2. 北京航空航天大学 可靠性与环境工程重点实验室, 北京 10008;
3. 中国航空综合技术研究所, 北京 100028;
4. 中国空间技术研究院 总体部, 北京 100029;
5. 北京卫星制造厂, 北京 100094
2. Science and Technology on Reliability and Engineering Laboratory, Beijing University of Aeronautics and Astronautics, Beijing 10008;
3. AVIC China Aero-poly Technology Establishment, Beijing 100028, China;
4. Institute of Spacecraft System Engineering CAST, Beijing 100029, China;
5. Beijing Satellite Manufacturer, Beijing 100094, China
谐波减速器是一种建立在弹性变形理论上的新型机械传动方式,在航空航天、机器人等精密定位领域广泛应用。相比于其他传统的传动方式,谐波传动具有精度高、传动比大、传动平稳和传动效率高等优点。为了提高谐波减速器的传动精度,国内外学者对谐波传动误差进行了深入的研究[1-4],主要从谐波减速器的加工、装配因素等传动机理方面考虑,并没有考虑谐波减速器的柔性和非线性摩擦等动态特性对其传动误差的影响。另外,在对谐波减速器进行静态误差分析时,并没有考虑相关参数不确定性的影响以及静态误差和动态误差之间的耦合关系。Hsia[5]提出了谐波传动误差主要是由波发生器带动柔轮运动时柔轮变形引起的,从设计的角度研究柔轮柔性变形对动态误差的影响,但未考虑其他非线性因素的影响。Tuttle和Seering[6]考虑了谐波减速器柔性以及静态误差等非线性因素,同时深入研究了刚柔轮齿啮合的非线性机理,建立了非线性动力学模型,并通过试验验证了所建模型的正确性。游斌弟和赵阳[7]研究了谐波减速器动态误差在考虑齿轮啮合摩擦和扭转刚度非线性因素时的动力学响应,利用拉格朗日方程建立谐波减速器的误差动力学方程,研究了不同频率谐波强迫激励下动态误差影响。Preissner等[8]建立了综合考虑柔轮非线性扭转、滞回特性和运动误差的谐波传动模型,着重研究了柔轮的滞回特性对动态精度的影响。这些研究大多只考虑谐波减速器的柔轮动力学特性对动态精度的影响,没有考虑动力学模型参数的不确定性。
本文从静态误差和动态误差产生的机理方面着手,综合考虑了加工、装配误差、柔轮柔性以及刚柔轮摩擦综合作用下的谐波减速器的动态精度问题;在此基础上考虑相关参数的不确定性对动态精度的影响,利用多项式混沌展开(Polynomial Chaos Expansion,PCE)方法进行动态精度不确定性分析,构建动态精度可靠度模型计算动态精度可靠度。结果表明,利用PCE方法能够很好地处理谐波减速器非线性以及多变量的动态精度不确定性问题。
1 谐波减速器动态精度分析谐波减速器主要由3部分组成:波发生器、柔轮和刚轮,其结构简图如图 1所示。电机输出一定转速驱动波发生器旋转,波发生器带动柔轮发生柔性变形与刚轮啮合传递运动。其中波发生器是输入端,柔轮是输出端。
谐波减速器的动态误差定义为在一定转速条件下柔轮理想输出位置和实际输出位置之差
(1) |
式中:θm为波发生器输入角位置;θl为柔轮输出角位置;N为谐波传动比。
在一定转速下谐波减速器柔轮会发生柔性变形,当其做正反往复运动时,就会出现如图 2所示的滞回误差。
谐波减速器动态误差由2部分组成[9]:一部分是由谐波减速器基本部件加工、装配引起的静态误差
(2) |
静态误差[10]主要由谐波减速器各个部件的加工、装配误差所产生。其中刚轮加工误差产生的运动误差Δc1为
(3) |
式中:ΔFp2为刚轮周节累积误差;Δf′f2为刚轮切向相邻齿综合误差;ωb为波发生器角频率;t为时间;Z2为刚轮齿数。
柔轮加工误差引起的运动误差Δr1为
(4) |
式中:ΔFp1为柔轮周节累积误差;Δf′f1为柔轮切向相邻齿综合误差;Z1为柔轮齿数。
刚轮、柔轮装配误差,不考虑装配误差初相角的影响,由刚轮安装偏心误差Ec产生的运动误差为
(5) |
式中:αn为啮合角。由柔轮安装偏心误差Ef产生的运动误差为
(6) |
由波发生器安装偏心误差Eb产生的运动误差Δb为
(7) |
综上,可以得到总的静态误差为
(8) |
式中:d为刚轮分度圆直径。
动态精度包括机构的静态误差和柔轮柔性和机构摩擦引起的滞回误差,分析谐波减速器的动态精度,应该在静态误差模型基础之上建立机构的动力学模型。本文研究的谐波减速器物理简化模型如图 3所示,图中:Jm和Jl分别为输入和输出端转动惯量;Bm、Bl和Bsp分别为谐波减速器输入端、输出端阻尼和刚柔啮合处阻尼。
将谐波减速器的柔轮等效为一个非线性扭转弹簧,柔轮扭转力矩特性用滞回误差角
(9) |
式中:Tk为柔轮的非线性扭矩;k1和k2为谐波减速器柔轮等效扭转刚度。由图 3所示的谐波减速器各部件之间的动力学关系,系统的动能T为
(10) |
系统的势能V为
(11) |
系统的瑞利耗散函数为
(12) |
因此,谐波减速器拉格朗日动力学方程[11]为
(13) |
(14) |
(15) |
式中:Tm为输入转矩。
根据式(1)和式(13)~式(15)可以得到谐波减速器动力学模型,动力学模型含有非线性静态误差项
PCE基本理论是用一个属于某个对应分析类型的正交多项式混沌之和(含有一个或多个随机变量)来近似地表示一个随机过程。
假设{ξi}i=1n表示机构的随机参数,根据Wiener[12]提出的高斯随机过程的均匀Hermite多项式展开模型,随机过程输出响应可以表示为
(16) |
c=(c0, c1, …)为待定系数矢量;ξ=[ξ1, ξ2, …, ξn]为服从标准正态分布的随机变量矢量;Πd(ξi1, ξi2, …, ξid)为d次多维Hermite多项式。将式(16)截断用s项来近似精度则可简化[12]为
(17) |
式中:η为随机事件;cj为待求解的确定性系数;Πj(ξ1, ξ2, …, ξn)为广义Wiener-Askey多项式混沌。对于一个n维Hermite多项式,可表示为
(18) |
Hermite多项式的随机变量如果是标准正态分布,则满足
(19) |
式中:〈··〉为希尔伯特空间上的内积,此处定义为
(20) |
式中:多项式基对应的权重函数W(ξ)为
(21) |
如果机构输入的随机变量的个数为n,又机构输出响应量的多项式展开式最高阶次为p,则待定系数的个数P可以用式(22)求得:
(22) |
谐波减速器的动态精度不仅跟机构部件的制造公差、装配间隙等参数有关,而且还受机构动力学参数(如等效刚度、阻尼和转动惯量)的影响。上述参数在制造、测量过程中必然会存在不确定性而不是一个固定的数值,参数的不确定性通过机构的动力学模型影响动态精度响应的不确定性。
设谐波减速器的动力学模型用M表示,机构部件的制造公差、间隙、等效刚度和转动惯量等不确定参数可以用向量ψ表示,即
(23) |
式中:ψi为各个具体的不确定参数,n为不确定参数总的个数。谐波减速器动态精度
(24) |
设ψi~N(μi, σi),通过标准变换
(25) |
(26) |
(27) |
式中:ci, m、cij, m、cijj, m和cijk, m为展开式中的待定系数。
本文采用随机响应面配点法[14]计算多项式混沌的展开系数。在随机向量展成的空间中,每一组样本{ξ1, ξ2,…}都会对应一个点,这些点称为配点,对于Wiener的多项式混沌,其展开式基函数为Hermite多项式,若Hermite多项式的最高阶数为p,则相应的配点通常取p+1阶Hermite多项式的根。未知系数可以通过式(28)和式(29)计算:
(28) |
(29) |
式中:ξ1, ξ2, …, ξk为采样点;k为采样点数。
根据Hermite多项式正交性,动态精度的均值可以通过式(30)求得[15]:
(30) |
由式(30)可以看出,动态精度均值是多项式混沌展开式的0阶项,同样的方法得到动态精度
(31) |
根据谐波减速器动力学模型,本文针对XB1-50型号谐波减速器采用Modelica语言在Dymola[16]编译环境下建立谐波减速器的动力学仿真模型(见图 4), 主要包括直PID控制模块、电机模块和谐波减速器模块。此模型能够较好处理谐波减速器的机、电耦合问题,在求解非线性微分方程时精度较高。
采用如图 5所示的谐波减速器精度测试平台对所建立的非线性谐波减速动力学仿真模型进行验证。采用恒温箱封闭,电机输入端和谐波减速器输出端分别连接有角速度和角度传感器,测定输入、输出端角位置以及角速度。实验条件下,电机输入转速为100 r/min。
根据实测和设计参数数据确定谐波减速器动力学模型各参数如表 1所示,静态误差模型参数如表 2所示。
参数 | 数值 |
Jm/(kg·m2) | 3.2×10-4 |
Jl/(kg·m2) | 8.5×10-4 |
Bm/(N·m·s·rad-1) | 1.7×10-4 |
Bl/(N·m·s·rad-1) | 5.0×10-4 |
Bsp/(N·m·s·rad-1) | 2.8×10-4 |
k1/(N·m·rad-1) | 7 160 |
k2/(N·m·rad-3) | 21 576 |
N | 90 |
参数 | 数值 |
Z2 | 182 |
Z1 | 180 |
αn/(°) | 20 |
Ef/m | 3.1×10-5 |
Ec/m | 1.78×10-4 |
Eb/m | 2.05×10-4 |
ΔFp2/m | 3.6×10-5 |
Δf′f2/m | 1.14×10-4 |
ΔFp1/m | 3.6×10-5 |
Δf′f1/m | 1.14×10-4 |
在Dymola仿真模型中通过PID调节器控制电机输出转速为100 r/min,仿真时间为10 s,待电机输入转速及误差波动曲线稳定后,通过采样输入角在3个周期均匀变化下动态误差值,实验条件下,通过角位移传感测得同样周期内谐波减速器输入端角位移θm以及输出端角位移θl,求得动态误差的测量值。将动态误差仿真值与测量值对比如表 3所示。
电机 输入角/(°) |
动态误差 实验值/(°) |
动态误差 仿真值/(°) |
真实误差/ (°) |
相对 误差/% |
50 | 0.023 473 | 0.021 908 | 0.001 565 | 4.7 |
100 | 0.018 379 | 0.017 472 | 0.000 907 | 4.9 |
200 | 0.016 727 | 0.015 486 | 0.001 241 | 7.4 |
300 | 0.018 503 | 0.017 256 | 0.001 247 | 6.7 |
400 | 0.016 205 | 0.017 369 | 0.001 164 | 7.1 |
500 | 0.019 181 | 0.018 688 | 0.000 493 | 2.5 |
700 | 0.017 800 | 0.019 320 | 0.002 820 | 8.5 |
800 | 0.019 794 | 0.017 913 | 0.001 881 | 9.5 |
1 000 | 0.020 950 | 0.022 243 | 0.001 293 | 4.1 |
通过表 3中数据可以看出,实验得到的动态误差值与仿真值比较接近,相对误差在2.5%~9.5%范围内,考虑到相关参数不确定性,相对误差在接受范围之内,由此可见所建立的谐波减速器Dymola仿真模型能够很好地模拟考虑间隙和柔性非线性因素综合作用下的动态精度。
仿真得到谐波减速器输入、输出转速,动态误差和静态误差曲线如图 6~图 9所示,图 6为谐波减速器输入转速曲线,转速经过短暂波动后稳定后在10.499 rad/s。由于谐波减速器柔性和控制器惯性的作用,图 7所示谐波减速器输出转速先出现大范围波动,然后稳定在0.1 rad/s附近波动。图 8为动态误差随输入角θm波动曲线,开始波动程度较大,稳定后,呈周期性波动。由图 9可以看出,综合考虑静态误差、柔性和摩擦作用的动态误差曲线比只考虑静态误差曲线波动幅值要大0.005°,相比增加25.4%,因此在对谐波减速器进行动态误差分析时,有必要考虑柔性和摩擦的影响。
3.2 动态精度灵敏度及不确定性分析Sobol敏感度分析方法[17]是一种基于方差分解的Monte Carlo方法,Sobol方法考虑了随机输入参量在整个取值范围内对输出响应的贡献以及随机参数的交互作用,Sobol敏感度指标可以直接从PCE式的系数得到。
动态精度PCE式(26)和式(27)可以改写为
(32) |
式中:cα为PCE多项式系数;ψα为PCE多项式的项;
(33) |
(34) |
由式(34)可以看出分解项表征不同随机变量及其相互作用对谐波减速器动态精度输出响应方差的贡献,因此可以定义Sobol敏感度指标为
(35) |
其满足
(36) |
式中:Si为主效应敏感度指标,表征各个随机变量对动态精度响应方差的贡献。
考虑到当不确定参数较多时,PCE式的待定系数呈指数增长,为了简化分析,对谐波减速器动力学模型参数进行灵敏度分析。考虑到输出轴转速相对输出轴转速小很多,输出轴转动惯量相对输入轴转动惯量大很多,故不考虑输出端阻尼和输入端转动惯量不确定参数的影响。静态误差模型中柔轮周节累积误差与柔轮安装偏心向量、刚轮周节累积误差与刚轮安装偏心向量、刚柔轮切向相邻齿综合误差在静态误差模型中形式上相同,为了简化不确定性分析,故只考虑刚轮、柔轮、波发生器安装偏心向量和柔轮切向相邻齿综合误差不确定性参数的影响。取谐波减速器不确定随机向量
参数 | 分布 | 均值 | 标准差 |
Bm | 正态分布 | 1.7×10-4 N·m·s·rad-1 | 0.000 03 |
Jl | 正态分布 | 8.5×10-4kg·m2 | 0.000 05 |
k1 | 正态分布 | 7.16×103 N·m·rad-1 | 200 |
k2 | 正态分布 | 2.157 6×104 N·m·rad-3 | 500 |
Bsp | 正态分布 | 2.8 ×10-4 N·m·s·rad-1 | 0.000 04 |
Ef | 正态分布 | 3.1 ×10-5 m | 0.000 005 |
Ec | 正态分布 | 1.78 ×10-4 m | 0.000 04 |
Eb | 正态分布 | 2.05 ×10-4 m | 0.000 05 |
Δf′f1 | 正态分布 | 1.14×10-4 m | 0.000 04 |
参数 | 多项式混沌展开配点 |
Bm | (0.000 187,0.000 17,0.000 152) |
Jl | (0.000 879,0.000 85,0.000 821) |
k1 | (0.000 721 8,0.000 716 0,0.000 710 2) |
k2 | (0.020 999,0.021 576,0.022 153) |
Bsp | (0.000 303,0.000 28,0.000 257) |
Ef | (0.000 033 8,0.000 031,0.000 028 1) |
Ec | (0.000 201,0.000 178,0.000 154) |
Eb | (0.000 176 3,0.000 205,0.000 233 7) |
Δf′f1 | (0.000 116 9,0.000 14,0.000 163 1) |
在随机向量空间ξ=[ξ1,ξ2,…,ξ9]关于中心点对称进行110次抽样仿真,求得动态精度
项 | 系数 |
1 | 0.012 5 |
ξ1 | 0.032 8 |
ξ2 | -0.154 0 |
ξ3 | 0.357 6 |
ξ4 | -0.009 3 |
ξ5 | 0.000 9 |
ξ6 | 0.014 0 |
ξ7 | 0.074 9 |
ξ8 | 0.007 2 |
ξ9 | 0.102 3 |
根据展开式一次项系数可以得到各随机变量对动态精度
可以得出动态误差主要与谐波减速器的动力学参数有关。在仿真模型中,分别将扭转刚度k1,输出轴转动惯量Jl,柔轮切向相邻齿综合误差Δf′f1的取值增加20%,得到动态误差曲线如图 10~图 12所示,对应的
为了减少参数不确定性分析的计算量,选取主要的不确定性随机参数向量ζ=[Jl,Δf′f1,k1,Ec,Bm],PCE响应为动态精度在t=5~6 s内等间距50个时间点的仿真值,每个时间点对称抽取42个样本点,图 13给出了PCE方法和1 000次Monte Carlo方法所得到的谐波减速器动态误差的均值,可以看出这2种方法所求的结果十分吻合。图 14给出了2种方法的动态误差的均方差,这2种方法求得的均方差之间的误差很小,PCE方法精度接近于Monte Carlo方法的精度,而PCE方法仅通过42次仿真计算就得到了Monte Carlo方法1 000次仿真的计算结果,效率明显优于Monte Carlo方法。
3.3 动态精度可靠性分析本文定义谐波减速器动态精度可靠度为谐波减速器在一定初始条件和时域内,动态误差绝对值的最大值不超过一定的阈值的概率,则动态精度可靠性功能函数可表示为
(37) |
式中:
项 | 系数 |
1 | 0.024 151 |
ζ1 | -0.236 1 |
ζ2 | 0.036 6 |
ζ3 | 0.720 9 |
ζ4 | -0.023 6 |
ζ5 | -0.021 8 |
ζ1-1 | 0.0365 |
ζ2-1 | -0.000 8 |
ζ3-1 | 0.049 0 |
ζ4-1 | -0.000 2 |
ζ5-1 | -0.001 5 |
ζ1ζ2 | 0.003 7 |
ζ1ζ3 | -0.036 5 |
ζ1ζ4 | 0.000 6 |
ζ1ζ5 | -0.005 3 |
ζ2ζ3 | -0.005 4 |
ζ2ζ4 | 0.000 1 |
ζ2ζ5 | 0.001 5 |
ζ3ζ4 | 0.003 9 |
ζ3ζ5 | 0.004 2 |
ζ4ζ5 | -0.000 5 |
设动态误差的阈值为
(38) |
式中:ζi*为设计验算点;
(39) |
通过迭代计算得到动态精度可靠度计算结果如表 8所示。
在Dymola中进行500次Monte Carlo抽样仿真实验,得到动态误差仿真曲线如图 15所示。
根据式Pf=nf/nt求动态精度在一定时域内的失效概率,nf为动态误差绝对值最大值超过阈值的次数,nt为抽样仿真次数,表 9给出了nt=500~10 000次的模拟结果。
nt | nf | Pf |
500 | 14 | 0.028 0 |
2 000 | 61 | 0.030 5 |
4 000 | 126 | 0.031 5 |
6 000 | 205 | 0.034 2 |
7 000 | 254 | 0.036 3 |
8 000 | 305 | 0.038 1 |
10 000 | 381 | 0.038 1 |
从表 9中可以看出, 当n=8 000时,Pf收敛,Monte Carlo仿真实验得到动态精度可靠度为0.961 9,利用PCE方法和FORM法得到的可靠度0.958 4,两者求得的结果相近,而PCE方法仅用到42次抽样仿真,效率明显优于Monte Carlo方法。
4 结论本文考虑了静态误差、柔性和摩擦等因素综合作用下的谐波减速器动态精度问题,建立了谐波减速器非线性动力学模型,利用多项式混沌展开(PCE)方法处理动态精度的随机不确定性。解决了以下问题:
1) 将静态误差与动态误差结合起来,建立了含有静态误差项的非线性动力学模型,通过实验验证了考虑动力学因素和静态误差综合作用下的动态精度更接近实际情况。
2) 基于PCE方法,通过灵敏度分析得到影响动态精度的主要参数是柔轮等效刚度、输出轴转动惯量和柔轮切向相邻齿综合误差。可以看出,柔轮的加工精度对动态精度影响较大。通过灵敏度分析选取影响较大的参数进行不确定性分析,可以减少不确定性分析的计算量。
3) 利用PCE方法得到动态精度的随机统计特性,并将PCE方法与传统的Monte Carlo方法比较,效率更高。利用PCE方法和FORM法计算一定时域内的动态精度可靠度,比传统的Monte Carlo方法计算可靠度效率更高。
[1] | NYE T, KRAML R. Harmonic drive gear error: Characterization and compensation for precision pointing and tracking[C]//Processing of the 25th Aerospace Mechanics, Symposium. Washington, D. C. : NASA, 1991: 237-252. |
[2] | EMELYANOV A F. Calculation of the kinematic error of a harmonic gear transmission taking into account the compliance of elements[J]. Soviet Engineering Research, 1983, 3 (7): 7–10. |
[3] | GRAVAGNO F, MUCINO V H, PENNESSTRI E. Influence of wave generator profile on pure kinematic error and centrodes of harmonic drive[J]. Mechanism and Machine Theory, 2016, 104 : 100–107. DOI:10.1016/j.mechmachtheory.2016.05.005 |
[4] |
沙晓晨, 范元勋. 谐波减速器传动误差的研究[J].
机械制造, 2015, 44 (5): 50–54.
SHA X C, FAN Y X. Study of transmission error of harmonic drive reducer[J]. Machine Building Automation, 2015, 44 (5): 50–54. (in Chinese) |
[5] | HSIA L M. The analysis and design of harmonic gear drives[C]//Processing of the 1988 IEEE International Conference on Systems, Man, and Cybernetics. Piscataway, NJ: IEEE Press, 1988: 616-620. |
[6] | TUTTLE T, SEERING W. Kinematic error, compliance, and friction in a harmonic drive gear transmission[C]//ASME Design Technical Conferences 19th Design Automation. New York: ASME, 1993: 319-324. |
[7] |
游斌弟, 赵阳. 考虑非线性因素的谐波齿轮传动动态误差研究[J].
宇航学报, 2010, 31 (5): 1297–1282.
YOU B D, ZHAO Y. Study on dynamic error of harmonic drive with nonlinear factors[J]. Journal of Astronautics, 2010, 31 (5): 1297–1282. (in Chinese) |
[8] | PREISSNER C, ROYSTON T J, SHU D. A high-fidelity harmonic drive model[J]. Journal of Dynamic Systems, Measurement, and Control, 2012, 134 (1): 011002. DOI:10.1115/1.4005041 |
[9] | KIRCANSKI N, GOLDENBERG A, ANGELS J. Nonlinear modeling and parameter identification of harmonic drive gear transmissions[C]//Processing of the 32nd IEEE Conference on Robotics and Automatic. Piscataway, NJ: IEEE Press, 1995: 3027-3032. |
[10] |
王爱东. 机器人用谐波齿轮传动装置的运动精度分析[D]. 北京: 中国科学院, 2001: 11-15.
WANG A D. Kinematic accuracy analysis of gear transmission for harmonic reducer of robot[D]. Beijing: Chinese Academy of Sciences, 2001: 11-15(in Chinese). |
[11] | FATHI H, PRASANNA S, FRIEDHELM A. On the kinematic error in harmonic drive gears[J]. Journal of Mechanical Design, 2001, 123 (1): 90–97. DOI:10.1115/1.1334379 |
[12] | WIENER N. The homogeneous chaos[J]. American Journal of Mathematics, 1938, 60 (4): 897–936. DOI:10.2307/2371268 |
[13] |
赵珂, 高正红, 黄江涛, 等. 基于PCE方法的翼型不确定性分析及稳健设计[J].
力学学报, 2014, 46 (1): 11–19.
ZHAO K, GAO Z H, HUANG J T, et al. Uncertainty quantification and robust design of airfoil based on polynomial chaos technique[J]. Chinese Journal of Theoretical and Applied Mechanics, 2014, 46 (1): 11–19. (in Chinese) |
[14] | LOEVETT T, PONCI F, MONTI A. A polynomial chaos approach to measurement uncertainty[J]. IEEE Transations on Instrumentation and Measurement, 2006, 55 (3): 729–736. DOI:10.1109/TIM.2006.873807 |
[15] | BENJAMIN L, HOSAM K, JEFFREY L. Recursive maximum likelihood parameter estimation for state space systems using polynomial chaos theory[J]. Automatica, 2011, 47 (11): 2420–2424. DOI:10.1016/j.automatica.2011.08.014 |
[16] |
陶海川, 来新民. 基于Dymola的无刷直流电机仿真模型[J].
计算机仿真, 2005, 22 (5): 63–65.
TAO H C, LAI X M. Computer simulation of brushless DC motor system based on Dymola[J]. Computer Simulation, 2005, 22 (5): 63–65. (in Chinese) |
[17] | SUDRET B. Global sensitivity analysis using polynomial chaos expansion[J]. Reliability Engineering and System Safety, 2008, 93 (7): 964–979. DOI:10.1016/j.ress.2007.04.002 |
[18] | PENG W S, ZHANG J G, ZHU D T. ABCLS methods for high-reliability aerospace mechanism with truncated random uncertainties[J]. Chinese Journal of Aeronautics, 2015, 28 (4): 1066–1075. DOI:10.1016/j.cja.2015.06.012 |