禹言芳 李中根 孟輝波 劉振江 吳劍華

(沈陽(yáng)化工大學(xué) 遼寧省高效化工混合技術(shù)重點(diǎn)實(shí)驗(yàn)室， 沈陽(yáng) 110142)

引 言

靜態(tài)混合器是一種依靠管道內(nèi)周期排列的螺旋元件對(duì)流體不斷的分割、變形、移位及混合的在線、高效混合設(shè)備，與攪拌混合器相比，具有結(jié)構(gòu)簡(jiǎn)單、空間占有率小、停留時(shí)間短、耗能低及維修成本低等優(yōu)點(diǎn)[1-2]，因此廣泛應(yīng)用于精細(xì)化工、水處理、制藥工業(yè)、生物化工等工業(yè)過(guò)程[3]。

掌握靜態(tài)混合器內(nèi)流場(chǎng)的湍流特性規(guī)律對(duì)現(xiàn)有靜態(tài)混合器的結(jié)構(gòu)優(yōu)化及新型靜態(tài)混合元件的開(kāi)發(fā)具有重要意義[4-6]。Tajima等[7]采用Kenics型靜態(tài)混合器(KSM)對(duì)水和液態(tài)CO2生成CO2水合物的混合過(guò)程進(jìn)行實(shí)驗(yàn)研究，結(jié)果表明湍流狀態(tài)下隨著KSM元件數(shù)量增加，液滴尺寸逐漸變小。Jaworski等[8]采用大渦模擬(LES)和非定常雷諾平均NS方程(RANS)兩種方法，對(duì)KSM內(nèi)兩相流的湍動(dòng)能和耗散率等參數(shù)估計(jì)進(jìn)行了有效性分析。Regner等[9]利用計(jì)算流體動(dòng)力學(xué)(CFD)軟件分析了層流狀態(tài)下不同長(zhǎng)徑比的Lightnin靜態(tài)混合器(LSM)內(nèi)的流場(chǎng)結(jié)構(gòu)，研究發(fā)現(xiàn)帶有曲率的混合元件交錯(cuò)處會(huì)產(chǎn)生徑向二次流旋渦。本課題組前期采用大渦模擬分析KSM內(nèi)二次流徑向渦的形成誘因及演化過(guò)程，從多縱向渦的角度探討了扭旋葉片對(duì)流體力學(xué)性能及傳熱性能的影響機(jī)理[10]；并對(duì)LSM中高黏度流體的層流流動(dòng)和混合特性進(jìn)行數(shù)值模擬，發(fā)現(xiàn)LSM內(nèi)高黏度流體的混合性能優(yōu)于KSM[11]；隨后分析了LSM內(nèi)湍流的傳熱特性和旋渦演化過(guò)程，結(jié)果表明半徑尺度的強(qiáng)制渦受到LSM混合元件的影響逐漸耗散且被誘導(dǎo)形成兩對(duì)迪恩渦，強(qiáng)化了LSM的傳熱性能[6]。

本征正交分解(POD)方法被許多學(xué)者用于對(duì)湍流流場(chǎng)的測(cè)量結(jié)果進(jìn)行分析，Rempfer等[12]將POD應(yīng)用于平板邊界層內(nèi)的相干結(jié)構(gòu)和動(dòng)力學(xué)行為的研究中，結(jié)果表明高階的流動(dòng)相干結(jié)構(gòu)對(duì)轉(zhuǎn)捩的產(chǎn)生有著十分重要的作用；余穎等[13]利用POD方法提取內(nèi)置蝸桿傳熱管流場(chǎng)中的擬序結(jié)構(gòu)，發(fā)現(xiàn)多對(duì)大尺度縱向渦結(jié)構(gòu)的產(chǎn)生會(huì)提高強(qiáng)化傳熱效率；張建偉等[14]對(duì)雙組分層對(duì)置撞擊流的湍流速度場(chǎng)進(jìn)行POD分解，發(fā)現(xiàn)低階模態(tài)能夠準(zhǔn)確還原流場(chǎng)中的渦結(jié)構(gòu)。

