机构运动时变可靠性是现有机构可靠性理论的重要内容，机构运动时变可靠度定义为在规定的运动区间上，规定的运动环境中，机构完成规定动作的能力，数学中表示为机构的实际运动输出满足期望运动输出的概率，即机构运动误差位于最大允许误差范围内的概率，其依赖于机构的运动误差建模和分析。机构运动误差问题一直受到机构学研究者的广泛关注，并提出了诸如最坏情况分析、概率分析和模糊分析等方法。Tuo等^[1-2]考虑机构参数的不确定性，建立了机构动态可靠性分析模型，并对连杆机构等典型传动机构进行了动态精度可靠性分析。董玉革等^[3]提出将机构中不确定性因素处理为模糊变量进而构造了机构模糊可靠度分析算法。张义民等^[4-5]采用Edgeworth级数方法研究了随机变量不完全概率信息和随机变量为任意分布的机构可靠性问题。但上述方法均为静态可靠性分析问题，并未衡量机构在整个工作范围内的可靠度。Zhang和Du^[6]结合结构可靠性分析中上穿越率和下穿越率概念，推导了计算机构时变运动可靠性解析算法。除此以外，传统的时变可靠性方法，也可以解决一般的机构运动精度可靠性求解问题，例如PHI2法^[7-8]，极值响应面法^[9-10]等。但这些研究只限于两态假设，即非成功即失败，这并不适用于实际工程应用中由认知不确定性带来的非两态问题分析。

基于包络思想的解决机构运动精度可靠性问题的方法由Du^[11]首次提出，其通过误差带包络和一次二阶矩方法相结合，解决机构运动误差时变可靠度的求解问题。Zhang和Du^[12]运用包络函数的方法，对铰链间隙尺寸和结构尺寸同时具有随机特性的机构进行了时变可靠性分析。Wei等^[13]提出了基于包络函数法的参数可靠性灵敏度分析和全局可靠性灵敏度分析方法，求出了参数可靠性灵敏度指标和全局可靠性灵敏度指标。Wei等^[14]利用包络函数与一次二阶矩相结合的方法，针对输入变量具随机过程特性的杆系，对整个运动行程进行了可靠度求解。但这些研究并未考虑参数的模糊性，仅考虑了随机特性，这与实际工程应用是不相符的。

综上所述，传统的包络函数方法虽然可以很好地求解机构在整个工作范围内的可靠度，但其并没有考虑普遍存在的认知不确定性问题，特别是一些参数除随机特性外，往往伴随着模糊性问题；同时，其并未考虑失效准则的模糊性，即认为产品只有“成功”和“失效”2个状态，但实际工况中，由于对故障机理、失效模式等的认知不确定性，使得常规“两态”假设无法满足实际要求。因此, 对于同时具有模糊和随机特性的机构系统的运动可靠度的研究十分必要。本文在关于机构模糊-随机时变可靠性研究的基础上，针对机构运动误差存在模糊判据和不确定参数具有模糊-随机混合特性的情况进行可靠性建模分析，对模糊判据和参数进行等效转化的同时，应用包络思想对运动误差带建立模糊-随机时变可靠性模型，将时变可靠性问题转化为时不变可靠性问题，并最终在每一截集下求出相应可靠度并进行平均加权，得到了机构系统运动在全行程内的可靠度，并最终应用于四连杆机构。

1 模糊失效判据的等效

组成机构的构件存在着加工和装配误差，这就使得构件的尺寸具有不确定性，从而造成机构运动误差的不确定性^[12]，设g(X, θ)为机构运动误差，其中X=(X₁, X₂, …, X_n)为机构构件尺寸随机变量，θ=θ(t)为机构运动角度，运动误差定义为机构实际输出ψ(X, θ)和理想输出ψ_d(θ)之差，即

由机构运动可靠度定义^[3]可知，若运动误差不大于机构允许的运动误差，则认为机构可靠；反之，则认为机构处于失效状态。机构运动误差功能函数G(X)可表示为

式中：ε为运动误差阈值。相应地，在输入角θ的范围为[θ₀, θ_e]上，其时变运动可靠度为

但上述处理方法会带来一些困难，如设机构运动误差阈值为0.025 mm，则当机构运动误差为0.025 mm时，机构是可靠的；而当机构运动误差为0.025 01 mm时，则认为机构不可靠，但这2种情况之间并无本质差别。造成这种矛盾的原因是将完好与失效状态截然分开，而未考虑中间过渡状态，即失效和完好这2个概念的外延是模糊的。为解决上述矛盾，更加准确、真实地反映机构的运动可靠性，有必要将模糊数学方法引入机构的运动可靠性分析中。

