接觸網(wǎng)周圍流體流動(dòng)及傳熱特性的數(shù)值模擬

楊曉彬，林國(guó)松,2，李原，張碩，潘飛

(1. 西南交通大學(xué)電氣工程學(xué)院，四川 成都 611756；2. 西南交通大學(xué)國(guó)家軌道交通電氣化與自動(dòng)化工程技術(shù)研究中心，四川 成都 611756)

接觸網(wǎng)作為電氣化鐵路牽引供電系統(tǒng)的重要組成部分，承擔(dān)著將從牽引變電所獲得的電能供給電力機(jī)車的重要任務(wù)，其運(yùn)行性能和工作狀態(tài)是否良好直接影響著弓網(wǎng)系統(tǒng)的受流質(zhì)量，進(jìn)而關(guān)系到鐵路運(yùn)輸?shù)陌踩院涂煽啃?。接觸網(wǎng)常年裸露于空氣中，大風(fēng)、懸掛結(jié)構(gòu)覆冰等隨機(jī)環(huán)境荷載都會(huì)造成接觸線偏移或振動(dòng)，嚴(yán)重時(shí)甚至?xí)霈F(xiàn)大幅舞動(dòng)現(xiàn)象，導(dǎo)致線索、金具磨損嚴(yán)重，使用壽命縮短[1]。因此，研究接觸網(wǎng)導(dǎo)線的氣動(dòng)參數(shù)具有重要意義。另外，隨著高寒、高海拔電氣化鐵路投入運(yùn)營(yíng)，因地形、氣候等環(huán)境因素影響，接觸網(wǎng)覆冰現(xiàn)象頻發(fā)，在線防冰技術(shù)應(yīng)運(yùn)而生，其技術(shù)核心之一是根據(jù)接觸線表面的熱平衡方程來(lái)確定防冰電流[2-3]。但高速鐵路接觸線的截面形狀非標(biāo)準(zhǔn)圓形，而是兩側(cè)帶有懸吊溝槽的規(guī)則圓形，其傳熱特性與圓形導(dǎo)線不同，這就需要研究分析不同雷諾數(shù)下接觸線表面的傳熱參數(shù)。目前，對(duì)于氣動(dòng)參數(shù)的研究，主要有風(fēng)洞試驗(yàn)和數(shù)值模擬2 種方法。謝強(qiáng)等[4]采用風(fēng)洞試驗(yàn)方法，測(cè)量了不同湍流強(qiáng)度下AC-150 型接觸線的氣動(dòng)參數(shù)。然而，風(fēng)洞試驗(yàn)存在費(fèi)用高、周期長(zhǎng)、測(cè)量精度有限以及試驗(yàn)?zāi)Ｐ秃蛯?shí)際問(wèn)題的等效性難以確定等缺點(diǎn)[5]。因此，眾多學(xué)者通過(guò)計(jì)算流體力學(xué)軟件數(shù)值計(jì)算了多種不同情況下導(dǎo)線的氣動(dòng)參數(shù)。張友鵬等[6]通過(guò)建立正饋線的外流場(chǎng)模型，研究了大風(fēng)區(qū)段擋風(fēng)墻對(duì)正饋線氣動(dòng)特性的影響；宋洋[7]基于CFD 仿真模型計(jì)算了不同磨損程度接觸線的氣動(dòng)參數(shù)，進(jìn)而確定了接觸網(wǎng)舞動(dòng)的發(fā)生區(qū)域；汪宏睿等[8]則利用CFD 軟件(Fluent)模擬了接觸線周圍的空氣流場(chǎng)，并研究了不同迎風(fēng)攻角、不同風(fēng)速下接觸線氣動(dòng)參數(shù)的變化規(guī)律。但上述文獻(xiàn)均未系統(tǒng)地分析不同雷諾數(shù)下承力索和接觸線氣動(dòng)參數(shù)的差異，且流動(dòng)模型的模擬精度有待考量。對(duì)于導(dǎo)線傳熱參數(shù)的研究，多以電力系統(tǒng)中的架空輸電線為主，而與接觸線相關(guān)的研究則集中在接觸線溫度場(chǎng)方面。MATTERA 等[9]通過(guò)建立接觸線的二維徑向熱路網(wǎng)絡(luò)，對(duì)因徑向溫度梯度過(guò)高導(dǎo)致的機(jī)械故障進(jìn)行了熱分析；郭鳳儀等[10]則利用COMSOL 軟件研究了電弧位置和電弧能量對(duì)弓網(wǎng)系統(tǒng)接觸線溫度的影響。因此，為探究不同雷諾數(shù)下接觸網(wǎng)導(dǎo)線的氣動(dòng)參數(shù)和傳熱參數(shù)，本文以承力索和接觸線為研究對(duì)象，通過(guò)建立接觸網(wǎng)導(dǎo)線的流體-熱耦合數(shù)值模型，分別對(duì)承力索和接觸線的流場(chǎng)和溫度場(chǎng)進(jìn)行數(shù)值模擬，并將數(shù)值模擬結(jié)果與SSTk-ω模型以及公開(kāi)的實(shí)驗(yàn)數(shù)據(jù)進(jìn)行比較，以分析不同流動(dòng)模型的模擬精度。

