谐波减速器是一种建立在弹性变形理论上的新型机械传动方式，在航空航天、机器人等精密定位领域广泛应用。相比于其他传统的传动方式，谐波传动具有精度高、传动比大、传动平稳和传动效率高等优点。为了提高谐波减速器的传动精度，国内外学者对谐波传动误差进行了深入的研究^[1-4]，主要从谐波减速器的加工、装配因素等传动机理方面考虑，并没有考虑谐波减速器的柔性和非线性摩擦等动态特性对其传动误差的影响。另外，在对谐波减速器进行静态误差分析时，并没有考虑相关参数不确定性的影响以及静态误差和动态误差之间的耦合关系。Hsia^[5]提出了谐波传动误差主要是由波发生器带动柔轮运动时柔轮变形引起的，从设计的角度研究柔轮柔性变形对动态误差的影响，但未考虑其他非线性因素的影响。Tuttle和Seering^[6]考虑了谐波减速器柔性以及静态误差等非线性因素，同时深入研究了刚柔轮齿啮合的非线性机理，建立了非线性动力学模型，并通过试验验证了所建模型的正确性。游斌弟和赵阳^[7]研究了谐波减速器动态误差在考虑齿轮啮合摩擦和扭转刚度非线性因素时的动力学响应，利用拉格朗日方程建立谐波减速器的误差动力学方程，研究了不同频率谐波强迫激励下动态误差影响。Preissner等^[8]建立了综合考虑柔轮非线性扭转、滞回特性和运动误差的谐波传动模型，着重研究了柔轮的滞回特性对动态精度的影响。这些研究大多只考虑谐波减速器的柔轮动力学特性对动态精度的影响，没有考虑动力学模型参数的不确定性。

本文从静态误差和动态误差产生的机理方面着手，综合考虑了加工、装配误差、柔轮柔性以及刚柔轮摩擦综合作用下的谐波减速器的动态精度问题；在此基础上考虑相关参数的不确定性对动态精度的影响，利用多项式混沌展开(Polynomial Chaos Expansion，PCE)方法进行动态精度不确定性分析，构建动态精度可靠度模型计算动态精度可靠度。结果表明，利用PCE方法能够很好地处理谐波减速器非线性以及多变量的动态精度不确定性问题。

1.   谐波减速器动态精度分析

谐波减速器主要由3部分组成：波发生器、柔轮和刚轮，其结构简图如图 1所示。电机输出一定转速驱动波发生器旋转，波发生器带动柔轮发生柔性变形与刚轮啮合传递运动。其中波发生器是输入端，柔轮是输出端。

图  1  谐波减速器结构简图

Figure  1.  Structure diagram of harmonic reducer

下载: 全尺寸图片 幻灯片

谐波减速器的动态误差定义为在一定转速条件下柔轮理想输出位置和实际输出位置之差为

式中：θ_m为波发生器输入角位置；θ_l为柔轮输出角位置；N为谐波传动比。

在一定转速下谐波减速器柔轮会发生柔性变形，当其做正反往复运动时，就会出现如图 2所示的滞回误差。

图  2  谐波减速器往复运动动态误差曲线

Figure  2.  Dynamic error curve of harmonic reducer reciprocation

下载: 全尺寸图片 幻灯片

谐波减速器动态误差由2部分组成^[9]：一部分是由谐波减速器基本部件加工、装配引起的静态误差；另一部分是由柔轮的柔性和摩擦引起的波发生器与柔轮产生相对运动产生的滞回误差。因此动态误差也可以表示为

静态误差^[10]主要由谐波减速器各个部件的加工、装配误差所产生。其中刚轮加工误差产生的运动误差Δ_c1为

式中：ΔF_p2为刚轮周节累积误差；Δf′_f2为刚轮切向相邻齿综合误差；ω_b为波发生器角频率；t为时间；Z₂为刚轮齿数。

柔轮加工误差引起的运动误差Δ_r1为

式中：ΔF_p1为柔轮周节累积误差；Δf′_f1为柔轮切向相邻齿综合误差；Z₁为柔轮齿数。

刚轮、柔轮装配误差，不考虑装配误差初相角的影响，由刚轮安装偏心误差E_c产生的运动误差为

式中：α_n为啮合角。由柔轮安装偏心误差E_f产生的运动误差为

由波发生器安装偏心误差E_b产生的运动误差Δ_b为

综上，可以得到总的静态误差为

式中：d为刚轮分度圆直径。

动态精度包括机构的静态误差和柔轮柔性和机构摩擦引起的滞回误差，分析谐波减速器的动态精度，应该在静态误差模型基础之上建立机构的动力学模型。本文研究的谐波减速器物理简化模型如图 3所示，图中：J_m和J_l分别为输入和输出端转动惯量；B_m、B_l和B_sp分别为谐波减速器输入端、输出端阻尼和刚柔啮合处阻尼。

