劉 剛 胡萬君 郝世緣 劉云鵬 李 琳

油浸式變壓器繞組瞬態(tài)溫升降階快速計算方法

劉 剛1胡萬君1郝世緣1劉云鵬1李 琳2

（1. 華北電力大學河北省輸變電設(shè)備安全防御重點實驗室 保定 071003 2. 新能源電力系統(tǒng)全國重點實驗室（華北電力大學） 北京 102206）

油浸式電力變壓器繞組溫升監(jiān)測是保證其安全穩(wěn)定運行的重要手段。為了改善采用有限元方法計算油浸式電力變壓器繞組瞬態(tài)溫升時存在效率不高的問題，提出一種結(jié)構(gòu)保留的本征正交分解（SPOD）與離散經(jīng)驗插值方法（DEIM）相結(jié)合的計算策略。首先，該文采用最小二乘有限元法（LSFEM）與迎風有限元法（UFEM）構(gòu)建變壓器繞組瞬態(tài)溫升計算控制方程；其次，針對控制方程的特點，引入SPOD方法，通過將采樣時間內(nèi)的計算結(jié)果構(gòu)成快照矩陣，建立降階模型，降低有限元剛度矩陣的計算階數(shù)，提高求解有限元方程的效率；然后，為了改善本征正交分解方法對于非線性問題效率提升不高的缺陷，結(jié)合DEIM算法，對有限元方程中的非線性項進行插值處理，從而減少每一時步形成總體剛度矩陣的時間，進一步提高總體計算效率。為了驗證文章所提算法的精確性及高效性，根據(jù)油浸式電力變壓器繞組的基本特點，建立了單分區(qū)分匝繞組傳熱模型，對其瞬態(tài)傳熱過程進行計算，結(jié)果表明：基于SPOD-DEIM的有限元降階計算能夠在保證精度的前提下有效提高計算效率，與全階計算結(jié)果相比，流場與溫度場的計算誤差均不超過1.5%，且計算效率提升5.1倍。同時，為了充分說明SPOD-DEIM算法在工程應(yīng)用中的價值，該文基于110 kV變壓器繞組搭建了溫升實驗平臺，建立了八分區(qū)分匝繞組數(shù)值計算模型，對算法的精度、效率及工程應(yīng)用價值進行了驗證及討論，計算及實驗結(jié)果表明：精度方面，降階計算較全階計算的瞬態(tài)全過程計算誤差小于2.5%，且與實驗結(jié)果相比，誤差不超過5.41 K；效率方面，降階計算的全過程計算時間為54.28 h，與全階計算相比，計算效率提升至10.57倍，與商業(yè)仿真軟件Fluent相比，效率提升至6.37倍，充分說明所提算法的高效性及工程應(yīng)用價值，為大型電力設(shè)備快速仿真提供新思路。

繞組瞬態(tài)溫升 結(jié)構(gòu)保留 本征正交分解 離散經(jīng)驗插值 降階計算 溫升實驗

0 引言

電力變壓器作為電力系統(tǒng)當中的主要電氣設(shè)備，其穩(wěn)定性和安全性需要得到切實的保證。在實際運行過程中，變壓器繞組產(chǎn)生的歐姆損耗與渦流損耗使得繞組溫度上升，一旦繞組溫度過高，就會面臨損壞甚至燒毀的風險，因此對于變壓器繞組溫升的監(jiān)測至關(guān)重要。針對油浸式電力變壓器，繞組溫升的監(jiān)測手段目前主要包括三種：直接測量法[1-3]、代理模型法[4-6]及數(shù)值計算法[7-9]。其中，數(shù)值計算方法是較為主要的方法。

數(shù)值計算方法通過建立流熱耦合場的控制方程，對溫升問題進行數(shù)值計算，從而得到繞組整體的溫度分布或者變化曲線。針對變壓器繞組溫升問題的計算目前采用較多的數(shù)值計算方法有有限元法、有限體積法等。相較于有限體積法，有限元法具有靈活性高、直觀性強、求解精度高等諸多優(yōu)勢，應(yīng)用前景更為廣泛。同時，隨著有限元方法的不斷發(fā)展，迎風有限元（Upwind Finite Element Method, UFEM）[10]、最小二乘有限元（Least Squares Finite Element Method, LSFEM）等[11-12]方法克服了傳統(tǒng)有限元算法在數(shù)值上的不穩(wěn)定性，極大地提高了數(shù)值模型的計算精度，使其更加逼近實際結(jié)果。然而，計算時間成本過高是阻礙該類方法應(yīng)用于工程實際的關(guān)鍵問題。

有限元法的瓶頸在于計算效率與模型剖分的精細程度相關(guān)。在實際工程問題當中，為了精確計算，需要對模型剖分的網(wǎng)格數(shù)量非常多，從而導(dǎo)致計算中的有限元剛度矩陣規(guī)模特別大，這對計算機內(nèi)存和計算效率提出了考驗。為了提高計算繞組溫升計算的速度，謝裕清[13-14]、榮世昌等[15]提出采用本征正交分解（Proper Orthogonal Decomposition,POD）的方法得到降階正交基，降低有限元方程階數(shù)，進而提高計算效率，該方法效果顯著且計算成本較低。然而，當有限元方程中涉及的物理量較多且彼此之間量級相差較大時，傳統(tǒng)的POD方法可能需要更多的正交基才可以保證精度，因此，Shuai Yan等[16-17]提出采用結(jié)構(gòu)保留的本征正交分解方法（Structure-Preserved Proper Orthogonal Decomposition, SPOD），這種方法的好處在于簡化問題的同時能夠保證完整問題的塊結(jié)構(gòu)，更適用于物理量較多的有限元計算。

