自由匯流旋渦多相耦合輸運(yùn)演變機(jī)理*

李霖 陸斌 許煒鑫 顧則恒 楊遠(yuǎn)山 譚大鵬2)?

1)(浙江工業(yè)大學(xué)機(jī)械工程學(xué)院,特種裝備制造與先進(jìn)加工技術(shù)教育部/浙江省重點(diǎn)實(shí)驗(yàn)室,杭州 310014)

2)(浙江大學(xué)機(jī)械工程學(xué)院,流體動(dòng)力與機(jī)電系統(tǒng)國家重點(diǎn)實(shí)驗(yàn)室,杭州 310027)

含自由液面的匯流旋渦抽吸演變中存在多相耦合、物質(zhì)傳輸、能量劇烈交換等物理過程,其中所涉及的多相流體耦合輸運(yùn)機(jī)理是具有高度非線性特征的復(fù)雜動(dòng)力學(xué)問題,多相黏滯耦合輸運(yùn)動(dòng)力學(xué)建模與數(shù)值求解具有較高難度.針對(duì)上述問題,提出一種含自由液面的匯流旋渦多相耦合輸運(yùn)建模與求解方法.基于水平集-流體體積耦合(CLSVOF)計(jì)算方法,結(jié)合連續(xù)表面張力模型和可實(shí)現(xiàn)(k-ε)湍流模型,建立含自由液面的匯流旋渦多相耦合輸運(yùn)動(dòng)力學(xué)模型;利用一種有效的體積修正方案來計(jì)算高速旋轉(zhuǎn)多相流,保證流場質(zhì)量守恒和無散度的速度場;結(jié)合相間耦合求解策略對(duì)多相流體分布與多相界面進(jìn)行精確追蹤.基于旋流場多特征物理變量,得到多相耦合界面動(dòng)態(tài)演變與跨尺度渦團(tuán)輸運(yùn)規(guī)律,揭示了多相耦合輸運(yùn)過程與壓力脈動(dòng)特性之間的相互作用機(jī)理.研究結(jié)果表明:多相耦合輸運(yùn)過程是流體介質(zhì)過渡的關(guān)鍵狀態(tài),旋渦微團(tuán)受到不同時(shí)空擾動(dòng)模式在界面處形成層層螺紋波形;旋渦多相耦合輸運(yùn)過程隨著水口尺度增大而增強(qiáng),且耦合能量激波引起非線性壓力脈動(dòng)現(xiàn)象.研究結(jié)果可為旋渦輸運(yùn)機(jī)理、渦團(tuán)跨尺度求解、流型追蹤等方面的研究提供有益借鑒.

1 引言

含自由液面的匯流旋渦是自由液面處氣-液相互作用所形成的一種動(dòng)態(tài)耦合流動(dòng)狀態(tài),且形成過程中伴隨相間傳質(zhì)、界面破碎、對(duì)流傳熱、能量交換等復(fù)雜物理現(xiàn)象.在旋渦輸運(yùn)過程中,強(qiáng)烈的渦芯抽吸力會(huì)卷吸表層流體介質(zhì)和固體顆粒物,是一個(gè)復(fù)雜的多相、多場耦合動(dòng)力學(xué)問題,其較強(qiáng)的氣-液抽吸耦合作用會(huì)引起非線性壓力脈動(dòng)現(xiàn)象,給工業(yè)生產(chǎn)造成不利影響[1-3].因此,研究含自由液面的匯流旋渦多相耦合輸運(yùn)過程,探索其臨界過渡狀態(tài)的多相黏滯耦合輸運(yùn)規(guī)律,揭示跨尺度湍流渦團(tuán)輸運(yùn)過程與非線性壓力脈動(dòng)特性之間的相互作用機(jī)理,具有重要的科學(xué)研究價(jià)值和廣闊的工程應(yīng)用前景.

自由匯流旋渦雖然是一種常見的自然現(xiàn)象,但卻是一個(gè)復(fù)雜的湍流力學(xué)問題,其規(guī)模參數(shù)與初始擾動(dòng)條件、約束物理空間構(gòu)造、表面粗糙度等因素有關(guān)[4].目前尚未有成熟的理論模型對(duì)其多相黏滯耦合輸運(yùn)規(guī)律進(jìn)行精確的定量分析,只能在某種理想假設(shè)下或結(jié)合實(shí)驗(yàn)研究,得到其有限的局部特征[5-7].由于旋渦具有三維非定常、紊流非線性、時(shí)空多尺度耦合特性,多相耦合跨尺度輸運(yùn)過程與壓力脈動(dòng)特性之間的相互作用機(jī)制尚未揭示.因此,研究匯流旋渦的多相耦合輸運(yùn)機(jī)理,對(duì)界面演變、渦團(tuán)結(jié)構(gòu)、多相耦合流型等多維動(dòng)力學(xué)特性進(jìn)行精確測定,分析旋流表觀物理特征與壓力脈動(dòng)演變之間的內(nèi)在聯(lián)系具有重要意義.

針對(duì)上述關(guān)鍵問題,國內(nèi)外學(xué)者進(jìn)行了大量的研究工作.Tahershamsi等[8]通過實(shí)驗(yàn)研究了旋渦的流動(dòng)過程,發(fā)現(xiàn)在較低雷諾數(shù)時(shí)氣核是一個(gè)空氣管.Tan等[9]利用水平集方法對(duì)兩相旋渦的自由界面進(jìn)行了追蹤,但臨界過渡狀態(tài)的界面演變特征不明顯.Morales等[10]采用剪切應(yīng)力輸運(yùn)湍流模型對(duì)旋渦的傳熱過程進(jìn)行了研究,發(fā)現(xiàn)溫度梯度可以提供足夠的浮力.Yang等[11]利用流體體積(volume of fluid,VOF)模型得到了旋渦流場,分析了旋渦的流動(dòng)結(jié)構(gòu)和動(dòng)態(tài)演化過程.?kerlavaj等[12]提出了一種帶有曲率修正的尺度自適應(yīng)模擬模型對(duì)單相旋渦的形成過程進(jìn)行研究,驗(yàn)證了該模型的可行性.Ann等[13]采用簡化的儲(chǔ)層模型對(duì)渦輪機(jī)運(yùn)行過程中產(chǎn)生的匯流旋渦進(jìn)行數(shù)值模擬,發(fā)現(xiàn)由于強(qiáng)渦量和水體中空氣含量的存在,導(dǎo)致管道流動(dòng)的不均勻性,所預(yù)測的旋渦與各種分析模型吻合較好.Li等[14]采用計(jì)算流體-離散單元耦合方法研究了氣液固三相匯流旋渦的輸運(yùn)效應(yīng),得到了顆粒流場的流型演變規(guī)律.

