高Re數(shù)層流管道中顆粒聚集特性的數(shù)值研究*

劉唐京，王企鯤，鄒 赫

（上海理工大學(xué) 能源與動(dòng)力工程學(xué)院，上海 200093）

引言

均勻的稀釋懸浮液以層流方式通過直圓管時(shí)，流體中的剛性球形顆粒最終會(huì)遷移到半徑約為0.6R（R為管道半徑）的圓環(huán)上，這個(gè)環(huán)也被稱為Segre-Silberberg 環(huán)[1].這表明顆粒在隨流體運(yùn)動(dòng)時(shí)，除受到主流方向的驅(qū)動(dòng)力外，在徑向上還受到一個(gè)升力使得顆粒發(fā)生遷移.這個(gè)徑向上的升力是由于運(yùn)動(dòng)流體的慣性引起的，因此又被稱為“慣性升力”，由其“慣性升力”而引發(fā)的顆粒聚集現(xiàn)象稱為“慣性聚集”[2-3].

這種現(xiàn)象是在低Re數(shù)（Re為Reynolds 數(shù)）下管道流中發(fā)現(xiàn)的，之后有學(xué)者在高Re數(shù)流中也發(fā)現(xiàn)了顆粒的慣性聚集現(xiàn)象.Matas 等[4]對(duì)高Re數(shù)下懸浮粒子在Poiseuille 流中的“慣性聚集”進(jìn)行實(shí)驗(yàn)研究，得到了不同Re數(shù)下管道直徑與顆粒直徑比在8～42 范圍內(nèi)顆粒的聚集位置.Morita 等[5]通過實(shí)驗(yàn)發(fā)現(xiàn)，Re數(shù)在1 000范圍內(nèi)的稀釋懸浮液通過直圓管時(shí)，只有一個(gè)穩(wěn)定的聚集點(diǎn).由于實(shí)驗(yàn)很難獲取流場(chǎng)中的各項(xiàng)力學(xué)參數(shù)，隨著數(shù)值模擬（CFD）在流體力學(xué)中的廣泛運(yùn)用[6-8]，一些學(xué)者為了揭示顆粒所受慣性升力的形成特性，通過數(shù)值計(jì)算方法對(duì)顆粒的“慣性聚集”進(jìn)行了研究[9-10].但由于流動(dòng)是非定常的，計(jì)算模型的數(shù)學(xué)描述比較困難，Carlo 等[11]首次提出“相對(duì)運(yùn)動(dòng)模型”，將原本的非定常問題轉(zhuǎn)化為準(zhǔn)定常問題.基于相對(duì)運(yùn)動(dòng)模型，王企鯤[12-13]對(duì)方管中顆粒的受力特性進(jìn)行了數(shù)值研究，揭示了顆粒聚集的力學(xué)成因并歸納出顆粒穩(wěn)定聚集點(diǎn)的水動(dòng)力學(xué)判據(jù).

相對(duì)運(yùn)動(dòng)模型的提出，極大地簡(jiǎn)化了計(jì)算模型.但采用相對(duì)運(yùn)動(dòng)模型進(jìn)行數(shù)值模擬時(shí)，為了使流動(dòng)達(dá)到充分發(fā)展?fàn)顟B(tài)，一般需要預(yù)留L=0.058Re·D（D為管徑）[14]的管長(zhǎng)以消除入口段的影響.這在計(jì)算高Re數(shù)工況時(shí)就需要非常長(zhǎng)的管道，造成網(wǎng)格數(shù)量大，使得數(shù)值計(jì)算變得困難.考慮到管內(nèi)流動(dòng)屬于周期性流動(dòng)，因此本文對(duì)管道進(jìn)、出口設(shè)置周期性邊界條件，解決了高Re數(shù)工況下管道較長(zhǎng)的問題，同時(shí)計(jì)算高Re數(shù)管流中顆粒慣性聚集的受力特性.

1 計(jì)算模型與計(jì)算方法

1.1 計(jì)算模型

