基于粒子群優化算法的Jiles-Atherton磁滯模型參數計算*1

郝曉亮， 葉美盈

(浙江師范大學 數理與信息工程學院，浙江 金華 321004)

基于粒子群優化算法的Jiles-Atherton磁滯模型參數計算*1

郝曉亮， 葉美盈

(浙江師范大學 數理與信息工程學院，浙江 金華 321004)

粒子群優化算法是一種模擬鳥群捕食行為的群體智能算法，該算法具有簡潔、易于實現、沒有太多調整參數及不需要梯度信息等特點，且在大多數情況下可快速收斂于最優解．為了描述材料的磁滯特性，提出了一種粒子群優化算法結合MATLAB/Simulink動態仿真集成環境的Jiles-Atherton磁滯回線模型參數計算方法，并分別以無噪及加噪的仿真數據對2組參數值不同的Jiles-Atherton磁滯回線模型進行了數值實驗．結果表明，將粒子群優化算法及MATLAB/Simulink動態仿真集成環境應用于Jiles-Atherton磁滯模型的參數計算是有效的．

磁滯回線；粒子群優化算法；磁滯模型；參數計算

0 引 言

磁滯現象是物理系統和電磁設備中常見的一種非線性特性，它對設備是否能穩定運行和系統的安全運轉起著重要的作用．例如硅鋼片的磁滯特性對變壓器偏磁飽和時的運行性能有很大的影響．因此，了解和分析磁滯特性具有重要意義．

根據實驗數據建立磁滯模型，并以此進行數值模擬是研究磁滯特性的重要內容之一．目前,研究者已提出了多種描述磁滯現象的模型，常見的有：Bouc-Wen磁滯模型[1-3]、Preisach磁滯模型[4]、Jiles-Atherton(JA)磁滯模型[5]．其中，Bouc-Wen磁滯模型的缺點是各參數之間的關系不直觀，物理意義不明確，且表達式以導數或微分形式出現，增大了系統的維數，使分析的難度加大[6]．Preisach磁滯模型則是一種偏數學意義的模型，Preisach權重函數的確定及相關參數的實驗測量難度較大，以致求解困難[7-9]，且數學表達式也太過復雜．Preisach磁滯模型的另一缺陷是:若考慮各向異性、頻率諸因素，Preisach模型將變得極為復雜，從而給數值計算和計算機編程等帶來諸多麻煩．JA磁滯模型是依據物理機理推導而得的一種模型，該模型具有清晰的物理意義，能夠真實地描述磁場強度H與磁感應強度B之間(即B-H)的非線性關系，通過求解JA磁滯模型方程能得到準確的B-H磁滯回線，且參數也相對較少，特別適用電磁器件的分析[10-11]．因此，JA磁滯模型是目前最為常用的一種磁滯模型[5]．

值得注意的是，在建立JA磁滯模型的過程中，雖然模型參數的近似值只需要借助飽和磁滯回線的測量數據即可獲得，但因磁滯非線性的復雜性，如何精確計算這些參數一直是個不易解決的關鍵問題．在過去幾十年里，許多傳統的非線性優化方法已被用于解決JA模型參數計算問題,然而這些優化方法易受初始值選擇的影響，尋優過程易陷入局部極小．為了克服這一缺點，近年來迅速發展的智能優化算法吸引了人們的注意，例如，遺傳算法(Genetic Algorithm，GA)因其良好的全局收斂特性，易于得到比傳統優化方法更好的參數計算結果．文獻[12-19]已成功地將遺傳算法應用于JA模型的參數計算，但GA算法的“早熟”現象及接近最優點收斂較慢的問題限制了其進一步推廣應用[3]．

