張?zhí)祚?,梁志堅,朱?/p>
(1. 廣西大學(xué) 電氣工程學(xué)院,廣西 南寧 530004;2. 國網(wǎng)能源研究院有限公司,北京 102209)
光伏系統(tǒng)及變電站作為新能源發(fā)電、系統(tǒng)能量轉(zhuǎn)換與傳輸?shù)闹匾M成部分,受到了極大關(guān)注,其安全及可靠性顯得尤為重要[1-2]。然而,光伏變電站通常安裝在露天區(qū)域,尤其是在雷電災(zāi)害頻發(fā)的地區(qū),雷擊已經(jīng)成為威脅系統(tǒng)設(shè)備的主要安全隱患之一[3],因此需要對光伏變電站采取防雷措施加以保護。近年來,氧化鋅避雷器以其殘壓低、無續(xù)流、動作時延小、通流容量大等優(yōu)點,在高壓及特高壓輸變電系統(tǒng)中得到廣泛應(yīng)用[4]。但是,隨著直流輸電工程的快速發(fā)展與線路電壓等級的提高,避雷器出現(xiàn)了損耗較大、溫升較高等一系列問題[5]。而溫升的大小和分布對氧化鋅閥片及避雷器的使用壽命有重大影響,關(guān)乎整個電力系統(tǒng)安全運行[6]。因此,研究避雷器溫度場分布,特別是在老化等故障情況下,就顯得尤為重要[7]。
目前,很多學(xué)者對電力設(shè)備多物理場進行分析,主要的方法總結(jié)為:試驗法、解析法和有限元法等[8]。其中,試驗方法需要加工一個實際模型,導(dǎo)致成本過高。此外,一些特殊工況的實驗條件也過于苛刻,因此,不適用于模型的初步研究。文獻[9]提出場路耦合法對500 kV避雷器進行溫度場計算。光伏變電站配備的220 kV避雷器采用的是兩節(jié)元件結(jié)構(gòu),在高、低壓段間的電場分布發(fā)生了躍變現(xiàn)象。因此采用常規(guī)電磁解析方法計算損耗,再結(jié)合傳熱學(xué)方法對模型溫度場進行分析的方法,并不是溫度場分析的最優(yōu)選擇。此外,場路耦合方法雖然從物理概念上描述比較清晰,但是無法考慮流熱耦合計算中的集膚效應(yīng)、渦流效應(yīng)和流體流動等影響,難以全面而準確地描述各個場的具體分布[10]。目前,有限元法采用對稱正定的稀疏矩陣,可以在綜合考慮對流及輻射換熱效應(yīng)的基礎(chǔ)上,直接求解流固耦合方程[11],因此能夠比較全面地對電力設(shè)備的多物理場進行耦合分析。文獻[12]采用有限元法計算10 kV氧化鋅避雷器的溫升,但是計算的是二維簡化模型,忽略了熱量在垂直方向的傳導(dǎo)過程。文獻[13]采用商業(yè)軟件計算動車組車頂避雷器的溫升,但是忽略了熱輻射和空氣對流換熱對溫度計算的影響。雖然有限元法得到了成功的應(yīng)用,但是該方法是建立在區(qū)域變分原理和泛函計算的基礎(chǔ)上,一旦網(wǎng)格剖分結(jié)束后就不易進行修正,因此不具備較強的自適應(yīng)能力[14]。特別是在前處理過程中的網(wǎng)格劃分十分繁瑣,會占用大量時間。同時,網(wǎng)格扭曲或畸變現(xiàn)象導(dǎo)致計算精度差、收斂困難等問題也亟待解決[15]。因此迫切需要一種能保障計算收斂速度,滿足計算精度,并且不依賴網(wǎng)格劃分的新方法。
近年來,一種只需節(jié)點信息,無須劃分網(wǎng)格即可構(gòu)造高階連續(xù)函數(shù),具有較高計算精度的無網(wǎng)格法得到廣泛應(yīng)用[16-17]。但是,無網(wǎng)格法采用徑向基函數(shù)法構(gòu)造函數(shù)的近似表達式,導(dǎo)致計算得到的一般是節(jié)點參數(shù)而非真實節(jié)點值。因此函數(shù)不滿足Kronecker特性,不能像有限元法一樣便捷地施加邊界條件。這也使得針對溫度場計算邊界的處理十分麻煩[18]。而有限元法處理邊界問題比較方便,因此將有限元法與無網(wǎng)格法結(jié)合起來,可以較好地解決電力設(shè)備的多物理場耦合計算問題。但是,在有限元區(qū)和無網(wǎng)格交界面處的節(jié)點值是真實的,而無網(wǎng)格法的節(jié)點值是虛假的,無法直接匹配。因此,文獻[19]引入由有限元與無網(wǎng)格法形函數(shù)組成的界面單元,雖然保證了方程求解的連續(xù)性,但忽略了交界面函數(shù)導(dǎo)數(shù)的連續(xù)性,而對函數(shù)及其空間導(dǎo)數(shù)的近似計算恰恰是無網(wǎng)格法求解偏微分方程的關(guān)鍵所在。目前,考慮魯棒性的針對電力設(shè)備有限元與無網(wǎng)格法耦合計算的文獻還未見公開報道。
本論文研究如何將無網(wǎng)格方法應(yīng)用到電力設(shè)備的流-熱耦合計算中,重點研究有限元求解的流體場與無網(wǎng)格求解的溫度場在交接界面的數(shù)據(jù)傳遞。首先建立1臺經(jīng)過簡化處理的某光伏變電站220 kV氧化鋅避雷器三維流體場模型。采用商業(yè)有限元軟件Ansys分別對其額定工作狀態(tài)和老化故障下的流速分布進行計算,得到各壁面的對流換熱系數(shù)。其次,通過在有限元法與無網(wǎng)格法交界面之間構(gòu)造過渡區(qū)域,采用徑向基函數(shù)法構(gòu)造無網(wǎng)格法的形函數(shù),推導(dǎo)出基于無網(wǎng)格法的傳熱學(xué)控制方程。在此基礎(chǔ)上,建立了2種工況下的無網(wǎng)格溫度場模型,得到了相應(yīng)的溫度分布。然后,將無網(wǎng)格-有限元法間接耦合計算結(jié)果與純粹有限元流體-溫度直接耦合仿真結(jié)果進行對比。最后通過實驗,在證明了本文耦合計算方法準確性的同時,還能發(fā)現(xiàn)無網(wǎng)格計算的速度遠高于有限元法。本文所述方法也可以為其他電力設(shè)備多物理場耦合分析提供一定的參考。
圖1 避雷器三維仿真模型Fig. 1 Three dimensional simulation model of surge arrester
表1 避雷器主要參數(shù)Table 1 Major parameters of surge arrester
表2 材料主要參數(shù)Table 2 Major parameters of materials
為了便于分析,對避雷器做以下假設(shè)。
(1)在絕緣套管與空氣之間為大空間輻射換熱,因此將氣體和固體交界面設(shè)置為無滑移邊界[20];
(2)流體被視為不可壓縮氣體,忽略浮力影響,同時忽略溫度對空氣粘度的影響,計算得到的避雷器壁面的對流換熱系數(shù)均為平均值;
(3)假設(shè)初始溫度為室溫,忽略接觸熱阻,所有介質(zhì)均可視為各向同性介質(zhì)。
無網(wǎng)格法與有限元法結(jié)合的氧化鋅避雷器多物理場耦合分析流程如圖2所示。首先,建立避雷器三維流體場有限元模型,并有如下簡化處理:忽略避雷器中的一些細微結(jié)構(gòu),如不影響溫升計算的退刀槽、螺紋孔、螺栓等;另外將曲率較小弧形的導(dǎo)角均按照直角處理;其次,對模型進行網(wǎng)格劃分,并進行求解;然后,將流體場求解出的對流換熱系數(shù)作為邊界,采用有限元法對避雷器進行共軛傳熱計算;最后,通過在有限元與無網(wǎng)格計算區(qū)域中間建立過渡函數(shù),實現(xiàn)流體場與溫度場的數(shù)據(jù)傳遞。在同樣的激勵與邊界條件的基礎(chǔ)上,利用無網(wǎng)格法再次進行溫度場計算,得出不同工況下避雷器溫升的分布。
圖2 避雷器多物理場耦合計算流程Fig. 2 Process of fluid-thermal coupling analysis
在避雷器流體場有限元計算中,主要的散熱方式為自然冷卻,因此氣體流通模式主要為層流,滿足的流體場控制方程[21]如下。
在溫度場計算中,同時考慮了熱傳導(dǎo)、對流傳熱和熱輻射的影響,其控制方程與邊界條件[23]為
借助準滑移面方法研究RT位移模式非極限狀態(tài)被動土壓力,結(jié)果表明:繞墻頂轉(zhuǎn)動擋土墻非極限狀態(tài)被動土壓力為凹曲線分布,越靠近墻底,被動土壓力增大越快;轉(zhuǎn)動越大,被動土壓力越大,其合力作用點越低,合力作用點在墻高的下三分點以下。
氧化鋅避雷器運行時會通過套管向周圍不斷輻射熱量,計算公式[24]表示為
式中:εt為表面輻射率,對于絕緣套管表面取為0.95,而對于氧化鋅閥片則為0.71;σ為斯忒藩-玻耳茲曼(Stefan-Boltzmann)常數(shù),文中取為5.7×10-8W/(m2?K4);Sto為壁面輻射面積。
最終溫度場有限元離散形式為
式中:KFEA為剛度矩陣;N為由四面體單元組成的形函數(shù)。
圖3為避雷器的網(wǎng)格剖分,可看出網(wǎng)格質(zhì)量良好,符合計算精度要求??偟墓?jié)點數(shù)為1 815 403,網(wǎng)格數(shù)為8 841 828。因此,網(wǎng)格劃分需要占用大量的計算時間和資源,這也是本文所述無網(wǎng)格法優(yōu)勢所在。
圖3 流體計算網(wǎng)格Fig. 3 Meshes of fluid calculation
額定運行與上節(jié)閥片老化故障時避雷器流體場粒子分布如圖4所示??梢钥闯觯?種工況下,流體的分布區(qū)別不大,最大值分別為0.18 m/s和0.181 m/s,都出現(xiàn)在套管外部正上方。
圖4 避雷器流速分布Fig. 4 Flow velocity distribution of surge arrester
求出的對流換熱系數(shù)可以用以下方程表示。
對于避雷器側(cè)壁面,可表示為
對于避雷器頂面和底面,則分別為
式中:Cf為阻力系數(shù);Ra為瑞利數(shù)。
圖5為避雷器額定運行與上節(jié)閥片老化故障時的絕緣套管表面溫度分布情況。從圖5 a)可以看出,正常情況下套管表面最大溫升為21.65 °C,主要位于中間法蘭與上節(jié)底部外套連接處。由于受到重力加速度對氣體流動的影響,溫度整體呈現(xiàn)出在豎直方向的梯度變化。圖5 b)為上節(jié)老化故障情況下的溫度分布,與正常工況對比可以發(fā)現(xiàn),下節(jié)套管表面溫升要低一些,溫升最大為22.25 °C。
圖5 基于有限元法的避雷器絕緣套管溫升分布Fig. 5 Thermal distribution of arrester insulating sleeve based on FEA
圖6為氧化鋅閥片額定運行及老化2種狀態(tài)下的溫度分布情況??梢钥闯?,在正常工況下,閥片各節(jié)均呈現(xiàn)出中間溫度高而兩端溫度低的分布趨勢,最高溫度為31.42 °C。主要是因為避雷器內(nèi)部閥片柱兩端的金屬連接件散熱效果較好,所以溫度值較低。當(dāng)閥片老化故障時,最大溫升為38.45 °C??梢园l(fā)現(xiàn)老化的閥片溫升要明顯降低,故需要特別關(guān)注下節(jié)閥片的溫升變化。
圖6 基于有限元法的氧化鋅閥片溫升分布Fig. 6 Thermal distribution of ZnO valve plate based on FEA
出現(xiàn)這種現(xiàn)象的主要原因是當(dāng)氧化鋅閥片在老化時,電阻片吸收了大量的能量。在較大的熱載荷作用下,內(nèi)部晶界層中的離子活動非?;钴S,增大了離子遷移率,導(dǎo)致閥片內(nèi)部肖特基勢壘壁高度降低,因此對過電壓產(chǎn)生了抑制效果,有效降低了損耗[27]。此外,在長時間的大電流作用下,閥片會受到較大電磁力和熱應(yīng)力的影響,導(dǎo)致壓敏電阻片內(nèi)部結(jié)構(gòu)發(fā)生變化,產(chǎn)生明顯的極性效應(yīng)。因此,損耗分布發(fā)生較大變化引起溫升分布不均勻。
圖7為避雷器額定運行與上節(jié)閥片老化故障時絕緣套管表面溫度分布情況。對比圖5,可以看出二者溫升分布基本相同,最大溫升分別為22.22 °C 和 23.25 °C,相應(yīng)的誤差分別為 2.6% 和4.5%。
圖7 基于無網(wǎng)格法的絕緣套管溫升分布Fig. 7 Thermal distribution of insulating sleeve based on meshless method
圖8為氧化鋅閥片正常及老化2種狀態(tài)下的溫度分布情況。對比圖6,可以看出計算出的溫升分布基本相同。其中,無網(wǎng)格法計算出的閥片正常狀態(tài)工況的最大溫升為31.83 °C,誤差為1.3%。而老化后的最大溫升為38.69 °C,誤差為0.6%,可以證明所用無網(wǎng)格計算方法有足夠的精度。
圖8 基于無網(wǎng)格法的氧化鋅閥片溫升分布Fig. 8 Thermal distribution of ZnO valve plate based on meshless method
表3為流體-溫度耦合計算過程中各個環(huán)節(jié)的用時統(tǒng)計,可以發(fā)現(xiàn),有限元法由于產(chǎn)生大量網(wǎng)格,計算時間遠超無網(wǎng)格法。對比可以發(fā)現(xiàn)采用有限元法計算溫度場消耗的時間大約是采用無網(wǎng)格法的210倍。
表3 計算時間統(tǒng)計Table 3 3 Calculation time
為了驗證避雷器多物理場耦合計算的準確性,按照符合國家標準GB18802.1—2011《低壓電涌保護器(SPD)第1部分:低壓配電系統(tǒng)的電涌保護器性能要求和試驗方法》的規(guī)定來確定裝備的溫升。溫度由mv2000型溫度傳感器來監(jiān)測,檢測點的安裝位置如圖1所示。避雷器持續(xù)工作20 h,采樣時間為10 min,測量的結(jié)果如圖9所示。由圖9可知,各仿真結(jié)果與試驗結(jié)果基本吻合,滿足精度要求。structure surge arrester[J]. High Voltage Engineering, 2012, 38(8):2129-2136.
圖9 溫升實驗結(jié)果Fig. 9 Experimental results of temperature rise
本文提出了一種針對光伏變電站用氧化鋅避雷器溫升計算的新方法,構(gòu)造了斜坡函數(shù),將有限元與無網(wǎng)格法結(jié)合,綜合了2種方法的優(yōu)點。得出如下主要結(jié)論。
(1)描述了無網(wǎng)格-有限元耦合法的基本原理,根據(jù)溫度場的控制方程,推導(dǎo)出了綜合無網(wǎng)格法與有限元的溫度場矩陣方程;
(2)用有限元法分別計算避雷器額定和老化工況下的流速和溫升分布情況,結(jié)果顯示,額定運行時,避雷器溫度場分布較為均勻,而在老化故障下,氧化鋅閥片溫升要明顯降低;
(3)采用無網(wǎng)格-有限元耦合法分析溫度場問題時,不僅能夠保障計算精度,還能夠極大地提高計算效率。此種方法對配置在光伏變電站內(nèi)的電氣設(shè)備數(shù)值計算具有一定的參考價值及意義。