水下自激吸气式脉冲射流装置瞬时冲击力分解

于晓龙， 刘新阳， 徐洪增， 陈 豪， 汪顺生， 高传昌

(1.山东黄河河务局 工程建设中心，济南 250011； 2. 华北水利水电大学 水利学院，郑州 450045)

自激脉冲射流按射流频率可分为高频脉冲和低频脉冲两种，高频脉冲的射流频率一般在几百赫兹甚至上千赫兹[1-3]，而低频脉冲的射流频率一般在10 Hz以内甚至不到1 Hz[4-6]，而且两种脉冲射流产生的机理也不相同[7-8]。前期研究表明自激脉冲射流瞬时冲击力时均值、脉冲峰值和脉冲幅值明显优于连续射流[9-10]，同时为了改善自激脉冲射流装置瞬时冲击力脉冲效果，国内外学者提出了很多方法，如引入介质[11-12]、结构尺寸优化[13-14]和形式改变[15]等。论文以具有低频脉冲特点的自激脉冲射流装置为研究对象，装置在非淹没条件和合适的结构参数下瞬时冲击力具有较好的脉冲效果，但在淹没条件下，由于围压的改变，其脉冲效果变差，因此笔者及团队提出了一种低压大流量水下自激吸气式脉冲射流装置，通过在装置边壁上开设吸气孔，气体通过装置内外压力差进入装置，一方面使内部流动结构由水相流动变为水气两相流动；另一方面达到了改善瞬时冲击力脉冲效果的目的[16-17]。水下自激吸气式脉冲射流装置产生的瞬时冲击力具有明显的非平稳性，根据冲击射流特征，瞬时冲击力可分解为时变冲击力和脉动冲击力[18]，其中时变冲击力是由装置冲击靶盘后动量减小并转化为作用力引起的，所含频率成分较低，可以反映瞬时冲击力的脉冲幅值、脉冲峰值和脉冲周期，而脉动冲击力起源于射流中时均速度场畸变和大尺度紊动涡体引起的脉动速度场[19]，所含频率成分较高，可以反映瞬时冲击力的脉动强度，因此，如何将瞬时冲击力分解为时变冲击力和脉动冲击力，对于分析瞬时冲击力脉冲效果和时频特性至关重要。

目前关于非平稳信号时变部分的提取方法包括短时间平均法、拟合法、滤波法和模型法等[20]，这些方法有的存在时间段选择问题，有的需要假定时间序列趋势项类型，不适宜处理具有复杂变化趋势的非平稳信号，因此论文拟采用两种新的提取方法：小波变换和经验模态分解(Empirical Mode Decomposition,EMD)[21-22]，其中小波变换存在小波基函数和分解层次选择的问题，EMD可根据信号自身时间尺度特征自动分解，但存在数学理论依据不完善、端点处理和固有模态函数个数确定的问题。本文首先利用离散正交小波变换对水下自激吸气式脉冲射流装置瞬时冲击力进行时变冲击力和脉动冲击力的分解并进行消噪处理，确定小波基函数和分解层次，再利用EMD对消噪后的瞬时冲击力进行分解，确定计算时变冲击力所需的固有模态函数(Intrinsic Mode Functions,IMF)个数，并与离散正交小波变换得到的时变冲击力进行对比分析，最后得到适用于水下自激吸气式脉冲射流装置瞬时冲击力的分解方法，在此基础上，研究围压和吸气对时变冲击力脉冲幅值和瞬时冲击力时均值以及脉动冲击力脉动强度的影响规律。

1 水下自激吸气式脉冲射流冲击试验

水下自激吸气式脉冲射流装置如图1所示。由上喷嘴、腔体、碰撞体、下喷嘴和吸气孔组成。射流试验系统如图2所示。其工作过程简述如下：采用注水泵给压力容器罐注水，开启离心泵为水流提供工作压力，水流经过水下自激吸气式脉冲射流装置后吸气并形成脉冲射流，冲击靶盘产生冲击力。水流工作压力由泵系统提供，如图3所示。由多级离心泵和单级离心泵组成；压力容器罐可用于模拟水下环境，形成不同的围压，如图4所示。设计承受压力为1 MPa，射流装置安装在压力容器罐内部；电磁流量计用来测量流量；闸阀和压力表用于工作压力和围压的调节；气体涡轮流量计用于测量装置吸气量；靶盘上布置测压孔用于瞬时冲击力测量，在测压孔上布置气动接头且通过测压管与压力传感器连接，测压孔内径4 mm，深度2 mm，压力传感器测压范围在0～4 MPa。数据采集系统包括数据采集仪和计算机组成，如图5所示。数据采集仪用来接收压力传感器的信号，计算机用于对信号进行转换、分析和存贮，获得装置瞬时冲击力。

