鄭馨瑤 楊 剛 張會(huì)臣
(大連海事大學(xué)船舶與海洋工程學(xué)院 遼寧大連 116026)
微通道是指水力直徑在微米量級(jí)的流動(dòng)通道,具有高的表面體積比、高的傳熱效率及較少的流體使用量等特點(diǎn),廣泛應(yīng)用于微流體系統(tǒng)[1-3]。由于尺度的限制,微通道內(nèi)流體的流動(dòng)特性很難通過(guò)實(shí)驗(yàn)完全揭示清楚,并且有些實(shí)驗(yàn)現(xiàn)象的分析也亟待深入。對(duì)微通道內(nèi)流體流動(dòng)進(jìn)行數(shù)值模擬,是對(duì)實(shí)驗(yàn)中難以觀測(cè)記錄現(xiàn)象和對(duì)實(shí)驗(yàn)現(xiàn)象加以分析的有效手段[4]。這些數(shù)值模擬方法主要包括計(jì)算流體動(dòng)力學(xué)、格子玻爾茲曼方法以及分子動(dòng)力學(xué)。相比于宏觀尺度下的計(jì)算流體動(dòng)力學(xué)方法與微觀尺度下的分子動(dòng)力學(xué)方法,格子玻爾茲曼方法是一種介于流體連續(xù)性假設(shè)與分子動(dòng)力學(xué)之間的方法,節(jié)省計(jì)算資源的同時(shí),可揭示微通道內(nèi)流體流動(dòng)的機(jī)制[5]。
早期一些學(xué)者借助格子Boltzmann方法模擬微通道內(nèi)氣體的流動(dòng)。TANG等[6]采用格子玻爾茲曼方法模擬了微通道中的等溫氣體流動(dòng),將反彈邊界和鏡面反射格式組合模擬滑移邊界條件[7],推導(dǎo)出了松弛時(shí)間與努森數(shù)的關(guān)系。LIU和GUO[8]采用格子玻爾茲曼方法研究了長(zhǎng)矩形微通道內(nèi)壓力驅(qū)動(dòng)的氣體流動(dòng),并與其他方法所得數(shù)據(jù)進(jìn)行了對(duì)比。此外,由于格子玻爾茲曼方法易于模擬復(fù)雜邊界,還被廣泛應(yīng)用于液滴動(dòng)力學(xué)特性的模擬。梁宏等人[9]采用格子 Boltzmann 方法模擬了液滴在分叉微通道中的遷移過(guò)程,發(fā)現(xiàn)當(dāng)毛細(xì)數(shù)足夠大且壁面具有親水性時(shí),液滴會(huì)發(fā)生二次分裂;而毛細(xì)數(shù)足夠小時(shí),液滴在分叉處滯留并且不發(fā)生分裂,通過(guò)調(diào)節(jié)出口流量比可以使液滴實(shí)現(xiàn)非對(duì)稱分裂。LIM、KEIVAN等[10-11]采用格子玻爾茲曼方法模擬了微通道中液滴的形成,研究表明液滴的大小以及液滴之間的距離與微通道的幾何形狀關(guān)系密切。
應(yīng)用格子Boltzmann方法,通過(guò)局部計(jì)算剪切應(yīng)力可模擬非牛頓流體的流動(dòng)。DONG等[12]將流變方程引入格子Boltzmann方法中,建立了非牛頓流體兩相流的格子Boltzmann模型,模擬了存在矩形粗糙單元的微通道內(nèi)冪律流體非牛頓流體的位移,發(fā)現(xiàn)冪律指數(shù)的降低導(dǎo)致位移速度的增加,相對(duì)粗糙度的減小和粗糙單元間距的增大可加大位移速度,從而減少位移時(shí)間。謝馳宇等[13]將非牛頓流體模型引入格子Boltzmann方法中,研究多孔復(fù)雜結(jié)構(gòu)內(nèi)多相非牛頓流體的驅(qū)替過(guò)程,發(fā)現(xiàn)了界面張力在驅(qū)替過(guò)程中起阻礙作用,并分析了壓差、潤(rùn)濕性對(duì)驅(qū)替過(guò)程的影響。SHI和TANG[14]采用格子Boltzmann方法模擬黏性指進(jìn)現(xiàn)象,發(fā)現(xiàn)剪切增稠流體中的黏性指進(jìn)增強(qiáng),剪切稀化流體中的黏性指進(jìn)被抑制,重力加速度導(dǎo)致“手指”形狀不對(duì)稱,高毛細(xì)數(shù)與高黏度比導(dǎo)致“手指”狹窄。該研究從介觀的角度解釋了黏性指進(jìn)現(xiàn)象。
盡管通過(guò)實(shí)驗(yàn)研究與數(shù)值模擬,對(duì)微通道內(nèi)流體流動(dòng)的過(guò)程有了一定程度的了解,但對(duì)于流動(dòng)過(guò)程中的細(xì)節(jié)仍不夠清晰。本文作者采用格子Boltzmann方法,模擬不同角度的分叉微通道內(nèi)剪切稀化流體的流動(dòng),通過(guò)流動(dòng)過(guò)程中密度隨時(shí)間的變化趨勢(shì)以及穩(wěn)態(tài)流動(dòng)下的密度,得到微通道內(nèi)壓力的分布以及流動(dòng)區(qū)間的壓力降,著重分析溶液質(zhì)量分?jǐn)?shù)、入口速度與分叉角度對(duì)非牛頓流體流動(dòng)特性的影響,闡明流體特性和微通道幾何構(gòu)型對(duì)非牛頓流體流動(dòng)行為的影響機(jī)制。
文中模擬的是高度和寬度均為300 μm的矩形截面多角度分叉微通道,分叉角度分別為60°、75°、90°、105°和120°。圖1所示為分叉角度為120°的微通道結(jié)構(gòu)示意圖。微通道具有2個(gè)出口,非牛頓流體從左側(cè)入口以固定速度流入,流經(jīng)分叉點(diǎn)處后由2個(gè)出口流出。
圖1 微通道結(jié)構(gòu)示意Fig 1 Structure of microchannel
選取3種不同質(zhì)量分?jǐn)?shù)的聚丙烯酰胺(PAM)水溶液進(jìn)行模擬,探討流變特性對(duì)其流動(dòng)特性的影響。聚丙烯酰胺水溶液為剪切稀化流體,其流變特性采用冪律模型描述,稠度系數(shù)K和冪律指數(shù)n2個(gè)參數(shù)根據(jù)流變儀測(cè)得的實(shí)驗(yàn)數(shù)據(jù)擬合得到,K和n及其他特性參數(shù)如表1所示。
表1 流體的特性參數(shù)
首先采用格子Boltzmann方法求解流體速度場(chǎng),在引入BGK碰撞算子[15]后,具有二次精度的格子Boltzmann方程離散為
(1)
平衡態(tài)分布函數(shù)的具體形式一般與離散速度模型的構(gòu)造有關(guān),因此離散速度模型的構(gòu)造至關(guān)重要。文中采用DnQm模型[16],D2Q9是最具代表性的格子BGK模型,二維D2Q9模型如圖2所示。
圖2 D2Q9模型Fig 2 D2Q9 model
對(duì)于二維D2Q9,其局部平衡分布函數(shù)的表達(dá)式為
(2)
式中:wk表示k方向的權(quán)重因子;cs為格子聲速;ck為沿著格子遷移方向的速度矢量。
流體的宏觀密度經(jīng)平衡分布函數(shù)求和得到
(3)
壓力由宏觀密度得到[17]
(4)
廣義牛頓流體中的冪律模型形式簡(jiǎn)單,廣泛使用,且與剪切稀化流體的真實(shí)特性基本相符,所以選用冪律模型描述非牛頓流體的動(dòng)力黏度
(5)
式中:K為稠度系數(shù);n為冪律指數(shù);γ為剪切速率。
采用根據(jù)剪切黏度調(diào)整松弛時(shí)間的直接方法,進(jìn)行格子玻爾茲曼方法中非牛頓流體特性的模擬[18-20]。首先計(jì)算格子的黏度,再根據(jù)黏度調(diào)整松弛時(shí)間。其中,通過(guò)式(6)迭代計(jì)算形變速率張量,根據(jù)式(7)對(duì)剪切速率進(jìn)行計(jì)算。
(6)
(7)
最后通過(guò)格子玻爾茲曼方法中松弛時(shí)間與黏度的關(guān)系調(diào)整松弛時(shí)間,實(shí)現(xiàn)剪切稀化特性的模擬。
(8)
在格子玻爾茲曼方法的模擬中,一般采用格子單位,與實(shí)際有關(guān)的問(wèn)題模擬時(shí)需要進(jìn)行物理單位與格子單位的轉(zhuǎn)化。對(duì)于不可壓縮流體流動(dòng)問(wèn)題,需要匹配宏觀雷諾數(shù)和格子雷諾數(shù)。宏觀雷諾數(shù)經(jīng)驗(yàn)公式為
(9)
通常,格子速度Ulattice的數(shù)量級(jí)保持在0.1或者0.2,取值與宏觀速度無(wú)關(guān)。文中Ulattice選擇0.1,微通道截面的實(shí)際尺寸為300 μm,其格子單位取值為14,匹配宏觀雷諾數(shù)和格子雷諾數(shù)后得到密度ρ的初始值。其中,根據(jù)微通道截面形狀,a的取值為0.212 1,b的取值為0.676 6。
圖3所示是模擬趨于穩(wěn)定后入口截面處與分叉點(diǎn)處的壓力差值。數(shù)據(jù)整體分布在負(fù)值區(qū)域,表明從入口處流經(jīng)到分叉點(diǎn)處時(shí)壓力增大。從圖中可以看出,90°分叉角度壓力差值最小;當(dāng)分叉角度大于90°時(shí),壓力差值隨著分叉角度的增大而減??;當(dāng)分叉角度小于90°時(shí),壓力差值隨著分叉角度的增大而增大;壓力差值最大值出現(xiàn)在中線偏向分叉點(diǎn)壁面一側(cè),越靠近壁面壓力差值越小。從圖中還可以發(fā)現(xiàn),當(dāng)分叉角度互補(bǔ)時(shí),其壓力差值結(jié)果相近,相差分別為10.72%和12.51%,且分叉角度小于90°的微通道壓力差值大于分叉角度大于90°的微通道。
圖3 不同分叉角度下分叉點(diǎn)截面處的 壓力降(w=0.1%,v=2 m/s)Fig 3 Pressure drop at the cross section of the bifurcation point in different bifurcation angles(w=0.1%,v=2 m/s)
出口1截面處的壓力降如圖4所示,整體分布趨勢(shì)呈拋物線形,壓力降最大值出現(xiàn)在微通道中間截面處,越靠近壁面,其壓力降值越小;90°分叉角度時(shí)壓力降最小,壓力降最大值出現(xiàn)在分叉角度為75°和105°處,最大值相差4.96%。從圖4(a)、(b)中可以看出,入口速度為2 m/s時(shí),壓力降的最大值相差4.96%,入口速度為4 m/s時(shí),壓力降的最大值相差8.94%。從圖4(a)、(c)中可以看出,溶液質(zhì)量分?jǐn)?shù)為0.1%時(shí),壓力降的最大值相差4.96%,溶液質(zhì)量分?jǐn)?shù)為0.3%時(shí),壓力降的最大值相差22.57%,說(shuō)明隨著溶液質(zhì)量分?jǐn)?shù)的增大,分叉角度對(duì)出口1處的壓力降影響更大。因此,相比于入口速度,溶液質(zhì)量分?jǐn)?shù)使得分叉角度對(duì)出口1壓力降的影響更為顯著。
出口2截面處的壓力降如圖5所示,不論分叉角度如何,其壓力降趨勢(shì)一致。
圖6所示是模擬趨于穩(wěn)定后入口截面處與分叉點(diǎn)處的壓力差值,可以看出,入口速度為10 m/s時(shí)壓力差值最小,為2 m/s時(shí)壓力差值最大,即入口速度越大,壓力差值越小;而截面處壓力差值最大值出現(xiàn)在中線偏向分叉點(diǎn)壁面一側(cè),越靠近壁面壓力差值越小。
圖6 不同入口速度下分叉點(diǎn)截面處的 壓力降(w=0.1%,θ=120°)Fig 6 Pressure drop at the cross section of the bifurcation point at different inlet velocities(w=0.1%,θ=120°)
出口1截面處的壓力降如圖7所示。從圖7(a)、(b)中可以看出,當(dāng)溶液質(zhì)量分?jǐn)?shù)一定時(shí),分叉角度為120°時(shí),壓力降的最大值相差2.32%,分叉角度為105°時(shí),壓力降的最大值相差5.92%。從圖7(a)、(c)中可以看出,當(dāng)分叉角度一定時(shí),溶液質(zhì)量分?jǐn)?shù)為0.1%的情況下,壓力降的最大值相差2.32%,溶液質(zhì)量分?jǐn)?shù)為0.3%時(shí),壓力降的最大值相差18.42%??梢?jiàn),隨著溶液質(zhì)量分?jǐn)?shù)的增大,使入口速度對(duì)出口1處的壓力降影響更大。所以相比于分叉角度,溶液質(zhì)量分?jǐn)?shù)使入口速度對(duì)出口1壓力降的影響更為顯著。
圖7 不同入口速度下出口1截面處的壓力降Fig 7 Pressure drop at the cross section of the outlet 1 at different inlet velocities (a) w=0.1%,θ=120°; (b)w=0.1%,θ=105°;(c)w=0.3%,θ=120°
出口2截面處的壓力降如圖8所示。從圖8(a)、(b)中可以看出,不同入口速度條件下出口2截面處壓力降基本相同。從圖8(a)、(c)中可以看出,當(dāng)分叉角度一定時(shí),溶液質(zhì)量分?jǐn)?shù)為0.1%時(shí),不同入口速度下的微通道出口2截面處的壓力降趨勢(shì)一致,而溶液質(zhì)量分?jǐn)?shù)為0.3%時(shí),不同入口速度下壓力降的最大值相差2.69%,且壓力降隨著入口速度的增大而減小。可見(jiàn)隨著溶液質(zhì)量分?jǐn)?shù)的增大使入口速度對(duì)出口1處的壓力降影響更大,所以相比于分叉角度,溶液質(zhì)量分?jǐn)?shù)使入口速度對(duì)出口2壓力降的影響更為顯著。
圖8 不同入口速度下出口2截面處的壓力降Fig 8 Pressure drop at the cross section of the outlet 2 at different inlet velocities (a) w=0.1%,θ=120°; (b)w=0.1%,θ=105°;(c)w=0.3%,θ=120°
圖9所示是模擬趨于穩(wěn)定后入口截面處與分叉點(diǎn)處壓力差值,可以看出,水的壓力差值比任一質(zhì)量分?jǐn)?shù)剪切稀化流體的壓力差值都小;溶液質(zhì)量分?jǐn)?shù)為0.5%時(shí)壓力差值最大,隨著溶液質(zhì)量分?jǐn)?shù)增大,壓力差值增大。溶液質(zhì)量分?jǐn)?shù)為0.3%和0.5%時(shí),微通道分叉點(diǎn)截面處的壓力差值基本相同,分析其原因是在這2種質(zhì)量分?jǐn)?shù)下,描述其非牛頓特性的稠度系數(shù)K和冪律指數(shù)n的值相近,其流變特性的相似導(dǎo)致了流動(dòng)過(guò)程中產(chǎn)生的壓力差值相近。
圖9 不同溶液質(zhì)量分?jǐn)?shù)下分叉點(diǎn)截面 處的壓力降(v=2 m/s,θ=120°)Fig 9 Pressure drop at the cross section of the bifurcation point at different solution concentrations(v=2 m/s,θ=120°)
出口1截面處的壓力降如圖10所示,壓力降的最大值相差3.01%;溶液質(zhì)量分?jǐn)?shù)為0.5%時(shí),壓力降最大,隨著溶液質(zhì)量分?jǐn)?shù)的降低,壓力降減小。由于水是非牛頓流體,沒(méi)有流變特性,導(dǎo)致水的壓力降遠(yuǎn)遠(yuǎn)小于剪切稀化流體的壓力降。
圖10 不同溶液質(zhì)量分?jǐn)?shù)下出口1截面處的 壓力降(v=2 m/s,θ=120°)Fig 10 Pressure drop at the cross section of the outlet 1 at different solution concentrations(v=2 m/s,θ=120°)
出口1處壓力值隨時(shí)間的變化趨勢(shì)如圖11所示。根據(jù)圖11將模擬過(guò)程分為3個(gè)階段,在模擬初期(0~1 000),壓力值隨著溶液質(zhì)量分?jǐn)?shù)的增大而增大,且溶液質(zhì)量分?jǐn)?shù)為0.3%與0.5%的微通道內(nèi)壓力峰值差值為27.5%,遠(yuǎn)遠(yuǎn)小于與溶液質(zhì)量分?jǐn)?shù)為0.1%的壓力峰值差值55.21%。在這一階段,水的壓力值已經(jīng)達(dá)到穩(wěn)定。隨著模擬的進(jìn)行(1 000~3 000),溶液質(zhì)量分?jǐn)?shù)大的微通道內(nèi),壓力值迅速降低,但是在水的模擬過(guò)程中,模擬開(kāi)始不久就達(dá)到穩(wěn)定,沒(méi)有壓力值降低的階段。所以剪切稀化流體在流動(dòng)過(guò)程中產(chǎn)生了壓力值大幅度降低的現(xiàn)象,分析其原因是由于剪切稀化流體的流變特性導(dǎo)致的,因?yàn)槠漯ざ入S著剪切速率的增加而減小,導(dǎo)致壓力也隨之減小。第三階段即模擬穩(wěn)定后(3 000~5 000),壓力值的降低幅度隨著溶液質(zhì)量分?jǐn)?shù)的增大而增大,與出口截面處壓力降結(jié)果一致。
圖11 不同溶液質(zhì)量分?jǐn)?shù)下出口1處壓力隨時(shí)間的 變化趨勢(shì)(v=2 m/s,θ=120°)Fig 11 The trend of pressure over time at the outlet 1 at different solution concentrations(v=2 m/s,θ=120°)
出口2截面處的壓力降如圖12所示。不同溶液質(zhì)量分?jǐn)?shù)條件下,壓力降的最大值相差3.55%;溶液質(zhì)量分?jǐn)?shù)為0.1%時(shí),壓力降最大;隨著溶液質(zhì)量分?jǐn)?shù)的增大,壓力降下降。
圖12 不同溶液質(zhì)量分?jǐn)?shù)下出口2截面處 的壓力降(v=2 m/s,θ=120°)Fig 12 Pressure drop at the cross section of the outlet 2 at different solution concentrations(v=2 m/s,θ=120°)
通過(guò)將非牛頓流體冪律模型引入牛頓流體格子Boltzmann模型的方法,數(shù)值模擬了矩形截面微通道內(nèi)剪切稀化流體的流動(dòng)特性,探討了流體流變特性和微通道幾何構(gòu)型對(duì)非牛頓流體流動(dòng)行為的影響,得到如下結(jié)論:
(1)分叉角度為90°的微通道壓力降最小,當(dāng)分叉角度大于90°時(shí),壓力降隨著分叉角度的增大而減小,當(dāng)分叉角度小于90°時(shí),壓力降隨著分叉角度的增大而增大,且當(dāng)微通道分叉角度在數(shù)學(xué)概念上互補(bǔ)時(shí),其壓力降相近。
(2)微通道內(nèi)各截面處壓力降呈拋物線形,最大值出現(xiàn)在微通道中間截面處,越靠近壁面壓力降越小;在溶液質(zhì)量分?jǐn)?shù)與分叉角度一定時(shí),隨著入口速度的增大,壓力降增大;相比于分叉角度,溶液質(zhì)量分?jǐn)?shù)使得入口速度對(duì)出口處壓力降的影響更為顯著。
(3)在入口速度與分叉角度一定時(shí),隨著溶液質(zhì)量分?jǐn)?shù)的增大,壓力降增大,且當(dāng)流體的流變特性相近時(shí),其流動(dòng)特性相近;而相比于入口速度,溶液質(zhì)量分?jǐn)?shù)使得分叉角度對(duì)出口處壓力降的影響更為顯著。