1 流體?熱耦合數(shù)值模型

在流體力學(xué)中，任何流動(dòng)問(wèn)題都必須遵循3個(gè)基本的物理學(xué)原理：質(zhì)量守恒定律、牛頓第二定律和能量守恒定律。因此，首先，基于流體流動(dòng)的基本控制方程給出流體-熱耦合模型的數(shù)學(xué)描述；然后，針對(duì)接觸網(wǎng)導(dǎo)線繞流確定流動(dòng)模型；最后，給定其物理模型和參數(shù)定義。

1.1 模型的數(shù)學(xué)描述

由于本文研究的流體為空氣，流場(chǎng)中流體的流速范圍為0.001～30 m/s，且接觸網(wǎng)導(dǎo)線與流體的溫差不大，故為簡(jiǎn)化分析，做出如下假設(shè)：

1)流體為不可壓縮牛頓流體；

2)流體的物性參數(shù)為常數(shù)且無(wú)內(nèi)熱源；

3) 流體的體積力、體積加熱、黏性耗散熱可忽略不計(jì)。

在上述假設(shè)的基礎(chǔ)上，結(jié)合氣體的狀態(tài)方程，流體-熱耦合模型可描述為如下形式：

式中：ux，uy，uz分別表示坐標(biāo)x，y，z方向上流體的速度分量；ρ，μ，k，cv分別表示流體的密度、動(dòng)力黏度、導(dǎo)熱系數(shù)和比定壓熱容，p，e，T，t分別表示壓力、單位質(zhì)量流體的內(nèi)能、溫度和時(shí)間。該方程組包含了10個(gè)未知量，而列舉的方程僅有6個(gè)，再給定4個(gè)物性參數(shù)：ρ，μ，k，cv，便可使方程組封閉。

1.2 流動(dòng)模型

接觸網(wǎng)導(dǎo)線繞流的數(shù)值模擬，本質(zhì)上是鈍體繞流問(wèn)題，而圓柱繞流作為鈍體繞流中的典型代表，直到現(xiàn)在仍是CFD 領(lǐng)域的一大研究熱點(diǎn)。對(duì)于低雷諾數(shù)流動(dòng)，文獻(xiàn)[11]使用二維非定常層流求解器數(shù)值求解了40≤Re≤140內(nèi)圓柱的強(qiáng)制對(duì)流換熱特性，且模擬結(jié)果與試驗(yàn)結(jié)果非常吻合。對(duì)于高雷諾數(shù)流動(dòng)，CFD 軟件中引入了2 類湍流模擬方法：雷諾平均模擬方法(Reynolds Averaged Navier-Stokes Simulation, RANS)和大渦模擬方法(Large Eddy Simulation, LES)。目前，F(xiàn)luent 提供了多種基于RANS 的湍流模型，主要包括：?jiǎn)畏匠棠Ｐ?Spalart-Allmaras 模型)、雙方程模型(k-ε模型、k-ω模型及其變種模型)及雷諾應(yīng)力模型等，但文獻(xiàn)[12-13]指出，使用非穩(wěn)態(tài)RANS 模擬亞臨界區(qū)的圓柱繞流時(shí)，模型的適用性下降，不能準(zhǔn)確預(yù)測(cè)邊界層從層流到湍流的轉(zhuǎn)變以及邊界層的分離，進(jìn)而導(dǎo)致氣動(dòng)參數(shù)過(guò)度預(yù)測(cè)。與RANS相比，雖然LES對(duì)計(jì)算機(jī)的內(nèi)存和計(jì)算速度的要求較高，但文獻(xiàn)[14-15]使用LES 成功地模擬了亞臨界區(qū)的圓柱繞流，且其預(yù)測(cè)精度較RANS有很大提高。

基于上述描述，本文的流動(dòng)模型選擇層流模型和大渦模型，即低雷諾數(shù)的流動(dòng)模擬采用層流模型，高雷諾數(shù)的流動(dòng)模擬采用大渦模型。