图1 水下自激吸气式脉冲射流装置Fig.1 Underwater self-excitation inspiration pulsed jet equipment

图2 射流试验系统Fig.2 Jet test system

图3 泵系统Fig.3 Pump system

图4 压力容器罐Fig.4 Pressure vessel tank

图5 数据采集系统Fig.5 Data acquisition system

关于装置结构参数和数据采集系统，即上喷嘴直径、下喷嘴直径、腔径、腔长、最优靶距、靶盘压力测点布置和压力传感器精度在笔者发表的论文中有详细的描述，这里不再赘述。在前期研究基础上，选择瞬时冲击力脉冲效果较好的装置结构参数，即上喷嘴直径10 mm、下喷嘴直径16 mm、腔径85 mm和腔长55 mm，靶距取100 mm。装置工作压力取2.2 MPa，围压取0.1 MPa，0.2 MPa，0.3 MPa，0.4 MPa，0.5 MPa和0.6 MPa。具有低频脉冲的低压大流量自激脉冲射流装置工作压力范围一般在0～4 MPa，而水下自激吸气式脉冲射流装置在不同围压下吸气的起始工作压力也不相同，围压越高，起始工作压力越大[23]，在试验过程中发现工作压力2.2 MPa可以保证装置在不同围压下均有较好的吸气量，同时也在工作压力要求的范围之内。结合装置工作压力和压力容器罐设计承受压力，试验过程中围压最高取0.6 MPa，为了分析围压变化对时变冲击力和脉动冲击力的影响效果，因此围压按0.1 MPa间隔取6个等级。与不吸气条件相比，装置吸气后瞬时冲击力时均值产生了聚能效应且由靶心向外迅速衰减，因此取靶盘靶心测点的瞬时冲击力作为代表进行分解研究。

2 瞬时冲击力分解方法

瞬时冲击力可分解为

(1)

2.1 离散正交小波变换

小波变换是在克服短时傅里叶变换窗函数大小不变缺点的基础上发展起来的，提供了一个在时频平面上可调的分析窗口，具有自适应的特点，对信号中的快变成分，需要时域分辨率高，可提供窄时间窗，而慢变成分，需要频域分辨率高，可提供宽时间窗，解决了时域分辨率和频域分辨率之间的矛盾[24]。离散正交小波变换如式(2)所示

(2)

2.2 经验模态分解(EMD)

EMD具有时变幅度与时变频率信号的瞬时频率为理论基础，依据信号自身的时间尺度特征将信号分解为多个IMF，与小波变换不同，无需设定基函数，因此，EMD具有自适应的特点，其分解步骤可参见文献[25]，其中停止准则采用限制标准差法，端点处理采用镜像延拓法。

对瞬时冲击力x(t)进行经验模态分解后，可得到有限个IMF及余量，如式(3)所示。其中第1个IMF分量c1(t)瞬时频率最高，随着阶数增加，频率成分逐渐降低，余量rn(t)为频率最低的成分。

(3)

式中：ci(t)为第i个IMF；rn(t)为余量。

时变冲击力一般是变化缓慢的低频部分，不仅仅包括余量，还可能是余量和相邻的几个IMF求和的结果，如式(4)所示。其中m值的确定是提取时变冲击力的关键，表示从第m个IMF分量开始计算求和得到时变冲击力。

(4)

3 离散正交小波变换和EMD分解瞬时冲击力的对比

3.1 噪声分析

为了确定瞬时冲击力噪声的来源，对数字采集仪进行了空采并获取空采信号，由于空采信号具有平稳性，可采用经典谱估计-多锥体法和现代谱估计-修正协方差法对空采信号进行功率谱估计[26]，如图6所示。现代谱估计得到功率谱曲线较经典谱估计平滑，但两种方法得到的噪声频率区间基本一致，噪声主要集中在0～50 Hz，100 Hz，200 Hz和300 Hz，0～50 Hz噪声是由于外界干扰引起的，而100 Hz,200 Hz和300 Hz噪声则是由数学采集仪50 Hz倍频成分引起的，因此需对瞬时冲击力进行消噪处理。

