夏振炎，靳秀青，姜 楠

(1. 天津大學力學系，天津 300072；2. 天津市現(xiàn)代工程力學重點實驗室，天津 300072)

基于小波及Hilbert-Huang變換提取壁湍流相干結構

夏振炎1,2，靳秀青1,2，姜 楠1,2

(1. 天津大學力學系，天津 300072；2. 天津市現(xiàn)代工程力學重點實驗室，天津 300072)

為了研究壁湍流中的相干結構，在直流抽吸式低(變)湍流度風洞的實驗段中，應用恒溫熱線風速儀 IFA-300精確測量平板湍流邊界層不同法向位置的流向速度信號．分析壁湍流脈動速度的概率密度分布，并與高斯分布比較，驗證了壁湍流中相干結構引起的猝發(fā)現(xiàn)象的統(tǒng)計規(guī)律．提出了一種利用小波變換和 Hilbert-Huang(HHT)變換相結合來提取相干結構的方法．對壁湍流脈動速度進行 HHT變換，通過變換后各模態(tài)的能量分析，分離出主要含能模態(tài)對應的相干結構．利用小波重構法繼續(xù)提取未分離的小尺度相干結構．對剩余湍流信號進行概率統(tǒng)計，并比較提取的全部相干結構與原始湍流信號的波形和 HHT邊際譜，驗證了湍流信號既包含非隨機的相干結構信號也包含隨機的非相干信號．

Hilbert-Huang變換(HHT)；小波變換；壁湍流；相干結構；間歇性

壁湍流一直是實際流場中備受關注的熱點問題，壁湍流脈動速度的概率密度分布是力學的基本問題之一．Brodkey等[1]依據(jù)實驗資料分析了剪切流中湍流信號的分布，得到的脈動速度分布與高斯分布有一定偏離．Kim 等[2-3]應用直接數(shù)值模擬和大渦模擬方法對槽道湍流進行了細致的研究，得到脈動速度的偏斜度和平坦度也與正態(tài)分布不符．脈動速度概率密度分布、平坦度和偏斜度非正態(tài)分布的原因與壁面處的相干結構有關[4]．而相干結構出現(xiàn)在湍流信號的各種頻段上，其在大尺度結構(低頻部分)中表現(xiàn)為較典型的猝發(fā)現(xiàn)象，而在小尺度結構(高頻部分)中同樣也存在相干結構[5-6]，表現(xiàn)為強間歇事件．

應用小波變換中的小波系數(shù)重構法對相干結構的提取取得了較好的效果[7]．1998年，黃鍔提出了經(jīng)驗模態(tài)分解方法(EMD)．該方法不需預設基函數(shù)就能將信號中不同尺度的波動逐級分解，產(chǎn)生不同特征尺度的數(shù)據(jù)序列．該方法通過波內頻率調制和波間頻率調制實現(xiàn)其非線性的分析過程[8]．研究表明，湍流信號中含能渦呈現(xiàn)明顯的多尺度特征[9]．用 HHT分析湍流的分尺度能量時，能量更加集中于低頻(高模態(tài))部分，更有利于大尺度相干結構的提取．

筆者首先對實驗數(shù)據(jù)進行概率密度統(tǒng)計，并與高斯分布比較，驗證了由相干結構引起的猝發(fā)現(xiàn)象的統(tǒng)計規(guī)律．然后將湍流脈動速度進行 HHT變換分析，將大尺度(高模態(tài))相干結構提取出來．應用小波變換中小波系數(shù)重構的方法繼續(xù)處理剩余信號，將小尺度相干結構提取出來．最后將提取出的全部相干結構信號與對應壁面距離處的壁湍流脈動速度進行相關分析并進行比較．

1 實驗裝置與測量方法

實驗在直流抽吸式低(變)湍流度的閉口風洞中進行，在該風洞實驗段中放置一塊長 4.5,m，寬0.35,m，厚 5,mm，前緣為楔形的平板．平板上游粘貼砂紙以加速轉捩，在距平板前緣 2.4,m的實驗測量位置處放置一根標定好的熱線探針，探針由三維坐標架固定并與恒溫熱線風速儀IFA-300相連接，采集的數(shù)據(jù)由計算機進行處理．實驗裝置如圖1所示．

