繞圓盤空化流動的拉格朗日擬序結構分析

白澤宇，王國玉，吳欽，黃彪

(北京理工大學 機械與車輛學院，北京100081)

0 引言

超空泡技術是一種高效的水下航行器減阻方式，其通過將航行器包裹在空泡中以降低粘性阻力從而實現水下航行器的超高速運動。隨著俄羅斯暴風雪高速魚雷的出現，超空泡航行器及超空泡武器逐漸成為各國研究重點。空化器是水下超空泡航行器的重要組成部分之一，它是空泡起始的地方，不僅影響航行器阻力，也為航行器提供升力和姿態控制［1］。隗喜斌等［2］對錐形空化器的空泡形態進行了分析，發現空泡變化具有延遲性和波動性。孟慶昌等［3-4］考慮了水的壓縮性對圓盤空化器超空泡流場的影響，基于勢流理論計算了亞聲速和超聲速圓盤空化器的超空化流動。

空化是一種復雜的氣液兩相流動，涉及兩相之間的相互作用和相界面的非定常變化等問題，多相流模型的選取對空化流動的數值計算至關重要。目前，空化流動的計算多采用基于均質平衡流理論的均相流模型，即將氣液兩相視為一種介質求解控制方程。如袁緒龍等［5］基于均相流假設，計算并分析了3 種典型空化器的空化流動;時素果等［6］建立了一種基于狀態方程的當地均相介質模型計算了繞圓盤空化器的通氣超空化流動。采用均相流模型無法準確預測流動的非定常特性，往往需要對湍流模型加以修正，如Johansen 等［7］、黃彪等［8］的研究。近年來，隨著計算機技術的發展，非均相流模型逐漸應用到多相流數值計算中，該方法分別對各相求解控制方程，把相間界面看作一個移動的邊界，并考慮相間的相互作用，能更好地反映流動細節和非定常特性。Boisson 等［9］運用非均相流模型對二維方腔的氣泡上升過程進行了數值模擬，部分結果與實驗結果吻合較好;Chahed 等［10］運用非均相流模型對泡狀流進行了計算，考慮了湍流輸運引起的附加質量力變化，得到的速度場與實驗結果吻合較好。目前對于空化流動的數值計算，非均相流模型的應用還不廣泛，國內外文獻中相關報道很少。

渦旋結構的產生可以強化流場內部的湍流強度，尤其是在空化現象中，從而影響空化現象的發展［11］;在高雷諾數下，空化現象的產生，又會使得湍流場內部產生復雜的渦旋結構，增加湍流場和空泡之間的交互作用，使流動變得更加復雜［12］。因此，渦旋結構的研究對空化現象的理解至關重要。渦旋結構的辨識與分析手段較多，使用最為廣泛的是基于歐拉體系的渦量判據、Q 判據［13］、Δ 判據［14］、λ2判據［15］等準則，這些判據均由速度梯度張量中的分量構成，需要對速度梯度張量進行精確計算，因而對速度場的精度要求較高;同時，這些判據需要人為選定閾值才能對渦結構進行較好捕捉，這給未知流場的渦結構識別帶來了困難。近年來，基于拉格朗日體系的流動結構辨識方法逐漸發展起來。Haller等［16］將流場看做離散粒子組成的動態系統，識別出了流場中的“拉格朗日擬序結構(LCS)”.LCS 的識別基于“有限時間Lyapunov 指數(FTLE)”的求解，當相鄰粒子以不同運動特性運動時，會導致FTLE場中出現突出的“脊”的結構——LCS，Haller［17］及Shadden 等［18］對此方法進行了詳細描述。LCS 將具有不同運動特性的流場區域分離開來，因此能很好地捕捉渦旋結構。相比傳統的歐拉判據，該方法無需求解速度梯度張量，對速度場的空間分辨率要求低，對其中的異常值不敏感;同時不需要人為設定閾值，因而能反映真實的流動結構。Green 等［19］使用FTLE 方法分析了三維槽道湍流流動，成功辨識了馬蹄渦結構。Farrell 等［20］基于LCS 提出了在軸對稱射流中的渦結構判據準則。目前，FTLE 方法及LCS在空化流動分析中的應用才剛剛起步，Tang 等［21］用該方法分析了繞水翼非定常空化的演化過程，但對渦結構捕捉效果不理想。其他文獻中暫未見到該方法在空化流動研究中的應用。

