王曉偉，武 越，杜春林，韓 峰，張志軍*，張世偉

（1. 東北大學(xué) 機(jī)械工程與自動化學(xué)院，沈陽 110819； 2. 北京衛(wèi)星環(huán)境工程研究所，北京 100094）

0 引言

在沒有任何初始壓力梯度的稀薄氣體中，非均勻溫度場會誘導(dǎo)氣體流動，且浸沒在氣體中的結(jié)構(gòu)會受到輻射力[1-3]。該現(xiàn)象最早發(fā)現(xiàn)于19 世紀(jì)后期William Crookes 爵士發(fā)明的輻射計(jì)中[4]。該輻射計(jì)由安裝在抽真空的玻璃燈泡內(nèi)部的一組薄葉片組成，葉片一面涂成黑色，另一面為光亮的白色。將該裝置放置于陽光下，葉片開始旋轉(zhuǎn)，旋轉(zhuǎn)方向?yàn)閺暮谏幻嬷赶虬咨幻妗?/p>

隨著航空航天科學(xué)技術(shù)和微機(jī)電系統(tǒng)（MEMS）的發(fā)展，研究人員對應(yīng)用輻射力的興趣日益濃厚。Jones 異想天開地提出基于輻射力的直升機(jī)[5]和微型巡航導(dǎo)彈[6]。基于輻射力的推進(jìn)系統(tǒng)還可用在臨近空間飛行器上[7]，由安裝在航天器表面的多個(gè)方形葉片陣列為航天器提供升力或者阻力補(bǔ)償力。研究發(fā)現(xiàn)，在40～80 km 的高度上，這種推進(jìn)系統(tǒng)葉片上的輻射力能夠達(dá)到25 mN/m2。此外，由于對環(huán)境壓力、溫度以及氣體成分較為敏感，輻射力也被用于MEMS 傳感器中的壓力測量[8]和氣體檢測[9-20]。

最近，Strongrich 等[9-10]提出一種基于輻射力驅(qū)動的設(shè)備（Microelectromechanical In-plane Knudsen Radiometric Actuator, MIKRA），可用于氣體檢測。研究表明，物理和幾何參數(shù)顯著影響該設(shè)備的性能[21-22]：當(dāng)冷熱微梁溫差從10 K 增大到100 K 時(shí)，最大輻射力將增大7 倍以上；輻射力隨環(huán)境壓力的增大呈現(xiàn)出“鐘形”分布規(guī)律；微梁寬度的增大將使設(shè)備整體所受輻射力減小，且最大輻射力在更小的環(huán)境壓力下取得；微梁高度增大3 倍時(shí)，最大輻射力將增大3.66 倍；當(dāng)梁的高度與寬度的比值增加9 倍時(shí)，最大輻射力能增大5.33 倍。需要指出的是，在這些研究中氣體與所有梁表面之間的相互作用規(guī)律都是基于完全漫反射假設(shè)；然而，在實(shí)際應(yīng)用中這種假設(shè)并非總是成立的。有關(guān)研究表明，氣-面相互作用規(guī)律（表面調(diào)節(jié)系數(shù)）對流場結(jié)構(gòu)、傳熱特性、壓力分布以及輻射力都有顯著的影響[23-25]。Zhou 等的研究表明[23]，法向能量調(diào)節(jié)系數(shù)（NEAC）和切向動量調(diào)節(jié)系數(shù)（TMAC）對微通道內(nèi)的壓力分布和流場結(jié)構(gòu)會產(chǎn)生很大的影響；并且，相比于NEAC，改變TMAC 對流動特性的影響更為顯著。此外，Kalempa 和Wang 等[24-25]研究了氣-面相互作用規(guī)律對稀薄氣體中球體和微梁周圍流體特性的影響，發(fā)現(xiàn)對于輻射力，改變NEAC 和TMAC 會產(chǎn)生相反的作用。值得注意的是，在文獻(xiàn)[25]的研究中，微梁（一對冷熱梁）上的NEAC 和TMAC 是同步變化的，并沒有考慮其中一根梁上表面調(diào)節(jié)系數(shù)單獨(dú)變化的影響情況。

為了補(bǔ)充氣-面相互作用規(guī)律對輻射力影響的理論研究，本文采用直接模擬蒙特卡羅（DSMC）方法研究微梁周圍的稀薄氣體流動特性，主要考查過渡流態(tài)下單梁表面調(diào)節(jié)系數(shù)變化對流場結(jié)構(gòu)、傳熱特性、壓力分布以及輻射力等的影響，以期為相關(guān)器件的設(shè)計(jì)制造和應(yīng)用提供理論參考與指導(dǎo)。

1 氣-面相互作用模型