圖1 實驗裝置Fig.1 Experimental setup

對風洞內壁湍流脈動速度進行測量時，熱線探針的初始位置(見圖 1)為距平板前緣 2.4,m處，位于平板上方 0.3,mm的高度，考慮到邊界層流場的特點，以不同的增長間距共設置了84個不同法向位置的測量點．將測得的數(shù)據(jù)文件應用 Matlab編程處理，進行相關湍流分析．

2 實驗數(shù)據(jù)分析

采用概率統(tǒng)計方法分析脈動速度分布特征與壁湍流猝發(fā)的關系；利用 HHT和小波重構法提取含能的全尺度壁湍流相干結構，并分析其時頻分布特性．

2.1 脈動速度的統(tǒng)計特性分析

將實驗測量的湍流脈動速度做概率統(tǒng)計分析，取2個不同的壁面距離y+=10.0和y+=610.3，分別做概率密度分布(PDF)曲線，將其歸一化后與標準高斯曲線進行比較，如圖2所示．

圖2 湍流信號概率密度曲線與標準高斯線的比較Fig.2 Comparison of turbulent probability density curve with standard Gaussian line

從圖2中可以看出，在近壁區(qū)域y+=10.0處概率密度曲線較高斯線正偏，而在較遠壁區(qū)域 y+=610.3處概率密度曲線較高斯線負偏．流向脈動速度的這種分布特征與壁面附近湍流猝發(fā)的掃掠和噴發(fā)事件有關．遠壁區(qū)脈動速度主要受低速流體噴發(fā)事件的影響，因而流向脈動速度分布負偏；而在壁面附近區(qū)域，脈動速度主要受制于高速流體的掃掠運動，因而流向脈動速度概率分布正偏．這是由湍流邊界層中有組織的相干結構運動引起的典型現(xiàn)象．

湍流隨機變量的分布特征還常用偏斜度(skewness)來描述．偏斜度用 S表示，表征隨機變量概率密度函數(shù)的對稱性，即S=0表示概率密度函數(shù)f關于均值中心對稱；若 S＜0，則概率密度函數(shù)左偏或負偏；若S＞0，則概率密度函數(shù)右偏或正偏．

為了更深入地了解實驗測量數(shù)據(jù)的統(tǒng)計特性，研究了偏斜因子隨壁面距離的變化規(guī)律，如圖 3所示．容易看出，在壁面距離由小變大的過程中，偏斜因子的值由正值變?yōu)樨撝担@與圖 2近壁區(qū)脈動速度概率密度曲線由正偏過渡為遠壁的負偏相一致，說明近壁掃掠事件占優(yōu)勢，雖然從局部掃掠階段的后期到局部猝發(fā)的開始之間時間很短，但是大部分雷諾應力是在這段時間產(chǎn)生的；而到了遠壁，低速流體向后上方噴射進入邊界層外層區(qū)域，這導致反流向的脈動速度成分增加，流向脈動速度被削弱．猝發(fā)現(xiàn)象說明湍流脈動不是隨機過程，湍流脈動信號中存在具有一定內在規(guī)律的擬序結構，亦即相干結構．

圖3 由近壁到遠壁偏斜因子的變化曲線Fig.3 Variation of skew factor with wall distance

2.2 應用HHT方法提取大尺度相干結構

1998年，黃鍔提出了經(jīng)驗模態(tài)分解方法，Huang等認為，實際信號基本上都是多成分信號，即同一時刻信號包含多種頻率成分，某一時刻對應的瞬時頻率沒有明顯的物理意義．該方法將湍流信號進行平穩(wěn)化處理，對信號中不同尺度的波動或趨勢項逐級分解開來，產(chǎn)生一系列不同特征尺度的數(shù)據(jù)序列，每個序列稱為一個固有模態(tài)函數(shù)(IMF)分量．對分解得到的分量進行 Hilbert變換得到信號在頻域或時域上的分布規(guī)律．為了使得瞬時頻率具有明顯的物理意義，Huang等[10]在原有特征尺度分解的基礎上，提出了內模函數(shù)的概念和經(jīng)驗模態(tài)分解方法．Huang用內模函數(shù)來近似單一成分信號，并做出如下定義：

