喻 樂(lè) 和敬涵 王小君 薄志謙

（1.北京交通大學(xué)電氣工程學(xué)院 北京 100044 2.ALSTOM電網(wǎng)公司 英國(guó)斯塔福德 ST17 4LX）

1 引言

城市軌道交通牽引網(wǎng)大量采用第三軌（接觸軌道）供電方式。第三軌容易因絕緣老化、受損發(fā)生對(duì)走行軌或接地短路故障，嚴(yán)重威脅著軌道交通安全運(yùn)行，造成不同程度的停運(yùn)事故。

目前，國(guó)內(nèi)城市軌道交通直流饋線保護(hù)設(shè)備主保護(hù)原理為電流上升率與電流增量保護(hù)，該方法能區(qū)分大部分故障情況，但隨著城市軌道交通的快速發(fā)展，該原理對(duì)于運(yùn)行區(qū)段內(nèi)多機(jī)車(chē)用電與高阻短路等情況存在難以區(qū)分的問(wèn)題。因此，大量研究針對(duì)該問(wèn)題展開(kāi)：文獻(xiàn)[1-5]利用小波等方法以區(qū)分短路電流趨膚效應(yīng)特征與機(jī)車(chē)正常用電電流時(shí)間常數(shù)特征；文獻(xiàn)[6-8]采用穩(wěn)態(tài)分析方法對(duì)城市軌道交通供電系統(tǒng)饋線短路現(xiàn)象進(jìn)行了分析。

當(dāng)鋼軌發(fā)生短路故障時(shí)，暫態(tài)初始時(shí)刻鋼軌中電流會(huì)產(chǎn)生趨膚效應(yīng)，使得暫態(tài)初始時(shí)刻電阻大、內(nèi)電感小，但隨著時(shí)間的推移電阻逐漸減小而內(nèi)電感逐漸增大，由趨膚效應(yīng)引起暫態(tài)過(guò)程時(shí)間常數(shù)變化規(guī)律與MTRC地鐵短路試驗(yàn)結(jié)果基本一致[9]。對(duì)鋼軌短路趨膚效應(yīng)的準(zhǔn)確計(jì)算是研究針對(duì)城市軌道交通供電系統(tǒng)的配套直流饋線保護(hù)、監(jiān)測(cè)及故障定位等技術(shù)開(kāi)發(fā)的基礎(chǔ)，且對(duì)直流饋線保護(hù)定值整定有著重要作用。

目前鋼軌短路電流趨膚效應(yīng)的仿真計(jì)算所采用的方法主要有以下三類：①時(shí)域計(jì)算方法：由H. S.Carslaw與J. C. Jaeger[10]為解決導(dǎo)體熱傳導(dǎo)問(wèn)題最先提出，并由Edward J. Tuohy等人[11]引入趨膚效應(yīng)計(jì)算。文獻(xiàn)[12]將該算法應(yīng)用于鋼軌趨膚效應(yīng)計(jì)算，取得了良好的效果；②頻域計(jì)算方法：由R. J. Hill等人提出[13,14]，該方法將貝塞爾方程轉(zhuǎn)換至頻域，將電路整體進(jìn)行頻域計(jì)算，最后將結(jié)果采用傅里葉反變換法轉(zhuǎn)至?xí)r域，獲得其時(shí)域解。Lawrence J.Giacoletto[15]采用Laplace變換方式進(jìn)行頻域與時(shí)域的轉(zhuǎn)換，其思路與R. J. Hill較為類似，對(duì)于交流暫態(tài)計(jì)算較為實(shí)用；③導(dǎo)體細(xì)分法：由P.Silvester[16-18]提出，該方法的基本思路就是將整個(gè)導(dǎo)體劃分為有限個(gè)子部分，通過(guò)分析細(xì)分導(dǎo)體間的電磁關(guān)系，建立電感、電阻計(jì)算矩陣，獲得導(dǎo)體整體阻抗表達(dá)式從而進(jìn)行計(jì)算，該方法計(jì)算結(jié)果準(zhǔn)確，在短路計(jì)算中應(yīng)用較多。

本文采用階躍級(jí)數(shù)逼近的方法對(duì)鋼軌短路趨膚效應(yīng)建立了數(shù)學(xué)模型，并將該方法采用 S-Function建立鋼軌模型與牽引變電站模型在Matlab/Simulink環(huán)境下聯(lián)合仿真，提高了仿真計(jì)算速度及結(jié)果的準(zhǔn)確性。

