真空管道列車動態運行氣動特性研究

宋嘉源，李田，張繼業

西南交通大學 牽引動力國家重點實驗室，成都 610031

0 引言

真空管道運輸系統（ETT）通過低真空與磁懸浮技術顯著降低列車運行阻力，受到國內外研究人員的廣泛關注。2005年，沈志云院士[1]提出并論證了中國發展600km/h真空管道高速交通的設想。2013年，Musk[2]提出了Hyperloop 真空管道膠囊車的概念。2020年5月，西南交通大學啟動了“多態耦合軌道交通動模試驗平臺”項目，計劃建設超導磁浮技術與低真空技術結合的綜合性試驗平臺。與目前的高速列車[3]相比，真空管道列車高速運行時流場中存在一系列復雜的氣動現象[4]，如前驅激波、尾部激波、壅塞流動等。

數值計算是研究真空管道列車的主要方法之一。Niu 等[5-6]采用軸對稱的簡化模型作為研究對象，利用數值計算模擬了不同運行速度和阻塞比的真空管道列車氣動及氣動熱特性，擬合了壓力系數與列車運行速度的關系。Kim[7]、Choi[8]和Zhou[9]等研究了真空管道列車的氣動阻力與阻塞比、真空度及列車運行速度的關系。針對二維管道內復雜的激波現象，周鵬等[10-11]采用自適應網格技術精確捕捉激波位置，研究了1 250 km/h 二維管道列車激波的產生與反射傳播規律。三維磁浮列車需要考慮懸浮間隙和非軸對稱的幾何外形，管道內會出現更復雜的激波及氣動熱現象[12]。胡嘯等[13-14]采用IDDES 方法研究了三維磁浮列車的熱壓載荷分布及其瞬態特性，得到了熱壓載荷的頻率和分布規律。此外，為減少數值模擬計算的時間，Hou 等[15]指出準一維的管道列車模型和數值方法適用于真空管道流場研究，并推導了基于無黏與黏性氣體的典型流動模型和流型之間的過渡極限，其中等熵極限與激波脫離極限分別對應于壅塞流動和激波脫離現象。張曉涵等[16]結合進氣道理論與數值仿真闡明了管道氣流壅塞現象的形成機理。Li 等[17]采用一維流動理論快速預測了管道列車周圍的氣流參數。

目前多數學者的研究主要針對勻速運行的管道列車，對列車加/減速運行時的氣動阻力及特性的關注較少[18]。本文建立三維列車計算模型，采用剪切應力輸運（Shear Stress Transport,SST）k–ω模型求解管道內的流場，通過對比勻速工況和不同加速度工況的列車氣動阻力時程曲線以及不同時刻管道內的壓力分布，揭示加速度對列車氣動阻力和管道流場特性的影響規律，為真空管道試驗平臺的建設提供一定參考。

1 數值模型

1.1 控制方程

真空管道列車高速運行時的馬赫數大于0.3，會引起一系列復雜的氣動現象（如前驅激波、壅塞流動和尾部激波），因此數值模擬需要考慮空氣的可壓縮性。根據氣體連續性假設、動量守恒和能量守恒關系得到連續性方程、動量方程和能量方程[19]：

式中：ρ為氣體密度；t 為時間；xi、xj為直角坐標分量，ui、uj分別為流體速度u 在xi、xj上的分量；p 為壓力；τij為黏性應力張量；e 為內能；h 為焓；K 為熱傳導系數；T 為溫度；i、j 取1、2、3 分別代表 x、y、z方向。

管道內空氣劇烈壓縮和膨脹會引起較大的溫度變化，造成不同區域氣體黏度不同，采用Sutherland公式可以更準確地描述氣體溫度與黏度的關系[20]：

式中：μ為氣體動力黏度；μ0= 1.71 × 10?5Pa·s，為標準氣壓環境（101 325 Pa、273.15 K）下的氣體動力黏度；T0為溫度參考值，取273.15 K；S 為Sutherland常數，取110.4 K。