(1) 極值點的個數(shù)和過零點的數(shù)目相等，或者最多相差1；

(2) 在任一點，由極大值構成的包絡和由極小值構成的包絡的平均值為零．

其中條件(1)將信號限制為窄帶信號，條件(2)可以保證求得的 IMF瞬時頻率有意義．這2個條件保證了后續(xù) Hilbert變換的可行性．凡滿足以上2個條件的分量就看作是一個內模函數(shù)．為了從原始數(shù)據(jù)中分解出內模函數(shù)，Huang也給出了具體的經(jīng)驗模態(tài)分解方法[10]．

應用該方法將湍流脈動速度信號做 Huang變換分解成多個模態(tài)，再做Hilbert變換，求得不同壁面距離下各模態(tài)的能量變化．需要說明的是 Hilbert-Huang變換過程中，由于內模函數(shù)由極值包絡的對稱性來判定，第一模態(tài)對應的分量包絡波動性最大，相應頻率范圍也最大．模態(tài)越高對應的頻率范圍越窄，值也越小，即高模態(tài)對應低頻率段，低模態(tài)對應高頻率段．本文選取 4個不同壁面距離的湍流脈動速度，分別為y+=8.0、y+=30.0、y+=61.4和y+=162.5．圖4所示為HHT分解各模態(tài)能量變化分布曲線．

圖4 HHT分解模態(tài)的能量變化分布曲線Fig.4 Modal energy of HHT

圖 4中，橫軸為 HHT分解的模態(tài)，縱軸為各模態(tài)能量的歸一化結果．可以看出，高模態(tài)(低頻)部分含能較高，這首先與湍流脈動信號本身有關，因為低頻率部分代表著大尺度湍渦結構，而大尺度渦結構構成相干結構的主要部分，含能也最高．其次，流場的不穩(wěn)定性和剪切作用，使來自平均流場的湍能沿低頻(高模態(tài))部分逐級向高頻(低模態(tài))部分傳遞[11]，湍能依靠湍渦之間的相互作用進行傳輸，同時該傳輸過程存在著反向能量級串[9]，這使得湍能在傳遞過程中逐級減少，呈現(xiàn)在圖中則表現(xiàn)為隨著模態(tài)的增加，能量也逐級增加．再次，HHT分解方法以 IMF當作展開的基底，分解完全基于信號本身的性質．將IMF分量作 Hilbert變換對應的模態(tài)能量更能代表真實頻段湍流信號的能量．研究表明，HHT分解信號時，能量更加集中于低頻部分[8，12]，因此 HHT分解更適合提取大尺度的相干結構．

統(tǒng)計主要含能的大尺度低頻模態(tài)，這些高能量模態(tài)所含能量占該壁面距離處湍流信號總能量的比例如表1所示．

表1 高模態(tài)能量eh占總能量et百分比Tab.1 Percentage of high modal energy share the total energy %

由表 1看出，高模態(tài)能量所占總能量比例可達到 80%左右，占總能量的絕大部分．分析高含能模態(tài)并將其與傳統(tǒng)方法提取出的相干結構進行比較，結果表明上述高含能模態(tài)與相干結構的特性近似相同，是相干結構的主要部分．基于上述結論，將高模態(tài)的IMF分量從原信號中分離出來，將剩余模態(tài)分量進行重構．重構后的剩余信號已去除主要含能的湍流相干結構，為近似隨機的湍流信號．分析重構后剩余信號的概率密度分布曲線，并與標準高斯曲線進行比較．圖 5所示為 y+=30.0壁面距離處剩余信號的概率密度分布與標準高斯分布的比較．

圖5 重構信號概率密度曲線與高斯線比較Fig.5 Comparison of the probability density curve reconstructed with the Gaussian line