常見的氣-面相互作用模型包括鏡面反射、漫反射、Maxwell 以及CL（Cercignani-Lampis）模型。對于不同的氣-面相互作用模型，氣體分子撞擊壁面并發(fā)生反射的速度分布規(guī)律[26]如圖1 所示。鏡面反射模型是最簡單的氣-面相互作用模型，它假設(shè)撞擊到表面的氣體分子類似于完全彈性的球體一樣被反射，即，氣體分子的法向速度分量反向而切向速度分量保持不變，如圖1(a)所示，這意味著在這種氣-面相互作用過程中，氣體分子與壁面之間不發(fā)生熱傳遞和摩擦行為。對于漫反射模型，撞擊到壁面的氣體分子首先與表面進(jìn)行能量交換，使氣體分子溫度和壁面溫度達(dá)到熱平衡；然后根據(jù)由壁面溫度確定的半范圍Maxwell 分布，隨機(jī)分配反射分子的速度，如圖1(b)所示。前面這2 種氣-面相互作用模型都不太切合實(shí)際，因此出現(xiàn)了Maxwell模型。實(shí)際上，Maxwell 模型是漫反射與鏡面反射模型的組合[27]，在模型中引入調(diào)節(jié)系數(shù)α（0≤α≤1）來指定入射分子在壁面上的反射形式，即，入射分子的α部分為漫反射，剩余的(1-α)部分為鏡面反射。相比于Maxwell 模型，CL 模型能更真實(shí)地反映物理現(xiàn)實(shí)中反射分子的速度和能量分布[28]，如圖1(c)所示。該模型包括2 個(gè)表面調(diào)節(jié)系數(shù)（NEAC和TMAC），二者相互獨(dú)立，對于一般的光滑表面，這2 個(gè)調(diào)節(jié)系數(shù)的數(shù)值均在0～1 之間變化。本研究的所有DSMC 模擬中，氣-面相互作用模型均選用CL 模型。

圖1 不同氣-面相互作用模型中反射分子的速度分布示意Fig. 1 Schematics diagram of velocity distributions of reflected molecules based on different gas-surface interaction models

2 問題陳述與數(shù)值方法

由于MIKRA 設(shè)備中冷熱梁的長度（1.4 mm）遠(yuǎn)大于其橫截面尺寸（50 μm×50 μm），所以只需進(jìn)行二維模擬以近似表征實(shí)際的氣體流動情況，模擬域如圖2 所示。左右兩表面為對稱邊界，上方是自由流邊界以表征H2-N2混合氣體的流入或流出，下方為基底壁面邊界。冷熱梁懸浮于距基底上方4 μm處。冷梁和熱梁具有相同的尺寸，且在模擬過程中假設(shè)兩者的間距g保持不變，為20 μm。此間距也用于計(jì)算Knudsen 數(shù)（Kn）。

圖2 模擬域和邊界條件Fig. 2 Simulation domain and boundary conditions

DSMC 方法基于粒子求解Boltzmann 方程[29]，使用大量的模擬粒子再現(xiàn)實(shí)際氣體分子的運(yùn)動情況，跟蹤并記錄這些粒子的位置、速度和內(nèi)能等信息，最后對模擬粒子的微觀值進(jìn)行平均從而獲得宏觀量。該方法為了降低統(tǒng)計(jì)誤差，需要較長的模擬時(shí)間。本文使用開源的DSMC 求解器dsmcFoam-Plus[30]來模擬微梁周圍的氣體流動；采用無時(shí)間計(jì)數(shù)器（NTC）方案選取粒子碰撞對；另外，因?yàn)槟M域氣體介質(zhì)是50%H2+50%N2的雙原子混合氣體[13,18,25]，所以采用變徑軟球（VSS）模型和Larsen-Borgnakke 模型分別模擬粒子間的碰撞和粒子本身動能與內(nèi)能的交換過程。

所有模擬均在環(huán)境壓力P=387 Pa（Kn=0.74）的氣氛中進(jìn)行；冷梁和均勻加熱的熱梁的溫度分別為Tc=304 K 和Th=341.5 K。值得一提的是，這2 個(gè)溫度是通過在純N2氣氛中將相應(yīng)壓力條件下測得的實(shí)驗(yàn)數(shù)據(jù)[9-10]取平均值得到的。根據(jù)這些參數(shù)，使用Δx=Δy=2 μm 的網(wǎng)格均勻地劃分模擬域[10,13,18,25]，每個(gè)網(wǎng)格中的模擬粒子數(shù)不少于20 個(gè)。模擬計(jì)算的時(shí)間步長Δt=1 ns，總模擬時(shí)間不短于5 ms[10,13,18,25]；另外，對于表面調(diào)節(jié)系數(shù) NEAC（αn）和TMAC（σt）均獨(dú)立設(shè)為0.2、0.4、0.6、0.8 和1。本研究中只單獨(dú)改變冷梁的表面調(diào)節(jié)系數(shù)，而熱梁的表面調(diào)節(jié)系數(shù)恒為1（對于冷、熱梁同步變化的情況，請參見文獻(xiàn)[25]）。在物理實(shí)際中，可通過改變制造冷梁和熱梁所用的材料來實(shí)現(xiàn)其表面調(diào)節(jié)系數(shù)的調(diào)整。本文為區(qū)分冷梁和熱梁的表面調(diào)節(jié)系數(shù)，分別以上標(biāo)Cold 和Hot 予以標(biāo)記。