粒子群優化(Particle Swarm Optimization，PSO)算法是基于觀察動物的社會行為，例如模擬鳥類、魚群、蟻群捕食行為的智能優化算法，該算法最早由美國電氣工程師Eberhart和社會心理學家Kennedy于1995年提出[20]．相比于GA算法，PSO算法有一些更吸引人的特性[21]，其主要特點是具有記憶，因而好的解決方案信息易被保留至下一代.而在GA算法中，先前的信息會隨著種群的改變而破壞．此外，PSO算法得益于粒子間的合作，種群中的粒子之間可相互共享信息．目前，PSO算法在眾多鄰域都引起了很大的關注，已在醫療診斷[22]、神經網絡訓練[23]、太陽能電池參數提取[24]、電氣設備的功率反饋和電壓控制[25]、機器人路徑規劃[26]等方面獲得了廣泛應用．為了得到更有效的JA磁滯模型建模方法，本文嘗試在MATLAB/Simulink動態仿真集成環境中，將PSO算法應用于JA磁滯模型的參數計算,通過Simulink動態仿真JA模型，并應用PSO算法最小化磁感應強度的實驗值和計算值之間的差別,指導JA模型的參數向真值逼近．

1 JA模型

JA模型考慮了實際材料受雜質、非磁性包含物、不均勻壓力及顆粒、孔洞等因素的影響，將磁化強度分解為不可逆磁化分量和可逆磁化分量2部分．其中:不可逆磁化分量緣于疇壁取代，即釘扎的阻礙作用；而可逆磁化分量則緣于疇壁彎曲，從而得到了關于磁滯現象的一個簡單微分表達式．JA模型參數較少,僅有一個一階微分方程,計算量低，并且包含了磁滯現象的主要特征：初始磁化曲線、飽和磁化曲線、矯頑力、剩磁和磁滯損耗等．

設JA模型的不可逆磁化分量為Mirr，可逆磁化分量為Mrev，則總磁化強度為

M=Mirr+Mrev. (1)

式(1)中的不可逆磁化分量Mirr滿足微分方程

式(2)中：H為磁場強度；Man為非磁滯磁化強度；α為平均場參數；k為釘扎阻礙作用引起的磁滯損失參數；μ0為真空磁導率；δ為方向系數．式(2)中的δ與dH/dt有關，當dH/dt>0時(H向正方向增加)，δ=1，而當dH/dt<0時(H向負方向增加)，δ=-1．式(2)中的Man可表示為

式(3)中：a為形狀參數；Ms為飽和磁化強度；He為有效磁場強度且滿足

由于非磁滯磁化位移而導致2個釘扎之間疇壁的膨脹，式(1)中的可逆磁化分量Mrev可表示為

Mrev=c(Man-M). (5)

式(5)中，c為磁疇壁彎曲常數.對式(5)求微分，可得

由式(1)、式(2)、式(6)可得到JA模型磁化強度的微分表達式

式(7)中共涉及5個參數：c,k,α,a,Ms．這些參數共同決定了初始磁化曲線及磁滯回線的形狀特性．

為了得到B-H磁滯回線，需算得磁感應強度B．若已知磁場強度H及磁化強度M，則磁感應強度B為

B=μ0(H+M). (8)

根據實驗數據計算式(7)中的5個參數，本質上是求解優化問題．本文采用粒子群優化算法計算這些參數．計算時，衡量優化性能優劣的目標函數定義為

式(9)中：SV為磁感應強度B的實驗值；OV為磁感應強度B的計算值；N為實驗數據的個數．計算時通過最小化磁感應強度B的實驗值和計算值之間的差別,指導JA模型的5個參數向真值逼近．在本文的數值實驗中，將以仿真值替代實驗值．

2 PSO算法基本原理

基本PSO算法初始化時，首先隨機產生一群粒子(隨機解)，迭代過程中每個粒子在搜索空間中“飛翔”．在每個具體的時刻，每個粒子都有一個位置矢量和速度矢量，其維數代表了問題解空間的維數．若第i個粒子的位置和速度分別記為Xi，Vi，并假設問題的搜索空間是一個d維空間，則粒子的位置和速度可以表示為：Xi=[xi1,xi2,…,xid]，Vi=[vi1,vi2,…,vid]．

