任意層導電板上方線圈阻抗的級數計算模型*

張思全，胡盛斌，陸文華

(上海工程技術大學航空運輸學院，上海 201620)

在石油、化工、航空等領域，許多關鍵設備和部件由多層導電材料制成，對其內部缺陷和結構完整性的檢測是一個較困難的問題。超聲檢測法雖具有可檢測較厚材料的優勢，但不適于檢測相鄰層間含有空氣的多層導電結構。電渦流檢測技術利用材料中缺陷與渦電流、線圈激勵電磁場之間相互作用引起線圈阻抗變化的原理來檢測缺陷并對結構完整性進行評價，并且具有不需耦合劑、可檢測多層導電結構隱藏缺陷的優勢［1－4］。

如何快速、精確求解線圈阻抗是渦流檢測的重要問題。數十年來，許多學者對此進行了研究。Dodd和Deeds研究了正弦電流激勵線圈位于有均勻厚度導電材料覆層的半無限大導電平板上方的情況，得到了積分形式的矢量磁位閉合解析解［5］;其后，Luquire等人又研究了線圈置于任意層導電結構中間時的阻抗變化模型［6］。在國內，雷銀照研究了半無限大導體上方覆蓋一定厚度金屬涂層構成的兩層平板導電結構的軸對稱渦流解析解［7］;黃平捷研究了置于多層導電結構上方線圈響應與多層導電結構厚度變化之間的關系模型［8］;范孟豹推導了位于任意層導電結構上方渦流探頭的阻抗解析模型，在阻抗計算中引入符號運算法求解矢量磁位表達式系數，減小了程序計算量［9］。但這些研究采用徑向無限大導體，所得解為積分形式，存在計算困難、時間花費大的問題。

本文在前人研究基礎上，采用解域截取法分析矩形截面圓柱線圈置于任意層導電平板上方的正向問題，推導渦流檢測探頭線圈與任意層導電結構相互作用的級數形式阻抗變化計算模型。為驗證模型的正確性，針對單層導電結構的特例，計算了線圈置于其上方時的阻抗變化，分析了影響計算精度的主要因素，并與有限元計算結果進行了對比。

1 任意層導電板上方線圈阻抗計算

1.1 任意層導電板上方單匝線圈阻抗計算

1.1.1 任意層導電板上方單匝線圈的矢量磁位

如圖1a)所示，一個單匝通電圓形線圈置于任意層(設為n層)導電平板上方，其軸線與任意層平板表面垂直;設各層平板的厚度、電導率、介電常數都是任意的，各層介質都是線性、各向同性、均勻的。設線圈中為正弦激勵電流。采用圓柱坐標系O－(r，φ，z)，各單位矢量分別為er，eφ，ez。設線圈軸與Z軸重合，與任意導電平板表面的交點為原點O，線圈的半徑為r0，線圈在平板上方的提離距離為z0。任意層導電結構最下方第1層，記為L1，設其相對磁導率為μr1，電導率為σ1，在實際運算中可假設其厚度為無限大或為0。自下往上，各層結構編號、電磁參數與厚度表達式依次增加，即第i層記為Li，相對磁導率為μri，電導率為σi。最上方第n層導電平板記為Ln。通電圓環線圈下方空氣層(0＜z＜z0)記為Ln+1，線圈上方(z＞z0)空氣層記為Ln+2。

圖1 線圈位于徑向被圓柱邊界截取的多層導電結構上方

許多文獻在研究多層導電結構時，研究區域為無窮大，電磁場以傅立葉－貝塞爾積分項的形式表示，在數值計算時收斂很慢。將積分變為級數形式，并通過調整級數的求和項數來控制計算的精度和速度可達到提高計算速度的目的。為此，在距離線圈適當距離處施加一個圓柱形邊界并令邊界處的解為零，就得到一個被截取的圓柱形解域，這樣雖然會增加解的誤差，但只要包含足夠多項，用級數解代替積分解就不會影響計算精度［10－11］。

1.1.2 矢量磁位的邊值問題

麥克斯韋方程組是渦流檢測的理論基礎。為了簡化電磁場的計算及分析，引入矢量磁位A。

圖1中所有區域被與線圈同軸、半徑為b的圓柱面截取，從下至上，共有n+2個區域，有n+1個相鄰層之間界面zk，其中k=1，2，…，n，n+1。zk也表示各界面的z軸坐標。

半徑為r0的單線圈位于內邊界面z=z0平面內，其電流密度可以用δ函數寫為

式(1)中I為單匝線圈激勵電流的幅值。由于介質和場源關于對稱軸r=0呈軸對稱分布，所以圖1中場的矢量磁位僅含有周向分量Aφ。

忽略位移電流，在一個各向同性、線性均勻媒質中，單匝線圈在(r0，z0)處產生的矢量磁位A(r，z)滿足下列約束方程及其邊界條件［5］:

(1)約束方程

式(2)中i=1，2，3，…，n+2，表示相應的層數;=－jωμriμ0σi;當I≠0 時，各區域矢量磁位非零且有界。

(2)邊界條件

在所有界面zi，i=1，2，…，n，n+1，滿足:

在界面zi，i=1，2，…，n，滿足:

圓形單線圈的上下區域，即第Ln+1區域和Ln+2區域之間界面zn+1滿足:

1.1.3 矢量磁位約束方程的通解

在圓柱坐標系下，約束方程式(2)可寫為:

用分離變量法求解式(6)得:

式(7)中J1(λjr)和Y1()分別是第一類和第二類一階貝塞爾函數。λj稱為分離常數，它與r和z無關，由于在r=b處施加了截取邊界，所以引入離散本征值λj來代替Dodd和Deeds模型中的連續值，由J1(λjb)=0可求得本征值λj。式(7)中

第L1層

第Li層，i=2，…，n+1

第Ln+2層

根據邊界條件式(3)，在界面z=zi，i=1，2，…，n，n+1，滿足

根據邊界條件式(4)，在界面z=zi，i=1，2，…，n，滿足

根據式(5)、式(10)和式(11)，并利用貝塞爾函數的正交性

線圈上下區域交界面z=zn+1=z0，滿足

這樣，共得到2n+2個方程，需要求解2n+2個未知系數。

1.1.4 矢量磁位的求解

在導電結構層數較多時，各層矢量磁位的系數較復雜，將式(13)重寫為

式(16)中定義

利用式(12)、式(16)，對界面z=zi，i=1，2，…，n，因為C12=0，可以求解得到以下關系

在界面z=zn+1，因為C(n+2)1=0，并且第Ln+1、Ln+2層為空氣，αn+1=αn+2=λj，由式(12)可得

由式(15)、式(20)可得

式(18)、式(19)可寫成矩陣形式

令Ti為式(23)中的轉換矩陣，則各系數之間關系可寫為

令

因此

將式(21)代入式(29)得

由式(22)、式(30)可得

由式(21)、式(28)得

這樣可以得到任意區域矢量磁位的系數，進而可以寫出任意區域的矢量磁位。其中解析計算中要利用線圈上、下方區域的矢量磁位。

線圈上方區域的矢量磁位

線圈下方區域的矢量磁位

1.2 任意層導電平板上方多匝線圈阻抗計算

確定了由單匝線圈在Li(i=1，2，…，n+1，n+2)區域的矢量磁位Ai(r，z)后，對于圖1b)所示具有矩形截面的多匝線圈，可以通過將單匝線圈導出的矢量磁位沿線圈橫截面積分獲得各區域總的矢量磁位:

在獲得矩形截面多匝線圈上方Ln+2區域總的矢量磁位和下方Ln+1區域總的矢量磁位后，令中的l2為z，令中的l1為z，并將兩式相加，可以得到線圈區域任一點的矢量磁位Atotal。利用公式V(r，z)=jω2πrAtotal(r，z)可以獲得單線圈上的感應電動勢，通過疊加可以得到m匝線圈的自感電壓為

從而得出任意層導電結構上方m匝線圈阻抗為:

式(37)中，I為單匝線圈內的電流幅值。

2 單層導電結構上方線圈阻抗計算

為了驗證所導出任意層導電結構上方m匝線圈阻抗表達式(37)的正確性，下面計算單層導電結構上方線圈阻抗的特例。

2.1 單層導電結構上方單匝線圈的矢量磁位

則置于單層導電結構上方線圈的阻抗變化為

令式(23)中的轉換矩陣Ti中的i取1，可得單層導電結構中轉換矩陣T1;令式(33)、式(34)中的n=1，可分別得到單層導電結構線圈上部區域L3和下方區域L2的矢量磁位A3、A2

圖2 多線圈位于單層導電結構上方

2.2 單層導電結構上方多匝線圈的阻抗計算

確定了單匝線圈上下區域的矢量磁位A3、A2后，對于圖2所示具有矩形截面的圓柱線圈，可以利用式(35)通過將單匝線圈的矢量磁位沿線圈橫截面積分獲得總的矢量磁位。

計算中采用線圈電流密度

可分別得到m匝線圈上方區域矢量磁位和m匝線圈下方區域矢量磁位。令中的l2=z中的l1=z，并將兩式相加，可得m匝線圈區域的矢量磁位A(2－3)為

利用式(36)，通過疊加可以得到探頭線圈的自感電壓為

從而得出單層導電結構上方圓柱形線圈的阻抗表達式為:

式(42)中

m為矩形線圈的匝數，I為單匝線圈內的電流幅值。m匝線圈的電流密度如式(39)所示。

式(38)中的單層導體上方線圈阻抗Zw即為式(42)。令Zw中α1=λj、μr1=1即得到無導體時的線圈阻抗Zh

利用式(38)即可計算出單層導電結構中渦電流引起的線圈阻抗變化為

2.3 級數展開計算及有限元仿真驗證

為了檢驗級數計算法的正確性，采用Mathematica軟件，根據式(45)，在不同線圈激勵頻率下，對線圈阻抗變化情況進行了計算。計算中所采用線圈及導電結構參數如表1所示。

表1 線圈及導電結構相關參數