本文考慮單個(gè)剛性球形顆粒在直圓管Poiseuille 流中運(yùn)動(dòng)，管道直徑為D，長(zhǎng)度為L(zhǎng)，小球直徑為d，如圖1所示.目前比較常用的模型是6 自由度模型，然而這種計(jì)算模型是非定常的，模型的數(shù)學(xué)描述比較復(fù)雜，并且計(jì)算過程繁瑣也很耗時(shí)，此外6 自由度模型雖然能獲得顆粒的運(yùn)動(dòng)軌跡，但無(wú)法獲取顆粒升力在通道內(nèi)的空間分布規(guī)律.本文采用的相對(duì)運(yùn)動(dòng)模型是在6 自由度模型的基礎(chǔ)上進(jìn)行簡(jiǎn)化，它雖然無(wú)法獲得顆粒的運(yùn)動(dòng)軌跡，但能計(jì)算出顆粒在通道內(nèi)的受力特性.通過分析受力特性來(lái)確定顆粒的聚集點(diǎn)，這在文獻(xiàn)[11-13]中給出了比較詳細(xì)的討論.

圖1 計(jì)算模型示意圖Fig.1 The sketch for the numerical model

相對(duì)運(yùn)動(dòng)模型是通過創(chuàng)建一個(gè)慣性坐標(biāo)系將顆粒的平移速度轉(zhuǎn)移給管壁，使得計(jì)算時(shí)的流場(chǎng)為準(zhǔn)定常，從而簡(jiǎn)化計(jì)算.即將計(jì)算坐標(biāo)系設(shè)置在顆粒的中心（如圖1所示）并隨顆粒一同平移，則在此慣性坐標(biāo)系下，顆粒的平移速度為零，僅存在以原點(diǎn)為中心的旋轉(zhuǎn)速度，管壁則以顆粒的速率反向移動(dòng).

1.2 計(jì)算方法

本文基于有限體積法進(jìn)行三維數(shù)值模擬，求解穩(wěn)態(tài)Navier-Stokes 方程，其控制方程如式（1）所示.流動(dòng)介質(zhì)為常溫常壓液態(tài)水，考慮到顆粒是懸浮顆粒，取其密度與水相同，為ρ=1 000 kg/m3.在CFD 計(jì)算中，采用雙精度來(lái)進(jìn)行計(jì)算，壓力與速度的耦合采用coupled 算法，壓力方程的離散采用二階格式，動(dòng)量方程的離散采用三階精度的QUICK 格式[14-15]：

式中，u為平移坐標(biāo)系Oxyz中流體的相對(duì)速度，p為壓強(qiáng)，ρ 和υ 分別為流體的密度和運(yùn)動(dòng)黏度.

在CFD 計(jì)算中，顆粒的平移速度UP被轉(zhuǎn)化為管壁運(yùn)動(dòng)的邊界條件，這里需要通過試湊的方式來(lái)獲得顆粒在該位置的最終恒定運(yùn)動(dòng)速度.首先假定一個(gè)速度UP0進(jìn)行試算，輸出顆粒在沿流動(dòng)方向（x方向）上的合力Fx0，在一定精度內(nèi)，判斷合力是否為零（零指的是零量階，不是數(shù)值上為零）；若不為零則需不斷地迭代更新壁面速度，直到顆粒在流動(dòng)方向上合力為零.同時(shí)，顆粒的旋轉(zhuǎn)速度也用同樣的方式進(jìn)行試湊，直到顆粒在所有方向上轉(zhuǎn)矩為零.此時(shí)，顆粒在y方向上所受的力便為徑向升力.

為了提高計(jì)算效率，平移速度的初始參考值為UP=2U[1?(r+)2]，旋轉(zhuǎn)速度的初始參考值為對(duì)于UP和ω 的更新，本文利用具有超線性收斂性的“割線法”更新下一步迭代數(shù)據(jù)，分別由式（2）、（3）計(jì)算，采用這種方法通常試湊3 至4 次就能得到結(jié)果，而試湊一次也只需迭代150 次左右：

式中，UPn和ωPn分別為平移速度和旋轉(zhuǎn)速度，F(xiàn)xn和MPn分別為阻力和轉(zhuǎn)矩，n=0,1,2.

顆粒的升力Fy和阻力Fx分別按式（4）、（5）計(jì)算：

式中，i和j分別為x方向和y方向的單位向量，P為應(yīng)力張量，n為顆粒表面外法線單位向量，dS為面積微元，Σ 為顆粒表面.

顆粒的轉(zhuǎn)矩按式（6）計(jì)算：

