做滚动轴承故障诊断这些年我遇到的第一个“劝退”瞬间是看着采回来的振动信号不知道从哪里下手。传感器贴在轴承座上采到的从来都不是单纯的轴承振动而是转频、保持架游动、结构共振、工频干扰和噪声搅在一起的混合物。早期故障的冲击能量本来就弱这些成分再一叠加频谱上连个像样的故障特征频率尖峰都找不到。于是拆信号就成了绕不开的环节。VMD变分模态分解是我这几年用下来最顺手的信号分解工具它能把原始振动信号自适应地分解成多个模态每个模态对应一个相对窄带的调幅调频分量然后再对选中的模态做功率分解和包络解调故障位置和严重程度就清楚多了。这篇内容适合两类人一是刚接触VMD看着K、alpha这些参数不知道从哪下手的同学二是被EMD的模态混叠、端点效应反复折磨想换个稳定方案的老手。我会从轴承信号特征讲起把VMD的参数机制、分解流程、功率分析、故障判别以及调试中踩过的坑全部梳理一遍。看完这套流程你自己也能搭出一个十分钟出结论的轴承诊断模块。1. 滚动轴承故障信号为什么难诊断1.1 轴承振动信号的主要成分与故障特征频率滚动轴承出现局部损伤后每当损伤点经过承载区或者与其他表面接触就会产生一次短时冲击。这个冲击会激励起轴承内外圈和传感器安装结构的高频共振表现为一段衰减振荡同时这种冲击按旋转周期重复出现又会形成对低频转频信号的调制。从信号上看就是一组高频衰减脉冲串压在一个低频周期性包络下面典型的调幅结构。这类信号最麻烦的地方在于频带重叠。转频分量、故障特征频率分量、共振频带成分、随机噪声堆在一起用普通带通滤波器很难一刀切开。而且早期故障的冲击幅度非常微弱很可能比噪声低一个数量级直接看时域波形就是一团毛刺。要判断故障位置必须用到轴承的特征频率公式。假设滚动体数量为N转频为fr滚动体直径为d节圆直径为D接触角为α那么频率公式如下外圈故障特征频率BPFO (N/2) × fr × (1 - d/D × cosα)内圈故障特征频率BPFI (N/2) × fr × (1 d/D × cosα)滚动体故障特征频率BSF (D/d) × fr × (1 - (d/D × cosα)²)保持架故障特征频率FTF (fr/2) × (1 - d/D × cosα)举一个具体例子转频fr30Hz滚动体数量N10d/D0.2接触角接近0°那么BPFO约等于120HzBPFI约等于180HzBSF约等于144HzFTF约等于12Hz。注意BPFI和BSF只差了36Hz在频谱上如果分辨率不够或者有边带干扰很容易混淆。这就是为什么不能只看一个峰值而要把分解、功率谱、包络谱结合起来判断。1.2 为什么不用EMD而选VMD很多初学者最初接触的分解工具是EMD经验模态分解。不能说它没用但我在实操中被它坑过几次主要是模态混悬和端点效应。举一个我遇到的真实场景一段含噪声的轴承信号EMD分解出来的第一个IMF有时候是高频冲击有时候是共振分量换个数据段结果又变了频率相近的成分会在不同IMF之间跳来跳去这就是模态混叠。另一个问题是信号两端的包络会突然甩出一条大尾巴也就是端点效应导致分解出的模态在边界上不可信只能裁掉一大段数据。VMD是Dragomiretskiy和Zosso在2014年提出的变分模态分解方法它把分解过程变成了一个约束优化问题要求所有模态之和等于原始信号同时每个模态都是围绕某个中心频率的有限带宽分量并且让所有模态的带宽总和最小。我认为VMD在实际应用中有三个看得见的优势模态数量是事先指定的K不会像EMD那样飘忽不定每个模态有明确的中心频率和带宽概念可解释性强对噪声相对稳健不会因为一点扰动就让模态面目全非。对比项EMDVMD数学框架经验筛选无统一优化目标变分约束优化问题模态混叠常见噪声下尤其严重可控靠参数调节缓解端点效应明显需要大量延拓较轻但仍需留意模态数量自动生成可能不稳定需预先指定K抗噪能力一般信号稍脏就全乱通过alpha约束带宽相对稳定参数干预几乎无法干预K、alpha、tau等均可调整当然VMD也有自己的短板最典型的是K必须人为指定而且参数设置不当会得到完全不同的分解结果。不过这个问题有规律可寻下面两章我详细讲。2. VMD的核心原理把信号拆成多个模态的变分框架2.1 变分问题怎么理解VMD的数学过程听起来有点劝退但理解起来并不难。你可以把原始信号想象成一桌菜VMD要做的是把荤素、冷热、咸甜按一定规则分到K个盘子里每个盘子的口味范围尽量窄所有盘子合在一起正好等于原来的整套菜。对应到信号上每个“盘子”就是一个本征模态IMF它有自己相对集中的频率范围并且是个调幅调频分量。用公式描述的话VMD要解决的问题是对每一个模态信号先通过Hilbert变换构造解析信号再乘一个和中心频率相关的指数项把它移到基带然后对时间求导并计算L2范数得到该模态带宽的估计。算法目标是最小化所有模态带宽之和同时保证所有模态加回去等于原始信号。我摘两个核心更新公式方便你理解alpha的作用。模态在频域的更新形式可以写成ûk(ω) (f̂(ω) - Σ_{i≠k} ûi(ω) λ̂(ω)/2) / (1 2α(ω - ωk)²)中心频率的更新形式是ωk ∫ ω|ûk(ω)|²dω / ∫ |ûk(ω)|²dω第一个公式的分母里出现了alpha和频率差alpha越大分母越大意味着偏离中心频率的分量被压得越厉害模态带宽就越窄。所以alpha在中文资料里常被称为“惩罚因子”或“带宽约束参数”它的实质是决定每个模态可以占据多宽的频率范围。如果alpha取太大模态会窄成一根细线有些真实成分直接被切掉alpha取太小模态带宽过宽又容易把噪声和相邻频段一起卷进来。实际求解VMD用的是ADMM交替方向乘子法简单说就是轮流更新模态、中心频率和拉格朗日乘子在迭代中让这三组变量逐步收敛到最优解。这个过程不需要用户干预收敛精度由tol控制。我刚接触VMD时也想过要不要自己推一遍证明后来发现工程上先把参数行为摸清楚应用反而跑得更快。2.2 每个参数到底在控制什么VMD的输入参数一共就几个K、alpha、tau、DC、init、tol。每个参数都有明确的物理含义参数表如下参数作用常见取值调整思路K模态个数决定分几盘3-8从5起步过小欠分解过大会把单个分量切碎alpha带宽惩罚因子决定每盘多窄2000噪声大时调大到3000-10000共振带宽宽时调小tau更新步长和噪声容忍相关0一般保持0大噪声环境可试0.1但要观察重构误差DC是否强制第一模态为直流分量0交流信号一般设0init中心频率初始化方式1频率均匀初始化收敛稳定tol迭代收敛容差1e-7追求速度可放宽到1e-6DC这个参数可能有人忽略。它设成1时第一个模态会被强制保留为信号的直流分量。轴承振动分析关注的是动态冲击成分直流不是重点所以我通常设0让直流成分自然被分到某个低频模态里或者直接在预处理时把均值去掉。init参数我建议固定用1也就是让中心频率在频带内均匀初始化。有些同学喜欢试init0想看看会不会有“新发现”实测下来差别不大但收敛速度可能变慢。既然轴承故障诊断追求的是可复现就不需要在这个参数上折腾。2.3 K值选择中心频率收敛判断法K到底取多少这是VMD使用中最头疼的问题。我给一个我最常用的实操方法叫“中心频率收敛判断法”。步骤是这样的固定alpha2000从K3开始逐步增加到K7或者8每次分解完都把各模态最终收敛的中心频率记录下来。然后看两件事第一不同K值下中心频率是否稳定。如果K从4增加到5前4个中心频率基本没动新增的第5个模态中心频率离第4个也不是特别近说明K5是合适的。第二如果某次增加K之后出现两个模态的中心频率非常接近比如都挤在1800Hz附近甚至频谱大幅重叠说明K已经取大了应该退回上一个值。另一个辅助判断是看残差。VMD分解完把原始信号减去所有模态之和得到残差信号。对残差做FFT如果还能看到明显的窄带尖峰说明K取少了这台机器里确实还有成分没被拆出来如果残差基本是平坦噪声那K就基本合适。对于滚动轴承信号我一般从K5起步。大多数轴承故障信号里低频转频、故障调制成分、结构共振、噪声背景四五个主导频带就能覆盖。K取得过大不仅浪费算力还容易出现虚假模态后患比K取小更麻烦。3. VMD滚动轴承故障诊断实操全流程3.1 数据采集和预处理的基本准则再好的分解算法也经不住烂数据。轴承振动采集这块我有几条固定遵循的规则。采样率要留足余量。滚动轴承故障特征频率通常不高但冲击激励出的高频共振往往在几千赫兹甚至更高因此采样率建议不低于12.8kHz常规情况下用25.6kHz。Nyquist定理只保证不混叠工程上要给共振频带和滤波器过渡带都留位置。采样时长要覆盖足够多的故障周期。比如转速3000转/分保持架故障频率只有十几赫兹采样太短根本看不出来。我一般要求至少采集10秒以上确保最低特征频率也有几十个周期可以统计。传感器安装位置也很关键尽量靠近轴承承载区贴在轴承座或壳体刚度较高的位置。安装面要打磨平整用黏合剂或磁座固定避免松动引入额外的非线性失真。预处理环节我比较克制去除均值、去除趋势项就够了。不建议在VMD之前做太激进的带通滤波因为VMD本来就是用来做自适应频带划分的提前把频带卡死反而失去了它的意义。如果信号里存在明显的工频干扰可以在分解前用陷波滤波器先滤掉随后再让VMD处理剩余成分。3.2 VMD分解的主流程与代码完整流程概括成四步读入信号去均值设定参数调用VMD分解检查残差和中心频率。下面是我最常用的一段代码基于Python的vmdpy库安装方式很简单pip install vmdpy分解主流程import numpy as np from vmdpy import VMD fs 25600 # 采样率 data np.loadtxt(vibration.txt) # 读入振动信号 data data - np.mean(data) # 去均值 K, alpha 5, 2000 # 参数起点 tau, DC, init, tol 0, 0, 1, 1e-7 u, u_hat, omega VMD(data, alpha, tau, K, DC, init, tol) # u: K行N列的数组每一行是一个模态信号 # omega: 每个模态收敛后的中心频率Hz或rad/s注意单位 residual data - np.sum(u, axis0) # 残差信号 residual_power np.std(residual) / np.std(data) print(残差占比:, residual_power)这里有个单位陷阱vmdpy返回的中心频率和采样率有关不同版本、不同调用方式可能返回归一化频率或角频率。我习惯在分解后把中心频率打印出来和原信号的FFT频谱峰值对照一下如果中心频率落在明显的主峰位置附近说明分解基本合理。残差占比控制在0.05以内也就是残差标准差不到原始信号标准差的5%这个分解就基本完整。参数初始值方面我建议K先取5alpha取2000tau取0这也是文献和库默认值。第一轮跑完再根据第2.3节的方法微调不追求一次到位。3.3 模态功率分解从时域模态到故障特征标题里提到的“功率分解”在实际操作中就是基于VMD分解得到的多个模态再做功率谱分析把信号的能量按不同频率模态拆开找出故障信息集中的频带。这一步很有必要因为VMD输出的是时域波形直接观察很难看出名堂转到频率域才能定位特征。我通常两步走。第一步对每个模态用Welch方法计算功率谱密度第二步计算每个模态的能量占比用于判断哪些模态值得继续分析。from scipy.signal import welch def modal_power(imf_list, fs): n_modes len(imf_list) powers [] ratios [] total 0.0 for imf in imf_list: freq, pxx welch(imf, fsfs, npersegmin(4096, len(imf))) total_power np.sum(pxx) powers.append((freq, pxx, total_power)) total total_power for freq, pxx, p in powers: ratios.append(p / total) return powers, np.array(ratios)跑完打印各模态能量占比你会看到一种典型分布低频模态占一部分能量背景噪声模态占一部分而包含故障冲击的模态往往表现为能量占比不高、但频谱上有清晰的特征频率尖峰。所以不能只选能量占比最高的模态还要看频谱成分。我举一个实际算例某外圈故障信号分解成5个模态能量占比分别是8.5%、31.2%、27.4%、22.6%、10.3%。前三个模态频谱都含明显的120Hz或240Hz峰值但第一个模态还混着大量转频能量第三模态在2800-3200Hz共振频带有突出峰包络谱反而最干净。所以我的习惯是先把所有模态功率谱都画出来挑出那些在特征频率附近出现独立尖峰的模态再做包络解调而不是机械地看谁能量大。3.4 包络解调与故障频率判别包络解调是轴承故障诊断的标准动作。故障冲击产生的是调制信号直接FFT只能看到共振频带附近的高频簇很难直接看到调制频率。Hilbert变换可以把信号的瞬时幅值包络解出来对这个包络做FFT才能在低频区看到特征频率。对选中的模态做包络谱的代码很简单from scipy.signal import hilbert def envelope_spectrum(imf, fs): env np.abs(hilbert(imf)) env env - np.mean(env) freq np.fft.rfftfreq(len(env), 1 / fs) amp np.abs(np.fft.rfft(env)) / (len(env) / 2) return freq, amp拿到包络谱之后对着预先算好的特征频率表做核查注意看基频、谐波和边带故障位置包络谱核心特征判别要点外圈故障BPFO处有明显峰且其2倍、3倍频逐渐衰减峰值稳定转频边带不明显内圈故障BPFI处有明显峰峰两侧出现fr边带边带间隔约等于转频fr滚动体故障BSF处有明显峰可能伴随FTF谐波丰富波形幅值波动大保持架故障FTF处低频峰幅值相对小需要较长数据确认回到前面那个数值例子fr30HzBPFO120Hz。如果包络谱在120Hz、240Hz、360Hz都出现峰且谐波间隔正好是120Hz基本可以判定为外圈故障。如果包络谱在180Hz附近有峰旁边还有150Hz和210Hz两个小峰边带间隔是30Hz那就是内圈故障的典型特征因为内圈故障的冲击要随旋转周期性改变承载区位置转频形成了调制边带。这里要提醒一句在包络谱上看到峰值不代表就确诊了还要看谐波结构、边带间隔和时域波形中的冲击间隔是否一致。我见过有人仅凭一个单峰就报故障最后发现是其他机械部件的固有频率教训很值得记住。4. 调试VMD时最常踩的坑4.1 K与alpha的联动欠分解还是过分解参数调试时最常见的问题是K和alpha互相拉扯。K取小了模态带宽被迫变宽一个模态里可能塞进去两个不同来源的频率成分导致后面包络谱上一堆假峰K取大了单个真实分量会被切成两个相邻模态两个模态中心频率几乎贴在一起输出高度相关。alpha过大会怎样模态带宽被压缩得过窄冲击共振频带被截断能量漏到残差里残差占比明显升高。alpha过小呢模态带宽宽松噪声也被一起包进来包络谱基线抬高小峰值被淹没。我的调节经验是先用alpha2000跑一遍根据中心频率冲突情况调K如果调完K之后残差仍然偏大再回头调alpha。具体来说强噪声环境下alpha可以调到3000到10000但要时刻盯着残差残差超过5%就说明alpha压得太狠了。如果信号带内存在宽频共振alpha反而要往1000以下调让模态把整个共振频带覆盖住。实际测试中我遇到过alpha10000加上K6的组合五个模态都在窄频带内打转故障冲击的主要成分全跑到了残差里包络谱干净得过分反而贻误了诊断。所以调参的终点不是“好看”而是物理意义清晰。4.2 端点效应和模态混叠的现场处理VMD的端点效应虽然没有EMD那么严重但在长数据两端仍然会出现可见的振荡偏离。处理方式我推荐“镜像延拓”在信号两端各自延拓一段数据分解完成后把延拓部分裁掉。具体长度视模态波长而定我一般延拓200到500个点延拓太多会引入额外计算延拓太少又起不到稳定作用。模态混叠在VMD里多半是参数设置问题。如果发现两个模态在频谱上有明显交叠先不要怀疑算法按三步排查K是不是取大了alpha是不是太小数据里是不是存在超出常规的动态范围。如果两步参数调整都无效再考虑用逐次分解策略也就是先对原信号做一次VMD选出感兴趣的模态对该模态再做一次VMD相当于把宽频带逐级拆细。还有一个容易被忽略的问题分段处理时的边界效应。如果数据很长从中间截出几段分别做VMD每段分解结果可能不完全一致。为了保持一致性我会让段与段之间有50%的重叠分析时优先采用中间段的结果。这个方法比单靠一端数据要稳得多。4.3 判断分解质量的三条硬指标被问最多的一个问题是“我怎么知道这次分解是好是坏”。我总结了三条硬指标分享出来作为参考。第一中心频率是否稳定收敛。跑完分解把omega打出来正常的模态应该有一个明确的收敛值并且重复运行结果一致。如果同一段数据连续跑两次中心频率漂移明显说明初始化或参数设置有问题。第二频谱分离度。画出每个模态的功率谱中心频率位置应该有清晰的单峰相邻模态的-3dB频带边缘尽量不重合。重合度越高分解质量越差。第三残差占比是否低于5%。把原始信号减去所有模态之和求残差标准差与原始信号标准差之比。小于5%算正常大于10%基本可以判定分解失效。但要注意残差占比也不是越低越好降到接近0意味着算法把噪声也强行拆成了“模态”反而会带来大量虚假分量。我每次跑完VMD都会顺手打印这三项指标。这组数据积累到几十条之后你会慢慢建立一种直觉看到中心频率和残差占比就知道该往哪个方向调参数。5. 参考代码、参数速查与个人建议5.1 可直接改用的Python调用代码最后给一套能直接改用的完整脚本。它在前面代码的基础上加入了能量占比排序和包络谱输出基本覆盖了从原始数据到故障频率判别的全过程import numpy as np from scipy.signal import hilbert, welch from vmdpy import VMD def vmd_analysis(data, fs, K5, alpha2000): data data - np.mean(data) u, u_hat, omega VMD(data, alpha, 0, K, 0, 1, 1e-7) residual data - np.sum(u, axis0) residual_ratio np.std(residual) / np.std(data) print(中心频率(Hz):, omega * fs / (2 * np.pi) if omega.max() 6.28 else omega) print(残差占比:, residual_ratio) total_power 0.0 powers [] for k in range(K): freq, pxx welch(u[k], fsfs, npersegmin(4096, len(u[k]))) p np.sum(pxx) powers.append((freq, pxx, p)) total_power p for k in range(K): ratio powers[k][2] / total_power peak_freq powers[k][0][np.argmax(powers[k][1])] print(f模态{k1}: 能量占比 {ratio:.3f}, 功率谱峰值频率 {peak_freq:.1f} Hz) return u, omega, residual调用时把你的震动数据和采样率传进去即可。如果主频分辨率不够可以适当调大nperseg代价是计算变慢。5.2 参数速查表与故障特征对照表这里整理两个速查表方便实际干活时快速对照。参数速查表场景K建议alpha建议备注常规滚动轴承信号无强噪声52000先跑基准结果强噪声环境4-63000-5000监控残差占比防过压宽频带共振明显5-71000-1500让模态覆盖完整共振峰高速小型轴承3-42000-3000频带范围窄K不宜大低速大型轴承6-82000-2500低频成分多适当增加K故障特征对照表特征频率包络谱表现补充判断BPFO基频和谐波间隔均为BPFO边带少幅值平稳BPFI基频两侧出现fr边带边带间隔约等于转频BSF基频和谐波丰富伴随FTF幅值波动大FTF低频峰值需要长数据验证用的时候先按公式计算理论特征频率再对照包络谱峰值不要反着来。最后说点我调VMD的真实体会。刚上手时我总想把参数一次调对后来发现这个工具更像一个需要反复对话的过程先设一组起点参数分解完看中心频率、能量占比和残差再决定加模态还是减模态。记住一句话VMD只负责把信号拆开拆得对不对、哪一堆才是故障信息最终还是要靠你对故障机理的熟悉程度。我习惯把每次实验的K、alpha、中心频率、残差能量记在同一张表里样本多了之后自然会有一种手感。这套流程跑顺之后从原始振动数据到给出诊断结论基本上十分钟以内就能完成。
