祁影霞，張 燦，車(chē)閆瑾，王禹賀，張 華

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

0 引言

隨著科學(xué)技術(shù)的發(fā)展，航天航空技術(shù)、超導(dǎo)技術(shù)、低溫電子學(xué)、低溫醫(yī)學(xué)等對(duì)低溫制冷技術(shù)的需求與日俱增[1-2]。脈管制冷機(jī)作為一種小型制冷機(jī)，在低溫制冷領(lǐng)域有著獨(dú)特的優(yōu)勢(shì)，具有無(wú)運(yùn)動(dòng)部件、結(jié)構(gòu)簡(jiǎn)單、機(jī)械振動(dòng)小、無(wú)噪音等優(yōu)點(diǎn)[3]，在上述領(lǐng)域應(yīng)用廣泛。

脈管是脈管制冷機(jī)的關(guān)鍵部位連接著冷端換熱器和熱端換熱器，氣體在其中交變流動(dòng)來(lái)實(shí)現(xiàn)熱量和冷量的傳遞。雖然過(guò)程簡(jiǎn)單，但內(nèi)部機(jī)理至今也沒(méi)有統(tǒng)一的定論。Gifford等[4]研究了基本型脈管制冷機(jī)，并提出了“泵熱理論”?；拘兔}管缺乏調(diào)相能力，效率低下，Radebaugh等[5]在此基礎(chǔ)上研究了調(diào)相型脈管制冷機(jī)中焓流和壓力、速度之間的關(guān)系并提出了“焓流調(diào)相理論”。后來(lái)又出現(xiàn)了如向量分析法[6]、熱聲理論[7-8]、特征線法[9]等新的理論來(lái)解釋脈管的制冷機(jī)理，使得脈管制冷機(jī)得到了很好的發(fā)展。

近年來(lái)，脈管制冷機(jī)的基礎(chǔ)理論不斷發(fā)展。De Waele等[10-11]在假設(shè)氣庫(kù)無(wú)限大，且壓力為正弦變化時(shí)，獲得了小孔型和雙向進(jìn)氣型脈管制冷機(jī)的冷熱端的動(dòng)態(tài)溫度變化。Jafarian等[12]運(yùn)用數(shù)值模擬的手段，通過(guò)控制容積方法預(yù)測(cè)了大冷量脈管制冷機(jī)的流體動(dòng)力學(xué)和熱力學(xué)特性，并且表明脈管制冷機(jī)的運(yùn)行參數(shù)是影響其制冷效率的主要因素之一。Sexl等[13]研究了脈管制冷機(jī)脈管內(nèi)部氣體振蕩流動(dòng)的壓力梯度。

脈管制冷機(jī)在運(yùn)行時(shí)涉及熱力學(xué)、流體動(dòng)力學(xué)過(guò)程，同時(shí)伴有復(fù)雜的流動(dòng)過(guò)程，因此其機(jī)理研究至今仍不完善。為了有效地提升脈管制冷機(jī)的制冷效率，需要進(jìn)一步研究其機(jī)理及影響因素?；诜肿觿?dòng)力學(xué)方法研究了脈管形狀、粗糙度以及壓力對(duì)脈管內(nèi)部溫度場(chǎng)的影響，并給出了脈管內(nèi)部交變流動(dòng)過(guò)程溫度、壓力、密度、速度隨時(shí)間的變化趨勢(shì)。

1 模擬方法及物理模型

1.1 分子動(dòng)力學(xué)

分子動(dòng)力學(xué)模擬方法（Molecular Dynamics Simulation Method）是當(dāng)前學(xué)者們用來(lái)計(jì)算規(guī)模龐大且復(fù)雜，但系統(tǒng)本身尺度微小的一種方法。隨著力學(xué)迅速發(fā)展，研究人員建立了適用于生化分子體系、聚合物、金屬與非金屬材料等不同粒子領(lǐng)域的力場(chǎng)，這些力場(chǎng)的應(yīng)用使得計(jì)算結(jié)果精準(zhǔn)度大幅提升。分子動(dòng)力學(xué)模擬方法是將這些力場(chǎng)和牛頓力學(xué)合起來(lái)的一種新型的計(jì)算方法。

