吳海龍，聶萬勝，鄭 直，劉 瑜

（航天工程大學，北京，101416）

0 引 言

近年來，液氧煤油火箭發(fā)動機由于其無污染、可靠性高、成本低、可重復使用等優(yōu)點得到各航天大國的青睞。為了提高液氧煤油火箭發(fā)動機的性能，有必要對其工作過程進行詳細的研究。煤油作為推進劑，一般要經歷噴注、霧化、蒸發(fā)、混合、燃燒等過程，其中煤油液滴的蒸發(fā)與燃燒過程尤為重要，直接影響著發(fā)動機的性能。由于液氧煤油火箭發(fā)動機燃燒室的溫度和壓力高于煤油推進劑的臨界溫度和臨界壓力[1，2]，因此有必要對超臨界環(huán)境下的煤油液滴的蒸發(fā)和燃燒過程進行研究。

由于煤油是一種十分復雜的混合物，仿真計算難以實現對這種復雜混合物各項性質的仿真，因此對于煤油液滴在超臨界環(huán)境下的蒸發(fā)的研究目前以實驗為主，主要研究了超臨界環(huán)境下各個環(huán)境參數對液滴蒸發(fā)過程的影響[3，4]及活塞驅動器技術用于超臨界液滴蒸發(fā)與燃燒實驗[5]；在數值仿真方面，多以單組分液滴的超臨界仿真計算為主，Umemura等[6]研究了超臨界氮氣環(huán)境下的正丁烷液滴的蒸發(fā)規(guī)律；Zhang Hongtao等[7]研究了亞臨界和超臨界無重力環(huán)境下處于運動狀態(tài)的正庚烷液滴的蒸發(fā)常數的變化規(guī)律；Aggarwal等[8]研究了超臨界環(huán)境下液滴跨臨界蒸發(fā)特性，得到達到臨界點所需的最小壓力的規(guī)律；Balaji等[9]研究了正十二烷液滴在亞臨界和超臨界環(huán)境下液滴生存時間與溫度和壓力的關系；王宏楠等[10]研究了環(huán)境壓力、環(huán)境溫度及液滴初始溫度對處于超臨界氮氣環(huán)境下的壬烷液滴的蒸發(fā)特性的影響；Singh等[11]對比研究了液滴高壓和低壓蒸發(fā)模型，得到了壓力和溫度對正庚烷液滴蒸發(fā)的影響；Xiao Guowei等[12]用分子動力學的方法研究了正十二烷和氮氣的超臨界混合過程；鄭麗[13]通過多組分高溫高壓蒸發(fā)模型仿真研究了二甲醚與液化石油氣混合液的蒸發(fā)過程，但其超臨界模型十分簡化。

本文建立了一個球對稱的全瞬態(tài)液滴蒸發(fā)模型，利用Fortran語言編寫可以計算單組分和雙組分燃料液滴的仿真程序，并驗證了程序的有效性；同時尋找合適的兩組分的煤油替代物，通過計算研究了超臨界環(huán)境下兩組分的煤油替代物的蒸發(fā)特性，得到環(huán)境壓力和環(huán)境溫度對其蒸發(fā)特性的影響規(guī)律。

1 液滴蒸發(fā)模型的建立

1.1 基本假設

由于實際的液體火箭發(fā)動機中經過霧化形成的液滴形狀并不完全為球形，且液滴所處的環(huán)境較復雜，為了便于仿真研究，對實際狀態(tài)下的液滴模型進行如下簡化：

a）假設液滴處于靜止環(huán)境中，忽略液滴周圍氣流速度及壓力波動；

b）假設液滴處于微重力環(huán)境下，忽略重力影響及自然對流；

c）假設液滴形狀在整個蒸發(fā)過程中保持為規(guī)則的球形；

d）忽略氣相的Dufour及Soret效應和熱輻射效應。