綜上可以推斷,當(dāng)前對(duì)自由匯流旋渦的研究主要集中在兩相動(dòng)態(tài)建模、介入擾動(dòng)、界面演變和流型追蹤等方面.針對(duì)多相旋渦演化過程研究,多層界面耦合求解、多相耦合輸運(yùn)規(guī)律、非線性壓力脈動(dòng)特性等問題尚不明確.由于旋渦抽吸輸運(yùn)過程的復(fù)雜性,相間界面的黏滯耦合輸運(yùn)與跨尺度渦團(tuán)流型追蹤具有高度非線性特性,這無疑增加了旋渦多相耦合數(shù)值計(jì)算的難度.因此,利用一種含自由液面的匯流旋渦多相耦合輸運(yùn)建模與求解方法,得到多相流黏滯耦合過程的相間輸運(yùn)規(guī)律,揭示跨尺度湍流渦團(tuán)輸運(yùn)過程與非線性壓力脈動(dòng)特性之間的相互作用機(jī)理是非常有必要的.

本文首先建立了基于水平集-流體體積耦合(coupled level set and volume-of-fluid,CLSVOF)耦合的旋渦數(shù)學(xué)模型,并結(jié)合體積修正方法和相間耦合求解方法.基于上述模型,建立了含自由液面的匯流旋渦多相耦合輸運(yùn)動(dòng)力學(xué)模型,討論了流量與旋渦多相耦合輸運(yùn)演變過程的關(guān)聯(lián)性,揭示了旋渦多相黏滯耦合輸運(yùn)規(guī)律.最后,進(jìn)一步討論了不同流量條件下的旋渦多相耦合輸運(yùn)過程與非線性壓力脈動(dòng)特性之間的相互作用機(jī)理.

2 匯流旋渦數(shù)學(xué)模型

如前所述,含自由液面的匯流旋渦是一個(gè)復(fù)雜的多相耦合流動(dòng)現(xiàn)象,具有高度非線性力學(xué)特性.常用的多相流體模型有兩種:VOF 模型[15]和CLSVOF 模型[16].后者結(jié)合了VOF 模型和水平集方法的優(yōu)點(diǎn),通過耦合計(jì)算來跟蹤自由界面的運(yùn)動(dòng),在多相流模擬中得到了廣泛的應(yīng)用.因此,基于CLSVOF 耦合方法可以更準(zhǔn)確地計(jì)算法向量和曲率,同時(shí)保持質(zhì)量守恒,更適合模擬具有自由界面的多相匯流旋渦.

2.1 基于CLSVOF 耦合的旋渦數(shù)學(xué)模型

本文研究有限物理空間中具有自由表面的非定常等溫旋轉(zhuǎn)多相流,其中流體域包含水、油、空氣的流體相以及氣-液界面和液-液界面.液體和氣體被視為不可壓縮的牛頓流體,多相流體的連續(xù)方程和動(dòng)量方程為

式中,u為流體速度,ρ和μ分別為流體相的密度和黏度,p為流體壓力,g為重力加速度,F為單位體積的表面張力.圖1 表示流體相的分布示意圖,由兩種不可壓縮牛頓流體和一個(gè)界面組成,故動(dòng)量方程應(yīng)包括界面處的表面張力.本節(jié)考慮連續(xù)表面張力模型[17,18]耦合計(jì)算表面張力,則修正的動(dòng)量方程為式中,σ為表面張力系數(shù)且假定為常數(shù),該函數(shù)在界面上為1,在其他區(qū)域?yàn)?;n為界面的單位法向量,κ為界面的平均曲率,δs為表面函數(shù).

圖1 含自由界面的多相流體示意圖Fig.1.Schematic diagram of multiphase fluid with a free interface.

為了實(shí)時(shí)追蹤旋渦多相耦合的動(dòng)態(tài)演變過程,本文利用了一種考慮體積修正的CLSVOF 耦合方法對(duì)多相流體分布與多相界面進(jìn)行精確追蹤.該耦合方法是通過VOF 函數(shù)追蹤對(duì)流界面,利用光滑的水平集(level set,LS)函數(shù)計(jì)算界面的法曲率,然后采用光滑的Heaviside 函數(shù)更新界面的物理性質(zhì).

在計(jì)算區(qū)域中,界面跟蹤通過控制體積中定義的體積分?jǐn)?shù)α來實(shí)現(xiàn).在位置矢量r和時(shí)間t條件下,流體相對(duì)應(yīng)的體積分?jǐn)?shù)α(r,t)表達(dá)式為[19]

在流體域中,相鄰層流體(如水-油和油-空氣)被定義為下層流體和上層流體.在水平集方法中,LS 函數(shù)φ作為一個(gè)平滑的函數(shù)表示到界面的符號(hào)距離.函數(shù)φ(r,t)可以表示為

這里d=d(r)表示在t時(shí)刻位置r到界面的最短距離.該值在界面處為零,在上層流體中為正值,在下層流體中為負(fù)值,如圖2 所示.根據(jù)(5)式,LS函數(shù)在每一個(gè)時(shí)間步都需要重新初始化.體積分?jǐn)?shù)α和LS 函數(shù)φ的平流方程可以表示為