1) 層流模型

在笛卡爾坐標(biāo)系下，給出非定常流動(dòng)層流模型的控制方程：

式中：i,j=x,y,z；τij表示作用在垂直于i方向的平面上沿j方向的應(yīng)力；fi表示作用在單位質(zhì)量流體上的體積力沿i方向的分量。

2) 大渦模型

LES 通過(guò)在Fourier 空間或物理空間對(duì)非穩(wěn)態(tài)N-S 方程進(jìn)行過(guò)濾，以濾除小于過(guò)濾器寬度或網(wǎng)格尺寸的旋渦，從而得到大渦的控制方程。以不可壓流動(dòng)為例，大渦模型的控制方程如下：

式中：帶有上劃線的變量為濾波后的場(chǎng)變量，σij，τij分別表示作用在垂直于i方向的平面上沿j方向的黏性應(yīng)力張量和亞格子應(yīng)力。

1.3 物理模型和參數(shù)定義

電氣化鐵路中，承力索多采用銅合金絞線，其截面近似為圓形，而接觸線截面形狀則為兩側(cè)帶有懸吊溝槽的規(guī)則圓形。雖然它們的截面形狀不同，但外流場(chǎng)模型具有相似的物理形式，因此，這里僅給出承力索外流場(chǎng)的三維物理模型，具體見(jiàn)圖1。外流場(chǎng)的三維物理模型采用笛卡爾坐標(biāo)系作為參考坐標(biāo)系，坐標(biāo)原點(diǎn)位于承力索截面中心處，x，y，z方向分別為順流向、橫流向和展向。以承力索直徑D為特征尺寸，三維物理模型的計(jì)算域尺寸為60D×40D×πD，承力索中心距上、下游邊界的距離分別為20D和40D，距側(cè)面邊界的距離都為20D，這種配置可使阻塞效應(yīng)忽略不計(jì)。

圖1 承力索外流場(chǎng)的三維物理模型Fig.1 Three-dimensional physical model of the external flow field of carrier wire

對(duì)于鈍體繞流，氣動(dòng)參數(shù)主要包括阻力系數(shù)CD，升力系數(shù)CL以及斯特勞哈爾數(shù)St，都屬于無(wú)量綱量，具體定義如下：

式中：FD和FL分別表示鈍體受到的阻力和升力；A⊥表示鈍體垂直于流向的投影面積；f表示旋渦脫落頻率。

傳熱參數(shù)中，鈍體表面與周圍流體的對(duì)流換熱可通過(guò)下式描述：

式中：h表示對(duì)流換熱系數(shù)；Tw，T∞分別表示鈍體壁面溫度和流體溫度；n表示鈍體壁面的法向方向。

鈍體壁面的局部努賽爾數(shù)Nuw和面積平均努賽爾數(shù)Nua可分別通過(guò)式(6)和式(7)計(jì)算。

式中：hw表示局部換熱系數(shù)；A表示鈍體壁面面積。

氣動(dòng)參數(shù)和傳熱參數(shù)的時(shí)均值均可通過(guò)下式表示：

式中：t1–t0為統(tǒng)計(jì)時(shí)間，需至少為15 倍的渦脫周期。

2 數(shù)值模擬

目前，國(guó)內(nèi)主要的高速鐵路運(yùn)營(yíng)線路，如京滬線、哈大線均采用截面積為150 mm2的接觸線和120 mm2的承力索。因此，以這2 種導(dǎo)線為研究對(duì)象，進(jìn)行數(shù)值模擬。承力索截面近乎圓形，一般可將其視為圓柱體。接觸線截面為帶有懸吊溝槽的規(guī)則圓形，具體見(jiàn)圖2。因合金種類識(shí)別溝槽和圓角的半徑很小，故圖2忽略了它們對(duì)接觸線截面的影響。

圖2 接觸線截面Fig.2 Cross section of contact wire

由于鈍體壁面處存在速度邊界層和熱邊界層，邊界層內(nèi)速度梯度和溫度梯度較大，為精確求解近壁面處流場(chǎng)和溫度場(chǎng)，就需要在邊界層內(nèi)布置較多的節(jié)點(diǎn)。因此，文中接觸網(wǎng)導(dǎo)線外流場(chǎng)采用分塊劃分方式，近壁面采用160×30 的O 型網(wǎng)格，且邊界層網(wǎng)格的設(shè)置需符合流動(dòng)模型的要求，其他區(qū)域的網(wǎng)格則為四邊形網(wǎng)格(二維)或六面體網(wǎng)格(三維)，展向方向的網(wǎng)格尺寸約為0.1D，具體的網(wǎng)格劃分見(jiàn)圖3和圖4。