本文采用非均相流模型計算了繞圓盤空化器的非定常空化流動，通過求解FTLE 場識別了流場中的LCS，并結合傳統的歐拉分析方法對繞圓盤空化器空化流場渦旋結構的發展規律及流動機理進行了分析。

1 非均相流模型

非均相流模型對每一相單獨求解，需要考慮相間傳輸和相間作用力。第n 相的控制方程為

連續性方程:

式中:αn為第n 相的體積分數為第n 相中的質量源;Smn'為單位體積流體中從n'相到n 相的質量傳輸流量，本文中的質量傳輸由空化引起。

動量方程:

式中:μn為第n 相的分子粘性系數;μt為湍流粘性系數;αnρnfi為第n 相中外部體積力;Mni為其他相對第n 相的作用力;為由于相間質量輸運引起的動量輸運項。

以上控制方程采用標準k-ε 湍流模型封閉［22］，空化模型采用Kubota 空化模型［23］。

2 拉格朗日擬序結構

2.1 Lyapunov 指數

Lyapunov 指數是描述混沌現象的一個重要參數，被用來度量相空間中兩條相鄰軌跡隨時間按指數率收斂或發散的程度［21］。對于初始位置無限接近的兩個質點x0和x0+δx0(δx0趨于無窮小)，他們的運動軌跡間發散的平均指數率，即Lyapunov 指數，可通過(3)式計算得到:

式中:δx(x0，t)表示質點的離散程度，其值與質點初始位置x0和時間t 有關。σ ＜0 時，表示相鄰點最終靠攏為一點，對應于不動點和周期運動;σ=0 時，系統維持在穩定分離狀態;σ ＞0 時，表示相鄰點最終將分離，對應于軌跡的局部不穩定性。

2.2 FTLE 和LCS

FTLE 是描述初始時刻系統空間中兩個無限接近的流體質點的相對距離在有限時間范圍內的平均變動率。x0處的有限時間柯西-格林變形張量［18］定義為

式中:tLE為運動時間;x(t0+tLE;t0，x0)為t0時刻位于x0位置的質點運動 tLE時間后所處位置;為變形梯度張量。FILE 定義為

式中:λmax((x0))定義為(x0)的最大特征值。FTLE 場是一個標量場，它反映了當前時刻流體的特性。當相鄰粒子以不同運動特性運動時，會導致FTLE 場中出現突出的“脊”的結構，這種結構稱之為LCS.利用LCS 可以分割不同運動特性的流場，從而捕捉流場潛在的動態力學和幾何學特性。

3 結果與討論

本文所取的計算流域和實驗段［24］完全相同，為190 mm×70 mm×700 mm，計算采用的圓盤與實驗［24］中的幾何尺寸相同，其直徑為D=15 mm，圓盤后部采用倒錐結構，倒錐角為120°;圓盤后有一支桿，直徑d =6 mm.圖1給出了計算區域其網格劃分。圓盤前端的區域采用C 型結構化網格劃分，以較好匹配圓盤頭部的形狀，在圓盤鈍體周圍近壁區域進行了網格加密，近壁面y +值為20 ～80 之間，滿足壁面函數要求，網格總數為230 000.文獻［6］中對本文計算工況做了網格無關性驗證，得出網格總數大于180 000 時，空化器阻力系數差別較小。考慮到求解FTLE 場對速度場樣本點數量的要求，本文選取的網格總數大于文獻［6］.采用速度入口作為邊界條件，所取速度ui=10.5 m/s 也和實驗工況一致，此時以圓盤直徑為特征長度的雷諾數為1.79×105.出口采用壓力出口，圓盤表面采用無滑移固壁條件。

圖1 計算域及網格劃分Fig.1 Computational domain and gir d

數值計算的時間步長取5×10-4s，計算FTLE場時采取的積分時間tLE為0.037 5 s(75 倍時間步長)。

3.1 計算模型驗證