式中，r為由顆粒中心指向顆粒表面的矢徑.

在CFD 計(jì)算中，當(dāng)對(duì)管道進(jìn)、出口采用周期性邊界條件時(shí)，計(jì)算區(qū)域內(nèi)的每一處壓力分為周期性壓力和線性壓力，而在CFD 實(shí)際的計(jì)算過程中只顯現(xiàn)出周期性壓力少了線性壓力，因此真實(shí)的壓力應(yīng)為：preal=p+β?x（p為周期性壓力，β為一個(gè)周期內(nèi)的平均壓力梯度，?x=x1?x2）[14].那么，式(4)～ (6)在計(jì)算顆粒的升力、阻力以及轉(zhuǎn)矩時(shí)未包含線性壓力提供的那部分力.線性壓力對(duì)顆粒貢獻(xiàn)的升力、阻力及轉(zhuǎn)矩由式(7)、(8)計(jì)算，從中可以發(fā)現(xiàn)，線性壓力對(duì)顆粒提供的力只在流動(dòng)方向上，并且對(duì)轉(zhuǎn)矩的貢獻(xiàn)為零，也就是說(shuō)線性壓力只對(duì)顆粒的阻力有影響：

式中，I為單位二階張量.因此顆粒在流向上的實(shí)際阻力為

為了方便下文的討論與分析，本文定義如下無(wú)量綱參數(shù)：

升力系數(shù)CFL為

式中，d為顆粒直徑，D為管道直徑.顆粒的無(wú)量綱直徑為

顆粒無(wú)量綱徑向位置為

式中，r為管中心與顆粒球心的距離，R為管道半徑.Re數(shù)為

式中U為管內(nèi)流體的平均速度，μ為流體的動(dòng)力黏度.擾動(dòng)強(qiáng)度為

式中，v和w分別為y方向和z方向的流速分量.

2 結(jié)果與討論

2.1 網(wǎng)格無(wú)關(guān)性驗(yàn)證

本文采用結(jié)構(gòu)化網(wǎng)格對(duì)計(jì)算域進(jìn)行網(wǎng)格劃分，如圖2（a）所示.為了提高計(jì)算的穩(wěn)定性和精度，在流體與顆粒間的邊界層網(wǎng)格進(jìn)行了加密處理，如圖2（b）所示.為了保證計(jì)算的準(zhǔn)確性和經(jīng)濟(jì)性，以及獲得顆粒的真實(shí)受力，本文先對(duì)網(wǎng)格無(wú)關(guān)性進(jìn)行了驗(yàn)證.本次驗(yàn)證的計(jì)算工況為：a+=1/9，顆粒的徑向位置r+=0.1，管長(zhǎng)L=4D.

圖2 網(wǎng)格示意圖：（a）管道網(wǎng)格；（b）顆粒周圍網(wǎng)格Fig.2 Schematic diagrams of the grids:a）pipeline grid;b）grids around particle

從圖3 中可看出，當(dāng)網(wǎng)格數(shù)量達(dá)到50 萬(wàn)左右時(shí)，顆粒的升力系數(shù)基本保持不變，考慮到在滿足計(jì)算精度的同時(shí)盡可能節(jié)約計(jì)算時(shí)間，本文最終確定計(jì)算域的網(wǎng)格數(shù)量控制在60 萬(wàn)左右.對(duì)于下文采用不同管長(zhǎng)計(jì)算的工況，網(wǎng)格劃分將以本次驗(yàn)證的管道網(wǎng)格為基準(zhǔn)，網(wǎng)格尺寸保持一致.

圖3 網(wǎng)格無(wú)關(guān)性驗(yàn)證Fig.3 Grid independence verifications

2.2 周期性邊條的可行性分析

在低Re數(shù)下，文獻(xiàn)[2-13]基于相對(duì)運(yùn)動(dòng)模型，邊界條件為進(jìn)口給定均勻的相對(duì)速度，出口為壓力出口.而在高Re數(shù)下，這種邊界條件通常需要很長(zhǎng)的管道才能計(jì)算出可靠的結(jié)果，因此考慮采用周期性邊界條件來(lái)減小管長(zhǎng).為了驗(yàn)證周期性邊界條件的可行性，本文對(duì)Re=350，a+=1/9 的工況進(jìn)行了數(shù)值模擬.本次模擬分別選取管長(zhǎng)L=4D和L=50D的管道，L=4D的管道對(duì)進(jìn)、出口采用周期性邊界條件，L=50D的管道則與文獻(xiàn)[13]采用相同的邊界條件，模擬結(jié)果如圖4所示.