2 鋼軌趨膚效應(yīng)計(jì)算

2.1 Matlab/Simulink仿真模型

本文基于 Matlab/Simulink建立城市軌道交通牽引供電系統(tǒng)模型，包含直流牽引變電站（10kV交流三相電源、二十四脈波牽引整流機(jī)組）模塊及第三軌及走行軌模塊[19]，如圖1所示。由于牽引變電站可由Matlab/Simulink自帶模塊組成，對(duì)此不再贅述。本文采用S-Function對(duì)鋼軌進(jìn)行建模，將鋼軌等效為可控電流源，通過(guò)對(duì)每個(gè)計(jì)算步長(zhǎng)獲得鋼軌電壓并計(jì)算其電流返回模型，以實(shí)現(xiàn)對(duì)鋼軌趨膚效應(yīng)特性仿真計(jì)算。

圖1 牽引供電系統(tǒng)Matlab/Simulink模型圖Fig.1 Matlab/Simulink model of DC railway power supply system

由于鋼軌具有較大的不規(guī)則橫截面，其截面示意圖如圖2所示。趨膚效應(yīng)使得鋼軌電阻和內(nèi)電感的計(jì)算變得非常復(fù)雜。本文采用等效截面積法將鋼軌等效為等截面積的圓柱形導(dǎo)體，截面半徑為

式中，A為鋼軌截面積；R為鋼軌等效半徑。

圖2 鋼軌等效圓柱體截面示意圖Fig.2 Equivalent circular cylinder conductor of rail

分析和計(jì)算任何電磁場(chǎng)問(wèn)題的基本出發(fā)點(diǎn)都是麥克斯韋的電磁場(chǎng)方程組。將麥克斯韋電磁場(chǎng)基本方程轉(zhuǎn)換為相量形式

式中，μ 為磁導(dǎo)率；σ 為電導(dǎo)率；J 為導(dǎo)體電流密度。

將式（2）在如圖3所示單位圓柱體模型下展開(kāi)可得

圖3 單位圓柱形導(dǎo)體模型Fig.3 Model of unit circular cylinder conductor

式（3）為圓柱型導(dǎo)體在時(shí)變電磁場(chǎng)中考慮趨膚效應(yīng)時(shí)的求解方程，是典型的貝塞爾方程。對(duì)于此類方程的求解過(guò)程需要對(duì)實(shí)際應(yīng)用中的邊界條件展開(kāi)討論并采用貝塞爾函數(shù)求解獲得。

2.2 導(dǎo)體表面疊加階躍電壓信號(hào)的趨膚效應(yīng)計(jì)算

考慮以下情況：階躍電壓信號(hào)疊加在圓柱形導(dǎo)體表面：

（1）在 t=0時(shí)刻：根據(jù)階躍信號(hào)定義可知：dU/dt→∞，即信號(hào)角頻率趨于無(wú)窮大。根據(jù)趨膚效應(yīng)定義

式中，Js為導(dǎo)體表面電流密度；ρ為導(dǎo)體電阻率；ω為電流角頻率，ω=2πf；d為電流與導(dǎo)體表面間的距離。

此時(shí)導(dǎo)體除表面電流密度以外，內(nèi)部電流密度為零。因此，此時(shí)導(dǎo)體表面電流密度可采用Js=σEs=σU；計(jì)算獲得由于導(dǎo)體內(nèi)部電流趨于無(wú)窮小，因此導(dǎo)體電流為無(wú)窮小。在實(shí)際計(jì)算中，可以近似認(rèn)為It=0=0。

在 t→∞時(shí)刻：導(dǎo)體內(nèi)部任意位置電流密度等于Js。因此，此時(shí)電流大小為：It=∞=πR2Js。

以上兩種情況可以視為電流計(jì)算的邊界條件。

計(jì)算導(dǎo)體表面疊加單位階躍電場(chǎng)強(qiáng)度E時(shí)導(dǎo)體電流密度J(r，t)關(guān)于時(shí)間t與半徑r的解為

式中，J0為零階第一類貝塞爾函數(shù)；J1為一階第一類貝塞爾函數(shù)；αn為零階第一類貝塞爾函數(shù)的根。