1.2 等熵極限與激波脫離極限

根據一維等熵流動理論，當來流馬赫數超過臨界值時，列車周圍氣流可通過的最大質量流量受到限制，氣流發生壅塞并產生前驅激波。壅塞現象對應的來流馬赫數Ma 為等熵極限的臨界速度[16]：

式中：阻塞比β=Strain/Stube，其中Strain為列車在運行方向（x 方向）的投影面積，Stube為管道在列車運行方向的投影面積；γ為空氣絕熱指數，取1.4。

當列車運行速度不低于等熵極限的臨界速度時，列車前方與尾部分別產生前驅激波和尾部激波。在管道參考系下，當列車速度繼續增大至大于尾部激波速度，尾部激波發生脫離并向列車后方傳播。根據一維流動理論[3]可以確定尾部激波脫離極限的臨界速度，如圖1所示。

圖1 等熵極限與激波脫離極限Fig.1 Isentropic limit and shock detachment limit

1.3 計算工況

為模擬列車速度達到不同極限（等熵極限和激波脫離極限）時的氣動現象，并對比不同加速度對氣動阻力及流場的影響，設置了3 種不同的計算工況。圖2 為不同工況下列車的運行速度及距離。工況1、2 的加速度a 分別為385 和278 m/s2，當列車加速至333 m/s 時，停止加速并保持勻速行駛；工況3 列車一直保持333 m/s 勻速運動。

圖2 不同計算工況的速度和距離時程曲線Fig.2 Curve of speed and distance under different cases

2 數值模擬方法

2.1 幾何模型

如圖3所示，采用三維真空管道列車模型進行數值模擬。為保證頭尾車流線型的細長比與高速磁浮列車接近，設置列車中心截面的輪廓曲線為列車采用3 編組結構，總車長度Lt為40 m，頭、尾車長度都為10 m，中間車長度為20 m，直徑D 為3.5 m。管道總長為550 m，在運動初始位置，列車頭車鼻尖距離隧道入口400 m。管道直徑為11.07 m，不同工況下列車的阻塞比均為0.1。

圖3 幾何模型Fig.3 Geometric model

2.2 計算域與邊界條件

圖4 為模型計算域與邊界條件，環境壓力設為20 265 Pa，溫度設為288.15 K。為保證出入口邊界不影響管道內流場分布，設置出入口邊界為壓力出口，相對壓力為0，管道壁面設置為固定絕熱壁面。采用重疊網格技術模擬列車運行過程，計算域分為重疊區域與背景區域，重疊區域直徑為4.7 m，重疊區域前、后邊界與頭、尾車距離分別為10 和30 m。

圖4 計算域與邊界條件Fig.4 Calculation domain and boundary conditions

2.3 計算模型

根據流體控制方程建立了三維可壓的數值計算模型，雷諾平均方法被廣泛應用于高速列車與真空管道列車數值模擬中[5]，本文采用了SST k–ω湍流模型。采用隱式非穩態耦合流算法求解控制方程，時間離散項為二階。計算的時間步長為1 × 10?4s，內部迭代步數為20。初始壓力和溫度分別為20 265 Pa和288.15 K，參考壓力與環境壓力相同。

2.4 網格劃分方法

由于真空管道列車運動計算域較大，為提高計算效率，采用重疊網格技術模擬列車運動過程。重疊網格技術通過重疊區域網格單元與背景區域網格單元的數據交換進行耦合，可以有效降低總網格數量。設置重疊區域與背景區域的網格單元尺寸一致可以降低計算誤差。此外，重疊網格的收斂速率與相同分辨率的單個網格相同。