到目前為止，關(guān)于LSM的研究主要集中在湍流傳熱、層流混合等方面，而針對(duì)LSM的瞬態(tài)湍流特性及分布混合特性的研究尚未見(jiàn)文獻(xiàn)報(bào)道。因此，本文利用POD方法提取流場(chǎng)的湍流結(jié)構(gòu)，以期明晰LSM內(nèi)流體的流動(dòng)特性。采用平均渦量、拉伸率等參數(shù)表征LSM內(nèi)流體的混合特性，為豐富LSM混合機(jī)理提供理論依據(jù)，并擴(kuò)大其在工業(yè)上的應(yīng)用范圍。

1 數(shù)值模型

1.1 控制方程

大渦模擬是流體力學(xué)中一種重要的數(shù)值模擬方法，其基本思想是直接計(jì)算大尺度渦，對(duì)于小渦則通過(guò)亞格子模型進(jìn)行模擬。區(qū)別于直接數(shù)值模擬(DNS)和RANS方法，LES在對(duì)大渦解析的同時(shí)還對(duì)小渦進(jìn)行模化，通過(guò)空間濾波操作過(guò)濾掉小渦，濾波后大渦控制方程為[15]

(1)

(2)

式中，x為笛卡爾坐標(biāo)系；ui和uj為流體速度；ρ、p和ν分別為流體的密度、壓力和運(yùn)動(dòng)黏度。

亞格子應(yīng)力τij定義為

(3)

其中，

(4)

式中，Sij為應(yīng)變速率張量，τkk為亞格子應(yīng)力的各向同性部分。Smagorinsky- Lilly模型定義的亞格子湍流渦黏系數(shù)為[16]

(5)

式中，Ls為亞網(wǎng)格尺度的混合長(zhǎng)度。

1.2 本征正交分解

本征正交分解方法由Lumley[17]首次提出，經(jīng)Sirovich[18]改進(jìn)后發(fā)展為“快照POD”，逐漸被廣大學(xué)者用于復(fù)雜流場(chǎng)的降階及特征提取。假設(shè)在混合器橫截面內(nèi)有m個(gè)離散采樣點(diǎn)，每次采樣同時(shí)獲取每個(gè)采樣點(diǎn)的三維速度值u(xm)。因此，獲取N張快照脈動(dòng)速度場(chǎng)數(shù)據(jù)，以矩陣U的形式表示為[19]

(6)

式中，上標(biāo)N表示快照數(shù)目，下標(biāo)m代表三維速度的矢量位置。自相關(guān)矩陣為

C=UTU

(7)

速度相關(guān)矩陣的特征向量和特征值為

CA=λA

(8)

用式(8)計(jì)算得到的特征值及對(duì)應(yīng)的特征向量去構(gòu)建POD模態(tài)φi的基函數(shù)

(9)

POD系數(shù)

an=ψTun

(10)

線性重構(gòu)任意時(shí)刻的速度場(chǎng)

(11)

1.3 物理模型

假設(shè)圓管直徑D=40 mm，總長(zhǎng)度為560 mm，對(duì)含有元件長(zhǎng)徑比Ar分別為1.0、1.5、2.0和3.0的LSM流場(chǎng)進(jìn)行比較分析。LSM內(nèi)混合元件由兩個(gè)半橢圓形板和兩個(gè)三角形板組成，相鄰兩元件首尾相連且相互交錯(cuò)90°布置安裝，幾何模型如圖1所示。對(duì)于長(zhǎng)徑比Ar=1.5的LSM，兩個(gè)交叉半橢圓板分別與軸截面成34°和146°，其他模型的具體幾何參數(shù)見(jiàn)表1。

表1 幾何參數(shù)Table 1 Geometrical parameters

圖1 不同長(zhǎng)徑比的LSM幾何模型Fig.1 Geometric models of LSM with different aspect ratios

1.4 模擬策略

