钱佳诚，陈换过，陈文华，李剑敏，黄春邵，沈建洋，许小芬

(1.浙江省机电产品可靠性技术研究重点实验室，浙江 杭州 310018;2.杭州机床股份有限公司，浙江 杭州 311305)

磨削过程中颤振是机床故障的主要表现形式之一.磨削颤振大致可以分成3个阶段：磨削稳定阶段、颤振孕育阶段和颤振稳定阶段[1].颤振孕育阶段是抑制颤振的关键时机，该阶段颤振特征信息开始增多，并且振动幅值在可控范围，如果不能及时预测并采取抑制措施，颤振就会迅速增大，对磨削加工造成严重影响.因此，寻求有效的颤振预测方法，对实际生产具有重大的现实意义和经济价值.

国内外许多学者曾对颤振预测方法进行了研究：Bediaga I等采用支持向量机的方法对颤振进行了预测[2]；吴石等提出一种应用小波系数特征和多类超球支持向量机进行铣削颤振预报的方法[3]；Karpuschewski B等提出了基于声发射信号的功率谱，采用神经网络方法对磨削颤振进行预测[4]；Tansel I N等对机床尾座采集到的加速度信号采用拉普拉斯变换进行处理，以阻尼系数作为预测颤振发生与发展的特征量，并采用模糊逻辑方法对车削加工中的颤振进行了预测[5].上述方法均通过监督学习实现了对颤振的预测，但是并未给出一个明确的颤振发生阈值.王民等提出基于控制图理论的磨削颤振预测方法，成功地给出了颤振发生的阈值.但是，该方法中直接监测的是原始数据，并未去除信号中无效成分的干扰[6].沈建洋借助二维经验模态分解(BEMD)方法，将磨床振动信号逐级分解为若干二维的特征模式函数后，找到了一些更能全面反映颤振信息的特征模式函数[7].

本文采用BEMD和控制图相结合的方法对磨削颤振进行预测，解决传统预测方法难以给出阈值的难题.首先，通过BEMD方法对信号进行分解，得到有效的二维固有模态函数(BIMF)，提取颤振预测特征量-实时方差；然后，选取一段平稳阶段的特征量序列进行正态性和相关性检验；接着，采用拉依达准则(±3σ)得到控制图边界，建立控制图；最终将后续特征值按时间序列输入模型进行颤振预测.

1 BEMD与控制图

1.1 BEMD基本理论

BEMD适用于分析非线性、非平稳信号，是针对二维复值信号的分解方法.

(1)输入原始二维复值信号s(x,y)，计算投影方向φm，有

φm=2mπ/N，1≤m≤N

(1)

式中：N代表信号投影方向数；m代表迭代次数.

(2)将信号s(x,y)投影到φm的各个方向上，得到投影pφm(x,y)，有

pφm(x,y)=Re[e-iφms(x,y)]

(2)

(3)

(4)

判断分量h(x,y)是否满足固有模态函数的筛选停止条件，如果不满足，则将h(x,y)替换成新的原始信号s(x,y)，并返回步骤(2)，直至满足固有模态函数的条件为止；否则，将分量h(x,y)作为第一个固有模态函数分量c1(x,y)，用原始信号s(x,y)减去第一个固有模态函数分量c1(x,y)，获得残余量r1(x,y)，有

(5)

(6)将r1(x,y)作为新的输入序列，重复上述步骤直至获得所有的固有模态函数和一个最终的残余量rn(x,y)为止.因此，原始信号s(x,y)可以表示为：

(6)

式中：ck(x,y)代表第k阶BIMF；n代表BIMF的总个数.

1.2 控制图基本理论

控制图的基本原理是利用预测对象的历史观察值来确定控制中心线(Central Line,CL)、上控制限(Upper Control Limit,UCL)和下控制限(Lower Control Limit,LCL)，并按时间顺序输入待监测的数据，观察数据是否在控制限范围内，若出现异常现象，则采取相应措施进行维护.

控制图形成的前提是数据为正态分布.正态分布有两个很重要的参数:平均值μ和标准差σ.不论两者取何值，特性值落在(μ-3σ,μ+3σ)内的概率为99.73%，而特性值落在(μ-3σ,μ+3σ)之外的概率为0.27%.本文采用±3σ法进行计算，即以被控制产品质量特性值分布中心μ为CL,则UCL=μ+3σ,LCL=μ-3σ.