另一方面，SPOD方法雖然能夠降低方程的階數(shù)，使得有限元方程能夠在一個階數(shù)較低的子空間內(nèi)求解，但是，由于非線性項的存在，其在形成剛度矩陣的過程中仍然需要遍歷所有單元，這是另一個限制計算效率的瓶頸問題。為了克服這一問題，本文結(jié)合采用了S. Chaturantabut和D. S. Sorensen[18]提出的離散經(jīng)驗插值方法（Discrete Empirical Interpolation Method, DEIM），該方法不需要將非線性項中的所有元素進行回代計算，其僅需通過貪婪算法選取方程非線性項中的部分元素，進而用少量元素對全體進行插值，從而達到提高計算效率的目的。

綜上所述，本文首先結(jié)合LSFEM以及UFEM推導(dǎo)油浸式電力變壓器流固耦合傳熱過程的瞬態(tài)控制方程；然后詳細介紹SPOD與DEIM的算法原理，并將其引入所推導(dǎo)的瞬態(tài)控制方程中，得到降階計算模型；最后根據(jù)110 kV電力變壓器結(jié)構(gòu)特點，搭建了油浸式電力變壓器繞組溫升實驗平臺，并建立二維數(shù)值傳熱模型, 結(jié)合相關(guān)實驗，提出并討論了SPOD-DEIM算法的相關(guān)性能及實際應(yīng)用價值。

1 瞬態(tài)流熱耦合控制方程

油浸式電力變壓器繞組的溫升過程顯然是一個流熱耦合計算問題。對于其中反映油流變化的流場，基于流體力學中的質(zhì)量守恒與動量守恒定律，結(jié)合渦量，同時對時間項采用Crank-Nicolson方法進行離散，其二維瞬態(tài)控制方程可表示為[19]

其中

采用LSFEM方法對式（1）進行空間離散，得到其有限元求解格式為

其中

式中，e為有限單元總數(shù)；e為求解場域；為有限元形函數(shù)；為各節(jié)點的自由度矩陣。

對于溫度場，基于能量守恒定律，采用UFEM方法得到其瞬態(tài)溫度場控制方程為[19]

其中

同樣采用Crank-Nicolson方法對其中的時間項進行離散，可以得到式（3）的最終求解格式為

其中

將式（2）與式（4）聯(lián)立，結(jié)合相應(yīng)的邊界條件與初始條件，即可以得到油浸式電力變壓器繞組瞬態(tài)溫升計算方程。根據(jù)LSFEM-UFEM結(jié)合的有限元方法即可較為精確地獲得繞組瞬態(tài)溫升變化過程中的流場及溫度場分布。然而，在實際工程應(yīng)用中，模型的網(wǎng)格規(guī)模和節(jié)點數(shù)量一般較大，采用有限元方法進行求解面臨時間成本過高的問題，因此，需要考慮引入降階技術(shù)，降低計算規(guī)模，從而提高計算效率。本文采用結(jié)構(gòu)保留的本征正交分解方法與離散經(jīng)驗插值方法改善該問題。

2 基于SPOD-DEIM的瞬態(tài)流熱計算

2.1 本征正交分解

基于變分方法與拉格朗日乘子法[20]，可以將上述最優(yōu)問題轉(zhuǎn)化為式（6）的矩陣特征值計算問題。

式中，為特征值矩陣，=[12×××λ]，且1＞2＞×××＞λ。該算法的截斷誤差可以表示為[21]

可以看出，在每個迭代步內(nèi)，系統(tǒng)方程由全階計算時的階降低為了階，矩陣計算規(guī)模得以減少，能夠很好地提高求解效率。然而，對于非線性問題而言，剛度矩陣的形成與前一時間步的物理量有關(guān)，因此每時刻計算都需要遍歷所有單元重新形成，這是降低計算效率的另一瓶頸問題，結(jié)合離散經(jīng)驗插值算法可以對其做出改善。

2.2 離散經(jīng)驗插值法

離散經(jīng)驗插值源于經(jīng)驗插值，其主要針對非線性問題提出了一種離散點插值思路，從而減小非線性計算引起的效率降低。

對于非線性項()∈R×1，通過奇異值分解選取一組正交基，使其可以表示為

上述方程對于求解，其未知數(shù)有個，方程有組，屬于超定方程組，因此需要剔除中的行。采用的方法是引入布爾矩陣，有

式中，1為第1個元素為1的單位向量。通過引入布爾矩陣，可以使得矩陣可逆，從而求得為

再將的求解結(jié)果代入式（10），則有

式中，為DEIM插值矩陣，在求解過程中僅需計算一次。同時，由于布爾矩陣的選擇作用，實際運算中僅僅需要中的個離散點，而其他個點的值則由插值得到。布爾矩陣的形成與個插值點的選取通過基于貪婪算法的搜索完成，其搜索思路為：將第一個插值點的位置選擇為中的第一個基向量1絕對值最大元素的位置，然后依次用當前的基向量組表示下一基向量，選擇殘差向量的最大元素位置作為下一個插值點的位置。程序偽代碼見表1。

特別地，在有限元方法中，由于插值點表示的是有限元剖分區(qū)域內(nèi)的節(jié)點，其值的大小與該節(jié)點所屬的所有單元均有關(guān)，因此在形成總體剛度矩陣時需要注意，將選擇的插值點屬于的所有單元剛度矩陣均進行計算并組裝，如圖1所示。

綜上所述，僅需提前計算一次POD正交基及DEIM插值矩陣和插值點，即可從有限元方程降階和非線性項插值兩個方面提高原系統(tǒng)的計算效率。

