飞机自动化装配设备有效工时分析

（西北工业大学机电学院，西安 710072）

目前，国内各大飞机制造厂所生产机型众多，装配任务繁重，但总体上自动化装配设备的利用率不足40%。导致设备利用率普遍偏低，除了某些关键瓶颈技术及专业人员培训周期长等因素制约外，另一个重要原因就是对设备的生产能力预估不足，生产计划安排不合理，盲目追求以高技术带来高生产率。具体表现为对设备的可靠性掌控不足，缺乏对装配设备有效工时的统计、分析及编排计划的经验。现有的装配设备可靠性评价指标主要包括设备可用度、故障间隔时间、平均修复时间等，其中设备可用度描述了可维修的装配设备在规定的时间、装配工艺流程和使用条件下，在规定的寿命周期内运行无故障的概率[1]。因此，研究飞机自动化装备的可用度，了解设备在确定的时刻、环境下能正常工作的概率以及实际工时消耗情况具有重要意义。

设备可用度分析理论中，威布尔分布模型广泛应用于拟合故障数据、计算设备的可靠率[2]，该模型包含3个参数，分别控制模型函数的位置、形状和尺度[3]。由于飞机自动化装配设备通常包含多类子系统，各类子系统的失效形式及规律多种多样，且各类失效近似为相互独立事件。因此，理论上可利用混合威布尔分布模型分析飞机自动化装配装备的可用度，主要原因包括：（1）混合威布尔模型的参数多样化，模型更有弹性，能够描述包含多种形状的可用度分布曲线，适用性较强；（2）混合威布尔模型能描述具有多个子系统的复杂设备可用度分布，分析总体中各子系统的可用度及其对总体的影响[4]。

文献[5]～[8]中，现有的国内外研究中，针对复杂机械设备的可靠性机理分析及工程应用方面已有大量文献，但研究成果难以直接应用于本文环境中，主要原因包括：（1）现有文献分析设备可用度的威布尔分布函数多为静态模型，难以反映设备各系统可用度的动态变化；（2）针对飞机自动化装配设备的可用度分析并计算其有效工时输出方面鲜有研究。综上所述，本文重点考虑了飞机自动化装配设备的可用度对设备工时的影响，运用动态混合威布尔分布模型求解设备的动态可用度，并根据可用度结果计算设备的有效工时。

动态混合威布尔分布建模

飞机装配线整体由包含大量工序的装配站群构成，每道工序又包含了多个工步，工步又可进一步拆分为基本操作单元[9]。这些数量庞大、不可进一步分割的操作单元构成了飞机装配作业的基本组成元素。为提高飞机自动化装备的有效工时解算精度，有必要将计算单位由工序进一步分解为更小的统计元素，分解成“元素”形式后作业时间以更容易被测量及度量[10-11]。基于以上考虑，提出“装配作业元素”（Assembly Operating Elements，AOEs），即为产品装配过程中无法进一步分解的最小操作单位和任务实体[12]，每个装配作业元素都对应一个确定的标准时间。

首先，结合本文应用环境定义假设条件为：（1）由自动化装配设备构成的工作站，每个站位仅有单台设备；（2）每台设备仅完成一类特定的装配作业，不考虑多种作业模式混合装配；（3）由于单台设备作业，装配站中的工序均为串行作业模式；（4）不考虑装配对象的多样性引起的工时差异。

然后，设自动化装配设备的定检周期为TD，表示一台设备在一个维护周期内可完成的最大装配作业元素总数。在任意一个定检周期内，以装配作业元素的序号x作为随机变量，设在装配第x个AOE时设备的可用度概率密度函数服从三参数威布尔分布，记做f（x）=W（m，η，γ，x），且各参数间相互独立，其概率密度函数和可用度函数分别为[3]：

式中，m为形状参数，决定分布曲线的形状；η为尺度参数，起缩小和放大曲线横坐标尺度的作用，控制了分布曲线在横坐标轴方向上的离散程度；γ为位置参数，决定分布曲线的起始位置，不改变曲线的形状，当γ=0时f（x）简化为双参数威布尔分布。