本文采用ANSYS FLUENT V16.1軟件中的LES湍流模型，對(duì)圖1中4種不同長(zhǎng)徑比的LSM進(jìn)行非穩(wěn)態(tài)數(shù)值模擬。以水(密度ρ=998.2 kg/m3，黏度μ=0.001 003 Pa·s)作為工質(zhì)。采用課題組前期對(duì)LSM湍流流動(dòng)傳熱預(yù)測(cè)中的網(wǎng)格無(wú)關(guān)性驗(yàn)證結(jié)果[6]，利用多面體網(wǎng)格對(duì)計(jì)算域進(jìn)行離散化處理，采用10層邊界層和1.2的網(wǎng)格增長(zhǎng)率對(duì)壁面進(jìn)行細(xì)化，總網(wǎng)格數(shù)量為1 597 947。

采用三維壓力基求解器，并使用SIMPLEC算法對(duì)壓力和速度耦合；梯度項(xiàng)求解采用Least Squares Cell Based離散法；壓力離散采用標(biāo)準(zhǔn)形式，動(dòng)量離散采用Bounded Central Differencing；采用用戶自定義函數(shù)(user define function，UDF)設(shè)置速度入口邊界[6]，對(duì)應(yīng)的雷諾數(shù)分別為1 990、3 981、5 971和7 962。出口邊界條件設(shè)置為outflow。連續(xù)性方程收斂殘差小于10-4，其余變量的收斂殘差為10-6。

1.5 模型有效性

大渦模擬在降低空間和時(shí)間分辨率的前提下可以對(duì)具有較高雷諾數(shù)和較復(fù)雜的湍流流動(dòng)進(jìn)行模擬，計(jì)算精度高。課題組前期采用LES方法對(duì)KSM內(nèi)的流場(chǎng)進(jìn)行研究，數(shù)值計(jì)算的速度分布和實(shí)驗(yàn)一致性較好[10]。LSM和KSM混合器都屬于螺旋流靜態(tài)混合器[9]，具有相似的工作原理。因此，本文選用大渦模擬對(duì)Lightnin型靜態(tài)混合器進(jìn)行研究。

2 POD分析與結(jié)果討論

2.1 采樣數(shù)量選擇

為探究采樣數(shù)量對(duì)POD分析精確度的影響，分析了不同采樣數(shù)量下的模態(tài)能量分布，如圖2所示。從圖中可以看出，隨著模態(tài)階數(shù)s的增加，POD模態(tài)的特征值λ迅速衰減，當(dāng)模態(tài)階數(shù)超過(guò)10之后，高階模態(tài)所攜帶的能量較低。其中，前四階模態(tài)攜帶的能量最多；在采樣數(shù)量超過(guò)2 000時(shí)，前四階的模態(tài)能量占整體的29%以上，且隨著樣本數(shù)量的增加低階模態(tài)所包含的能量不再發(fā)生變化。因此，綜合考慮計(jì)算精度及成本，選取采樣數(shù)量為2 000。

圖2 Re=3 981、Ar=1.5條件下，z/l=2處不同樣本數(shù)量的模態(tài)能量分布Fig.2 Modal energy distributions for different numbers of samples at z/l=2，Re=3 981 and Ar=1.5

2.2 POD模態(tài)收斂性驗(yàn)證

對(duì)LSM數(shù)值模擬所得的樣本進(jìn)行計(jì)算，評(píng)估不同模態(tài)下的二范數(shù)εL2(L2-norm)[20]，具體公式如下。

εL2=‖δi,1-δi,2‖2

(12)

式中，δi,1、δi,2分別為兩次不同采樣頻率樣本中POD分解之后的第i階模態(tài)。

樣本數(shù)量與模態(tài)收斂性之間的關(guān)系如圖3所示。由圖可知，模態(tài)之間的差異隨著樣本數(shù)量的增加逐漸減小。當(dāng)樣本數(shù)量n小于500時(shí)，差異值變化明顯，在樣本數(shù)量很小的情況下，L2-norm可達(dá)到1.3左右。當(dāng)n>1 000時(shí)，L2-norm小于0.3，并且逐漸趨于0，表明經(jīng)POD得到的模態(tài)差異性穩(wěn)定，即POD模態(tài)是收斂的。