表1給出了采用非均相流模型計算得到的空穴形態及截面含氣量分布的非定常變化過程，并與實驗結果進行了對比，T 為數值計算得到的空穴發展周期。由表1可以發現，采用非均相流模型計算得到的結果與實驗結果吻合較好，二者表現出相似的非定常特性。在計算工況下，繞圓盤的空化流動處于典型的云狀空化階段，在圓盤后部產生了一個環狀空穴，空穴前端為透明氣相區，中后部為水氣混合區，如表1中圓圈所示。兩個區域之間存在強烈的相互作用，在反向射流的作用下，透明氣相區的范圍周期性增大和減小;整個云狀空泡團表現出較為規律的非定常特性，空泡繞支桿呈螺旋狀向后發展，支桿四周的環狀空穴交替性的膨脹收縮，收縮時發生斷裂脫落，表現出明顯的周期性，數值計算得到的周期為T=21 ms，與文獻［24］實驗值19 ms 較為接近。

表1 實驗結果與數值計算非定常空穴形態結果對比Tab.1 Experimental and calculated unsteady cavity forms

3.2 FTLE 分布與LCS 分析

表2給出了截面含氣量αg分布、Q 分布和FTLE分布。Q 分布為歐拉體系下的渦判據準則，其定義為Q=1/2(|Ω|2- |S|2)［21］，式中:Ω 為渦量張量，反映旋轉效應;S 為應變速率張量，反映變形效應。Q ＞0 的區域定義為渦，即在渦旋的區域內流體的旋轉比較應變率大小而言起主導作用。

表2中黑色實線為流線，用以直觀表現在當前時刻下的流動狀態和渦旋結構。由表2可以發現，繞圓盤的空化流動屬于典型的渦空化，在圓盤后部的空化區域內存在一個大尺度渦旋結構，圍繞在支桿四周。由含氣量分布圖可以發現，渦心的位置處于環狀空穴的閉合端，隨著支桿四周空穴交替性的膨脹收縮，渦心的位置也交替前后運動，且渦心的運動相對于空穴的運動有所延遲，t0+2/6T 時刻，支桿上側的空穴膨脹到最大，下側的空穴收縮到最小，而渦心在t0+3/6T 時刻才運動到上下游的極限位置。由Q 分布圖可以發現，在圓盤前端，有一個Q ＜0 區域，說明該區域內，流體的變形效應占據主導地位，從流線圖中可以看出該區域沒有形成渦結構;圓盤后部，存在Q ＞0 的區域，說明該區域內，流體的旋轉效應占據主導地位，存在渦旋結構，該區域與空穴覆蓋區域重疊較好。根據Q 判據的理論，渦心處應該具有最大的Q 值，渦邊界處的Q 值較小。而對于繞圓盤的空化流動，存在若干個Q 極大值，可以發現Q 的極大值區域并沒有與渦心重疊，且Q 最大值出現在圓盤圓錐段后部，這說明Q 判據不能準確捕捉渦心的位置，只能大致表現存在渦旋結構的區域;同時也說明對于繞圓盤空化流動中的大尺度渦旋結構，Q 判據的效果不好，沒有準確識別出渦的邊界。文獻［25］中指出，對于大尺度渦的識別，基于歐拉體系的判據準則存在較大的困難，這與本文得到的結論一致。

觀察表2中的FTLE 分布可以發現，繞圓盤的空化流場中出現了明顯的LCS，即圖中紅色的部分。LCS 分為兩種:1)在圓盤前端出現了兩條平直的LCS;2)在圓盤后部，支桿的兩側出現了環狀的LCS.需要注意的是，求解FTLE 場時會過度預測壁面的FTLE 值，文獻［21］中也提到了這一點，因此只關注遠離壁面區域的LCS.

1)對于圓盤前端的平直LCS，觀察不同時刻的FTLE 分布可以發現，隨著時間推移，兩條LCS 始終保持平行，在圓盤前端劃分出一個窄流道。觀察流線圖可以發現，在LCS 形成的流道內，流線延伸至圓盤端面后貼壁向兩側發展，而流道外的流線直接繞過圓盤，延伸至圓盤后部流場區域。這說明流道內的質點流動受到圓盤阻礙，而流道外的質點能暢通地往下游流動，二者具有不同的運動特性，因而被LCS 分割開來。