同时，设装配设备具有N类子系统，每类子系统具有特定的失效模式，每种失效模式的规律均服从三参数威布尔分布，每种模式的动态权重因子为pi（i=1，2，…，N），则该设备整体的失效概率密度及可用度函数可表示为：

式中，共有4N-1个未知参数，、其中包含3N个威布尔分布参数组合mi、ηi、γi和N-1个权重系数pi，如何准确估计这些未知参数是建立动态混合威布尔分布模型的核心和关键[13]。下文中将详细推导动态权重系数pi的计算流程；利用最大似然估计法估算威布尔分布参数组合，求解设备可用度后分析有效工时。

装配设备子系统动态权重值计算

为描述自动化装配设备各类子系统对系统总体可用度的影响，本节以飞机装配线中的自动钻铆设备为对象，分析了此类自动化装配设备的系统构成及特性，确定各子系统的动态权重运算规则。一般地，一套功能完整的自动钻铆设备主要包括3大系统：控制系统、定位及执行系统、测量及反馈系统。

1 控制系统

控制系统作为自动化装配设备的大脑，负责设备所有机构运动执行的逻辑控制，其可靠性与稳定性直接决定了设备能否按装配指令完成任务[14]。在设备使用的初始阶段，对控制系统的开发及调试任务繁重，所占工作量大，通常需要专业人员频繁进行仿真模拟、在线/离线编程以及程序校准等工作，容易积累错误并致使可用度降低。因此，控制系统的可用度权重在装配前期应具有较高的影响度。随着设备进入稳态运行阶段，其可用度权重影响度开始下降，并趋于常量。为定量化描述该动态权重的非线性衰减规律，设其权重为p1（x），服从双参数威布尔分布，并预设初始值为0.6，则p1（x）满足：

当形状参数mp1＜1时，威布尔曲线变换为单调递减区间的伽玛分布。在满足初始值的条件下，可解得mp1由0.5增加至1过程中对应的尺度参数值ηp1，并绘制出动态权重p1（x）随x变化的函数曲线，如图1所示。其中，横坐标TD的数量级为万个，定检周期TD设为10万个。由图1可知，随着mp1的增加，权重曲线的衰减过程更加缓和。当mp1=1时，p1（x）退化为指数分布。为避免控制系统可用度权重的过速或过缓衰减，取mp1=0.9、ηp1=2.846 （图1中蓝色曲线）。

图1 形状参数对控制系统动态权重的影响Fig.1 Effect of shape parameters on dynamic weight of control system

2 定位及执行系统

定位及执行系统是自动化装配设备的硬件主体及具体执行机构，主要包括定位装夹机构、运动驱动机构和末端执行机构。该系统的可靠性主要取决于各硬件的精度、质量及稳定性等。由于该系统经安装调试定型后结构不再变化，因此其可用度权重在定检周期的前、中期相对稳定；进入中后期时由于硬件磨损及老化，可用度逐渐降低，此时权重影响度开始上升。同理，设该系统动态权重为p2（x）并服从双参数威布尔平移曲线，权重初始值设为0.2，则p2（x）满足：

当形状参数mp2=4时，威布尔曲线接近于正态分布，其单调递增区间恰好能描述定位及执行系统可用度权重的动态增长规律。当mp2确定后，尺度参数ηp2由5增加至15时，权重曲线沿横轴方向伸展并在纵向压缩，如图2所示。此时设定当x的数量到达定检周期时，权重值达到上限0.3，即p2（x）∈[0.2，0.3]，可解得ηp2=12.931（图2中蓝色曲线）。

3 测量及反馈系统

测量及反馈系统用于从装配现场测量相关数据，通过实时修正误差指导控制系统使设备正常运行，保证装配准确度并提高自动化程度[15]。由于运行过程中随着设备的磨损、疲劳，对测量及反馈系统的依赖程度不断提高，因此该系统的权重影响度在定检周期内应逐步上升并趋于稳定。在满足∑pi=1的条件下，p1、p2的动态权重表达式已确定，因此可直接求出p3（x）。最终，3种系统可用度权重的动态变化关系如图3所示。由图3可知，p2（x）在总体上随x的递增而上升，并在后段稳定，因此符合预设规律。