劃分了如圖5所示的重疊網格及背景網格用于數值計算，為保證2 套網格過渡良好，背景網格的基礎尺寸設置為200 mm，重疊區域添加了2 個加密區，基礎尺寸分別為50 和100 mm，重疊網格最外層尺寸為200 mm。真空管道列車尾部流場更加復雜，為更好地捕捉尾部激波及后方尾流，設置尾車鼻尖與重疊網格交界面后側的距離為30 m。如圖5 右上框所示，列車近壁面劃分了總厚度為10 mm 的15 層邊界層網格，列車表面網格的基礎尺寸為50 mm。

圖5 計算域網格劃分Fig.5 Mesh of computational domain

2.5 激波管驗證

1955年，Trimpi 和Cohen[21]提出了預測激波管內氣流流動的理論公式，并與激波管試驗數據進行了對比。圖6 為該試驗的裝置示意圖，激波管左側為高壓段，右側為低壓段，交界面處有1 或2 片0.003 175或0.003 81 cm 的軟黃銅墊片作為膜片，將高低壓段隔開。初始狀態下，低壓段氣體壓力等于環境壓力（101 325 Pa），高壓段氣體壓力分別為45、70 和95 psi（即411 481、584 037 和755 378 Pa），氣體溫度均為室溫300 K。該試驗開始時，通過裝置戳破薄膜，交界面處將產生激波向低壓段傳播并衰減，在距膜片3.23、5.23、7.23、9.18 和11.18 英尺（即0.984 5、1.594 1、2.203 7、2.798 和3.407 7 m）處監測壓力，得到激波通過前后壓力比與膜片距離的關系。在數值模擬時使用了SST k–ω湍流模型進行求解，考慮氣體傳熱，求解精度為隱式二階。時間步長設為10?6s，最大迭代次數為20 次。

圖6 激波管裝置示意圖[21]Fig.6 Schematic diagram of shock tube[21]

圖7 為0.001 s 時刻激波管內的壓力分布和波系結構。圖8 為激波通過前后壓力比的試驗和數值模擬結果對比，高壓段壓力為45、70 和95 psi 條件下對應最大誤差分別為1.08%，1.16%和1.38%。結果表明，當前模型能準確捕捉激波管內激波的瞬態傳播特征。

圖7 0.001 s 時刻激波管內壓力分布Fig.7 Pressure distribution in shock tube at 0.001 s

圖8 Trimpi 實驗[21]與數值計算結果對比Fig.8 Comparison of Trimpi experiment and simulation[21]

3 數值模擬結果分析

3.1 氣動阻力瞬態特性分析

根據1.2 小節關于等熵極限與激波脫離極限的討論，在管道列車阻塞比為0.1 時，壅塞極限與激波脫離極限對應的臨界馬赫數分別為0.678 和0.741，在工況2 中對應的時刻分別為0.83 和0.91 s。圖9為工況2 頭、尾車氣動阻力的時程曲線，呈現先增加后穩定的趨勢；在0～0.5 s，管道內流場處于非壅塞狀態，尾車阻力增長緩慢，而頭車阻力隨列車加速基本不變。圖10 和11 分別為列車周圍壓力分布和列車表面壓力分布情況。0～0.1 s，列車啟動形成了流場初始擾動，正壓主要集中在頭車周圍；0.1～0.4 s，氣流在頭車鼻尖位置被劇烈壓縮產生正壓聚集區域，經過鼻尖時壓力快速降低；當氣流通過頭車周圍的流通區域時，頭車表面壓力由正變負，負壓在頭車與中間車連接處達到最大值。0.1～0.4 s，頭車阻力受到正負壓區域共同影響，正負壓力同時增大，因此阻力變化不明顯，而尾車阻力隨表面負壓增大逐漸增大。如圖11所示，在0.1、0.3 和0.5 s 時刻，頭車鼻尖的最大正壓分別為2.756、1.071 和0.214 kPa，與運行距離基本呈線性關系。列車勻加速過程中運行距離x ∝t2，因此列車表面最大壓力pmax∝t2。

圖9 工況2 氣動阻力時程曲線Fig.9 Curve of aerodynamic resistance under case 2