圖3 不同樣本數(shù)量下不同模態(tài)的二范數(shù)Fig.3 L2-norm of the modes for differentnumbers of samples

2.3 POD模態(tài)能量分析

圖4為不同軸截面LSM內(nèi)湍流流體前四階模態(tài)的能量分布。從圖4可以看出，不同階的模態(tài)能量在z/l=0～2范圍內(nèi)呈減小趨勢(shì)，表明流場(chǎng)復(fù)雜度不斷增加。在流場(chǎng)復(fù)雜度穩(wěn)定后，1階模態(tài)能量所占的湍動(dòng)能分?jǐn)?shù)在8%～10%，2階模態(tài)能量所占湍動(dòng)能分?jǐn)?shù)在6%～8%，3階模態(tài)能量所占湍動(dòng)能分?jǐn)?shù)在5%～7%，4階模態(tài)能量所占湍動(dòng)能分?jǐn)?shù)在5%～6%。前四階模態(tài)能量占總湍動(dòng)能的24%～31%，包含了主要的大尺度旋渦結(jié)構(gòu)信息。

圖4 Re=3 981和Ar=1.5條件下不同軸向位置的模態(tài)能量分布Fig.4 Modal energy distributions for different axial positions at Re=3 981 and Ar=1.5

圖5為前四階POD模態(tài)的徑向速度矢量圖。由于1階模態(tài)包含著最多的流場(chǎng)信息，從圖5(a)可以清晰地看出在元件壁面存在兩對(duì)大尺度旋渦結(jié)構(gòu)且對(duì)稱(chēng)分布，周?chē)鸁o(wú)小尺度旋渦結(jié)構(gòu)。隨著模態(tài)能量的降低，在第2、3階模態(tài)中小尺度旋渦數(shù)量逐漸增加，所提取的流場(chǎng)信息逐漸減少，這主要是由流場(chǎng)的復(fù)雜程度增加所致。在第4階模態(tài)中旋渦結(jié)構(gòu)開(kāi)始散亂，小尺度的旋渦結(jié)構(gòu)尺寸逐漸增大。

圖5 前四階POD模態(tài)的徑向速度矢量圖Fig.5 Radial velocity vectors of the first four POD modes

對(duì)LSM內(nèi)的原始瞬態(tài)流場(chǎng)(圖6(a))進(jìn)行40階本征正交分解，利用前四階模態(tài)重構(gòu)流場(chǎng)流線圖，如圖6(b)所示。對(duì)比發(fā)現(xiàn)，原始流場(chǎng)在管壁邊緣處存在小尺度渦對(duì)，干擾流場(chǎng)大尺度相干結(jié)構(gòu)的提取與分析。POD方法在時(shí)均速度場(chǎng)的基礎(chǔ)上對(duì)任意時(shí)刻的原始流場(chǎng)進(jìn)行線性重構(gòu)，能夠過(guò)濾小尺度的旋渦結(jié)構(gòu)，同時(shí)保留半圓形流道內(nèi)原始流場(chǎng)的特征，使核心流區(qū)域的徑向二次流大尺度相干結(jié)構(gòu)的空間分布更加明晰。

圖6 Re=3 981時(shí)原始和重構(gòu)流場(chǎng)流線圖對(duì)比Fig.6 Comparison of original and reconstructed flow field streamlines at Re=3 981

2.4 回流系數(shù)

二次流能夠強(qiáng)化靜態(tài)混合器的混合性能。為了更好地了解LSM內(nèi)的二次流分布，通過(guò)計(jì)算整個(gè)采樣過(guò)程流場(chǎng)中選定位置處流速方向與主流方向相反的時(shí)間與總采樣時(shí)間的比值得到回流系數(shù)γt。圖7為在z/l=4.1～6.0處LSM內(nèi)的回流系數(shù)云圖。由圖中的高亮區(qū)域可以看出，回流主要發(fā)生在管壁與元件連接處附近以及元件中心截面處。