装配设备有效工时求解

1 动态混合威布尔分布的参数估计

针对参数估计问题，现广泛使用的方法主要包括：图解法、回归分析法、最大似然估计法以及Bayes估计等。其中，图解法及回归分析法的计算模型相对简单，能快速获取各参数的预估值，适合在没有先验信息的条件下进行点估计或区间估计，但随未知参数的增加后求解准确度下降，适应度较差。Bayes估计利用后验分布统计推断结果，其后验分布集中了样本信息、总体信息和先验信息，排除了与未知量无关的信息，通常比先验分布更接近实际情况，估计结果相对理想[16]。但由于混合威布尔分布的参数众多且形式多样，导致Bayes估计中的后验分布更加复杂，因积分没有显式表达式而难以准确推断[17]，因此Bayes估计难以应用于本文环境中。而最大似然估计利用总体的概率密度分布及其子样本提供的信息求解未知参数的估计量，兼具较好的求解精度和效率，适合混合威布尔分布的多参数估计。

因此，建立最大似然估计模型：设自动化装配设备各子系统i的概率密度函数的参数为θi（mi，ηi，γi）。在一个定检周期内，测量各子系统的AOE失效次序，其独立同分布的子样本试验观测结果为xij（j=1，2，…，nSi），即有nSi个AOE装配失效。则设备在[xij，xi（j+1）]区间内失效的概率是f（xij）dxij，建立最大似然函数：

对式（5）左右两边同时取对数：

求解似然函数的极值，令∂L/∂θi=0，获得方程组：

图2 尺度参数对定位及执行系统动态权重的影响Fig.2 Effect of scale parameter on dynamic weight of position and execution system

图3 3种子系统可信度动态权重Fig.3 Dynamic weights of reliability of three sub-systems

由于上式为复杂非线性方程组，难以直接求出精确解，本文采用Newton-Raphson法进行迭代求解[18]，具体过程为：首先选取待估计参数初始值，现有方法多采用图解法，防止因初始值选择不当引起的迭代发散。在此基础上，本文将该问题先简化为双参数威布尔分布的参数估计，可快速得到形状、尺度参数的估计值，并将其作为混合三参数威布尔分布的初始值mi0、ηi0，再结合图解法经验预估位置参数的初始值γi0。然后，在初始值θi0（mi0，ηi0，γi0）处将似然方程组的左侧各项级数展开，一阶近似后将方程组转化为线性方程组以解算Δmi、Δηi、Δγi。最后定义误差界限并判断Δmi、Δηi、Δγi是否小于该界限，若小于则取Δmi、Δηi、Δγi为估计值；若大于则用θi0+Δθi替换Δmi、Δηi、Δγi重新计算，经迭代循环后获得各参数的估计值

2 设备有效工时求解

设有自动化装配设备的工作站内，共包括NE个AOE，每个AOE的首次测量值为hE，站内的生产准备时间、清理时间等辅助工时数设为HR，则该站位的理论工时数HT为：

由于理论工时为理想状态下的静态标定工时，不能反映实际生产条件下随装配架次的增加设备的可靠性下降所引起的工时延长作用。因此，考虑设备可用度后，利用获取的可用度值将设备完成任意AOE消耗的理论工时进行等效放大，再累加后即可获得并模拟设备在实际环境中所耗的有效工时。具体过程为：结合公式（2）～（7），可求出设备在第x个AOE的可用度为R（x），失效率为F（x）=1-R（x），则该设备站位装配首架次飞机的有效工时数HE为：

同理，该站位装配第X架次的有效工时数RE（X）及有效装配周期TE（X）为：

式中，x=k+（X-1）NE，hD表示自动化装配设备单工作日的制度工时数，表示上取整函数。