然而，對于本文所探討的流場問題而言，其節(jié)點自由度包含四個物理量：橫向及縱向速度、壓力、渦量，且四個物理量之間量綱不同，數(shù)量級不同，四者均為構(gòu)成流場的重要特征。因此，在降階模型中，需要保持這四個物理量在有限元矩陣中原本的塊結(jié)構(gòu)，本文采用了SPOD-DEIM方案。

表1 DEIM插值點搜索算法

Tab.1 DEIM interpolation point search algorithm

注：max()表示max函數(shù)，其輸出結(jié)果依次表示為最大值以及最大值所在位置序號。

圖1 DEIM算法應(yīng)用于有限元

2.3 基于SPOD-DEIM方法的繞組溫升降階計算模型

油浸式電力變壓器繞組溫升計算涉及流場與溫度場的耦合，因此，引入降階算法需要同時考慮溫升計算中的流場控制方程與溫度場控制方程。首先，對于反映變壓器油流分布的瞬態(tài)流場計算方程，根據(jù)式（2），其離散形式可以用矩陣表示為

式中，、為有限元剛度矩陣，矩陣中各元素根據(jù)油流場中各節(jié)點的物性參數(shù)、位置以及時間步長形成；為有限元右端項矩陣，與各節(jié)點的油流速度有關(guān)；+1與分別為第+1時間步和第時間步的變壓器油狀態(tài)參數(shù)，包括橫向速度、縱向速度、壓強以及渦量。通過設(shè)置入口油流速度等邊界條件、場域內(nèi)油流初始速度等初始條件即可通過式（14）計算得到瞬態(tài)過程中油流的速度變化。

=α=α=α=α

這樣通過塊分解的方式即可合理選擇各部分降階正交基的階數(shù)，充分提取并保留各物理量的數(shù)據(jù)特征，同時也使得構(gòu)建的降階模型更具有物理意義。

獲得降階正交基后，為了降低計算階數(shù)，將原系統(tǒng)投影至降階子空間中，采用伽遼金投影方法，將系統(tǒng)方程重構(gòu)為

其中

在此基礎(chǔ)上，基于DEIM算法思路，分別選取非線性項u、v、p、ω、u、v、p、ω前個時間步的數(shù)值計算結(jié)果形成快照矩陣，通過奇異值分解得到正交基，再利用搜索算法得到個插值點對應(yīng)的布爾矩陣以及插值矩陣。綜上所述，采用SPOD-DEIM算法的流場計算方程為

其中

式中，ROM、ROM、ROM、ωROM、ROM、ROM、ROM、ωROM分別為四個物理量在SPOD-DEIM降階模型中的塊結(jié)構(gòu)。各部分的數(shù)據(jù)特征及物理特征在降階計算方程中得到了充分保留。

同理，繞組溫度場計算方程參考式（4）可以表示為

式中，()、()為與材料物性參數(shù)（如導(dǎo)熱系數(shù)等）相關(guān)的剛度陣；(,)、(,)為與流體速度相關(guān)的剛度陣；為右端項。同樣地，采用SPOD-DEIM方法對其進行降階，將控制方程重構(gòu)為

其中

式中，T為應(yīng)用SPOD方法獲得的正交基矩陣，由于溫度場的解向量中不存在類似于流場的多種狀態(tài)量，因此此處的T與采用傳統(tǒng)POD方法得到的正交基相同；λ與U為采用DEIM方法獲得的插值矩陣。

綜上所述，結(jié)合SPOD和DEIM方法的油浸式電力變壓器繞組瞬態(tài)溫升降階計算流程如圖2所示。

圖2 SPOD-DEIM算法流程

顯然，隨著模型的節(jié)點與網(wǎng)格規(guī)模逐漸增大，降階模型相較于全階模型的效率提升也將逐漸增加。

3 SPOD-DEIM算法油浸式電力變壓器繞組溫升計算中的應(yīng)用

為了討論SPOD-DEIM算法應(yīng)用于油浸式電力變壓器溫升計算時的適用性及高效性，參考油浸式電力變壓器結(jié)構(gòu)特點，建立了單分區(qū)分匝繞組傳熱模型對其精度和效率進行了討論。同時，為了進一步驗證該算法的實用價值，文章基于110 kV變壓器繞組基本結(jié)構(gòu)，搭建了產(chǎn)品級油浸式電力變壓器繞組溫升實驗平臺，并根據(jù)其結(jié)構(gòu)參數(shù)建立了二維繞組傳熱模型，在此基礎(chǔ)上進行計算分析。

3.1 SPOD-DEIM算法有效性驗證

3.1.1 二維單分區(qū)分匝繞組模型建立

油浸式電力變壓器的基本結(jié)構(gòu)如圖3所示，繞組內(nèi)部的擋板將繞組分成了幾個結(jié)構(gòu)大致相同的分區(qū)，各分區(qū)內(nèi)油流均從入口進入，以S形線路與繞組充分接觸后從出口流出，因此在初步驗證算法有效性時，可以選取其中一個單分區(qū)繞組模型進行討論。

圖3 油浸式電力變壓器繞組一般結(jié)構(gòu)

二維單分區(qū)分匝繞組模型如圖4所示，該分區(qū)包含九個線餅，十條油道，將線餅從上至下標號為線餅1～9，將油道從上至下標號為油道1～10。模型的總體尺寸為155.5 mm × 154.2 mm，油流入口寬度為10 mm，出口寬度為8 mm，除1號油道為3 mm外，其余油道均為6 mm。

圖4 二維單分區(qū)分匝繞組模型

設(shè)置出口為壓力邊界條件，壓力為0 Pa；將絕緣筒壁與擋板設(shè)置為無滑移壁面邊界條件，即速度為0 m/s；該模型的熱源為繞組的歐姆損耗，考慮溫度對繞組電阻的影響，給出損耗計算公式[12]為