3 結(jié)果與討論

3.1 流場結(jié)構(gòu)

圖3 不同NEAC 下的速度流線和溫度分布云圖Fig. 3 Velocity streamlines and temperature contours for different NEACs (In each of the sub-figures, the left half shows the results of αnCold=αn, αnHot=1, the right half shows the results extracted from Ref. [25], αnCold=αnHot=αn)

TMAC 變化對微梁附近流場的影響如圖4 所示，同樣在每幅圖中均給出了本文和文獻(xiàn)[25]的模擬結(jié)果。不難看出，相比于NEAC 對流場結(jié)構(gòu)的影響，切向TMAC 的影響力較弱，且在僅改變冷梁TMAC 以及冷梁和熱梁TMAC 同步變化的情況中均是如此。這一結(jié)果對微梁制作材料的選取具有一定的指導(dǎo)意義。同樣由于熱梁TMAC 的降低導(dǎo)致氣體與壁面交換的能量減少，使得冷梁和熱梁TMAC 同步改變所引起的流場結(jié)構(gòu)的變化比僅改變冷梁TMAC 的情況更為明顯，且流場溫度更高。

圖4 不同TMAC 下的速度流線和溫度分布云圖Fig. 4 Velocity streamlines and temperature contours for different TMACs (In each of the sub-figures, the left half shows the results of =σt, =1, the right half shows the results extracted from Ref. [25], ==σt)

3.2 傳熱特性

不同NEAC 對冷梁各表面?zhèn)鳠崃康挠绊懭鐖D5所示，每幅圖中虛線數(shù)據(jù)為文獻(xiàn)[25]的模擬結(jié)果，下同。傳熱量數(shù)值為正表示該表面處于吸熱狀態(tài)，為負(fù)則表示其處于放熱狀態(tài)。可以看到，受到高溫?zé)崃旱挠绊懀淞旱纳媳砻婧陀冶砻嫣幱谖鼰釥顟B(tài)；對于下表面，由于其溫度比底部基底的更高且距離較近，故處于放熱狀態(tài)；左表面處于放熱狀態(tài)的原因是，其遠(yuǎn)離熱梁且表面溫度比其周圍氣體的溫度高。正如預(yù)期的那樣，NEAC 的減小導(dǎo)致表面?zhèn)鳠崃康慕^對值減小，即，表面的吸熱或放熱能力減弱。另外，對于僅改變冷梁的NEAC 的情況，表面的吸熱或放熱能力比冷梁和熱梁的NEAC 同步變化時(shí)更強(qiáng)。由于冷梁右表面距離熱梁較近，這種現(xiàn)象尤為明顯，如圖5(d)所示。

圖5 不同NEAC 對冷梁各表面?zhèn)鳠崃康挠绊慒ig. 5 Heat transfer distributions along surfaces of cold beam with different NEACs (In each of the sub-figures, dashed lines show the results extracted from Ref. [25])

不同TMAC 對冷梁各表面?zhèn)鳠崃康挠绊懭鐖D6所示。可以看出，各表面的吸熱或放熱狀態(tài)以及傳熱分布規(guī)律與改變NEAC 時(shí)的一致。從傳熱量上來看，改變TMAC 比改變NEAC 能夠獲得更大的絕對值。

圖6 不同TMAC 對冷梁各表面?zhèn)鳠崃康挠绊慒ig. 6 Heat transfer distributions along surfaces of cold beam with different TMACs (In each of the sub-figures, dashed lines show the results extracted from Ref. [25])

這對于物理應(yīng)用中的熱管理具有一定的參考價(jià)值，即，針對材料不同的微梁需要配備不同換熱能力的散熱系統(tǒng)。另外，僅改變冷梁的TMAC 獲得的結(jié)果與冷梁和熱梁的TMAC 同步變化所得到的結(jié)果之間的差異在TMAC 較大的情況下并不明顯，只在TMAC 較小時(shí)才能凸顯。

3.3 壓力分布特性

圖7(a)和圖7(b)分別為NEAC 和TMAC 變化對冷梁左右表面壓力分布的影響。