图  3  谐波减速器物理简化模型

Figure  3.  A physical simplified model of harmonic reducer

下载: 全尺寸图片 幻灯片

将谐波减速器的柔轮等效为一个非线性扭转弹簧，柔轮扭转力矩特性用滞回误差角的三次多项式表示为

式中：T_k为柔轮的非线性扭矩；k₁和k₂为谐波减速器柔轮等效扭转刚度。由图 3所示的谐波减速器各部件之间的动力学关系，系统的动能T为

系统的势能V为

系统的瑞利耗散函数为

因此，谐波减速器拉格朗日动力学方程^[11]为

式中：T_m为输入转矩。

根据式(1)和式(13)~式(15)可以得到谐波减速器动力学模型，动力学模型含有非线性静态误差项和非线性扭转刚度k₂，同时动力学模型含有参数较多，因此谐波减速器动力学模型具有高维非线性的特点，采用数值仿真方法求解。

2.   基于PCE动态精度不确定性模型

2.1   PCE基本理论

PCE基本理论是用一个属于某个对应分析类型的正交多项式混沌之和(含有一个或多个随机变量)来近似地表示一个随机过程。

假设{ξ_i}_i=1ⁿ表示机构的随机参数，根据Wiener^[12]提出的高斯随机过程的均匀Hermite多项式展开模型，随机过程输出响应可以表示为

c=(c₀, c₁, …)为待定系数矢量；ξ=[ξ₁, ξ₂, …, ξ_n]为服从标准正态分布的随机变量矢量；Π_d(ξ_i1, ξ_i2, …, ξ_id)为d次多维Hermite多项式。将式(16)截断用s项来近似精度则可简化^[12]为

式中：η为随机事件；c_j为待求解的确定性系数；Π_j(ξ₁, ξ₂, …, ξ_n)为广义Wiener-Askey多项式混沌。对于一个n维Hermite多项式，可表示为

Hermite多项式的随机变量如果是标准正态分布，则满足

式中：〈··〉为希尔伯特空间上的内积，此处定义为

式中：多项式基对应的权重函数W(ξ)为

如果机构输入的随机变量的个数为n，又机构输出响应量的多项式展开式最高阶次为p，则待定系数的个数P可以用式(22)求得：

2.2   动态精度不确定性模型

谐波减速器的动态精度不仅跟机构部件的制造公差、装配间隙等参数有关，而且还受机构动力学参数(如等效刚度、阻尼和转动惯量)的影响。上述参数在制造、测量过程中必然会存在不确定性而不是一个固定的数值，参数的不确定性通过机构的动力学模型影响动态精度响应的不确定性。

设谐波减速器的动力学模型用M表示，机构部件的制造公差、间隙、等效刚度和转动惯量等不确定参数可以用向量ψ表示，即

式中：ψ_i为各个具体的不确定参数，n为不确定参数总的个数。谐波减速器动态精度在任一时刻t可以表示为向量ψ的确定性函数：

设ψ_i~N(μ_i, σ_i)，通过标准变换，ξ=[ξ₁, ξ₂, …, ξ_n]，则动态精度用一阶、二阶、三阶Hermite多项式展开如式(25)~式(27)所示^[13]:

式中：c_{i, m}、c_{ij, m}、c_{ijj, m}和c_{ijk, m}为展开式中的待定系数。

本文采用随机响应面配点法^[14]计算多项式混沌的展开系数。在随机向量展成的空间中，每一组样本{ξ₁, ξ₂，…}都会对应一个点，这些点称为配点，对于Wiener的多项式混沌，其展开式基函数为Hermite多项式，若Hermite多项式的最高阶数为p，则相应的配点通常取p+1阶Hermite多项式的根。未知系数可以通过式(28)和式(29)计算：

式中：ξ¹, ξ², …, ξ^k为采样点；k为采样点数。

根据Hermite多项式正交性，动态精度的均值可以通过式(30)求得^[15]：

由式(30)可以看出，动态精度均值是多项式混沌展开式的0阶项，同样的方法得到动态精度的均方差为

3.   不确定性及可靠性求解分析

3.1   谐波减速器Dymola仿真模型

根据谐波减速器动力学模型，本文针对XB1-50型号谐波减速器采用Modelica语言在Dymola^[16]编译环境下建立谐波减速器的动力学仿真模型(见图 4), 主要包括直PID控制模块、电机模块和谐波减速器模块。此模型能够较好处理谐波减速器的机、电耦合问题，在求解非线性微分方程时精度较高。

