ARTICLE DETAIL

资讯详情

深耕网站建设与运营推广的一线实战洞察。

传递路径分析(TPA)在齿轮箱故障诊断中的原理与Matlab实现

传递路径分析(TPA)在齿轮箱故障诊断中的原理与Matlab实现 去年处理一台矿山破碎机齿轮箱的振动异常时我盯着测点屏幕上1300Hz左右一群密密麻麻的边带发了很久的呆。齿轮啮合频率的谐波和疑似轴承特征频率重叠在一起时域波形上还能看到明显冲击谁都不敢拍板说问题在齿轮还是轴承。拆机验证的结果是齿轮完好输出端轴承保持架疲劳断裂——测点信号狠狠误导了我们。问题出在哪出在我们把测点信号直接当成了源信号。齿轮系统的信号从啮合点传到传感器中间经过轴、轴承、箱体壁每走一段就被“过滤”一次最终测到的是多条传递路径叠加后的结果。想在这团混合信号里做故障诊断先要把“信号从哪来、走哪条路、贡献有多大”搞清楚这正是传递路径分析TPATransfer Path Analysis做的事。这篇文章围绕齿轮系统的TPA展开从原理拆解到Matlab代码实现再到实测中的坑和故障闭环验证都是我实际跑项目时沉淀下来的做法。适合做设备状态监测的工程师、机械故障诊断方向的研究生以及所有被“测点信号不好使”折磨过的人。1. 齿轮箱信号为什么测不准传递路径在背后“搅局”1.1 从“测点信号”到“源信号”中间隔着多个环节齿轮箱的振动信号从产生到被传感器记录下来通常经过这样一条物理链路齿轮啮合点激励 → 齿轮本体 → 轴系 → 轴承 → 轴承座 → 箱体壁板 → 测点加速度计这条链路里的每一级都是质量-刚度-阻尼构成的机械系统本质上就是一个滤波器。齿轮啮合产生的宽频激励经过这些环节后不同频段的能量被重新分配有的频段被放大有的频段被严重衰减。我见过不少工程师拿着测点频谱直接做齿轮故障判断结果低频段的异常能量被误判为“不平衡”高频段的齿轮故障特征被箱体衰减后几乎消失诊断结论自然就跑偏了。最典型的例子是齿轮齿面点蚀。点蚀冲击在源处是宽频信号但经过轴承滚子与滚道的高频阻尼、箱体壁板的弯曲模态过滤之后测点处能留下的往往只剩啮合频率及其低阶谐波冲击特征被压得几乎看不见。这时候你不管做包络谱还是小波分析得到的结果都会大打折扣。这就是传递路径的第一个作用它改变了信号的频域形态让“测点所见”不等于“源处所有”。1.2 多路径叠加边带不仅来自故障源齿轮箱内部并不只有齿轮一个源。电机电磁激励、轴系不平衡、不对中、轴承早期故障、齿轮啮合冲击甚至负载端的周期性波动这些源产生的振动同时沿着各自路径传到同一个测点最终测点信号是它们按矢量和叠加的结果。这个叠加不是功率相加而是带相位的矢量叠加。两个频率相同、相位相反的信号甚至可能相互抵消。现场最容易踩的坑就在这里某个测点处的边带结构看起来像齿轮故障调制实际可能是轴承故障信号沿着某条路径与齿轮啮合信号叠加后形成的“伪边带”。举个例子。二级齿轮箱输出端轴承出现内圈故障时故障特征频率会调制到该轴转频及其倍频上而如果齿轮正好也有一点轻微磨损啮合频率两侧也会出现转频调制。两种边带在频谱上看形态相似距离相近单靠测点频谱很难区分。只有搞清楚“输出端轴承这一路传递函数如何、齿轮啮合这一路传递函数如何、两者各自贡献了多少”才能把边带归属到具体源上。这就是TPA要解决的核心问题把测点响应拆解为“若干源×若干路径”的贡献之和让你能看到到底是哪条路径、哪个源在主导测点处的振动能量。2. TPA在齿轮系中的数学原理与载荷识别逻辑2.1 路径贡献的“标准公式”怎么理解TPA的理论基础是线性时不变系统下的叠加原理。假设有 p 个等效源力或力矩m 个目标测点那么在某一角频率 ω 下第 k 个测点的响应可以写成y_k(ω) Σ_i H_ik(ω) · F_i(ω)其中 H_ik 表示第 i 个源到第 k 个测点之间的频率响应函数FRFF_i 是第 i 个源处施加的等效载荷。这个公式看起来平淡但它包含了TPA最核心的思想每个测点响应不是某个源的“独奏”而是所有源经过各自“音箱”路径播放后的“混响”。如果你想评估某条路径对最终振动的贡献就取 H_ik·F_i 这一项所有路径贡献累加就重构出测点响应。理解这个公式可以用一个生活化的类比水池里装了六个水龙头每个水龙头的管路长度、粗细、弯头数量都不同。你站在池边测量水位变化想知道“是哪个水龙头往里注水注得最多”只测总水位是没用的必须知道每个管路的“水力特性”相当于 H_ik以及每个龙头当前的“开度”相当于 F_i。2.2 齿轮系中的“载荷”为什么必须反求理论上拿到 FRF 矩阵和源载荷 F_i就能算出所有路径贡献。但齿轮系统的特殊难点在于源载荷几乎测不到。齿轮啮合点内部的动态啮合力无法直接安装传感器测量轴承内部滚动体与滚道之间的接触力也测不了。我们能测的只有箱体表面、轴承座这些位置的响应。所以必须走逆问题的路线——先通过试验或仿真得到 FRF 矩阵再由实测响应反推等效源载荷。载荷识别的数学形式很直接F_est(ω) H(ω)^ · y_实测(ω)H^ 是 FRF 矩阵的伪逆。听起来简单但这一“反推”恰是整套TPA流程中最容易出问题的地方。我见过不少刚入坑的同行在第一步就翻车直接把 H 矩阵求逆然后拿实测信号一乘得到的载荷结果千奇百怪——今天算出来齿轮啮合力是几千牛明天同一工况下算出来只有几十牛。原因就是 H 矩阵严重病态。测点都布在同一个箱体上各测点 FRF 之间高度线性相关矩阵接近奇异微小噪声在求逆时被放大了几个数量级。这种“放大”会让载荷识别结果完全失真后续路径贡献自然也是错的。2.3 “把路线反过来走”等效载荷识别的思路这里需要明确一个概念TPA 里识别的“源载荷”并不是物理意义上真实的啮合力而是一个“等效载荷”。它是在给定 FRF 不变的前提下为了重构测点响应而反推出来的等效激励。这个等效载荷包含了真实源的信息也吸收了模型误差和路径误差所以不要对它做过度物理解读关键是看它结合路径后能多大程度还原目标响应。在齿轮箱实操中我一般这样处理载荷识别问题先布设足够测点测点数必须大于等于等效源数否则方程组欠定。获取 FRF 矩阵可以用锤击试验实测也可以用有限元模型计算甚至用运行状态下被动估算的 OPA 方法。规划阶段用实测 FRF 更稳。对 FRF 矩阵做奇异值分解观察奇异值谱判断有效秩确定“真正起作用”的等效源数量。采用正则化方法做逆运算而不是直接求伪逆。3. Matlab代码实现频响函数矩阵、逆运算与路径贡献可视化3.1 用集中参数模型生成系统FRF矩阵工程项目的头一版我不会上来就上有限元而是先用一个多自由度集中参数模型把 TPA 的逻辑跑通。这个模型把齿轮箱简化为若干惯性质量和弹性连接环节虽然不能精确模拟结构模态但用来验证载荷识别算法、训练团队对 TPA 的理解效率很高。假设简化系统有 6 个自由度对应输入端轴承座、齿轮啮合点、输出端轴承座等关键位置。质量矩阵、刚度矩阵、阻尼矩阵都是 6×6% 集中参数模型6自由度齿轮箱简化 % 自由度对应1输入端轴承 2输入轴齿轮 3中间轴齿轮 4输出轴齿轮 5输出端轴承 6箱体测点 n_dof 6; M diag([2.5, 3.2, 8.5, 4.0, 2.8, 15.0]); % 质量 kg % 刚度矩阵简化对角耦合项 K [8e6 -3e6 0 0 0 -1e5; -3e6 9e6 -4e6 0 0 0; 0 -4e6 1.2e7 -5e6 0 0; 0 0 -5e6 1.1e7 -4e6 0; 0 0 0 -4e6 8e6 -1e5; -1e5 0 0 0 -1e5 3e5]; % 比例阻尼 C alpha*M beta*K alpha 5; beta 2e-6; C alpha*M beta*K; % 频率轴 freq (0:0.5:2000); H zeros(n_dof, n_dof, length(freq)); for fi 1:length(freq) omega 2*pi*freq(fi); A -omega^2 * M 1i*omega*C K; H(:, :, fi) inv(A); % 每列 单位激励下各测点的响应 end循环里对每个频点做一个矩阵求逆就是整个 FRF 矩阵的生成过程。实际工程中这条 H 是通过锤击试验测出来的原理一模一样只是矩阵数据来自现场更可靠。3.2 载荷识别与路径贡献计算核心代码拿到 H 矩阵之后载荷识别是最关键的一步。直接伪逆的代码很简单但前面说过容易病态。我习惯用带 Tikhonov 正则化的伪逆替代普通伪逆% 载荷识别以某一频率点为例 fi 500; % 对应 500Hz 左右的频点 H_k squeeze(H(:, :, fi)); y_meas [0.12; 0.18; 0.25; 0.31; 0.42; 0.15]; % 实测响应的频域幅值示例 % 直接伪逆 F_direct pinv(H_k) * y_meas; % Tikhonov正则化伪逆 [U, S, V] svd(H_k); s diag(S); lambda 0.05 * max(s); % 正则化参数后面专门说怎么调 s_reg s ./ (s.^2 lambda^2); H_inv_reg V * diag(s_reg) * U; F_reg H_inv_reg * y_meas; % 计算各路径贡献 contrib zeros(n_dof, 1); for i 1:n_dof contrib(i) H_k(:, i) * F_reg(i); % 第i路径贡献的复数值 end % 各测点总响应重构 y_reconstruct sum(contrib, 1);正则化背后的物理直觉并不复杂当某条奇异值远小于主奇异值时它对应的特征方向放大噪声严重。直接求逆会让这个方向数值爆炸加上 lambda 之后相当于给“可疑方向”的放大倍数封了顶以微小拟合误差为代价换来载荷识别结果的稳定。我跑过大量仿真和实测数据lambda 取最大奇异值的 1%10% 通常能兼顾稳定性和拟合精度。这个比例不是固定不变的后面 4.2 节会展开讲如何用 L 曲线选得更准。3.3 贡献度可视化与报告输出TPA 算完不是只给自己看得让现场工程师和管理层一眼看懂。我最常用的是两类图路径贡献占比柱状图、路径贡献随频率变化的曲线图。% 路径贡献占比取幅值占比 contrib_amp abs(contrib); contrib_percent contrib_amp / sum(contrib_amp) * 100; figure; bar(contrib_percent, FaceColor, [0.3 0.6 0.9]); xlabel(路径编号); ylabel(贡献占比 (%)); title(sprintf(TPA路径贡献分析 %.1f Hz, freq(fi))); grid on; set(gca, XTickLabel, {输入轴承,输入齿轮,中间齿轮,输出齿轮,输出轴承,箱体测点}); % 保存报告图片 saveas(gcf, TPA_contrib.png);这段代码生成的柱状图信息非常直观哪条路径贡献大一眼就能看出来。实际做报告的时候我会再补一张“各路径贡献 vs 频率”的二维图横轴频率、纵轴贡献幅值多条曲线叠加这样能看出不同频段由哪条路径主导。诊断固件做久了你会明白TPA 结果不只是拿来定位单一故障的更多时候是拿来看“哪个频率段、哪条路径、在什么工况下占主导”这是后续做状态监测传感器布点和报警阈值设计的依据。4. 实测踩坑记录路径越多、矩阵越病态4.1 矩阵病态的本质与现场快速判断TPA 在实验室跑得漂亮一到现场就容易露馅。第一个坑就是 FRF 矩阵病态。原因并不神秘齿轮箱是个连续金属结构测点之间距离又不远各点 FRF 有很强的同源性——都是同一条结构路径在不同位置取点。反应在数学上就是 H 矩阵各列高度相关奇异值差距巨大。如果你设了 10 个测点、10 个源但实际有效奇异值只有 3 个那这个矩阵就是“名义满秩、实际缺秩”直接求逆必翻车。我在现场处理数据时会先做一个快速判断对每个频点的 H 矩阵求条件数并把条件数随频率画出来。条件数超过 1000 的频段就要万分小心超过 10000 的频段基本不能用普通伪逆。cond_list zeros(size(freq)); for fi 1:length(freq) cond_list(fi) cond(squeeze(H(:, :, fi))); end figure; semilogy(freq, cond_list); xlabel(频率 (Hz)); ylabel(条件数 (log)); title(FRF矩阵条件数随频率变化);如果某频段条件数长得像悬崖一样那这段的 TPA 结果就别写了。要么增加测点、要么去掉冗余路径、要么用正则化硬压实在不行就放弃该频段的路径分解只保留总响应重构。4.2 正则化参数怎么选L曲线与广义交叉验证正则化参数 lambda 选多大是个永远躲不开的问题。lambda 太小病态没有改善lambda 太大载荷被过度“压缩”拟合误差变大。我先后用过经验值法、L 曲线法、广义交叉验证GCV法最终常驻的是 L 曲线。L 曲线的思想很直接把解范数||F_est||作为横轴、残差范数||H·F_est - y||作为纵轴取双对数坐标lambda 从很小到很大扫一遍会画出一条“L”形曲线。拐点处曲率最大对应的 lambda就是拟合误差与解稳定性之间的最佳折中点。lambdas logspace(-4, 1, 40); residual_norm zeros(size(lambdas)); solution_norm zeros(size(lambdas)); for li 1:length(lambdas) lam lambdas(li) * max(s); % s是之前SVD分解得到奇异值 s_reg s ./ (s.^2 lam^2); H_inv_tmp V * diag(s_reg) * U; F_tmp H_inv_tmp * y_meas; residual_norm(li) norm(H_k * F_tmp - y_meas, 2); solution_norm(li) norm(F_tmp, 2); end figure; loglog(residual_norm, solution_norm, o-); xlabel(残差范数 ||H F - y||); ylabel(解范数 ||F||); title(L曲线正则化参数选择);选拐点有点靠眼力但跑熟之后基本一眼就能挑出来。还拿不准的话可以用 GCV 做交叉验证Matlab 里自带的tikhonov或者写到 CVX 包里也能很好解决。我的经验是L 曲线选出的 lambda 通常落在最大奇异值的 0.5%8% 区间和 3.2 节的经验范围吻合。4.3 相位同步、平均次数与现场采集纪律TPA 的运算全是复数域运算相位哪怕偏 20 度矢量叠加的结果都可能完全变样这是很多人在现场最容易忽略的硬约束。相位一致的前提是采样同步。多通道采集卡每个通道必须严格同步采样曾经有人图省事用几个独立采集器接同一个转速信号分别采集结果每条通道独立触发相位差完全是随机漂移做出来的贡献占比毫无参考价值。你回头去问现场工程师“触发同步了吗”答案是“同步了啊都是转速触发”但不同采集器之间的启动延迟根本没有校准。所以做 TPA 测试必须确认采集设备是同一个时钟源统一同步采样各通道之间不能有毫秒级的延迟差。第二个容易栽跟头的是平均次数。现场振动信号噪声大若平均次数太少载荷识别结果会剧烈跳动。我做锤击 FRF 测试一般至少平均 35 次做运行响应测 TPA 时每个工况下至少采集 30 秒以上分段加窗平均不少于 50 次窗函数首选汉宁窗、重叠率 50%这样能压掉大部分随机噪声和脉冲干扰。另外测试期间必须保证转速稳定。齿轮系统对转速敏感转速一会高一会低传递特性实测不变但激励特性已经变了测点响应的频率结构也随之变化。低速阶段和高速阶段的数据混在一起平均出来的 TPA 结果谁也不代表。4.4 源数估计与实际匹配测点数量和等效源数量之间必须满足“测点数 ≥ 源数量”否则逆问题欠定Matlab 直接报错或给出荒谬结果。但“等效源数量”到底取多少本身是个变量。我的做法是把 FRF 矩阵做奇异值分解画出奇异值谱看有效奇异值的数目同时结合齿轮箱的结构列出所有可能的源输入轴承、输出轴承、齿轮啮合点、电机脚等作为候选源。然后跑一轮 TPA 预分析把贡献占比始终低于 1% 的候选源剔除掉再重新做一次识别。这种方式在工程上叫“逐步约简”过程不复杂但很有效能显著改善矩阵条件数。5. 从路径贡献到故障诊断结论一个二级齿轮箱示例的闭环5.1 一个聚焦示例输出端轴承内圈故障如何被TPA锁定用一个实际场景串一下整个流程。有一台二级齿轮箱输入转速 1470rpm输出轴转速约 245rpm。现场布了三个测点输入端轴承座、箱体中部、输出端轴承座。运行状态下采集到的振动信号在 950Hz 附近出现明显的调制边带边带间隔约 4.1Hz这个频率既接近输入轴转频的谐波也与输出端轴承故障特征频率比较接近一时间两派观点争执不下。我按照 TPA 流程处理首先锤击三个测点得到 3×3 的 FRF 矩阵然后采集运行响应对每个频点做正则化载荷识别算出三个等效源输入端轴承路径、齿轮啮合路径、输出端轴承路径对输出端测点的贡献。结果在 950Hz 处输出端轴承路径的贡献占比达到 68%齿轮啮合路径 22%输入端轴承路径 10%。这就说明输出端测点看到的主要边带能量来自输出端轴承路径故障源大概率在输出端轴承。配合对响应做包络谱分析在包络谱中确认了输出端轴承内圈故障特征频率BPFI及其谐波双重证据锁定了故障位置。拆机后果然发现输出端轴承内圈剥落。TPA 在这个案例里的价值不是替代包络谱而是在包络谱有歧义时给出第二维证据即使边带频率识别困难路径贡献的能量分布也能指向正确方向。5.2 从TPA数据反推传感器布点优化完成一次 TPA 分析后留下的 FRF 矩阵和路径贡献数据不应该只用一次就丢。它可以反向指导设备日常状态监测系统的传感器布置。沿用上面的案例如果 TPA 结果显示输出端轴承对输出端轴承座测点的贡献在 950Hz 附近最高那么长期监测系统在输出端轴承座上安装加速度计时分析频带就应该重点保留 9001000Hz而不是动不动只看全程频谱。又比如输入端轴承故障的特征频率集中在 2000Hz 以上输入端轴承座测点在该频段的路径贡献最大那么监测该轴承时采样率就不能低于 10kHz以免高频特征被“过滤”掉。很多厂搞设备监测系统传感器装得很密采样率也拉得很高但诊断率却不理想很大一部分原因就是布点和分析频带与路径贡献不匹配。拿 TPA 结果做一次“传感器点位体检”往往比盲目加装传感器更有效。5.3 边界条件变速工况下TPA还灵不灵TPA 的理论基石是线性时不变假设。齿轮系统在工作时啮合刚度按转频周期性变化严格说是时变系统。但工程上只要在稳态工况恒转速、恒负载下测取数据系统可以近似认为是时不变的TPA 的结果具备参考价值。这台二级齿轮箱的案例就建立在稳态工况下问题在于如果你遇上升速、降速、变载工况直接用运行响应做 TPA 就会出错——转速一变激励频率跟着变H 矩阵在这个频点对应的传递特性已经不是当前状态下的特性了。我的处理对策是变转速工况先做阶次跟踪把时域信号转换到角度域再按转速段分段做准稳态 TPA每个转速段单独计算路径贡献。这是一个可行的工程化扩展但前提是采集系统必须同步记录转速信号没有转速度做不了。还有一些极端情况比如齿轮出现断齿这种强冲击非线性故障结构响应已经不再是线性叠加的结果TPA 的假设前提就失效了。这时候硬套 TPA得到的贡献占比会很奇怪甚至完全颠倒了故障源顺序。所以我会先把时域包络和阶次特征看一遍如果冲击过强、明显非线性就直接改用冲击能量分析和时频域诊断没必要在 TPA 上死磕。回到齿轮系统诊断本身TPA 不是万能的但它是诊断逻辑里重要的一环。我自己在配合保质期内的新设备验收时也会用 TPA 来验证不同传感器位置的频响匹配度避免后期监测数据互相打架。这套从原理到 Matlab 实现再到现场坑点的流程跑下来最大的心得是TPA 算法本身不复杂真正把结果从“能算”变成“可信”的是你对传递路径物理含义的理解以及每一次测试中采集纪律的严格执行。用熟了之后它帮你解决的不只是“故障在哪”的问题更是“为什么测点在这里看不出来”的疑惑——这一点比单纯得到一个诊断结论有价值得多。
返回列表