失效边界(误差阈值ε)定义为某一模糊区间[z_L, z_U]，其中z_L和z_U分别为模糊区间的下界和上界，把机构运动误差失效准则看作模糊事件，描述事件状态程度隶属函数μ_G(z)可以用来表示这一过渡情况，其数值越大，事件失效的倾向越大，数值越小，失效的倾向越小。机构运动失效概率可表示为

式中: z为机构模糊-随机混合空间Ω中的随机变量，其概率密度函数为f_G(z)；μ_G(z)为描述事件状态程度的隶属函数，其范围为0≤μ_G(z)≤1。

若μ_G(z)为递减函数，即使得机构运动的失效程度随z值的减小而增大，根据文献[15]，结合随机变量的概率分布函数的定义，可以把1－μ_G(z)看作一个新的随机变量(记为Z′)概率分布函数，模糊-随机失效域可以描述为{X|G(X)≤Z′}，等效的功能函数为G^e=G(X)－Z′。

若μ_G(z)为递增函数，则可以把μ_G(z)看作一个新的随机变量(记为Z″)概率分布函数，同理，这种情况下机构运动的失效域描述为{X|G(X)≥Z″}，对应的等效功能函数为G^e=Z″－G(X)。

2 运动可靠度的一般包络方法

对于机构输入角范围为[θ₀, θ_e]，由X的随机特性易知在误差阈值ε下，g(X, θ)=ε实为一族曲线，包络法^[11]定义了对[θ₀, θ_e]上运动误差曲线族g(X, θ)=ε中所有曲线的包络，定义包络方程G⁺(X)=0和G^－(X)=0分别对应上界包络和下界包络，即曲线族的上边界和下边界，机构运动误差的包络函数如图 1所示^[12]。由定义知在包络曲线上，函数值应处处等于ε，同时由于包络曲线包络了不同θ下的曲线，所以其对θ求偏导应为0。

在包络方程给出之后，考虑在整个运行周期中的机构运动误差时变可靠性问题，可表示为

式中：S⁺={G⁺(X) < 0}；S^－={G^－(X)>0}。

上边界包络方程G⁺(X)=0为

下边界包络方程G^－(X)=0为

式中:代表偏导数∂g(·)/∂θ。

由以上分析可知，可对每个包络方程的2个式子进行联立进而消除θ，时变可靠性问题则转变为时不变可靠性问题。

3 模糊-随机时变可靠性建模

考虑具有模糊-随机混合不确定性的机构，在其输入角范围[θ₀, θ_e]内，其运动误差可表示为

式中：为模糊-随机时变运动误差；θ)为实际运动输出；ψ_d(θ)为理想运动输出；为模糊-随机混合变量。

结合前文所述，其模糊-随机时变运动可靠性功能函数可表示为

根据第1节的方法，考虑模糊失效判据，则上边界失效判据μ_G^u(z)为递增函数，下边界失效判据μ_G^l(z)为递减函数，可用新的随机变量来描述失效模糊-随机事件，由于机构运动误差上下界的对称性，设Z′=Z″=Z，则运动误差的上界失效域描述为，同样，下界失效域描述为。运动误差的失效隶属度函数如图 2所示，本文认为其模糊失效状态的隶属函数在过渡区间为线性变化的，γ为由模糊失效判据造成的状态变化值。阈值模糊的可靠性失效事件描述如图 3所示。应当注意到，根据以上分析可知这里的Z为均值为0的某种对称随机分布。

对于功能函数中的模糊变量，本文全部取三角形隶属度函数，根据文献[16]的理论及文献[17]的方法，将模糊变量处理成为概率密度函数中含有参数α的随机变量，其中0≤α≤1。设实际运动输出中的模糊变量的隶属函数为，得到给定的α水平下的截集和对应的概率密度函数f_α(x)，相应地，θ)转换为α水平下的运动误差g(X_α, θ)。

至此，在含有模糊-随机混合变量和模糊判据情况下，α水平下的机构运动的可靠度为

后续通过对α的离散求和或数值积分可以求得机构运动的失效概率。

根据以上分析，在输入角范围[θ₀, θ_e]上定义α水平下的上界包络方程为G_α⁺(X)=0，下界包络方程G_α^-(X)=0。考虑模糊失效判据，α水平下机构运动可靠度可通过下列事件来表示:

式中：。

G_α⁺(X)=0表达式为

G_α^－(X)=0表达式为