2)對于圓盤后部的環狀LCS，圖2給出了其局部放大圖。Haller 等［16-17］通過嚴格的數學推導證明，流體微團將在渦旋邊界面上出現最大拉伸率或壓縮率，因而能被LCS 捕捉到。結合流線圖可以發現環狀LCS 很好地捕捉了圓盤后端空化區域的兩個渦旋結構，LCS 與渦旋的邊界吻合較好;渦心處于環狀LCS 的包裹之中，LCS 的運動規律與渦心保持一致。由此說明LCS 能準確識別大尺度渦旋結構，且能捕捉渦心的位置。

表2 含氣量分布、Q 分布和FTLE 分布Tab.2 Gas volume fraction，Q-criterion and FTLE

值得一提的是，LCS 是FTLE 場中的極大值構成的“脊”結構，其對渦旋結構的識別不需要人為設定閾值，而基于歐拉體系的渦判據的識別效果很大程度上依賴于設定的閾值，如前文所用到的Q 分布。

結合表2和圖2，可以說明空化發展與渦旋結構運動之間的關系。空穴尺度增大時，其內部含氣量升高，兩相界面更為清晰，由表2的流線圖和Q分布圖可以發現，流線分布更為密集且Q 值增大，說明位于空穴閉合端的渦旋結構得到加強;由FTLE分布圖可以發現，空穴尺度增大時，LCS 尺度變小，主要分布在空穴閉合端，說明此時的渦旋結構范圍變小，分布更為集中。

圖2 環狀LCS 局部放大圖Fig.2 Details of circular LCS

LCS 將流場劃分為運動特性不同的若干區域，為了進一步分析不同流場區域的流體運動特性，分別在5 個區域放置了一簇追蹤粒子，用不同顏色表示。表3給出了0.93T(T 為周期)內典型粒子的運動軌跡即跡線，給出了運動0.465T 后粒子所處的位置，初始位置用矩形符號表示，末了位置用倒三角符號表示。由表3可以發現，對于1 和2 區域的粒子，其跡線繞過圓盤和渦旋區向下游流動，由粒子的始末位置可以看出，初始時刻呈團狀排列的粒子在運動過程中被拉伸成直線排列，這是由于粒子間的運動速度存在差異，受到圓盤的阻礙作用，沿中軸線向兩側，粒子的運動速度逐漸增大。對于3 區域的粒子，其跡線到達圓盤后終止，由始末位置分布同樣可以看出，初始時刻呈團狀分布的粒子在運動0.465T后在圓盤前端面均勻分布，表示該區域的粒子到達圓盤端面后貼壁向兩側運動。對于4 和5 區域的粒子，其運動的軌跡較為復雜，粒子的運動受渦旋結構的大小以及渦心與粒子的相對位置影響，通過觀察不同時刻的FTLE 及LCS 分布，可以得出以下規律:1)由表3的跡線圖可知，環狀LCS 中的粒子呈螺旋狀向外運動，運動到LCS 邊界后，粒子隨主流向下游流動;部分粒子會受到反向射流的影響，運動到圓盤端面，如t0+4/6T 時刻4 區域的粒子。2)由表3的粒子始末位置圖可知，分布于大尺度LCS 中心的追蹤粒子，其運動規律較為一致，在運動0.465T 后，其相對位置變化不大，仍保持團狀分布，如t0時刻4區域和t0+4/6T 時刻5 區域的粒子，這是因為大尺度LCS 內部的FTLE 值較小，相鄰粒子的分離程度較小;對于小尺度LCS 中的粒子以及分布于環狀LCS 邊界上的粒子，其分離程度較大，如t0+2/6T時刻4、5 區域的粒子。

表3 追蹤粒子跡線圖及始末位置分布Tab.3 The tracks and position changes of particles

3.3 反向射流分析