圖2 界面附近的水平集函數(shù)等值線Fig.2.Level set function contours near the interface.

經(jīng)過幾個(gè)時(shí)間步后,LS 函數(shù)由于數(shù)值擴(kuò)散不再是距離函數(shù)|?φ|≠1.此時(shí)需要執(zhí)行界面重新初始化程序?α/?τS(φ0)(1-|?φ|)重新設(shè)置LS 函數(shù)的距離.利用LS 函數(shù)獲得界面法向量n,得到曲率和光滑的物理性質(zhì),然后計(jì)算界面的平均曲率.

為了評(píng)估含界面的控制單元的VOF 通量,采用分段線性界面重構(gòu)(piecewise linear interface construction,PLIC)方案實(shí)現(xiàn)復(fù)雜剪切流、旋渦流場自由液面的高精度重構(gòu)[20,21].因此,本節(jié)在CLSVOF 方法基礎(chǔ)上,通過PLIC 界面追蹤法對(duì)旋渦界面演化進(jìn)行捕捉,重構(gòu)Euler 網(wǎng)格下的自由界面形態(tài).針對(duì)平滑界面附近不連續(xù)的物理特性,利用平滑Heaviside 函數(shù)H(φ)確保物理性質(zhì)的連續(xù)變化[22]:

式中,Δ為網(wǎng)格的大小決定數(shù)值界面的厚度,且界面不會(huì)覆蓋超過兩個(gè)單元.基于平滑函數(shù)H(φ),流體性質(zhì)(如密度和黏度)將被分配到計(jì)算域中的每個(gè)控制體積.則每個(gè)控制體的流體相密度和黏度可表示為[23]

式中,下標(biāo)b 和a 分別表示上層流體和下層流體,則動(dòng)量方程可修改為

為了補(bǔ)償質(zhì)量變化并保持無散度的速度場,結(jié)合一種有效的體積修正方案來模擬高速旋轉(zhuǎn)多相流.將(4)式和(5)式重寫為變換后的定義域ξi=(ξ,η,ζ)如下:

式中,Ui=(U,V,W)為變換域內(nèi)的逆變速度分量,J為物理域到變換域變換的雅可比矩陣.界面重構(gòu)的關(guān)鍵是確定界面段在每個(gè)單元中的方向.在本研究中,每個(gè)界面段可以在變換域ξi=(ξ,η,ζ)內(nèi)重構(gòu):

式中,φ為從原點(diǎn)到界面最近距離的參數(shù),由(7)式中的LS 函數(shù)可以得到法向量n=(n1,n2,n3).

基于上述變換方法,該界面可適用于一般曲線網(wǎng)格,并通過體積修正方案保持質(zhì)量守恒.速度場可以等量修改為(0.5U,0.5V,0),(0,0.5V,0.5W)和(0.5U,0,0.5W).然后,采用混合分裂三維歐拉隱式-拉格朗日顯式方法將VOF 函數(shù)完全更新到新的時(shí)間級(jí)別.最后,基于PLIC 重構(gòu)界面,對(duì)LS函數(shù)進(jìn)行質(zhì)量守恒調(diào)整,采用三線性插值方法將LS 函數(shù)從變換域轉(zhuǎn)換到物理域[24].

2.2 湍流模型

匯流旋渦是一種復(fù)雜的氣-液-液多相耦合流動(dòng),由于初始擾動(dòng)的影響,旋渦流場在初始條件下可達(dá)到完全湍流狀態(tài).可實(shí)現(xiàn)k-ε模型作為一種經(jīng)典的湍流模型,在分析大時(shí)均應(yīng)變情況時(shí),滿足雷諾應(yīng)力的約束條件,保證流動(dòng)更符合湍流物理定律.同時(shí)該模型在湍流黏性計(jì)算公式引入曲率和旋度相關(guān)項(xiàng),依據(jù)渦量脈動(dòng)方程對(duì)耗散率方程進(jìn)行修正,保證模型在旋轉(zhuǎn)流場、邊界層及分離計(jì)算時(shí)更符合真實(shí)物理情況[25,26],其輸運(yùn)方程為

式中,ρ為流體的密度;u為流體的速度;σk,σε分別為k和ε對(duì)應(yīng)的普朗特?cái)?shù);E1為時(shí)均應(yīng)變率張量的模量;ν為流體運(yùn)動(dòng)黏度;μ為流體的動(dòng)力黏度;Gk,Gb和YM與標(biāo)準(zhǔn)k-ε模型的對(duì)應(yīng)項(xiàng)相同;μt為湍流黏度,并可通過公式μtρCμk2/ε求解.采用的其他經(jīng)驗(yàn)參數(shù)如下:C1=1.44,C2=1.92,σk=1.0,σε=1.3,A0=4.0.

在可實(shí)現(xiàn)k-ε模型中,湍流黏性系數(shù)Cμ是流體動(dòng)力黏度μt的一個(gè)關(guān)鍵變量,可看作是含旋轉(zhuǎn)速度、時(shí)均應(yīng)變和湍流強(qiáng)度等參數(shù)的函數(shù).顯然.可實(shí)現(xiàn)k-ε模型在模擬射流、邊界層流和旋轉(zhuǎn)剪切流場特性方面具有更好的優(yōu)勢[27,28].該方法可以得到更精確的雷諾應(yīng)力,符合湍流的實(shí)際物理特征.

2.3 CLSVOF 相間耦合求解方法

基于上述數(shù)學(xué)模型,本節(jié)利用了CLSVOF 相間耦合求解算法,如圖3 所示.該方法計(jì)算旋渦多相耦合輸運(yùn)的主要過程包括相函數(shù)初始化、流動(dòng)控制方程和湍動(dòng)能-耗散方程求解、相函數(shù)對(duì)流輸運(yùn)方程求解、相界面重構(gòu)以及相函數(shù)再次重新初始化等.在CLSVOF 方法中,首先調(diào)用自定義初始化函數(shù),求解多相流體的控制方程和湍動(dòng)能-耗散方程,得到多相流體的壓力和速度.

