周 蠡 彭 迎 李 靖 陳 然 張毅軍 但林烽 黃雄峰
(1. 國網湖北省電力公司經濟技術研究院, 武漢 430077; 2. 國網湖北省電力公司, 武漢 430077; 3. 武漢供電設計院有限公司, 武漢 430030; 4. 三峽大學 電氣與新能源學院,宜昌 443002)
近年來城市輸電網中地下電纜逐漸取代架空線路,其運行安全受到廣泛關注.電纜導體溫度過高將引起絕緣材料加速老化,影響電纜使用安全、縮短電纜使用壽命[1].對電纜導體溫度實時監測是提高運維效率、保障輸電安全的重要手段.
目前主要是通過在電纜內部埋設感溫光纖等來對導體溫度進行直測,此種手段會破壞電纜結構.可由較易測量的表面溫度間接反演出導體溫度.借鑒于有限元法在電纜溫度場計算中的有效應用[2-3],利用有限元法可實現多種環境下電纜溫度場的分析計算[4-5],從而獲取不同載荷、邊界條件下的表面溫度和導體溫度的對應關系,為利用表面溫度間接估計導體溫度提供了可能,從而避免了對電纜結構的破壞.該方法的不足之處是在建立電纜的溫度場有限元模型時需要知道電纜內部不同材料的熱參數值,而廠家往往不能提供這些參數值.即使可以提供,這些值在電纜的長期運行中也會由于老化等原因發生很大變化[6].
本文提出了一種電纜導體溫度間接測量方法,通過對電纜內部材料的熱參數值進行辨識從而校正了電纜的溫度場有限元模型,在校正后的溫度場模型基礎上得到正確的表面溫度與導體溫度的對應關系,再結合表面溫度的測量值間接得出導體溫度.
測量原理如圖1所示.針對某種類型電纜,根據其結構參數、材料熱參數、電纜表面的邊界條件、負荷參數、初始條件,利用有限元法對電纜溫度場進行計算,根據分析結果歸納出不同載荷、邊界條件下表面溫度與導體溫度值的對應曲線,再結合外表皮溫度的監測值反演出導體溫度.

圖1 測量過程原理圖
電纜溫度場的有限元計算是反演實現的關鍵,溫度場有限元模型中負荷參數、邊界條件、結構參數一般已知,但某些熱參數往往無法由電纜規格書上準確獲得,從而影響了模型的準確性.本文基于測量點溫度的實測值與計算值建立目標函數,通過對兩者之間的誤差進行處理反辨識出相關熱參數,將辨識后的熱參數代替原先的初始熱參數載入到溫度場有限元模型中,從而實現了對模型的校正.
本文以典型單芯110 kV高壓電力電纜為例,其截面如圖2所示,E(x0,y0)為測量點.由于電纜線路與其直徑相比可認為無限大,因此電纜溫度場可按二維溫度場進行分析和計算.

圖2 電纜截面圖
電纜中,含內部熱源區域的控制方程為:
(1)
式中,T為點(x,y)處的瞬態溫度(℃);t為過程進行的時間(s);λ為導熱系數[W/(m2·℃)];S為體積生熱率(W/m3);φ為材料的密度(kg/m3);c為材料的比熱[J/(kg·℃)].
無熱源區域的溫度控制方程為:
(2)
電纜外邊界處為第三類邊界條件,邊界方程為:
(3)
式中,α指電纜的對流換熱系數[W/(m2·K-1)];Tf指空氣溫度(K).
電纜的初始條件指傳熱過程開始時刻計算(t=0)電纜內部任一點的溫度,即
T(x,y,0)=T0(x,y)(4)
采用三角形單元計算有限元溫度場,結合Galerkin法,列出電纜平面溫度場有限元方程為:




l=i,j,m;e為平面電纜溫度場的定義域.
對上式進行合并得:
(9)
式中,通過三角形單元面積積分合并可以分別求出K、N和P.
利用Grank-Nicolson差分格式計算t時刻和t-Δt時刻電纜暫態溫度值,即得出電纜暫態溫度場計算公式:

運用加權法處理三類邊界條件可得出線積分方程.利用迭加法和消去法求解,即可求出不同時刻電纜溫度場內各點的溫度值.
電纜的瞬態熱特性取決于以下5個參數:1)負荷電流I;2)初始溫度T0;3)電纜表面的邊界條件B;4)電纜結構參數G;5)電纜材料的熱參數P.在以上參數中,除了熱參數P,I、T0、B、G都可以獲得.因此,電纜表面處點E(x0,y0)的有限元計算溫度值TE可用熱參數P的函數表達:
TE=TE(t,x1,x2,…,xj,…,xn)(11)
式中,x1,x2,…,xj,…,xn是熱參數集P中的元素,例如導熱系數λ、比熱c.
另一方面,E點的溫度可以用設備直接測量,測量值用Tm(t)表示.令x=(x1,x2,…,xn),參數辨識即通過測量值Tm(t)來擬合出熱參數集合x中的所有未知元素,本文選擇最優化方法對目標函數F進行最小化以辨識出x的最優解,其中F為:
(12)
式中,x為所有熱參數集合,x=(x1,x2,…,xn);N是實測值總數;Δt測量時間間隔,Δt=TD/N,其中TD為測量時間.
結合以上分析,本文是通過選擇最優化方法對目標函數F進行最小化來辨識熱參數,是把溫度場有限元計算程序作為最優化方法的子程序,從而將有限元程序與主程序優化算法相銜接,流程如圖3所示.

