朱 嵩,劉國華,程偉平,黃躍飛
(1.廣東省電力設(shè)計(jì)研究院,廣東 廣州 510663;2.浙江大學(xué)建筑工程學(xué)院,浙江 杭州 310058;3.清華大學(xué)水沙科學(xué)與水利水電工程國家重點(diǎn)實(shí)驗(yàn)室,北京 100084)
在湍流傳質(zhì)和傳熱過程中,濃度擴(kuò)散系數(shù)和溫度擴(kuò)散系數(shù)不僅和流體的屬性有關(guān),而且和流動(dòng)狀態(tài)密切相關(guān)。由于湍流在時(shí)間上的非定常性和空間上的各向異性,因而濃度或溫度的擴(kuò)散系數(shù)是時(shí)間和空間的函數(shù),所以相對于物性參數(shù),擴(kuò)散系數(shù)的取值相對較為困難。在湍流傳質(zhì)研究過程中,一般引入Schmidt數(shù)來表達(dá)湍流的擴(kuò)散能力。Schmidt數(shù)是一個(gè)無量綱參數(shù),定義為動(dòng)量擴(kuò)散和質(zhì)量擴(kuò)散之比。
由于湍流Schmidt數(shù)和具體的湍流場有較大關(guān)系,因而很多學(xué)者進(jìn)行了大量的調(diào)查、試驗(yàn)和數(shù)學(xué)模型研究。Tominaga等[1]調(diào)查總結(jié)了射流、濁流、邊界層羽流彌散和建筑物附近的流動(dòng)中的湍流Schmidt數(shù)的經(jīng)驗(yàn)取值,指出湍流Schmidt數(shù)對預(yù)測結(jié)果具有較大的影響,應(yīng)該根據(jù)不同類型的流動(dòng)具體研究。He等[2]采用RANS方法對在橫流中的射流的湍流Schmidt數(shù)對射流混合的影響進(jìn)行了研究,認(rèn)為此情況下湍流Schmidt數(shù)取0.2能較好地符合試驗(yàn)數(shù)據(jù)。張曉航等[3]對無剪切湍流混合層中被動(dòng)標(biāo)量的擴(kuò)散進(jìn)行了數(shù)值模擬,結(jié)合煙粒子擴(kuò)散風(fēng)洞試驗(yàn)研究了不同Schmidt數(shù)對被動(dòng)標(biāo)量擴(kuò)散的影響。Xu[4]研究了高Schmidt數(shù)弱擴(kuò)散條件下湍流中被動(dòng)標(biāo)量的混合,計(jì)算中Schmidt數(shù)達(dá)到了1024,泰勒尺度雷諾數(shù)約為8。Dudukovi?等[5]研究了湍流Schmidt數(shù)對降液膜的質(zhì)量傳輸率的影響,指出湍流雷諾數(shù)通過影響湍流譜來影響湍流Schmidt數(shù),其主要原因來源于高頻和低頻湍流脈動(dòng)對傳熱傳質(zhì)輸運(yùn)的不同影響。隨著雷諾數(shù)的增大,不僅湍流強(qiáng)度增大,而且湍流譜從低頻向高頻移動(dòng),結(jié)果使得雷諾數(shù)對質(zhì)量輸運(yùn)的影響不如對動(dòng)量輸運(yùn)的影響大,最終使得湍流的Schmidt數(shù)隨著湍流雷諾數(shù)的增大而增大。Flesch等[6]通過大范圍氣體示蹤試驗(yàn)獲得了大氣邊界層的Schmidt數(shù)(約為0.6)。Ojo等[7]采用水動(dòng)力觀測數(shù)據(jù)研究了湍流擴(kuò)散過程及擴(kuò)散系數(shù)的確定方法。在湍流Schmidt數(shù)識別方面的研究較少,Guo等[8]采用遺傳算法估計(jì)了橫流中的射流的變系數(shù)Schmidt數(shù),取得了較好的結(jié)果。
筆者將從概率統(tǒng)計(jì)的角度,采用Metropolis-Hastings算法對湍流標(biāo)量輸運(yùn)的Schmidt數(shù)進(jìn)行識別,其中湍流和濃度耦合場的求解采用有限單元法,湍流場求解采用Boussinesq渦黏性假設(shè)。采用Metropolis-Hastings算法識別流動(dòng)輸運(yùn)參數(shù)參見文獻(xiàn)[9-12]。
湍流計(jì)算采用穩(wěn)態(tài)標(biāo)準(zhǔn) k-ε模型,濃度場計(jì)算采用非穩(wěn)態(tài)對流擴(kuò)散方程,控制方程如下:

式中:ui為流速;xi,xj,xk分別為i,j,k方向的空間坐標(biāo);ρ為密度;p為壓強(qiáng);η為流體動(dòng)力黏度;ηt為湍流運(yùn)動(dòng)黏度;k為湍動(dòng)能;ε為耗散率;c為濃度;νt為湍流運(yùn)動(dòng)黏度;Dt為湍流擴(kuò)散系數(shù);Sc為湍流Schmidt數(shù);c1,c2,σk,σε為k-ε模型系數(shù)。
Metropoils-Hastings算法是由Metropolis提出并由Hastings發(fā)展完善起來的一種概率抽樣算法,屬于馬爾科夫鏈蒙特卡羅方法的一種。它是一種隨機(jī)優(yōu)化算法,通過隨機(jī)游走的方式在Metropolis-Hastings準(zhǔn)則下獲得未知參數(shù)的樣本,在對后驗(yàn)樣本統(tǒng)計(jì)的基礎(chǔ)上獲得未知參數(shù)的估計(jì)值。結(jié)合湍流Schmidt數(shù)識別,提出參數(shù)識別算法如下:
第1步:構(gòu)造正態(tài)分布似然函數(shù)作為湍流Schmidt數(shù)隨機(jī)優(yōu)化的目標(biāo)函數(shù)。
第2步:在先驗(yàn)范圍內(nèi)隨機(jī)產(chǎn)生湍流Schmidt數(shù)的初始值。
第3步:采用有限單元法和標(biāo)準(zhǔn) k-ε模型計(jì)算當(dāng)前湍流Schmidt數(shù)Scnow的似然函數(shù)值Lnow。
第4步:對當(dāng)前湍流Schmidt數(shù)進(jìn)行一個(gè)無偏擾動(dòng),獲得湍流Schmidt數(shù)新的測試值Sctest,計(jì)算在測試值時(shí)的似然函數(shù)值Ltest。
第5步:如果 L test大于 L now,那么用 Sc test取代Scnow。
第6步:如果 Ltest小于Lnow,那么產(chǎn)生1個(gè)[0,1]區(qū)間的均勻分布的隨機(jī)數(shù)m,如果 m小于Ltest和Lnow的比值,那么亦用Sctest取代 Scnow。
第7步:返回第 4步,反復(fù)迭代直到預(yù)定的次數(shù)。
第8步:對所有隨機(jī)產(chǎn)生的樣本進(jìn)行統(tǒng)計(jì)以獲得湍流Schmidt數(shù)的估計(jì)值。
取一個(gè)擴(kuò)散器內(nèi)的二維湍流擴(kuò)散為算例。如圖1所示,具有一定濃度的流體從左邊界流入、右邊界流出,邊界條件設(shè)置見表1。進(jìn)口流速為0.1m/s,進(jìn)口濃度為1 mol/m3,瞬態(tài)計(jì)算總時(shí)間為20 s,時(shí)間步長為 0.1 s,運(yùn)動(dòng)黏度為 10-6m2/s。共剖分了2307個(gè)二次拉格朗日型三角形單元。

圖1 二維湍流擴(kuò)散流動(dòng)計(jì)算域(單位:m)

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


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

圖3 湍流Schmidt數(shù)的后驗(yàn)概率直方圖
從圖2中可以看出,馬爾科夫鏈隨機(jī)游走經(jīng)歷過一個(gè)初始化階段(迭代大約100次)后進(jìn)入了統(tǒng)計(jì)收斂域,表明計(jì)算參數(shù)的選取是合理的。為了驗(yàn)證計(jì)算結(jié)果的正確性,對迭代100次以后的400個(gè)抽樣樣本進(jìn)行統(tǒng)計(jì),樣本均值為0.40498,均值的估計(jì)誤差為1.245%,標(biāo)準(zhǔn)差為0.02199,均值95%置信區(qū)間為[0.40282,0.40715]。由此可以看出基于Metropolis-Hastings算法和有限單元法的參數(shù)識別方法識別精度較高。此外,從圖3中可以看出,湍流Schmidt數(shù)后驗(yàn)分布具有較好的正態(tài)分布性質(zhì),其可識別性和識別精度都較高。
為了對湍流輸運(yùn)過程中的關(guān)鍵控制參數(shù)進(jìn)行識別,提出基于Metropolis-Hastings算法和有限單元法的參數(shù)識別方法,湍流計(jì)算采用標(biāo)準(zhǔn)k-ε模型,標(biāo)量場計(jì)算采用非穩(wěn)態(tài)對流擴(kuò)散方程。算例計(jì)算結(jié)果表明本文提出的算法能對湍流Schmidt數(shù)進(jìn)行可靠的識別。
[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ù)對湍流被動(dòng)標(biāo)量擴(kuò)散的影響[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]朱嵩,毛根海,劉國華,等.利用貝葉斯推理估計(jì)二維含源對流擴(kuò)散方程參數(shù)[J].四川大學(xué)學(xué)報(bào):工程科學(xué)版,2008,40(2):38-43.
[10]朱嵩.基于貝葉斯推理的環(huán)境水力學(xué)反問題研究[D].杭州:浙江大學(xué),2008.
[11]朱嵩,劉國華,王立忠,等.水動(dòng)力-水質(zhì)耦合模型污染源識別的貝葉斯方法[J].四川大學(xué)學(xué)報(bào):工程科學(xué)版,2009,41(5):1-6.
[12]朱嵩,毛根海,劉國華,等.改進(jìn)MCMC方法及應(yīng)用[J].水利學(xué)報(bào),2009,40(8):1019-1023.