式中，S0為繞組在0溫度下的損耗密度，在本算例中參考了文獻[12]取值為227 000 W/m3，并在Fluent軟件中通過udf設(shè)置與自編程程序相同的熱源；為導(dǎo)體溫度系數(shù)，對于銅導(dǎo)體一般可以取為0.003 93 K-1[22]。相關(guān)材料的物性參數(shù)見表2。

表2 模型物性參數(shù)

Tab.2 Model physical parameters

3.1.2 全階算法正確性驗證

由于本文所采用的降階算法均以式（2）與式（4）為基礎(chǔ)，因此在進行降階算法的驗證之前，有必要對基于LSFEM-UFEM的全階算法正確性進行說明。設(shè)置單分區(qū)繞組入口油流速度為0.05 m/s，初始溫度為290 K，分別采用仿真軟件Fluent和基于LSFEM-UFEM的自編程程序?qū)υ摾@組模型的瞬態(tài)溫升過程進行計算，模型網(wǎng)格數(shù)量為46 827，節(jié)點數(shù)量為187 783。本文所采用的算法編程與計算軟件為：Matlab R2021a，計算機配置為：CPU Intel Core i9-12900KF，內(nèi)存128 GB。當繞組的傳熱過程達到穩(wěn)態(tài)時，流場與溫度場的云圖對比如圖5所示。

從圖5中可以看出，本文所采用的全階計算方法與Fluent的計算結(jié)果大致相同，各線餅及油道的最高溫升和最大流速如圖6所示。

圖5 穩(wěn)態(tài)云圖對比

圖6 各線餅（油道)最高溫升（最大油速）

可以看出，二者的計算結(jié)果幾乎一致，流場最大誤差為0.002 3 m/s；溫度場的最大誤差為0.443 0 K。選取油速最快的第2油道以及溫度最高的第7線餅的中心位置為監(jiān)測點，將油流流速與繞組溫度的瞬態(tài)變化過程繪圖如圖7所示。

圖7 參考點油流流速及溫度變化

在監(jiān)測點的瞬態(tài)變化過程中，第2油道的流速最大誤差出現(xiàn)在第0.04 s，為0.000 7 m/s，相對誤差為3.6%；第7線餅的溫度最大誤差出現(xiàn)在第870 s，為0.470 1 K，攝氏溫標相對誤差為1.75%。由于仿真軟件Fluent采用的數(shù)值計算方法為有限體積法，而本文的自編程方法為LSFEM-UFEM結(jié)合的有限元法，因此二者的計算過程中必然會存在一定誤差，但根據(jù)前文的分析可以看出，這個誤差處在可以接受的范圍內(nèi)，因此，可驗證本文所采用的基于LSFEM-UFEM全階算法的正確性，在該算法的基礎(chǔ)上進行降階計算具有一定可信度。

3.1.3 降階算法正確性驗證

基于3.1.2節(jié)中采用的全階計算模型，以繞組溫升從0～100 s的計算數(shù)據(jù)為樣本數(shù)據(jù)，各物理量根據(jù)式（8）選取模態(tài)能量大于99.9%的SPOD正交基數(shù)量，并選取不會使得計算發(fā)散的最小DEIM插值點數(shù)量，基于式（16）與式（18）構(gòu)造流場與溫度場的降階模型并進行降階計算，將計算結(jié)果與全階計算結(jié)果進行對比。

計算達到穩(wěn)態(tài)時，全階計算與SPOD-DEIM降階計算的流場和溫度場對比如圖8所示。

圖8 穩(wěn)態(tài)結(jié)果對比

圖8a與圖8b表示分區(qū)各線餅和油道中心線上的溫度與速度大小。從圖8可以看出，全階計算模型與降階計算模型達到穩(wěn)態(tài)時的計算結(jié)果幾乎完全相同，其中，速度場的最大誤差出現(xiàn)在第2油道，為0.000 01 m/s；溫度場的最大誤差出現(xiàn)在第4線餅，為0.067 K。

定義第時間步的降階誤差為

表3 計算效率對比

Tab.3 Calculation efficiency comparison

式中，為求解變量，∈(，)；下標Full表示全階計算結(jié)果，ROM表示降階計算結(jié)果；node為流場域或溫度場域所有節(jié)點數(shù)量。計算瞬態(tài)計算過程中各個時間步下的降階誤差如圖9所示。

圖9　流場及溫度場降階計算誤差

可以看出，在該過程中，流場降階計算誤差最大為1.26%，溫度場降階計算誤差最大為0.56%。兩場的降階計算精度均處于可接受范圍內(nèi)，由此可初步說明降階算法的正確性。同時，為了反映SPOD-DEIM算法的計算效率，將全階計算、SPOD降階計算、SPOD-DEIM降階計算以及Fluent單線程計算的平均單步計算時間記錄于表3中。表中“單步總時間”列計算時間由流場和溫度場兩部分計算時間構(gòu)成，如全階模型“單步總時間”中流場計算時間為33.09 s，溫度場計算時間為22.34 s，總計算時間為55.43 s。

從表3中可以看出，通過引入SPOD算法，能夠有效降低有限元方程的階數(shù)，極大地提高有限元方程求解時間，但從單步計算時間來看，流場計算時間從33.09 s降低至12.07 s，后者為前者的36.47%；溫度場計算時間從22.34 s降低到19.94 s，降階后計算時間為前者的89.26%。效率提升并不明顯，這是由于流固耦合問題具有的強非線性，使得每一時步的有限元剛度矩陣需要重新形成，耗費較多時間。