基于以上基本假設，即可建立一維的球對稱全瞬態(tài)液滴蒸發(fā)模型。

1.2 控制方程

氣相控制方程：

液相控制方程：

式中 r為物理平面中距離液滴中心的距離；t為時間；Yi為氣相和液相中各組分的質量分數；ρg，ρl分別為氣相和液相密度；Tg，Tl分別為氣相和液相溫度；vg，vl分別為氣相和液相速度；Dg，Dl分別為氣相和液相自擴散系數；Cpg，Cpl分別為氣相和液相定壓比熱容；λg，λl分別為氣相和液相熱傳導系數。

1.3 初始條件及邊界條件

1.3.1 初始條件

在t=0時刻，設置液滴內部溫度均勻分布且為Tl0，環(huán)境溫度均勻分布且為 Tg0=T∞，液滴內部只含燃料組分，不含環(huán)境氣體組分，壓力在整個計算域內均勻分布且保持為一恒值，即Pl0=Pg0=P∞。

1.3.2 亞臨界邊界

在氣相的無窮遠處，滿足的邊界條件為：T=T∞、Yi=0或1（當i組分為燃料組分時為0，當i組分為氮氣時為1）；

1.3.3 超臨界邊界條件

當液滴表面氣相側溫度和壓力均超過混合物的臨界值時，即認為達到超臨界狀態(tài)，此時液滴表面張力減小為零，液滴表面消失，氣相與液相之間的熱物性趨于連續(xù)[14]。因此，當液滴表面達到超臨界狀態(tài)時，不再區(qū)分氣液相，只需定義液滴中心處和液滴表面處的邊界條件，其定義與亞臨界相同。

1.4 數值計算方法

由于在液滴蒸發(fā)過程中液滴表面隨著蒸發(fā)過程不斷變化，如果直接對控制方程進行離散和求解，就需要在每個時間步長計算完成之后重新劃分計算網格。為了避免該情況的發(fā)生，采用式（7）和式（8）進行坐標轉換[14]，將物理平面轉換到計算平面中，使液滴表面始終固定在計算坐標為1的位置處，如圖1所示。

兩坐標系的空間位置轉換為

兩坐標系中的速度轉換為

式中 ζ為計算平面中距離液滴中心的距離；Rs為液滴半徑；SR˙為液滴半徑的變化率。

由于液滴表面附近的流場參數變化較大，因此計算網格的劃分采用非均勻網格，對液滴表面附近的網格進行加密處理，從而提高液滴表面附近的計算精度。

圖1 網格節(jié)點及其控制容積分布示意Fig.1 Schematic Diagram of the Grid Nodes and the Control Volume W，E—控制容積N兩側的網格節(jié)點；w，e—控制容積 N的控制容積邊界點

經過坐標轉換后的控制方程的離散方法采用有限體積法，其控制容積邊界取為相鄰網格節(jié)點的中心處，如圖1所示。

以氣相連續(xù)性方程為例，經過坐標轉換得：

對式（9）在控制容積N上，根據上一時刻(t - Δt)控制容積N的密度，對式（9）的時間及空間導數進行積分，可得計算平面上的離散方程：

2 液滴蒸發(fā)模型的驗證

為了檢驗計算模型及程序的有效性，利用程序計算了直徑為50 μm的正庚烷液滴在不同壓力和溫度的氮氣環(huán)境下的蒸發(fā)過程，并將本文的計算結果與文獻[15]的結果進行對比，如圖2所示。

圖2 仿真結果對比驗證曲線Fig.2 Simulation Result Comparison and Verification Pr—環(huán)境壓力與正庚烷臨界壓力的比值；Tr—環(huán)境溫度與正庚烷臨界溫度的比值

從圖 2中可以看出，本文所建立的計算模型與文獻的計算結果在部分位置存在偏差，出現這種偏差的原因是本文中采用的計算氣液相熱物性參數的估算方法及其在高壓下的修正方法與文獻中采用的方法不同，但總體來說，本文的仿真結果與文獻的仿真結果基本吻合，可以認為本文所建立的計算模型是正確的。