在迭代過程中，粒子通過跟蹤2個“極值”不斷調整自己的位置，進行更新．第一個就是粒子本身所找到的最優解，這個解被稱為“個體極值”或個體最優解，第i個粒子的個體最優解記為Pi=[pi1,pi2,…,pid]；另一個極值是整個種群目前所找到的最優解，稱為“全局極值”或“全局最優解”，記為Pg=[pg1,pg2,…,pgd]．粒子根據式(10)來更新自己的速度和位置：

(10)

式(10)中：t為當前迭代次數；c1,c2為常數，稱為加速因子；r1,r2為[0,1]的隨機數；w稱為慣性因子．第j維的位置和速度變化范圍分別為[-xj,max,xj,max]和[-vj,max,vj,max]，迭代中若某一維的xi,j或vi,j超過邊界則可取邊界值．

3 JA磁滯模型的參數計算

為了計算JA模型的參數，本文建立了基于MATLAB/Simulink的JA模型參數計算仿真平臺．在該平臺中，JA磁滯模型參數與B-H磁滯回線之間的關系由Simulink動態仿真集成環境構成的仿真模型確定，即在PSO方法優化程序運行過程中，只要提供磁滯模型的5個參數(c,k,α,a,Ms)，就可以由Simulink仿真模型精確算得B-H磁滯回線數據．考慮式(7)是對磁場強度H求微分的表達式，因此,對其兩邊同乘dH/dt，得

由此可將式(7)轉化為對時間的微分．所以根據式(11),JA模型可由Simulink進行與時間相關的動態仿真．在Simulink動態仿真JA模型時，采用計算精度較高的4階龍格-庫塔(Runge-Kutta)法求解，時間步長為0.005 s，激勵信號頻率為5 rad/s．

在利用本文提出的PSO方法計算JA模型參數時，種群中的每一個粒子被視為JA模型參數的一個潛在解，因此,第i個粒子的位置可表示為xi(ci,ki,αi,ai,MSi)，即每個粒子的位置和速度的維數都是5維的．在計算過程中，以目標函數式(9)作為PSO方法的適應度函數，通過PSO方法自動調節5個參數的值，使磁感應強度B的實驗值(或仿真值)與計算值的偏差最小，從而實現對JA模型參數的計算，得到磁滯回線．

以下是應用PSO及MATLAB/Simulink計算JA模型參數的步驟：

1)設定PSO算法的控制變量，包括種群規模、慣性因子、最大迭代次數、模型參數的維數及搜索空間等;

2)隨機產生模型參數的初始值，即初始化群體中各粒子的位置和速度；

3)調用Simulink形式的JA模型，按式(9)計算粒子的適應度，若得到更優的pij和pgj，則更新并保存當前的最優值；

4)依照式(10)更新每個粒子的速度和位置；

5)若未能達到預先設定的最大迭代次數，則返回3)；否則停止運算，并輸出最優的適應度值(即目標函數值)及對應的優化模型參數值.

4 數值實驗結果與討論

本文數據實驗中所用的JA模型參數真值分別來源于文獻[27]與[28]．其中，文獻[27]給出的是一組低碳鋼SAE1010環形材料(內徑60 mm，外徑79.04 mm，高度9.52 mm)的JA模型參數集(以下稱參數集1)；而由文獻[28]給出的則是一組碳化鐵(Fe3C)材料對應的JA模型參數集(以下稱參數集2)．在本文的JA模型參數計算數值實驗中，參數集1的采樣點數為261個，磁場強度H的最大值為2 600 A/m；參數集2的采樣點數為361個，磁場強度H的最大值為6 100 A/m．本文分別在仿真數據無噪和加噪2種情況下，驗證了PSO方法計算JA磁滯回線模型參數的可行性．

對于PSO的控制變量，本文取種群規模PS=50，最大迭代次數tmax=500，5個模型參數的搜索空間均設定為其真值附近±100%以內的區域．仿真數據無噪時，對于參數集1，取w=1，c1=c2=2；而參數集2，則取w=0.6，c1=c2=1.7．仿真數據加噪時，w在迭代過程中均由0.9到0.4呈線性下降，并取c1=c2=2．