圖7 Ar=1.5，Re=3 981時(shí)不同軸截面的回流系數(shù)云圖Fig.7 Cross sectional contours of the reverse flow coefficient under the conditions Re=3 981 and Ar=1.5

2.5 渦量分布

圖8描述了不同雷諾數(shù)下平均渦量〈ω〉隨軸向位置的變化。從圖中可知，受到元件旋向影響，不同軸截面內(nèi)平均渦量分布具有正弦狀的周期性變化規(guī)律。在混合器入口段無(wú)元件分布，沒(méi)有旋渦產(chǎn)生，〈ω〉為0；在z/l=4.0～5.0范圍內(nèi)，〈ω〉逐漸增大，表明在一個(gè)混合元件周期內(nèi)旋渦從正向渦逐漸向反向渦轉(zhuǎn)變，峰值出現(xiàn)在左旋和右旋元件連接處；在z/l=5.0～6.0范圍內(nèi)，〈ω〉逐漸減小，表明旋渦從反向渦逐漸向正向渦轉(zhuǎn)變。隨著雷諾數(shù)的增加，渦量峰值出現(xiàn)不同程度的變化，當(dāng)Re=1 990時(shí)，正向渦值的峰值平均值為6.37，而Re=3 981、5 971和7 962對(duì)應(yīng)的正向渦值峰值平均值分別為15.83、25.33和36.39，渦量分別增大了148.5%、297.6%和322.7%，表明增大Re可提高LSM內(nèi)流體的湍流程度。

圖8 不同雷諾數(shù)下LSM內(nèi)不同軸向截面的平均渦量(Ar=1.5)Fig.8 Mean vorticity of different coaxial sections for different Reynolds numbers(Ar=1.5)

2.6 混合特性

前期的研究成果表明，靜態(tài)混合器內(nèi)的流動(dòng)具有混沌混合特性[21]。采用拉伸率Λn這一指標(biāo)對(duì)不同長(zhǎng)徑比和不同雷諾數(shù)下的LSM分布混合性能進(jìn)行定量分析，拉伸率數(shù)值越大，代表強(qiáng)化分布混合能力越強(qiáng)。根據(jù)Ottino混合理論[22]，每個(gè)粒子所經(jīng)歷的拉伸量可以通過(guò)小矢量的變形求出。

圖9為不同長(zhǎng)徑比和雷諾數(shù)下LSM內(nèi)粒子拉伸率的幾何平均值。由圖9可以發(fā)現(xiàn)，LSM內(nèi)粒子的拉伸率均隨軸向長(zhǎng)度的增加而線性增長(zhǎng)；在同一軸向位置，Ar=1.0對(duì)應(yīng)的拉伸率最大，Ar=3.0對(duì)應(yīng)的拉伸率最小，說(shuō)明長(zhǎng)徑比越小，拉伸效果越明顯。這是由于相同軸向距離內(nèi)小長(zhǎng)徑比的混合葉片扭曲程度更大，徑向二次流更明顯，并且在相同長(zhǎng)度的管路內(nèi)，可以布置的小長(zhǎng)徑比的元件個(gè)數(shù)多于大長(zhǎng)徑比的元件個(gè)數(shù)，因此，流體微元在經(jīng)過(guò)布置有小長(zhǎng)徑比元件的管路中時(shí)被拉伸的概率更大。從圖9(b)可以看出，隨著雷諾數(shù)的增加，LSM的最大拉伸率和橫截面的平均拉伸率均增大，Re=3 981、5 971和7 962下的拉伸率與Re=1 990時(shí)的相比，分別提高了11.5%、15.1%和18.9%。

圖9 不同長(zhǎng)徑比和雷諾數(shù)下的平均拉伸率變化Fig.9 The average stretching rate for different aspect ratios and Reynolds numbers

