張健 陳紹露 張紹華 郭濤 劉超波 肖琦

（1 北京宇航系統(tǒng)工程研究所，北京 100076；2 北京衛(wèi)星環(huán)境工程研究所，北京 100094）

0 引言

磁場模擬裝置在航天器磁性控制[1-4]、艦船消磁等領(lǐng)域有著重要的應(yīng)用[5,6]。磁場線圈是磁場模擬裝置的主要組成部分，中心區(qū)域的磁場均勻度是關(guān)鍵技術(shù)指標(biāo)之一[7-9]，對其性能有直接影響。方形赫姆霍茲線圈是最簡單的磁場線圈形式，加工制作方便，并且可以通過增加幾何尺寸來滿足均勻區(qū)指標(biāo)要求，但是其代價是線圈常數(shù)降低，電源功耗增加，尤其在線圈空間尺寸受到限制時，均勻性指標(biāo)很難滿足設(shè)計要求。所以，在磁場均勻性要求較高，工程安裝空間有限，勵磁電源功耗受限等情況下，可以采用共軸多組線圈組合的方式加以解決。

隨著線圈組數(shù)增多，各組線圈尺寸、安裝位置、線圈匝數(shù)、電纜截面、電感電阻等參數(shù)，存在眾多的組合形式，如果選擇不合理會使得磁場均勻性降低，工程代價增加，甚至不能滿足設(shè)計要求。所以，需要結(jié)合工程可行性約束條件進行優(yōu)化設(shè)計與計算分析，使其綜合指標(biāo)性能達到最優(yōu)。

線圈參數(shù)的計算有多種方法，常用的有高階求導(dǎo)法[10,11]、智能優(yōu)化法[12-16]等。隨著線圈組數(shù)的增加，待優(yōu)化參數(shù)也相應(yīng)增加。張偉等利用磁位的拉普拉斯方程在球坐標(biāo)系下的解，給出了典型的多段環(huán)形線圈的設(shè)計方法和過程[7]；對于高階求導(dǎo)法，4對線圈參數(shù)優(yōu)化需要8階求導(dǎo)[10]，9對參數(shù)優(yōu)化需要16階求導(dǎo)[11]，高階求導(dǎo)法十分復(fù)雜。呂志峰等針對多線圈均勻磁場優(yōu)化設(shè)計中的高階求導(dǎo)及優(yōu)化結(jié)果可信度評估問題，提出一種基于粒子群優(yōu)化算法和有限元法相結(jié)合的多線圈均勻磁場優(yōu)化設(shè)計方法[12]；王志峰等利用差分進化算法對線圈的均勻磁場進行優(yōu)化設(shè)計[13]。

多組線圈參數(shù)的計算問題可以歸結(jié)為求解由物理指標(biāo)（如磁場均勻度指標(biāo)等）和工程約束條件（如線圈安裝空間等）耦合的目標(biāo)函數(shù)。線圈優(yōu)化的目標(biāo)函數(shù)中包含非線性函數(shù)，是非線性規(guī)劃問題。非線性規(guī)劃目前還沒有適用于各種問題的一般算法，已有的各種方法都有其特定的適用范圍。序列二次規(guī)劃（Sequence quadratic programming，SQP）算法是目前公認的求解約束非線性優(yōu)化問題最有效的方法之一[17-19]。與其他算法相比，序列二次規(guī)劃法的優(yōu)點是收斂性好、計算效率高、邊界搜索能力強，因此受到了廣泛的重視及應(yīng)用。本文使用序列二次規(guī)劃方法進行共軸多組充退磁線圈系統(tǒng)優(yōu)化設(shè)計，通過優(yōu)化線圈位置和驅(qū)動電流等參數(shù)來實現(xiàn)所需的磁場均勻度指標(biāo)，并且給出了典型的多組矩形線圈在不同工程邊界約束條件下的線圈優(yōu)化參數(shù)、磁場均勻區(qū)及電參數(shù)的計算結(jié)果。

