基于改进VMD的风电齿轮箱不平衡故障特征提取_周福成.pdf

 振动与冲击 第 39 卷第 5 期JOURNAL OF VIBRATION AND SHOCKVol． 39 No． 5 2020 基金项 目 国 家 自 然 科 学 基 金 51377074 ; 中 央 高 校 基 金 项 目 2017MS188 收稿日期2018 －08 －24修改稿收到日期2018 －11 －26 第一作者 周福成 男， 博士， 副教授， 1977 年生 基于改进 VMD 的风电齿轮箱不平衡故障特征提取 周福成1，唐贵基2，何玉灵2 1． 华北电力大学 科技学院， 河北 保定071003; 2． 华北电力大学 机械工程系， 河北 保定 071003 摘要变分模态分解 Variational Mode Decomposition，VMD 是一种不同于递归式模态分解新方法， 具有优良的 频率剖分特性， 但其在处理信号时受分量个数影响严重， 通过主观经验难以合理设置该参数。针对该问题， 利用奇异值分 解清晰的信噪分辨能力， 根据奇异值最佳有效秩阶次自动搜寻 VMD 的分量个数， 提出了一种改进变分模态分解的风电齿 轮箱不平衡故障特征提取方法。通过仿真信号及轴不平衡实验信号对该方法进行了验证， 并将其应用于风电齿轮箱稳定 工况下的现场故障诊断中， 均成功提取出微弱特征频率信息， 实现对齿轮箱不平衡故障的有效判别， 具有一定可靠性。 关键词改进变分模态分解; 风力发电机; 不平衡故障; 故障特征提取; 现场数据 中图分类号TH165． 3文献标志码ADOI 10． 13465/j． cnki． jvs． 2020． 05． 023 Unbalanced fault feature extraction for wind power gearbox based on improved VMD ZHOU Fucheng1，TANG Guiji2，HE Yuling2 1． School of Science and Technology，North China Electric Power University，Baoding 071003，China; 2． Department of Mechanical Engineering，North China Electric Power University，Baoding 071003，China Abstract Variational mode decomposition VMDis a new modal decomposition different from recursive one，and has a good frequency partition characteristics． However，it is seriously affected by number of components in signal processing，so it is difficult to set up its parameters rationally with subjective experience． Here，to solve this problem，the singular value decomposition SVDwith clear signal- to- noise resolution ability was used to automatically search component number of VMD according to the optimal effective rank order of singular value， and propose an improved VMD for unbalanced fault feature extraction of a wind power gearbox． Simulated signals and shaft unbalance test ones were used to verify the proposed ． The proposed was applied in field fault diagnosis of a wind power gearbox under stable working condition to successfully extract weak fault feature frequency ination，and realize effective judgement for wind power gearbox’ s unbalanced fault with certain reliability． Key wordsimproved variational modal decomposition;wind turbine;unbalanced fault;fault feature extraction; field data 经验模态分解 Empirical Mode Decomposition， EMD ［1 ］是一种自适应信号处理方法， 由于它可将复杂 非平稳信号分解为几个本征模态函数 Intrinsic Mode Function， IMF 分量， 故在信号处理中应用较广， 但 EMD 方法缺乏严格的数学理论基础导致迭代计算量 大， 在包络拟合时易出现端点效应及模态混叠等诸多 问 题 ［2 ］。Dragomiretskiy 等［3 ］ 提 出 变 分 模 态 分 解 Variational Modal Decomposition，VMD 方法， VMD 不 同于 EMD 这种递归式分解算法， 其中心思想是约束性 变分问题， 具有坚实的理论基础 ［4 ］。VMD 算法的前提 是设定每个分解模态分量都围绕在某个中心频率周 围， 将模态带宽的确定过程转化成了约束变分问题， 通 过对其求解实现了各个模态分量的分离。在合理调整 VMD 算法参数基础上， 能够实现目标信号频域的灵活 分割， 进而有效挖掘信号中掩盖的特征信息。目前部 分学者将该方法应用到机械故障诊断领域， 取得非常 不错的效果。文献［ 5］ 利用参数优化变分模态分解方 法对轴承故障进行分析， 成功提取出微弱特征频率信 息， 对滚动轴承早期故障进行有效判别。文献［ 6］ 提出 了基于变分模态分解和奇异值分解的特征提取方法， 采用标准模糊 C 均值聚类进行故障识别。文献［ 7］ 针 对滚动轴承早期故障振动信号信噪比低、 故障特征提 取困难的问题， 提出了基于变分模态分解和能量算子 的滚动轴承故障特征提取方法。文献［ 8］提出了依据 ChaoXing 功率谱密度极值点自适应确定模态数量与中心频率的 参数自适应变分模态分解方法， 对行星齿轮箱第 2 级 太阳轮裂纹的故障进行了诊断。文献［ 9］提出一种基 于变分模态分解的多特征量风电机组轴承故障诊断方 法， 优化核参数和惩罚参数基础上能对轴承故障进行 最佳预测。 VMD 处理信号前需要设定分量个数 K， K 值的准 确判别对信号分解后的模态影响至关重要， K 值偏小 导致某个分量模态无法识别， K 值偏大导致模态混叠。 文献［ 10- 11］ 依靠经验进行人工设定 K 值， 因缺乏判据 而不够严谨; 文献［ 12- 14］根据 EMD 分解后的自适应 分量进行确定 K 值， 但 EMD 在进行包络拟合时出现模 态混叠导致分量个数也不精准; 唐贵基等通过智能算 法进行多次迭代优化选取参数， 但计算时间长导致效 率低， 而当中心频率误差较大时， 一些算法无法正确实 现信号分解。本文利用奇异值最佳有效秩阶次自动适 应 VMD 的分量个数， 提出了一种改进变分模态分解的 风电齿轮箱不平衡故障特征提取方法， 通过仿真信号 与实测及现场数据分析表明该方法在故障特征提取上 的有效性和准确性， 并通过对比 EMD 及小波方法， 验 证了该方法的优越性。 1变分模态分解基本原理 在 VMD 算法中， 将“本征模态函数 Intrinsic Mode Function， IMF ” 定义为一个调幅调频信号， 其表达式为 uk t Ak t cos φk t 1 式中 φk t 为信号的相位; Ak t 为信号瞬时幅值; ωk t 为信号瞬时频率， ωk t φk t dφ k t dt ， 同相 位相比， 瞬时幅值和瞬时频率均都是缓变量。 VMD 约束变分模型实质为搜寻 K 个具有特定稀 疏性的 IMF 分量 uk t ， 使各分量的估计带宽之和达到 最小， 这里限定约束条件满足各分量之和等于原始信 号 f t ， 其模型的构造步骤如下 步骤 1首先对原信号进行 Hilbert 变换， 得到每 个 IMF 分量 uk t 的解析信号， 获取其单边频谱为 δ t j π t * uk t 2 步骤2其次对每个 IMF 分量的解析信号， 逐一估 计其相应的中心频率 ωk， 将它与指数信号 e － jωkt 相乘， 将频谱转移到基频带 δ t j π t * uk t [] e －jωkt 3 步骤 3最后计算该调制信号梯度的 L2 范数的平 方， 估算出各个 IMF 分量的带宽， 结果构造如下形式的 约束变分模型 min { uk} ， { ωk} ∑ K k 1 tδ t j π t * uk t [] e －jωkt {} 2 2 s． t．∑ K k 1 uk t f t { 4 式中 δ t 为单位脉冲函数; j 为虚数单位; * 为卷积运 算; t为对函数求偏导; { uk} 为分解后的 k 个 IMF 分 量， { uk} { u1， u2， ， uK } ; { ωk} 为各分量中心频率， { ω k} { ω1 ， ω 2， ， ωK} 。 为了求解上述变分问题， 在式 4 中引入惩罚因子 α 与拉格朗日乘子 λ， 这样就将约束变分问题化为非约 束变分问题， 得到增广拉格朗日表达式 L { uk} ， { ωk} ， λ α∑ K k 1 tδ t j π t * uk t [] e －jωkt 2 2 f t－∑ K k 1 uk t 2 2 〈λ t ， f t－∑ K k 1 uk t 〉 5 利用交替方向乘子算法 Alternate Direction of Multipliers， ADMM 不断迭代更新 un 1 k 、 ω n 1 k 和 λn 1 来寻找式 5 的 “鞍点” ， 最终完成式 4 中约束变分问 题求解。在搜索“鞍点” 的过程中， 三个变量更新表达 式如下 u n1 k ω f ω－ ∑ k－1 i 1 u n1 i ω－∑ K i k1u n i ω λ n ω 2 1 2α ω － ωnk 2 6 ωn1 k ∫ ∞ 0 ω u n1 k ω2 dω ∫ ∞ 0 u n1 k ω2 dω 7 λ n1 ω λn ω τf ω－ ∑ K k 1 u n1 k ω 8 式中 为傅里叶变换; n 为迭代次数; τ 为保真系数。 当迭代求解变分模型时， 每个 IMF 分量的频率中 心和 带 宽 不 断 更 新，直 到 满 足 迭 代 停 止 条 件 ∑ K k 1 un1 k ω－ u n k ω 2 2 u n k ω 2 2＜ ε 后， 再结束 整个循环， 最后根据实际信号的频域特性进行信号频 带的自适应分割， 将所得的{ uk ω } 通过傅里叶逆变 换转变为时域的 IMF 分量。 2改进变分模态分解方法 根据奇异值分解原理 其具体算法由于篇幅所限， 不再赘述 ， 其去噪过程是保留前 K 个较大奇异值， 将 K 以后的奇异值置零， 重构以后信号基本能较准确反 映原信号。从其构造原理可知， 前 K 个较大的奇异值 反映了信号的主要成分， K 以后较小的奇异值反映了 噪声成分， K 值的搜索可通过奇异值分布曲线来实现， 171第 5 期周福成等基于改进 VMD 的风电齿轮箱不平衡故障特征提取 ChaoXing 当奇异值由极大值快速下降到稳定的极小值时， 其转 折点就是奇异值阶数 K 值。而 VMD 算法在处理复合 信号时， 也是将主要信号成分分成 K 个， 剩下 1 个残次 分量是干扰噪声， 并被滤波处理。因此奇异值分解后 的突变点 K 值和 VMD 的分量个数 K 在信号处理过程 中起的作用是一致的， 所以 VMD 的分量个数 K 可以根 据奇异值分解的有效阶次来定。 目前较多学者依靠经验来寻找降噪阶次， 由于实 际信号奇异值变化成曲线趋势， 奇异值突变点不易观 察， 直接寻有效秩阶次不够严谨。如果根据估算在有 效秩附近取值， 有可能就会导致噪声过多或者损失有 效信号， 从而影响故障特征的提取。针对这种情况， 为 顺利搜寻奇异值突变点阶次， 本文提出一种自适应的 方法， 具体实施如下 1根据奇异值分布曲线， 计算各阶次对应奇异 值的斜率 gm gm ds dm 9 式中 s 为奇异值; m 为奇异值对应阶次。 2计算相邻两个奇异值斜率的差值 dgm dgm gm1－ gm 10 3如果相邻两个奇异值斜率的差值满足下列 条件 if dgm＞ thr1＆dgm1＜ thr2 ， then k m 11 式中 thr1为斜率差值的上阈值; thr2 为斜率差值的下 阈值， 其具体数值得根据信号强度大小而定， 通常取 thr150， thr25 便能求出奇异值突变点阶次 k。改进 变分模态分解方法其具体流程如图 1 所示。 图 1改进变分模态分解方法 Fig． 1 of improved variational mode decomposition 3仿真信号分析 为验证上述问题， 设置一个仿真信号进行说明。 仿真信号 x t 由间隔为 0． 01 s 的冲击信号 x1 t 5e －300πtsin 6 000πt 、 正弦信号 x 2 t sin 100πt 及 随机噪声 x3 t 组成， 设置采样频率为 20 000 Hz， 采样 点数 4 096 点。该合成信号及分量波形如图 2 所示， 其 表达式如下 x t x1 t x2 t x3 t 12 图 2合成信号及分量波形 Fig． 2Wave of composite signal and its components 对合成信号进行奇异值分解后得到图 3 所示奇异 值分布曲线， 根据式 9～ 11 计算得到奇异值突变点 为 3， 这跟仿真信号由三个分量组成完全一致， 因此确 定 VMD 分量个数 K 3。 图 3合成信号奇异值分布曲线 Fig． 3Singular value distribution curve of the synthetic signal 对该复合信号进行 VMD 处理， 设置惩罚因子 α 2 000 该参数影响各模态带宽和信号幅值， 根据原文 建议取 2 000 具有普遍适应性 ， 分量个数 K 3， 得到 如图 4 所示波形。从图上看出， VMD 将合成信号分解 成 3 个分量后， 按照随机噪声、 冲击信号、 正弦信号顺 序排列， 分解后的各分量波形相对原信号几乎没有变 化， 说明 3 个分量被完整剥离， 且不存在互相混叠或相 互干扰现象。 图 4合成信号 VMD 分解结果 Fig． 4Decomposition result of synthesis signal by VMD 271振 动 与 冲 击2020 年第 39 卷 ChaoXing 4实验信号分析 转子实验台及传感器的布置如图 5 所示， 为使转 盘的重心偏离中心位移， 在转盘 45角旋入一颗 1 g 的 螺钉， 模拟轴不平衡故障， 圆盘两侧各布置两个涡流传 感器， 成90夹角， 对转子系统 X、 Y 向摆度进行测量， 在 靠近端点轴的上部布置一个键相传感器， 用来测量实 时转速。实验时， 电机转速为 1 550 r/min， 经计算转轴 频率 fl 25． 8 Hz， 采样频率为 1 280 Hz， 采样点数为 4 096。 图 5转子实验台及传感器布置 Fig． 5Rotor test bed and sensor layout 先分析电机运转时采集的转轴信号， 其波形及频 谱如图 6 所示， 从时域波形上看， 信号由于受噪声干扰 严重， 振幅沿正弦曲线“抖动” ， 且幅值较大， 频谱图中 从低频到高频分布一系列幅值不等的信号， 在低频段 电机转速基频 fl很突出， 二倍频 2fl也较明显， 其他高 次谐波被强烈噪声淹没。因此该轴不平衡故障的特征 并不突出， 据此判断效果欠佳。 a实测信号波形 b实测信号频谱 图 6转轴信号及频谱 Fig． 6Signal and spectrum of rotating shaft 为提取该信号故障特征， 并作出准确判别， 下面采 用本文所提方法进行验证。将实测信号先进行奇异值 分解， 得到其奇异值分布曲线如图 7 所示。 根据前文算法寻找奇异值突变点为 3， 以此为 VMD 分量个数， 自适应处理后得到图 8 所示分解结果。 图 7转轴信号奇异值分布曲线 Fig． 7Singular value distribution curve of shaft a分量波形 b分量频谱 图 8低速轴信号 VMD 处理结果 Fig． 8Decomposition result of low- speed shaft signal by VMD 分解出来的 3 个分量中， 第 1 个分量 IMF1 信号振 幅小、 波动快、 频率高， 满足随机噪声的特点， 第 3 个分 量 IMF3 由一部分波动慢、 振幅更小的低频噪声和转轴 的高次谐波组成， 第 2 个分量 IMF2 的振动时域波形近 似为正弦波， 频谱分布中， 谐波能量集中在基频， 振幅 突出， 并且出现了较小振幅的高次谐波， 整个频谱呈现 了所谓的 “枞树形” ， 其它干扰少， 转轴不平衡故障的主 要振动特征已经呈现， 证明故障特征提取成功。 为对比说明 VMD 方法的优势， 利用 EMD 方法对 上述转轴实测故障信号分别进行分析， 得到如图 9 所 示波形和频谱。信号经 EMD 处理后共获得 10 个 C 分 量， 且前 4 个 C 分量能量比较大， 后 6 个 C 分量能量较 小， 并且与原始信号相关度不大， 属于缺乏实际物理意 义的虚假分量， 因此分析前 4 个能量较大、 相关性较强 分量 C1 ～ C4 的波形及频谱即可。从图 9 的分解信号 上可观察到， C4 分量波形接近正弦波， 但是中间段被 噪声干扰打断， 频谱低频信号较多， 分解不彻底， 但可 见转频及其二三次谐波， 前 4 个 C 分量频谱边频带重 371第 5 期周福成等基于改进 VMD 的风电齿轮箱不平衡故障特征提取 ChaoXing 叠， 没有很好区分所在频谱段信号的中心频率， 造成模 态混叠现象， 相比 VMD 处理结果， 时域波形和频谱的 主要振动特征体现不够充分。 a分量波形 b分量频谱 图 9转轴信号 EMD 处理结果 Fig． 9Decomposition result of shaft signal by EMD 5风电现场故障信号分析 本文现场数据采集于某风电场 S50/750 kW 机组 齿轮箱， 该齿轮箱采取一级行星两级平行结构， 采集信 号时机组处于额定工况附近， 传感器安装在齿轮箱的 高速、 中速、 低速轴承处及发电机轴承处， 测试位置分 垂直、 水平和轴向三个方向。齿轮箱结构及测点位置 如图 10 所示。 图 10齿轮箱结构及测点位置 Fig． 10Structure and measurement point location of gear box 采集数据时， 齿轮箱振动异常， 发电机转速为1 520 r/min， 经计算转频约等于 25． 3 Hz， 采样频率为 16 384 Hz。对各测点数据一一进行分析与比对， 发现测点 4 振动数据存在异常， 其时域波形和频谱如图 11 所示。 a信号波形 b信号频谱 图 11故障信号及频谱 Fig． 11Signal and spectrum of fault 图中信号波形除了振幅增大不存在其他明显规 律， 频谱受噪声干扰严重， 复杂的传递路径和恶劣的周 边环境带来的振动淹没了故障特征信号， 显然从频谱 中无从判断故障。用本文所提改进 VMD 方法对该故 障信号进行处理。奇异值分解后得到图 12 所示的奇 异值分布曲线， 依据奇异值突变点所求算法， 计算得出 奇异值阶数为 7， 据此确定 VMD 分量个数 K 7 后， 对 原信号进行 VMD 处理， 得到图 13 所示的波形和频谱。 图 12奇异值分布曲线 Fig． 12Singular value distribution curve 471振 动 与 冲 击2020 年第 39 卷 ChaoXing a分量波形 b分量频谱 图 13测点 4 信号 VMD 处理结果 Fig． 13Decomposition result of high- speed shaft signal by VMD 从分解的 7 个 IMF 分量来看， IMF7 分量波形具有 明显的正弦波特征， 频谱中可见 24 Hz 的一条谱线， 对 IMF7 分量进行细化谱分析， 得到图 14 所示频谱图， 图 中频谱出现了“枞树形” 分布特征， 且 24 Hz 的基频突 出， 跟发电机转频 25． 3 Hz 接近， 故可判断测点 4 处发 电机转轴出现不平衡故障。后经维修人员检修， 认证 了这一判断， 经过不平衡故障修复， 风力发电机运转 正常。 图 14IMF7 分量细化谱 Fig． 14Thinning spectrum of IMF7 component 为对比说明， 对该信号进行小波分析， 设小波分解 层数为 7， 得到如图 15 所示的波形和频谱。从分解波 形上看， S7 分量出现正弦波特征， 但是频谱上前 6 个分 量出现了明显的模态混叠现象， 说明该复合信号并没 有被完整分离， 为进一步观察 S7 分量的频谱构成， 对 其细化分析， 得到如图 16 所示频谱图。图中频谱虽然 也出现了 “枞树形” 排列特征， 但是基频幅值明显减少， 而且 29 Hz 的基频偏离发电机转频较远， 在其前 3 Hz 处出现干扰频率， 据此不能判断故障。反观本文所提 改进 VMD 方法， 处理信号后各分量模态清晰， 故障信 号特征明显， 表明了该方法在处理该类信号时的可 靠性。 a分量波形 b分量频谱 图 15测点 4 信号小波处理结果 Fig． 15Decomposition result of high- speed shaft signal by wavelet 图 16S7 分量细化谱 Fig． 16Thinning spectrum of S7 component 571第 5 期周福成等基于改进 VMD 的风电齿轮箱不平衡故障特征提取 ChaoXing 6结论 VMD 的频域剖分特性与 EMD 完全不同， 具有类似 但不同于小波包变换的带通滤波性质， 能够对信号区 域进行精细分析。本文利用奇异值分解有效降噪阶次 来确定 VMD 的分量个数， 通过仿真信号、 轴不平衡实 测信号及风电现场数据， 证明了该方法的正确性， 通过 和 EMD 方法、 小波方法对比分析， 验证了该方法的优 越性。 参 考 文 献 ［1］ HUANG N E，SHEN Z，LONG S R，et al． The empirical mode decomposition and the Hilbert spectrum for nonlinear and non- stationary time series analysis［ J］ ． Proceeding of the Royal Society A， 1998， 454 1971 903- 995． ［2］ LIU H H，HAN M H． A fault diagnosis based on local mean decomposition and multi- scale entropy for roller bearings ［J］ ． Mechanism and Machine Theory，2014，75 67- 78． ［3］ DRAGOMIRETSKIYK， ZOSSOD．Variationalmode decomposition［ J］ ． IEEE Transactions on Signal Processing， 2014， 62 3 531- 544． ［4］ 王晓龙． 基于振动信号处理的滚动轴承故障诊断方法研 究［ D］ ． 北京 华北电力大学 北京 ， 2017． ［5］ 唐贵基， 王晓龙． 参数优化变分模态分解方法在滚动轴承 早期故障诊断中的应用［ J］ ． 西安交通大学学报， 2015， 49 5 73- 81． TANGGuiji， WANGXiaolong．Parameteroptimized variational mode decomposition with application to incipient fault diagnosis of rolling bearing［J］ ．Journal of Xi’ an Jiaotong University， 2015， 49 5 73- 81． ［6］ 刘长良， 武英杰， 甄成刚． 基于变分模态分解和模糊 C 均值 聚类的滚动轴承故障诊断［ J］ ． 中国电机工程学报， 2015， 35 13 3358- 3365． LIU Changliang，WU Yingjie，ZHEN Chenggang．Rolling bearingfaultdiagnosisbasedonvariationalmode decomposition and fuzzy c means clustering［J］ ． Proceedings of the CSEE， 2015， 35 13 3358- 3365． ［7］ 马增强， 李亚超， 刘政， 等． 基于变分模态分解和 Teager 能 量算子的滚动轴承故障特征提取［ J］ ． 振动与冲击， 2016， 35 13 134- 139． MA Zengqiang，LIYachao，LIUZheng，etal．Rolling bearings’fault feature extraction based on variational mode decomposition and Teager energy operator［J］ ．Journal of Vibration and Shock， 2016， 35 13 134- 139． ［8］ 孙灿飞， 王友仁， 沈勇， 等． 基于参数自适应变分模态分解 的行星齿轮箱故障诊断［J］ ． 航空动力学报， 2018 11 2756- 2765． SUN Canfei，WANG Youren，SHEN Yong，et al．Fault diagnosis of planetary gearbox based on adaptive parameter variational mode decomposition［J］ ．Journal of Aerospace Power， 2018 11 2756- 2765． ［9］ 张瑶， 张宏立． 基于 VMD 多特征量风电机组轴承故障诊断 法［ J］ ． 计算机仿真， 2018， 35 9 98- 102． ZHANG Yao， ZHANG Hongli． Bearing fault diagnosis for wind turbine based on VMD［J］ ． Journal of Computer Simulation， 2018， 35 9 98- 102． ［ 10］ ABDOOS A A． Detection of current transer saturation based on variational mode decomposition analysis ［J］ ． Transmission ＆ Distribution， 2016， 10 11 2658- 2669． ［ 11］ SIVAVARAPRASAD G，SREE PADMAJA R，VENKATA RATNAM D． Mitigation of ionospheric scintillation effects on gnss signals using variational mode decomposition ［J］ ． Geoscience and Remote Sensing Letters，2017，14 3 389- 393． ［ 12］ LIU W，CAO S，WANG Z，et al． Spectral decomposition for hydrocarbon detection based on VMD and teager- kaiser energy ［ J］ ． Geoscience and Remote Sensing Letters， 2017， 14 4 539- 543． ［ 13］ LAHMIRI S．Comparing variational and empirical mode decomposition in forecasting day- ahead energy prices ［J］ ． Systems Journal， 2017， 11 3 1907- 1910． ［ 14］ YANG W， PENG Z， WEI K， et al．Superiorities of variationalmodedecompositionoverempiricalmode decomposition particularly in time- frequency feature extraction and wind turbine condition monitoring［ J］ ． Renewable Power Generation， 2017， 11 4 檿檿檿檿檿檿檿檿檿檿檿檿檿檿檿檿檿檿檿檿檿檿檿檿檿檿檿檿檿檿檿檿檿檿檿檿檿檿檿檿檿檿檿檿檿檿檿檿檿 443- 452． 上接第 169 页 ［ 15］ SHI Y， LI Zhongxian， HONG H．Mesh size effect in numerical simulation of blast wave propagation and interaction with structures［ J］ ． Transactions of Tianjin University， 2008， 14 6 396- 402． ［ 16］ LEI S，DU X L，XIN F． A study on the Mesh generation for numerical simulation of blast wave［ J］ ． Journal of Beijing University of Technology， 2010， 36 11 1465- 1470． ［ 17］ KRAUTHAMMER T． Modern protective structures civil and environmental engineering ［M］ ． Boca RatonCRC Press， 2008． 671振 动 与 冲 击2020 年第 39 卷 ChaoXing