实例验证

为检验本文方法的有效性，取某型飞机部装生产线中的机翼壁板组件装配站位内的自动钻铆设备为验证对象。首先，根据历史统计记录，提取该设备的3大系统在一个定检周期（10万个AOE）内发生失效的AOE次序，如表1所示。

根据表1的分布数据，对动态混合三参数威布尔分布参数估计，利用Matlab迭代计算各参数值，包括9个待估计参数和9个初始值，计算结果见表2。根据参数估计结果，可求出自动钻铆设备3大系统的可用度函数R1（x）、R2（x）、R3（x），如图4所示。同时，结合公式（3）、（4）解出的动态权重函数p1（x）、p2（x）、p3（x），最终可解得该设备整体的运行可用度分布函数，如图5所示。

表1 某自动钻铆设备加工失效序列

表2 混合威布尔分布参数估计值

图4 3种系统可用度分布曲线Fig.4 Distribution curves of reliability of three sub-systems

图5 自动钻铆设备的可用度示意图Fig.5 Reliability of the automatic drilling & riveting equipment

在解得设备可用度的基础上，现场统计该设备所在的站位内包含的AOE总数为NE=1852，单个AOE的标准测量时间hE=35s，辅助工时数HR=12h。利用公式（10）可求出该设备在其站位内装配第X架产品所耗的有效工时。同时，根据某型机装配履历系统的数据记录，获取了该站位在前30架次的生产中各架次所耗的工时数，作为实测参考值，验证本文方法的有效性，两者的数据对比如图6所示，两组数据的偏差分析如图7所示。

根据图6，本文方法解出的设备有效工时曲线为随架次X的提高单调递增的凸函数，模拟了设备随着装配周期的增长，实际消耗工时的缓慢递增过程，能够反映出实测工时的基本走势，两组数据的相关系数达到91.6%。根据图7，两组数据的最大偏差绝对值为4.195h，最大偏差百分比小于11.2%；偏差的累计标准差为1.717h，波动峰值小于4.6%，说明本文方法求解的自动化装配设备有效工时能够较好的模拟设备在实际环境中的生产耗时，求解偏差较小，满足产能分析的准确度要求。

由图6、图7可见，前10架次的实测工时均高于计算获得的有效工时。10架次以后实测工时增长缓慢，略低于有效工时，主要原因在于：在装配的初始阶段，人员对新设备的使用需要一定的熟练周期或磨合期，需要经过大量的试验、调试后使设备逐渐进入稳态，因此在该阶段实测工时高于有效工时。随着产量的增加，人员的学习作用对辅助工时产生影响，使得辅助工时数逐渐减小并部分抵消了设备可用度引起的工时增长，因此该阶段的实测工时低于有效工时。

结论

（1）本文方法求解飞机自动化装配设备有效工时准确度较高，最大偏差百分比小于11.2%，偏差的波动峰值小于4.6%，能够较好的模拟此类设备在实际飞机装配环境中的生产耗时情况。

图6 设备有效工时解算结果与实测值的对比Fig.6 Comparison of working hours between the proposed method and the measured values

图7 设备有效工时与实测值偏差比较Fig.7 Comparison of deviation values between the proposed method and the measured working hours

（2）在自动化装备生产的初始阶段，由于人员同设备的磨合及调试校准会导致实测工时高于有效工时；在批生产阶段，由于人员的学习作用使辅助工时数缩减，致使该阶段的实测工时低于有效工时。

[1]BLANCHARD B S. Logistics engineering and management[M]. London:Prentice Hall, 2004.

[2]高萍. 基于可靠性分析的复杂设备预防性维修决策研究 [D]. 北京：清华大学，2008.GAO Ping. The research on preventive maintenance decision of complex equipment based on reliability analysis [D]. Beijing:Tsinghua University, 2008.

[3]WEIBULL W. A statistical distribution function of wide applicability[J]. Journal of Applied Mechanics, 1951,18(3): 293-298.