4.1仿真數據無噪時的參數計算

此時，通過直接將參數集1或參數集2加入Simulink的JA模型，并求解該模型即可得到無噪的仿真數據序列．

在用無噪仿真數據計算JA模型參數時，理論上目標函數最小值應該為0，即式(9)中磁感應強度B的仿真值與計算值應相等．這種情況下可認為計算所得的最優JA模型參數與產生仿真數據所用的參數集1或參數集2相同(實際因數值計算誤差會略有不同)．因此，可以通過無噪仿真數據的數值實驗觀測PSO方法的計算精度．

圖1(a)與(b)分別是以參數集1、2利用無噪仿真數據時PSO方法JA模型參數計算所得的磁滯特性．比較圖1中的磁滯回線與原始無噪仿真數據，可見無論是參數集1還是參數集2，對于無噪仿真數據，兩者均能很好地吻合．圖2(a)與(b)則分別展現了PSO方法在參數集1、2的JA模型參數計算過程中，目標函數值F隨迭代次數t的變化情況，可見PSO方法能快速收斂．相關的JA模型參數真值與PSO計算值的比較列于表1，顯然模型參數的計算值與真值幾乎一致．特別地，表1所列參數集1對應的計算值與真值完全相同(實際上若數據位數取足夠多，兩者仍會略有差異)．

圖1 無噪仿真數據時PSO方法JA模型參數計算的磁滯特性

應用無噪仿真數據以4種PSO程序運行模式計算JA模型參數的結果列于表2．其中程序運行的時間tp與軟硬件環境有關，本文數值實驗所用的計算機配置：CPU為INTEL雙核酷睿2，內存為2 GB，操作系統為Windows 7/64位，MATLAB為R2010a版本．由表2可知：對于無噪仿真數據，不同情況的參數計算耗時tp相差并不大，均不需花太多的時間，且計算結束時，目標函數值F除參數集2在PSO程序Model=0和3時略差些(說明正確選擇PSO算法的控制參數是有必要的)，其余情況均很小．總之，應用PSO方法在大多數情況下都可得到誤差較小的計算結果．

圖2 無噪仿真數據時PSO方法目標函數值隨迭代次數的變化

表1 JA模型參數真值與無噪仿真數據時PSO計算值的比較

表2 無噪仿真數據時PSO在不同模式下對JA模型參數集1、2計算結果的比較

4.2仿真數據加噪時的參數計算

考慮到在測量B-H磁滯數據的物理實驗時不可避免地會受到來自測量儀器本身的誤差及外部環境干擾等因素的影響，以致實際得到的B-H實驗數據并不“純凈”，其中必定包含測量噪聲．因此，筆者在此增加了仿真數據加噪后的數值實驗，以進一步驗證PSO方法的有效性．

設在無噪數據序列上所疊加的測量噪聲為

式(12)中：r為一組數據個數與無噪數據序列相同的、均值為0、方差與標準差均為1的正態分布隨機數;Ne則表示所疊加的噪聲水平．在以下數值實驗中分別取Ne=0.015及0.030，以檢驗不同測量噪聲水平的影響．

仿真數據加噪后應用PSO方法計算JA模型參數獲得的磁滯特性見圖3．由圖3(a)～(d)可看出：參數集1、2的磁滯回線與加噪后的仿真數據吻合程度有所降低，顯然噪聲水平越高,兩者的差異就越大，但即使Ne=0.030,其誤差在工程上也仍可接受．可見無論是參數集1還是參數集2，即使對于加噪仿真數據，兩者均能很好地吻合.