圖3 參數辨識過程
很多最優化方法可以用來求解式(12),本文選擇最速下降法來對目標函數F進行最小化,優化過程如下.
步驟1:給定初點x(1)、允許誤差ε=0.1,置k=1;
步驟2:計算搜索方向d(k),其中
d(k)=-(13)
式中,目標函數F相對x中的第j個參數xj的偏微分可由式(14)求得:
(14)
其中TE(jΔt,x)相對xj的偏微分由式(15)求得:
(15)
式中,xe=(x1,x2,…,xj+ε,…,xn);xd=(x1,x2,…,xj-ε,…,xn);TE(jΔt,xe)和TE(jΔt,xd)的值是由有限元溫度場計算中得到的;
步驟3:若‖d(k)‖≤ε,停止計算;否則,從x(k)出發,沿d(k)進行一維搜索,求λk,使得
(16)
式中,λk是從x(k)出發沿方向d(k)進行一維搜索的步長;
步驟4:令x(k+1)=x(k)+λkd(k),置k=k+1,轉步驟2.
利用有限元法對實驗室中空氣敷設狀態中的標稱截面為1×400 mm2的YLJW02 64/110 kV電纜進行溫度場計算,在電纜表面選定監測點E,將此點的溫度計算值與利用熱電偶實測得到的溫度值相結合,采用最速下降法對相關熱參數進行辨識.利用辨識后的熱參數對電纜溫度場模型進行校正,計算后得出表面溫度和導體溫度的對應曲線.
如圖4所示,電纜溫度測量實驗裝置由以下部分組成:開關柜、調壓器、升流器、電流互感器、熱電偶測溫儀、數據接收轉換器、電腦等.其中開關柜控制調壓器的開斷,調壓器用于升流器的加壓;升流器選用穿心變壓器,用于為電纜加載電流,其輸入端的三相電壓為380 V,單臺升流器提供的最大電流為600 A,實驗將兩臺升流器并聯使用.

圖4 電纜模型
提供長5 m、截面積400 mm2的110 kV電纜.在電纜表面放置熱電偶用以測量表面溫度,另外在電纜本體鉆孔,孔洞直徑1~2 mm,深度約5 mm,在電纜導體處放置熱電偶用以測量導體溫度.熱電偶通過光纖連接到數據接收轉換器,轉換器將模擬信號轉換為數值信號后傳送到電腦用以顯示.
1)電纜的敷設條件
樣品電纜敷設于3×3 m2的實驗室內,實驗室溫度在23℃左右小范圍波動,門窗關閉后無明顯空氣流動,電纜與室內空氣自然對流來進行散熱.電纜表面對流換熱系數h可通過經驗公式(17)計算得到[7-8],為7.371.電纜的初始溫度T0與環境溫度相同,為20℃.
h=7.371+6.43v0.75(17)
式中,v為空氣流動速度(m/s).
2)電纜的參數
實驗室樣品為標稱截面1×400 mm2的YLJW02 64/110 kV電纜,電纜結構參數見表1.

表1 電纜結構參數
樣品電纜中有5種材料:銅(導體)、鋁(皺紋鋁護套)、交聯聚乙烯(絕緣層)、半導體材料(導體屏蔽層與絕緣屏蔽層)、橡膠與聚乙烯的混合材料(外護套),如圖5所示.對應這5種材料,結合式(1)、式(2)可知溫度場模型中需要15種熱參數:λcon(導體熱導率)、φcon(導體密度)、ccon(導體比熱);λAl(鋁護套熱導率)、φAl(鋁護套密度)、cAl(鋁護套比熱);λins(絕緣層熱導率)、φins(絕緣層密度)、cins(絕緣層比熱);λjac(外護套熱導率)、φjac(外護套密度)、cjac(外護套比熱);λscr(屏蔽層熱導率)、φscr(屏蔽層密度)、cscr(屏蔽層比熱).

圖5 電纜材料示意圖
在對電纜進行溫度場計算中,導體、鋁護套、屏蔽層的熱參數在電纜運行中變化較小,通過產品說明書得到,見表2.至于其它的6個熱參數(λins、φins、cins;λjac、φjac、cjac),其在運行中受老化等因素影響較大因此不能準確獲得,可先根據經驗數據設定它們的初始值[9],再通過后續步驟對它們進行修正,設定初始值見表3.

表2 已知電纜熱參數