1 線圈參數(shù)及磁場計算

方形線圈和圓形線圈是工程中常用的兩種線圈形式。方形線圈具有結(jié)構(gòu)簡單，加工制作方便的優(yōu)點。不失一般性，本文以共軸正方形多組線圈為例，對磁場分布的均勻性指標(biāo)進行參數(shù)優(yōu)化研究。

圖1 給出4個方形赫姆霍茲線圈系統(tǒng)及線圈離散單位電流元產(chǎn)生磁場示意圖，坐標(biāo)原點位于充退磁線圈系統(tǒng)的中心點。4個線圈大小相同，邊長均為a。在實際應(yīng)用中為保證勵磁電流的同步性，4個線圈串聯(lián)，通過電流為I，并假定逆著x軸正方向看電流是逆時針方向。需要說明的是，雖然每個矩形線圈的尺寸可以不同，但是考慮實際加工制作及工程應(yīng)用，通常設(shè)置為相同的大小尺寸。為了滿足軸向整體磁場的對稱性，線圈繞制匝數(shù)也成對相等。因此，內(nèi)側(cè)線圈距離d1，線圈匝數(shù)相同，均為N1匝；外側(cè)線圈距離d2，線圈匝數(shù)相同，均為N2匝。預(yù)先固定線圈邊長a，需要優(yōu)化的線圈參數(shù)有三組，即d1，d2和內(nèi)外線圈電流比例（即，ck12=N1/N2）。隨著線圈匝數(shù)增多，相應(yīng)的待優(yōu)化參數(shù)增加。

圖1 4個方形赫姆霍茲線圈系統(tǒng)及線圈離散單位電流元產(chǎn)生磁場示意圖Fig.1 Schematic diagram of 4 square Helmholtz coil systems and coil discrete unit current elements to generate magnetic field

優(yōu)化時需要計算線圈系統(tǒng)在優(yōu)化區(qū)域內(nèi)多個空間點的磁場。對于方形線圈的磁場計算有很多種方法。解析法在中軸線上的磁場計算相對簡單，但是偏離中軸線的磁場矢量解析公式非常復(fù)雜。有限元方法在計算大尺寸空心線圈的三維空間磁場時，對于細導(dǎo)線的網(wǎng)格劃分與空間點的網(wǎng)格劃分差別巨大，導(dǎo)致離散網(wǎng)格點過多、計算耗時甚至解不收斂。本文采用離散積分的方法來計算消磁線圈在空間點的磁場分布。把線圈電流等效為理想線電流模型，線圈離散成微小的電流矢量單元，根據(jù)畢奧薩伐爾公式，每個電流矢量單元在某個空間點P的磁場矢量計算公式如下：

2 優(yōu)化過程

共軸多組充退磁線圈均勻磁場優(yōu)化設(shè)計過程如下：1）確定線圈組數(shù)及待優(yōu)化參數(shù)；2）給出磁場均勻度的目標(biāo)函數(shù)，并結(jié)合工程實際給出物理參數(shù)邊界條件；3）采用序列二次規(guī)劃優(yōu)化算法對目標(biāo)函數(shù)進行迭代尋優(yōu)；4）線圈匝數(shù)取整化，再次尋優(yōu)計算；5）數(shù)值仿真驗證指標(biāo)。

2.1 確定線圈組數(shù)及待優(yōu)化參數(shù)

優(yōu)化設(shè)計首先需要確定充退磁線圈所采用的線圈組數(shù)及對應(yīng)的待優(yōu)化參數(shù)，通常待優(yōu)化參數(shù)隨著線圈匝數(shù)增多而增多。例如：4線圈系統(tǒng)待優(yōu)化參數(shù)有3個，即內(nèi)線圈距離d1、外線圈距離d2及內(nèi)外線圈電流比例ck12；6線圈系統(tǒng)的優(yōu)化參數(shù)增加為5個，即內(nèi)中外三對線圈的距離d1、d2、d3，中線圈與最內(nèi)線圈電流比ck12，最外線圈與最內(nèi)線圈電流比ck13。

