劉向東, 戴含暉, 陳永平,,3

(1.揚(yáng)州大學(xué) 電氣與能源動力工程學(xué)院，江蘇 揚(yáng)州 225127；2.東南大學(xué) 能源與環(huán)境學(xué)院，南京 210096；3.蘇州科技大學(xué) 環(huán)境科學(xué)與工程學(xué)院，江蘇 蘇州 215009)

0 引 言

流體力學(xué)是力學(xué)的一個重要分支，主要研究流體內(nèi)的力學(xué)規(guī)律，是水利工程、化學(xué)工程、動力工程等諸多學(xué)科的理論基礎(chǔ)，在航空航天、船舶制造、能源電力等領(lǐng)域有著廣泛應(yīng)用[1-2]。流體力學(xué)有著理論性強(qiáng)、抽象難懂的特點(diǎn)。特別是，在多數(shù)高校的流體力學(xué)教學(xué)中，各類理論性極強(qiáng)的概念、案例往往只能通過陳述/推導(dǎo)公式的方式傳授給學(xué)生，學(xué)生在學(xué)習(xí)過程中難以獲得直觀體驗，這大大削弱了他們的接受度和學(xué)習(xí)興趣，進(jìn)而明顯影響到了教學(xué)效果[3-5]。

隨著計算機(jī)技術(shù)和計算流體力學(xué)理論的發(fā)展，對流體力學(xué)教學(xué)及實(shí)驗案例進(jìn)行數(shù)值虛擬仿真能夠直觀清晰地展現(xiàn)流體力學(xué)相關(guān)規(guī)律，從而為改觀傳統(tǒng)教學(xué)方法形式單一、教學(xué)成效有限的現(xiàn)狀提供了一種有效手段[6-8]。傳統(tǒng)計算流體力學(xué)方法通過數(shù)值化求解以非線性偏微分納維爾-斯托克斯方程組(Navier-Stokes equations)為核心的宏觀連續(xù)性流體力學(xué)控制方程來實(shí)現(xiàn)對流體力學(xué)問題的虛擬仿真，其數(shù)值化編程和計算求解過程都相當(dāng)復(fù)雜，并且求解計算量大、計算效率一般，從而限制了其在流體力學(xué)教學(xué)與實(shí)驗領(lǐng)域的普及。為解決這一問題，本文將使用有別于傳統(tǒng)計算流體力學(xué)方法的格子玻爾茲曼方法(Lattice Boltzmann method，LBM)[9-10]來進(jìn)行典型流體力學(xué)問題的數(shù)值模擬仿真。與傳統(tǒng)計算流體力學(xué)模擬方法不同，LBM基于分子動理論，通過對微觀粒子行為的統(tǒng)計平均來反應(yīng)流體宏觀動力學(xué)特征，是一種典型的介觀模擬方法。相較于傳統(tǒng)計算流體力學(xué)方法，LBM具備先天的并行特性，且物理背景清晰、邊界條件處理簡單、程序易于實(shí)現(xiàn)、計算效率高[11]，已經(jīng)在許多傳統(tǒng)計算流體力學(xué)模擬方法難以勝任的領(lǐng)域都獲得了成功應(yīng)用[12-13]，近年來逐漸受到了國內(nèi)外研究者的關(guān)注。

值得注意的是，Matlab軟件中的圖形用戶界面(Graphic User Interface，GUI)工具可通過添加控件設(shè)計人機(jī)交互界面，便于搭建可視化仿真實(shí)驗平臺，已應(yīng)用于數(shù)學(xué)、物理、力學(xué)等課程教學(xué)中[14-16]。鑒于LBM方法在流體力學(xué)高效模擬計算和Matlab GUI工具在用戶操控與可視化方面所具有的優(yōu)勢，結(jié)合LBM方法和Matlab GUI工具設(shè)計搭建流體力學(xué)虛擬實(shí)驗仿真平臺，用于高校流體力學(xué)教學(xué)實(shí)驗中，為進(jìn)一步豐富流體力學(xué)課堂與實(shí)驗教學(xué)手段、提高課堂與實(shí)驗的教學(xué)效果提供一種嶄新工具。

