一種變壓器時域場路耦合計算的自適應算法

潘 超 何雯婷 郝文波 徐冰亮

（1. 東北電力大學電氣工程學院 吉林 132012 2. 黑龍江省電力科學研究院 哈爾濱 150030）

1 引言

變壓器直流偏磁條件下的時域磁場計算以穩態電磁場分析為基礎[1-3]，通常直接求解有限元微分方程，數據運算規模龐大，精確性較低[4,5]。文獻[6,7]采用諧波平衡有限元法對單元內磁矢量位的各次諧波同時求解再疊加，結果具有較高的精確性；這種方法原理復雜，占用資源較多[8]。文獻[9]通過計算分解的方法來優化資源，雖然減少了系統資源的占用率，但是計算時間較傳統方法更長。基于電路-磁路耦合的方法能夠反映出變壓器直流偏磁時的勵磁特性，滿足工程分析需要，但是沒有充分考慮漏磁情況，導致準確性降低[10]。場路耦合方法分為直接耦合與間接耦合兩種方式，目前的大部分研究都采用直接方法實現場與路的耦合，將磁場方程與電路方程聯立計算，效率較低[11]。上述方法大都建立二維模型模擬變壓器鐵心磁場分布，具有較高的精度，但二維場或對稱場僅可作為局部場域的近似處理，只能用于特定運行方式或模型特性的分析，無法模擬變壓器直流偏磁條件下的內部磁場分布[12]。這些問題在變壓器等電磁設備的時域電磁計算中尤為突出，并亟需解決。

另一方面，傳統時域計算通常采用單步低階方法，其截斷誤差較大，且由于時間迭代造成誤差累積，易引起求解振蕩。若細化空間網格或縮短步長，雖然可以在一定程度上減小振蕩，但計算規模會大幅度增加[13,14]。研究表明，變壓器瞬態磁場和時域場路計算中均存在計算穩定性、效率及精確性的問題[15,16]。文獻[17]討論了變壓器時域場路計算過程中非線性勵磁、步長選擇及時間離散在直流擾動下的變化情況。但是關于變壓器時域場路耦合計算的穩定性、精度及效率等優化問題，目前罕有文獻研究。

本文結合變壓器直流偏磁問題，推導計算收斂時的穩定域邊界，確定時間步長與計算穩定性的關系，研究變壓器三維時域場路耦合計算效率與精確性的優化方法。根據穩定性分析提出一種自適應變步長算法，在計算效率和精確性方面與定步長算法對比，并通過實驗測量對其計算結果進行驗證。

2 時域場路耦合

利用時域場路耦合方法計算變壓器磁場和等效電路，動態電感和時域電流為關鍵耦合參數，分別由三維棱邊元磁場模型和等效電路模型計算。

建立變壓器磁場模型，不考慮磁滯效應，棱邊有限元法采用矢量磁位A，得到非線性磁場方程

變壓器磁鏈方程

式中，ψ 為磁鏈向量；i為繞組電流向量；LS為靜態電感矩陣，表示激勵電流與磁鏈的穩態關系。

實際上，變壓器勵磁過程中電感與電流均體現出時變非線性的特點[18]。由u=dψ /dt，推導變壓器等效電路的時域響應方程

式中，u為交流電壓向量；UDC為直流電壓；LD為考慮勵磁變化的動態電感矩陣，需要通過磁場模型計算獲得。

基于能量擾動的思想，根據系統能量變化計算動態電感參數。當線圈電流增加δip（0≤δ≤1，單相變壓器中 1≤p≤2）時，將電源總能量與動態電感和電流關聯，得到

若由電流增量δi引起的場量變化為δH、δB，變壓器內部系統的磁場能量增量為

由能量平衡原理，方程（5）和（6）的能量相等，則可計算動態電感LD。

方程（4）的求解可采用四階龍格庫塔法，由tk時刻的線圈電流ik計算tk+1時刻ik+1

式中，h為步長；s1～s4為步長內的分段計算斜率。

利用繞組電流i與動態電感LD即可實現時域場路耦合，其迭代計算原理如圖1所示。

式中，μ 為磁導率；J為電流密度，需要通過電路模型計算得到。

應用格林定理，得伽遼金加權余量方程

式中，Mm為權函數序列，且權函數與基函數相同。

圖1 時域場路耦合Fig.1 Time-domain magnetic field and electrical circuit coupling

3 自適應變步長算法

