陳雨成,秦斌,梁習鋒
(中南大學 交通運輸工程學院,湖南 長沙 410075)
近年來,隨著高速列車技術的飛速發展,國內外研究人員通過風洞試驗和數值模擬等手段對高速列車的空氣動力學特性進行了大量研究[1-2]。然而,隨著傳統高速列車的運行速度越來越高,可能會出現相當大的問題:空氣阻力和氣動噪聲功率迅速增加[3-4],導致列車能耗急劇上升。此外,橫風等外部天氣條件會增加列車運行穩定性的安全隱患。諸如此類問題的出現制約了列車的進一步提速,因此構建真空管道交通系統的想法應運而生。真空管道交通系統是磁懸浮列車在抽成一定真空度的管道中運行。與在標準大氣壓下運行的高速列車相比,真空管道交通系統具有以下優點:1) 低壓管道環境支持更高的運行速度,同時顯著降低運行能耗;2) 更低的氣動噪聲;3) 更安全的運行環境。2013 年,SpaceX 公司的MUSK[5]提出了被稱為Hyperloop 的一種類似真空管道磁懸浮列車的概念性運輸方式。隨后TAYLOR 等[6]從可行性和經濟性的角度對Hyperloop 進行了分析。BRAUN 等[7-8]對hyperloop 的吊艙進行了最小化阻力、最大化升力的外形設計,并利用數值模擬以提高其空氣動力學性能。雖然對Hyperloop 的研究針對的是真空管道中長度較短的飛行器,但對真空管道高速列車的研究有一定參考價值。對于真空管道列車,研究人員對其外部流場的氣動特性進行了大量的理論和數值研究。總的來說,影響真空管道列車氣動性能的主要參數包括:列車速度v,管道壓力P和阻塞比β(=Atrain/Atube,其中Atrain和Atube分別為列車和管道截面積)等。KIM 等[9]使用二維數值模型研究了真空管道列車在v=500~700 km/h,P=0.01~0.05 atm(1 atm=101 325 Pa)和β=0.25~0.75條件下的氣動特性,得到在不同條件下v-β,β-P和P-v的關系曲線,并找到了在不同壓力和阻塞比的管道內可能引起激波的臨界速度。ZHANG[10]建立了一個直徑為3 m 的真空管道運輸系統二維數值模型,在v=50~300 km/h,P=0.001~0.1 atm 和β=0.21~0.69 的條件下對真空管道列車進行了數值模擬,并根據實用性和經濟性對各種參數給出了合理的取值范圍。MA 等[11]利用高溫超導磁懸浮實驗系統證實了列車速度、管壓和阻塞比是影響真空管道列車能量損失的3個重要因素,并得出當管壓和阻塞比恒定時,空氣阻力是速度的二次函數。對于真空管道列車的外形設計,CHEN 等[12-14]在二維尺度上對不同頭型和尾型的列車進行了數值模擬。在β=0.25~0.4 和P=1 000 Pa,CHEN 等[12-14]得出鈍頭形或橢圓形(2:1)的頭、尾型能達到最佳的空氣動力學性能。在β=0.21~0.69 和P=0.000 1~0.1 atm 時,BIBIN 等[13]得出橢圓和三角形(2:1)的頭、尾型為減阻的最佳選擇。對于真空管道中氣流特性的研究,JIA 等[15]建立了帶有壓力循環管道(PRD)的真空管道列車二維數值模型。PRD 的存在可以有效降低列車頭尾壓差,在P=10 000 Pa,β=0.30,v=400 m/s 時可降低高達30%。此外,還討論了管中出現的車頭壅塞現象。NIU 等[16]采用二維數值模型對v=400~2 000 km/h的真空管道列車進行模擬,再現了管內沖擊波在亞聲速、跨聲速和超聲速下的產生和發展。ZHOU 等[17-18]用二維數值模型展示了真空管道列車在v=1 250 km/h時管道內弓形激波、正激波、反射激波、Lambda 激波和菱形激波等激波簇結構的產生、發展與消失的過程,還對真空管道中臨界阻塞比和臨界來流馬赫數之間的關系進行了研究。隨后,ZHOU 等[19]建立真空管道磁懸浮列車三維模型,探究了P=0.01 atm,β=0.3 的真空管道中馬赫數Ma=1 的磁懸浮列車懸浮間隙對列車氣動特性的影響,得出間隙減小有利于縮短正激波的形成時間。此外,張曉涵等[20-21]采用二維數值模型,胡嘯等[22]采用三維數值模型,分別研究了不同氣壓、阻塞比和車速條件下真空管道中壅塞區的長度和激波的強度及結構變化的規律。綜上所述,首先,以往對真空管道列車氣動特性的研究較多使用二維數值模型,而采用三維模型對管道內部激波的產生及分布的研究較少。使用三維數值模型可以更準確地捕捉流場結構,獲得比二維模型更可信的結果。其次,對于較高的來流馬赫數對真空管道列車的影響研究較少。第三,較為缺乏對激波、渦流和邊界層的相互作用的詳細研究。因此,本文基于三維模型處理較高來流馬赫數下的真空管道列車問題,重點關注尾流流場結構。本文采用數值模擬方法,再現了不同壓力或來流條件下真空管道磁懸浮列車的流場分布。管道阻塞比β=0.15,壓力為P∞=0.1 和0.3 atm。來流馬赫數Ma∞(=v/c,其中v和c分別為流速和當地聲速)設置為0.490,0.654,0.817 和0.980。本文的結構如下:首先介紹算法和真空管道磁懸浮列車模型,包括幾何模型、網格及邊界條件及其驗證。然后分析了不同壓力和來流馬赫數下的數值模擬結果,包括整個流場的分布和激波在尾流中的傳播。最后給出結論。
真空度管道磁懸浮列車幾何模型采用基于高速磁懸浮列車EMU1的5 車編組模型,幾何拓撲模型如圖1 所示。磁懸浮列車由1 節頭車、3 節中間車以及1 節尾車構成,頭尾具有相同的流線型結構。坐標軸原點設置在頭車末端中心點,x,y,z軸分別沿列車的長、寬、高方向。計算域長830 m,軌道在計算域正中間,列車頭車鼻尖點距離入口200 m,尾車鼻尖點距離出口500 m。管道橫截面為扇形,面積為80 m2,管道阻塞比為0.15。以上計算域尺寸為全尺寸情況,實際計算以1:8 縮比尺度進行。