由圖4 可知，兩種邊界條件得到的升力系數(shù)曲線完全重合，升力系數(shù)為零且該點(diǎn)一階導(dǎo)數(shù)小于零的點(diǎn)即為顆粒的穩(wěn)定聚集點(diǎn)，因此Re=350 時(shí)a+=1/9 的顆粒主要聚集在r+≈0.76 處.而文獻(xiàn)[4]的實(shí)驗(yàn)結(jié)果為r+≈0.77，兩者比較吻合，這說(shuō)明周期性邊界條件是可行的.當(dāng)Re=350 時(shí)，對(duì)進(jìn)、出口采用周期性邊界條件只需4 倍管徑長(zhǎng)度的管道便可計(jì)算出結(jié)果，而用文獻(xiàn)[13]中的邊界條件卻需要50 倍管徑長(zhǎng)度.因此在求解高Re數(shù)流中顆粒的慣性聚集時(shí)，采用周期性邊界條件可以有效地減小管長(zhǎng)，降低計(jì)算量.

圖4 不同邊界條件的模擬結(jié)果Fig.4 Simulation results under different boundary conditions

2.3 周期長(zhǎng)度的確定

當(dāng)計(jì)算高Re數(shù)工況采用周期性邊界條件時(shí)，需知道管長(zhǎng)L（周期長(zhǎng)度）為多少時(shí)，獲得的計(jì)算結(jié)果才真實(shí)可靠.針對(duì)這個(gè)問題，本文采用不同周期長(zhǎng)度的管道對(duì)無(wú)量綱直徑a+=1/9 的顆粒進(jìn)行了模擬分析，Re=350.

從圖5 中可看出，當(dāng)L≥3D時(shí)，隨著周期長(zhǎng)度的增加，計(jì)算結(jié)果將保持不變，并且與實(shí)驗(yàn)結(jié)果是吻合的（2.2 小節(jié)中已做對(duì)比），這從力的角度來(lái)看L=3D的管道便可以得到穩(wěn)定的計(jì)算結(jié)果.在相對(duì)運(yùn)動(dòng)模型中，顆粒是相對(duì)靜止的，這會(huì)對(duì)管中的流體產(chǎn)生擾動(dòng)，而較短的周期長(zhǎng)度有可能會(huì)使這種擾動(dòng)延伸到邊界上，影響計(jì)算結(jié)果.而且這有可能會(huì)出現(xiàn)不同周期長(zhǎng)度得到相同的升力分布，但流場(chǎng)卻不一定是相同的情況.為了確保計(jì)算結(jié)果的可靠性，本文對(duì)不同周期長(zhǎng)度顆粒附近以及靠近進(jìn)、出口處的擾動(dòng)強(qiáng)度進(jìn)行了對(duì)比，如圖6所示.考慮到越靠近壁面，顆粒對(duì)流體造成的擾動(dòng)越強(qiáng)，因此本文選取顆?？拷诿妫╮+=0.8）的工況進(jìn)行對(duì)比.

圖5 不同周期長(zhǎng)度下的升力分布Fig.5 The lift distribution under different period lengths

從圖6 中可以發(fā)現(xiàn)，隨著周期長(zhǎng)度的增加，靠近管道進(jìn)、出口處的擾動(dòng)強(qiáng)度是不斷減小的，當(dāng)L≥4D時(shí)，靠近進(jìn)、出口處的擾動(dòng)基本可以忽略不計(jì).這說(shuō)明當(dāng)周期長(zhǎng)度大于4D時(shí)，流場(chǎng)將不會(huì)在發(fā)生變化，結(jié)合上文升力分布結(jié)果，對(duì)于Re=350 的工況，L=4D是可行的計(jì)算周期.而對(duì)更高的Re數(shù)，本文也進(jìn)行了驗(yàn)證，結(jié)果如圖7所示.當(dāng)Re數(shù)達(dá)到800 時(shí)，計(jì)算結(jié)果也符合：1）靠近管道進(jìn)、出口處擾動(dòng)強(qiáng)度為零；2）顆粒的聚集點(diǎn)與文獻(xiàn)[4]的實(shí)驗(yàn)結(jié)果相符.因此，對(duì)于Re<1 000、a+≤1/9 的工況，本文確定4D為計(jì)算周期.