图  4  谐波减速器Dymola仿真模型

Figure  4.  Dymola simulation model of harmonic reducer

下载: 全尺寸图片 幻灯片

采用如图 5所示的谐波减速器精度测试平台对所建立的非线性谐波减速动力学仿真模型进行验证。采用恒温箱封闭，电机输入端和谐波减速器输出端分别连接有角速度和角度传感器，测定输入、输出端角位置以及角速度。实验条件下，电机输入转速为100 r/min。

图  5  谐波减速器精度测试平台

Figure  5.  Precision test platform of harmonic reducer

下载: 全尺寸图片 幻灯片

根据实测和设计参数数据确定谐波减速器动力学模型各参数如表 1所示，静态误差模型参数如表 2所示。

表  1  谐波减速器动力学模型参数

Table  1.  Parameters of harmonic reducer dynamic model

参数	数值
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

下载: 导出CSV 

| 显示表格

表  2  静态误差模型参数

Table  2.  Parameters of static error model

参数	数值
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

下载: 导出CSV 

| 显示表格

在Dymola仿真模型中通过PID调节器控制电机输出转速为100 r/min，仿真时间为10 s，待电机输入转速及误差波动曲线稳定后，通过采样输入角在3个周期均匀变化下动态误差值，实验条件下，通过角位移传感测得同样周期内谐波减速器输入端角位移θ_m以及输出端角位移θ_l，求得动态误差的测量值。将动态误差仿真值与测量值对比如表 3所示。

表  3  不同电机输入角下动态误差实验值与仿真值

Table  3.  Experimental values and simulation values of dynamic errors under different motor input corners

电机 输入角/(°)	动态误差 实验值/(°)	动态误差 仿真值/(°)	真实误差/ (°)	相对 误差/%
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

下载: 导出CSV 

| 显示表格

通过表 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%，因此在对谐波减速器进行动态误差分析时，有必要考虑柔性和摩擦的影响。

图  6  谐波减速器输入转速

Figure  6.  Speed of harmonic reducer input

下载: 全尺寸图片 幻灯片

图  7  谐波减速器输出转速

Figure  7.  Speed of harmonic reducer output

下载: 全尺寸图片 幻灯片

图  8  谐波减速器动态误差

Figure  8.  Dynamic error of harmonic reducer

下载: 全尺寸图片 幻灯片

图  9  静态误差和动态误差曲线

Figure  9.  Static error and dynamic error curves

下载: 全尺寸图片 幻灯片

3.2   动态精度灵敏度及不确定性分析

Sobol敏感度分析方法^[17]是一种基于方差分解的Monte Carlo方法，Sobol方法考虑了随机输入参量在整个取值范围内对输出响应的贡献以及随机参数的交互作用，Sobol敏感度指标可以直接从PCE式的系数得到。

动态精度PCE式(26)和式(27)可以改写为

式中:c_α为PCE多项式系数；ψ_α为PCE多项式的项；, 根据多项式混沌的正交性，可得

由式(34)可以看出分解项表征不同随机变量及其相互作用对谐波减速器动态精度输出响应方差的贡献，因此可以定义Sobol敏感度指标为

其满足

式中：S_i为主效应敏感度指标，表征各个随机变量对动态精度响应方差的贡献。

考虑到当不确定参数较多时，PCE式的待定系数呈指数增长，为了简化分析，对谐波减速器动力学模型参数进行灵敏度分析。考虑到输出轴转速相对输出轴转速小很多，输出轴转动惯量相对输入轴转动惯量大很多，故不考虑输出端阻尼和输入端转动惯量不确定参数的影响。静态误差模型中柔轮周节累积误差与柔轮安装偏心向量、刚轮周节累积误差与刚轮安装偏心向量、刚柔轮切向相邻齿综合误差在静态误差模型中形式上相同，为了简化不确定性分析，故只考虑刚轮、柔轮、波发生器安装偏心向量和柔轮切向相邻齿综合误差不确定性参数的影响。取谐波减速器不确定随机向量，对应的分布如表 4所示。多项式混沌展开的响应为时刻t₀=5.7 s谐波减速器动态误差值，取二阶混沌多项式为最高阶多项式基，n=9，p=2多项式混沌展开的待定系数为55个，三阶Hermite多项式的根为，则PCE各个随机变量的配点如表 5所示。

表  4  谐波减速器不确定性参数及分布

Table  4.  Uncertainty parameters and distribution of harmonic reducer

参数	分布	均值	标准差
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

下载: 导出CSV 

| 显示表格

表  5  动态精度多项式混沌展开式配点

Table  5.  Collocation of dynamic accuracy polynomial chaos expansion

参数	多项式混沌展开配点
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)

下载: 导出CSV 