1 LBM簡介

LBM基于分子動理論，將流體抽象為大量只有質(zhì)量沒有體積的微觀粒子，這些粒子可以向空間的若干方向任意移動。LBM將粒子的演化過程分為碰撞和遷移兩個階段，因此格子玻爾茲曼方程可以分解為兩個部分。

(1)

在數(shù)值模擬過程中根據(jù)流場空間維數(shù)的不同，有不同的離散速度集和與之對應(yīng)的平衡態(tài)分布函數(shù)。針對二維的計算區(qū)域，選擇D2Q9模型(2維空間，9個離散速度，見圖1)，其流場被一系列正方形劃分。其速度、局部平衡分布函數(shù)、權(quán)函數(shù)及聲速分別如下：

(2)

圖1 D2Q9模型示意圖

在獲得速度分布函數(shù)后，流體的宏觀密度與速度通過以下公式計算得到：

(5)

計算過程中固體界面處的邊界節(jié)點(diǎn)使用非平衡外推的方法，該節(jié)點(diǎn)上分布函數(shù)計算方法如下

(6)

式中,x+為壁面法向的相鄰節(jié)點(diǎn)的位置。

LBM的演化過程可分解為碰撞和遷移兩個過程，其計算流程圖如圖2所示，具體的計算步驟如下：

步驟1給定初始密度ρ和速度場u；

步驟2計算平衡態(tài)分布函數(shù)feq；

步驟3在時刻t執(zhí)行碰撞過程，發(fā)生改變的節(jié)點(diǎn)x上的粒子分布函數(shù)為原有粒子分布和碰撞項之和；

步驟4執(zhí)行遷移，節(jié)點(diǎn)x上的粒子以速度ei運(yùn)動到與節(jié)點(diǎn)x相鄰的x+ei，粒子分布函數(shù)也隨之遷移；

步驟5邊界處理；

步驟6計算時刻t+1的密度ρ與速度u；

步驟7重復(fù)步驟3～6直到滿足終止條件。

圖2 LBM方法計算流程圖

2 功能模塊

2.1 整體設(shè)計思路

采用Matlab GUI工具對LBM模擬計算方法進(jìn)行封裝，便可建立起流體力學(xué)虛擬實(shí)驗仿真平臺，其用戶界面如圖3所示。由圖可知，本文所建立的仿真平臺主要包含頂蓋驅(qū)動流、泊肅葉流和卡門渦街模塊。3個模塊對應(yīng)著流體力學(xué)中的幾個經(jīng)典問題。該仿真平臺的邏輯結(jié)構(gòu)如圖4所示。

圖3 流體力學(xué)虛擬實(shí)驗仿真平臺用戶界面

圖4 流體力學(xué)虛擬實(shí)驗仿真平臺模塊結(jié)構(gòu)圖

2.2 頂蓋驅(qū)動流模塊

頂蓋驅(qū)動流是通過腔體頂蓋水平滑移運(yùn)動驅(qū)動腔體內(nèi)流體產(chǎn)生流動的典型流體力學(xué)問題。不同工況條件下的頂蓋驅(qū)動流流場內(nèi)蘊(yùn)含著復(fù)雜流體力學(xué)現(xiàn)象(諸如不穩(wěn)定層流、多級渦及湍流等)，從而使之成為流體力學(xué)教學(xué)中的經(jīng)典物理模型之一。為此，本文建立的流體力學(xué)虛擬實(shí)驗仿真平臺特包含了頂蓋驅(qū)動流模塊，其用戶界面如圖5(a)所示。該模塊為用戶提供了自定義計算區(qū)域尺寸與流動雷諾數(shù)的功能。用戶可通過文本框分別鍵入計算區(qū)域水平方向網(wǎng)格數(shù)、豎直方向網(wǎng)格數(shù)、流動雷諾數(shù)，然后點(diǎn)擊“計算”按鈕完成模擬操作。頂蓋驅(qū)動流的不同流動狀態(tài)主要受流動雷諾數(shù)影響。在輸入計算區(qū)域幾何尺寸后，在顯示計算結(jié)果的區(qū)域內(nèi)會給出網(wǎng)格劃分示意圖(圖5(b))。圖6顯示了不同雷諾數(shù)條件下，同一個方腔內(nèi)的流體流動的LBM模擬計算結(jié)果。