圖6 不同周期長(zhǎng)度下的擾動(dòng)強(qiáng)度對(duì)比：（a）L=2D；（b）L=3D；（c）L=4D；（d）L=5DFig.6 Comparisons of disturbance intensities under different period lengths:a）L=2D;b）L=3D;c）L=4D;d）L=5D

圖7 Re=800 的計(jì)算結(jié)果：（a）升力分布；（b）擾動(dòng)強(qiáng)度Fig.7 Calculation results for Re=800:a）lift distribution;b）disturbance intensity

2.4 周期性邊條的應(yīng)用

2.4.1Re<1 000

前文對(duì)周期性邊條的可行性進(jìn)行了驗(yàn)證，并對(duì)Re<1 000 工況給出了一個(gè)可行的計(jì)算周期為L(zhǎng)=4D.在此基礎(chǔ)上，本文在這里對(duì)不同粒徑的顆粒進(jìn)行了數(shù)值模擬，研究不同Re數(shù)下顆粒的受力特性以及對(duì)顆粒聚集點(diǎn)的影響.

如圖8所示，在低Re數(shù)下，顆粒在徑向上升力分布是類拋物線，不同大小的顆粒主要聚集在r+=0.6 ～0.7 之間，與管壁有一定的距離.而隨著Re數(shù)的增大，顆粒的升力分布以及聚集點(diǎn)都出現(xiàn)了明顯的變化.主要表現(xiàn)如下：高Re數(shù)下顆粒的升力不再呈類拋物線分布，顆粒受到的升力具有一定的波動(dòng)，在r+=0.5 ～ 0.7 之間出現(xiàn)一段升力相對(duì)較小區(qū)域，而在這個(gè)區(qū)域內(nèi)有出現(xiàn)第二個(gè)聚集點(diǎn)（內(nèi)環(huán)）的趨勢(shì).此外，隨著Re數(shù)的增大，顆粒主要聚集在r+=0.75 ～ 0.85 之間（外環(huán)），向著壁面靠近.即Re數(shù)越大，顆粒的聚集位置越靠近壁面；而與顆粒粒徑的關(guān)系則與之相反，粒徑越大，聚集位置越向著管中心遷移.

圖8 不同Re 下顆粒的升力分布：（a）Re=50；（b）Re=350；（c）Re=500；（d）Re=800Fig.8 The lift distribution of particles under different values of Re:a）Re=50;b）Re=350;c）Re=500;d）Re=800

在本次的計(jì)算結(jié)果中，只發(fā)現(xiàn)了一個(gè)穩(wěn)定聚集點(diǎn)，這與Morita 等[5]的實(shí)驗(yàn)結(jié)果是一致的，當(dāng)Re<1 000 時(shí)，如果管道足夠長(zhǎng)，內(nèi)環(huán)將消失，所有的粒子都將聚集在外環(huán)上.而在Matas 等[4]的實(shí)驗(yàn)中，他們?cè)诠艿赖纳嫌螀^(qū)域發(fā)現(xiàn)了內(nèi)環(huán)的存在，但在下游區(qū)域觀察到內(nèi)環(huán)上的顆粒有向外環(huán)遷移的趨勢(shì).本文認(rèn)為這可能是Matas 等[4]實(shí)驗(yàn)的管道不夠長(zhǎng)，顆粒的遷移未完全發(fā)展，顆粒要脫離r+=0.5 ～ 0.7 這個(gè)升力相對(duì)較小的區(qū)域需要較長(zhǎng)的時(shí)間.

2.4.2Re≥1 000

從圖8 的升力分布來(lái)看，在更高Re數(shù)下有可能出現(xiàn)第二個(gè)穩(wěn)定聚集點(diǎn).為了探究Re≥1 000 時(shí)是否有第二個(gè)穩(wěn)定聚集點(diǎn)的出現(xiàn)，本文選用了a+=1/17 的顆粒進(jìn)行模擬.