此时，通过改进包络方程，功能函数已转变为仅含参数α的随机变量的表达式，但继续通过联立消除θ的方法，一般很难获得其精确的解析表达式，下面用近似方法对包络函数进行求解。

4 基于改进包络函数的时变可靠度计算 4.1 模糊-随机时变可靠度的近似求解

由于包络函数一般具有较强的非线性特性^[15]，若通过传统的可靠度求解方法，即在一个展开点展开构造其近似表达式，会导致与原式误差较大。为充分利用包络函数特性，将包络函数在多个展开点上分段线性化，构造近似包络函数并求出每段的展开点，最终通过求这些展开点的联合分布函数最终求得可靠度。

假设X_α=[X_1α, X_2α, …, X_nα]为相互独立的正态分布随机变量，对g(X_α, θ)－Z在X_α的均值μ_Xα处用一阶泰勒级数展开来近似运动误差，并把X_α化成标准正态随机变量U，展开后表达式记为L_α(U, θ)，则由式(12)与式(13)知：

G_α⁺(U)=0表达式为

G_α^－(U)=0表达式为

式中：U=[U₁, U₂, …, U_n]，U_i=[(X_iα－μ_Xiα)/σ_Xiα]_{i=1, 2, …, n}，σ_Xiα为X_iα标准差；a_α(θ)=ψ(μ_X, θ)－ψ_d(θ)－S(G)μ_Z；b_α(θ)=((∂ψ/∂X_i|_μXiα·σ_Xiα)_{i=1, 2, …, n-1}, S(G)·σ_Z)；S(G)=。

设待求展开点为U_α(θ_i)(i=1, 2, …, m)，由展开点特性易知原点到展开点的向量U_α(θ)和L_α(U, θ)在U_α(θ)处的梯度共线，即U_α(θ)=c·b_α(θ)/，其中c为常数。与式(14)联立，并根据G_α⁺(U)特性，去除解中运动误差为负对应的角度，即可求出上边界包络G_α⁺(U)的展开点U_α(θ_i⁺)与展开点对应的角度θ_i⁺(i=1, 2, …, m⁺)。同理，也可以得到上边界包络G_α^－(U)的展开点U_α(θ_i^－)与展开点对应的角度θ_i^－(i=1, 2, …, m^－)。至此，上下边界的近似包络展开点及对应角度已全部求出。

另一方面，由式(14)与式(15)可知，近似运动误差S(θ_i)L_α(U, θ_i)为正态分布，由此可以用展开点的高维正态分布函数(均值μ_α和协方差阵Σ_α)对失效概率进行数值求解：

式中：为保守计算，θ_i包括θ_i⁺，θ_i^-和两端点(θ₀和θ_e)，m=m⁺+m^－+2，ε_i=s(θ_i)ε。易知，Σ_α应是正定阵，如果条件不满足，那么并不是所有角度都需要。由于失效概率小的展开点必定被包含在失效概率大的展开点所线性展开的近似包络内，即在算可靠度时失效概率小的点是多余的，设其秩为r，则应算出所有θ_i对应的失效概率，并去除那些m-r个最小失效概率点所对应的时刻，失效概率的算法如下

从而得到新的均值向量μ′_α和新的协方差阵Σ′_α，则通过式(20)的数值算法(通过积分)可算出α水平下的可靠度累积分布函数(CDF)R_α(θ₀, θ_e):

将α在[0, 1]上离散n等份，则阈值为ε下的可靠度可表示为

特别地，若α连续且所有水平中去除的多余展开点均相同，则阈值为ε下的可靠度可表示为

4.2 算法流程

步骤1   将模糊功能函数通过第1节的方法，用新的随机变量来描述失效判据。

步骤2   将α在[0, 1]上离散n等份，α初始化为0。将步骤一的结果转化成α水平下的仅含随机变量的功能函数G_α。

步骤3   构建α水平下的包络函数G_α⁺(X)和G_α^－(X)。

步骤4  对包络函数G_α⁺(X)和G_α^－(X)进行分段线性化近似，得到G_α⁺(U)和G_α^－(U)。

步骤5  根据式(14)、式(15)与向量U_α(θ)分别联立求解G_α⁺(U)和G_α^－(U)的解θ_i⁺和θ_i^-，此过程中要用到通过G_α求解得到的a_α、、b_α(θ)和。

步骤6   通过式(16)和式(17)来计算均值μ_α和协方差阵Σ_α，并通过协方差阵的秩判断其是否正定，若正定，则进行步骤7；若非正定，则结合式(19)求每一展开点的失效概率决定去除失效概率最小的m-r个时刻，得到μ′_α和Σ′_α。