图6 空采信号功率谱Fig.6 Empty signal power spectrum

3.2 离散正交小波变换时变冲击力的提取

3.2.1小波基函数和分解层数的选择

在进行消噪之前，先进行小波基函数和分解层数的选择，目前从非平稳信号提取时变部分的小波基函数主要sym10小波[27]和db10小波，论文将主要对这两种小波基函数在瞬时冲击力分解过程中的适用性进行讨论和分析，运用离散正交小波变换对围压0.4 MPa下瞬时冲击力进行分解，图7和图8给出了分解后8层、9层和10层的时变冲击力提取效果。从图中可知，无论是sym10小波还是db10小波，当分解层数为8层时，时变冲击力还有较快的波动情况，如11～12 s时段，说明低频部分还需要进一步分解；而分解10层时，则时变冲击力在有些时段波动过于平缓，如sym10小波在4～8 s时段，db10小波在12～16 s时段，说明低频部分分解层数过多；当分解层数为9层时，时变冲击力随时间变化的平滑性较好，脉冲周期明显，说明时变冲击力提取效果较好。

图7 sym10小波时变冲击力提取Fig.7 Extraction of time-varying impact force by sym10 wavelet

图8 db10小波时变冲击力提取Fig.8 Extraction of time-varying impact force by db10 wavelet

图9给出了sym10小波和db10小波对时变冲击力提取效果。从图9可知，从整体上db10小波时变冲击力提取效果要略优于sym10小波，如在0～1 s，9～10 s以及17～18 s时段，db10小波得到时变冲击力峰值比sym10小波略微偏左，与瞬时冲击力拟合较好，但在6～7 s和12～13 s时段，则是sym10小波提取效果较好，原因在于，瞬时冲击力随时间变化波形不一样，而两种小波基函数的适用特点又不相同，sym10小波与db10小波具有同样的支撑，但sym10小波波形对称性比db10小波要好，所以当瞬时冲击力随时间变化具有较好的对称性时，sym10小波提取精度较高。另外，当小波变换为9层分解时，采样频率1 500 Hz，则奈斯奎特频率为750 Hz，瞬时冲击力频率分析区间在0～750 Hz，第9层分解得到时变冲击力频率区间在0～2 Hz。时变冲击力的周期约在0.5～1 s，因此频率在1～2 Hz(见图9)。

图9 sym10和db10小波时变冲击力提取效果比较Fig.9 Comparison of time-varying impact force extraction effect by sym10 and db10 wavelet

3.2.2瞬时冲击力小波去噪

在“3.2.1”小节的基础上，将瞬时冲击力减去时变冲击力得到脉动冲击力，仍采用修正协方差法计算脉动冲击力功率谱，修正协方差法中的协方差函数具有两个变量，因此在非平稳信号中也适用[28]，图10(a)给出了脉动冲击力消噪前的功率谱。从图10(a)可知，脉动冲击力频率约为0～50 Hz。结合图6可知，高频部分的100 Hz，200 Hz和300 Hz频率为数字采集系统的电噪声及其倍频，因此，需将瞬时冲击力进行消噪处理，对于100 Hz，200 Hz和300 Hz的频率可用小波强制消噪的方法，即将包括两个频率部分频带的小波系数置为0，而对于低频部分噪声采用软阈值消噪的方法进行处理，对消噪后的脉动冲击力再进行功率谱估计，如图10(b)所示。从图10(b)可知，高频噪声全部都被消掉。图11给出了消噪前后瞬时冲击力对比图。其中黑色代表消噪后瞬时冲击力。

图10 脉动冲击力功率谱Fig.10 Fluctuation impact power spectrum

图11 瞬间冲击力消噪前后对比Fig.11 Comparison of instantaneous impact force before and after denoising

3.3 EMD提取时变冲击力

运用EMD对消噪后的瞬时冲击力分解，得到11个IMF和1个余量，为了得到时变冲击力，需取有限个IMF值和余量之和，图12给出了不同m值下时变冲击力提取效果，此时瞬时冲击力为消噪后的值。从图12可知，与小波分解类似，当m=6时，时变冲击力随时间变化还有较快的波动情况，如约为2 s和10 s的时段，说明m值取值太小；m=8时，时变冲击力随时间变化波动过于平缓，说明m取值太大；当m=7时，此时时变冲击力与瞬时冲击力拟合较好，脉冲周期明显，脉冲效果较好。