圖1 真空管道磁懸浮列車EMU1幾何拓撲模型Fig. 1 Geometric topology model of the EMU1 vacuum tube maglev train
本文采用ANSYS ICEM CFD 軟件的拓撲優化、多層網格加密技術、附面層網格技術與網格拉伸技術展開精細化網格劃分。本文全局采用三角形網格劃分,網格模型如圖2所示。頭、尾車網格40 mm 尺度,中車網格60 mm 尺度,流線型部分面的尺度為30 mm。考慮到尾流區復雜的流動現象,設置尾車流線型部分區域網格尺寸為10 mm。整體的加密區有3 層,尺度分別為200,400 和600 mm。共設置8 層附面層,0.3 mm 尺度。以上網格尺寸均為全尺寸情況,計算時采用1:8 縮比,模型總體網格規模約為2.1億。

圖2 計算網格示意圖Fig. 2 Grid model of tube and train
本文采用ANSYS FLUENT 軟件進行網格模型的流場仿真。穩態計算基于顯式壓力基求解方法,湍流模型選用SSTk-ω模型。針對馬赫數小于1 的真空管道列車氣動特性研究,采用SIMPLE 算法,王海明等[23]研究了列車首尾壓差變化規律(Ma=0.35-0.71),黃尊地等[24]分析了其氣動阻力(Ma=0.16-0.82)。本文為保證運算效率,穩態計算壓力-速度耦合選用SIMPLE 算法,壓力采用Standard 離散格式,動量、湍動能、比耗散率均采用2階迎風離散格式。非穩態計算以穩態流場作為初始流場,湍流模型選用LES,亞格子模型為Smagorinsky模型。壓力-速度耦合選用Coupled 方法,非穩態計算時間步長取5×10-5s,共計算10 000個時間步。
由于計算來流速度較高,空氣壓縮效應明顯,湍流模型選擇Density-Based,氣體采用理想氣體。為模擬無窮遠處的自由可壓縮流動,計算域入口處采用壓力遠場(給定來流馬赫數)邊界條件,流域出口處為壓力出口邊界條件。隧道壁面及軌道設置為滑移壁面邊界條件,其速度和來流速度一致,列車模型表面設定為固定壁面邊界條件。管道內氣壓P∞設置為0.1 atm 和0.3 atm。設置來流馬赫數Ma∞=0.490,0.654,0.817,0.980。
1.3.1 數值算法驗證
為了驗證數值仿真算法的合理性,以翼型Rae 2822 為驗證實例進行建模[25]并根據本文使用的湍流模型及邊界條件進行計算,壓力遠場邊界用于模擬無限流入的流動。圖3為翼型Rae 2822的數值驗證結果,其中顯示了壓力系數的分布以及數值模擬與實驗所得壓力系數的比較。選擇0.74 的來流馬赫數和25 000 Pa 的壓力條件,可以在翼型的上表面觀察到激波,如圖3(a)所示。如圖3(b)所示,模擬獲得的結果與實驗測量的結果非常吻合。

