朱 嵩,劉國華,程偉平,黃躍飛

(1.廣東省電力設(shè)計研究院,廣東 廣州 510663;2.浙江大學(xué)建筑工程學(xué)院,浙江 杭州 310058;3.清華大學(xué)水沙科學(xué)與水利水電工程國家重點實驗室,北京 100084)

在湍流傳質(zhì)和傳熱過程中,濃度擴散系數(shù)和溫度擴散系數(shù)不僅和流體的屬性有關(guān),而且和流動狀態(tài)密切相關(guān)。由于湍流在時間上的非定常性和空間上的各向異性,因而濃度或溫度的擴散系數(shù)是時間和空間的函數(shù),所以相對于物性參數(shù),擴散系數(shù)的取值相對較為困難。在湍流傳質(zhì)研究過程中,一般引入Schmidt數(shù)來表達湍流的擴散能力。Schmidt數(shù)是一個無量綱參數(shù),定義為動量擴散和質(zhì)量擴散之比。

由于湍流Schmidt數(shù)和具體的湍流場有較大關(guān)系,因而很多學(xué)者進行了大量的調(diào)查、試驗和數(shù)學(xué)模型研究。Tominaga等[1]調(diào)查總結(jié)了射流、濁流、邊界層羽流彌散和建筑物附近的流動中的湍流Schmidt數(shù)的經(jīng)驗取值,指出湍流Schmidt數(shù)對預(yù)測結(jié)果具有較大的影響,應(yīng)該根據(jù)不同類型的流動具體研究。He等[2]采用RANS方法對在橫流中的射流的湍流Schmidt數(shù)對射流混合的影響進行了研究,認為此情況下湍流Schmidt數(shù)取0.2能較好地符合試驗數(shù)據(jù)。張曉航等[3]對無剪切湍流混合層中被動標量的擴散進行了數(shù)值模擬,結(jié)合煙粒子擴散風(fēng)洞試驗研究了不同Schmidt數(shù)對被動標量擴散的影響。Xu[4]研究了高Schmidt數(shù)弱擴散條件下湍流中被動標量的混合,計算中Schmidt數(shù)達到了1024,泰勒尺度雷諾數(shù)約為8。Dudukovi?等[5]研究了湍流Schmidt數(shù)對降液膜的質(zhì)量傳輸率的影響,指出湍流雷諾數(shù)通過影響湍流譜來影響湍流Schmidt數(shù),其主要原因來源于高頻和低頻湍流脈動對傳熱傳質(zhì)輸運的不同影響。隨著雷諾數(shù)的增大,不僅湍流強度增大,而且湍流譜從低頻向高頻移動,結(jié)果使得雷諾數(shù)對質(zhì)量輸運的影響不如對動量輸運的影響大,最終使得湍流的Schmidt數(shù)隨著湍流雷諾數(shù)的增大而增大。Flesch等[6]通過大范圍氣體示蹤試驗獲得了大氣邊界層的Schmidt數(shù)(約為0.6)。Ojo等[7]采用水動力觀測數(shù)據(jù)研究了湍流擴散過程及擴散系數(shù)的確定方法。在湍流Schmidt數(shù)識別方面的研究較少,Guo等[8]采用遺傳算法估計了橫流中的射流的變系數(shù)Schmidt數(shù),取得了較好的結(jié)果。

筆者將從概率統(tǒng)計的角度,采用Metropolis-Hastings算法對湍流標量輸運的Schmidt數(shù)進行識別,其中湍流和濃度耦合場的求解采用有限單元法,湍流場求解采用Boussinesq渦黏性假設(shè)。采用Metropolis-Hastings算法識別流動輸運參數(shù)參見文獻[9-12]。

1 數(shù)學(xué)模型

湍流計算采用穩(wěn)態(tài)標準 k-ε模型,濃度場計算采用非穩(wěn)態(tài)對流擴散方程,控制方程如下:

式中:ui為流速;xi,xj,xk分別為i,j,k方向的空間坐標;ρ為密度;p為壓強;η為流體動力黏度;ηt為湍流運動黏度;k為湍動能;ε為耗散率;c為濃度;νt為湍流運動黏度;Dt為湍流擴散系數(shù);Sc為湍流Schmidt數(shù);c1,c2,σk,σε為k-ε模型系數(shù)。

2 參數(shù)識別算法

Metropoils-Hastings算法是由Metropolis提出并由Hastings發(fā)展完善起來的一種概率抽樣算法,屬于馬爾科夫鏈蒙特卡羅方法的一種。它是一種隨機優(yōu)化算法,通過隨機游走的方式在Metropolis-Hastings準則下獲得未知參數(shù)的樣本,在對后驗樣本統(tǒng)計的基礎(chǔ)上獲得未知參數(shù)的估計值。結(jié)合湍流Schmidt數(shù)識別,提出參數(shù)識別算法如下:

第1步:構(gòu)造正態(tài)分布似然函數(shù)作為湍流Schmidt數(shù)隨機優(yōu)化的目標函數(shù)。

第2步:在先驗范圍內(nèi)隨機產(chǎn)生湍流Schmidt數(shù)的初始值。

第3步:采用有限單元法和標準 k-ε模型計算當(dāng)前湍流Schmidt數(shù)Scnow的似然函數(shù)值Lnow。

第4步:對當(dāng)前湍流Schmidt數(shù)進行一個無偏擾動,獲得湍流Schmidt數(shù)新的測試值Sctest,計算在測試值時的似然函數(shù)值Ltest。