在結(jié)合了DEIM插值方法之后，SPOD-DEIM降階計算較SPOD降階計算在形成有限元剛度矩陣方面效率得到了進一步提高。以流場為例，剛度矩陣的形成時間從7 s左右減少到0.056 s；同時，對于單步計算時間，較全階計算模型，流場計算時間從33.09 s降低至5.89 s，減少為原先的17.78%，溫度場計算時間從22.34 s降低至4.99 s，減少為原先的22.33%，SPOD-DEIM降階總體計算時間為全階計算的19.63%。

大幅度減少了計算時間。且與仿真軟件Fluent的串行計算速度相比，降階模型的效率也提升近2倍。這說明了本文所采用降階算法的高效性。

可以看到，在單分區(qū)分匝模型中，降階的瞬態(tài)計算與全階計算效率提升并不是很高，這是由于該模型的規(guī)模較小，降階計算的效率提升并不明顯。從2.3節(jié)中可以得到，降階模型效率的提升顯然與降階之后節(jié)點與網(wǎng)格數(shù)量減少的倍數(shù)有關(guān)，因此，隨著全階模型的網(wǎng)格與節(jié)點規(guī)模逐漸擴大，基于SPOD-DEIM算法的降階計算相較于全階計算將有更為明顯的效率優(yōu)勢，這一點將在3.2節(jié)的溫升實驗中加以驗證。

3.2 溫升實驗驗證

為了充分討論SPOD-DEIM算法應(yīng)用于實際油浸式電力變壓器溫升計算時的適用性及高效性，參考110 kV電力變壓器結(jié)構(gòu)特點，本文搭建了產(chǎn)品級油浸式電力變壓器繞組溫升實驗平臺，并根據(jù)其結(jié)構(gòu)參數(shù)建立了二維繞組傳熱模型，在此基礎(chǔ)上進行計算分析。

該實驗平臺由空心餅式無感線圈、器身絕緣、隔熱油箱、片式散熱器、油泵、風機、導(dǎo)油管路和熱電偶等部件組成。圖10與圖11為實驗?zāi)Ｐ驼w布置外觀圖與示意圖。

圖10　溫升實驗平臺

繞組整體分為8個分區(qū)，從頂端到底端依次編號為分區(qū)1～8，其中，分區(qū)1～3各有7個線餅，分區(qū)4～8各有9個線餅，為了便于比較與說明，將構(gòu)成繞組的線餅依次編號為1～66。每個線餅由15匝導(dǎo)線組成，每匝導(dǎo)線由2根扁銅導(dǎo)線并繞而成，由內(nèi)軸向到外軸向?qū)⒚總€線餅中的導(dǎo)線依次編號為線匝1～30。

圖11　實驗平臺示意圖

根據(jù)實驗采用的餅式繞組軸對稱的結(jié)構(gòu)特點，構(gòu)造二維等效數(shù)值模型如圖12所示。為精確數(shù)值計算的結(jié)果，該模型考慮了餅式繞組的分匝結(jié)構(gòu)以及絕緣紙、擋板等對繞組溫升有顯著影響的結(jié)構(gòu)件。變壓器油、擋板以及絕緣紙的物性參數(shù)同表2。

圖12 二維八分區(qū)分匝繞組模型

實驗條件設(shè)置如下：

1）本實驗的測溫裝置為熱電偶，分別布置于12、20、30與38號線餅，同時，每個測溫線餅上設(shè)置11個測點，分別為1、4、7、10、13、16、19、22、25、28、30號線匝。

2）本實驗的目的在于測量油浸式電力變壓器繞組在強油散熱條件下的溫升變化趨勢，記錄進口油量為14.4m3/s，進口油溫隨時間變化，根據(jù)入口處的熱電偶記錄每一時刻的油溫用以仿真計算，初始油溫為16.24℃，功率分析儀得到繞組的功率為25.000 4 kW ，其功率因數(shù)測量為0.999 82，可以假定作為熱源的繞組損耗中只有歐姆損耗，計算過程中熱源密度同式（19）。電源和進口流量保持恒定。

3）溫度記錄儀每10 s采集一次各測點的繞組溫度。當頂層油溫變化率為1℃/h并維持2 h，認為油溫達到穩(wěn)定狀態(tài)。實驗共進行了7.66 h，收集各測點數(shù)據(jù)2 758組。

進行溫升實驗的同時，采用圖12所示的八分區(qū)二維軸對稱模型進行仿真分析，模型網(wǎng)格數(shù)量為211 613，節(jié)點數(shù)量為846 211，設(shè)置時間間隔Δ=10 s，并以0～1 000 s的全階計算結(jié)果作為樣本數(shù)據(jù)形成降階子空間，基于SPOD-DEIM方法對繞組的溫升全過程進行降階計算。將Fluent軟件的仿真結(jié)果、基于SPOD-DEIM算法的自編程結(jié)果、實驗結(jié)果進行對比，做如下討論。

達到穩(wěn)態(tài)時，12號、20號、30號、37號線餅上繞組溫度的分布如圖13所示。

圖13　穩(wěn)態(tài)時各線餅監(jiān)測點溫度分布

從圖13可以看出，三種方法達到穩(wěn)態(tài)時的溫度差異不大，本文所采用的SPOD-DEIM降階算法與實驗結(jié)果的最大誤差出現(xiàn)在20號線餅28號測點上，為4.4 K，在工程上屬于可接受的范圍。同時，與Fluent算法的計算結(jié)果相比，本文自編程結(jié)果在這幾個測點上更接近于實驗結(jié)果，初步說明了本文所提出算法的正確性。

選擇12號繞組4號測點、20號繞組22號測點、30號繞組4號測點與38號繞組28號測點，導(dǎo)出其運行過程中的溫度變化數(shù)據(jù)，并繪制溫度變化曲線如圖14所示。