圖3 CLSVOF 相間耦合求解流程圖Fig.3.Flow chart of interphase coupling solution with CLSVOF.

隨后,利用VOF 函數(shù)和LS 函數(shù)計(jì)算界面法向量,并通過PLIC 方案重構(gòu)界面.通過幾何程序重新距離LS 函數(shù),結(jié)合體積修正方法以實(shí)現(xiàn)質(zhì)量守恒.LS 方法和VOF 方法的耦合發(fā)生在界面重構(gòu)和LS 重新距離化過程中(如圖3 虛線框所示).但當(dāng)耦合計(jì)算不滿足收斂條件時(shí),需要對(duì)修正的動(dòng)量方程再次求解,得到更新后的流體壓力和速度數(shù)據(jù).上述求解過程不斷重復(fù),直到流場耦合求解過程滿足收斂條件.

3 多相匯流旋渦輸運(yùn)動(dòng)力學(xué)模型

3.1 物理對(duì)象模型

為了研究含自由液面的匯流旋渦多相耦合輸運(yùn)過程,建立有限物理空間的幾何尺寸模型,如圖4所示.在初始狀態(tài)下,多相流體(水、油和空氣)由于密度不同且互不相容,均勻分布在上部容器中.旋渦的形成演化過程發(fā)生在底部帶有排流管道的容器內(nèi)部,其中容器的幾何參數(shù)主要包括:容器高度H=0.55 m,容器的直徑D=0.5 m,管道長度L=0.15 m,管道直徑d=0.030 m.容器的頂部與外界大氣壓相連,參考?jí)毫?.01×105Pa,且流體重力方向與排流方向保持一致.

圖4 物理對(duì)象模型Fig.4.Physical objective model.

3.2 流體動(dòng)力學(xué)模型

基于上述幾何尺寸模型,建立了含自由液面的匯流旋渦輸運(yùn)動(dòng)力學(xué)模型,如圖5 所示.網(wǎng)格劃分是采取離散化處理有限空間的容器模型,其劃分尺寸直接關(guān)系到數(shù)值計(jì)算的時(shí)間和精度,具有高質(zhì)量合適尺度的網(wǎng)格劃分,也可提高數(shù)值精度和計(jì)算效率.鑒于六面體結(jié)構(gòu)化網(wǎng)格具有網(wǎng)格質(zhì)量高、數(shù)值擴(kuò)散小以及精度高等方面的優(yōu)勢,適用于多相匯流旋渦數(shù)值模型的網(wǎng)格劃分.

本文基于ICEM 網(wǎng)格劃分軟件進(jìn)行模型網(wǎng)格劃分,并對(duì)位于圓柱中心區(qū)域和排流管道的網(wǎng)格進(jìn)行加密處理,加密尺度為0.003 m,精確捕捉臨界過渡狀態(tài)的多相耦合輸運(yùn)動(dòng)力學(xué)特征.在保證網(wǎng)格質(zhì)量的前提下,采用較大尺度的網(wǎng)格(尺度為0.035 m)對(duì)容器模型的剩余部分進(jìn)行劃分以便提高計(jì)算收斂速度,網(wǎng)格總數(shù)為335460,且網(wǎng)格質(zhì)量在計(jì)算允許情況下調(diào)整至最高(網(wǎng)格質(zhì)量0.8 以上).從圖5 中的網(wǎng)格劃分可以看出,模型網(wǎng)格劃分較為均勻,排流管道附近的網(wǎng)格相對(duì)較密,采用較高質(zhì)量的網(wǎng)格保證數(shù)值計(jì)算的精度要求.

圖5 流體力學(xué)模型與邊界條件(a)流體域模型;(b)區(qū)域A 的局部視角;(c)區(qū)域B 的局部視角;(d)容器頂部視角.1-入口壁面;2-固定壁面;3-出口壁面Fig.5.Fluid dynamic mechanic model and boundary conditions:(a)Fluid-domain model;(b)local view of region A;(c)local view of region B;(d)top view of the container.1-inlet wall;2-fixed wall;3-out wall.

3.3 模型邊界條件和初始條件

基于上述數(shù)值模型,研究匯流旋渦的多相耦合輸運(yùn)過程,其中容器的幾何特征和邊界條件如表1所列.多相流體介質(zhì)的基本參數(shù)如表2 所列,為了建立水-油界面和油-空氣界面,通過區(qū)域分割方法將油相和水相分別從整個(gè)數(shù)值模型中分離,其余區(qū)域設(shè)為氣相.上述3 個(gè)流體區(qū)域相互聯(lián)系,可以進(jìn)行質(zhì)量、動(dòng)量和能量的相互傳輸,且初始切向速度和重力加速度作用于整個(gè)流場區(qū)域.相關(guān)研究結(jié)果[1-3]表明,當(dāng)初始切向速度達(dá)到臨界值時(shí),旋轉(zhuǎn)方向與初始切向方向一致,本文選取初始切向速度(ω0=1.0π rad/s)來研究旋渦多相耦合輸運(yùn)形成過程,流場在該條件下可以達(dá)到完全發(fā)展的湍流狀態(tài)[29-32].

表1 流體動(dòng)力學(xué)模型的邊界條件.Table 1.Boundary conditions of the fluid dynamic model.

表2 流體介質(zhì)的物理參數(shù).Table 2.Physical parameters of fluid mediums.