第5步:如果 L test大于 L now,那么用 Sc test取代Scnow。

第6步:如果 Ltest小于Lnow,那么產(chǎn)生1個[0,1]區(qū)間的均勻分布的隨機數(shù)m,如果 m小于Ltest和Lnow的比值,那么亦用Sctest取代 Scnow。

第7步:返回第 4步,反復(fù)迭代直到預(yù)定的次數(shù)。

第8步:對所有隨機產(chǎn)生的樣本進行統(tǒng)計以獲得湍流Schmidt數(shù)的估計值。

3 算例研究

3.1 流動正問題的參數(shù)設(shè)置

取一個擴散器內(nèi)的二維湍流擴散為算例。如圖1所示,具有一定濃度的流體從左邊界流入、右邊界流出,邊界條件設(shè)置見表1。進口流速為0.1m/s,進口濃度為1 mol/m3,瞬態(tài)計算總時間為20 s,時間步長為 0.1 s,運動黏度為 10-6m2/s。共剖分了2307個二次拉格朗日型三角形單元。

圖1 二維湍流擴散流動計算域(單位:m)

表1 邊界條件設(shè)置

以下取不同的湍流Schmidt數(shù)進行非穩(wěn)態(tài)的濃度場計算。計算發(fā)現(xiàn),不同湍流Schmidt數(shù)對標量場的時空分布有較大的影響。在相同的射流場內(nèi),當(dāng)湍流Schmidt數(shù)較小時,射流濃度分散性較大,對流作用較弱;當(dāng)湍流Schmidt數(shù)較大時,射流濃度較為集中,擴散作用較弱。由此可以看出標量場對湍流Schmidt數(shù)的敏感性較強,因而湍流Schmidt數(shù)具備良好的可識別性條件。

3.2 參數(shù)識別

圖2 湍流Schmidt數(shù)的迭代過程

圖3 湍流Schmidt數(shù)的后驗概率直方圖

從圖2中可以看出,馬爾科夫鏈隨機游走經(jīng)歷過一個初始化階段(迭代大約100次)后進入了統(tǒng)計收斂域,表明計算參數(shù)的選取是合理的。為了驗證計算結(jié)果的正確性,對迭代100次以后的400個抽樣樣本進行統(tǒng)計,樣本均值為0.40498,均值的估計誤差為1.245%,標準差為0.02199,均值95%置信區(qū)間為[0.40282,0.40715]。由此可以看出基于Metropolis-Hastings算法和有限單元法的參數(shù)識別方法識別精度較高。此外,從圖3中可以看出,湍流Schmidt數(shù)后驗分布具有較好的正態(tài)分布性質(zhì),其可識別性和識別精度都較高。

4 結(jié) 語

為了對湍流輸運過程中的關(guān)鍵控制參數(shù)進行識別,提出基于Metropolis-Hastings算法和有限單元法的參數(shù)識別方法,湍流計算采用標準k-ε模型,標量場計算采用非穩(wěn)態(tài)對流擴散方程。算例計算結(jié)果表明本文提出的算法能對湍流Schmidt數(shù)進行可靠的識別。

[1]TOMINAGA Y,STATHOPOULOS T.Turbulent Schmidt numbers for CFD analysis with various types offlow field[J].Atmospheric Environment,2007,41(37):8091-8099.

[2]HE Guang-bin,GUO Yan-hu,HSU A T.The effectof Schmidt number on turbulent scalar mixing in a jet-in-crossflow[J].International Journal of Heat andmAss Transfer,1999,42:3727-3738.

[3]張曉航,許春曉,張兆順.Sc數(shù)對湍流被動標量擴散的影響[J].工程力學(xué),2004,21(3):36-39.

[4]XU Shu-yi.Turbulent mixing of passive scalars at high Schmidtnumber[D].Georgia:Georgia Institute of Technology,2005.

[5]DUDUKOVI?A,PJANOVIC R.Effect of turbulent Schmidt number onmAss-transfer rates to falling liquid films[J].Ind Eng Chem Res,1999,38:2503-2504.

[6]FLESCH T K,PRUEGER J H,HATFIELD J L.Turbulent Schmidtnumber from a tracer experiment[J].Agricultural and Forest Meteorology,2002,111(4):299-307.

[7]OJO T O,BONNER J S,PAGE C.Studies on turbulent diffusion processes and evaluation of diffusivity values from hydrodynamic observations in Corpus Christi Bay[J].Continental Shelf Research,2006,26(20):2629-2644.

[8]GUO Yan-hu,HE Guang-bin,HSU A T.Application ofgenetic algorithms to the development of a variable Schmidt number model for jet-in-crossflows[J].International Journal of Numerical Methodsfor Heat and Fluid Flow,2001,11(8):744-761.

[9]朱嵩,毛根海,劉國華,等.利用貝葉斯推理估計二維含源對流擴散方程參數(shù)[J].四川大學(xué)學(xué)報:工程科學(xué)版,2008,40(2):38-43.

[10]朱嵩.基于貝葉斯推理的環(huán)境水力學(xué)反問題研究[D].杭州:浙江大學(xué),2008.

[11]朱嵩,劉國華,王立忠,等.水動力-水質(zhì)耦合模型污染源識別的貝葉斯方法[J].四川大學(xué)學(xué)報:工程科學(xué)版,2009,41(5):1-6.

[12]朱嵩,毛根海,劉國華,等.改進MCMC方法及應(yīng)用[J].水利學(xué)報,2009,40(8):1019-1023.