同軸圓環(huán)永磁體磁力計算與實驗對比研究

程書晗，高 琳，蘇宇鋒，王天航

（鄭州大學機械與動力工程學院，河南 鄭州 450001）

1 引言

隨著稀土永磁材料的發(fā)展，得益于磁力機械結(jié)構(gòu)簡單、環(huán)保、非接觸、損耗小等優(yōu)點，磁力在機械系統(tǒng)中有了越來越廣泛的應用。磁力齒輪［1］、磁力泵［2］和磁力彈簧［3-4］等磁力機械都是通過永磁鐵或電磁鐵的磁場耦合實現(xiàn)力或力矩的傳遞，通過無接觸傳動的方式可以有效解決傳統(tǒng)機械結(jié)構(gòu)存在的問題。電磁式能量采集裝置［4-10］也得到了廣泛的研究，這種裝置一般都通過不同形狀、大小的永磁體布置使磁場相互耦合，在外界激勵的作用下產(chǎn)生相應的響應，利用法拉第電磁感應原理進行能量的采集。為了對磁力機械進行結(jié)構(gòu)設計和響應分析需要引入磁力的計算，磁力的計算一般抽象為兩個永磁體之間力的計算。文獻［11］基于磁路分析，使用虛位移方法得出了軸向混合磁軸承的徑向承載力解析數(shù)學模型，得到了永磁體徑向承載力與徑向位移的關(guān)系。文獻［12-13］采用等效磁荷法模型對同軸或偏心時兩圓環(huán)磁體之間的磁力進行了建模，得到含有四重積分的磁力計算公式。文獻［14］通過磁偶極子之間的力模型計算兩磁體之間的磁力。文獻［15］對文獻［12］的計算公式進行簡化并增加校正項，得到相對精確的同軸圓柱磁體之間磁力計算結(jié)果。文獻［16-17］通過等效電流法對同軸圓柱磁體之間磁力進行求解，得到含有四重積分的解析表達式。

引入了文獻［18］的磁場求解方法，對其進行補充分析，使用等效電流法得到圓環(huán)永磁體含有廣義完全橢圓積分的磁場求解公式，分別使用自適應辛普森法則和文獻［19］中的數(shù)值解法進行Matlab 編程計算，得到圓環(huán)永磁體在空間中任意位置產(chǎn)生的磁場。進一步對同軸布置的兩個圓環(huán)磁體進行建模求解磁力，通過與COMSOL中仿真結(jié)果進行對比驗證模型的正確性，最后通過該模型在一種抗磁懸浮能量采集器［8］中的應用，在實驗中驗證這種磁力計算方法的準確性。

2 等效電流法進行磁場計算

2.1 圓環(huán)永磁體磁場計算

根據(jù)等效電流模型，均勻磁化的圓環(huán)磁體內(nèi)部等效圓電流相互抵消，只剩余圓環(huán)邊界上的電流。圓環(huán)厚度為H，外圓半徑為R，內(nèi)圓半徑為R’，磁場強度矢量為M，磁化方向為軸向，如圖1所示。

圖1 磁場求解示意圖Fig.1 Schematic Diagram of Magnetic Field Solution

把坐標系Oxyz建立在圓環(huán)磁體中心，圓環(huán)磁體等效為一個個邊界電流環(huán)，每個邊界電流環(huán)的寬度為dz'，電流環(huán)L的參數(shù)方程可以表示為：