本次模擬仍是用L=4D的管道進(jìn)行計(jì)算，首先對(duì)周期長(zhǎng)度的可靠性進(jìn)行驗(yàn)證，結(jié)果如圖9（a）所示.從其擾動(dòng)強(qiáng)度來(lái)看，4 倍管徑周期長(zhǎng)度符合計(jì)算精度，且顆粒的升力分布與上文中高Re數(shù)的分布特征一樣.因此本文認(rèn)為對(duì)于Re≤1 600，a+=1/17 的工況，L=4D的管道依然是可行的.

從圖9（b）顆粒的升力分布可以發(fā)現(xiàn)，當(dāng)Re數(shù)達(dá)到1 200 時(shí)，a+=1/17 的顆粒在徑向上有三個(gè)聚集點(diǎn)，其中有兩個(gè)是穩(wěn)定聚集點(diǎn)，分別在r+≈0.63，0.87 處.這說(shuō)明當(dāng)Re>1 000 時(shí)，對(duì)于小粒徑的顆粒是有可能存在兩個(gè)穩(wěn)定聚集點(diǎn)的.由于大粒徑的顆粒在Re數(shù)達(dá)到1 000 時(shí)計(jì)算不穩(wěn)定，所以對(duì)于更大粒徑的顆粒本文沒有繼續(xù)進(jìn)行深入研究.

圖9 a + =1/17，Re≥1 000 的計(jì)算結(jié)果：（a）Re=1 600 時(shí)的擾動(dòng)強(qiáng)度；（b）升力分布Fig.9 The calculation results of a+ =1/17,Re≥1 000:a）the disturbance intensity at Re=1 600;b）the lift distribution

2.5 流場(chǎng)分析

為了探究低Re數(shù)和高Re數(shù)下管內(nèi)顆粒慣性升力分布不同的原因，本文對(duì)顆粒所在橫截面（即x=0 截面）的流場(chǎng)進(jìn)行了分析，以a+=1/9 的顆粒為例，如圖10 ～ 12所示.圖10 ～ 12 分別為r+=0.4，0.6，0.8 時(shí)不同Re數(shù)下z方向的速度云圖和該截面上的速度矢量圖，對(duì)z方向的速度無(wú)量綱化為.

圖10 x=0 截面的速度云圖和矢量圖（r+ =0.4）：（a）Re=50；（b）Re=350；（c）Re=500；（d）Re=800Fig.10 Velocity contours and velocity vectors of section x = 0r+ =0.4）:a）Re=50;b）Re=350;c）Re=500;d）Re=800

圖12 x=0 截面的速度云圖和矢量圖（r+ =0.8）：（a）Re=50；（b）Re=350；（c）Re=500；（d）Re=800Fig.12 Velocity contours and velocity vectors of section x = 0r+ =0.8）:a）Re=50;b）Re=350;c）Re=500;d）Re=800

從圖11、12 可以明顯地看出，在顆粒的周圍有二次流產(chǎn)生，而二次流可能會(huì)對(duì)顆粒的升力造成影響.從速度云圖及矢量圖來(lái)看，當(dāng)Re=50 時(shí)，顆粒周圍并沒有明顯的二次流動(dòng)，尤其是顆粒更靠近通道中心時(shí)，但隨著顆粒向壁面靠近，可以發(fā)現(xiàn)有微弱的二次流產(chǎn)生.然而當(dāng)Re≥350 時(shí)，即使顆粒更靠近通道中心也會(huì)有微弱的二次流產(chǎn)生，且隨著Re數(shù)的增大以及顆粒向壁面靠近，二次流變得越來(lái)越強(qiáng)烈.此外從速度矢量圖可以發(fā)現(xiàn)：在低Re數(shù)時(shí)，二次流主要向著顆粒的下方流動(dòng)；而在高Re數(shù)時(shí)，二次流在顆?？拷ǖ乐行臅r(shí)先是向顆粒下方流動(dòng)，而后隨著顆粒向壁面靠近以及Re數(shù)的增大，其逐漸向顆粒左右兩側(cè)流動(dòng).

圖11 x=0 截面的速度云圖和矢量圖（r+ =0.6）：（a）Re=50；（b）Re=350；（c）Re=500；（d）Re=800Fig.11 Velocity contours and velocity vectors of section x = 0（r+ =0.6）:a）Re=50;b）Re=350;c）Re=500;d）Re=800