圓導(dǎo)體內(nèi)電流可由式（5）積分獲得

將式（5）代入式（6）可得圓柱形導(dǎo)體中的時(shí)域電流解為

在t =0與t =∞間的時(shí)間電流大小可以通過(guò)暫態(tài)內(nèi)電阻與內(nèi)電感計(jì)算獲得。

一種暫態(tài)內(nèi)電阻定義方法為通過(guò)功率損失定義[11]

另一種適用于時(shí)間較短時(shí)的電阻定義為

導(dǎo)體內(nèi)電感可以采用內(nèi)部磁鏈計(jì)算

式中，ψint為導(dǎo)體內(nèi)部磁鏈；Ik為半徑為r的圓形積分路徑內(nèi)的電流；IR為圓導(dǎo)體內(nèi)電流。

將式（7）~式（10）代入式（11）計(jì)算獲得導(dǎo)體電壓，并令 Js=Uσ =1·σ，即穩(wěn)態(tài)時(shí)電壓為 1V時(shí)所產(chǎn)生導(dǎo)體表面電流密度；為使計(jì)算結(jié)果準(zhǔn)確，本文在0～0.000 1s內(nèi)，采用Rint2作為內(nèi)電感計(jì)算，部分結(jié)果見(jiàn)表1；在0.000 1～1s內(nèi)，采用Rint1作為內(nèi)電感計(jì)算，部分結(jié)果見(jiàn)表2。

表1 時(shí)間較短時(shí)電壓計(jì)算結(jié)果Tab.1 Voltage calculation results with short time

表2 時(shí)間較長(zhǎng)時(shí)電壓計(jì)算結(jié)果Tab.2 Voltage calculation results with long time

由表1與表2可知：計(jì)算結(jié)果除在初始時(shí)刻與1V有一定誤差，其他時(shí)刻基本符合。該誤差可能由內(nèi)電感定義造成。并由式（11）可知：導(dǎo)體疊加電壓信號(hào)幅值與電流響應(yīng)幅值為正比關(guān)系。因此，導(dǎo)體疊加階躍電壓信號(hào)所產(chǎn)生的電流響應(yīng)可以采用式（7）進(jìn)行計(jì)算。

2.3 階躍電壓級(jí)數(shù)分析方法

根據(jù)上節(jié)分析，可知對(duì)于鋼軌疊加單位階躍電壓信號(hào)所產(chǎn)生電流響應(yīng)可用式（7）表示。因此，對(duì)電壓信號(hào)的電流響應(yīng)，本文采用信號(hào)的階躍電壓級(jí)數(shù)逼近方法，其步驟如下：

（1）將對(duì)鋼軌施加的電壓信號(hào)采用階躍函數(shù)的級(jí)數(shù)逼近方法進(jìn)行表示。設(shè)在離散采樣的形式下，電壓信號(hào) f(t)可以用下式離散表示：f[n]=f(nT)，T為采樣時(shí)間間隔。

將電壓離散信號(hào) f[n]用單位階躍信號(hào) U(t)級(jí)數(shù)進(jìn)行逼近

式中，當(dāng)n=0時(shí)，f[(n-1)T]=0。

由f[n]的級(jí)數(shù)表達(dá)式可知n-1，f[n]可表示為一系列階躍電壓分量的疊加形式，而單位階躍電壓U(nT)的系數(shù)為前后兩次采樣電壓值的差，即{f(nT)-f[(n-1)T]}。

（2）對(duì)鋼軌施加單位階躍電壓信號(hào)U(t)時(shí)所產(chǎn)生的電流響應(yīng)為 I(t)，即式（7），所以單位階躍電壓U(nT)所對(duì)應(yīng)的電流響應(yīng)為I(t-nT)。因此，階躍電壓分量{f(nT)-f((n-1)T)}U(nT)所對(duì)應(yīng)的電流響應(yīng)可表示為{f(nT)-f((n-1)T)}I(t-nΤ)。

（3）計(jì)算電壓信號(hào)用階躍級(jí)數(shù)表示的逼近式中的每個(gè)階躍電壓分量所對(duì)應(yīng)的電流響應(yīng)，并將這些電流響應(yīng)進(jìn)行時(shí)域內(nèi)的組合，即電流響應(yīng)的時(shí)域解：

該方法具體計(jì)算過(guò)程如圖4所示。

