简介这份资源面向从事风电结构强度校核与疲劳寿命评估的工程师及研究生聚焦风力发电机塔筒筒体在风载、自重与扭转载荷下的疲劳分析问题。包内以MATLAB脚本为主体共12个文件含11个m文件与1份readme说明压缩包约11KB涵盖雨流计数、螺栓校核、屈曲分析、疲劳寿命计算及数据采集等模块便于按流程调用与二次修改。雨流计数法将应力历程转化为等效疲劳载荷序列是金属结构疲劳寿命预测的关键手段配合有限元思路可评估塔筒在复杂工况下的安全裕度。目前已有900人学习下载适合希望快速搭建塔筒校核计算框架、理解雨流矩阵构建与疲劳寿命估算逻辑的读者参考也可作为课程设计或工程复核的脚本模板。1. 风力发电机塔筒校核里雨流计数法到底在算什么一台兆瓦级风机立在旷野里塔筒看着是个静止的铁疙瘩实际它每天都在被反复揉搓。风载方向随机、大小随机叶轮旋转还叠加一个周期分量塔筒底部焊缝处的应力时程是一条又乱又长的曲线。做疲劳校核时真正决定寿命的不是最大应力而是这条曲线里藏了多少个幅值—均值组合的循环。雨流计数法就是干这件事的把一段看似无规律的应力时程拆成一个个可累加的疲劳循环。这个标题里的活儿本质是一条完整的链路拿到塔筒危险点的应力时程 → 用雨流计数法统计循环 → 配合 S-N 曲线和 Miner 线性累积损伤算出总损伤 → 判断校核是否通过。MATLAB 在这里承担的是数据处理和计数算法实现.rar说明它大概率是一个打包好的脚本集合。适合谁看做风电结构、做疲劳寿命评估、或者手头有一堆应力时序数据不知道怎么变成损伤值的工程师。下面我按自己实际做过的顺序把这条链路拆开讲。2. 从应力时程到循环计数雨流法的实现骨架2.1 为什么塔筒疲劳必须用雨流而不是峰值计数很多人第一反应是数峰值不就行了把应力时程里的极大值极小值挑出来配对。这个做法在塔筒这种宽带随机载荷下会严重高估或低估损伤因为峰值计数丢掉了循环之间的嵌套关系。雨流法的核心思想来自一个直观的物理图像把应力时程旋转 90 度当成屋顶雨水从每个峰值内侧往下流遇到更大的峰值就停流过的路径构成一个闭合循环。塔筒载荷的典型特征是低频大风载叠加高频涡激和叶轮通过频率形成大量小循环嵌套在大循环里。雨流法能正确识别这种嵌套把每个闭合滞回环计为一个完整循环剩下的不闭合部分计为半循环。这是它和峰值计数、量程计数的本质区别也是疲劳规范里普遍推荐它的原因。实现上有两个主流路线一是基于四点法的在线计数适合流式数据二是先做峰谷提取再配对适合离线批处理。塔筒校核通常是离线做我一般用后者逻辑清晰、好调试。2.2 峰谷提取先把冗余点砍掉原始应力时程采样率可能上千赫兹直接送进雨流算法又慢又容易受数值噪声干扰。第一步是提取峰谷序列把单调段中间的过渡点全部删掉只保留局部极值。function [pv, idx] extractPeaksValleys(sig) % sig: 原始应力时程列向量 % pv: 峰谷序列 % idx: 峰谷在原序列中的位置 n length(sig); if n 3 pv sig; idx (1:n); return; end % 一阶差分判断单调性 d diff(sig); % 去掉差分为 0 的平段避免误判 keep [true; d ~ 0]; sig sig(keep); idx0 find(keep); d diff(sig); % 极值点前后差分符号相反 isExt [false; (d(1:end-1) .* d(2:end)) 0; false]; pv sig(isExt); idx idx0(isExt); % 保证首尾纳入便于后续闭合判断 if idx(1) ~ 1 pv [sig(1); pv]; idx [1; idx]; end if idx(end) ~ n pv [pv; sig(end)]; idx [idx; n]; end end逻辑说明先做差分把连续相等的平段剔除否则平段会产生虚假的符号变化。然后用相邻差分乘积小于零判断极值点。最后强制把首尾点纳入因为雨流配对需要序列首尾参与闭合判断。参数说明sig是应力列向量单位建议统一成 MPa。采样率不在这里处理但采样率太低会漏掉高频小循环塔筒分析里一般要求至少覆盖到 10 倍叶轮通过频率。如果数据里噪声大建议先做低通滤波再提取峰谷否则会数出一堆幅值极小但数量巨大的伪循环把损伤算飞。2.3 四点法配对雨流计数的核心循环峰谷序列拿到后用四点法反复扫描。规则是取连续四点 A、B、C、D如果 B 到 C 的幅值不超过 A 到 B 和 C 到 D就认为 B-C 构成一个完整循环把它计出来并从序列里删掉然后回退重新检查。function [ranges, means, counts] rainflowCounting(pv) % pv: 峰谷序列 % ranges: 每个循环的应力幅2*幅值 % means: 每个循环的均值 % counts: 每个循环的计数1 为全循环0.5 为半循环 s pv(:); ranges []; means []; counts []; i 1; while i length(s) - 3 A s(i); B s(i1); C s(i2); D s(i3); rAB abs(B - A); rBC abs(C - B); rCD abs(D - C); if rBC rAB rBC rCD % B-C 构成完整循环 ranges(end1) rBC; % 应力范围 means(end1) (B C) / 2; % 循环均值 counts(end1) 1.0; % 删除 B、C 两点 s(i1:i2) []; % 回退重新检查前面的点 i max(1, i - 1); else i i 1; end end % 剩余序列按半循环处理 for k 1:length(s)-1 ranges(end1) abs(s(k1) - s(k)); means(end1) (s(k) s(k1)) / 2; counts(end1) 0.5; end end逻辑说明while循环里每次只看连续四点满足条件就计一个全循环并删除中间两点然后回退一步。回退是关键因为删除后前面可能形成新的可配对四点。循环结束后剩下的点无法闭合按半循环各计 0.5。参数说明ranges是应力范围不是幅值幅值是范围的一半后续查 S-N 曲线时要注意用的是范围还是幅值这是最容易翻车的地方。means用于考虑平均应力修正塔筒焊缝对平均应力敏感不能忽略。counts区分全循环和半循环Miner 累积时直接乘进去。2.4 用 Goodman 修正把均值折进等效幅值雨流计出来的循环均值各不相同但 S-N 曲线通常是在对称循环均值 0下测的。塔筒焊缝存在残余应力和自重引起的静应力均值修正必须做。常用的是 Goodman 修正function Sa_eq goodmanCorrection(Sa, Sm, Su) % Sa: 应力幅值 % Sm: 平均应力 % Su: 材料抗拉强度 % Sa_eq: 等效对称应力幅值 Sa_eq Sa ./ (1 - Sm ./ Su); % 防止分母过小导致数值爆炸 Sa_eq(Sm Su) Inf; end逻辑说明Goodman 把拉伸平均应力对疲劳的削弱作用折算成等效幅值放大。Sm为正拉应力时等效幅值变大损伤加重Sm为负时等效幅值变小偏保守可截断。参数说明Su取塔筒钢材的抗拉强度常见 Q345 系列在 470~630 MPa 区间具体按材料证书取值。注意Sm/Su接近 1 时公式发散工程上一般限制Sm不超过0.8*Su超出部分按截断处理。这一步做完所有循环都变成均值 0 的等效幅值可以统一查 S-N 曲线。3. 损伤累积与校核把循环数变成寿命结论3.1 S-N 曲线取值塔筒焊缝为什么不能用母材曲线塔筒筒体的疲劳危险点几乎都在焊缝和热影响区母材的 S-N 曲线在这里完全不适用。焊缝的疲劳强度受构造细节控制规范里按细节类别detail category给曲线常见做法是取 90 MPa 或 80 MPa 级别的曲线具体取决于焊缝形式和检测等级。S-N 曲线在双对数坐标下是两段直线高周区斜率为 -1/m到某个循环数后转为水平截止线。塔筒分析里循环数通常在 10^6~10^8落在斜线段和截止段交界附近截止线取不取对结果影响很大。function N snCurve(Sa_eq, detail, Ncut) % Sa_eq: 等效应力幅值MPa % detail: 细节类别对应的 2e6 循环疲劳强度MPa % Ncut: 截止循环数低于此幅值不再计损伤 m 3; % 斜率倒数焊缝常用 3 Sref detail; % 2e6 循环对应的幅值 Nref 2e6; N zeros(size(Sa_eq)); for k 1:numel(Sa_eq) if Sa_eq(k) 0 N(k) Inf; else N(k) Nref * (Sref / Sa_eq(k))^m; end end % 超过截止循环数的按无限寿命处理 N(N Ncut) Inf; end逻辑说明按N Nref * (Sref/Sa)^m计算每个循环的许用次数。m3是焊缝常用值母材有时取 4 或 5用错斜率会让寿命差出好几倍。参数说明detail按规范查表常见焊接细节在 80~100 MPa。Ncut一般取 5e6 或 1e7取决于采用的规范。注意Sa_eq单位必须和detail一致都是 MPa混用 Pa 会让结果差 10^6 倍这是血泪级的低级错误。3.2 Miner 线性累积把每个循环的损伤加起来Miner 准则假设每个循环独立贡献损伤总损伤是各循环损伤之和损伤等于循环数除以该幅值下的许用次数。function D minerDamage(ranges, counts, detail, Ncut, Su) % ranges: 应力范围 % counts: 循环计数 Sa ranges / 2; % 范围转幅值 Sm zeros(size(Sa)); % 若已有均值替换此处 Sa_eq goodmanCorrection(Sa, Sm, Su); N snCurve(Sa_eq, detail, Ncut); d counts ./ N; d(isinf(N)) 0; % 无限寿命循环不贡献损伤 D sum(d); end逻辑说明先把范围折半成幅值做均值修正查 S-N 得到许用次数再按counts/N累加。N为无穷的循环损伤置零。参数说明D是总损伤D 1认为校核通过D 1不通过。工程上通常还留安全系数把D控制在 0.5 甚至更低。注意这里假设各循环顺序无关Miner 准则的局限就在这但塔筒这种随机载荷下它是被广泛接受的工程近似。3.3 一个最小可跑通的端到端脚本把上面的函数串起来用一段模拟应力时程验证整条链路。% 生成模拟应力时程低频大风载 高频叶轮通过分量 噪声 fs 50; % 采样率 Hz t 0:1/fs:600; % 10 分钟 low 30 * sin(2*pi*0.01*t); % 低频大风载 high 8 * sin(2*pi*1.2*t); % 叶轮通过频率 noise 1.5 * randn(size(t)); % 测量噪声 sig low high noise; % 链路 [pv, ~] extractPeaksValleys(sig); [ranges, means, counts] rainflowCounting(pv); D minerDamage(ranges, counts, 90, 5e6, 490); fprintf(循环总数: %d\n, length(ranges)); fprintf(总损伤 D %.4f\n, D); if D 1 fprintf(校核通过\n); else fprintf(校核不通过需优化\n); end逻辑说明模拟信号刻意叠加了低频和高频两个尺度用来检验雨流法能否正确分离嵌套循环。跑完后看循环总数是否合理再看损伤值。参数说明fs50对 1.2 Hz 分量足够实际塔筒数据采样率更高。detail90、Ncut5e6、Su490都是示例值实际按规范和材料替换。这段脚本的价值在于验证算法链路通不通不用于真实校核。4. 避坑与排查雨流计数最容易翻车的五个地方4.1 循环数出来是几千个损伤却大得离谱现象跑完发现D是几百甚至上千明显不合理。原因多半是噪声没滤峰谷提取把每个采样点的抖动都当成极值数出大量幅值极小但数量巨大的循环。虽然小循环单个损伤小但数量一多counts/N累加就爆了。解决进雨流之前先做低通滤波截止频率取到关心的最高载荷频率的 1.5~2 倍。或者设一个最小应力范围阈值低于阈值的循环直接丢弃。阈值一般取疲劳截止限对应的幅值。4.2 范围当幅值用寿命差一倍现象损伤值比预期大一倍或小一半反复检查算法没发现问题。原因雨流计出来的是应力范围峰谷差S-N 曲线用的是应力幅值范围的一半。这两个概念在代码里如果没统一Sa ranges和Sa ranges/2会让损伤差 2^m 倍m3 时就是 8 倍。解决在函数命名上强制区分变量名用range和amplitude不要都用S。查 S-N 前打印几个样本值人工核对量级。4.3 均值修正把压缩应力也放大了现象某些循环均值为负Goodman 修正后等效幅值反而变小损伤被低估。原因Goodman 公式对压缩平均应力给出小于 1 的分母修正理论上压缩应力对疲劳有利但塔筒焊缝存在残余拉应力实际不一定有利。解决工程保守做法是均值修正只对拉伸均值生效压缩均值不折减即Sm 0时按Sm 0处理。这样偏保守但安全。4.4 首尾点没闭合半循环数量对不上现象全循环数加半循环数不等于预期或者半循环特别多。原因峰谷序列首尾如果不是极值或者提取时没强制纳入首尾雨流配对会剩下一大堆无法闭合的点全按半循环计。解决提取峰谷时强制把序列首尾纳入如 2.2 节代码所示。另外检查数据长度太短的时程本身就无法形成完整统计。4.5 采样率不够高频循环直接丢失现象和别人的结果对比自己的循环数明显偏少损伤偏小。原因采样率低于载荷最高频率的 2 倍高频分量被混叠或直接漏掉。塔筒的叶轮通过频率和涡激频率如果没覆盖到小循环全丢。解决分析前先确认数据采样率必要时做频谱分析看能量分布。采样率至少取到关心最高频率的 5~10 倍。如果原始数据采样率不够任何后处理都补不回来。5. 把校核做成可复用的批处理流程单次跑通只是起点真实项目里塔筒有多个危险点、多个工况、多台机组靠手工一段段跑不现实。我一般会把整条链路封装成一个批处理函数输入是各测点各工况的应力时程文件输出是一张损伤汇总表。function summary batchFatigueCheck(fileList, cfg) % fileList: 应力时程文件路径列表 % cfg: 配置结构体含 detail、Ncut、Su、fs、fc 等 summary table(); for k 1:numel(fileList) raw load(fileList{k}); sig raw.sig(:); % 低通滤波截止频率 cfg.fc [b, a] butter(4, cfg.fc / (cfg.fs/2), low); sigF filtfilt(b, a, sig); [pv, ~] extractPeaksValleys(sigF); [ranges, ~, counts] rainflowCounting(pv); D minerDamage(ranges, counts, cfg.detail, cfg.Ncut, cfg.Su); summary [summary; table(string(fileList{k}), ... length(ranges), D, D cfg.Dlimit, ... VariableNames, {file,nCycles,damage,pass})]; end end逻辑说明filtfilt做零相位滤波避免滤波引入的相位偏移影响峰谷位置。每个文件独立算损伤最后汇总成表pass列直接给出是否通过。参数说明cfg.fc是低通截止频率取关心最高载荷频率的 1.5 倍左右。cfg.Dlimit是损伤限值考虑安全系数后一般取 0.5。butter的阶数 4 是常用折中阶数太高会振铃太低滚降不够。批处理跑起来后真正要盯的是损伤分布如果某个测点损伤远高于其他先别急着下结论回去看它的应力时程是不是有异常尖峰或者数据缺失段。我踩过的一次坑是某段数据采集时传感器松动时程里全是高频毛刺滤波后仍然残留导致该点损伤虚高一个量级。后来加了一步数据质量检查统计时程的峰度和过零率异常段直接标记出来人工确认。另一个习惯是每次改参数细节类别、截止循环数、均值修正策略都保留一份结果快照因为这几个参数对最终损伤的影响是乘性的改一个可能让结论从通过变不通过。把参数和结果一起存档后面评审或者复核时能说清楚每个数字是怎么来的。这套流程不复杂但能省掉大量返工。希望帮到你。本文还有配套的精品资源点击获取 
长河君 · 2026/10/20 20:35:20 · 阅读 1.2k 风力发电机塔筒疲劳校核：MATLAB雨流计数法实现与损伤累积
简介这份资源面向从事风电结构强度校核与疲劳寿命评估的工程师及研究生聚焦风力发电机塔筒筒体在风载、自重与扭转载荷下的疲劳分析问题。包内以MATLAB脚本为主体共12个文件含11个m文件与1份readme说明压缩包约11KB涵盖雨流计数、螺栓校核、屈曲分析、疲劳寿命计算及数据采集等模块便于按流程调用与二次修改。雨流计数法将应力历程转化为等效疲劳载荷序列是金属结构疲劳寿命预测的关键手段配合有限元思路可评估塔筒在复杂工况下的安全裕度。目前已有900人学习下载适合希望快速搭建塔筒校核计算框架、理解雨流矩阵构建与疲劳寿命估算逻辑的读者参考也可作为课程设计或工程复核的脚本模板。1. 风力发电机塔筒校核里雨流计数法到底在算什么一台兆瓦级风机立在旷野里塔筒看着是个静止的铁疙瘩实际它每天都在被反复揉搓。风载方向随机、大小随机叶轮旋转还叠加一个周期分量塔筒底部焊缝处的应力时程是一条又乱又长的曲线。做疲劳校核时真正决定寿命的不是最大应力而是这条曲线里藏了多少个幅值—均值组合的循环。雨流计数法就是干这件事的把一段看似无规律的应力时程拆成一个个可累加的疲劳循环。这个标题里的活儿本质是一条完整的链路拿到塔筒危险点的应力时程 → 用雨流计数法统计循环 → 配合 S-N 曲线和 Miner 线性累积损伤算出总损伤 → 判断校核是否通过。MATLAB 在这里承担的是数据处理和计数算法实现.rar说明它大概率是一个打包好的脚本集合。适合谁看做风电结构、做疲劳寿命评估、或者手头有一堆应力时序数据不知道怎么变成损伤值的工程师。下面我按自己实际做过的顺序把这条链路拆开讲。2. 从应力时程到循环计数雨流法的实现骨架2.1 为什么塔筒疲劳必须用雨流而不是峰值计数很多人第一反应是数峰值不就行了把应力时程里的极大值极小值挑出来配对。这个做法在塔筒这种宽带随机载荷下会严重高估或低估损伤因为峰值计数丢掉了循环之间的嵌套关系。雨流法的核心思想来自一个直观的物理图像把应力时程旋转 90 度当成屋顶雨水从每个峰值内侧往下流遇到更大的峰值就停流过的路径构成一个闭合循环。塔筒载荷的典型特征是低频大风载叠加高频涡激和叶轮通过频率形成大量小循环嵌套在大循环里。雨流法能正确识别这种嵌套把每个闭合滞回环计为一个完整循环剩下的不闭合部分计为半循环。这是它和峰值计数、量程计数的本质区别也是疲劳规范里普遍推荐它的原因。实现上有两个主流路线一是基于四点法的在线计数适合流式数据二是先做峰谷提取再配对适合离线批处理。塔筒校核通常是离线做我一般用后者逻辑清晰、好调试。2.2 峰谷提取先把冗余点砍掉原始应力时程采样率可能上千赫兹直接送进雨流算法又慢又容易受数值噪声干扰。第一步是提取峰谷序列把单调段中间的过渡点全部删掉只保留局部极值。function [pv, idx] extractPeaksValleys(sig) % sig: 原始应力时程列向量 % pv: 峰谷序列 % idx: 峰谷在原序列中的位置 n length(sig); if n 3 pv sig; idx (1:n); return; end % 一阶差分判断单调性 d diff(sig); % 去掉差分为 0 的平段避免误判 keep [true; d ~ 0]; sig sig(keep); idx0 find(keep); d diff(sig); % 极值点前后差分符号相反 isExt [false; (d(1:end-1) .* d(2:end)) 0; false]; pv sig(isExt); idx idx0(isExt); % 保证首尾纳入便于后续闭合判断 if idx(1) ~ 1 pv [sig(1); pv]; idx [1; idx]; end if idx(end) ~ n pv [pv; sig(end)]; idx [idx; n]; end end逻辑说明先做差分把连续相等的平段剔除否则平段会产生虚假的符号变化。然后用相邻差分乘积小于零判断极值点。最后强制把首尾点纳入因为雨流配对需要序列首尾参与闭合判断。参数说明sig是应力列向量单位建议统一成 MPa。采样率不在这里处理但采样率太低会漏掉高频小循环塔筒分析里一般要求至少覆盖到 10 倍叶轮通过频率。如果数据里噪声大建议先做低通滤波再提取峰谷否则会数出一堆幅值极小但数量巨大的伪循环把损伤算飞。2.3 四点法配对雨流计数的核心循环峰谷序列拿到后用四点法反复扫描。规则是取连续四点 A、B、C、D如果 B 到 C 的幅值不超过 A 到 B 和 C 到 D就认为 B-C 构成一个完整循环把它计出来并从序列里删掉然后回退重新检查。function [ranges, means, counts] rainflowCounting(pv) % pv: 峰谷序列 % ranges: 每个循环的应力幅2*幅值 % means: 每个循环的均值 % counts: 每个循环的计数1 为全循环0.5 为半循环 s pv(:); ranges []; means []; counts []; i 1; while i length(s) - 3 A s(i); B s(i1); C s(i2); D s(i3); rAB abs(B - A); rBC abs(C - B); rCD abs(D - C); if rBC rAB rBC rCD % B-C 构成完整循环 ranges(end1) rBC; % 应力范围 means(end1) (B C) / 2; % 循环均值 counts(end1) 1.0; % 删除 B、C 两点 s(i1:i2) []; % 回退重新检查前面的点 i max(1, i - 1); else i i 1; end end % 剩余序列按半循环处理 for k 1:length(s)-1 ranges(end1) abs(s(k1) - s(k)); means(end1) (s(k) s(k1)) / 2; counts(end1) 0.5; end end逻辑说明while循环里每次只看连续四点满足条件就计一个全循环并删除中间两点然后回退一步。回退是关键因为删除后前面可能形成新的可配对四点。循环结束后剩下的点无法闭合按半循环各计 0.5。参数说明ranges是应力范围不是幅值幅值是范围的一半后续查 S-N 曲线时要注意用的是范围还是幅值这是最容易翻车的地方。means用于考虑平均应力修正塔筒焊缝对平均应力敏感不能忽略。counts区分全循环和半循环Miner 累积时直接乘进去。2.4 用 Goodman 修正把均值折进等效幅值雨流计出来的循环均值各不相同但 S-N 曲线通常是在对称循环均值 0下测的。塔筒焊缝存在残余应力和自重引起的静应力均值修正必须做。常用的是 Goodman 修正function Sa_eq goodmanCorrection(Sa, Sm, Su) % Sa: 应力幅值 % Sm: 平均应力 % Su: 材料抗拉强度 % Sa_eq: 等效对称应力幅值 Sa_eq Sa ./ (1 - Sm ./ Su); % 防止分母过小导致数值爆炸 Sa_eq(Sm Su) Inf; end逻辑说明Goodman 把拉伸平均应力对疲劳的削弱作用折算成等效幅值放大。Sm为正拉应力时等效幅值变大损伤加重Sm为负时等效幅值变小偏保守可截断。参数说明Su取塔筒钢材的抗拉强度常见 Q345 系列在 470~630 MPa 区间具体按材料证书取值。注意Sm/Su接近 1 时公式发散工程上一般限制Sm不超过0.8*Su超出部分按截断处理。这一步做完所有循环都变成均值 0 的等效幅值可以统一查 S-N 曲线。3. 损伤累积与校核把循环数变成寿命结论3.1 S-N 曲线取值塔筒焊缝为什么不能用母材曲线塔筒筒体的疲劳危险点几乎都在焊缝和热影响区母材的 S-N 曲线在这里完全不适用。焊缝的疲劳强度受构造细节控制规范里按细节类别detail category给曲线常见做法是取 90 MPa 或 80 MPa 级别的曲线具体取决于焊缝形式和检测等级。S-N 曲线在双对数坐标下是两段直线高周区斜率为 -1/m到某个循环数后转为水平截止线。塔筒分析里循环数通常在 10^6~10^8落在斜线段和截止段交界附近截止线取不取对结果影响很大。function N snCurve(Sa_eq, detail, Ncut) % Sa_eq: 等效应力幅值MPa % detail: 细节类别对应的 2e6 循环疲劳强度MPa % Ncut: 截止循环数低于此幅值不再计损伤 m 3; % 斜率倒数焊缝常用 3 Sref detail; % 2e6 循环对应的幅值 Nref 2e6; N zeros(size(Sa_eq)); for k 1:numel(Sa_eq) if Sa_eq(k) 0 N(k) Inf; else N(k) Nref * (Sref / Sa_eq(k))^m; end end % 超过截止循环数的按无限寿命处理 N(N Ncut) Inf; end逻辑说明按N Nref * (Sref/Sa)^m计算每个循环的许用次数。m3是焊缝常用值母材有时取 4 或 5用错斜率会让寿命差出好几倍。参数说明detail按规范查表常见焊接细节在 80~100 MPa。Ncut一般取 5e6 或 1e7取决于采用的规范。注意Sa_eq单位必须和detail一致都是 MPa混用 Pa 会让结果差 10^6 倍这是血泪级的低级错误。3.2 Miner 线性累积把每个循环的损伤加起来Miner 准则假设每个循环独立贡献损伤总损伤是各循环损伤之和损伤等于循环数除以该幅值下的许用次数。function D minerDamage(ranges, counts, detail, Ncut, Su) % ranges: 应力范围 % counts: 循环计数 Sa ranges / 2; % 范围转幅值 Sm zeros(size(Sa)); % 若已有均值替换此处 Sa_eq goodmanCorrection(Sa, Sm, Su); N snCurve(Sa_eq, detail, Ncut); d counts ./ N; d(isinf(N)) 0; % 无限寿命循环不贡献损伤 D sum(d); end逻辑说明先把范围折半成幅值做均值修正查 S-N 得到许用次数再按counts/N累加。N为无穷的循环损伤置零。参数说明D是总损伤D 1认为校核通过D 1不通过。工程上通常还留安全系数把D控制在 0.5 甚至更低。注意这里假设各循环顺序无关Miner 准则的局限就在这但塔筒这种随机载荷下它是被广泛接受的工程近似。3.3 一个最小可跑通的端到端脚本把上面的函数串起来用一段模拟应力时程验证整条链路。% 生成模拟应力时程低频大风载 高频叶轮通过分量 噪声 fs 50; % 采样率 Hz t 0:1/fs:600; % 10 分钟 low 30 * sin(2*pi*0.01*t); % 低频大风载 high 8 * sin(2*pi*1.2*t); % 叶轮通过频率 noise 1.5 * randn(size(t)); % 测量噪声 sig low high noise; % 链路 [pv, ~] extractPeaksValleys(sig); [ranges, means, counts] rainflowCounting(pv); D minerDamage(ranges, counts, 90, 5e6, 490); fprintf(循环总数: %d\n, length(ranges)); fprintf(总损伤 D %.4f\n, D); if D 1 fprintf(校核通过\n); else fprintf(校核不通过需优化\n); end逻辑说明模拟信号刻意叠加了低频和高频两个尺度用来检验雨流法能否正确分离嵌套循环。跑完后看循环总数是否合理再看损伤值。参数说明fs50对 1.2 Hz 分量足够实际塔筒数据采样率更高。detail90、Ncut5e6、Su490都是示例值实际按规范和材料替换。这段脚本的价值在于验证算法链路通不通不用于真实校核。4. 避坑与排查雨流计数最容易翻车的五个地方4.1 循环数出来是几千个损伤却大得离谱现象跑完发现D是几百甚至上千明显不合理。原因多半是噪声没滤峰谷提取把每个采样点的抖动都当成极值数出大量幅值极小但数量巨大的循环。虽然小循环单个损伤小但数量一多counts/N累加就爆了。解决进雨流之前先做低通滤波截止频率取到关心的最高载荷频率的 1.5~2 倍。或者设一个最小应力范围阈值低于阈值的循环直接丢弃。阈值一般取疲劳截止限对应的幅值。4.2 范围当幅值用寿命差一倍现象损伤值比预期大一倍或小一半反复检查算法没发现问题。原因雨流计出来的是应力范围峰谷差S-N 曲线用的是应力幅值范围的一半。这两个概念在代码里如果没统一Sa ranges和Sa ranges/2会让损伤差 2^m 倍m3 时就是 8 倍。解决在函数命名上强制区分变量名用range和amplitude不要都用S。查 S-N 前打印几个样本值人工核对量级。4.3 均值修正把压缩应力也放大了现象某些循环均值为负Goodman 修正后等效幅值反而变小损伤被低估。原因Goodman 公式对压缩平均应力给出小于 1 的分母修正理论上压缩应力对疲劳有利但塔筒焊缝存在残余拉应力实际不一定有利。解决工程保守做法是均值修正只对拉伸均值生效压缩均值不折减即Sm 0时按Sm 0处理。这样偏保守但安全。4.4 首尾点没闭合半循环数量对不上现象全循环数加半循环数不等于预期或者半循环特别多。原因峰谷序列首尾如果不是极值或者提取时没强制纳入首尾雨流配对会剩下一大堆无法闭合的点全按半循环计。解决提取峰谷时强制把序列首尾纳入如 2.2 节代码所示。另外检查数据长度太短的时程本身就无法形成完整统计。4.5 采样率不够高频循环直接丢失现象和别人的结果对比自己的循环数明显偏少损伤偏小。原因采样率低于载荷最高频率的 2 倍高频分量被混叠或直接漏掉。塔筒的叶轮通过频率和涡激频率如果没覆盖到小循环全丢。解决分析前先确认数据采样率必要时做频谱分析看能量分布。采样率至少取到关心最高频率的 5~10 倍。如果原始数据采样率不够任何后处理都补不回来。5. 把校核做成可复用的批处理流程单次跑通只是起点真实项目里塔筒有多个危险点、多个工况、多台机组靠手工一段段跑不现实。我一般会把整条链路封装成一个批处理函数输入是各测点各工况的应力时程文件输出是一张损伤汇总表。function summary batchFatigueCheck(fileList, cfg) % fileList: 应力时程文件路径列表 % cfg: 配置结构体含 detail、Ncut、Su、fs、fc 等 summary table(); for k 1:numel(fileList) raw load(fileList{k}); sig raw.sig(:); % 低通滤波截止频率 cfg.fc [b, a] butter(4, cfg.fc / (cfg.fs/2), low); sigF filtfilt(b, a, sig); [pv, ~] extractPeaksValleys(sigF); [ranges, ~, counts] rainflowCounting(pv); D minerDamage(ranges, counts, cfg.detail, cfg.Ncut, cfg.Su); summary [summary; table(string(fileList{k}), ... length(ranges), D, D cfg.Dlimit, ... VariableNames, {file,nCycles,damage,pass})]; end end逻辑说明filtfilt做零相位滤波避免滤波引入的相位偏移影响峰谷位置。每个文件独立算损伤最后汇总成表pass列直接给出是否通过。参数说明cfg.fc是低通截止频率取关心最高载荷频率的 1.5 倍左右。cfg.Dlimit是损伤限值考虑安全系数后一般取 0.5。butter的阶数 4 是常用折中阶数太高会振铃太低滚降不够。批处理跑起来后真正要盯的是损伤分布如果某个测点损伤远高于其他先别急着下结论回去看它的应力时程是不是有异常尖峰或者数据缺失段。我踩过的一次坑是某段数据采集时传感器松动时程里全是高频毛刺滤波后仍然残留导致该点损伤虚高一个量级。后来加了一步数据质量检查统计时程的峰度和过零率异常段直接标记出来人工确认。另一个习惯是每次改参数细节类别、截止循环数、均值修正策略都保留一份结果快照因为这几个参数对最终损伤的影响是乘性的改一个可能让结论从通过变不通过。把参数和结果一起存档后面评审或者复核时能说清楚每个数字是怎么来的。这套流程不复杂但能省掉大量返工。希望帮到你。本文还有配套的精品资源点击获取