采用有限體積法離散控制方程保證嚴(yán)格的流體質(zhì)量守恒.由于旋渦模擬屬于非穩(wěn)態(tài)過程,相間耦合演化非常復(fù)雜,屬于典型的過渡流動(dòng)計(jì)算,旋渦界面形態(tài)以及速度和壓力場演變復(fù)雜,采用壓力的隱式算子分割(pressure implicit with split of operator,PISO)算法處理壓力、速度耦合以保證收斂效率[1,4].為了解決壓力離散插值問題,防止內(nèi)壓急劇變化和高旋流,采用壓力交錯(cuò)方式(pressure staggering option,PRESTO)[33].LS 函數(shù)屬于哈密頓-雅可比方程,空間上采用五階加權(quán)本質(zhì)無振蕩(weighted essentially non-oscillatory,WENO)格式離散,時(shí)間上采用三階總變差減小(total variation diminishing,TVD)龍格-庫塔格式離散[34,35].此外,采用二階迎風(fēng)格式離散動(dòng)量、湍動(dòng)能和耗散率來獲得精確解.

3.4 旋渦動(dòng)力學(xué)模型驗(yàn)證

為了驗(yàn)證所提方法的有效性,研究了二維駐波算例的數(shù)值精度,得到了解析解和VOF 方法的數(shù)值測試[36],如圖6 所示.測試案例與Deshpande等[36]的相似,界面的初始輪廓是一個(gè)正弦曲面形狀,計(jì)算的初始邊界條件如圖6(a)所示.重力水平設(shè)為零,沿著左右兩側(cè)的邊界是循環(huán)的,沿著頂部和底部的邊界是壁面.圖6(b)表示表面振蕩的衰減,且存在一個(gè)最大振蕩點(diǎn).計(jì)算得到的衰減高度與解析解吻合較好,兩次振蕩后計(jì)算得到的高度略大于解析值,VOF 和CLSVOF 方法的誤差分別為4.03%和1.01%.該耦合方法在測試算例中得到了更準(zhǔn)確的結(jié)果.

圖6 二維駐波算例驗(yàn)證(a)駐波;(b)表面高度的衰減.Fig.6.Validation of the two-dimensional standing wave:(a)Standing wave;(b)decay of the surface elevation.

網(wǎng)格數(shù)量對(duì)旋渦流場瞬態(tài)模擬結(jié)果的精度具有很大影響,有必要進(jìn)行網(wǎng)格無關(guān)性研究以滿足仿真結(jié)果的精度和可重復(fù)性[37-39].數(shù)值算例生成了5 種不同密度網(wǎng)格比較檢測點(diǎn)的總壓強(qiáng)變化曲線,如圖7(a)所示,得到了壓力數(shù)值計(jì)算結(jié)果.可以看出,當(dāng)網(wǎng)格密度較低(S1—S3)時(shí),總壓p的結(jié)果偏差較大,數(shù)值誤差分別為0.51%,0.34%和0.77%.當(dāng)網(wǎng)格密度達(dá)到一定值時(shí),曲線S4的壓力p與曲線S5的值呈現(xiàn)一致的均勻分布,誤差分別為0.26%和0.23%.因此,網(wǎng)格分辨率(S4和S5)已經(jīng)滿足了網(wǎng)格無關(guān)性的要求,保證了數(shù)值計(jì)算的準(zhǔn)確性和可重復(fù)性.

圖7 旋渦動(dòng)力學(xué)模型驗(yàn)證(a)網(wǎng)格無關(guān)性;(b)CLSVOF 與VOF 模型的液位高度驗(yàn)證Fig.7.Validation of the vortex dynamic model:(a)Mesh independence;(b)liquid level height validation of the CLSVOF and VOF models.

針對(duì)自由匯流旋渦的形成過程,與Park等[40]的實(shí)驗(yàn)和理論結(jié)果對(duì)比分析旋渦的液面高度變化.圖7(b)采用所提出的建模方法建立了氣-液兩相界面的數(shù)值模型,并保證相同的邊界條件.通過VOF方法和CLSVOF 方法得到旋渦的液面高度.可以看出,數(shù)值結(jié)果與理論值和實(shí)驗(yàn)結(jié)果吻合較好,且液面高度呈現(xiàn)連續(xù)下降趨勢.但隨著液體表面的下降和自由匯流旋渦的形成,VOF 方法與S3網(wǎng)格的數(shù)值結(jié)果略大于實(shí)驗(yàn)值,這與氣液耦合界面的高曲率變形有關(guān).因此,采用基于CLSVOF(S4)的旋渦動(dòng)力學(xué)建模方法,可以獲得更精確的自由界面變形特性.

4 數(shù)值計(jì)算結(jié)果和分析

4.1 旋渦多相耦合輸運(yùn)演變過程

旋渦輸運(yùn)過程中,由于抽吸力變化與界面黏性摩擦影響,多相耦合輸運(yùn)規(guī)律具有高度非線性特性.本節(jié)選取數(shù)值案例(ω0=1.0π rad/s,d=30 mm)研究旋渦多相耦合輸運(yùn)演變規(guī)律.流體體積分?jǐn)?shù)剖面可以反映相間界面的演化規(guī)律和流場過渡特征,如圖8 所示,其中紅色區(qū)域?yàn)樗?綠色區(qū)域?yàn)橛拖?藍(lán)色區(qū)域?yàn)闅庀?

在旋渦形成前期,流體的流動(dòng)模式表現(xiàn)為重力引起的軸向運(yùn)動(dòng)和初始擾動(dòng)誘導(dǎo)的切向運(yùn)動(dòng),如圖8(a)—(c)所示.初始切向速度增加了流動(dòng)雷諾數(shù),湍流強(qiáng)度增大,這不僅增強(qiáng)了流體的切向運(yùn)動(dòng),而且提高了從壁面附近向中心運(yùn)動(dòng)的收斂速率.在圖8(b)和圖8(c)中,旋渦中心的液面高度明顯低于兩端,油層底部中心形成凸點(diǎn).當(dāng)抽吸力達(dá)到一定強(qiáng)度時(shí),凸點(diǎn)逐漸向排流孔發(fā)展并形成細(xì)長油柱.隨著油滴和水從排流管流出,旋渦的尺度不斷增大,說明出口抽吸力促進(jìn)了流體向排流孔收斂,并能以能量的形式向上輸送.但由于水相和油相的黏性阻力不同,能量的衰減使得剩余能量不足以突破水相表面的能量壁壘.