圖3給出圓盤后端面FTLE 分布以及速度矢量的局部放大圖，由圖3可以發現，圓盤后端面附近的位置，近壁面存在低FTLE 值的區域，見圓圈標示的位置。t0+2/6T 時刻，該區域出現在支桿下側，t0+5/6T 時刻出現在上側，由速度矢量圖可以發現，該區域的反向射流強度較大。圖4給出了4 條監測曲線上X 方向的速度分布，監測曲線位置見圖3(a)所示。由圖4可以發現，遠離壁面的區域，速度大小與來流一致;在氣液交界面具有很大的速度梯度，速度迅速減小;近壁面附近，由于反向射流的存在，速度減小為負值。t0+2/6T 時刻，4 個位置的負向速度均為支桿上側較小，支桿下側較大，說明下側低FTLE 區域反向射流強度較大;t0+5/6T 時刻，速度分布規律相同，與t0+2/6T 時刻呈上下對稱狀分布。對比4 個位置的速度分布圖發現，沿圓盤端面往后，反向射流的強度逐漸增大，這是因為圓盤限制了反向射流的發展。

圖3 圓盤后端面FTLE 分布與速度矢量圖Fig.3 FTLE and velocity vector distribution on the back face of disc

圖5給出了相應時刻的空穴形態及含氣量分布圖。由圖5可以發現，在低FTLE 區域的一側，透明氣相區范圍較小，且空穴收縮至最小，空泡出現斷裂和脫落，這是由于該時刻反向射流強度較大，沿支桿推進的距離較大造成的。支桿四周的反向射流的強度交替增大和減弱，呈現周期性變化，故引起了環狀空穴周期性的膨脹收縮，造成了空穴的交替斷裂和空泡脫落。

圖4 X 方向速度分布Fig.4 Velocity distribution in X direction

反向射流與主流共同構成了圓盤下游的大尺度渦旋結構，因此反向射流強度的變化必然會引起渦旋結構的變化和運動。結合圖2、表3可以發現，反向射流強度較大的一側，環狀LCS 結構較大，且渦心位置距圓盤較近。由于主流強度與來流速度有關，在流動中較為穩定，因此反向射流強度的增大會引起大尺度渦旋結構的增大，從而引起渦心向上游運動。3.2 節中提到渦心的運動與空穴的膨脹收縮正是由于二者的運動均由反向射流強度變化引起。

圖5 空穴形態及含汽量分布圖Fig.5 Cavity forms and gas volume fraction

圖6 截面流線分布圖Fig.6 Streamlines on the axial cross section

為了分析LCS 與壁面分離的現象以及反向射流強度變化的原因，在距離圓盤后端面1 mm 處取了一個截面，如圖5(a)所示，圖6給出了截面上的流線圖。由圖6可以發現，支桿周圍的流動存在明顯的三維效應，流線存在同一個起始點，且起始點位置并不在支桿壁面上，而是與支桿存在一定距離，如圖中圓圈所示。起始點引出的流線分為不對稱的兩股，從支桿的兩側形成三維繞流，并在起始點逆時針方向90°的位置匯聚成一股。隨著時間發展，起始點的位置逆時針繞支桿轉動，這與空泡繞支桿螺旋狀發展規律一致。觀察t0+2/6T 和t0+5/6T 時刻的流線圖，可以發現起始點正好處于支桿四周低FTLE 區域邊界的位置。結合圖3和圖4可以發現，流線起始點實際上是反向射流最集中、強度最大的位置，反向射流沿軸向發展到圓盤端面附近后，受到圓盤的阻礙作用，從支桿兩側繞過，形成三維繞流。近壁面低FTLE 區域位于起始點內側，與遠離壁面的區域具有不同的流動狀態，從流線圖可以發現，起始點靠近壁面的一側形成繞支桿的三維流動，遠離壁面的一側流線向流場中發散。

4 結論

本文采用非均相流模型計算了繞圓盤的非定常空化流動，通過求解FTLE 場識別了流場中的LCS，并結合傳統的歐拉分析方法對繞圓盤空化流場渦旋結構的發展過程及流動機理進行了分析，得到了以下結論:

1)非均相流模型能準確預測繞圓盤空化流動的非定常特性，表現出云狀空泡團繞支桿的螺旋狀發展以及支桿四周空穴的交替膨脹收縮和空泡脫落，空穴形態呈周期性變化，與實驗結果一致。

