趙曉賦,王 林,楊文韜,張健欣
(內蒙古工業大學自動化系,呼和浩特 010080)
微精餾系統的優點包括使用安全性高、體積小、傳質傳熱的效率高以及易于操控,已經廣泛應用于很多領域,其中包括生物醫學、有機合成和食品加工等。在微尺度效應的影響下,微通道與常規尺寸的通道相比,兩相流型存在很大不同[1]。在微精餾系統中研究最多的、最重要的是微通道內氣-液兩相流的流動。受不同熱邊界條件的影響,在微通道內,傳熱特性對流質的流動有著顯著影響[2-3];微通道內局部熱點對流動換熱、流型又有很明顯的影響[4-5]。在微通道中,改變混合流質的流量對流動及傳熱特性有顯著的影響[6-7]。
顧娟利用微元玻爾茲曼方法,研究了在恒壁溫與恒熱流2種熱邊界條件下,稀薄效應對微通道內流質流動、傳熱特性的影響。微通道內氣體受流動特性與稀薄的變化相同,傳熱特性與稀薄的變化截然相反[2]。KEEPAIBOON 等在水力直徑0.68 mm矩形微通道中,在傳熱方面發現傳熱系數與流型有一定關系。在低熱流范圍情況下,傳熱系數隨著飽和溫度的變化具有一致性[3]。影響通道內流質的氣化率的因素有很多,其中包括氣相和液相表觀流速的變化和通道壁面加熱溫度的變化。申利梅等采用VOF(volume of fluid)方法,對局部熱點對微通道內壓降、換熱系數、流型和溫度的影響進行了分析[4]。王鑫寶等利用Fluent 軟件研究了微通道在滑移區的流動換熱情況[5]。劉冬稔和劉肖研究了當微通道內流質的流型為彈狀流時,當氣相的流量不變時,液相的流量發生改變對流型的影響很小[6]。
本研究建立微通道模型,并利用計算機仿真軟件Fluent對其進行求解。在入口流質比例不變的情況下,改變入口流速、管壁加熱溫度,通過Fluent的后處理計算得到當前工況下流質的氣化率。利用最小二乘法,建立流質氣化率與流質入口流速、管壁加熱溫度的多元回歸預測模型,并應用其進行不同流速下、不同加熱溫度下的氣化率預測。
KANDLIKAR 針對通道內沸騰問題,基于水力直徑得出了通道尺度的劃分[8]:常規通道(Conventional channels)→3 mm;細通道(Minichannels)→200 μm~3 mm;微通道(Microchannels)→10 μm~200 μm。
本研究建立模型的截面積為800 μm×100 μm,長度為5 mm。水力直徑的計算:

式中,R為水力直徑,A為過水斷面面積,X為過水斷面上水流所濕潤的邊界長度。
將上述條件代入公式(1)得R=44.4 μm,符合微通道的水力直徑范圍。
利用仿真軟件Fluent對模型進行求解,模型的入口流速范圍是由現場實驗設備的流量范圍確定的。
實驗設備是型號為2PB-10005Ⅲ的平流泵,其體積流量為0.1~100 mL/min。模型的流速范圍:

其中:x為換算后模型的流速,S為模型的截面積,y為平流泵的體積流量。
將平流泵的體積流量0.1~100 mL/min、模型的截面積S=0.08 mm2代入式(2),求得模型的流速為20.84 mm/s~20.83 m/s。因此Fluent 求解器對模型進行求解時,流速范圍應當在此范圍內。
對模型進行流體仿真時,選擇液態酒精、氣態酒精、液態水和氣態水完成氣液混合相的確定。流體屬性參數如表1所示。

表1 流體屬性參數Tab 1 Fluid property parameters
仿真采用VOF 模型,在模型入口處液態酒精、氣態酒精、液態水、氣態水的體積比為0.7:0:0.3:0,在出口處測量混合流質的氣化率,因此將主相phase-1設置為液態酒精,將phase-2設置為氣態酒精,phase-3 設置為液態水,phase-4 設置為氣態水。在模型的計算中會發生相變,液態酒精→氣態酒精相變溫度78.4 ℃,氣態水→氣態酒精相變溫度100 ℃。
2.3.1 模型的建立
微通道模型的建立采用與現場實驗所定制的微流控芯片尺寸一致。因此微通道模型是截面積為800 μm×100 μm,長度為5 mm的管道。
2.3.2 網格的劃分
劃分后的網格按照網格點之間的相接關系可分為非結構化網格劃分和結構化網格劃分。結構化網格的網格點之間都是有順序的并且規則的,此類網格內部各層的網格都要相同的網格數;當幾何模型的形狀不是很簡單時,結構化網格很難劃分出符合要求的網格[9]。非結構化網格之間的網格點是無順序的、不規則的,其劃分復雜的結構更為準確,因此選用非結構化網格對微通道模型進行網格劃分。
邊界層是在高Re繞流中緊貼壁面的粘性力必須考慮的流動薄層。在任何流體域的壁面上都需要設置邊界層。在此微通道模型內含有多種氣液混合,并且通道內發生液態酒精到氣態酒精的相變,為保證結果收斂和計算精度,需要設置邊界層,進行局部網格加密[10]。邊界層設置情況見圖1,邊界層數為5。