2.2 定義磁場均勻度指標(biāo)目標(biāo)函數(shù)

一般采用磁場偏差來衡量磁場均勻度指標(biāo)大小，磁場偏差值越小，表明磁場均勻度越好。在考察空間內(nèi)取m個計算點，以線圈系統(tǒng)中心的磁場強度B0為參考，第i個考察空間點的磁場偏差定義為：

本文采用磁場平均相對誤差作為磁場均勻度指標(biāo)的目標(biāo)函數(shù)，平均相對磁場偏差定義為：

式（3）中，m為網(wǎng)格點總數(shù)。根據(jù)式（3），可以引入平均相對磁場均勻度指標(biāo)越大代表磁場均勻度越高。根據(jù)2.1知道，對于4線圈系統(tǒng)是關(guān)于參數(shù)(d1,d2,ck12)的函數(shù)，對于6線圈系統(tǒng)是關(guān)于參數(shù)(d1,d2,d3,ck12,ck13)的函數(shù)。實際工程對線圈距離有限制，體現(xiàn)為待優(yōu)化參數(shù)的邊界條件，一般要求d2>d1，d3>d2，d3

2.3 序列二次規(guī)劃算法迭代尋優(yōu)

線圈優(yōu)化的目標(biāo)函數(shù) minf是非線性函數(shù)，并且通常包含多個約束條件，本文采用的序列二次規(guī)劃算法是公認的求解約束非線性優(yōu)化問題最有效的方法之一。序列二次規(guī)劃把問題的求解主要分為兩個步驟：第一步為通過泰勒展開，將原非線性優(yōu)化問題轉(zhuǎn)換為簡單的二次優(yōu)化問題；第二步為對簡單的二次優(yōu)化（QP）問題進行迭代搜索求解。在序列二次規(guī)劃（SQP）的迭代過程中，每一步都需要求解一個或多個二次規(guī)劃（QP）子問題。序列二次規(guī)劃（SQP）算法流程見圖2，具體算法參考文獻[16-18]。

圖2 序列二次規(guī)劃（SQP）優(yōu)化算法流程圖Fig.2 Sequential Quadratic Programming (SQP) Optimization Algorithm Flowchart

2.4 線圈匝數(shù)整數(shù)化及二次尋優(yōu)

線圈匝數(shù)需要整數(shù)化處理。步驟2.3對于磁場均勻度目標(biāo)函數(shù) minf計算中以內(nèi)線圈電流I0=1A為基準(zhǔn)，利用2.3中優(yōu)化參數(shù)計算出線圈系統(tǒng)中心磁場強度B0。如果指標(biāo)要求的中心點最大磁場強度為Bmax，勵磁電源系統(tǒng)允許的最大電流為Imax，那么Bmax/B0就是線圈系統(tǒng)的安匝數(shù)，N1=Bmax/（B0*Imax）就是系統(tǒng)最內(nèi)線圈組單個線圈的匝數(shù)，N2=N1/ck12=Bmax/（B0*Imax*ck12）是第2個線圈組單個線圈的匝數(shù)，N3=N1/ck13=Bmax/（B0*Imax*ck13）是3個線圈組單個線圈的匝數(shù)。這里用迭代尋優(yōu)線圈電流比例系數(shù)計算得到的線圈匝數(shù)（N1、N2、N3）通常不是整數(shù)，而實際加工繞制都是以整數(shù)繞制的，所以需要對線圈匝數(shù)進行整數(shù)化處理。為保證線圈中心最大磁場強度，一般需要對小數(shù)向上取整。線圈匝數(shù)整數(shù)化之后，其它線圈參數(shù)就不是最優(yōu)了，因此需要固定線圈匝數(shù)整數(shù)化后的電流比例，對其它的可變參數(shù)以步驟2.3的結(jié)果為初始值，并在較小區(qū)域范圍內(nèi)進行二次尋優(yōu)計算，得到最終線圈系統(tǒng)最優(yōu)參數(shù)。