可以看出，在這幾個測點上，三種方式所得到的溫度變化趨勢近乎一致，且差異不大。SPOD-DEIM算法與實驗結(jié)果的最大誤差出現(xiàn)于2 500 s左右，在12號線餅4號測點上，為5.41 K，F(xiàn)luent仿真結(jié)果與實驗結(jié)果的最大誤差出現(xiàn)在12號線餅4號測點上，為7.13 K，與Fluent相比，本文算法更加接近于實際結(jié)果，考慮實驗過程中的不可控因素，該誤差在可接受范圍內(nèi)。為進一步說明算法的精度，以12號線餅各測溫點為例，將實驗過程中，降階計算與全階計算在各個時間點根據(jù)式（20）計算得到的計算誤差記錄如圖15所示。

由于采用前1 000 s的計算結(jié)果作為樣本數(shù)據(jù)，因此，前100個時間步的計算誤差為0。可以看出，SPOD-DEIM降階計算過程中，與全階計算相比最大計算誤差不超過2.5%，這進一步說明了本文所采用的基于SPOD-DEIM降階方法的混合有限元算法的高精確性。為了分析采用SPOD-DEIM降階算法之后，模型計算的效率提升，將計算時間記錄見表4。

圖14　瞬態(tài)時各線餅監(jiān)測點溫度變化

圖15　12線餅各線匝監(jiān)測點降階計算誤差

表4 計算時間對比

Tab.4 Calculation time comparison（單位：h）

表4中，預(yù)處理時間表示形成降階正交基與插值矩陣的時間，在降階計算中，二者均只需形成一次便可用于全過程瞬態(tài)計算；SPOD-DEIM降階模型中的全階計算時間表示形成快照矩陣的時間，該算例中將0～100時間步的全階計算結(jié)果組成快照矩陣。由表4可知，對實驗的全過程采用全階模型計算，計算時間約為573.59 h，而若采用基于SPOD-DEIM的降階算法計算，則包括預(yù)處理時間在內(nèi)的計算時間僅為54.28 h，較全階計算模型有大幅度提高，計算效率為后者的10.57倍。且相較于商業(yè)仿真軟件Fluent的單核串行運算，該算法的效率提升也在6.37倍以上。由此可說明，基于SPOD-DEIM的降階算法針對實際工程模型也有較好的應(yīng)用價值，能夠在可接受的誤差范圍內(nèi)有效地提高溫升計算效率。同時與3.1.3節(jié)中得到的結(jié)論相比也進一步說明隨著實際工程模型規(guī)模的逐漸增大，降階計算模型的效率提升也將更加明顯。

4 結(jié)論

本文針對油浸式電力變壓器繞組瞬態(tài)溫升過程，采用LSFEM-UFEM混合有限元方法建立了二維繞組溫升數(shù)值模型，研究并討論了將SPOD-DEIM算法應(yīng)用于其中的潛在能力，相關(guān)計算及實驗結(jié)果表明：

1） 精度和效率方面，本文根據(jù)油浸式電力變壓器繞組基本結(jié)構(gòu)，建立了二維繞組單分區(qū)分匝繞組傳熱模型，基于SPOD-DEIM混合有限元方法對繞組的瞬態(tài)溫升過程進行計算，結(jié)果表明：在1 000 s的仿真時間中，降階計算的結(jié)果與全階計算幾乎完全相同，其流場最大誤差不超過1.26%，溫度場最大誤差不超過0.56%，同時，效率方面降階計算效率較全階計算效率提升了5倍以上，且相較于仿真軟件Fluent的單核串行運算，該算法的計算效率也提升了近2倍。這初步說明了SPOD-DEIM降階算法的有效性。

2）工程應(yīng)用方面，本文基于110 kV油浸式電力變壓器繞組基本結(jié)構(gòu)，搭建了溫升實驗平臺，并基于平臺繞組建立了二維八分區(qū)分匝繞組傳熱模型進行數(shù)值計算，將實驗結(jié)果、SPOD-DEIM混合有限元算法計算結(jié)果、Fluent仿真結(jié)果進行對比，結(jié)果表明：精度方面，在繞組的瞬態(tài)溫升實驗過程中，與實驗結(jié)果相比，基于SPOD-DEIM混合有限元算法在監(jiān)測點的最大誤差僅為5.41 K，考慮實驗過程中的不可控因素，該誤差在可接受范圍內(nèi)；效率方面，全階計算模型的全過程計算時間為573.59 h，而SPOD-DEIM降階算法包括預(yù)處理時間在內(nèi)的計算時間僅為54.29 h，后者較前者效率提升約為10.57倍。同時，相較于同等網(wǎng)格規(guī)模的Fluent單核串行計算，本文所采用的算法效率提升也在6.37倍以上。由此可充分說明本文所提算法具有較高的工程應(yīng)用價值。

[1] 劉云鵬, 李歡, 高樹國, 等. 分布式光纖傳感在大型變壓器溫度和繞組變形監(jiān)測中的應(yīng)用研究[J]. 中國電機工程學報, 2022, 42(16): 6126-6135.

Liu Yunpeng, Li Huan, Gao Shuguo, et al. Research on application of distributed optical fiber sensing in monitoring of temperature and winding deformation of large transformer[J]. Proceedings of the CSEE, 2022, 42(16): 6126-6135.

[2] 徐征宇, 張書琦, 廖和安, 等. 傳感光纖與變壓器電磁線一體化技術(shù)[J]. 中國電機工程學報, 2021, 41(19): 6816-6826.

Xu Zhengyu, Zhang Shuqi, Liao Hean, et al. Integration technology of sensing optical fiber and transformer electromagnetic wire[J]. Proceedings of the CSEE, 2021, 41(19): 6816-6826.

[3] Liu Yunpeng, Li Xinye, Li Huan, et al. Experimental and numerical investigation of the internal temperature of an oil-immersed power transformer with DOFS[J]. Applied Sciences, 2020, 10(16): 5718.