2)LCS 將流場劃分為運動特性不同的區域:①平直段LCS 在圓盤前端劃分出窄流道，流道外的流體質點繞過圓盤向下游流動，流道內的流體質點流動受到圓盤阻礙，貼壁向兩側運動;②環狀LCS 準確捕捉了圓盤后部大尺度渦旋結構的邊界及渦心的位置，且能捕捉渦的運動規律，對渦結構的識別效果明顯好于Q 判據;LCS 內部的流體質點呈螺旋狀往外運動，達到LCS 邊界后隨主流往下游流動，環狀LCS 結構越大，其內部質點的運動規律越接近。

3)隨著空穴尺度的增大，渦旋結構得到加強，且渦旋結構的尺度減小，分布范圍更加集中。

4)反向射流強度的周期性變化是造成空穴和渦旋結構周期性變化和運動的原因;圓盤后端面出現低FTLE 值區域，這是由于反向射流強度增大造成的;支桿附近繞流存在明顯的三維效應，其軸截面上存在流線起始點，對應低FTLE 區域的邊界，起始點兩側流動狀態不同。

References)

［1］ 黃海龍，王聰，黃文虎 等.變攻角圓盤空化器生成自然超空泡數值模擬［J］.北京理工大學學報，2008，28(2):100 -103.HUANG Hai-long，WANG Cong，HUANG Wen-hu，et al.Numerical simulation on natural supercavity based on disk-cavitator with variable attack angle［J］.Transactions of Beijing Institute of Technology，2008，28(2):100 -103.(in Chinese)

［2］ 隗喜斌，王聰，榮吉利 等.錐體空化器非定常超空泡形態分析［J］.兵工學報，2007，28(7):863 -866.WEI Xi-bin，WANG Cong，RONG Ji-li，et al.Unsteadysupercavitating flow on cone cavitator［J］.Acta Armamentarii，2007，28(7):863 -866.(in Chinese)

［3］ 孟慶昌，張志宏，劉巨斌 等.壓縮性對亞聲速圓盤空化器超空泡流動的影響［J］.華中科技大學學報:自然科學版，2011，39(3):51 -54.MENG Qing-chang，ZHANG Zhi-hong，LIU Ju-bin，et al.Effect of the compressibility over supercavitating flow on subsonic disk cavitators［J］.Journal of Huazhong University of Science and Technology:Natural Science Edition，2011，39(3):51 - 54.(in Chinese)

［4］ 孟慶昌，張志宏，劉巨斌 等.超聲速圓盤空化器超空泡流動數值計算方法［J］.上海交通大學學報，2011，45(10):1435-1439.MENG Qing-chang，ZHANG Zhi-hong，LIU Ju-bin，et al.Numerical method for supercavitating flow over disk cavitator of underwater supersonic projectile［J］.Journal of Shanghai Jiaotong University，2011，45(10):1435 -1439.(in Chinese)

［5］ 袁緒龍，張宇文，楊武剛.高速超空化航行體典型空化器多相流CFD 分析［J］.彈箭與制導學報，2005，25(1):53 -59.YUAN Xu-long，ZHANG Yu-wen，YANG Wu-gang.Multi-phase CFD analysis of typical cavitator for high-speed supercavitating underwater vehicle［J］.Journal of Projectiles，Rockets.Missiles and Guidance，2005，25(1):53 -59.(in Chinese)

［6］ 時素果，王國玉，權曉波，等.當地均相介質模型在通氣超空化流動計算中的應用［J］.兵工學報，2011，32(2):147 -154.SHI Su-guo，WANG Guo-yu，QUAN Xiao-bo，et al.The application of a local homogenous medium model in the ventilated-supercavitation flow computations［J］.Acta Armamentarii，2011，32(2):147 -154.(in Chinese)

［7］ Johansen S T，Wu J，Shyy W.Filter-based unsteady RANS computations［J］.International Journal of Heat and Fluid Flow，2004，25(1):10 -21.