圖4是仿真數據加噪后目標函數值隨迭代次數的變化情況，表明此時PSO算法仍能比較快地收斂．表3列出了仿真數據加噪后PSO方法的JA模型參數計算結果．由表3可知，加噪后磁滯模型中的4個參數k,α,a,Ms的計算值與真值仍然接近，但參數c對計算誤差比較敏感，計算值與真值偏離較大．比較表1～表3，總體來說仍是噪聲水平越高,誤差越大，其原因不難理解，仿真數據加噪后目標函數值必然大于0．根據表3可知，加噪后對計算耗時并無影響，與未加噪聲前依然相近．

圖3 不同加噪水平仿真數據時PSO方法JA模型參數計算的磁滯特性

圖4 不同加噪水平仿真數據時PSO方法目標函數值隨迭代次數的變化

表3 仿真數據加噪后PSO方法的JA模型參數計算結果

5 結 論

本文嘗試將PSO算法結合MATLAB/Simulink動態仿真集成環境應用于JA磁滯模型參數計算．數字實驗結果表明，當仿真數據未加噪聲時，應用本文提出的方法可以得到幾乎與模型參數真值一致的計算結果；即使在仿真數據加噪情況下，該方法仍能得到較高的參數計算精度，且由計算所得參數繪制的磁滯回線與仿真數據相吻合，說明該方法有較高的穩定性．因此，本文提出的方法是有效的．

[1]Wen Y K.Method for random vibration of hysteretic systems[J].Journal of the Engineering Mechanics Division,1976,102(2):249-263.

[2]Bouc R.Forced vibration of mechanical systems with hysteresis[C]//Proceedings of the Fourth Conference on non-linear Oscillation.Prague:1967:315.

[3]Ye Meiying,Wang Xiaodong.Parameter estimation of the Bouc-Wen hysteresis model using particle swarm optimization[J].Smart Materials and Structures,2007,16(6):2341.

[4]Mayergoyz I D.Mathematical models of hysteresis[J].IEEE Transactions on Magnetics,1986,22(5):603-608.

[5]Jiles D C,Atherton D L.Theory of ferromagnetic hysteresis[J].Journal of Applied Physics,1984,55(6):2115-2120.

[6]李韶華,楊紹普.滯后非線性模型的研究進展[J].動力學與控制學報,2006,4(1):8-15.

[7]Mayergoyz I D,Adly A A.Numerical implementation of the feedback Preisach model[J].IEEE Transactions on Magnetics,1992,28(5):2605-2607.

[8]Mayergoyz I D,Friedman G.Generalized Preisach model of hysteresis[J].IEEE Transactions on Magnetics,1988,24(1):212-217.

[9]Ge P,Jouaneh M.Generalized preisach model for hysteresis nonlinearity of piezoceramic actuators[J].Precision Engineering,1997,20(2):99-111.

[10]Annakkage U D,McLaren P G,Dirks E,et al.A current transformer model based on the Jiles-Atherton theory of ferromagnetic hysteresis[J].IEEE Transactions on Power Delivery,2000,15(1):57-61.

[11]Leite J V,Sadowski N,Kuo-Peng P,et al.The inverse Jiles-Atherton model parameters identification[J].IEEE Transactions on Magnetics,2003,39(3):1397-1400.

[12]Wilson P R,Ross J N,Brown A D.Optimizing the Jiles-Atherton model of hysteresis by a genetic algorithm[J].IEEE Transactions on Magnetics,2001,37(2):989-993.

[13]Chan Chehang,Liu Guangjun.Hysteresis identification and compensation using a genetic algorithm with adaptive search space[J].Mechatronics,2007,17(7):391-402.

[14]Chwastek K,Szczyglowski J.Identification of a hysteresis model parameters with genetic algorithms[J].Mathematics and Computers in Simulation,2006,71(3):206-211.

[15]Lu Hailiang,Wen Xishan,Lan Lei,et al.A self-adaptive genetic algorithm to estimate JA model parameters considering minor loops[J].Journal of Magnetism and Magnetic Materials,2015,374:502-507.

[16]Leite J V,Avila S L,Batistela N J,et al.Real coded genetic algorithm for Jiles-Atherton model parameters identification[J].IEEE Transactions on Magnetics,2004,40(2):888-891.