[4] 唐釗, 劉軒東, 陳銘. 考慮流體動力學的干式變壓器熱網(wǎng)絡(luò)模型仿真分析[J]. 電工技術(shù)學報, 2022, 37(18): 4777-4787.

Tang Zhao, Liu Xuandong, Chen Ming. Simulation analysis of dry-type transformer thermal network model considering fluid dynamics[J]. Transactions of China Electrotechnical Society, 2022, 37(18): 4777-4787.

[5] 駱仁松, 汪濤, 文繼峰, 等. 大功率高頻變壓器三維溫升計算及優(yōu)化設(shè)計方法[J]. 電工技術(shù)學報, 2023, 38(18): 4994-5005, 5016.

Luo Rensong, Wang Tao, Wen Jifeng, et al. Three dimensional temperature rise calculation and optimization design method for high-power high-frequency transformers[J]. Transactions of China Electrotechnical Society, 2023, 38(18): 4994-5005, 5016.

[6] 彭道剛, 陳躍偉, 錢玉良, 等. 基于粒子群優(yōu)化-支持向量回歸的變壓器繞組溫度軟測量模型[J]. 電工技術(shù)學報, 2018, 33(8): 1742-1749, 1761.

Peng Daogang, Chen Yuewei, Qian Yuliang, et al. Transformer winding temperature soft measurement model based on particle swarm optimization-support vector regression[J]. Transactions of China Electrote-chnical Society, 2018, 33(8): 1742-1749, 1761.

[7] 劉剛, 郝世緣, 胡萬君,等. 基于SCAS時間匹配算法油浸式變壓器繞組瞬態(tài)溫升計算方法[J/OL]. 電工技術(shù)學報, 2023: 1-14. https://doi.org/10.19595/ j.cnki.1000-6753.tces.222137.

Liu Gang, Hao Shiyuan, Hu Wanjun, et al. Transient temperature rise calculation of oil immersed transformer winding based on SCAS time matching algorithm[J]. Transactions of China Electrotechnical Society, 2023: 1-14. https://doi.org/10.19595/j.cnki.1000-6753.tces.222137.

[8] 李永建, 閆鑫笑, 張長庚, 等. 基于磁-熱-流耦合模型的變壓器損耗計算和熱點預(yù)測[J]. 電工技術(shù)學報, 2020, 35(21): 4483-4491.

Li Yongjian, Yan Xinxiao, Zhang Changgeng, et al. Numerical prediction of losses and local overheating in transformer windings based on magnetic-thermal-fluid model[J]. Transactions of China Electrotechnical Society, 2020, 35(21): 4483-4491.

[9] 王澤忠, 李明洋, 宣夢真, 等. 單相四柱式變壓器直流偏磁下的溫升試驗及仿真分析[J]. 電工技術(shù)學報, 2021, 36(5): 1006-1013.

Wang Zezhong, Li Mingyang, Xuan Mengzhen, et al. Temperature rise test and simulation of single-phase four-column transformer under DC-bias[J]. Transactions of China Electrotechnical Society, 2021, 36(5):1006-1013.

[10] Hendriana D, Bathe K J. On upwind methods for parabolic finite elements in incompressible flows[J]. International Journal for Numerical Methods in Engineering, 2000, 47(1/2/3): 317-340.

[11] 靳立鵬, 劉剛, 任增強, 等. 無量綱最小二乘有限元法GPU實現(xiàn)及其在變壓器繞組流場仿真中的應(yīng)用研究[J/OL]. 華北電力大學學報(自然科學版), 2022: 1-9. http://kns.cnki.net/kcms/detail/13.1212.TM. 20220715.0954.002.

Jin Lipeng，Liu Gang，Ren Zengqiang，et al. Research on dimensionless least square finite element method GPU implementation and its application in transformer winding flow field simulation[J/OL]. Journal of North China Electric Power University (Natural Science Edition): 2022: 1-9. http://kns.cnki. net/kcms/detail/13.1212.TM.20220715.0954.002.

[12] 謝裕清, 李琳, 宋雅吾, 等. 油浸式電力變壓器繞組溫升的多物理場耦合計算方法[J]. 中國電機工程學報, 2016, 36(21): 5957-5965.

Xie Yuqing, Li Lin, Song Yawu, et al. Multi-physical field coupled method for temperature rise of winding in oil-immersed power transformer[J]. Proceedings of the CSEE, 2016, 36(21): 5957-5965.

[13] 謝裕清, 李琳, 王帥兵. 瞬態(tài)Navier-Stokes方程的最小二乘有限元降階計算方法[J]. 華北電力大學學報(自然科學版), 2017, 44(3): 74-80.

Xie Yuqing, Li Lin, Wang Shuaibing. Reduced-order method for the transient Navier-Stokes equation based on least-squares finite element method[J]. Journal of North China Electric Power University, 2017, 44(3): 74-80.

[14] 謝裕清, 李琳. 基于本征正交分解的換流變壓器極性反轉(zhuǎn)電場降階模型[J]. 中國電機工程學報, 2016, 36(23): 6544-6551.

Xie Yuqing, Li Lin. A reduced order model via proper orthogonal decomposition for polarity reverse electric fields in converter transformers[J]. Proceedings of the CSEE, 2016, 36(23): 6544-6551.

[15] 劉剛, 榮世昌, 武衛(wèi)革, 等. 基于混合有限元法和降階技術(shù)的油浸式變壓器繞組2維瞬態(tài)流-熱耦合場分析[J]. 高電壓技術(shù), 2022, 48(5): 1695-1704.