表3 未知電纜熱參數的設定初始值x(1)
3)熱源計算
對樣品電纜加載400 A、50 Hz的零相位電流作為載荷,利用貝塞爾函數計算熱源值[10].經計算求得:導體損耗Q1=48.351 W/m;絕緣損耗Q2=0.145 W/m;金屬護套損耗Q3=0.614 W/m.
4)溫度場計算
根據實際參數建立電纜幾何模型,如圖6(a)所示.對其進行網格剖分后得到有限元模型,如圖6(b)所示.

圖6 電纜模型
在溫度場計算部分,設置求解域的熱參數,包括空氣溫度、電纜各層的熱導系數.在空氣邊界處設置溫度22.5℃作為邊界條件,將之前計算的損耗值作為載荷施加到熱-流體模型中來,仿真時長為6 h,求解后得到溫度場分布云圖,如圖7所示.

圖7 溫度分布云圖
5)參數辨識
提取監測點E(0,34.6 mm)的計算溫度值繪制曲線;同時利用熱電偶測溫儀對E點的溫度以及空氣溫度進行實時監測并繪制溫度曲線,監測時間也為6 h,溫度采集間隔為2 min.曲線在圖8中標示.由圖8可知在初始參數x(1)的條件下計算值與實測值有一定的偏差,需要對初始參數進行辨識得到優化值.

圖8 表面溫度實測值與計算值的比較
根據E點的溫度測量值和相應有限元計算結果建立目標函數,將相關熱參數的辨識問題化為該函數的最小值問題.基于最速下降法編寫程序對該問題進行求解,目標函數在大約40 min后完成收斂,迭代34次,如圖9所示.6個相關熱參數最后收斂于:λins=0.312、φins=896.706、cins=2 703.12;λjac=0.268、φjac=948.3、cjac=1 315.8.

圖9 目標函數的收斂過程
為對優化結果進行驗證,將優化過的參數代入有限元程序進行計算,得出的結果提取后示于圖8中.
另外,為驗證優化后參數的有效性,本文針對負載電流800 A、環境溫度29℃的工況,將優化參數代入有限元程序進行計算并提取E點的溫度值繪制曲線,如圖10所示.實測值與計算值的最大溫差為1.1℃,誤差在可接受范圍內,證明了優化后參數的有效性.

圖10 800 A時表面溫度實測值與計算值的比較

圖11 表面溫度與導體溫度的對應曲線
6)對應關系的確立
通過以上分析,可知校正后的有限元模型具有較高的準確性,在此基礎上進行溫度場計算得到的表面溫度與導體溫度的對應曲線可用于導體溫度的間接測量.本文提取出了負載電流400 A、環境溫度23℃時表面溫度與導體溫度的對應曲線,如圖11所示.測量到表面溫度后即可在曲線上搜尋對應的導體溫度.
針對電力電纜導體溫度的測量,以通過表面溫度測量間接得到纜芯溫度為思路,基于有限元法對電纜溫度場進行分析計算,將分析結果與電纜表面溫度測量值結合從而反演出纜芯溫度值.針對電纜中某些材料熱參數難以標定從而影響溫度場模型準確性的問題,本文通過對電纜內部材料的熱參數值進行辨識從而校正了電纜的溫度場有限元模型,在校正后的溫度場模型基礎上得到正確的表面溫度與導體溫度的對應關系,再結合表面溫度的測量值間接得出導體溫度.
[1] Vu K, Begovic M, Novosel D. Use of Local Measurements to estimate Voltage Stability Margin[J]. IEEE Transactions on Power System, 1999, 24(3).
[2] Rachek, M. Nait Larbi. S. Magnetic Eddy-Current and Thermal Coupled Models for the Finite- Element Models for the Finite-Element Behavior Analysis of Underground Power Cables [J]. IEEE Trans on Magnetics, 2008,44(12):4739-4746.
[3] 王有元,陳仁剛,陳偉根,等.有限元法計算地下電纜穩態溫度場及其影響因素[J]. 高電壓技術,2008,34(12):3086-3092.
[4] 王東濤,高 丹.基于組態王的動力電纜溫度在線監測系統[J]. 中國電力,2006,39(4):79-82.
[5] 黃詩雅,吳 勇,李 磊.土壤直埋敷設單芯電力電纜溫度場與載流量計算[J]. 武漢大學學報,2014,47(4):502-505.
[6] H.J. Li. Estimation of thermal parameters and prediction of temperature rise in crane cables [J]. IEEE Proc-Gener. Transm Distib, 2004, 151(3): 355-360.
[7] Mitchell J K. Temperature Distribution Around Buried Cables [J]. IEEE Trans. Power Apparat Syst, 1979,98(4):1158-1166.
[8] Tarasiewicz E, Kuffel E, Grzybowski S. Calcultion of Temperature Distribution Within Cable Trench Backfill and Surrounding Soil [J]. IEEE Trans. Power Apparat Syst, 1985,3(8):1973-1978.
[9] Buonanno G. Effect of Radiative and Convective Heat Transfer on Thermal Transients in Power Cables[J]. IEE Proc. Gener. Transm. Distrib, 1995, 142, pp. 436-444.
[10] 梁永春,柴進愛,李彥明,等.有限元法計算交流電纜渦流損耗 [J]. 高電壓技術,2007,33(9):196-199.