在數(shù)量龐大的原子運(yùn)動(dòng)系統(tǒng)中，系統(tǒng)能量是系統(tǒng)中原子的動(dòng)能與系統(tǒng)總勢(shì)能之和。其表達(dá)式為：

式中：U為原子間總勢(shì)能；UVDW為非鍵結(jié)合范德瓦耳斯作用力；Uint為分子內(nèi)部勢(shì)能。

根據(jù)經(jīng)典力學(xué)可知，系統(tǒng)中原子i所受的力為勢(shì)能梯度：

由牛頓運(yùn)動(dòng)定律可得i原子的加速度為：

4、創(chuàng)新考核機(jī)制，保證任務(wù)的落實(shí)到位。對(duì)于各自決策和負(fù)責(zé)的事權(quán)，各自權(quán)責(zé)內(nèi)恰，管理相對(duì)容易。對(duì)于中央決策、地方執(zhí)行的事權(quán)，與自主權(quán)下放相適應(yīng)，需要?jiǎng)?chuàng)新監(jiān)督機(jī)制，既要保障地方自主權(quán)，又要保障地方執(zhí)行與中央大政方針和戰(zhàn)略目標(biāo)一致。

開(kāi)始計(jì)算時(shí)，固定立方體中的分子數(shù)和工質(zhì)氣體原子的初始位置、初始速度是確定的。系統(tǒng)中全部原子運(yùn)動(dòng)動(dòng)能總和為：

式中：KE為系統(tǒng)的總動(dòng)能；N為總原子數(shù)；kB為玻爾茲曼常數(shù)；T為熱力學(xué)溫度。

原子運(yùn)動(dòng)的初始速度可由式（2）計(jì)算出，得到原子運(yùn)動(dòng)的初始位置和初始速度后，之后每一步會(huì)產(chǎn)生當(dāng)前速度和當(dāng)前位置，通過(guò)當(dāng)前速度可計(jì)算固定格子的溫度。計(jì)算公式為：

系統(tǒng)運(yùn)行過(guò)程中，在不同時(shí)間段對(duì)格子內(nèi)工質(zhì)氣體原子數(shù)目進(jìn)行統(tǒng)計(jì)，可以得到不同時(shí)間段內(nèi)軸向壓力分布。聯(lián)合上面的公式，可以計(jì)算出各個(gè)時(shí)間段內(nèi)脈管內(nèi)部軸向速度和溫度場(chǎng)的分布。

1.2 物理模型

建立三種不同壁面的脈管通道，模擬帶有氣庫(kù)的脈管制冷機(jī)在膨脹時(shí)系統(tǒng)內(nèi)部工質(zhì)氣體的流動(dòng)過(guò)程。圖1為幾種通道的模擬物理模型，分別為平直光滑通道、彎曲通道、粗糙通道。所有通道高度為200 nm，長(zhǎng)度為730 nm，總長(zhǎng)（含熱端）1 000 nm。其中，彎曲通道曲率半徑為100 nm；粗糙通道有四種參數(shù)，分別如表1所列。左側(cè)氣體壓力分別取4.05×105Pa、6.08×105Pa、8.11×105Pa，工質(zhì)為Ar氣，模擬過(guò)程中，壁面保持固定，氣體與壁面之間沒(méi)有熱交換。首先用正則（NVT）系綜（系統(tǒng)中粒子數(shù)N，體積V，溫度T均為定值），使系統(tǒng)平均溫度保持在298 K。之后采用微正則（NVE）系綜（系統(tǒng)中粒子數(shù)N，體積V，總能E均為定值），使系統(tǒng)與外界沒(méi)有能量交換。當(dāng)系統(tǒng)運(yùn)行時(shí)間在3 000 ps以后，通道內(nèi)氣體基本達(dá)到平衡，如圖2所示。

表1 粗糙通道粗糙元參數(shù)Table1 Roughness parameter of rough channels

圖2為Ar氣在運(yùn)行到3 000 ps時(shí)，原子位置分布圖?？梢钥闯?，直通道的原子分布沿軸向已經(jīng)較均勻，而彎曲通道次之，粗糙通道更次之。