用控制图对加工过程是否失控进行判别的一般标准为：如果有数据点超过控制界限或者排列不随机，则加工过程处于失控状态.对于故障状态的判断准则有以下4条:①点落在控制限外;②连续8个点落在中心线的同一侧;③连续6个点递增或者递减;④连续14个相邻点交替上下.

控制图用于颤振预测的主要步骤为：①获取正常磨削状态的特征值时间序列，对其进行单一性、正态性和相关性分析；②计算统计特性值的均值和标准差，使用±3σ法计算出控制图的上、下控制界限，从而确定控制图边界；③利用建立的控制图对后续特征值时间序列进行判断，以达到颤振预测的目的.

2 仿真信号与BEMD结果

2.1 仿真信号时域图

通过Matlab Simulink模块创建磨削颤振信号发生系统，来模拟颤振的状态.

仿真总时间为2.5 s，采样频率为800 Hz.在0～1 s内，模拟平稳磨削阶段，得到的信号振动幅值在较小范围内波动；在1～1.6 s内，模拟颤振孕育阶段，得到的信号振动幅值随时间逐渐增大；在1.6～2.5 s内，模拟颤振稳定阶段，得到的信号振动趋于稳定，而且振动的幅值明显增大.仿真信号的时域如图1所示.

图1 仿真信号时域图

2.2 BEMD结果

设定投影方向数目为32，迭代次数为10，对上述仿真信号用BEMD方法进行分解，共得到5阶BIMF，然后根据文献[7]的分量筛选原则，计算出各阶分量与原信号的相关系数(表1).将相关性较低的分量(BIMF3、BIMF4、BIMF5)归入趋势项，得到的最终分解结果如图2所示.

表1 各阶BIMF分量与原始信号的相关系数

图2 仿真信号及筛选后的BIMF

3 颤振特征值序列及其控制图

3.1 颤振特征值序列的提取

方差作为描述随机变量偏离期望程度的标志，能有效地反映时间序列变化程度.在磨床振动信号诊断中，实时方差可以反映离散随机信号在期望值附近的振荡及其总体分布情况.因此，实时方差可作为一个理想的特征向量来监测颤振的发生情况.

将筛选出来的两阶BIMF的实部和虚部每隔8个点计算一次实时方差，并记实部为Rs，虚部为Rx，然后，通过式(7)分别求解两阶有效BIMF的实时方差特征量Rv1、Rv2和总的特征量RV，再由式(8)得到总的实时方差特征值序列(图3).

(7)

(8)

图3 仿真信号颤振总的实时方差特征值序列

3.2 控制图的建立

以磨削稳定阶段的实时方差作为特征值，构造一个时间序列{R(t)}.采用稳定磨削状态下25～74个数据点，即对应0.25～0.74 s时间段的50个特征值，作为统计特性值进行分析.

首先，对时间序列进行单一性、正态性和相关性检验.这些特征值都是针对同一信号的，能确保单一性.其次，采用常用的正态性检验方法中的Q-Q图(Quantile-Quantile Plot)对数据点进行正态性说明(图4).图4中纵坐标代表标准正态分布分位数，横坐标代表样本分位数，直线斜率代表标准差.若数据点靠近图中的直线，则满足正态分布.

由图4可以看出，数据点靠近直线分布，故该特征值序列满足正态性.采用自相关图观察的结果如图5所示.由图5可以看出，自相关系数(ACF)大致处于±0.5范围内，相关性一般，符合控制图对统计特性值相关性的要求.

图4 正态性检验

图5 自相关图

完成上述检验后，可以计算出平均值μ和标准差σ，得出控制图的上、下控制限(表2).

表2 仿真信号控制图参数

注：表中的量均无量纲.

已知均值、标准差和上下控制限，就可以得到磨削稳定阶段特征值序列的控制图(图6).

图6 磨削稳定阶段特征值序列的控制图

4 预测结果处理

应用上述控制图进行预测时，需将实时方差特征值序列数据按时间顺序输入控制图中，并采用控制图判异准则实时对数据进行异常波动检验.检验结果如图7所示.

图7 仿真信号预测检验结果

从图7可以看出，仿真信号大致在第1.11 s时出现了异常.于是，将该时间点前后的数据列入了表3.

表3 特征值数据

注：表中特征值对应的上控制限UCL=4.634 50；特征值为无量纲量.