隨著旋渦輸運(yùn)過程的動(dòng)態(tài)演變,流場中心有少量空氣被排流孔抽吸并懸浮在油相中,如圖8(d)所示.強(qiáng)烈的抽吸力克服了流體的黏滯阻力,打破了水相和油相混合表面中心的能量壁壘,此時(shí)為吸氣狀態(tài).在圖8(e)中,液面中心的最低位置與容器底部對(duì)齊,由于油相吸入大量氣泡,流體在管道的流動(dòng)不穩(wěn)定,具有高度非線性特征,流體成分為水、油、空氣的混合物,此時(shí)為抽氣狀態(tài).當(dāng)流場與出口大氣壓相連時(shí),水油兩相的流速基本為零,表明旋渦到達(dá)貫穿狀態(tài),如圖8(f)所示,可以看出吸氣-抽氣-貫穿過程是多相流體介質(zhì)過渡的關(guān)鍵狀態(tài).

針對(duì)旋渦的多相耦合輸運(yùn)過程,得到了旋渦流場的油相三維形貌圖,如圖9 所示,其中油相流體表面為空氣和油兩相自由液面.從圖9 可觀測到旋渦底部中心油層凸點(diǎn)和油柱形成的現(xiàn)象,與圖8(b)和圖8(c)所示的體積分?jǐn)?shù)剖面相對(duì)應(yīng),自由液面形成層層螺紋波形,如圖9(a)所示.當(dāng)流場受到初始擾動(dòng)速度的影響時(shí),旋渦表面的流體微團(tuán)在不同的時(shí)間和空間上存在特性的擾動(dòng)模式,導(dǎo)致流體微團(tuán)沿徑向方向具有不同的速度梯度.該現(xiàn)象在油層底部的油水界面處同樣存在,水油兩相間的速度存在差值,如圖9(b)所示.從圖9(c)可以看到,隨著旋渦多相耦合輸運(yùn)過程,排流管道存在明顯的氣泡抽吸現(xiàn)象,氣泡的尺度和數(shù)量具有高度非線性特性,管道內(nèi)演變?yōu)閺?fù)雜的水-油-氣多相耦合流動(dòng)模式.

圖8 含自由液面的匯流旋渦體積分?jǐn)?shù)云圖(a)t=1.00 s;(b)t=3.50 s;(c)t=6.00 s;(d)t=21.50 s;(e)t=23.00 s;(f)t=25.00 sFig.8.Volume fraction cloud chart of the sink vortex with a free-liquid surface(ω0=1.0π rad/s):(a)t=1.00 s;(b)t=3.50 s;(c)t=6.00 s;(d)t=21.50 s;(e)t=23.00 s;(f)t=25.00 s.

圖9 旋渦流場的油相三維形貌圖(a)凸點(diǎn)形成;(b)油相抽吸;(c)多相耦合輸運(yùn);(d)油相體積分?jǐn)?shù)Fig.9.Three-dimensional morphology of oil-phase in the vortex flow field:(a)Salient point formation;(b)oil-phase suction;(c)multiphase coupling transport;(d)volume fraction of the oil phase.

圖10 為旋渦輸運(yùn)過程的流線剖面圖,對(duì)其內(nèi)部流動(dòng)模式進(jìn)行分析.在圖10(a)—(c)中,流場中心受到排流口的抽吸力向管口流動(dòng),但在排流口左右兩側(cè),流線出現(xiàn)了Ekman 層上升流動(dòng)結(jié)構(gòu).油層內(nèi)流體沿著水油兩相交界面流動(dòng)時(shí),受到交界面束縛與凸點(diǎn)區(qū)域結(jié)構(gòu)影響,流動(dòng)模式發(fā)生輕微的改變,如圖10(c)所示.在圖10(d)中,隨著旋渦達(dá)到吸氣階段,在旋渦中心形成較大尺度的局部渦旋.根據(jù)伯努利方程可知,旋渦中心區(qū)域由于較高的抽吸旋轉(zhuǎn)速度具有較低的壓力梯度,導(dǎo)致氣相流動(dòng)模式演化為局部渦旋.此時(shí),容器底部壁面兩側(cè)流動(dòng)模式受到Ekman 邊界層的螺旋耦合作用演變?yōu)閷?duì)稱的局部渦旋,流動(dòng)模式受到水口抽吸出現(xiàn)水平輸運(yùn),表現(xiàn)為先上升后向下流動(dòng),其輸運(yùn)過程呈現(xiàn)較強(qiáng)的非線性流動(dòng)特征,在x=0.06 m 出現(xiàn)先上升流動(dòng)后被管口抽吸流動(dòng)的現(xiàn)象.在圖10(e)和圖10(f)中,隨著旋渦多相耦合輸運(yùn)過程的演化,中心氣相渦旋逐漸縮小,且隨著貫穿階段到達(dá)而消失.

圖10 旋渦多相耦合輸運(yùn)流線圖(a)t=1.00 s;(b)t=3.50 s;(c)t=6.00 s;(d)t=21.50 s;(e)t=23.00 s;(f)t=25.00 sFig.10.Transport streamline chart of the vortex multiphase coupling:(a)t=1.00 s;(b)t=3.50 s;(c)t=6.00 s;(d)t=21.50 s;(e)t=23.00 s;(f)t=25.00 s.

4.2 水口尺度對(duì)旋渦多相耦合輸運(yùn)的影響