[17] Salvini A,Fulginei F R.Genetic algorithms and neural networks generalizing the Jiles-Atherton model of static hysteresis for dynamic loops[J].IEEE Transactions on Magnetics,2002,38(2):873-876.

[18]Cao Shuying,Wang Boweng,Yan Rongge,et al.Optimization of hysteresis parameters for the Jiles-Atherton model using a genetic algorithm[J].IEEE Transactions on Applied Superconductivity,2004,14(2):1157-1160.

[19]Fulginei F R,Salvini A.Softcomputing for the identification of the Jiles-Atherton model parameters[J].IEEE Transactions on Magnetics,2005,41(3):1100-1108.

[20]Kennedy J,Eberhart R.Particle swarm optimization[C]//Purdue school of engineering and technology.Proceedings of the IEEE International Conference on Neural Networks.Perth:IEEE Press,1995:1942-1948.

[21]Slade W H,Ressom H W,Musavi M T,et al.Inversion of ocean color observations using particle swarm optimization[J].IEEE Transactions on Geoscience and Remote Sensing,2004,42(9):1915-1923.

[22]Hsieh Y Z,Su M C,Wang P C.A PSO-based rule extractor for medical diagnosis[J].Journal of Biomedical Informatics,2014,49:53-60.

[23]Kennedy J,Eberhart R C,Shi Y.Swarm Intelligence[M].San Francisco:Morgan Kaufmann,2001.

[24]Ye Meiying,Wang Xiaodong,Xu Yousheng.Parameter extraction of solar cells using particle swarm optimization[J].Journal of Applied Physics,2009,105(9):094502.

[25]Yoshida H,Kawata K,Fukuyama Y,et al.A particle swarm optimization for reactive power and voltage control considering voltage security assessment[J].IEEE Transactions on Power Systems,2000,15(4):1232-1239.

[26]孫波,陳衛東,席裕庚.基于粒子群優化算法的移動機器人全局路徑規劃[J].控制與決策,2005,20(9):1052-1060.

[27]Hernandez E D M,Muranaka C S,Cardoso J R.Identification of the Jiles-Atherton model parameters using random and deterministic searches[J].Physica B:Condensed Matter,2000,275(1/3):212-215.

[28]Jiles D C,Atherton D L.Theory of ferromagnetic hysteresis[J].Journal of Magnetism and Magnetic Materials,1986,61(1/2):48-60.

(責任編輯 杜利民)

ParametercalculationoftheJiles-Athertonhysteresismodelusingparticleswarmoptimizationalgorithm

HAO Xiaoliang, YE Meiying

(CollegeofMathematics,PhysicsandInformationEngineering,ZhejiangNormalUniversity,JinhuaZhejiang321004,China)

Particles swarm optimization algorithm was a kind of swarm intelligence algorithm simulating birds feeding behavior. The algorithm had the advantages of concision, easy implement, few control parameters, did not need the gradient information, and fast convergence to optimal solution in most cases．In order to describe the hysteresis characteristics of the material， it was proposed a method of the particle swarm optimization algorithm combination with MATLAB/Simulink dynamic simulation integration environment to calculate the Jiles-Atherton hysteresis loop model parameters. By means of noise-free and noisy simulation data, the numerical experiments of Jiles-Atherton hysteresis loop model with the two groups of different parameter values were carried out. The results indicated that the particle swarm optimization algorithm combination with MATLAB/Simulink dynamic simulation integration environment was an effective technique for the parameters calculation of Jiles-Atherton hysteresis model．

hysteresis loop; particle swarm optimization algorithm; Jiles-Atherton hysteresis model; parameter calculation

10.16218/j.issn.1001-5051.2015.02.003

2015-01-12

國家自然科學基金資助項目(51305407)

郝曉亮(1988-)，男，山西晉中人，碩士研究生．研究方向：理論物理．

葉美盈．E-mail： ymy@zjnu.cn

TM936.3

A

1001-5051(2015)02-0133-09