圖1 三種通道模擬物理模型圖Fig.1 Simulation physical model of three channels

圖2 模擬結(jié)束時(shí)原子位置分布圖Fig.2 Atomic position distribution status at the end of simulation

2 分子動(dòng)力學(xué)仿真結(jié)果

2.1 軸向速度分布

圖3、圖4和圖5分別是直通道、彎曲通道、粗糙通道原子軸向運(yùn)動(dòng)速度分布圖。可以看出，通道內(nèi)氣體流動(dòng)速度在1 260 ps時(shí)達(dá)到最大，這是由于Ar原子剛開(kāi)始運(yùn)動(dòng)時(shí)，通道兩端壓差較大，壓差驅(qū)動(dòng)其快速運(yùn)動(dòng)。彎曲通道和粗糙通道軸向速度分布達(dá)到最大值在1 380 ps，可知軸向速度分布受壁面的光滑度和彎曲程度影響。在運(yùn)動(dòng)一段時(shí)間后，原子基本分布在整個(gè)通道內(nèi)，兩端壓差減小，原子運(yùn)動(dòng)速度降低。從圖中看出，最終在3 180 ps時(shí)原子分布均勻且速度接近0，說(shuō)明系統(tǒng)已達(dá)到穩(wěn)定狀態(tài)。

圖3 直通道原子軸向運(yùn)動(dòng)速度分布曲線Fig.3 Axial velocity distribution of straight channel atom

圖4 彎曲通道原子軸向運(yùn)動(dòng)速度分布曲線Fig.4 Axial velocity distribution of bent channel atom

圖5 粗糙通道原子軸向運(yùn)動(dòng)速度分布曲線Fig.5 Axial velocity distribution of rough channel atom

2.2 通道形狀對(duì)溫度場(chǎng)的影響

圖6為三種不同形式的通道，在通道內(nèi)冷端溫度達(dá)到最低時(shí)的軸向溫度分布。圖中虛線左側(cè)為通道，右側(cè)為通道外（熱端）?？梢钥闯?，在冷、熱端之間存在一個(gè)溫差，這是因?yàn)橥ǖ纼?nèi)氣體初始?jí)毫Ω?，通道外氣體初始?jí)毫Φ?，因而通道?nèi)氣體向外膨脹溫度降低，通道外氣體被壓縮因而溫度上升。這就是脈管內(nèi)氣體絕熱放氣制冷的熱力學(xué)機(jī)理。同時(shí)，在冷端和熱端都有一個(gè)溫度分布相對(duì)平緩的部分。

圖6 三種類(lèi)型通道沿軸向溫度分布對(duì)比圖Fig.6 Comparison of axial temperature distribution of three channels

從圖6還可以看出，和直通道相比，彎曲和粗糙通道冷端溫度大幅提升，且整體軸向溫度分布上移，這是由于在彎曲和粗糙通道內(nèi)，氣體與壁面的摩擦轉(zhuǎn)換為摩擦熱而使氣體溫度上升所致。粗糙度和彎曲度的影響程度幾乎一樣。直通道冷端溫度由初始值的298 K降低到138 K，彎曲通道、粗糙通道分別降低到190 K、215 K。可見(jiàn)彎曲通道、粗糙通道對(duì)絕熱放氣的制冷溫度影響較大。圖7是與圖6相對(duì)應(yīng)的軸向壓力分布圖。可以看出，直通道的壓力分布較平緩，彎曲通道、粗糙通道的壓力分布呈上升趨勢(shì)，直通道的最大壓差大約是20 kPa，彎曲、粗糙通道的壓差分別為110 kPa、140 kPa，說(shuō)明彎曲、粗糙通道的流動(dòng)阻力較大，因而壓力損失也較大。

圖7 三種類(lèi)型通道沿軸向壓力分布對(duì)比圖Fig.7 Comparison of axial pressure distribution of three channels

2.3 壓力對(duì)制冷溫度的影響

圖8 三種類(lèi)型通道冷端最低溫度隨系統(tǒng)壓力的變化圖Fig.8 The change of the minimum temperature of the cold end with the system pressure in three channels