由圖6可見，整個流場中的流線由藍(lán)色線條表示，計算區(qū)域中的速度矢量場用紅色箭頭表示。值得一提的是，為了讓圖中表示速度的矢量箭頭顯示得更加清晰，在繪制過程中對整個速度場做了標(biāo)準(zhǔn)化處理，即所有速度矢量的模都設(shè)為相同。這樣，圖中紅色箭頭便可直觀展示流場內(nèi)各處的速度方向信息。從圖中可以發(fā)現(xiàn)，隨著雷諾數(shù)Re的增大，在中央1級渦的右下部、左下部和左上部會依次逐漸出現(xiàn)2級渦。Re=200時，方腔右下部的2級渦已有形成的趨勢，到Re=1 000時右下部的2級渦已完全成型。當(dāng)Re提升到4 000時方腔左下部會出現(xiàn)2級渦。當(dāng)Re=700時，方腔左上部的2級渦也已初步形成?；谠撃K中流線與流體速度矢量的可視化展示，學(xué)生可以在圖中清晰的獲得方腔內(nèi)流體流動的詳細(xì)流場分布，并能夠直觀地感受到流動狀態(tài)在很大程度上與流體流動的雷諾數(shù)直接相關(guān)。

2.3 泊肅葉流模塊

深入了解泊肅葉流的流動規(guī)律是研究更為復(fù)雜流動問題的基礎(chǔ)。本文建立的流體力學(xué)虛擬實(shí)驗仿真平臺也納入了泊肅葉流模塊，其用戶界面如圖7所示，計算區(qū)域的網(wǎng)格尺寸與入口流速可通過對應(yīng)的文本框鍵入，再點(diǎn)擊“計算”按鈕便可獲得模擬結(jié)果。圖中顯示的是計算區(qū)域內(nèi)水平速度的分布云圖。從圖右側(cè)的顏色圖例可以看出，顏色由淺入深表示的是速度由大到小的變化過程。水平方向的速度分布總體呈現(xiàn)從流道中央向兩側(cè)遞減的規(guī)律。圖中還顯示了中央某界面上的速度分布圖，用紅色箭頭與輪廓線標(biāo)出，這使得計算結(jié)果更加直觀清晰。值得注意的是，在計算區(qū)域左側(cè)流體流入部分存在一定的過渡段。這是由于在給定入口流速時是平均分布的，由速度平均分布演變?yōu)槌浞职l(fā)展段內(nèi)拋物線狀的速度分布需要一定的過渡段。同時，考慮到不同入口流速下泊肅葉流速度云圖分布相似，該模塊還設(shè)置了計算模擬數(shù)據(jù)導(dǎo)出功能，為用戶獲得不同工況下通道內(nèi)充分發(fā)展段主流速度分布的具體計算模擬結(jié)果提供了有效接口。利用該功能，圖8給出了不同入口流速條件下充分發(fā)展段的水平方向速度分布的模擬結(jié)果。如圖所示，學(xué)生可以觀察到兩平板間的流動大致呈現(xiàn)二次方函數(shù)的拋物線分布，垂直于流體流動方向的流動速度由流道中央向通道兩側(cè)壁面遞減，直至壁面上速度降為0。該模塊可以為學(xué)生直觀形象地認(rèn)識流體力學(xué)中典型的泊肅葉流主流速度分布特征提供一種有效手段。