步骤7   根据式(20)算出α水平下的可靠度并记录，返回步骤2；若α=1时，进行步骤8。

步骤8   通过对R_α(θ₀, θ_e)的离散求和或数值积分，得到阈值ε下的可靠度，算法流程图见图 4。

5 四连杆机构运动可靠度计算 5.1 四连杆机构运动学建模与分析

本例用本文方法对四连杆机构进行运动可靠度分析，首先，对四连杆机构进行运动学建模。在二维空间F(x, y)中，机构如图 5所示，输入角为θ，输出角为ψ。

四连杆机构的变量为，其中为模糊变量，R₁、R₃和R₄为随机变量，其运动输出方程组为

首先考虑模糊判据，根据第1节的分析，将模糊阈值等效为随机变量，本例γ取0.01，对比式(10)及分析可知，应引入随机变量Z，它的概率累积分布函数为

另一方面，针对模糊变量，根据实际测量数据和相关的专家经验，认为其隶属函数为

则α水平下的概率密度函数为

则α水平下R₂的均值和标准差分别为

其运动输出方程组在α水平下的表达式为

由式(27)和式(28)可得α水平下尺寸随机变量分布的数字特征，如表 1所示。

通过以上对模糊变量的等效，由式(23)解得α水平下的ψ(R_α, θ)为

式中：

进而由式(15)可得求得b=(b₁, b₂, b₃, b₄)。

5.2 四连杆机构模糊-随机时变可靠度求解

如图 5所示的四连杆机构理想的运动方程为

输入角范围为[θ₀, θ_e]，其中θ₀为0°，θ_e为90°，则其理想运动输出ψ_d(θ)为输入角从0°~90°间输出端的角度。

在可变阈值ε上，用本文提出的计算方法和蒙特卡罗仿真(Monte Carlo Simulation，MCS)法的结果进行比较，其中，由于蒙特卡罗法的仿真次数规模为10⁷级，认为其计算结果是趋于真实值的。

运动误差ψ(R_α, θ)－ψ_d(θ)在变量取均值时的变化如图 6所示。这里以阈值ε=0.5为例展示计算结果，本例将截集水平α在[0, 1]上分为50等份，α=1时，θ₀=0°，θ₁^- = 20.196 5°，θ₁⁺=53.812 1°，θ_e=90°，结合式(17)和式(18)得μ_α和Σ_α并知其秩r=3，即多余一个时刻，通过式(19)得到P_f(θ₀)=0.303 0，P_f(θ₁^-)=0.382 1，P_f(θ₁⁺)=0.074 6，P_f(θ_e)=0.248 1，去除失效概率最小的时刻θ₁⁺，则更新的μ_α和Σ_α分别为

由式(20)可得可靠度R₁(θ₀, θ_e)=0.313 3。其余截集水平下的可靠度结果重复以上步骤即可，计算结果如图 7所示。将其平均加权即为ε=0.5时的运动可靠度。

当ε在[0.2, 1.4]上变化时计算结果如图 8所示，由失效概率曲线图表明本文提出的模糊-随机时变可靠度计算方法与MCS法在失效概率的计算上差别较小，计算结果准确度可以认为满足机构运动误差的分析要求，对工程实际应用有一定的参考价值。部分误差阈值对应的结果见表 2。

从计算结果可以看出，基于本文方法所得的结果均略大于MCS法的结果，这是因为MCS法是对所有时刻(角度)进行遍历并进行求解，而基于本文提出的方法是对有限的、失效概率较大的展开点进行计算，所以可靠度计算结果均比MCS法略大一些。

6 结论

本文针对机构产品的时变可靠性建模分析问题，综合考虑了模糊性和随机性，针对机构运动误差建立了机构运动模糊-随机时变可靠性模型，提出了模糊-随机时变可靠度求解方法，经案例验证表明：

1) 本文提出的模糊-随机时变可靠度求解方法相较于传统时变可靠度求解方法，不但考虑了参数的模糊性，同时还考虑了判据的模糊性，解决了模糊-随机混合参数下的机构时变可靠性问题，更符合实际工程应用，对类似的机构运动可靠性分析具有一定的指导意义。

2) 本文方法求解方便，相比于MCS法计算误差较小，本案例中2种方法的计算误差最大不超过0.000 8，贴合度较高；同时计算效率大大提高，每α水平下MCS法计算次数为10⁷次，而本文方法为50次左右。

3) 本文提出的机构运动模糊-随机时变可靠性分析方法适用于随机样本不完善，或样本数量不够的情况，可计算出较精确的结果。