圖3 翼型 Rae 2822的數值驗證結果Fig. 3 Numerical verification result of Airfoil Rae 2822
此外,為了驗證本文所采用的數值算法捕捉激波間斷的能力,建立如圖4(a)所示的數值模型,基于LI 等[26]以縮放噴管原理設計的超聲速風洞試驗進行數值仿真算法驗證,并捕捉到了超聲速來流下管道內由楔形塊引起的初始斜激波與激波串。入口來流馬赫數為1.85,溫度為182 K。風洞試驗與數值算法捕捉到的管道激波串對比如圖4(b)~4(c)所示,兩者流場結構基本一致,數值算法不僅捕捉到了斜激波的反射與干涉,包括激波誘導管道壁面邊界層分離同樣也被精確捕捉。圖5比較了沿管道壁面壓力系數的實驗結果與數值仿真結果,兩者吻合較好。因此,本文所采用的數值算法與邊界條件可以精確捕捉真空管道中激波結構以及激波與邊界層之間的相互干擾,具有一定的可信度。

圖4 超聲速風洞數值驗證Fig. 4 Numerical verification of supersonic wind tunnel

圖5 頂板中心線壓力系數實驗數據與數值模擬結果的比較Fig. 5 Comparison between experimental data and numerical simulation results of roof centerline pressure coefficient
1.3.2 網格獨立性驗證
為驗證網格無關性,車體面網格尺寸和附面層設置保持不變,僅改變包裹列車的3層空間加密區尺寸,得到粗網格和細網格模型。從內到外,細網格空間加密區尺寸分別為100,200 和400 mm,粗網格空間加密區尺寸分別為300,600和900 mm。圖6給出3套網格在10 000至20 000步的整車氣動阻力系數,可知粗網格、中網格和細網格氣動阻力系數分別穩定在0.295,0.273 和0.268,中網格和細網格整車阻力系數幾乎重合,可認為將中網格作為計算網格模型滿足網格無關性要求。另外LES 模型要求列車表面y+在1 的量級,圖7 給出列車表面y+分布云圖,其值均滿足小于1 的條件,可知本文網格能滿足LES 模型計算要求。

圖6 不同網格模型整車阻力系數Fig. 6 CD of the entire train with different grid models

圖7 磁懸浮列車表面y+分布云圖Fig. 7 y+ distribution on the surface of the train
圖8 展示了壓力P∞=0.3 atm,來流馬赫數為0.490,0.654,0.817 和0.980 對應的馬赫數分布云圖。取列車上方氣流x=-10~25 m,y=0,z=0.7 m繪制馬赫數曲線,并將列車附近流場分為5 個區域:I 區(x<0,頭車上游區域)、II 區(x=0~1.64 m,頭車流線型區域、收縮段)、III 區(x=1.64~14.62 m,中車段區域)、IV 區(x=14.62~16.34 m,尾車流線型區域、擴張段)和V 區(x=16.34~25 m,尾流區域)。

圖8 0.3 atm壓力下不同來流條件的流場馬赫數云圖Fig. 8 Variation of Ma of the flow field in tube at different Ma∞, P∞=0.3 atm
一般來說,流場結構可分為3 種類型,如圖10(a)所示。對于類型1(Ma∞=0.490),整個外部流場為亞聲速,流速在I 區幾乎保持不變,在II 區域增加,在III 區馬赫數從0.591 加速到0.641,在IV 區下降,然后在V區經歷小幅波動。對于類型2(Ma∞=0.654),流場馬赫數在II 區和III 區中比類型1 增長更為迅速,在IV區馬赫數達到1。在尾車肩部開始出現激波,之后流場馬赫數在1附近波動,然后由于膨脹波而下降。在V 區氣流馬赫數的振蕩比類型1更劇烈,這可能導致更強的渦旋脫落。對于類型3(Ma∞>0.654),流場在III 區中變為超聲速,然后在IV 區繼續加速。在尾流中,可以觀察到斜激波和“X”型激波,干擾了邊界層和渦流。

圖10 0.3 atm壓力下不同來流條件列車上方氣流變化曲線Fig. 10 Changes in airflow above the train under different Ma∞, P∞ = 0.3 atm
圖9 和圖10(b)展示了壓力系數CP分布云圖及列車上方氣流的壓力系數曲線。對于類型1,頭車流線型的前端有一個高壓區。頭車前部附近區域為正壓區,中車和尾車及下游的流場為負壓區,中車附近形成高負壓區。對于類型2,頭車前部附近的壓力增加,正壓區逐漸向后移動。對于類型3,高負壓區向尾車后方移動,尾車后方形成顯著的壓力波動。