2.5 線圈系統(tǒng)均勻性指標(biāo)驗證

通過以上過程計算得到線圈系統(tǒng)的最優(yōu)參數(shù)后，需要進行數(shù)值仿真計算，檢驗在該參數(shù)條件下磁場的均勻性指標(biāo)是否滿足技術(shù)指標(biāo)要求。此外，需要綜合考慮其它工程約束條件是否滿足，如線圈系統(tǒng)電阻、電源最大功耗、電纜總重量等。如果不滿足，則需要調(diào)整線圈系統(tǒng)的約束條件，如線圈個數(shù)或者線圈組合形式等，重新進行優(yōu)化計算。在線圈加工過程中受到加工精度、加工公差及線圈繞制精度的影響，實際的磁場均勻度指標(biāo)可能會與理論仿真計算結(jié)果有一定偏差，所以對于目標(biāo)函數(shù)的優(yōu)化計算一般需要留有一定余量，例如，如果某區(qū)域的目標(biāo)函數(shù)是平均相對磁場均勻度 %100)1(min×-f不低于90%，也就是平均相對磁場偏差 minf<10%，在理論計算中可設(shè)置稍微嚴(yán)苛一些，如要求目標(biāo)函數(shù) minf必須不大于8%，為工程加工留有2%的余量。

3 線圈設(shè)計實例

結(jié)合工程應(yīng)用，給出一個共軸6線圈充退磁系統(tǒng)的設(shè)計實例。為了給某系列包絡(luò)尺寸不大于2m的衛(wèi)星充退磁，要求設(shè)計一個空間尺寸長寬高不大于3m的方形共軸赫姆霍茲充退磁線圈系統(tǒng)，最大線圈間距3m，中心待考察衛(wèi)星工作區(qū)域2m×2m×2m立方體均勻區(qū)的磁場均勻度不低于95%，中心點最大磁場強度不低于5mT，電源峰值電流不大于100A，充退磁勵磁電源最大功耗不大于40kW。

下面給出該設(shè)計實例的優(yōu)化計算過程及結(jié)果。根據(jù)普遍采用的方式，優(yōu)化前的線圈寬度a=3m，6個線圈等間距（線圈間距3/5m），電流相等。設(shè)計優(yōu)化計算中采用平均相對磁場偏差函數(shù)作為優(yōu)化目標(biāo)函數(shù)，以待考察中心區(qū)域的總體磁場偏差作為優(yōu)化目標(biāo)，而不是只考慮某一個點的磁場偏差。

表1給出了優(yōu)化前、初次優(yōu)化及線圈匝數(shù)取整二次優(yōu)化的參數(shù)結(jié)果。優(yōu)化前采用最常用的線圈參數(shù)，即6線圈等間距分布，匝數(shù)相等，d1=3/5m，d2=9/5m，d3=3m，ck12=1，ck13=1，中心待考察的2m×2m×2m立方體區(qū)域內(nèi)的平均相對磁場偏差值為0.073428，最低均勻度為88.22%。初次優(yōu)化后得到的最優(yōu)線圈參數(shù)，d1=0.4426m，d2=1.3610m，d3=3m，ck12=0.8827，ck13=0.2630，中心待考察的2m×2m×2m立方體區(qū)域內(nèi)的平均相對磁場偏差值為0.003259，最低均勻度為95.95%。優(yōu)化前與初次優(yōu)化對比，看出預(yù)定區(qū)域內(nèi)的平均相對磁場偏差降低到了優(yōu)化前的1/23，表明區(qū)域內(nèi)的整體磁場均勻度都得到了很大提高。同時，預(yù)定區(qū)域內(nèi)的最低磁場均勻度也由88.22%提高到95.95%。

表1 參數(shù)計算結(jié)果對比Table 1 Comparison of parameter calculation results

3.1 線圈匝數(shù)整數(shù)化處理

要求的中心磁場強度為5mT（即5000000nT），優(yōu)化前單位電流（1A）得到的線圈常數(shù)為1536.58nT，最大允許電流為100A，則單個線圈匝數(shù)5000000/1536.58/100=32.53，取整為33匝，6個線圈共198匝；初次優(yōu)化，以最內(nèi)線圈通電1A為基準(zhǔn)，計算得到的線圈常數(shù)是2579.58，所以計算得到充退磁線圈系統(tǒng)最內(nèi)、中間、最外三組線圈的匝數(shù)分別為N1=5000000/2579.58/100=19.38，N2=N1/ck12=19.38/0.8827=21.95，N3=N1/ ck13= 19.38/0.2630=73.69。對稱線圈匝數(shù)向上取整數(shù)之后，N1=20匝，N2=22匝，N3=74匝，線圈系統(tǒng)的總匝數(shù)為2×(N1+N2+N3)=232匝。

線圈匝數(shù)整數(shù)化之后，對于線圈距離參數(shù)進行二次尋優(yōu)，得到d1=0.4272m，d2=1.3753m，d3=3m，ck12=20/22=0.86956，ck13=20/74=0.27027，中心待考察的2m×2m×2m立方體區(qū)域內(nèi)的平均相對磁場偏差值為0.004858，最低均勻度為96.26%。線圈匝數(shù)取整并二次優(yōu)化之后與初次優(yōu)化相比較，最內(nèi)與中間線圈組的間距d1與d2有稍微調(diào)整，平均相對磁場偏差值從0.003259增加到0.004858，符合預(yù)期。

圖3是沿著線圈中心軸x=[-1m，1m]范圍內(nèi)，優(yōu)化前、初次優(yōu)化及線圈匝數(shù)取整二次優(yōu)化的磁場均勻度對比圖，虛線是優(yōu)化前線圈等間距、等電流（匝數(shù)）的計算結(jié)果，實線是初次優(yōu)化之后的結(jié)果，虛點線是匝數(shù)整數(shù)化之后二次優(yōu)化的結(jié)果。由圖2可知，優(yōu)化前線圈系統(tǒng)在中心軸的磁場均勻度分布極其不均勻，不能滿足技術(shù)指標(biāo)要求，只有在x=[-0.6m，0.6m]范圍內(nèi)大于95%。初次優(yōu)化后，在線圈中心軸x=[-1m，1m]區(qū)域范圍內(nèi)磁場均勻度整體提高到97.5%以上。在線圈匝數(shù)取整二次優(yōu)化之后，均勻度稍微降低，但整體均勻度變化不大，仍然在97.5%以上。

圖4中的子圖（a）（b）（c）分別是優(yōu)化前、初次優(yōu)化及線圈匝數(shù)取整二次優(yōu)化在中心區(qū)2m×2m×2m立方體表面的磁場均勻度分布圖。圖中的方框是6個方形線圈的位置，給出了中心區(qū)立方體6個表面的磁場均勻度等值面分布圖。由圖3（a）看出，待考察的立方體區(qū)域表面最低磁場均勻度88.22%，并且大部分表面磁場均勻度值低于90%，不滿足磁場均勻度95%的指標(biāo)要求。從圖4（b）、（c）看出，全部區(qū)域滿足磁場均勻度大于95%的指標(biāo)要求；磁場均勻度在立方體的8個角點附近相對略差，其值大約在96%附近，這也符合矩形線圈在角點附近的磁場均勻度偏差的普遍規(guī)律。在4個側(cè)面的大部分區(qū)域磁場均勻度在98.9%以上。所以，線圈系統(tǒng)優(yōu)化之后整體磁場均勻度分布有顯著提高，更能滿足充退磁線圈對于中心工作區(qū)域的均勻度指標(biāo)要求。

圖3 沿著線圈中心軸x=[-1m，1m]范圍內(nèi)，優(yōu)化前、初次優(yōu)化及線圈匝數(shù)取整二次優(yōu)化的磁場均勻度對比Fig.3 Along the coil central axis x=[-1m, 1m], the comparison of the magnetic field uniformity before optimization, the first optimization and the second optimization of the rounded number of coil turns

圖4 均勻區(qū)磁場均勻度分布圖，a）b）c）分別是優(yōu)化前、初次優(yōu)化及線圈匝數(shù)取整二次優(yōu)化的結(jié)果Fig.4: The distribution diagram of the uniformity of the magnetic field in the uniform area, a) b) c) are the results of the optimization before optimization, the first optimization and the second optimization of coil turns