為了研究流量對(duì)旋渦多相耦合輸運(yùn)過程的影響,通過定義無量綱水口尺度D*=d/D(d=0.003 m,D=0.5 m),研究不同水口尺度條件下的結(jié)果,如圖11 所示.流場徑向速度流線可以描述多相耦合過程的旋渦強(qiáng)度和尺度特征.從圖11 可以看出,水口尺度不同時(shí),速度流線分布較為均勻,徑向速度在壁面處最大,在旋渦中心處較小.隨著水口尺度的增大,旋渦徑向速度不斷增大,周向動(dòng)能和流體軸向擾動(dòng)具有較快的累積速率,流體表面在短時(shí)間內(nèi)不穩(wěn)定.說明了多相耦合輸運(yùn)過程的高速旋轉(zhuǎn)使湍流運(yùn)動(dòng)增強(qiáng),流體黏滯力的影響較小,慣性力占主導(dǎo)地位.在高速旋轉(zhuǎn)時(shí),沙漏區(qū)域使表面流體和雜質(zhì)容易流入排流孔,加快多相耦合輸運(yùn)過程.

圖11 不同水口尺度條件下的旋渦多相耦合輸運(yùn)流線圖(a)0.73D*;(b)0.87D*;(c)1.00D*;(d)1.13D*Fig.11.Transport streamline chart of the vortex multiphase coupling under different nozzle diameter:(a)0.73D*;(b)0.87D*;(c)1.00D*;(d)1.13D*.

圖12 描述了不同水口直徑條件的液位高度.在旋渦多相耦合輸運(yùn)過程中,隨著水口直徑的增大,液位高度在吸氣狀態(tài)的液高有明顯增加.隨著旋渦的動(dòng)態(tài)演變,抽氣狀態(tài)和貫穿狀態(tài)的液位高度差明顯減小,說明在初始擾動(dòng)速度一定的情況下,水口直徑越大,旋渦的強(qiáng)度和規(guī)模越大,對(duì)氣相的抽吸力顯著增強(qiáng),這使得在較高的液位下旋渦可以快速完成貫穿過程.顯然,在實(shí)際工業(yè)應(yīng)用中,要想減小旋渦卷吸表層液體的過程,降低排流量是一種有效的控制措施.

圖12 不同水口直徑的液位高度.Fig.12.Liquid-level height of different nozzle diameter.

為了進(jìn)一步探索不同水口尺度對(duì)旋渦多相耦合輸運(yùn)過程的影響,得到了截面z=0.18 m 時(shí)旋渦流場沿徑向坐標(biāo)分布的關(guān)鍵流動(dòng)特性,如圖13所示.圖13(a)表示切向速度沿半徑方向的演變規(guī)律.可以發(fā)現(xiàn),旋渦周向速度隨管徑增大而增大,且極大值出現(xiàn)在半徑x=0.032 m,呈現(xiàn)出較好的相似性,這與Ahn等[13]旋渦抽吸數(shù)值結(jié)果相吻合.水口尺度1.13D*的周向速度最大值為0.32 m/s 是水口尺度0.73D*周向速度最大值的1.9 倍.這說明流量的增加促進(jìn)了旋渦周向流動(dòng)的收斂速度,提高了旋渦多相耦合輸運(yùn)的規(guī)模和強(qiáng)度.

在圖13(b)中,湍動(dòng)能是反映流場湍流流動(dòng)的重要指標(biāo).在湍流狀態(tài)下,不同尺度和渦量的渦旋不斷破碎、合并,流體微團(tuán)作復(fù)雜的、無規(guī)則的、隨機(jī)的不定常運(yùn)動(dòng).可以看出,湍動(dòng)能在旋渦渦芯處具有最大值,且沿著徑向方向單調(diào)遞減.這說明在多相耦合輸運(yùn)過程中,渦芯處受到管口抽吸力作用以軸向流動(dòng)為主,管口處具有最大的流動(dòng)速度,此時(shí)湍動(dòng)能在x＞ 0.15 m 時(shí)趨于零值.上述現(xiàn)象說明旋渦多相耦合過程主要以管口的軸向輸運(yùn)為主,沿徑向方向的流場湍動(dòng)能相對(duì)較弱,遠(yuǎn)離管口區(qū)域的旋渦多相耦合輸運(yùn)能量較低.因此,在旋渦主動(dòng)控制中,應(yīng)盡可能降低旋渦管口附近區(qū)域受到強(qiáng)烈抽吸力作用而導(dǎo)致的軸向輸運(yùn)過程,從而降低旋渦多相耦合輸運(yùn)的強(qiáng)度.

圖13 旋渦流場的關(guān)鍵流動(dòng)特性(a)切向速度分布;(b)湍動(dòng)能分布;(c)渦量;(d)動(dòng)壓.Fig.13.Key flow characteristics of the vortex flow field:(a)Tangential velocity distribution;(b)turbulent energy distribution;(c)vorticity;(d)dynamic pressure.

圖13(c)描述了旋渦流場中渦量沿徑向坐標(biāo)的變化曲線.可以看出,渦量值在x<0.05 m 范圍內(nèi)迅速減小,呈現(xiàn)出與湍動(dòng)能曲線相似的單調(diào)遞減趨勢.但是,在渦芯處的渦量在排流量改變時(shí)強(qiáng)度值相差并不大,此時(shí)的渦量主要是管口附近軸向主流動(dòng)以及底部受到Ekman 層作用產(chǎn)生的先上升后向下流動(dòng)共同作用形成的.因此,旋渦在管口的中心上方具有最大的渦量值.圖13(d)表示旋渦流場中動(dòng)態(tài)壓力沿徑向左邊的演變趨勢,流場的動(dòng)壓與流速有關(guān).可以看出,在當(dāng)前的流動(dòng)狀態(tài)下,不同排流量下的動(dòng)態(tài)壓力曲線差別不大,且沿著徑向坐標(biāo)方向不斷增大,最后動(dòng)態(tài)壓力曲線在容器壁面處具有相近值,這與譚大鵬等[3]兩相旋渦抽吸過程的數(shù)值結(jié)果相吻合.