圖9 0.3 atm壓力下不同來流條件的流場壓力系數云圖Fig. 9 Variation of CP in tube under different Ma∞, P∞ = 0.3 atm
圖11 顯示了不同壓力的管道中不同來流條件下列車上方氣流的變化曲線。可見不同壓力條件下,氣流的變化趨勢基本相同。從氣流馬赫數來看,亞聲速區I,II 和III 中P∞=0.1 atm 管道氣流馬赫數略高于P∞=0.3 atm 管道。而在區域IV 出現超聲速區的情況下,壓力P∞=0.3 atm 管道氣流馬赫數在區域IV 和V 中高于P∞=0.1 atm 管道,波動也更為劇烈。從氣流壓力系數來看,總體上不同壓力條件下相同來流馬赫數的氣流壓力系數差異并不明顯。可以觀察到在更低真空度的P∞=0.1 atm 管道氣流在區域IV和V中波動更為平緩。

圖11 不同來流條件和壓力下列車上方氣流變化曲線Fig. 11 Changes in airflow above the train under different Ma∞ and P∞
在對真空管道列車空氣動力學特性的研究中,對不同來流條件下尾流區流場結構的分析具有極其重要的意義。隨著來流馬赫數的增加,激波對尾流區流場結構的影響逐漸增大。
圖12 展示了在壓力P∞=0.3 atm 的不同來流條件下尾車附近渦量的變化情況。尾車區域出現了較強烈的擾動渦量,分離流向后形成了典型的脫落尾渦結構。對于流場類型1 和2,并沒有形成尾部激波。而與類型1 相比,類型2 在氣流馬赫數增大后尾流渦的強度相應增大。對于類型3,尾渦的脫落相對于類型1 和2 受到抑制,渦脫過程在尾流中變得更加規律。這說明激波的出現會損耗附近的擾動能量,進而干擾和削弱尾渦結構的發展和向后傳遞。

圖12 不同來流條件的尾流渦量云圖Fig. 12 Vorticity contours in the wake at different Ma∞
在來流馬赫數接近1 之后,類型3 中不同來流馬赫數的跨聲速流場具有相似的流場馬赫數分布,而壓力系數呈現梯度分布,如圖13 所示。從圖中也可以觀察到當氣流通過激波位置時馬赫數的驟降和壓力系數的驟升。

圖13 不同來流條件的尾流流場結構,l為取值線Fig. 13 Wake structure of flow field with different Ma∞, l represents the line of the value
圖14 對來流馬赫數為0.980條件下列車尾流區域的流場結構進行了放大展示。當超聲速氣流經過尾車時,尾車流線型表面形成低壓區,邊界層在尾車流線型端部分離,誘導產生斜激波1。斜激波1遇到管道壁面后產生反射,并與壁面以及尾流中形成的渦相互干涉后,誘導形成一道“X”型激波結構(2,3,4和5)。與二維模型相比,三維數值模型對圓柱體管道內激波結構的模擬將導致管道壁面“菱形”激波結構的入射角逐漸增大,如圖15所示。

圖14 尾流流場結構圖Fig. 14 Wake structure behind TC

圖15 尾流下游流場結構Fig. 15 Wake structure of downstream extension area
圖16給出了Ma∞=0.980下的Q 判據[27]渦結構云圖,可以觀察到I,II,III 和IV 區管道壁面分布大量小尺度渦結構,且尾流區強度大于頭車前方區域。頭車流線型和尾車流線型上表面形成周期性條紋狀小尺度渦結構,直至氣流在頭車流線型肩部或尾車流線型端部發生分離。頭車底部與軌道之間的狹窄間隙受到高速氣流沖擊而形成了周期性的蠕蟲狀渦結構。尾流區兩側形成雙馬尾狀的大量發夾渦結構向下游發展。此外,在頭、尾車流線型兩側形成了平行于車體的長條狀渦,這種結構類似于機翼末端的脫落渦結構。

圖16 Q判據渦結構云圖Fig. 16 Vortices distribution based on Q-criterion
1) 在阻塞比β=0.15,真空度P∞=0.3 atm 的管道中:當來流馬赫數Ma∞=0.490時,整個流場區域均為亞聲速,沒有出現激波。當Ma∞=0.654時,氣流在尾車擴張段區域起點處馬赫數達到1,在擴張段中圍繞馬赫數為1上下波動并由多個激波的影響變為亞聲速。當Ma∞=0.817~0.980 時,氣流在尾車擴張段區域起點處馬赫數達到1,在擴張段區域繼續加速,并在尾流區形成多道斜激波。當管道真空度由0.3 atm 降為0.1 atm 后,氣流變化趨勢與真空度降低前基本保持一致,尾流區波動更為平緩。
2) 對于Ma∞>0.654的情況,不同來流馬赫數條件下流場中的氣流馬赫數分布相似,而壓力系數呈現梯度分布。隨著尾車附近出現超聲速區,將產生一系列斜激波和“X”型激波等結構。由于激波、渦流和邊界層的相互影響,激波和渦出現相互干擾與融合現象。研究成果能為真空管道列車不同來流速度和不同真空度情況尾流激波抑制以及氣動阻力優化設計提供工程指導。