图12 EMD时变冲击力提取Fig.12 Extraction of time-varying impact force by EMD

将EMD和db10小波提取的时变冲击力进行对比分析，如图13所示。从图13可知，与EMD相比，离散正交小波变换得到的时变冲击力与瞬时冲击力在整体上拟合地更好，如在0～2 s和16～18 s时段的波峰以及10～11 s时段的波谷，说明小波变换虽然存在小波基函数的选择问题，但如果选择合理，时变冲击力提取精度还是较高的。

图13 EMD和db10小波时变冲击力提取效果比较Fig.13 Comparison of time-varying impact force extraction effect by EMD and db10 wavelet

因此，在围压0.4 MPa条件下，db10小波时变冲击力提取效果要优于sym10小波和EMD。当围压改变时，瞬时冲击力也会发生改变，此时可能还需要重新选择小波基函数或采用EMD。根据上述研究思路，选择低围压0.1 MPa和高围压0.6 MPa条件下瞬时冲击力进行同样的分析，图14和图15给出了两种围压条件下时变冲击力提取效果。在围压0.1 MPa时，sym10与db10小波各有优势，如在0～2 s，2～3 s和9～10 s等时段的波峰，db10小波提取效果要好，但在4～7 s，15～16 s和17～18 s等时段的波峰以及0～2 s时段的波谷，sym10小波提取效果要好，综合考虑波峰和波谷的提取效果，sym10小波提取效果要略优于db10小波；将sym10小波与EMD对比分析，如图14(b)所示。在2～4 s时段的波峰和波谷，EMD提取效果较好，但6～8 s和12～14 s时段的波峰，则sym10小波提取效果较好，整体上看sym10小波提取效果较好。在围压0.6 MPa时，由图15可知，db10小波和sym10小波提取效果相差不大，db10小波略微优于sym10小波，但明显优于EMD法，如在4～6 s和约为18 s时段。从以上分析可知，对于小波基函数，围压0.1 MPa时， sym10小波提取效果较好，而对于0.4 MPa和0.6 MPa围压，则db10小波提取效果较好，三种典型围压下离散正交小波变换提取效果整体上均优于EMD。因此，下面主要分析在不同围压下sym10小波和db10小波对时变冲击力的提取效果。

图14 围压0.1 MPa下时变冲击力提取效果比较Fig.14 Comparison of extraction effect of time-varying impact force under confining pressure 0.1MPa

图15 围压0.6 MPa下时变冲击力提取效果比较Fig.15 Comparison of extraction effect of time-varying impact force under confining pressure 0.6 MPa

3.4 围压对时变冲击力和脉动冲击力的影响

图16给出了三个围压条件下时变冲击力提取效果。可知，与上述围压类似，在围压0.2 MPa、围压0.3 MPa和围压0.5 MPa条件下，db10小波时变冲击力提取效果在某些时段要优于sym10小波，但在另外一些时段提取效果差于sym10小波，但整体上db10小波要略优于sym10小波。

由以上分析可知，在不同围压条件下，除了在围压0.1 MPa条件下选择sym10小波，其余围压条件均选用db10小波来提取时变冲击力。

图16 不同围压下小波变换时变冲击力提取效果比较Fig.16 Comparison of time-varying impact force extraction effect by wavelet transform under different confining pressure

图17给出了不同围压下的时变冲击力，从上到下围压依次为0.1 MPa，0.2 MPa，0.3 MPa，0.4 MPa，0.5 MPa和0.6 MPa，将时变冲击力时间平均值作为瞬时冲击力时均值，即时均冲击力，并与不吸气条件下的瞬时冲击力时均值进行对比，如图18所示。对于不吸气条件，通过对消噪后的瞬时冲击力求时间平均即可得到瞬时冲击力时均值。由图17可知，随着围压的增加，时变冲击力脉冲周期变化不大，但脉冲幅值减小，围压在0.1 MPa条件下，脉冲幅值略大于0.1 MPa，围压在0.2 MPa条件下，脉冲幅值在0.1 MPa左右，当围压大于0.2 MPa时，脉冲幅值均小于0.1 MPa。由图18可知，随着围压的增加，在吸气和不吸气条件下，瞬时冲击力时均值均减小，但在吸气条件下，下降幅度较快，当围压从0.1 MPa变化到0.3 MPa时，下降幅度较大，为0.31 MPa，当围压从0.3 MPa变化到0.6 MPa时，则下降幅度较小，为0.14 MPa。与不吸气相比，装置吸气后瞬时冲击力时均值明显增加，且增加幅度随着围压的增加而减小。