由表1可知，y+=30.0處分離出的大尺度相干結構能量占總能量的 76.67%，剩余信號能量僅占總能量的23.33%．從圖5中可以看出，重構后的剩余信號概率密度曲線與標準高斯曲線相比，略微存在一些差別．原因在于相干結構分布于湍流信號的整個頻段內，由于 HHT分解具有“低頻能量集中”特性，利用 HHT將高模態(tài) IMF分量分離是將主要含能的大尺度(低頻段)湍渦結構去除，而對于全部相干結構包含的小尺度(高頻段)部分相干結構，HHT方法未對其進行提?。?/p>

2.3 小波變換提取小尺度間歇事件

相干結構具有多尺度分布特征，不僅存在于大尺度結構中，也存在于小尺度結構中[7]，其中大尺度(低頻率)部分中的相干結構已在第 2.2節(jié)中采用 HHT方法從原始信號中分離出來，剩余的非隨機成分為小尺度相干結構，這部分相干結構表現(xiàn)為可由間歇因子表征的含能較多的強間歇結構．所謂間歇因子是根據(jù)小波系數(shù)引入的表征間歇性的量化指標，其計算公式為

式中：I(r，x)為小波分析中各尺度的局部間歇強度；w(r，x)為不同壁面距離、不同尺度下的小波系數(shù)；r為壁面距離；x為尺度．

將第 2.2節(jié)中的重構信號采用小波變換進行多尺度分解．為了減少能量泄露，對比多個小波函數(shù)，最后選用“db20”作為小波基函數(shù)，使能量泄露最少．由于小尺度相干結構含有更多的能量，分析實驗數(shù)據(jù)的統(tǒng)計結果，選取小波分解后各尺度小波系數(shù)間歇因子的均值作為判定強間歇事件的門限值 L，其計算公式為

式中：I(r，x(i))為同一壁面距離同一尺度下的各小波系數(shù)對應的間歇因子；L(x)為第x尺度下由該尺度n個間歇因子決定的門限值．

對強間歇事件以外的各尺度小波系數(shù)均設為零，然后經(jīng)過小波重構后可以提取出具有強間歇性的小尺度相干結構．同樣地，將強間歇事件對應的各尺度小波系數(shù)置為零，各尺度小波系數(shù)信號重構得到去除相干結構之后的剩余信號．圖6所示為y+=30.0時，剩余信號的概率密度分布曲線與標準高斯曲線進行比較，2條曲線已非常接近，表明經(jīng)過HHT以及小波變換處理后，已將絕大部分湍流相干結構提取出來，剩余部分為近似高斯型的帶有隨機特征的湍流信號．比較圖2和圖5，可以看出經(jīng)過HHT后，已經(jīng)將大部分含能最多的大尺度相干結構從原湍流信號中分離出來，但仍留存有部分小尺度湍流信號．而原始湍流信號u(t)可表示為

圖6 小波處理后概率密度曲線與高斯線比較Fig.6 Comparison of probability density curve continuely treated by wavelet with the Gaussian line

式中：ulc(t)為大尺度湍流信號，該部分已由 HHT方法提取出來；usc(t)為小尺度間歇事件湍流信號，該部分由小波系數(shù)重構的方法提取完成；ur(t)為非相干湍流成分，即為圖6中的剩余信號．度(高頻)部分的相干結構后，經(jīng)小波系數(shù)重構獲得的湍流脈動速度的概率密度分布曲線與標準高斯曲線基本吻合，說明在經(jīng)過小波變換對小尺度強間歇事件進行提取后，剩余的湍流接近隨機信號，因此壁湍流相干結構與原始湍流信號完全分離．

2.4 提取的相干結構與原始湍流信號的比較

為了進一步研究提取出的相干結構的性質，將采用HHT方法提取出的大尺度相干結構信號與提取出的小尺度相干結構引起的間歇事件信號疊加，以y+= 30.0為例，獲得湍流脈動信號中全部相干結構的信號(見圖7)，并與原始湍流信號(圖8)比較．