圖3 承力索計(jì)算域網(wǎng)格Fig.3 Mesh of the computational domain of carrier wire

圖4 接觸線計(jì)算域網(wǎng)格Fig.4 Mesh of the computational domain of contact wire

根據(jù)1.2 節(jié)中的流動(dòng)模型，給出數(shù)值模擬中不同雷諾數(shù)下流場(chǎng)的求解設(shè)置及流動(dòng)模型，具體見(jiàn)表1。

表1 不同雷諾數(shù)下流場(chǎng)的求解設(shè)置及流動(dòng)模型Table 1 Solution setting and flow model of flow field under different Reynolds numbers

數(shù)值模擬過(guò)程中，入口設(shè)置為速度入口，來(lái)流速度u可根據(jù)不同的雷諾數(shù)計(jì)算得到，來(lái)流溫度T∞=293.15 K；導(dǎo)線表面設(shè)置為無(wú)滑移、恒溫壁面，壁面溫度Tw=298.15 K；出口設(shè)置為自由流出口，上下壁面都設(shè)置為無(wú)滑移、絕熱壁面，展向方向的2個(gè)對(duì)稱面設(shè)置為周期面，流體的物性參數(shù)則根據(jù)表2設(shè)置。

表2 293.15 K時(shí)流體的物性參數(shù)Table 2 Physical parameters of fluid at 293.15 K

至于方程的離散和求解，穩(wěn)態(tài)時(shí)，壓力與速度耦合的求解采用SIMPLE方法，非穩(wěn)態(tài)時(shí)，則采用PISO 方法；層流模型時(shí)，動(dòng)量方程和能量方程中的對(duì)流項(xiàng)采用2 階迎風(fēng)格式，時(shí)間項(xiàng)則采用2 階隱式推進(jìn)法；大渦模型時(shí)，動(dòng)量方程和能量方程中的對(duì)流項(xiàng)采用有界中心差分格式，時(shí)間項(xiàng)則采用有界2階隱式推進(jìn)法。

3 模型驗(yàn)證及結(jié)果分析

3.1 流動(dòng)模型的驗(yàn)證

為驗(yàn)證流體-熱耦合數(shù)值模型中流動(dòng)模型的模擬精度，采用層流(大渦)模型、SSTk-ω模型分別對(duì)雷諾數(shù)為300，500，1 000，3 000 和5 000 時(shí)的承力索流場(chǎng)和溫度場(chǎng)進(jìn)行數(shù)值模擬，提取其氣動(dòng)參數(shù)和傳熱參數(shù)，并將數(shù)值模擬結(jié)果與公開(kāi)的實(shí)驗(yàn)數(shù)據(jù)[16]進(jìn)行比較，以分析各流動(dòng)模型的模擬精度，具體見(jiàn)表3。

由表3可知，當(dāng)Re≤1 000時(shí)，層流模型的模擬誤差最小，最大相對(duì)誤差僅為7.0%，而SSTk-ω模型的最大相對(duì)誤差可達(dá)34.7%。當(dāng)Re≥3 000時(shí)，大渦模型的模擬誤差最小，最大相對(duì)誤差僅為5.1%，而SSTk-ω模型的最大相對(duì)誤差可達(dá)39.8%。因此，與SSTk-ω模型相比，低雷諾數(shù)選用層流模型、高雷諾數(shù)選用大渦模型的流體-熱耦合數(shù)值模型的模擬精度更高。

表3 數(shù)值模擬結(jié)果與實(shí)驗(yàn)數(shù)據(jù)的比較Table 3 Comparison of simulation results and experimental data

3.2 流場(chǎng)及氣動(dòng)參數(shù)

為探究不同雷諾數(shù)下接觸網(wǎng)導(dǎo)線的氣動(dòng)參數(shù)，利用求解器Fluent分別對(duì)承力索和接觸線的流場(chǎng)進(jìn)行數(shù)值模擬。圖5給出了幾個(gè)典型雷諾數(shù)下接觸網(wǎng)導(dǎo)線的流場(chǎng)分布。