图17 不同围压下时变冲击力比较Fig.17 Comparison of time-varying impact force under different confining pressure

图18 不同围压下瞬时冲击力时均值Fig.18 Time-varying value of instantaneous impact force under different confining pressure

按照“3.2”节所述方法对瞬时冲击力进行去噪处理并提取不同围压下脉动冲击力，如图19所示。以脉动冲击力均方根作为脉动冲击力脉动强度，不仅可以表示脉动冲击力平均能量，也可以反映脉动幅值，图20给出了不同围压下装置吸气和不吸气条件下的脉动冲击力脉动强度，对于不吸气条件，只需利用离散正交小波变换对瞬时冲击力进行消噪处理，减去瞬时冲击力时均值后，得到脉动冲击力，再计算其脉动强度。结合图19和图20可知，随着围压的增加，脉动冲击力脉动强度逐渐减小，在吸气条件下，脉动强度从0.07 MPa下降到0.02 MPa，当围压从0.1 MPa变化到0.3 MPa时，脉动强度下降幅度较大，为0.037 MPa，当围压从0.3 MPa变化到0.6 MPa时，则下降幅度较小，为0.015 MPa。与不吸气相比，装置吸气后脉动强度在围压0.1～0.4 MPa条件下增加，但在围压0.5～0.6 MPa条件下减小，原因在于脉动强度大小不仅取决于射流时均速度场和脉动速度场，而且与射流流体密度有关，装置吸气后时均速度增加，但是水气混合密度减小，在0.1～0.4 MPa围压下，时均速度和脉动速度增加的程度大于水气混合密度减小的程度，导致脉动强度变大，而在0.5～0.6 MPa围压下，时均速度和脉动速度增加影响程度小，脉动强度变小。

图19 不同围压下脉动冲击力Fig.19 Fluctuation impact force under different confining pressure

图20 不同围压下脉动冲击力脉动强度Fig.20 Fluctuation strength of fluctuation impact force under different confining pressure

4 结 论

论文提出了利用离散小波变换和EMD来提取水下自激吸气式脉冲射流装置瞬时冲击力中的时变冲击力和脉动冲击力，得到以下结论：

(1)离散正交小波变换结合功率谱估计可以有效地消除瞬时冲击力中所含的噪声。

(2)通过离散正交小波变换和EMD两种方法对时变冲击力提取效果的分析，确定了合适的小波基函数，证明了离散正交小波变换对时变冲击力提取效果的优越性，在围压0.1 MPa条件下，sym10小波提取效果较好，而在其它围压条件下，db10小波提取效果较好。

(3)在不同围压条件下，装置吸气后时变冲击力脉冲周期基本保持一致，当围压为0.1MPa和0.2MPa时，脉冲幅值在0.1 MPa以上，其它围压条件下则脉冲幅值在0.1 MPa以下；在吸气条件下，随着围压的增加瞬时冲击力时均值和脉动冲击力脉动强度减小的幅度越来越小，当围压小于0.3 MPa时减小幅度大，围压大于0.3 MPa时减小幅度小，而在不吸气条件，两者下降幅度较为平缓。在不同围压条件下，与不吸气条件相比，装置吸气后瞬时冲击力时均值均变大，在低围压条件下吸气后脉动冲击力脉动强度变大，高围压条件下则变小。

论文借助离散小波变换和EMD对水下自激吸气式脉冲射流装置的瞬时冲击力进行了有效地消噪和准确地分解，从时变冲击力和脉动冲击力可以得到瞬时冲击力的脉冲幅值、脉冲周期和脉动强度，下一步需采用时频分析技术对脉动冲击力的时频特性进行研究。另外，瞬时冲击力分解是水下自激吸气式脉冲射流装置外部冲击特性研究范畴，未来需结合装置内部水气两相流动特性研究来揭示装置脉冲射流的产生机理。