為驗證級數計算法的正確性，同時采用有限元法對置于單層導電板上方線圈的阻抗變化進行仿真計算。一般文獻通常采用ANSYS軟件［12］，這里采用Maxwell 2D電磁場計算軟件進行仿真。首先構造如圖3所示的2D模型，由于模型關于線圈軸對稱，所以只需畫出一半;然后指定模型的材料屬性，并設定邊界條件和激勵源;在模型左邊界施加奇對稱邊界條件，在線圈上施加電流激勵源;設定求解參數，將線圈上的阻抗與電感值包含在計算過程中。

如圖3所示進行網格劃分，并設定求解器殘差，最后啟動求解過程。

圖3 2D有限元模型及網格劃分

在不同激勵頻率下，將級數展開法與有限元仿真計算所得的阻抗變化結果進行了對比，其中電阻和電抗的對比分別如圖4、圖5所示。在級數法計算中，取解域截取半徑b=30r2，計算級數前150項的和。

圖4 線圈電阻變化與激勵頻率關系

圖5 線圈電抗變化與激勵頻率關系

由圖4、圖5可見，級數法計算所得線圈阻抗變化值與有限元仿真計算結果能夠很好吻合，驗證了級數展開計算模型的正確性。

3 討論

為進一步研究解域截取半徑與求和項數的不同選擇對級數法計算線圈阻抗變化值精度的影響，在激勵頻率為1 kHz、10 kHz、50 kHz和100 kHz下，分別取解域截取半徑b=10r2和b=30r2，分別計算級數前50項和前150項的和，計算中所采用線圈及導電結構參數如表1所示。計算各參數下級數法與仿真計算法所得線圈阻抗變化值之間的相對誤差，結果如表2所示。

表2 級數計算與有限元仿真結果之間的相對誤差單位:%

由表2中可見，在較大的截取半徑b=30r2情況下，如果級數相加項數為較少的30項，則所得計算結果誤差大于較小截取半徑b=10r2的情況;但如果將級數求和項數增大到150項時，計算精度就大于b=10r2的情況。這說明一個較大的截取半徑需要更多的求和項來保證較好的計算精度。

級數計算式(42)、式(45)中只有截取半徑b和求和項數需要確定，積分項Int(x1，x2)僅依賴于線圈內外半徑r1和r2，而與線圈的提離和材料特性無關，可以提前計算。因此級數計算法可獲得較快的計算速度。

4 結論

從Maxwell方程組出發，引入矢量磁位，采用解域截取和變量分離的渦流檢測線圈阻抗計算方法，將有限層導電平板推廣到任意層，導出了位于任意層有限區域導電平板上方線圈阻抗變化的模型，用易于計算的級數解代替積分解，拓寬了經典線圈阻抗計算公式的應用范圍，提高了計算速度，這在需要計算大量阻抗值的情況下非常有用。對單層導電結構的特例進行了計算并用有限元仿真方法進行了驗證，證明了級數展開計算方法的正確性，可以應用在多層導電結構缺陷渦流檢測法定量檢測正向模型中。

［1］Theodoulidis T，Poulakis N，Dragogias A.Rapid Computation of Eddy Current Signals from Narrow Cracks［J］.NDT＆E International，2010，43:13－9.

［2］Guohou Li，Pingjie Huang，Peihua Chen，et al.Quantitative Nondestructive Estimation of Deep Defects in Conductive Structures［J］.International Journal of Applied Electromagnetics and Mechanics，2010，33(3－4):1273－1278.

［3］Sun Y S，Ouyang T，Udpa S.Remote Field Eddy Current Testing:One of the Potential Solution for Detecting Deeply Embedded Discontinuities in Thick and Multilayer Metallic Structures［J］.Materials Evaluation，2001，59(5):632－637.

［4］張玉華，孫慧賢，羅飛路，等.一種用于盤孔裂紋檢測的差動式渦流探頭的設計與實現［J］.傳感技術學報，2008，21(6):1079－1083.

［5］Dodd C V，Deeds W E.Analytical Solutions to Eddy Current Probe-Coil Problems［J］.Journal of applied physics，1968，39(6):2829－2838.

［6］Luquire J W，Deeds W E，Dodd C V.Alternating Current Distribution Between Planar Conductors［J］.Journal of Applied Physics，1970，41(10):3983－3991.

［7］雷銀照，馬信山.渦流線圈的阻抗計算［J］.電工技術學報，1996，11(1):17－20.

［8］黃平捷，吳昭同，鄭建才，等.多層厚度電渦流檢測反演算法及實驗研究［J］.儀器儀表學報，2005，26(4):428－432.

［9］范孟豹，黃平捷，葉波，等.多層導電結構厚度電渦流檢測解析模型及實驗驗證［J］.中國電機工程學報，2008，28(27):142－147.

［10］Theodoulidis T，Kriezis E.Series Expansions in Eddy Current Nondestructive Evaluation Models［J］.Journal of Materials Processing Technology，2005，161(5):343－347.

［11］于亞婷，杜平安，李代生.電渦流傳感器線圈阻抗計算方法［J］.機械工程學報，2007，43(2):210－214.

［12］任吉林，刁海波，唐繼紅，等.渦流傳感器提離效應的ANSYS模擬［J］.傳感技術學報，2008，21(6):967－971.