設單位寬度上電流匝數(shù)為n，外圈每匝電流為I，所以外圈每個電流環(huán)的電流大小為nIdz'。內(nèi)圈每匝電流為I'，有I'=-I，內(nèi)圈每個電流環(huán)的電流大小為nI'dz'。設圓環(huán)頂面外圓面積為S，內(nèi)圓面積為S'，則有：MH(S-S')=ISnH+I'S'nH即：

因此外圈、內(nèi)圈電流環(huán)電流分別可以寫成Mdz'、-Mdz'。

對圓環(huán)外側(cè)面進行求解，P(x'，y'，z')為電流環(huán)上任意一點，則可用柱坐標表示向量OP為(Rcosφ'，Rsinφ'，z')。設待求解磁場點Q(x，y，z)，則根據(jù)畢奧-薩伐爾定律，首先對單個電流環(huán)每個位置進行積分，再對所有軸向上所有電流環(huán)累加磁場進行求解，得到如下等式：

式中：μ0—真空磁導率。

在電流環(huán)L上的P點位置，dL∕dφ'可以表示為：

為方便求解，調(diào)整坐標系的x軸位置，使待求點Q在xOz平面上，則向量OQ可用柱坐標表示為(ρ，0，z)。

將dL、OP和OQ代入方程（2）中，且將線積分轉(zhuǎn)換為對φ'的積分后可以得到：

由于電流環(huán)的對稱性，且點Q在xOz平面上，對于電流恒定的電流環(huán)來說不存在軸向的電流和位移電流，所以Q點在y軸方向的磁感應強度為0，即對于柱坐標來說，切向磁感應強度大小為0。因而可以舍掉等式（3）中的y方向磁場，可將等式（3）分解為徑向磁場Bρ(Q)和軸向磁場Bz(Q)。

將等式（4）和等式（5）對z'以及φ'積分，擴展到空間中任意一點可以得到：

廣義完全橢圓積分：

同理推出圓環(huán)內(nèi)側(cè)面產(chǎn)生的磁場。得到圓環(huán)產(chǎn)生的總磁場Bρtot和Bztot為：

2.2 廣義完全橢圓積分求解

對于圓環(huán)永磁體所產(chǎn)生的磁場為式（14）和式（15），從解析式中可以看出，其中包含廣義完全橢圓積分項，想要求解該方程，需要對積分項進行求解。對于廣義完全橢圓積分，一般采用數(shù)值法進行求解，這里分別采用文獻［19］中的算法與自適應辛普森法則［20］進行求解。

自適應辛普森法則是一種對常規(guī)積分求解的數(shù)值解法，也可應用于橢圓積分。與之相比，文獻［19］提出的算法是一種基于拉登變換，利用高斯和拉格朗日等人提出的算術(shù)幾何平均的快速收斂尺度，只適用于求解橢圓積分的快速收斂算法。通過對兩種方法在Matlab中編程對比發(fā)現(xiàn)，對于同樣的求解精度，快速收斂算法求解速度遠快于自適應辛普森法，為提升后續(xù)磁力求解速度與求解精度，采用快速收斂算法求解磁場。

3 磁力計算推導

根據(jù)圓環(huán)磁體的磁場計算式（14）與式（15），進行兩同軸圓環(huán)磁體之間磁力的推導。

如圖2所示，兩個圓環(huán)磁體同軸，磁體1外圈半徑為R1，內(nèi)圈半徑為R1'，厚度為H1，磁場強度矢量為M1；磁體2 外圈半徑為R2，內(nèi)圈半徑為R2'，厚度為H2，磁場強度矢量為M2；兩磁體磁化方向均為軸向，且磁體間距為d。待求磁體2對磁體1的磁力。

圖2 圓環(huán)磁體示意圖Fig.2 Schematic Diagram of Ring Magnet

由于對于永磁體的參數(shù)已知的經(jīng)常是磁體的剩余磁通密度Br，為簡化計算，當不考慮兩磁體之間磁場相互作用的對磁化的影響時，可以通過M=Br∕u0推導出磁體的磁化強度M。

由于磁場具有疊加性，磁體1所受的力為：

式中：Fij—磁體i所產(chǎn)生的磁場對磁體j的力。易得F11=0，所以磁體1所受的力為磁體2產(chǎn)生的磁場對其上等效電流的力。

將坐標軸建立在磁體2中心處將磁體1軸向劃分為一個個邊界電流環(huán)，對應的外環(huán)與內(nèi)環(huán)分別為L1、L'1，根據(jù)安培力公式，單個電流環(huán)受磁力為：

式中：n1—磁體1單位長度電流環(huán)匝數(shù)；dz1—電流環(huán)寬度；I1—每匝外圈電流大小；I2—每匝內(nèi)圈電流大小；B2—磁體2在磁體1處磁感應強度。

由等式（6）和等式（7）推導出磁體2在磁體1側(cè)面上各位置磁感應強度：

電流環(huán)在軸向和徑向磁場中的位置示意圖，如圖3所示。

圖3 電流環(huán)受力示意圖Fig.3 Schematic Diagram of Current Loop Force

由圖3（a）可知，徑向磁場對電流環(huán)產(chǎn)生的安培力合力方向垂直電流環(huán)平面。由圖3（b）可知，軸向磁場所產(chǎn)生的安培力方向都為徑向，且由于電流環(huán)的中心對稱性，所有軸向磁場所產(chǎn)生的安培力合力為0。所以將等式（18）帶入等式（17），取z軸方向為正方向，并將軸向磁場舍去，對圓環(huán)積分后可以得到：

將M1=n1I1=-n1帶入等式（19），對z1進行積分后可以獲得磁體2對磁體1的磁力F：

等式（20）是一個含有二重積分的解析式，需要使用數(shù)值法進行求解，內(nèi)部的廣義橢圓積分可以使用快速收斂法進行求解，外部對z1的積分可以使用自適應辛普森法則進行求解，這樣只相當于一重積分的求解，與含有四重積分的磁力求解公式相比，速度大有提高。采用Matlab數(shù)值計算編程語言利用等式（20）可以對同軸圓環(huán)磁體間的相互作用力計算求解。

4 COMSOL磁力仿真結(jié)果對比研究

本節(jié)分別在有限元軟件COMSOL中進行磁力仿真分析，與Matlab求解出來的磁力結(jié)果進行比較。

4.1 COMSOL中磁場及磁力計算方法

在COMSOL中進行永磁體之間磁力的仿真需要使用磁場、無電流接口。在COMSOL軟件中內(nèi)置的磁場求解方法為等效磁荷模型與有限元相結(jié)合。通過

式中：Vm—磁標勢；Br—剩余磁化強度；μrec—回復磁導率。

COMSOL中對磁力的求解采用Maxwell應力張量法。

式中：n1—表面的法向量；H—磁場強度；T2—空氣的應力張量。通過磁體所有邊界表面應力張量的積分可以求得磁力為：

式中：S—磁體邊界表面。

表1 圓柱磁體組參數(shù)Tab.1 Cylindrical Magnet Group Parameters

表2 圓環(huán)磁體組參數(shù)Tab.2 Ring Magnet Group Parameters

使用COMSOL進行仿真時，設置各組圓環(huán)磁體的尺寸，為保證求解精度，空氣域半徑統(tǒng)一設置為55mm，無限元域半徑設置為25mm，為保證網(wǎng)格質(zhì)量，對圓環(huán)磁體進行結(jié)構(gòu)化網(wǎng)格劃分，根據(jù)各組具體尺寸調(diào)節(jié)結(jié)構(gòu)化網(wǎng)格最大單元大小，保證收斂精度以及平衡運算時間，其他部分采用自由四面體網(wǎng)格，網(wǎng)格大小采用超細化，網(wǎng)格劃分后模型，如圖4所示。穩(wěn)態(tài)求解器的相對容差設置為1e-6。

圖4 COMSOL仿真模型Fig.4 COMSOL Simulation Model

由于使用COMSOL進行磁力仿真時，會有迭代誤差和空氣域和無限元域網(wǎng)格劃分所產(chǎn)生的誤差，使用文獻［21］中的研究方法進行校正，通過計算含誤差力與誤差力之差，得到兩磁體之間真實作用力，增加仿真結(jié)果的準確度。

4.2 兩種計算方法的對比

經(jīng)過仿真與計算后結(jié)果對比，如圖5所示。可以看出由等效電流法求解出的兩圓環(huán)永磁體之間的磁力與使用COMSOL仿真得到的磁力結(jié)果非常相近。使用等式（20）擬合COMSOL中求解的離散點，六組數(shù)據(jù)得到的R2值分別為：0.9958，0.9999，0.9998，0.99997，0.9991，0.9996。可以看出擬合效果非常好，關(guān)于圓環(huán)磁體之間的磁力求解，使用等式（20）完全可以代替COMSOL軟件建模仿真。在求解速度上，同樣精度的情況下，使用COMSOL求解一個位置磁力的平均時間約2min，使用Matlab 模型進行求解的平均時間小于1ms，極大地提高了計算效率。

圖5 各組磁體不同間距下磁力對比圖Fig.5 Comparison of the Magnetic Force of Each Group of Magnets at Different Distances

5 實驗驗證

為了對磁力求解公式進行進一步驗證，對一種抗磁懸浮振動能量采集結(jié)構(gòu)［22］進行分析，該結(jié)構(gòu)主要由一個圓環(huán)提升磁鐵、一個圓柱懸浮磁鐵以及一對抗磁體構(gòu)成，通過提升磁體對懸浮磁鐵的吸引力平衡懸浮磁鐵的重力判斷懸浮間距，將該懸浮間距求解問題抽象為磁力求解模型，分別進行模型求解和實驗求解的對比。實驗平臺，如圖6所示。實驗采用15組數(shù)據(jù)進行計算，為保證實驗數(shù)據(jù)的完備性，實驗中包含提升磁體和懸浮磁體分別為圓柱或圓環(huán)的四種情況。實驗使用的磁體具體尺寸，如表3所示。

表3 磁力驗證實驗中磁體參數(shù)Tab.3 Magnet Parameters in the Magnetic Verification Experiment

圖6 實驗平臺Fig.6 Experiment Platform

磁鐵的材料為不同牌號的燒結(jié)釹鐵硼，密度為7.5g∕cm2，提升磁體與懸浮磁體磁化方向相同且為軸向。尋找間距采用二分法，結(jié)合前面磁力計算模型自動尋找懸浮間距，最終懸浮間距的實驗測量結(jié)果和模型計算結(jié)果及誤差，如圖7所示。

圖7 磁體懸浮間距對比Fig.7 Comparison of Magnet Suspension Spacing

通過對比圖可以看出，實驗結(jié)果和模型預測結(jié)果之間有一定誤差，經(jīng)計算十五組數(shù)據(jù)使用公式求解出的懸浮間距與實驗得到的結(jié)果誤差都小于5%，由于實際磁體存在磁化強度不準確，磁化不均勻，實驗磁體布置不對中等問題，該計算結(jié)果在可接受范圍之內(nèi)。

6 結(jié)論

首先通過等效電流法求解同軸圓環(huán)永磁體之間的磁力，得到含有廣義完全橢圓積分的解析式，通過Matlab編程對橢圓積分部分采用快速數(shù)值解法進行求解。其次通過解析方程對不同參數(shù)下通過有限元仿真軟件獲得的仿真結(jié)果進行擬合，得到解析方程的R2值均大于0.995，驗證了解析式的準確性，并對比兩種求解方式的求解速度，分析結(jié)果表明此磁力計算模型速度遠優(yōu)于有限元軟件仿真計算速度。最后將磁力求解模型應用于一種抗磁懸浮能量采集裝置中，用于求解提升磁體與懸浮磁體之間的懸浮間隙，通過建立15組實驗的實驗測量結(jié)果與模型預測結(jié)果進行對比，得到的模型計算誤差均小于5%，進一步驗證了同軸磁體間相互作用力的有效便捷計算方法。這種同軸圓環(huán)磁體磁力計算方法，由于求解速度快，精度較高的特點，可用于基于磁力作用機理的驅(qū)動器和傳感器分析與設計。