由表3可见，第1.12 s时有特征值掉落至上控制限外，同时，从1.07～1.13 s这0.07 s内，特征值序列呈现持续增长趋势.由于数据时间间隔为0.01 s，故上述时间段内有7点连续增长，由上文准则①、③可得到结论：磨削在1.12 s时出现异常，可能发生颤振，应及时采取保护措施，使磨床正常工作.

5 试 验

本文参考文献[7]，以杭州杭机股份有限公司的KD4020X16数控动梁龙门导轨磨床为试验对象，搭建磨床振动检测试验平台，采集磨床颤振信号.试验中选取1.25～31.25 s之间30 s的数据，用BEMD方法对采集到的颤振数据进行处理(图8)，并求解采集信号的实时方差，作为特征量.

图8 试验信号及BEMD分解结果

选取正常磨削状态1～5 s间的实时方差数据，验证时间序列的正态性和相关性后，计算出平均值和标准差，得出控制图的上下控制限(表4)，根据控制图参数建立控制图模型.最后，将特征值按时间顺序输入模型中，得到图9所示的控制图预测结果.

表4 试验信号控制图参数

注：表中的量均无量纲.

图9 试验信号的控制图预测结果

从图9可以看出，试验信号特征值在5 s处开始持续增大，到5.35 s时超出上控制限，故认为磨床在此时发生颤振.

6 讨 论

尽管控制图在磨床颤振研究领域已经有过应用，但几乎都以一维振动信号为对象，并且未对信号进行合适的预处理，而实际生产中磨床振动信号具有多维性、非平稳性和非线性的特点.为解决实际振动信号的上述问题，本文采用BEMD与控制图相结合的方法，对二维非平稳仿真信号进行颤振预测，经过BEMD分解并筛选后得到有效的两阶BIMF和一个趋势项，对实时方差进行控制图建模，有效预测了颤振.

预测结果表明，1.12 s时出现颤振，仿真信号颤振发生时间为1 s，两者相对误差在12%，预测精度处于较高水平.通过验证，实测信号大致在6 s处振动幅值开始增大，进入颤振状态.这与控制图预测的结果5.35 s基本吻合.

本文方法有效去除了信号中的无效成分，利用控制图对颤振进行预测，并给出了颤振阈值.但是，控制图对颤振信号的正态性有较高的要求，在应用范围上存在一定局限，并且本文方法没有解决模型自适应调整的问题，预测精度受到了影响，因此有待于深入研究和改进.

7 结 论

(1)采用BEMD方法可以对二维颤振信号进行分解，得到全面反映颤振信息的特征模式函数，获取更合理的颤振特征量，克服了传统信号分析方法在二维非平稳信号上应用的局限性，去除了混杂在原始信号中的干扰成分，提高了预测的精度.

(2)控制图的制作过程简单，并且控制图可给出传统预测方法不能给出的颤振发生阈值.

(3)采用BEMD与控制图相结合的方法对二维非平稳信号进行颤振预测，比单纯运用控制图预测更科学，精度更高，对颤振预测研究具有一定的指导意义.

[1] 王立刚，刘习军，贾启芬. 机床颤振的若干研究和进展[J]. 机床与液压，2004(11):1-5.

[2] Bediaga I,Munoa J,Hernández J,et al. An automatic spindle speed selection strategy to obtain stability in high-speed milling[J]. International Journal of Machine Tools & Manufacture,2009,49(5):384-394.

[3] 吴 石,林连冬,肖 飞,等. 基于多类超球支持向量机的铣削颤振预测方法[J]. 仪器仪表学报,2012,33(11):2414-2421.

[4] Karpuschewski B,Wehmeier M,Inasaki I.Grinding monitoring system based on power and acoustic emission sensors[J].CIRP Annals Manufacturing Technology,2000,49(1):235-240.

[5] Tansel I N,Wang X,Chen P,et al. Transfor-mations in machining(Part2):Evaluation of machi-ning quality and detection of chatter in turning byusing stransformation[J]. International Journal of Machine Tools & Manufacture,2006,46(1):43-50.

[6] 王 民,刘国付,昝 涛,等. 基于控制图的磨削颤振预测方法[J]. 北京工业大学学报,2015,41(9):1296-1301.

[7] 沈建洋.基于BEMD和LSSVM的大型磨床磨削颤振在线检测方法研究[D].杭州：浙江理工大学,2015.