2.4 粗糙度對(duì)溫度場(chǎng)的影響

初始?jí)毫υ?.08×105Pa時(shí)，運(yùn)行到3 000 ps后，三種類(lèi)型通道及不同粗糙度通道冷端達(dá)到的最低溫度，如圖9所示?？梢钥闯鲭S著粗糙度的增加，冷端溫度也相應(yīng)提高，相對(duì)于平直光滑通道，case4的冷端溫度最大，提高了82 K。同時(shí)，也可以看出，粗糙度的影響比彎曲度對(duì)冷端溫度的影響要大。由此可見(jiàn)，通道壁面光潔度對(duì)脈管制冷的影響還是非常大的。

圖9 不同粗糙度通道冷端最低溫度比較圖Fig.9 Comparison of the minimum temperature of the cold end in different channels

2.5 通道各參數(shù)隨時(shí)間的變化

比較在冷端（150 nm）、熱端（790 nm）處，三種類(lèi)型通道內(nèi)溫度、壓力、密度、速度隨時(shí)間的變化。

圖10和圖11分別為150 nm、790 nm處，三種類(lèi)型通道內(nèi)溫度隨時(shí)間的變化。在圖10中，對(duì)于直通道，冷端溫度在2 500～2 800 ps時(shí)達(dá)到最低，之后有明顯回升；彎曲通道在1 500～2 800 ps、粗糙通道在1 700～2 800 ps時(shí)達(dá)到較低穩(wěn)定波動(dòng)狀態(tài)。直通道的冷端溫度要比彎曲和粗糙通道的溫度低很多，而彎曲通道又比粗糙通道的溫度低。在圖11中，位置為790 nm處，該位置處于通道出口處，對(duì)應(yīng)熱端溫度。直通道在2 400 ps達(dá)到最高，之后開(kāi)始下降；而彎曲、粗糙通道在2 000 ps之后達(dá)到較高穩(wěn)定波動(dòng)狀態(tài)。也就是說(shuō)，彎曲、粗糙通道的溫度梯度在2 000 ps之后基本不變了，彎曲和粗糙通道的熱端溫度要遠(yuǎn)高于直通道。

圖10 三種類(lèi)型通道在150 nm處溫度隨時(shí)間的變化圖Fig.10 The temperature varies of three channels with time at 150 nm

圖11 三種類(lèi)型通道在790 nm處溫度隨時(shí)間的變化圖Fig.11 Temperature varies of three channels with time at 790 nm

圖12～15分別為150 nm、790 nm處，三種類(lèi)型的通道內(nèi)壓力、密度隨時(shí)間的變化。

從圖12、13可知，在150 nm的位置，三種通道的密度和壓力隨著時(shí)間的推移都呈下降的趨勢(shì)，這是由于靠近冷端的位置，0 ps時(shí)，Ar氣密度較大、壓力較大；在壓差作用下，隨著Ar氣向熱端的流動(dòng)，靠近冷端處密度下降，壓力降低。對(duì)比三種通道，直通道的密度、壓力整體低于彎曲、粗糙通道，這是由于直通道壁面光滑，原子運(yùn)動(dòng)較快，密度和壓力下降也較快；彎曲、粗糙通道由于曲率半徑和粗糙度的存在，產(chǎn)生阻力，阻力抵消掉了一部分壓差驅(qū)動(dòng)力，所以彎曲、粗糙通道內(nèi)原子運(yùn)動(dòng)速度較慢，密度、壓力下降也較慢，導(dǎo)致密度和壓力明顯高于直通道的壓力。

圖12 三種類(lèi)型通道在150 nm處軸向密度分布隨時(shí)間的變化圖Fig.12 Axial density distribution varies of the three channels with time at 150 nm

圖13 三種類(lèi)型通道在150 nm處軸向壓力分布隨時(shí)間的變化圖Fig.13 Axial pressure distribution varies of the three channels with time at 150 nm