圖4 階躍電壓級(jí)數(shù)分解Fig.4 Analysis of step series of voltage method

2.4 頻域電流計(jì)算驗(yàn)證

本文采用穩(wěn)態(tài)頻域計(jì)算方法對(duì)該算法進(jìn)行驗(yàn)證。驗(yàn)證方法為假設(shè)圓導(dǎo)體上疊加高頻電壓，計(jì)算其頻域阻抗，從而獲得其電流響應(yīng)。圓柱形導(dǎo)體頻域阻抗[1,9,20,21]為

其中

因此導(dǎo)體電流可以通過(guò)下式計(jì)算獲得

采用不同頻率單位交流電壓對(duì)階躍級(jí)數(shù)算法進(jìn)行計(jì)算。測(cè)試中，采樣率為 10kHz，計(jì)算結(jié)果如圖5～圖8所示，與頻域電流計(jì)算結(jié)果對(duì)比見(jiàn)表3，可以發(fā)現(xiàn)：在高頻情況下，該算法與頻域計(jì)算方法基本符合，這證明了該算法的有效性。

圖5 100Hz電流響應(yīng)Fig.5 Current response of 100Hz excitement

圖6 300Hz電流響應(yīng)Fig.6 Current response of 300Hz excitement

圖7 700Hz電流響應(yīng)Fig.7 Current response of 700Hz excitement

圖8 1kHz電流響應(yīng)Fig.8 Current response of 1kHz excitement

表3 電壓階躍級(jí)數(shù)計(jì)算與頻域阻抗計(jì)算結(jié)果比較Tab.3 Comparison between calculation results of step series voltage method and frequency domain impedance method

3 地鐵供電系統(tǒng)仿真及分析

本文采用文獻(xiàn)[9]所提供實(shí)驗(yàn)參數(shù)對(duì) Matlab/Simulink模型進(jìn)行測(cè)試。該實(shí)驗(yàn)對(duì)倫敦地鐵進(jìn)行了短路測(cè)試，短路點(diǎn)距離牽引變電站 2.882km，采樣率為10kHz。計(jì)算結(jié)果與測(cè)試結(jié)果對(duì)比見(jiàn)表4。

表4 電壓階躍級(jí)數(shù)計(jì)算與實(shí)測(cè)結(jié)果比較Tab.4 Comparison between calculation results of step series voltage method and test results

從以上結(jié)果可知，該模型很好地反應(yīng)了鋼軌在短路時(shí)趨膚效應(yīng)的暫態(tài)特性，計(jì)算初始時(shí)刻存在一定誤差，為仿真計(jì)算與實(shí)際實(shí)驗(yàn)中開(kāi)關(guān)暫態(tài)特性有所不同導(dǎo)致。

4 結(jié)論

本文首先對(duì)鋼軌等效圓柱形導(dǎo)體疊加階躍電壓信號(hào)這一情況展開(kāi)分析，采用信號(hào)的階躍級(jí)數(shù)分析方法建立鋼軌趨膚效應(yīng)模型，并用頻域計(jì)算對(duì)其進(jìn)行了驗(yàn)證。最后，將該方法采用S-Function模塊在Matlab/Simulink仿真實(shí)現(xiàn)，并與城市軌道交通短路試驗(yàn)實(shí)測(cè)數(shù)據(jù)進(jìn)行對(duì)比，取得了良好的效果。

[1]衛(wèi)志農(nóng), 鄭玉平, 何樺. 基于小波變換的地鐵遠(yuǎn)方短路電流分析[J]. 電力系統(tǒng)及其自動(dòng)化學(xué)報(bào). 2003,15(6): 5-8.Wei Zhinong, Zheng Yuping, He Hua. The analyzing of distant short-circuit fault current in rapid transit system s based on wavelet transform[J]. Proceedings of the CSU-EPSA, 2003, 15(6): 5-8.

[2]Chang C S, Kumar S, Liu B, Khambadkone A,Real-time detection using wavelet transform and neural network of short-circuit faults within a train in DC transit systems[J]. IEE Procedings-Elects. Power Applications, 2001, 148(3): 251-256.

[3]Yu L, He J H, Yip H T, et al. A new DC feeder protection based on wavelet transform[C]. Powertech,2009 IEEE Bucharest, 2009.