通過(guò)對(duì)LSM內(nèi)的拉伸率數(shù)值進(jìn)行統(tǒng)計(jì)學(xué)分析，可以定量評(píng)價(jià)流體的微觀混合程度[21，23]。圖10為L(zhǎng)SM在不同長(zhǎng)徑比和不同雷諾數(shù)下出口截面的拉伸率概率密度分布。從圖中可以看出，對(duì)流體有切割分流作用的LSM，其拉伸率概率密度分布曲線形狀類(lèi)似高斯分布曲線，說(shuō)明拉伸率數(shù)值分布較為集中；徑向二次流隨著長(zhǎng)徑比的變小而增強(qiáng)，致使概率密度曲線向著高拉伸率方向移動(dòng)，這是由于流體經(jīng)過(guò)較大的拉伸變形，產(chǎn)生較高的拉伸率。Ar=1.0時(shí)曲線在尾部的震蕩相較于其他長(zhǎng)徑比的要明顯，這說(shuō)明其內(nèi)部的高拉伸率的粒子數(shù)量較多。隨著長(zhǎng)徑比的增加，各曲線峰值逐漸向左移動(dòng)，說(shuō)明小長(zhǎng)徑比的混合元件對(duì)流體的混合能力較強(qiáng)。隨著雷諾數(shù)的增加，拉伸率分布范圍沒(méi)有改變，但各個(gè)截面處平均拉伸率的概率密度增大，總體上提高了混合器的混合能力。

圖10 不同長(zhǎng)徑比和雷諾數(shù)下的拉伸率概率密度分布Fig.10 The probability distributions of stretching rate fordifferent aspect ratios and Reynolds numbers

采用壓降Δp和混合性能強(qiáng)化系數(shù)η來(lái)分析不同長(zhǎng)徑比的LSM的流動(dòng)混合強(qiáng)化行為[11]。不同長(zhǎng)徑比LSM內(nèi)的壓降及混合性能強(qiáng)化系數(shù)的變化如圖11所示。可以看出壓降Δp和混合性能強(qiáng)化系數(shù)η均隨長(zhǎng)徑比的增大逐漸減小。在Re=3 981條件下，Ar=1.0的LSM內(nèi)的壓降Δp=712.9 Pa，混合性能強(qiáng)化系數(shù)η=6.8。隨著長(zhǎng)徑比增加，Ar=1.5、2.0、3.0的LSM的壓降分別是Ar=1.0時(shí)的28.5%、14.2%和0.06%，而混合性能強(qiáng)化系數(shù)分別是Ar=1.0時(shí)的93.5%、85.2%和66.9%。綜上所述，在Ar=1.5時(shí)LSM的壓降明顯減小，同時(shí)具有較好的混合性能。

圖11 Re=3 981條件下不同長(zhǎng)徑比LSM的壓降及混合性能強(qiáng)化系數(shù)Fig.11 Pressure drop and mixing enhancement factor forLSM with different aspect ratios at Re=3 981

3 結(jié)論

(1)采用POD方法對(duì)LSM內(nèi)的瞬態(tài)湍流流場(chǎng)進(jìn)行分析，前四階模態(tài)能量占總湍動(dòng)能的24%～31%，蘊(yùn)含主要的大尺度旋渦結(jié)構(gòu)信息，可消除小尺度旋渦結(jié)構(gòu)對(duì)流場(chǎng)的干擾。

(2)對(duì)于同一長(zhǎng)徑比的LSM，渦量隨著元件的軸向位置變化呈正弦狀周期性變化，隨著雷諾數(shù)的增加，不同軸向截面的平均渦量均增大。

(3)對(duì)拉伸率的分析結(jié)果表明，LSM內(nèi)部的拉伸變形隨著雷諾數(shù)的增加而增加；流體微元在經(jīng)過(guò)布置有小長(zhǎng)徑比元件的管路時(shí)被拉伸的概率更大。

(4)Ar=1.0時(shí)LSM的混合強(qiáng)化性能最好，然而與Ar=1.5、2.0、3.0時(shí)相比壓降較高。綜合壓降和混合性能強(qiáng)化系數(shù)，Ar=1.5的LSM具有較好的混合性能，同時(shí)壓降降低顯著。