4.3 旋渦多相耦合輸運(yùn)的壓力脈動(dòng)特性

旋渦多相耦合輸運(yùn)過程中流體組分的不穩(wěn)定是壓力脈動(dòng)產(chǎn)生的主要原因.由于氣相密度較小且在多相耦合輸運(yùn)狀態(tài)不連續(xù),使得流場的壓力脈動(dòng)呈現(xiàn)高度非線性特征.因此,研究旋渦壓力脈動(dòng)特性與多相耦合流體組分過渡的內(nèi)在關(guān)聯(lián)具有重要的意義.選擇不同水口直徑條件下管口總壓分析旋渦多相耦合輸運(yùn)過程的壓力脈動(dòng)演變特征,如圖14所示.

從圖14 可以得到以下規(guī)律:1)管口總壓隨著排流過程不斷減小,壓力云圖的旋渦中心壓力也不斷減小,且在排流后期存在劇烈的壓力振蕩現(xiàn)象.當(dāng)流體成分為穩(wěn)定的水相時(shí),總壓相對(duì)穩(wěn)定,而當(dāng)流體演變?yōu)閺?fù)雜水-油-氣耦合狀態(tài)且流量不穩(wěn)定時(shí),壓力曲線呈現(xiàn)劇烈振蕩過程,且具有較強(qiáng)的隨機(jī)性.2)從不同流量角度來看,隨著流量的增大,水油兩相的體積分?jǐn)?shù)下降較快,總壓的變化速率更快.3)根據(jù)振蕩區(qū)間可以推斷,壓力振蕩的時(shí)間受流量控制,與流量成反比例關(guān)系.隨著流量增加,旋渦在排流后期具有較快的能量積累和釋放速率,導(dǎo)致旋渦壓力脈動(dòng)區(qū)間明顯減小,這主要是由于上層氣相可在極短時(shí)間內(nèi)貫穿管道.

圖14 排流過程的總壓變化曲線(a)總壓曲線;(b)總壓云圖Fig.14.Total pressure variation curves at the whole drain process:(a)Total pressure curve;(b)cloud diagram of the total pressure.

根據(jù)上述壓力脈動(dòng)現(xiàn)象可以發(fā)現(xiàn),在臨界貫穿狀態(tài)之前,總壓主要受到流體的靜壓的主導(dǎo).隨著排流過程的發(fā)展,流體靜壓迅速下降,以流體靜壓為主的總壓迅速下降.流體動(dòng)壓與流體速度有關(guān),流體流速隨旋渦的不斷發(fā)展而增大.因此,在臨界貫穿狀態(tài),流體速度達(dá)到最大值,流體動(dòng)壓占總壓的比例不斷增大.上述現(xiàn)象表明,在排流后期,流體介質(zhì)過渡轉(zhuǎn)變引起的非周期湍動(dòng)能擴(kuò)散導(dǎo)致總壓發(fā)生劇烈振蕩現(xiàn)象.當(dāng)流量較小時(shí),流體非周期湍動(dòng)能的擴(kuò)散持續(xù)時(shí)間較長,流體動(dòng)壓主導(dǎo)的總壓將產(chǎn)生持續(xù)振蕩.上述規(guī)律可為旋渦主動(dòng)抑制控制和實(shí)時(shí)檢測提供有益指導(dǎo).

5 結(jié)論

研究自由匯流旋渦的多相耦合輸運(yùn)機(jī)理具有重要的科研價(jià)值與工程應(yīng)用前景.本文提出了一種含自由液面的旋渦多相耦合輸運(yùn)建模與求解方法,研究不同流量條件下的旋渦多相黏滯耦合輸運(yùn)規(guī)律和非線性壓力脈動(dòng)機(jī)理.主要結(jié)論如下:

1)基于體積修正的CLSVOF 耦合方法和可實(shí)現(xiàn)湍動(dòng)能-耗散模型,建立了含自由液面的旋渦多相耦合輸運(yùn)動(dòng)力學(xué)模型,得到了多相黏滯耦合輸運(yùn)規(guī)律.多相耦合輸運(yùn)過程是流體介質(zhì)過渡的關(guān)鍵狀態(tài);容器底部的流動(dòng)模式受到Ekman 邊界層的螺旋耦合作用演變?yōu)榫植繙u旋,其輸運(yùn)模式呈現(xiàn)較強(qiáng)的非線性流動(dòng)特征.

2)自由液面受到擾動(dòng)速度影響形成了層層的螺紋狀流動(dòng)波形,表面流體微團(tuán)在不同的時(shí)空上具有特性的擾動(dòng)模式,使得流體微團(tuán)呈現(xiàn)不同的速度梯度;旋渦抽吸氣泡尺度和數(shù)量的演變導(dǎo)致管道內(nèi)部水-油-氣耦合流動(dòng)模式變得復(fù)雜.

3)水口尺度決定旋渦多相耦合輸運(yùn)過程.隨著水口尺度D*增大,旋渦多相耦合輸運(yùn)程度增強(qiáng),周向動(dòng)能和流體軸向擾動(dòng)具有較快的累積速率,周向速度的極大值出現(xiàn)在半徑x=0.032 m 處,提高了旋渦多相耦合輸運(yùn)的規(guī)模和強(qiáng)度.

4)多相耦合輸運(yùn)過程是壓力脈動(dòng)產(chǎn)生的主要原因.在排流后期,多相耦合輸運(yùn)的能量激波引起壓力振蕩現(xiàn)象,該振蕩時(shí)間在整個(gè)排流過程所占比例較小且受流量的控制.當(dāng)流量增大時(shí),旋渦輸運(yùn)能量聚合與耗散速率明顯提高,減小了壓力脈動(dòng)區(qū)間.上述規(guī)律可為旋渦多相耦合輸運(yùn)過程的主動(dòng)抑制控制提供有益指導(dǎo).