(a) 用戶界面圖

(b) 計算網(wǎng)格圖

圖5 頂蓋驅(qū)動流模塊用戶界面及計算網(wǎng)格圖

(a) Re=200

(b)Re=1 000

(c)Re=4 000

(d)Re=7 000

圖6 不同雷諾數(shù)條件下頂蓋驅(qū)動流的計算結(jié)果

圖7 泊肅葉流模塊用戶界面

圖8 不同入口流速條件下泊肅葉流主流速度分布

2.4 卡門渦街模塊

卡門渦街是當(dāng)定常來流繞過流場中的障礙物時，障礙物兩側(cè)周期性脫落形成旋轉(zhuǎn)方向相反、排列規(guī)則的雙列線渦的流體力學(xué)現(xiàn)象。在不同的來流速度與障礙物尺寸條件下，障礙物之后的卡門渦街特征各異。本文流體力學(xué)虛擬實(shí)驗仿真平臺所包含的卡門渦街模塊用戶界面如圖9所示，演示的是二維圓柱繞流問題。該模塊允許用戶調(diào)整入口流速u、障礙物直徑與流道寬度比pr，可以實(shí)現(xiàn)對卡門渦街特征的模擬與直觀展示。從模擬計算結(jié)果可見，該模塊可以直觀地展示出不同工況下計算區(qū)域內(nèi)流體速度大小連續(xù)分布的流場云圖，計算結(jié)果界面還給出了每個模擬計算工況所對應(yīng)的不同尺寸圓柱障礙物的輪廓示意。圖9(a)～(c)顯示的是圓柱直徑一定、入口流速不同的條件下，各個流場計算結(jié)果的對比。圖9(a)中的入口流速為30 m/s，該條件下整個計算區(qū)域內(nèi)的流動較為平穩(wěn)，在圓柱的后方基本沒有出現(xiàn)線渦脫落的現(xiàn)象。當(dāng)入口流速增大到250 m/s(見圖9(b))，圓柱之后出現(xiàn)了明顯的擾動，但隨著流動向下游發(fā)展，這一擾動也逐漸趨于緩和。而繼續(xù)增加入口流速到500 m/s(見圖9(c))，可以發(fā)現(xiàn)在圓柱之后的整個計算區(qū)域內(nèi)都出現(xiàn)了明顯的線渦脫落現(xiàn)象。當(dāng)在圖9(c)的基礎(chǔ)上將障礙物半徑占流道寬度的比例從0.3減小到0.1時(見圖9(d))，圓柱后的流場波動仍然顯著，但波動幅度相較于圖9(c)有所減小而波動頻率則有所增加，這就說明隨著障礙物半徑占流道寬度比例的減小，卡門渦街所產(chǎn)生的線渦強(qiáng)度有所減小，同時其脫落頻率則有所提升。通過以上算例，學(xué)生可以直觀地認(rèn)識到，在流體繞障礙物流動的過程中，流體流速較高時更易在障礙物后產(chǎn)生卡門渦街，并且渦街產(chǎn)生的頻率和強(qiáng)度與障礙物與流道間的尺寸比例直接相關(guān)。

3 虛擬仿真實(shí)驗平臺在教學(xué)活動中的應(yīng)用方式與優(yōu)勢所在

本文建立的流體力學(xué)虛擬實(shí)驗仿真平臺為學(xué)生提供了一個直觀、清晰認(rèn)知經(jīng)典流體力學(xué)的有效窗口，為改觀傳統(tǒng)陳述/推導(dǎo)公式型的教學(xué)方法提供了契機(jī)。任課教師在課堂教學(xué)過程中可以通過引導(dǎo)學(xué)生在計算機(jī)上進(jìn)行仿真平臺的實(shí)際操作來加深對經(jīng)典流體力學(xué)的理解認(rèn)知，并直觀體會到關(guān)鍵工況參數(shù)對流體流動特性的影響。同時，教師還可以就此向?qū)W生引申介紹后續(xù)《計算流體力學(xué)》《數(shù)值傳熱學(xué)》等相關(guān)課程的基本思想，進(jìn)一步拓展學(xué)生的知識面，并為其未來專業(yè)課的學(xué)習(xí)打下良好基礎(chǔ)。