對比圖7和圖8可以看出，提取前后脈動較強部分分布的時間比較一致，只是脈動速度的大小有一定差別，原始湍流脈動速度較大，這是由于湍流信號除了包含相干結構外還包含有隨機的非相干成分．

圖7 提取的相干結構信號(y+=30.0)Fig.7 Extracted coherent structure(y+=30.0)

圖8 原始湍流信號(y+=30.0)Fig.8 Original turbulent signal(y+=30.0)

分析相干結構信號的 HHT邊際譜，并與原始湍流脈動速度信號的 HHT邊際譜進行比較．所謂邊際譜是從統(tǒng)計意義上表征整組數(shù)據(jù)每個頻率點的累積幅值分布．用于處理湍流這樣的瞬態(tài)信號時，如果信號中出現(xiàn)某一頻率的能量，就表示一定存在該頻率的振動波，故邊際譜可以更準確地反映湍流信號的實際頻率成分，更適合湍流信號的分析．圖9所示為y+= 30.0處提取出的相干結構HHT邊際譜和原始湍流脈動速度HHT邊際譜的比較．

圖9中，橫軸為頻率，縱軸為對應頻域下信號能量(Abs.,magnitude)的歸一化值．可以看出，較低頻率(大尺度)部分提取相干結構前后的2條曲線基本重合；較高頻率(小尺度)部分提取相干結構后曲線快速衰減，說明主要含能的大尺度相干結構集中于低頻部分，而高頻部分含有較多非相干信號，這與Farge等[13]應用數(shù)值模擬方法提取相干結構獲得頻譜圖的趨勢相同．

圖9 相干結構與原信號的邊際譜比較(y+=30.0)Fig.9 Comparison of HHT marginal spectrum between coherent structure and the original signal (y+=30.0)

3 結 論

(1) 近壁區(qū)脈動速度概率密度曲線向正方向偏離高斯曲線，在稍遠壁區(qū)正好相反，說明由相干結構引起的典型“猝發(fā)”現(xiàn)象導致湍流脈動速度分布具有分層特征．

(2) 基于 HHT方法本身的優(yōu)點，其對低頻部分的相干結構提取更準確．且該提取部分能量約占湍流信號總能量的 80%．用小波變換中小波系數(shù)重構方法將剩余信號中的強間歇事件進行提取．最后剩余信號的概率密度分布曲線與標準高斯曲線吻合較好，相干結構得到了較完全的提?。?/p>

(3) 提取出的全部相干結構數(shù)據(jù)與原始湍流脈動速度比較，驗證了湍流信號既包含非隨機的相干結構信號也包含隨機的非相干信號．將相干結構數(shù)據(jù)的 HHT邊際譜與原始湍流信號 HHT邊際譜比較，在較低頻率部分提取相干結構前后的 2條曲線基本重合，而對于較高頻率部分提取相干結構后曲線快速衰減，2條曲線逐漸分離．

［1］ Brodkey R S，Wallace J M，Eckelmann H. Some properties of truncated turbulence signals in bounded shear flows[J]. Journal of Fluid Mechanics，1974，63(2)：209-224.

［2］ Kim J，Moin P，Moser R. Turbulence statistics in fully developed channel flow at low Reynolds number[J]. Journal of Fluid Mechanics，1987，177：133-166.

［3］ Moin P，Kim J. Numerical investigation of turbulent channel flow [J]. Journal of Fluid Mechanics，1982，118：341-377.

［4］ 周濟福，張 強，李家春. 近壁湍流脈動的概率分布函數(shù)[J]. 應用數(shù)學和力學，2005，26(10)：1135-1143.

Zhou Jifu，Zhang Qiang，Li Jiachun. Probability distribution function of near-wall turbulent velocity fluctuations [J]. Applied Mathematics and Mechanics，2005，26(10)：1135-1143(in Chinese).

［5］ Douady S，Couder Y，Brachet M E. Direct observation of the intermittency of intense vorticity filaments in turbulence[J]. Physics Review Letters，1991，67(8)：983-986.