由圖14、15可知，在790 nm的位置，直通道在900 ps時(shí)開(kāi)始出現(xiàn)密度和壓力波動(dòng)，彎曲通道和粗糙通道在1 400 ps時(shí)才出現(xiàn)壓力和密度波動(dòng)。這是由于彎曲通道和粗糙通道沿程阻力較大，原子運(yùn)動(dòng)進(jìn)程較緩慢，因而到達(dá)790 nm位置所花的時(shí)間也較長(zhǎng)，密度和壓力上升的也較慢。對(duì)比三種通道，直通道的密度、壓力整體高于彎曲、粗糙通道，這是由于直通道內(nèi)原子運(yùn)動(dòng)的速度更快，因而密度和壓力也上升的更快。

對(duì)比圖12～15可知，粗糙通道的阻力比彎曲通道的要大。

圖14 三種類(lèi)型通道在790 nm處軸向密度分布隨時(shí)間的變化圖Fig.14 Axial density distribution varies of the three channels with time at 790 nm

圖15 三種類(lèi)型通道在790 nm處軸向壓力分布隨時(shí)間的變化圖Fig.15 Axial pressure distribution varies of the three channels with time at 790 nm

圖16和圖17為三種類(lèi)型的通道在150 nm、790 nm處，軸向速度隨時(shí)間的變化。由圖16可知，在150 nm處，軸向速度隨時(shí)間逐漸下降，在3 000 ps時(shí)，基本處于穩(wěn)定，且接近0。彎曲、粗糙通道內(nèi)原子的運(yùn)動(dòng)速度整體低于直通道，說(shuō)明彎曲和粗糙通道對(duì)原子運(yùn)動(dòng)的阻力較大，原子會(huì)比在直通道中更大概率與壁面發(fā)生碰撞，導(dǎo)致速度下降。由圖17可知，在790 nm處，隨著運(yùn)行時(shí)間的增加，軸向速度也逐漸下降，在3 000 ps時(shí)，趨于穩(wěn)定，速度基本接近0，說(shuō)明系統(tǒng)此時(shí)達(dá)到平衡。在790 nm處，三種通道的軸向速度基本一致，說(shuō)明通道阻力對(duì)出口速度影響不大。

圖16 三種類(lèi)型通道在150 nm處軸向速度分布隨時(shí)間的變化圖Fig.16 Axial velocity distribution varies of the three channels with time at 150 nm

圖17 三種類(lèi)型通道在790 nm處軸向速度分布隨時(shí)間的變化圖Fig.17 Axial velocity distribution varies of the three channels with time at 790 nm

3 結(jié)論

利用分子動(dòng)力學(xué)模擬了脈管形狀、粗糙度以及壓力對(duì)脈管內(nèi)部溫度場(chǎng)的影響，給出了脈管內(nèi)部交變流動(dòng)溫度、壓力、密度、速度隨時(shí)間的變化趨勢(shì)。得到了以下結(jié)論：

（1）平均壓力的提高，冷端溫度呈現(xiàn)下降趨勢(shì)。故應(yīng)提高系統(tǒng)壓力，進(jìn)而達(dá)到降低冷端溫度的目的。

（2）彎曲、粗糙通道由于曲率半徑和粗糙元的存在，造成運(yùn)動(dòng)阻力，氣體無(wú)法充分膨脹，冷端壓力不能充分降低，制冷溫度有所上升；導(dǎo)致摩擦，造成氣體整體溫度上升。

（3）粗糙度對(duì)冷端溫度的影響比彎曲曲率半徑更大。所以在脈管加工過(guò)程中，要保證其表面盡可能光潔。

（4）冷端溫度隨運(yùn)行時(shí)間逐漸下降，直至穩(wěn)定；熱端溫度隨運(yùn)行時(shí)間逐漸上升，直至穩(wěn)定。直通道的冷、熱端溫度均顯著低于彎曲、粗糙通道。冷端壓力、密度隨運(yùn)行時(shí)間逐漸減小，直至達(dá)到平衡，直通道的冷端壓力、密度低于彎曲、粗糙通道；熱端壓力、密度隨運(yùn)行時(shí)間逐漸上升，直通道的熱端壓力、密度要顯著高于彎曲、粗糙通道。