Liu Gang, Rong Shichang, Wu Weige, et al. Two-dimensional transient flow-thermal coupling field analysis of oil-immersed transformer windings based on hybrid finite element method and reduced-order technology[J]. High Voltage Engineering, 2022, 48(5): 1695-1704.

[16] Yan Shuai, Tang Zuqi, Henneron T, et al. Structure-preserved reduced-order modeling for frequency-domain solution of the Darwin model with a gauged potential formulation[J]. IEEE Transactions on Magnetics, 2020, 56(1): 1-4.

[17] Lall S, Krysl P, Marsden J E. Structure-preserving model reduction for mechanical systems[J]. Physica D: Nonlinear Phenomena, 2003, 184(1/2/3/4): 304-318.

[18] Chaturantabut S, Sorensen D C. A state space error estimate for POD-DEIM nonlinear model reduction[J]. SIAM Journal on Numerical Analysis, 2012, 50(1): 46-63.

[19] Jiang Bonan. The least-squares finite element method: theory and applications in computational fluid dynamics and electromagnetics[M]. Berlin, Heidelberg: Springer Berlin Heidelberg, 1998.

[20] 陳剛, 呂計男, 龔春林, 等. 計算流固耦合力學[M].北京: 科學出版社, 2021.

[21] Sato Y, Igarashi H. Model reduction of three-dimensional eddy current problems based on the method of snapshots[J]. IEEE Transactions on Magnetics, 2013, 49(5): 1697-1700.

[22] 廖才波, 阮江軍, 劉超, 等. 油浸式變壓器三維電磁-流體-溫度場耦合分析方法[J]. 電力自動化設(shè)備, 2015, 35(9): 150-155.

Liao Caibo, Ruan Jiangjun, Liu Chao, et al. Comprehensive analysis of 3-D electromagnetic-fluid-thermal fields of oil-immersed transformer[J]. Electric Power Automation Equipment, 2015, 35(9): 150-155.

Reduced Order Calculation Method of Steady Temperature Rise of Oil Immersed Power Transformer

Liu Gang1Hu Wanjun1Hao Shiyuan1Liu Yunpeng1Li Lin2

（1. Hebei Provincial Key Laboratory of Power Transmission Equipment Security Defense North China Electric Power University Baoding 071003 China 2. State Key Laboratory of Alternate Electrical Power System with Renewable Energy Sources North China Electric Power University Beijing 102206 China）

In order to improve the efficiency problem of using finite element method to calculate the transient temperature rise of oil immersed power transformer windings, this paper proposes a calculation strategy that combines structure preserved property orthogonal decomposition (SPOD) with discrete empirical interpolation method (DEIM). Firstly, the article uses the least squares finite element method (LSFEM) and upwind finite element method (UFEM) to construct the control equation for calculating the transient temperature rise of transformer windings. The finite element method combined with LSFEM-UFEM can accurately obtain the flow field and temperature field distribution during the transient temperature rise change process of the winding; However, in practical engineering applications, the mesh size and number of nodes of the model are generally large, and using finite element method to solve faces the problem of high time cost. Therefore, in order to reduce the computational scale and improve computational efficiency, the article introduces a reduction method based on the characteristics of the control equation. The traditional POD method usually decomposes the entire dataset and only reflects the statistical characteristics of the dataset, However, the inability to preserve physical features results in low interpretability of the reduced order model. Therefore, this article adopts a reduced order calculation strategy of "data structure preservation" to improve the efficiency of solving finite element equations; At the same time, in order to improve the efficiency improvement of the intrinsic orthogonal decomposition method for nonlinear problems, the DEIM algorithm is combined to interpolate the nonlinear terms in the finite element equation, thereby reducing the time for forming the overall stiffness matrix at each step and further improving the overall computational efficiency. In order to verify the accuracy and efficiency of the algorithm proposed in the article, a single zone split turn winding heat transfer model was established and its transient heat transfer process was calculated. The results showed that the finite element reduced order calculation based on SPOD-DEIM can effectively improve the calculation efficiency while ensuring accuracy. Compared with the full order calculation results, the calculation errors of the flow field and temperature field are not more than 1.5%, and the calculation efficiency is increased by 5.1 times. At the same time, in order to fully demonstrate the value of SPOD-DEIM algorithm in engineering applications, the article built a temperature rise test platform based on 110 kV transformer windings, and established an eight zone split turn winding numerical calculation model to verify and discuss the accuracy, efficiency, and engineering application value of the algorithm. The calculation and experimental results show that in terms of accuracy, the transient whole process calculation error of reduced order calculation is less than 2.5% compared to full order calculation, and compared with experimental results, The error shall not exceed 5.41 K; In terms of efficiency, the entire process calculation time of reduced order calculation is 54.28 h, which increases the calculation efficiency to 10.57 times compared to full order calculation and 6.37 times compared to commercial simulation software Fluent. This fully demonstrates the efficiency and engineering application value of the algorithm proposed in this article.

Transient temperature rise of winding, structure-preserved, proper orthogonal decomposition (SPOD), discrete empirical interpolation method (DEIM), reduced order calculation, temperature rise experiment

TM411

10.19595/j.cnki.1000-6753.tces.222155

國家重點研發(fā)計劃（2021YFB2401703）和中央高?；究蒲袠I(yè)務(wù)費專項資金（2022MS073）資助項目。

2022-11-05

2023-04-03

劉 剛 男，1985年生，副教授，碩士生導(dǎo)師，研究方向為電氣設(shè)備多物理場建模及仿真、電力系統(tǒng)時域仿真和電磁場理論及其應(yīng)用等。E-mail：liugang_em@163.com（通信作者）

胡萬君 男，1999年生，碩士研究生，研究方向為多物理場快速計算等。E-mail：1366320104@qq.com

（編輯 郭麗軍）