圖1 邊界層設置Fig 1 Boundary layer setting
2.3.3 網格無關性檢驗
在進行數值模擬之前,要先對建立的模型進行網格劃分,網格質量和網格的大小對仿真的最終結果起著十分重要作用,后面數值計算結果的準確性直接被其影響。因此對網格的無關性的檢驗在利用Fluent 軟件對模型進行求解時顯得尤為重要。
按流質的流速為20.84 mm/s,壁面的加熱溫度分別設定為99 ℃,時間步長1 ms,對A、B、C、D 4 種網格數進行無關性檢驗,得到的結果如表2 所示。由表2可知,網格數量對模型氣化率結果具有極大影響,當網格超過10 萬個時,模擬結果變化差異不大。

表2 網格無關性檢驗結果Tab 2 Grid independence test results
結合模型的尺寸,選擇合理的網格大小,生成非結構性網格之后,最終得到網格單元總數為54 051。檢查劃分后所有網格偏度都小于0.8 以下,從而滿足計算精度要求。微通道模型的網格劃分如圖2所示。

圖2 模型網格劃分Fig 2 Model grid division
2.4.1 VOF模型
在液態酒精和水的加熱過程中,會產生氣態酒精,從而產生兩相流的混合問題,Fluent軟件中提供了多種多項流模型,其中流體體積(VOF)方法用于捕獲多相流交界處中的拓撲變化[11]。
模擬流體的非定常運動的總過程通常采用VOF 方法。在流體域中的每個網格,此函數為混合流質的體積與網格體積的比。對于氣-液兩相流,氣體的體積分數α定義為:
α=單位流體體積/單位體積。(3)
當α=1 時,區域內只有液體;當α=0 時,所定義的區域內只有氣體;當0<α<1時,區域內氣液均存在。
對于VOF 方法跟蹤相之間的界面是通過求解連續方程得到的。對于不可壓縮的牛頓液體,連續性方程和動量方程及流體運輸方程為:

式中,▽為梯度算子;u、ρ、μ、p和g分別為混合流質的速度矢量、密度、黏度、壓力和重力加速度。C為流體體積函數,表示網格內混合流質體積與網格體積的比,C=0表示網格內不含目標流體,0<C<1表示網格內存在不同的流體,C=1示網格內全部是目標流體。
2.4.2 層流模型
雷諾數(Re)是一個無因次數群,用來判別具有粘性力的流體流動狀態。圓管雷諾數的計算:
式中,ρ為混合流質密度,μ為混合流質的動力粘性系數,v為流體域的特征速度,L為特征長度。
當Re<2 300 時,流動是層流;當Re>2 300時,流動是湍流;當Re=2 300時,流動是臨界流。
對于非圓斷面管流,計算雷諾數時需要引入水力直徑R,其計算如式(1)所示。代替圓管雷諾數中的特征長度L,經計算本算例的Re<2 300,因此在模型中流體的流動是層流。
2.4.3 邊界條件
模擬中,求解器為壓力基求解器,時間為瞬態,由于是微通道模型,因此不考慮重力影響,模型需要對管壁提供連續的熱量,開啟能量守恒方程。Fluent 的主要邊界條件包括入口邊界條件、出口邊界條件和壁面條件,這些量是一些數學和物理變量。邊界條件的設定:采用標準的層流模型方程,速度入口作為入口邊界條件,混合流質的流速設置為7 個:2.084、6.25、12.5、18.756、25.0、31.26、37.512 cm/s,入口處溫度設置為30 ℃。出口邊界為壓力出口,出口處溫度同樣設置為30 ℃。壁面條件的設定:壁面的加熱溫度分別設定為99~79 ℃間每隔2 ℃。
由多個自變量組成的組合{x1,x2…,xn}(n≥2)和因變量y 構成隨機樣本數據集(x1,x2……,xn,y),并存在一種如式(8)的線性關系,這種關系通常叫做多元線性回歸。