3 煤油替代物的選取

目前，中國使用的航空煤油型號主要為RP-3，文獻[3]中通過氣相質譜色譜聯用儀對 RP-3航空煤油的組分進行了分析，結果表明RP-3由42.330%的鏈烷烴，21.348%的環(huán)烷烴、24.045%的苯系物、8.184%的萘類和4.093%的其它組成。可見RP-3航空煤油是一種組分十分復雜的混合物，目前還無法實現對這種復雜混合物液滴的超臨界蒸發(fā)過程直接進行仿真計算。因此必須選取一種或幾種物質作為其替代物進行仿真研究。最常見的是采用正十二烷作為煤油替代物進行相關的仿真研究。采用單一的一種組分替代煤油雖然簡單但忽略了煤油是一種混合物的重要特性，而采用過多組分又會造成仿真計算難以實現。因此本文采用文獻[16]和文獻[17]中提出的兩組分的 RP-3航空煤油的替代物：（質量分數）80%正癸烷和 20%1，2，4-三甲基苯。

4 計算結果與分析

4.1 兩組分煤油替代物液滴蒸發(fā)過程

本文在計算中采用了質量分數為 80%正癸烷和20% 1，2，4-三甲基苯組成的兩組分煤油替代物液滴，液滴直徑為100 μm，液滴初始溫度為300 K，環(huán)境氣體為氮氣。已知正癸烷的臨界壓力為2.1 MPa，臨界溫度為 617.7 K，1，2，4-三甲基苯的臨界壓力為 3.2 MPa，臨界溫度為649.1 K。

首先計算P=4 MPa、T=800 K，P=6 MPa、T=1000 K，P=8 MPa、T=1200 K，P=10 MPa、T=1400 K，P=12 MPa、T=1600 K 5種工況下的液滴的蒸發(fā)過程。

圖3和圖4分別為上述5種工況下液滴的相對直徑平方隨時間的變化曲線和液滴表面溫度隨時間的變化曲線?？梢钥闯?種工況中的環(huán)境壓力和溫度都已經超過了兩種組成物的臨界值，但只有P=10 MPa、T=1400 K和P=12 MPa、T=1600 K兩種工況在蒸發(fā)過程中出現了超臨界遷移，達到了臨界點。這主要是因為液滴初始溫度遠遠低于環(huán)境溫度，液滴在蒸發(fā)過程中首先要經歷液滴表面溫度的快速升溫過程，在這個過程中液滴表面的組分也發(fā)生著變化，液滴表面混合物的臨界溫度也在不斷變化，只有當液滴表面的溫度和壓力同時達到混合物的臨界溫度時才能發(fā)生超臨界遷移。因此，只有在溫度和壓力都遠遠大于組分的超臨界溫度和壓力時，才能在蒸發(fā)完畢之前出現超臨界遷移。

圖3 不同環(huán)境壓力和環(huán)境溫度下液滴直徑變化曲線Fig.3 Variation of Droplet Diameter under Differenten Environmental Pressure and Environmental Temperature

圖4 不同環(huán)境壓力和環(huán)境溫度下液滴表面溫度變化曲線Fig.4 Variation of Droplet Surface Temperature under Different Environmental Pressure and Environmental Temperature

在每個工況下液滴的蒸發(fā)過程可以分為3個階段：a）吸熱膨脹，該階段主要是液滴吸熱快速升溫，液滴半徑增大；b）緩慢蒸發(fā)，該階段液滴表面溫度已經達到一個較高溫度，液滴開始蒸發(fā)，液滴半徑緩慢減?。籧）快速蒸發(fā)，該階段液滴溫度已基本達到最大值且保持不變，所有吸收的熱量均用于液滴蒸發(fā)，液滴半徑快速減小。