圖10 工況2 下0～0.5 s 列車周圍壓力分布Fig.10 Pressure distribution around the train under case 2 of 0～0.5 s

圖11 工況2 不同時刻列車表面壓力分布Fig.11 Pressure distribution of surface at different times under case 2

圖12 為t=0.6～1.2 s 列車周圍壓力分布。t=0.83 s 時刻，真空管道列車達到壅塞極限，通過列車周圍流通區域的空氣質量流量受到限制，列車前方的壅塞現象造成頭尾車氣動阻力增長速度明顯變快。同時，氣流在中間車與尾車連接處膨脹加速，橫截面積的突變導致該位置產生膨脹波，隨著氣流繼續加速，在尾車流線型膨脹過度并產生尾部激波。t=0.91 s 時刻，真空管道列車達到激波脫離極限，隨著列車繼續加速，尾部激波和膨脹波相對于列車向后運動，波系在管道與列車表面之間反射傳播。當膨脹波反射并擊中尾車表面時，尾車氣動阻力產生波動，如圖9所示。t>1.2 s 時，列車保持速度為1 200 km/h 的勻速運動，壅塞氣流聚集產生的前驅激波會繼續向前方流場傳播，但對頭車周圍流場影響較小。由于尾部激波和膨脹波不斷向后傳播，尾車周圍流場繼續變化，尾車阻力緩慢減小并逐漸趨于穩定。

圖12 工況2 下0.6～1.2 s 列車周圍壓力分布Fig.12 Pressure distribution around the train under case 2 of 0.6～1.2 s

圖13 為管道內流場穩定即t=1.5 s 時刻管道內的速度矢量圖。由于前驅激波在管道參考系下的運行速度大于列車，隨運行時間增加，未受擾動的流場區域相對于頭車越來越遠。在頭車與中間車連接處，壅塞現象限制了氣流流動，如圖13 的局部放大圖所示，頭車帶動的一部分氣流沒有通過流通區域進入后方流場，而是堆積在列車前方區域，導致壅塞段變長。尾部流場的分布具有周期特性，尾部激波在管道壁面間反射，并在后方流場形成了多個馬赫桿。氣流通過間斷面[4]時速度發生突變，間斷面后的流場速度基本為0，但列車尾流的影響距離較遠，并在尾流兩側形成強度較低的渦旋。

圖13 工況2 下1.5 s 時刻管道內速度矢量圖Fig.13 Velocity vector diagram in the tube under case 2 at 1.5 s

3.2 加速工況與勻速工況對比

圖14 為工況3 頭、尾車氣動阻力的時程曲線。由圖9 和14 可知，當流場基本穩定時，工況2 和3 的頭車氣動阻力分別為64 929 和64 448 N，偏差為0.74%，尾車氣動阻力分別為57 778 和56 027 N，偏差為3.13%，因此，加速階段基本不影響流場穩定時的氣動阻力。與工況2 相比，在工況3 的初期階段（t=0～0.2 s），頭車氣動阻力的增長速率基本不變，而尾車氣動阻力發生了周期性波動，振蕩幅值隨運行時間增長不斷衰減，在t=0.09 和0.11 s 時刻，尾車氣動阻力明顯降低。

圖14 工況3 列車氣動阻力時程曲線Fig.14 Curve of aerodynamic resistance under case 3

為分析工況3 尾車氣動阻力大幅波動的原因，提取了t=0.09、0.1 和0.11 s 3 個時刻列車表面壓力分布曲線，如圖15所示。可以看出氣流通過頭尾車時壓力的變化趨勢和幅值基本相同，列車表面壓力分布差異主要集中在中間車。不同時刻中間車與尾車連接處的壓力不同，從而影響了尾車整體的壓力分布。

圖15 工況3 不同時刻列車表面壓力Fig.15 Pressure distribution of surface at different times under case 3