時域電路模型（4）的離散通解為式中，τ 為時間常數，取決于狀態矩陣的特征值λmax，且τ =1/λmax。

另一方面，在tk時刻式（4）中的泛函可化為常微分方程，得到式中，v為泰勒級數項的組合函數。

i(tk+1) 與 ik+1一致收斂到 h4階，聯立式（8）和式（10）有

式中，G為傳遞函數；χ =hλmax。

時域場路耦合計算收斂所對應的絕對穩定域判據為│G│＜1。確定時域場路耦合計算的穩定邊界為

式中，h0為穩定邊界對應的計算步長。

對時域場路耦合迭代計算進行改進，提出一種自適應變步長算法，其基本原理為

（1）時域場路耦合計算受初值、步長等因素的影響，在起始階段存在過渡過程，為了減少計算耗時，根據狀態矩陣特征值，在計算穩定的前提下增大步長，從而提高效率。

（2）經過某段計算時間直至計算穩定時，在每次磁場求解動態電感之后，由變壓器狀態矩陣計算特征值，根據穩定域邊界確定對應的h0，引入修正因子α(0＜α＜1)，確定下一時刻的計算步長，h=αh0，進行迭代計算。

（3）合理增大磁場計算的時間間隔，即降低磁場求解頻率，增加電路迭代計算次數，w=βh0（80≥w≥40，β取整數），從而保證計算效率與精確性。

定義tk時刻的動態電感變化函數g(LD)，表示在該時間段內動態電感參數的變化速率，用來評價磁場計算的精確性。

式中，max(LD) 為動態電感的最大值；ε 為變化率限值。

4 算例分析

采用時域場路耦合模型計算變壓器直流偏磁問題，編譯四階龍格庫塔法程序求解時域電路模型，利用Ansys軟件建立變壓器的八分之一模型（模型尺寸與實際比例為1∶1），對絕對穩定域邊界和自適應算法進行研究分析。實驗變壓器型號為 BK300，如圖 2所示，具體參數見表 1；鐵心硅鋼片型號為DW360—50，磁化曲線如圖3所示。

圖2 實驗變壓器Fig.2 Transformer for experiment

表1 變壓器參數Tab.1 Parameters of transformer

圖3 鐵心磁化曲線Fig.3 Magnetization curve of the core

4.1 時域場路耦合定步長算法

變壓器一次側空載電流峰值為Im，接入直流源UDC時的直流電流為 IDC。對空載運行 IDC=0、IDC=25%Im、IDC=50%Im、IDC=100%Im時的直流偏磁情況進行實驗和仿真計算。無直流時，穩定邊界對應步長 h0=6.0×10-4s，步長h取 5.0×10-4s時計算收斂；隨著直流電流的升高，鐵心勵磁飽和程度加深，動態電感變化加劇，狀態矩陣特征值增大，穩定性降低。當 IDC=50%Im時，h0=2.57×10-4s，h若取5.0×10-4s不滿足穩定性要求；選取步長 h=2.50×10-4s，計算收斂。

分別采用基于時域場路耦合模型的定步長算法和自適應變步長算法進行計算，結果如圖4所示。

圖4 電流仿真計算和實驗測量結果Fig.4 Results of currents of simulation and experiment

空載運行時一次電流 i1即為勵磁電流 ie，由于磁場模型未考慮磁滯效應，根據 u→dψ/dt→dΦ/dt→B-H→i1的電磁耦合關系，i1為對稱波，如圖4a、圖4b所示，兩種方法計算的電流誤差較小。當存在直流時，變壓器勵磁受直流水平影響，ie波形在正負半周不再對稱，隨著直流電流增大，變壓器勵磁飽和程度加深，ie波形畸變嚴重，正負半周的不對稱程度加劇。由圖4c可知，空載運行直流偏磁時電流的計算結果與實驗測量基本相同；兩者存在的誤差可能由磁滯所導致，在直流偏磁情況下兩者誤差更小，從而驗證了本文所提方法是正確可行的。

設計低通濾波和時域差分模塊測量變壓器瞬時電感波形，并與計算的動態電感對比，如圖5所示。

圖5 動態電感Fig.5 Results of dynamic inductances of experiment and calculation