3.2 線圈電參數(shù)計算

線圈尺寸優(yōu)化設(shè)計完成后，需要根據(jù)充退磁電源允許的峰值電流和最大功耗等參數(shù)，選取電纜的尺寸，估算線圈電纜總電阻、總重量、電源功耗等。本設(shè)計方案允許的最大電流為100A，根據(jù)經(jīng)驗選取直徑4mm的漆包銅線作為繞線電纜。銅線電阻率ρ通常取值1.75×10-8Ω/m，銅線電阻計算公式如下：

式（4）中a為方形線圈寬度，Nc為線圈系統(tǒng)總匝數(shù)，d為銅線直徑取值4×10-3m。

分別對優(yōu)化前及匝數(shù)取整二次優(yōu)化的線圈物理常量進行計算，表2給出線圈物理常量的計算結(jié)果對比。得到優(yōu)化前（Nc=198匝）的充退磁線圈參數(shù)：長度2376.00m，電阻3.31Ω，功耗33.09kW，重量265.73kg；線圈匝數(shù)取整二次優(yōu)化（Nc=232匝）的充退磁線圈參數(shù)：長度2784.00m，電阻3.87Ω，功耗38.77kW，重量311.36kg。雖然優(yōu)化設(shè)計后的線圈系統(tǒng)匝數(shù)、功耗與重量增加約17.2%，但是在中心待考察衛(wèi)星工作區(qū)（2m×2m×2m立方體區(qū)域）的平均相對磁場偏差降低到了優(yōu)化前的1/23，整體磁場均勻度得到了極大提升。