圖7 不同表面調(diào)節(jié)系數(shù)下冷梁左右表面的壓力分布Fig. 7 Pressure distributions along left and right surfaces of the cold beam with different surface accommodation coefficients (In each of the sub-figures, dashed lines show the results extracted from Ref. [25])

由圖7(a)可以看到，僅改變冷梁的NEAC 在引起梁表面壓力變化方面與冷梁和熱梁的NEAC 同步變化的情況顯示出完全不同的結(jié)果。更確切地說，如果僅僅改變冷梁的NEAC，冷梁右表面壓力始終大于左表面壓力，且隨著NEAC 的減小，右表面壓力明顯增大。這意味著冷梁與熱梁兩者之間相互排斥。而對于冷梁和熱梁NEAC 同步變化的情況，則會出現(xiàn)冷梁左表面壓力大于右表面壓力的現(xiàn)象，且隨著NEAC 的降低右表面壓力明顯減小，這種現(xiàn)象更為顯著[25]。這意味著冷梁與熱梁兩者之間在NEAC 較大時(shí)是相互排斥的，而當(dāng)NEAC 較小時(shí)則表現(xiàn)為相互吸引。因此，在物理實(shí)際中可通過選用不同的材料來實(shí)現(xiàn)對微梁移動方向的控制。

由圖7(b)可以看到，無論是僅改變冷梁的TMAC還是同步改變冷梁和熱梁的TMAC，冷梁右表面壓力均始終大于左表面壓力，即，冷梁與熱梁兩者之間相互排斥。此外，隨著TMAC 的減小，右表面壓力明顯增大，且同步降低冷梁和熱梁的TMAC 在增大右表面壓力方面比僅降低冷梁的TMAC 的情況表現(xiàn)出更顯著的效果。

3.4 輻射力

在DSMC 中，施加在冷梁上的輻射力F可通過撞擊在冷梁表面的入射粒子和反射粒子之間的動量變化來計(jì)算，

式中：Nequ為一個(gè)DSMC 模擬粒子所代表的真實(shí)氣體分子數(shù)量；N為與表面相互作用的總粒子數(shù)；m為粒子的質(zhì)量；c為粒子的速度，上標(biāo)i 和r 分別表示入射粒子和反射粒子。

圖8(a)和圖8(b)分別為NEAC 和TMAC 變化對冷梁所受輻射力的影響，其中虛線為文獻(xiàn)[25]的結(jié)果。

由圖8(a)可以看出，對于僅改變冷梁的NEAC與冷梁和熱梁的NEAC 同步變化這2 種情況，冷梁所受輻射力對NEAC 的依賴性顯示出完全不同的情形。更確切地說，如果同步減小冷梁和熱梁的NEAC，輻射力隨著NEAC 的減小呈非線性減小，當(dāng)αn＜0.2 之后出現(xiàn)輻射力方向反轉(zhuǎn)的現(xiàn)象，且其絕對值具有隨著NEAC 的減小而增大的趨勢（輻射力曲線的斜率隨著NEAC 的減小而增大）；而對于僅減小冷梁的NEAC 的情況，輻射力隨著NEAC的減小呈近似線性增大。

由圖8(b)可以看出，在僅改變冷梁表面TMAC的情況中，輻射力隨著TMAC 的減小呈非線性減?。欢鴮τ诶淞汉蜔崃旱腡MAC 同步變化的情況，輻射力卻隨著TMAC 的減小呈非線性增大。同時(shí)，我們注意到，在這2 種情況中，隨著TMAC 的減小輻射力曲線的斜率（絕對值）均增大。

圖8 不同表面調(diào)節(jié)系數(shù)下冷梁所受輻射力Fig. 8 Radiometric force of the cold beam with different surface accommodation coefficients (In each of the sub-figures, dashed lines show the results extracted from Ref. [25])

4 結(jié)束語

本文采用DSMC 方法并基于CL 氣-面相互作用模型研究法向能量調(diào)節(jié)系數(shù)（NEAC）和切向動量調(diào)節(jié)系數(shù)（TMAC）對一組微梁周圍稀薄氣體流動的影響。結(jié)果表明，表面調(diào)節(jié)系數(shù)對流場結(jié)構(gòu)、傳熱特性、壓力分布和輻射力的影響是顯著的。為了在實(shí)際應(yīng)用中盡可能地獲得最大的輻射力，在設(shè)計(jì)和制造MEMS 器件時(shí)，如果冷梁和熱梁使用相同的材料，則應(yīng)該選取NEAC 較大而TMAC 較小的材料；如果冷梁和熱梁使用不同的材料，則冷梁應(yīng)該選具有較小NEAC 和較大TMAC 的材料，而熱梁應(yīng)該選NEAC 和TMAC 均較大（最好接近于1）的材料。