相對于傳統(tǒng)流體力學(xué)教學(xué)方法，使用Matlab GUI設(shè)計的流體力學(xué)仿真模擬實(shí)驗平臺進(jìn)行輔助教學(xué)的優(yōu)勢明顯：

(1) 具有良好的可視化效果。虛擬實(shí)驗仿真平臺最顯著的優(yōu)點(diǎn)就是具備良好的可視化效果。通過矢量場圖、流線、速度云圖將各個算例清晰直觀地展現(xiàn)在學(xué)生面前，有利于強(qiáng)化學(xué)生對該方面知識的理解，為后續(xù)內(nèi)容的學(xué)習(xí)打好基礎(chǔ)。

(a) u=30 m/s，pr=0.3

(b) u=250 m/s，pr=0.3

(c) u=500 m/s，pr=0.3

(d) u=500 m/s，pr=0.1

(2) 操作簡單，運(yùn)行速度快。本文所建立的虛擬仿真實(shí)驗平臺界面清晰，操作簡便。同時，該平臺基于LBM方法進(jìn)行數(shù)值模擬仿真，且本文已根據(jù)Matlab GUI的運(yùn)行特點(diǎn)對程序內(nèi)的運(yùn)算邏輯做了優(yōu)化。因此，相較于基于傳統(tǒng)計算流體力學(xué)方法的數(shù)值仿真平臺，該平臺運(yùn)行速度得到了明顯提升，可在幾分鐘內(nèi)獲得合理計算結(jié)果。

(3) 其他優(yōu)勢。該平臺所基于的LBM方法，邊界條件處理簡單，程序易于實(shí)現(xiàn)且開源，學(xué)有余力的學(xué)生可以通過自行編譯程序，對幾個流體力學(xué)問題進(jìn)行自主模擬。這充分體現(xiàn)了“因材施教”的教育理念，并有助于進(jìn)一步提升學(xué)生的計算機(jī)程序編譯水平。虛擬仿真實(shí)驗平臺節(jié)約了實(shí)驗過程中的耗材耗能，成本低廉、節(jié)能環(huán)保，從而實(shí)現(xiàn)了“綠色教學(xué)”。

4 結(jié) 語

采用LBM對流體力學(xué)教學(xué)中的經(jīng)典案例進(jìn)行了數(shù)值模擬，并基于Matlab GUI工具對其進(jìn)行封裝，開發(fā)了流體力學(xué)虛擬實(shí)驗仿真平臺。該仿真平臺主要涵蓋了流體力學(xué)知識體系中頂蓋驅(qū)動流、泊肅葉流和卡門渦街3個經(jīng)典案例模塊。學(xué)生可以在平臺用戶界面上輸入可變工況參數(shù)進(jìn)行計算，并可以獲得流場中的速度矢量場、流線等重要信息的直觀展示。同時，學(xué)有余力的學(xué)生還可以通過自行編譯程序，對幾個流體力學(xué)問題進(jìn)行自主模擬。

相對于流體力學(xué)傳統(tǒng)課堂與實(shí)驗教學(xué)方法，本虛擬仿真實(shí)驗平臺具有可視化效果好，計算速度快，界面清晰操作簡單，教學(xué)成本低廉，方便交流等優(yōu)點(diǎn)，有助于激發(fā)學(xué)生的學(xué)習(xí)興趣，加深學(xué)生對抽象概念與復(fù)雜流動狀態(tài)的理解認(rèn)識，課堂與實(shí)驗教學(xué)效果的提高，值得在流體力學(xué)教學(xué)中進(jìn)行推廣。