圖16 為工況3 下不同時刻列車周圍壓力分布，t=0.09 s 時刻，頭車前方的壅塞氣流通過頭車與中間車連接處時產生斜激波，隨著時間增加，斜激波的角度逐漸變大，在管道壁面和列車表面的反射點位置發生移動，造成中間車表面壓力變化較大。當反射點位于中間車與尾車的連接處時，反射點附近的正壓區域會影響尾車表面壓力分布，使其出現周期性的壓力波動，從而影響尾車氣動阻力。隨著反射次數的增加，斜激波強度逐漸變小，尾車氣動阻力的振蕩幅值逐漸降低。

圖16 工況3 不同時刻列車周圍壓力分布Fig.16 Pressure distribution around the train at different times under case 3

3.3 不同加速工況對比

圖17～19 分別為工況3、2、1 不同時刻的管道壓力分布。在勻速工況下，頭尾車周圍壓力與運行距離基本呈線性關系，并在大約t=0.16 s 時刻頭尾車周圍壓力基本同時達到最值。在工況1 和2 下，尾車周圍壓力首先達到最小值，在加速階段完成后列車勻速運行，頭車周圍壓力繼續保持一段時間增長，在大約t=1.4 s（工況2）和t=1 s 時刻（工況1）達到最大值。這是由于在加速階段完成時刻，加速工況下頭車壓縮前方空氣過程較緩慢，對應的壓力梯度較小，前驅激波強度小于工況3，頭車將繼續壓縮前方空氣直到前驅激波后的氣流壓力和馬赫數穩定。總體看來，頭車氣動阻力和周圍壓力變化與運行速度變化相比表現出一定的滯后性，且加速度越小，滯后性越明顯。

圖17 工況3 不同時刻管道內壓力分布Fig.17 Pressure distribution in tube at different times under case 3

圖18 工況2 不同時刻管道內壓力分布Fig.18 Pressure distribution in tube at different times under case 2

圖19 工況1 不同時刻管道內壓力分布Fig.19 Pressure distribution in tube at different times under case 1

在頭尾車周圍壓力達到最值后，隨著列車勻速運動，前驅激波與頭車鼻尖間的壅塞段長度增加，尾部激波向后方流場的傳播距離也增加。以工況3 為例，當t=0.3、0.4 和0.5 s 時，壅塞段長度分別為12、18 和24 m，尾部激波段長度分別為11、16 和20 m，與運行時間成正比。對于尾部流場，勻速工況列車啟動速度較大，形成尾部激波、中心間斷面和尾部膨脹波等一系列氣動現象，氣流經過激波及膨脹波時壓力先降低后恢復。在加速度階段，由于啟動速度較低，尾部壓力波動不明顯。

4 結論

本文建立了三維真空管道列車模型，采用SST k–ω模型求解了管道內的流場，研究了勻速工況和兩個不同加速度工況下的瞬態列車氣動阻力、流場特性以及管道壓力分布，主要得到以下結論：

1）管道在非壅塞狀態下，尾車氣動阻力增長緩慢，頭車氣動阻力基本不變；達到壅塞狀態后，頭尾車氣動阻力快速增大；達到激波脫離極限和加速階段完成后，尾車與頭車氣動阻力都趨于穩定。

2）流場穩定時加速工況與勻速工況的氣動阻力相對誤差較小，由于勻速工況啟動速度較大，頭車與中間車連接處產生的斜激波影響列車表面壓力，導致尾車阻力大幅波動，波動幅值隨運行時間增長逐漸降低。此外，在勻速工況下，尾部流場還存在中心間斷面和尾部膨脹波等一系列氣動現象。

3）由于加速工況下頭車壓縮前方空氣過程較緩慢，前驅激波強度較低，頭車氣動阻力和周圍壓力的變化與運行速度變化相比表現出一定滯后性，并且加速度越小滯后性越明顯。

4）在加速階段，尾車周圍壓力比頭車先達到最值，勻速階段頭尾車周圍壓力基本不變；壅塞段長度和尾部激波段長度與列車運行時間成正比關系。