［6］ She Z S，Jackson E，Orszag S A. Intermittent vortex structures in homogeneous isotropic turbulence[J]. Nature，1990，344(15)：226-228.

［7］ 姜 楠，柴雅彬. 用子波系數(shù)概率密度函數(shù)研究湍流多尺度結構的間歇性[J]. 航空動力學報，2005，20(5)：718-724.

Jiang Nan，Chai Yabin. Experimental investigation of multi-scale eddy structures′ intermittency in turbulent flow using probability density function of wavelet coefficients[J]. Journal of Aerospace Power，2005，20(5)：718-724(in Chinese).

［8］ 黃永祥. 圓柱繞流和大氣湍流的Hilbert-Huang變換分析[D]. 上海：上海大學力學所，2006.

Huang Yongxiang. Circular Cylinder and Atmospheric Turbulence Analyzed by Hilbert-Huang Transform[D]. Shanghai： Institute of mechanics， Shanghai University，2006(in Chinese).

［9］ Xia Zhenyan，Tian Yan，Jiang Nan. Wavelet spectrum analysis on energy transfer of multi-scale structures in wall turbulence[J]. Applied Mathematics and Mechanics：English Edition，2009，30(4)：435-443.

［10］ Huang N E，Long S R，Shan Z. The mechanism for frequency downshift in nonlinear wave evolution[J]. Advances in Applied Mechanics，1996，32：59-117.

［11］ Xia Zhenyan，Jiang Nan，Tian Yan，et al. Experimental study on the local similarity scaling of the turbulence spectrum in the turbulent boundary layer[J]. Science in China Series：Physics，Mechanics & Astronomy，2009，39(3)：428-434.

［12］ 連海寧. 傅里葉、小波及希爾波特-黃變換在地震工程領域的應用研究[D]. 北京：中國地震局工程力學研究所，2004.

Lian Haining. Fourier，Wavelet and Hilbert-Huang Transform Applied in the Field of Earthquake Engineering Research[D]. Beijing：Institute of Engineering Mechanics of China Earthquake Administration，2004(in Chinese).

［13］ Farge M，Schneider K，Pellegrino G，et al. Coherent vortex extraction in three-dimensional homogeneous turbulence：Comparison between CVS-wavelet and PODFourier decompositions[J]. Physics of Fluids，2003，15(10)：2886-2896.

Extraction of Coherent Structures in Wall Turbulence Based on Wavelet and Hilbert-Huang Transforms

XIA Zhen-yan1,2，JIN Xiu-qing1,2，JIANG Nan1,2
(1. Department of Mechanics，Tianjin University，Tianjin 300072，China；2. Tianjin Key Laboratory of Modern Engineering Mechanics，Tianjin 300072，China)

In order to study the coherent structures in wall turbulence，the constant temperature hot-wire anemometer IFA-300 was used to measure the flow velocity signals at different normal locations of the turbulent boundary layer in the wind tunnel，which is of direct suction-type low(variable)turbulent intensity. The probability density distribution of wall turbulence fluctuating velocity was analyzed and compared with the Gaussian curves，and the results confirm the statistical laws of bursts caused by coherent structures in wall turbulence. A method for extracting coherent structures was proposed，which is hybridized by wavelet transform and Hilbert-Huang transform(HHT). The HHT was applied to analyze the turbulent fluctuation velocity，and the main coherent structures were extracted by analyzing the energy of each mode. The wavelet reconstruction method was used to extract small-scale coherent structures from the left data. The probability density distribution of the remained signals was analyzed. All of the coherent structures extracted by the above method were compared with the original turbulent signal，and their HHT marginal spectrums were also compared，which verifies that turbulence signal contains both non-random coherent signal and random noncoherent signal.

Hilbert-Huang transform；wavelet transform；wall turbulence；coherent structures；intermittency

O357.5

A

0493-2137(2012)04-0373-06

2010-11-08；

2011-04-12.

國家自然科學基金重點資助項目(10832001)；國家自然科學基金資助項目(10872145，11172205).

夏振炎（1969— ），男，副教授.

夏振炎，xia_zhy@tju.edu.cn.