從圖5 可看出，當(dāng)Re=1 時(shí)，流動(dòng)為蠕動(dòng)流，承力索上下表面的流場(chǎng)對(duì)稱，而接觸線的流場(chǎng)因溝槽喪失了對(duì)稱性；當(dāng)Re=30 時(shí)，流動(dòng)出現(xiàn)分離，承力索尾流中出現(xiàn)一對(duì)固定且對(duì)稱的旋渦，而接觸線尾流中的上旋渦稍比下旋渦大；當(dāng)Re=100時(shí)，尾流中都出現(xiàn)交替脫落的旋渦，并逐漸形成層流渦街；當(dāng)Re=5 000 時(shí)，自由剪切層中都存在層流到湍流的轉(zhuǎn)變，且尾流變?yōu)橥耆牧鳌?/p>

圖5 接觸網(wǎng)導(dǎo)線的流場(chǎng)分布Fig.5 Flow field distribution around catenary wires

圖6則給出了接觸網(wǎng)導(dǎo)線氣動(dòng)參數(shù)隨雷諾數(shù)的變化規(guī)律。

從圖6可看出，接觸線由于溝槽的存在，其氣動(dòng)參數(shù)與承力索的氣動(dòng)參數(shù)存在明顯的差異，主要表現(xiàn)在：1)當(dāng)Re≥300 時(shí)，接觸線的阻力系數(shù)高于承力索的阻力系數(shù)；2) 接觸線的升力系數(shù)明顯不同于承力索，其升力系數(shù)為負(fù)值，且當(dāng)Re≥1 000 時(shí)，其升力系數(shù)出現(xiàn)明顯的驟降；2) 當(dāng)Re≥300 時(shí)，接觸線的斯特勞哈爾數(shù)明顯低于承力索的斯特勞哈爾數(shù)。

圖6 不同雷諾數(shù)下接觸網(wǎng)導(dǎo)線的氣動(dòng)參數(shù)曲線Fig.6 Curves of aerodynamic parameters of catenary wires under different Reynolds numbers

3.3 溫度場(chǎng)及傳熱參數(shù)

為探究不同雷諾數(shù)下接觸網(wǎng)導(dǎo)線的傳熱參數(shù)，利用求解器Fluent分別對(duì)承力索和接觸線的溫度場(chǎng)進(jìn)行數(shù)值模擬。圖7給出了幾個(gè)典型雷諾數(shù)下接觸網(wǎng)導(dǎo)線的溫度場(chǎng)分布。

圖7 接觸網(wǎng)導(dǎo)線的溫度場(chǎng)分布Fig.7 Temperature field distribution around catenary wires

從圖7可看出，承力索和接觸線迎風(fēng)面的溫度邊界層厚度都比背風(fēng)面的薄，且溫度邊界層的厚度都隨雷諾數(shù)的增大而減小，這表明迎風(fēng)面的換熱性能要比背風(fēng)面的好，且導(dǎo)線整體的換熱性能隨雷諾數(shù)的增大而增強(qiáng)。但接觸線與承力索也存在不同之處，接觸線溝槽處的溫度梯度較小，這表明溝槽處的換熱性能較差。

工程中，對(duì)流換熱參數(shù)的計(jì)算一般是基于面積平均努賽爾數(shù)的，因此，圖8給出了接觸網(wǎng)導(dǎo)線面積平均努賽爾數(shù)——Nu a隨雷諾數(shù)的變化規(guī)律。

圖8 不同雷諾數(shù)下接觸網(wǎng)導(dǎo)線的傳熱參數(shù)曲線Fig.8 Curve of heat transfer parameters of catenary wires under different Reynolds numbers

從圖8可看出，承力索和接觸線的面積平均努賽爾數(shù)都隨雷諾數(shù)的增大而增大，同時(shí)，經(jīng)計(jì)算發(fā)現(xiàn)，當(dāng)Re≥50 時(shí)，接觸線的——Nu a比承力索的——Nu a低9%以上，這表明接觸線因溝槽的存在，其換熱性能一定程度上低于承力索的換熱性能。

4 結(jié)論

1)與SSTk-ω模型相比，低雷諾數(shù)選用層流模型、高雷諾數(shù)選用大渦模型的流體-熱耦合數(shù)值模型的模擬精度更高，氣動(dòng)參數(shù)的最大相對(duì)誤差僅為7.0%，傳熱參數(shù)的最大相對(duì)誤差僅為5.8%。

2) 在亞臨界區(qū)的氣動(dòng)參數(shù)中，接觸線的阻力系數(shù)高于承力索的阻力系數(shù)；接觸線的升力系數(shù)明顯不同于承力索，其升力系數(shù)為負(fù)值；接觸線的斯特勞哈爾數(shù)明顯低于承力索的斯特勞哈爾數(shù)。

3) 傳熱參數(shù)中，接觸線因溝槽的存在，其換熱性能一定程度上低于承力索的換熱性能(約9%以上)。