就可以確保模型的預測值最接近實際值[12]。在構建多元非線性回歸預測模型的過程中,提高模型預測精度的重要前提是選取契合具體的觀測數據的模型函數。
預測模型的精度可以通過回歸模型的檢驗指標實現,主要包括擬合優度檢驗和顯著性檢驗[13]。
使用判定系數R2作為擬合優度判定模型的擬合優度的主要依據。定義為:

式中,SSR為預測平方和,SST為總離差平方和,SSE為殘差平方和。
對預測模型擬合程度的高低,根據判定系數來確定。判定系數越大,預測模型與原始數據越擬合,判定系數越小,說明擬合程度越差,與原始數據偏差更大。R2的取0~1,因此R2越接近于1,模型的擬合效果就越理想。
采用F檢驗對建立的預測方程進行顯著性檢驗,F統計量定義為:

式中,m為自變量個數,n為樣本數,F為統計量服從第1自由度m、第2自由度為n-m-1的F分布,即F~(m,n-m-1)。
由F統計量的定義可以得到,當F偏大時,在對因變量造成的變化方面相比,自變量比隨機因素的影響大。由此可得,F統計量越大,回歸方程的擬合效果也就。
預測模型的建模過程包括數據預處理、建模、分析及預測4大部分,其具體步驟如圖3所示。

圖3 多元回歸流程Fig 3 Multiple regression process
在Fluent 仿真計算中,所得數據共7 組(7 種流速條件),總共77 個樣本點(7 種流速×11 個壁面溫度),其中1 組數據如表3 所示。其中氣化率為出口處氣相組分(氣態酒精)所占面積與出口截面面積的比。


表3 仿真結果(部分)Tab 3 Simulation results(partial)
測量值比預測值之間的差值,通過標準化殘差可以明確的表示出來,并且還能得到在絕對值上本次測量的樣本點的殘差與其他殘差之間的差值。當標準化殘差的絕對值大于等于某一數值時(在些選擇0.5%-線性&0.3%-非線性),表示此殘差對應的樣本點存在問題。在剔除某樣本點后,回歸方程的標準差有所下降,即認為此樣本點為異常,在擬合時,就對此樣本點進行刪減。線性和非線性預測模型的殘差如圖4和圖5所示。

圖4 線性回歸殘差Fig 4 Residual plot of linear regression

圖5 非線性回歸殘差Fig 5 Residual plot of nonlinear regression
將多元線性預測模型的樣本點中第1、8、13、50 個樣本點剔除,和多元非線性預測模型的樣本點中第1、6、8、13、50、71個樣本點剔除后的標準差如表4所示。

表4 剔除異常值前、后的標準差Tab 4 Standard deviation before and after excluding outliers
表5為回歸模型擬合檢驗結果。

表5 回歸模型的檢驗結果Tab 6 Estimated results of the regression model
由表5 可以看出,非線性模型比線性模型在R2統計量上提高了0.178,在F統計量上提高了32.90。
圖6為非線性預測方程的三維擬合。

圖6 非線性回歸模型擬合Fig 6 Fitting graph of nonlinear regression model
綜上所述,非線性回歸模型擬合結果達到了預期的目標。
為使氣化率預測模型對本課題的通用性得到驗證,利用預測模型分別在3 種工況下求出氣化率,然后與Fluent的計算結果進行比較,結果如表6 所示。考慮到微小截面氣相面積分的不精確性,設定的預測精度指標為20%。
由表6 可以看出,3 個仿真模型的氣化率誤差均低于預測模型精度指標,預測模型的有效性和可行性得到了充分的驗證。

表6 3組仿真結果Tab 6 Three sets of CFD results
通過利用仿真軟件Fluent對建立三維微通道模型在不同工況下,對模型出口處的截面混合流質的氣化率進行了求解。利用最小二乘法對77 組原始數據處理,建立多元回歸預測模型,通過線性預測模型和非線性預測模型的比較,確定了非線性模型結構。利用這個預測模型,對3種獨立的入口條件進行了計算。通過和Fluent 的計算結果比較,預測誤差均低于20%的預測精度指標。
通過建立入口流速、壁面加熱溫度與出口氣化率之間的關系,可以定量分析兩相流的傳質傳熱特性,為構建微精餾系統數學模型做好理論準備。
在下一步的工作中,還將在實驗平臺上,繼續驗證預測模型的精度,以期得到實驗-仿真-預測模型三者計算結果的統一。