表2 線圈物理常量計算結(jié)果對比Table2 Comparison of calculation results of coil physical constants

4 結(jié)論

本文基于序列二次規(guī)劃（SQP）方法，對共軸多組充退磁線圈進行系統(tǒng)優(yōu)化設(shè)計。用離散積分的方法計算消磁線圈在空間點的磁場，以待考察區(qū)域的平均相對磁場偏差作為優(yōu)化目標(biāo)函數(shù)，并考慮工程約束條件，優(yōu)化計算出不同線圈位置和驅(qū)動電流等參數(shù)。對線圈匝數(shù)整數(shù)化后進行二次優(yōu)化，提高了均勻區(qū)的指標(biāo)。結(jié)合一個實際工程中的共軸6線圈系統(tǒng)設(shè)計案例，給出了詳細的優(yōu)化步驟及結(jié)果，對于線圈匝數(shù)取整及二次尋優(yōu)有詳細闡述。仿真計算結(jié)果表明，經(jīng)過優(yōu)化設(shè)計后的共軸充退磁線圈系統(tǒng)的磁場均勻度得到了極大提升，在中心待考察衛(wèi)星工作區(qū)的平均相對磁場偏差降低到了優(yōu)化前的1/23，其它電參數(shù)均能滿足技術(shù)指標(biāo)要求。