一種變壓器時(shí)域場(chǎng)路耦合計(jì)算的自適應(yīng)算法

潘 超 何雯婷 郝文波 徐冰亮

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

1 引言

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

另一方面，傳統(tǒng)時(shí)域計(jì)算通常采用單步低階方法，其截?cái)嗾`差較大，且由于時(shí)間迭代造成誤差累積，易引起求解振蕩。若細(xì)化空間網(wǎng)格或縮短步長(zhǎng)，雖然可以在一定程度上減小振蕩，但計(jì)算規(guī)模會(huì)大幅度增加[13,14]。研究表明，變壓器瞬態(tài)磁場(chǎng)和時(shí)域場(chǎng)路計(jì)算中均存在計(jì)算穩(wěn)定性、效率及精確性的問(wèn)題[15,16]。文獻(xiàn)[17]討論了變壓器時(shí)域場(chǎng)路計(jì)算過(guò)程中非線性勵(lì)磁、步長(zhǎng)選擇及時(shí)間離散在直流擾動(dòng)下的變化情況。但是關(guān)于變壓器時(shí)域場(chǎng)路耦合計(jì)算的穩(wěn)定性、精度及效率等優(yōu)化問(wèn)題，目前罕有文獻(xiàn)研究。

本文結(jié)合變壓器直流偏磁問(wèn)題，推導(dǎo)計(jì)算收斂時(shí)的穩(wěn)定域邊界，確定時(shí)間步長(zhǎng)與計(jì)算穩(wěn)定性的關(guān)系，研究變壓器三維時(shí)域場(chǎng)路耦合計(jì)算效率與精確性的優(yōu)化方法。根據(jù)穩(wěn)定性分析提出一種自適應(yīng)變步長(zhǎng)算法，在計(jì)算效率和精確性方面與定步長(zhǎng)算法對(duì)比，并通過(guò)實(shí)驗(yàn)測(cè)量對(duì)其計(jì)算結(jié)果進(jìn)行驗(yàn)證。

2 時(shí)域場(chǎng)路耦合

利用時(shí)域場(chǎng)路耦合方法計(jì)算變壓器磁場(chǎng)和等效電路，動(dòng)態(tài)電感和時(shí)域電流為關(guān)鍵耦合參數(shù)，分別由三維棱邊元磁場(chǎng)模型和等效電路模型計(jì)算。

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

變壓器磁鏈方程

式中，ψ 為磁鏈向量；i為繞組電流向量；LS為靜態(tài)電感矩陣，表示激勵(lì)電流與磁鏈的穩(wěn)態(tài)關(guān)系。

實(shí)際上，變壓器勵(lì)磁過(guò)程中電感與電流均體現(xiàn)出時(shí)變非線性的特點(diǎn)[18]。由u=dψ /dt，推導(dǎo)變壓器等效電路的時(shí)域響應(yīng)方程

式中，u為交流電壓向量；UDC為直流電壓；LD為考慮勵(lì)磁變化的動(dòng)態(tài)電感矩陣，需要通過(guò)磁場(chǎng)模型計(jì)算獲得。

基于能量擾動(dòng)的思想，根據(jù)系統(tǒng)能量變化計(jì)算動(dòng)態(tài)電感參數(shù)。當(dāng)線圈電流增加δip（0≤δ≤1，單相變壓器中 1≤p≤2）時(shí)，將電源總能量與動(dòng)態(tài)電感和電流關(guān)聯(lián)，得到

若由電流增量δi引起的場(chǎng)量變化為δH、δB，變壓器內(nèi)部系統(tǒng)的磁場(chǎng)能量增量為

由能量平衡原理，方程（5）和（6）的能量相等，則可計(jì)算動(dòng)態(tài)電感LD。

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

式中，h為步長(zhǎng)；s1～s4為步長(zhǎng)內(nèi)的分段計(jì)算斜率。

利用繞組電流i與動(dòng)態(tài)電感LD即可實(shí)現(xiàn)時(shí)域場(chǎng)路耦合，其迭代計(jì)算原理如圖1所示。

式中，μ 為磁導(dǎo)率；J為電流密度，需要通過(guò)電路模型計(jì)算得到。

應(yīng)用格林定理，得伽遼金加權(quán)余量方程

式中，Mm為權(quán)函數(shù)序列，且權(quán)函數(shù)與基函數(shù)相同。

圖1 時(shí)域場(chǎng)路耦合Fig.1 Time-domain magnetic field and electrical circuit coupling

3 自適應(yīng)變步長(zhǎng)算法

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

另一方面，在tk時(shí)刻式（4）中的泛函可化為常微分方程，得到式中，v為泰勒級(jí)數(shù)項(xiàng)的組合函數(shù)。

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

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

時(shí)域場(chǎng)路耦合計(jì)算收斂所對(duì)應(yīng)的絕對(duì)穩(wěn)定域判據(jù)為│G│＜1。確定時(shí)域場(chǎng)路耦合計(jì)算的穩(wěn)定邊界為

式中，h0為穩(wěn)定邊界對(duì)應(yīng)的計(jì)算步長(zhǎng)。

對(duì)時(shí)域場(chǎng)路耦合迭代計(jì)算進(jìn)行改進(jìn)，提出一種自適應(yīng)變步長(zhǎng)算法，其基本原理為

（1）時(shí)域場(chǎng)路耦合計(jì)算受初值、步長(zhǎng)等因素的影響，在起始階段存在過(guò)渡過(guò)程，為了減少計(jì)算耗時(shí)，根據(jù)狀態(tài)矩陣特征值，在計(jì)算穩(wěn)定的前提下增大步長(zhǎng)，從而提高效率。

（2）經(jīng)過(guò)某段計(jì)算時(shí)間直至計(jì)算穩(wěn)定時(shí)，在每次磁場(chǎng)求解動(dòng)態(tài)電感之后，由變壓器狀態(tài)矩陣計(jì)算特征值，根據(jù)穩(wěn)定域邊界確定對(duì)應(yīng)的h0，引入修正因子α(0＜α＜1)，確定下一時(shí)刻的計(jì)算步長(zhǎng)，h=αh0，進(jìn)行迭代計(jì)算。

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

定義tk時(shí)刻的動(dòng)態(tài)電感變化函數(shù)g(LD)，表示在該時(shí)間段內(nèi)動(dòng)態(tài)電感參數(shù)的變化速率，用來(lái)評(píng)價(jià)磁場(chǎng)計(jì)算的精確性。

式中，max(LD) 為動(dòng)態(tài)電感的最大值；ε 為變化率限值。

4 算例分析

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

圖2 實(shí)驗(yàn)變壓器Fig.2 Transformer for experiment

表1 變壓器參數(shù)Tab.1 Parameters of transformer

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

4.1 時(shí)域場(chǎng)路耦合定步長(zhǎng)算法

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

分別采用基于時(shí)域場(chǎng)路耦合模型的定步長(zhǎng)算法和自適應(yīng)變步長(zhǎng)算法進(jìn)行計(jì)算，結(jié)果如圖4所示。

圖4 電流仿真計(jì)算和實(shí)驗(yàn)測(cè)量結(jié)果Fig.4 Results of currents of simulation and experiment

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

設(shè)計(jì)低通濾波和時(shí)域差分模塊測(cè)量變壓器瞬時(shí)電感波形，并與計(jì)算的動(dòng)態(tài)電感對(duì)比，如圖5所示。

圖5 動(dòng)態(tài)電感Fig.5 Results of dynamic inductances of experiment and calculation

由圖5可以確定動(dòng)態(tài)電感變化與勵(lì)磁非線性的對(duì)應(yīng)關(guān)系。無(wú)直流時(shí)，L1在勵(lì)磁的正負(fù)半周為對(duì)稱波形，L1的波峰、波谷分別對(duì)應(yīng)于勵(lì)磁的非飽和狀態(tài)與飽和狀態(tài)。ie接近零值時(shí)，鐵心勵(lì)磁處于不飽和區(qū)，L1數(shù)值趨于最大；當(dāng)ie趨于各半周內(nèi)的極值時(shí)，鐵心勵(lì)磁飽和程度逐漸加深，L1數(shù)值趨于最小。當(dāng)存在直流時(shí)，L1受直流水平影響，波形在正負(fù)半周不對(duì)稱，隨著直流電流增大，變壓器勵(lì)磁飽和程度加深，i1畸變嚴(yán)重，L1波形在正負(fù)半周的不對(duì)稱程度加劇。為了便于對(duì)比，計(jì)算時(shí)采用單調(diào)雙曲函數(shù)對(duì)磁化曲線進(jìn)行擬合，圖 5c的動(dòng)態(tài)電感波形為固定步長(zhǎng)計(jì)算結(jié)果，兩種算法計(jì)算的電感在整體波形上并沒(méi)有明顯差異，因此這里不再給出。

變壓器八分之一鐵心模型磁場(chǎng)計(jì)算結(jié)果如圖6所示。IDC=0時(shí)，鐵心最大磁感應(yīng)強(qiáng)度Bmax為1.24T，IDC=25%Im時(shí)，Bmax升至1.25T。通過(guò)計(jì)算可知，隨著直流水平升高，變壓器勵(lì)磁飽和程度加深，鐵心磁通密度增加。

圖6 鐵心磁場(chǎng)計(jì)算結(jié)果Fig.6 Magnetic distribution in iron core

采用定步長(zhǎng)算法計(jì)算變壓器正常運(yùn)行，研究h、w對(duì)電感變化率和計(jì)算速度的影響。分析不同方案的計(jì)算穩(wěn)定性、精確性及效率，結(jié)果見表2。

表2 定步長(zhǎng)不同方案的計(jì)算結(jié)果Tab.2 Results of cases with fixed-step algorithm

變壓器時(shí)域場(chǎng)路耦合計(jì)算過(guò)程中的電感參數(shù)具有時(shí)變性，導(dǎo)致勵(lì)磁飽和時(shí)對(duì)應(yīng)的穩(wěn)定域邊界較小。直流偏磁時(shí)隨著直流分量增大，變壓器勵(lì)磁飽和程度加劇，對(duì)計(jì)算穩(wěn)定性的要求也更為嚴(yán)格。當(dāng)?shù)?jì)算的時(shí)間步長(zhǎng)滿足穩(wěn)定域條件時(shí)，可以保證計(jì)算收斂，但是按該步長(zhǎng)計(jì)算將占用較多資源和時(shí)間，效率較低。總結(jié)其規(guī)律，計(jì)算的穩(wěn)定性由時(shí)間步長(zhǎng)確定，滿足絕對(duì)穩(wěn)定域時(shí)，步長(zhǎng)對(duì)計(jì)算效率和精確性的影響很小，此時(shí)效率和精確性主要取決于磁場(chǎng)的求解頻率。不難看出，定步長(zhǎng)算法在計(jì)算變壓器非線性磁場(chǎng)時(shí)存在效率和精確性的優(yōu)化問(wèn)題。

4.2 自適應(yīng)變步長(zhǎng)算法

采用變步長(zhǎng)算法進(jìn)行計(jì)算，并與定步長(zhǎng)算法的效率和精確性對(duì)比，結(jié)果見表3。

表3 不同方法的計(jì)算結(jié)果Tab.3 Results of different algorithms

結(jié)果表明，兩種方法都能模擬變壓器直流偏磁時(shí)的勵(lì)磁過(guò)程。當(dāng)時(shí)間步長(zhǎng)選定為0.25ms時(shí)，雖然能夠在一定程度上保證計(jì)算的穩(wěn)定性和精確性，但是整個(gè)過(guò)程將會(huì)耗費(fèi)大量時(shí)間。定步長(zhǎng)算法雖然能夠模擬變壓器的勵(lì)磁過(guò)程，但不能充分考慮其隨時(shí)間變化的勵(lì)磁特性，因此在計(jì)算時(shí)可能存在不收斂的問(wèn)題；同時(shí)由于該算法的h和w固定，導(dǎo)致其計(jì)算效率和精確性無(wú)法改善。自適應(yīng)變步長(zhǎng)方法能夠模擬變壓器直流偏磁時(shí)勵(lì)磁的時(shí)變非線性，時(shí)域計(jì)算每次迭代的h都能自動(dòng)滿足絕對(duì)穩(wěn)定域邊界，因此具有良好的計(jì)算穩(wěn)定性；在此基礎(chǔ)上合理選取 w能達(dá)到優(yōu)化計(jì)算效率和精確性的目的。

對(duì)于實(shí)際工程中普遍關(guān)注的大型電力變壓器直流偏磁問(wèn)題，定步長(zhǎng)算法可能受直流水平影響而不再適用，因而更能凸顯出自適應(yīng)優(yōu)化算法的計(jì)算效率和精確性。目前作者還未能獲取足夠的實(shí)際數(shù)據(jù)以進(jìn)行大型變壓器的直流偏磁計(jì)算，因此該內(nèi)容將是下一步研究的重點(diǎn)。

5 結(jié)論

結(jié)合變壓器直流偏磁問(wèn)題，基于時(shí)域場(chǎng)路耦合計(jì)算的穩(wěn)定性分析提出一種自適應(yīng)變步長(zhǎng)算法，通過(guò)研究得出以下結(jié)論：

（1）變壓器非線性勵(lì)磁具有時(shí)變性，直流偏磁時(shí)勵(lì)磁飽和程度加深，狀態(tài)矩陣特征值增大，對(duì)應(yīng)的穩(wěn)定域邊界變小。時(shí)域場(chǎng)路耦合計(jì)算時(shí)必須首先考慮穩(wěn)定性問(wèn)題，其計(jì)算步長(zhǎng)應(yīng)滿足收斂條件，在此基礎(chǔ)上合理選取磁場(chǎng)求解頻率能夠有效提高計(jì)算精確性和效率。

（2）與定步長(zhǎng)時(shí)域場(chǎng)路耦合算法相比，自適應(yīng)變步長(zhǎng)算法計(jì)算變壓器直流偏磁時(shí)能夠有效模擬內(nèi)部勵(lì)磁的時(shí)變非線性，具有良好的計(jì)算穩(wěn)定性，并能實(shí)現(xiàn)計(jì)算效率與精確性的合理優(yōu)化。空載直流偏磁的電流計(jì)算結(jié)果與實(shí)驗(yàn)基本相同，驗(yàn)證了該方法的正確性。

[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] 郭滿生, 梅佳華, 張喜樂(lè), 等. 直流偏磁條件下單相三柱電力變壓器的損耗計(jì)算[J]. 電工技術(shù)學(xué)報(bào),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] 趙志剛, 劉福貴, 張俊杰, 等. 直流偏磁條件下變壓器勵(lì)磁電流的實(shí)驗(yàn)與分析[J]. 電工技術(shù)學(xué)報(bào),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]. 保定: 華北電力大學(xué), 2010.

[10] 李泓志, 崔翔, 盧鐵兵, 等. 變壓器直流偏磁的電路-磁路模型[J]. 中國(guó)電機(jī)工程學(xué)報(bào), 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] 郭健, 林鶴云, 徐子宏, 等. 單螺旋繞組變壓器支路電流的場(chǎng)路耦合計(jì)算及分析[J]. 電工技術(shù)學(xué)報(bào),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] 潘超. 基于時(shí)域場(chǎng)路耦合模型的變壓器直流偏磁電磁特性研究[D]. 保定: 華北電力大學(xué), 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] 潘超, 王澤忠, 楊敬瑀, 等. 變壓器直流偏磁瞬態(tài)場(chǎng)路耦合計(jì)算穩(wěn)定性分析[J]. 電工技術(shù)學(xué)報(bào), 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] 王澤忠, 潘超, 周盛, 等. 基于棱邊有限元的變壓器場(chǎng)路耦合瞬態(tài)模型[J]. 電工技術(shù)學(xué)報(bào), 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.