| 显示表格

在随机向量空间ξ=[ξ₁，ξ₂，…，ξ₉]关于中心点对称进行110次抽样仿真，求得动态精度[t₀, ξ]二阶多项式混沌展开一次项系数如表 6所示。

表  6  动态精度多项式混沌展开式系数

Table  6.  Polynomial chaotic expansion coefficient of dynamic accuracy

项	系数
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

下载: 导出CSV 

| 显示表格

根据展开式一次项系数可以得到各随机变量对动态精度影响大小顺序依次为

可以得出动态误差主要与谐波减速器的动力学参数有关。在仿真模型中，分别将扭转刚度k₁，输出轴转动惯量J_l，柔轮切向相邻齿综合误差Δf′_f1的取值增加20%，得到动态误差曲线如图 10~图 12所示，对应的分别增加48.3%，减少41.7%，增加22.5%，这与PCE灵敏度分析结果基本一致。在动力学模型参数中，对动态精度影响最大的是柔轮的等效扭转刚度，当扭转刚度k₁增大时，柔轮刚性增大，动态误差值增大。其次是输出轴转动惯量，输出轴转动惯量增加，动态误差波动幅值减小。静态误差模型中对动态误差影响较大的因素是刚柔轮的切向相邻齿综合误差。可见刚柔轮的齿加工精度对谐波减速器的动态精度影响较大，是整个谐波减速器设计中的关键。

图  10  不同扭转刚度k₁时动态误差变化曲线

Figure  10.  Dynamic error curves at different torsional stiffness k₁

下载: 全尺寸图片 幻灯片

图  11  不同输出轴转动惯量J_l时动态误差变化曲线

Figure  11.  Dynamic error curves at different output shaft moment of inertia J_l

下载: 全尺寸图片 幻灯片

图  12  不同柔轮切向相邻齿综合误差Δf′_f1时动态误差变化曲线

Figure  12.  Dynamic error curves at different flexible wheel tangential adjacent gear comprehensive error Δf′_f1

下载: 全尺寸图片 幻灯片

为了减少参数不确定性分析的计算量，选取主要的不确定性随机参数向量ζ=[J_l，Δf′_f1，k₁，E_c，B_m]，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方法。

图  13  2种方法动态误差均值比较

Figure  13.  Dynamic error mean comparison between two methods

下载: 全尺寸图片 幻灯片

图  14  2种方法动态误差均方差比较

Figure  14.  Dynamic error's mean square error comparison between two methods

下载: 全尺寸图片 幻灯片

3.3   动态精度可靠性分析

本文定义谐波减速器动态精度可靠度为谐波减速器在一定初始条件和时域内，动态误差绝对值的最大值不超过一定的阈值的概率，则动态精度可靠性功能函数可表示为

式中：为动态精度的阈值；为一定时域内动态误差绝对值的最大值；t_m为最大值时间。不确定参数，取二阶多项式混沌展开，时域t₁~^{[5, 6]}，响应为，抽取42个样本点，求得各项系数如表 7所示。

表  7  动态误差多项式混沌展开式系数

Table  7.  Polynomial chaos expansion coefficients of dynamic error

项	系数
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

下载: 导出CSV 

| 显示表格

设动态误差的阈值为，根据一次二阶矩FORM法^[18]，通过式(38)求得动态精度可靠度系数β为

式中：ζ_i^*为设计验算点；，初始点取均值点。μ_ζi, σ_ζi分别为随机变量ζ_i的均值和标准差。ζ_i^*可通过式(39)求得：

通过迭代计算得到动态精度可靠度计算结果如表 8所示。

表  8  FORM法可靠度计算结果

Table  8.  Reliability calculation results based on FORM method

迭代次数	β	可靠度
1	2.783	0.997 3
2	1.782	0.962 6
3	1.737	0.958 8
4	1.732	0.958 4
5	1.732	0.958 4

下载: 导出CSV 

| 显示表格

在Dymola中进行500次Monte Carlo抽样仿真实验，得到动态误差仿真曲线如图 15所示。

图  15  动态误差的Monte Carlo仿真曲线

Figure  15.  Monte Carlo simulation curves of dynamic error

下载: 全尺寸图片 幻灯片

根据式P_f=n_f/n_t求动态精度在一定时域内的失效概率，n_f为动态误差绝对值最大值超过阈值的次数，n_t为抽样仿真次数，表 9给出了n_t=500~10 000次的模拟结果。

表  9  Monte Carlo仿真实验估算失效概率

Table  9.  Failure probability estimation under Monte Carlo simulation experiment

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

下载: 导出CSV 

| 显示表格

从表 9中可以看出, 当n=8 000时，P_f收敛，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方法计算可靠度效率更高。