［8］ 黃彪，王國玉，袁海濤，等.基于密度分域的混合湍流模型的應用與評價［J］.北京理工大學學報，2010，30(10):1170 -1174.HUANG Biao，WANG Guo-yu，YUAN Hai-tao，et al.Assessment of a density modify based turbulence closure for unsteady cavitating flows［J］.Transactions of Beijing Institute of Technology，2010，30(10):1170 -1174.(in Chinese)

［9］ Boisson N，Malin M R.Numerical prediction of two-phase flow in bubble columns［J］.International Journal for Numerical Methods in Fluids，1996，23(12):1289 -1310.

［10］ Chahed J，Roig V，Masbernat L.Eulerian-Eulerian two-fluid model for turbulent gas-liquid bubbly flows［J］.InternationalJournal of Multiphase Flow，2003，29(1):23 -49.

［11］ Arakeri V H，Acosta A J.Viscous effects in the inception of cavitation on axisymmetric bodies［J］.Journal of Fluids Engineering，1973，95(4):519 -527.

［12］ Belahadji B，Franc J P，Michel J M.Cavitation in the rotational structures of a turbulent wake［J］.Journal of Fluid Mechanics，1995，287(1):383 -403.

［13］ Hunt J C R，Wray A A，Moin P.Eddies，streams，and convergence zones in turbulent flows［C］//Studying Turbulence Using Numerical Simulation Databases.Washington:NASA，1988:193 -208.

［14］ Chong M S，Perry A E，Cantwell B J.A general classification of three-dimensional flow fields［C］// Topological Fluid Mechanics.Cambridge:IUTAM，1990:408 -420.

［15］ Jeong J，Hussain F.On the identification of a vortex［J］.Journal of Fluid Mechanics，1995，285(69):69 -94.

［16］ Haller G，Yuan G.Lagrangian coherent structures and mixing in two-dimensional turbulence［J］.Physica D:Nonlinear Phenomena，2000，147(3):352 -370.

［17］ Haller G.Distinguished material surfaces and coherent structures in three-dimensional fluid flows［J］.Physica D:Nonlinear Phenomena，2001，149(4):248 -277.

［18］ Shadden S C，Lekien F，Marsden J E.Definition and properties of Lagrangian coherent structures from finite-time Lyapunov exponents in two-dimensional aperiodic flows［J］.Physica D:Nonlinear Phenomena，2005，212(3):271 -304.

［19］ Green M A，Rowley C W，Haller G.Detection of Lagrangian coherent structures in three-dimensional turbulence［J］.Journal of Fluid Mechanics，2007，572(1):111 -120.

［20］ O’Farrell C，Dabiri J O.A lagrangian approach to identifying vortex pinch-off［J］.An Interdisciplinary Journal of Nonlinear Science，2010，20(1):1 -9.

［21］ Tang J N，Tseng C C，Wang N F.Lagrangian-based investigation of multiphase flows by finite-time Lyapunov exponents［J］.Acta Mechanica Sinica，2012，28(3):612 -624.

［22］ Launder B E，Spalding D B.The numerical computation of turbulent flows［J］.Computer Methods in Applied Mechanics and Engineering，1974，3(2):269 -289.

［23］ Kubota A，Kato H，Yamaguchi H，et al.Unsteady structure measurement of cloud cavitation on a foil section using conditional sampling technique［J］.Journal of Fluids Engineering，1989，111(2):204 -210.

［24］ 張敏弟，郭善剛，邵峰，等.繞空化器自然空化流動的研究［J］.工程熱物理學報，2010，31(8):1319 -1323.ZHANG Min-di，GUO Shan-gang，SHAO Feng，et al.Study on natural cavitation around a cavitator［J］.Journal of Engineering Thermophysics，2010，31(8):1319 -1323.(in Chinese)

［25］ 胡子俊，張楠，姚惠之，等.渦判據在孔腔渦旋流動拓撲結構分析中的應用［J］.船舶力學，2012，16(8):839 -846.HU Zi-jun，ZHANG Nan，YAO Hui-zhi，et al.Vortex identification in the analysis on the topology structure of vortical flow in cavity［J］.Journal of Ship Mechanics，2012，16(8):839 -846.(in Chinese)