[4]張秀峰, 茍紹波. 利用小波變換檢測(cè)地鐵饋線故障[J]. 西南交通大學(xué)學(xué)報(bào), 2002, 27(5): 548-552.Zhang Xiufeng, Gou Shaobo. Fault detection of subway feed line by means of wavelet transform[J].Journal of Southwest Jiaotong University, 2002, 27(5):548-552.

[5]孔瑋, 崔運(yùn)海, 秦立軍,等. 直流牽引系統(tǒng)饋線微機(jī)保護(hù)裝置[J]. 電力系統(tǒng)自動(dòng)化, 2005, 29(1): 92-96.Kong Wei, Cui Yunhai, Qin Lijun, et al. Research and development of microcomputer-based feeder protection for dc traction power supply system[J].Automation of Electric Power Systems, 2005, 29(1):92-96.

[6]張勛, 李夏青. 直流牽引供電系統(tǒng)的建模以及參數(shù)分析[J]. 北京石油化工學(xué)院學(xué)報(bào). 2009, 19 (1): 5-9.Zhang Xun, Li Xiaqing. Modeling and parameter analysis of DC traction power supply system[J].Journal of Beijing Istitute of Petro-Chemical Technology, 2009,19 (1): 5-9.

[7]張海波, 張穎. 直流牽引供電系統(tǒng)穩(wěn)態(tài)短路分析[J].城市軌道交通研究, 2010(3): 60-65, 68.Zhang Haibo, Zhang Ying. Analysis of steady-state short-circuit in DC traction power supply system[J].Urban Mass Transit, 2010(3): 60-65, 68.

[8]黃維軍, 周密, 王瀅. 直流牽引網(wǎng)穩(wěn)態(tài)短路電流計(jì)算[J]. 電氣化鐵道, 2010(1):39-41.Huang Weijun, Zhou Mi, Wang Ying. DC traction net steady state short circuit current calculation[J].Electric Railway, 2010(1): 39-41.

[9]Brown J C, Allan J, Mellitt B. Calculation of remote short circuit fault currents for DC railways[J]. IEE Proceedings-B. I992, 139(4): 289-294.

[10]Carslaw and Jaeger. Conduction of heat in solids[M].1st ed. London: Oxford University Press, 1947.

[11]Edward J T, Thomas T H, Fullerton H P, Transient resistance of conductors[J]. IEEE Transactions on Power Apparatus and Systems, 1968, 87 (2): 455-462.

[12]孔瑋. 城市軌道交通直流牽引系統(tǒng)故障分析及若干問(wèn)題的研究[D]. 北京: 華北電力大學(xué), 2005.

[13]Hill R J, Carpenter D C. Modelling of nonlinear rail impedance in ac traction power systems[J]. IEEE Transactions on Power Delivery, 1991, 6(4): 1755-1761.

[14]Hill R J, Carpenter D C, Determination of rail internal impedance for electric railway traction system simulation[J]. IEE Proceedings-B 1991, 138 (6): 311-321.

[15]Giacoletto L J, Frequency-and-time-domain analysis of skin effects[J]. IEEE Transactions on Magnetics,1996, 32(1): 220-229.

[16]Silvester P, modal theory of skin effect in single and multiple turn coils[J]. IEEE Transactions on Power Apparatus and Systems, 1970, 89(1): 29-34.

[17]Silvester P. AC resistance and reactance of isolated rectangular conductors[J]. IEEE Transactions on Power Apparatus and Systems, 1967, 86 (6):770-774.

[18]Silvester P, The accurate calculations of skin effect in conductors of complicated shape[J]. IEEE Transactions on Power Apparatus and Systems, 1968, 87 (3):735-741.

[19]Yu L, He J H, Hu J, et al. Accurate track modeling for fault current on DC railways based on MATLAB/Simulink[C]. IEEE Power and Energy Society General Meeting, 2010.

[20]諸葛向彬. 高頻載流圓柱形導(dǎo)體的電阻計(jì)算[J].杭州大學(xué)學(xué)報(bào), 1995(5): 74-77.Zhuge Xiangbing. The calculation of resistance of high frequency current cylinder conductor[J]. Journal of Hangzhou University, 1995(5): 74-77.

[21]崔運(yùn)海. 直流牽引饋線微機(jī)保護(hù)裝置的研制[D]. 北京：華北電力大學(xué), 2005.