由圖5可以確定動態電感變化與勵磁非線性的對應關系。無直流時，L1在勵磁的正負半周為對稱波形，L1的波峰、波谷分別對應于勵磁的非飽和狀態與飽和狀態。ie接近零值時，鐵心勵磁處于不飽和區，L1數值趨于最大；當ie趨于各半周內的極值時，鐵心勵磁飽和程度逐漸加深，L1數值趨于最小。當存在直流時，L1受直流水平影響，波形在正負半周不對稱，隨著直流電流增大，變壓器勵磁飽和程度加深，i1畸變嚴重，L1波形在正負半周的不對稱程度加劇。為了便于對比，計算時采用單調雙曲函數對磁化曲線進行擬合，圖 5c的動態電感波形為固定步長計算結果，兩種算法計算的電感在整體波形上并沒有明顯差異，因此這里不再給出。

變壓器八分之一鐵心模型磁場計算結果如圖6所示。IDC=0時，鐵心最大磁感應強度Bmax為1.24T，IDC=25%Im時，Bmax升至1.25T。通過計算可知，隨著直流水平升高，變壓器勵磁飽和程度加深，鐵心磁通密度增加。

圖6 鐵心磁場計算結果Fig.6 Magnetic distribution in iron core

采用定步長算法計算變壓器正常運行，研究h、w對電感變化率和計算速度的影響。分析不同方案的計算穩定性、精確性及效率，結果見表2。

表2 定步長不同方案的計算結果Tab.2 Results of cases with fixed-step algorithm

變壓器時域場路耦合計算過程中的電感參數具有時變性，導致勵磁飽和時對應的穩定域邊界較小。直流偏磁時隨著直流分量增大，變壓器勵磁飽和程度加劇，對計算穩定性的要求也更為嚴格。當迭代計算的時間步長滿足穩定域條件時，可以保證計算收斂，但是按該步長計算將占用較多資源和時間，效率較低。總結其規律，計算的穩定性由時間步長確定，滿足絕對穩定域時，步長對計算效率和精確性的影響很小，此時效率和精確性主要取決于磁場的求解頻率。不難看出，定步長算法在計算變壓器非線性磁場時存在效率和精確性的優化問題。

4.2 自適應變步長算法

采用變步長算法進行計算，并與定步長算法的效率和精確性對比，結果見表3。

表3 不同方法的計算結果Tab.3 Results of different algorithms

結果表明，兩種方法都能模擬變壓器直流偏磁時的勵磁過程。當時間步長選定為0.25ms時，雖然能夠在一定程度上保證計算的穩定性和精確性，但是整個過程將會耗費大量時間。定步長算法雖然能夠模擬變壓器的勵磁過程，但不能充分考慮其隨時間變化的勵磁特性，因此在計算時可能存在不收斂的問題；同時由于該算法的h和w固定，導致其計算效率和精確性無法改善。自適應變步長方法能夠模擬變壓器直流偏磁時勵磁的時變非線性，時域計算每次迭代的h都能自動滿足絕對穩定域邊界，因此具有良好的計算穩定性；在此基礎上合理選取 w能達到優化計算效率和精確性的目的。

對于實際工程中普遍關注的大型電力變壓器直流偏磁問題，定步長算法可能受直流水平影響而不再適用，因而更能凸顯出自適應優化算法的計算效率和精確性。目前作者還未能獲取足夠的實際數據以進行大型變壓器的直流偏磁計算，因此該內容將是下一步研究的重點。

5 結論

結合變壓器直流偏磁問題，基于時域場路耦合計算的穩定性分析提出一種自適應變步長算法，通過研究得出以下結論：

（1）變壓器非線性勵磁具有時變性，直流偏磁時勵磁飽和程度加深，狀態矩陣特征值增大，對應的穩定域邊界變小。時域場路耦合計算時必須首先考慮穩定性問題，其計算步長應滿足收斂條件，在此基礎上合理選取磁場求解頻率能夠有效提高計算精確性和效率。

（2）與定步長時域場路耦合算法相比，自適應變步長算法計算變壓器直流偏磁時能夠有效模擬內部勵磁的時變非線性，具有良好的計算穩定性，并能實現計算效率與精確性的合理優化。空載直流偏磁的電流計算結果與實驗基本相同，驗證了該方法的正確性。

[1] C M Liu, L G Liu, R Pirjola. Geomagnetically induced currents in the high-voltage power grid in China[J]. IEEE Transactions on Power Delivery, 2010,24(4): 2368-2374.

[2] 郭滿生, 梅佳華, 張喜樂, 等. 直流偏磁條件下單相三柱電力變壓器的損耗計算[J]. 電工技術學報,2010, 25(7): 67-71.Guo Mansheng, Mei Jiahua, Zhang Xile, et al.Calculation of losses in three phase limb transformer under DC-biasing[J]. Transactions on China Electrotechnical Society, 2010, 25(7): 67-71.