4.2 壓力對兩組分煤油替代物液滴蒸發(fā)的影響

圖5和圖6分別為環(huán)境溫度T=1400 K，環(huán)境壓力分別為P=4 MPa、P=6 MPa、P=8 MPa、P=10 MPa和 P=12 MPa的工況下液滴的相對直徑平方隨時間的變化曲線和液滴表面溫度隨時間的變化曲線。由圖5、圖6可以看出，隨著環(huán)境壓力的增加液滴的蒸發(fā)加快，這主要是因為環(huán)境壓力的增加使汽化焓減小和液滴蒸發(fā)率增加。同時，隨著環(huán)境壓力的增加液滴表面溫度在環(huán)境壓力較小時有所增加，但當環(huán)境壓力較大時基本保持不變。但隨著環(huán)境壓力的增加，液滴表面溫度的初始升溫速度增加，且更容易達到臨界點，這主要是因為隨著壓力的增加，導熱系數增大，加快了熱傳導過程，提高了氣相溫度向液相的傳輸速率，加快了液滴表面溫度的初始升溫速率。

圖5 T=1400 K時不同環(huán)境壓力下液滴直徑變化曲線Fig.5 Variation of Droplet Diameter under Different Environmental Pressure at T=1400K

圖6 T=1400 K時不同環(huán)境壓力下液滴表面溫度變化曲線Fig.6 Variation of Droplet Surface Temperature under Different Environmental Pressure at T=1400K

4.3 溫度對兩組分煤油替代物液滴蒸發(fā)的影響

圖7 和圖8分別為環(huán)境壓力P=10 MPa，環(huán)境溫度分別為T=800 K、T=1000 K、T=1200 K、T=1400 K、T=1600 K的工況下液滴的相對直徑平方隨時間的變化曲線和液滴表面溫度隨時間的變化曲線。

圖7 P=10 MPa時不同環(huán)境溫度下液滴直徑變化曲線Fig.7 Variation of Droplet Diameter under Different Environmental Temperature at P=10 MPa

圖8 P=10MPa時不同環(huán)境溫度下液滴表面溫度變化曲線Fig.8 Variation of Droplet Surface Temperature under Different Environmental Temperature at P=10MPa

由圖7、圖8可以看出，隨著環(huán)境溫度的增加，液滴的蒸發(fā)速率極大地增加，相比于環(huán)境壓力對蒸發(fā)速率的影響，環(huán)境溫度對蒸發(fā)速率的影響要大很多。同時，隨著環(huán)境溫度的增加，液滴表面溫度所能達到的最大值幾乎保持不變，但蒸發(fā)初期的液滴表面溫度的上升速度隨著環(huán)境溫度的增加而增加，與環(huán)境壓力的影響相似，隨著環(huán)境溫度的增加，液滴更容易達到臨界點。

5 結 論

本文首先建立了超臨界環(huán)境下液滴的非穩(wěn)態(tài)蒸發(fā)模型的仿真計算程序并驗證了計算程序的有效性，然后尋找了適合于RP-3航空煤油的兩組分替代物，最后通過仿真計算研究了超臨界環(huán)境下兩組分煤油替代物液滴的蒸發(fā)過程及環(huán)境壓力和溫度對其的影響。所得結論如下：

a）只有環(huán)境溫度和壓力大于組分的臨界溫度和壓力的2～3倍時，液滴在蒸發(fā)過程中才能發(fā)生超臨界遷移，達到臨界點；

b）環(huán)境壓力的增大能加快蒸發(fā)速率，提高液滴表面溫度的初始升溫速率；

c）相比于環(huán)境壓力，環(huán)境溫度的增大能極大地加快蒸發(fā)速率和提高液滴表面的升溫速率，且環(huán)境溫度和環(huán)境壓力較高時，其對液滴在蒸發(fā)過程中所能達到的最大溫度影響不大。