[4]凌丹. 威布尔分布模型及其在机械可靠性中的应用研究 [D]. 成都：电子科技大学，2011.LING Dan. Research on Weibull distribution and its applications in mechanical reliability engineering [D]. Chengdu: University of Electronic Science and Technology of China,2011.

[5]沈安慰, 郭基联, 王卓健. 航空装备现场数据可靠性评估方法有效性分析[J]. 航空学报, 2014, 35(5):1311-1318.SHEN Anwei, GUO Jilian, WANG Zhuojian.Validity analysis of reliability evaluation method in aviation equipment field data [J]. Acta Aeronautica et Astronautica Sinica, 2014,35(5):1311-1318.

[6]GEBIZLIOGLU O L, SENOGLUA B, KANTAR Y M. Comparison of certain value-at-risk estimation methods for the twoparameter Weibull loss distribution[J]. Journal of Computational and Applied Mathematics, 2011,235(11): 3304-3314.

[7]IWINSKA M, SZYMKOWIAK M.Characterizations of distributions through selected functions of reliability theory[J]. Communications in Statistics-Theory and Methods, 2017, 46(1):69-74.

[8]ALTUN M, COMERT S V. A changepoint based reliability prediction model using field return data [J]. Reliability Engineering &System Safety, 2016, 156: 175-184.

[9]XIN B, LI Y, YU J F, et al. A multifactor revision based analysis of the personnel operational capacity of aircraft assembly lines[J]. The International Journal of Advanced Manufacturing Technology, 2015, 78(1-4): 211-220.

[10]ORTIZ C A. Kaizen assembly:Designing, constructing, and managing a lean assembly line[M]. Boca Ratai: Taylor & Trancis Gxrp, 2006.

[11]NIEBEL B W. Motion and time study design and measurement of work[M]. IL: Irwin Professional Publishing, 1982.

[12]朱振杰. 机械产品装配线平衡问题优化研究[D]. 济南：山东大学，2010.ZHU Zhenjie. Study on mechanical products assembly line balancing optimal algorithms [D].Ji’nan: Shandong University, 2010.

[13]吴江. 飞机可靠性分析中的混合威布尔分布参数估计方法[J]. 工程设计学报,2015, 22(1): 26-29.WU Jiang. Parameter estimation method of mixed Weibull distribution in reliability analysis of aircraft[J]. Chinese Journal of Engineering Design, 2015, 22(1): 26-29.

[14]李冬磊. 飞机数字化柔性装配钻铆加工控制技术研究[D]. 南京：南京航空航天大学，2015.LI Donglei. Research on control technology of aircraft digital flexible assembly drilling &riveting [D]. Nanjing: Nanjing University of Aeronautics and Astronautics, 2015.

[15]贺美华. 基于多传感器集成的机器人自动钻铆在线检测技术研究[D]. 南京：南京航空航天大学，2013.HE Meihua. Research on online detection technology of robotic drilling and riveting system based on integration of multiple sensors[D].Nanjing: Nanjing University of Aeronautics and Astronautics, 2013.

[16]张月. 指数Weibull分布的贝叶斯估计与模拟 [D]. 武汉：华中科技大学，2012.ZHANG Yue. The Bayesian estimation for exponentiated-Weibull distribution and its simulation[D]. Wuhan: Huazhong University of Science & Technology, 2012.

[17]刘飞, 王祖尧, 窦毅芳, 等. 基于Gibbs抽样算法的三参数威布尔分布Bayes估计 [J]. 机械强度, 2007, 29(3):429-432.LIU Fei, WANG Zuyao, DOU Yifang, et al.Bayesian analysis of three-parameter Weibull distribution based on Gibbs sampling algorithm[J].Journal of Mechanical Strength, 2007, 29(3):429-432.

[18]李颖红. 定时截尾场合三参数Weibull分布的Bayes估计[D]. 北京：北京交通大学，2010.LI Yinghong. Bayes estimators of the threeparameter Weibull distribution based on type I censoring sample [D]. Beijing: Beijing Jiaotong University, 2010.