[3] 趙志剛, 劉福貴, 張俊杰, 等. 直流偏磁條件下變壓器勵磁電流的實驗與分析[J]. 電工技術學報,2010, 25(4): 71-76.Zhao Zhigang, Liu Fugui, Zhang Junjie, et al.Measurement and analysis of magnetizing current in DC-biased transformers[J]. Transactions on China Electrotechnical Society, 2010, 25(4): 71-76.

[4] Lu Shu, Liu Yilu. FEM analysis of DC saturation to assess transformer susceptibility to geomagnetically induced currents[J]. IEEE Transactions on Power Delivery, 1993, 8(3): 1367-1376.

[5] Mauricio Valencia Ferreira da Luz, Jean Vianei Leite,et al. Three-phase transformer modeling using a vector hysteresis model and including the eddy current and the anomalous losses[J]. IEEE Transactions on Magnetics, 2010, 46(8): 3201-3204.

[6] Zhao X J, Lu J W, et al. Analysis of DC bias phenomenon by the harmonic Balance Finite-Element method[J]. IEEE Transactions on Power Delivery,2011, 26(1): 475-485.

[7] Zhao X J, Li L, et al. Characteristic analysis of the square laminated core under dc-biased magnetization by the fixed-point harmonic-balanced FEM[J]. IEEE Transactions on Magnetics, 2012, 48(2): 747-750.

[8] Yamada R, Biringer P P, et al. Calculation of nonlinear eddy-current problems by the harmonic balance finite element method[J]. Transactions on Magnetics, 1991, 26(5): 4122-4125.

[9] 趙小軍. 基于諧波平衡有限元的變壓器直流偏磁特性研究[D]. 保定: 華北電力大學, 2010.

[10] 李泓志, 崔翔, 盧鐵兵, 等. 變壓器直流偏磁的電路-磁路模型[J]. 中國電機工程學報, 2009, 29(27):119-125.Li Hongzhi, Cui Xiang, Lu Tiebing, et al. Electric circuit and magnetic circuit combined model of DC biased power transformer[J]. Proceedings of the CSEE, 2009, 29(27): 119-125.

[11] 郭健, 林鶴云, 徐子宏, 等. 單螺旋繞組變壓器支路電流的場路耦合計算及分析[J]. 電工技術學報,2010, 25(4): 65-70.Guo Jian, Lin Heyun, Xu Zihong, et al. Calculation and analysis of branch currents of single spiral winding transformer based on field circuit coupled method[J].Transactions on China Electrotechnical Society, 2010,25(4): 65-70.

[12] 潘超. 基于時域場路耦合模型的變壓器直流偏磁電磁特性研究[D]. 保定: 華北電力大學, 2013.

[13] Shi Z W, Rajanathan C B. A method of approach to transient eddy current problems coupled with voltage sources[J]. IEEE Transactions on Magnetics, 1996,32(3): 1082-1085.

[14] D Rodger, N Allen. Calculation of transient 3D eddy currents in nonlinear media verification using a rotational test rig[J]. IEEE Transactions on Magnetics,1994, 30(5): 2988-2991.

[15] S L Ho, W N Fu. An incremental method for studying the steady state performance of induction motors using time stepping finite element model[J]. IEEE Transactions on Magnetics, 1997, 33(2): 1374-1377.

[16] S L Ho, S X Niu, et al. A power-balanced timestepping finite element method for transient magnetic field computation[J], IEEE Transactions on Magnetics,2012, 48(2): 291-294.

[17] 潘超, 王澤忠, 楊敬瑀, 等. 變壓器直流偏磁瞬態場路耦合計算穩定性分析[J]. 電工技術學報, 2013,28(5): 82-87.Pan Chao, Wang Zezhong, Yang Jingyu, et al. Stability analysis based on transient magnetic-circuit coupled method for DC-biased transformer[J]. Transactions on China Electrotechnical Society, 2013, 28(5): 82-87.

[18] 王澤忠, 潘超, 周盛, 等. 基于棱邊有限元的變壓器場路耦合瞬態模型[J]. 電工技術學報, 2012,27(9): 146-152.Wang Zezhong, Pan Chao, Zhou Sheng, et al. Transient magnetic-circuit coupled model of transformer based on edge finite element method[J]. Transactions on China Electrotechnical Society, 2012, 27(9): 146-152.