對(duì)此本文認(rèn)為，由于顆粒周圍二次流的影響，顆粒在徑向上的升力分布才出現(xiàn)圖8所示的變化.在低Re數(shù)時(shí)，二次流主要向顆粒下方流動(dòng)且其強(qiáng)度隨著顆?？拷诿娑龃螅@會(huì)給顆粒一個(gè)向下的力，而在圖8中也能明顯地看到在靠近壁面時(shí)升力下降得更快，這說(shuō)明二次流對(duì)顆粒的升力是有影響的.在高Re數(shù)時(shí)，顆粒所受升力在r+=0.4 ～ 0.7 之間相對(duì)平緩且升力系數(shù)較小，這與前文所說(shuō)的二次流強(qiáng)度隨顆粒徑向位置變化相對(duì)應(yīng).而在Re≥500，r+=0.75 時(shí)升力出現(xiàn)回升，本文認(rèn)為這是由于二次流流向變化所引起的，在r+=0.8 時(shí)，二次流主要向顆粒兩側(cè)成對(duì)稱流動(dòng)，這使得其對(duì)升力的影響減弱.

3 結(jié)論

Carlo[2]提出的相對(duì)運(yùn)動(dòng)模型在求解低Re數(shù)下顆粒的慣性聚集是比較成熟的.但在高Re數(shù)下，如果仍對(duì)進(jìn)口給定均勻流，為了消除入口段的影響需要很長(zhǎng)的管道，造成網(wǎng)格數(shù)量大，計(jì)算成本高，因此本文嘗試對(duì)管道進(jìn)、出口采用周期性邊界條件以減小計(jì)算域管長(zhǎng).本文主要對(duì)周期性邊界條件的可行性進(jìn)行了驗(yàn)證并求解了高Re數(shù)下顆粒受力特性，得到了以下結(jié)論：

1）在求解高Re數(shù)流的慣性聚集，周期性邊界條件的使用可以有效地減小管長(zhǎng)，這很大程度上提高了數(shù)值計(jì)算的效率以及經(jīng)濟(jì)性.當(dāng)Re<1 000 時(shí)，a+≤1/9 的顆粒用4D周期就可以計(jì)算出可靠的結(jié)果.對(duì)于粒徑細(xì)小的顆粒，如文中a+=1/17 的顆粒，4D周期可計(jì)算的Re數(shù)高達(dá)1 600.

2）在低Re數(shù)下，顆粒在徑向上的升力呈拋物線分布，且顆粒主要聚集在離壁面較遠(yuǎn)的區(qū)域.隨著Re數(shù)的不斷增大，顆粒的聚集位置向著壁面靠近，且其升力分布出現(xiàn)了較大的波動(dòng)，它將不再呈類拋物線分布，在r+=0.5 ～ 0.7 之間出現(xiàn)了一段升力相對(duì)較小的區(qū)域，而在這個(gè)區(qū)域內(nèi)有出現(xiàn)新聚集點(diǎn)的趨勢(shì).

3）當(dāng)Re≤1 000 時(shí)，本文只發(fā)現(xiàn)了一個(gè)聚集點(diǎn)，新的聚集點(diǎn)并沒有出現(xiàn).但當(dāng)Re>1 000 時(shí)，本文用a+=1/17 的顆粒進(jìn)行計(jì)算得到了兩個(gè)穩(wěn)定的聚集點(diǎn)，這說(shuō)明在高Re數(shù)流中小粒徑的顆粒有可能出現(xiàn)兩個(gè)穩(wěn)定的聚集點(diǎn).

4）顆粒周圍有二次流的產(chǎn)生，其強(qiáng)度隨著Re數(shù)的增大而增大，且隨著顆粒越靠近壁面，二次流的強(qiáng)度也會(huì)增加.在低Re數(shù)時(shí)，二次流主要向著顆粒的下方流動(dòng)；而在高Re數(shù)時(shí)，二次流在顆?？拷ǖ乐行臅r(shí)先是向顆粒下方流動(dòng)，而后隨著顆粒向壁面靠近以及Re數(shù)的增大，其逐漸向顆粒左右兩側(cè)流動(dòng).而受二次流的影響，顆粒所受升力在高Re數(shù)和低Re數(shù)呈現(xiàn)不同的空間分布規(guī)律.