表面波檢測的數(shù)值模擬與試驗對比

王 貝，竇曉蘭，韓曉華

（1.山東省科學(xué)院激光研究所，濟寧 272017；2.山東省濟寧市糧食局，濟寧 272019）

為獲得工件的超聲波檢測圖像，按照檢測的步驟，要先制定檢測工藝方案，再制備對比試塊，準備儀器和換能器，這一過程不可避免地耗費著大量的人力物力。另外，對于結(jié)構(gòu)復(fù)雜工件的檢測結(jié)果，還需要有經(jīng)驗的技術(shù)人員做出判定。而借助于有限元數(shù)值模擬可以彌補上述的不足，其能夠進行縱波、橫波、蘭姆波等超聲波聲場的仿真，可直觀形象地呈現(xiàn)超聲波在工件中的傳播過程，并根據(jù)需要進行再現(xiàn)。在聲場模擬計算方面，有學(xué)者做過一定的工作：宗侶等利用ANSYS軟件模擬了導(dǎo)波在管材中的傳播［1］；剛鐵等利用MARC 模擬計算了Lamb波在薄板中的傳播［2］；孫繼華 等［3］及劉長福等［4］進行了縱波及橫波的有限元模擬計算方面的工作。然而對于表面波檢測，尤其是遇到缺陷后所發(fā)生的波型轉(zhuǎn)換現(xiàn)象，相關(guān)的有限元分析工作較少。筆者以無損檢測中經(jīng)常使用的1.5 MHz表面波探頭為例，運用斜楔縱波入射在界面處獲得表面波的原理，對表面波的傳播過程進行了有限元分析。

1 理論模型

表面波在介質(zhì)表面?zhèn)鞑r，介質(zhì)表面質(zhì)點作橢圓運動，可等效為質(zhì)點縱向與橫向振動的合成?？v波折射法表面波換能器的原理如下：壓電晶片產(chǎn)生的縱波通過不同傾角的斜楔射向試件表面，經(jīng)波型轉(zhuǎn)換可在試件中相應(yīng)地產(chǎn)生表面波，斜楔一般用有機玻璃制作［5］。圖1所示為表面波探頭壓電晶片和斜楔的二維結(jié)構(gòu)簡圖，介質(zhì)1為有機玻璃，介質(zhì)2為被檢測工件。有機玻璃／工件界面上波的入射范圍為AB，縱波的入射角為αL，當(dāng)入射角大于第二臨界角時，工件中無縱波和橫波出現(xiàn)，經(jīng)波型轉(zhuǎn)換后僅在界面上形成表面波，斜楔中縱波的速度為CL，工件中表面波的波速為CR。

圖1 工件表面等效為表面波換能器的原理

設(shè)Di、Ei、Dj、Ej為介質(zhì)1 中任意兩條波線 上的點，壓電晶片產(chǎn)生的同相位縱波分別沿波線DiEi與波線DjEj傳播到界面上Fi與Fj，相對應(yīng)的時間分別為Ti與Tj，關(guān)系為：

式（2）說明同相位的縱波由壓電晶片發(fā)出，沿不同波線到達BC 上各點具有一定的時間差。于是，在進行數(shù)值模擬時，為了簡化模型，可以去除斜楔和介質(zhì)1中的壓電晶片，只考慮介質(zhì)2工件，將工件的AB區(qū)域等效看為表面波換能器壓電晶體激勵表面，那么界面上的各振點之間應(yīng)該存在著對應(yīng)于時間差ΔTi，j的相位差。

AB區(qū)域中各振點的激勵力方向也是一個很重要的因素，實際超聲波無損檢測過程中，由于介質(zhì)1中表面波探頭和介質(zhì)2之間存在著液體耦合劑，因為耦合劑屬于液態(tài)物質(zhì)，不能傳遞切向力，所以激勵力的方向應(yīng)該平行于介質(zhì)2表面法線的方向［6］。因此在有限元分析過程中，激勵信號施加的方向應(yīng)該與工件表面垂直，另外為了避免邊界反射造成影響，將兩側(cè)邊界設(shè)置為吸收邊界［7］。

2 數(shù)值計算與分析

2.1 材料屬性參數(shù)

假設(shè)選取的工件材料內(nèi)部的力學(xué)性能是各向同性的，工件材料為鋼。另外計算時有機玻璃中縱波的速度為2 700m／s。同時鋼中縱波的速度為5 940m／s，橫波速度為3 290m／s，表面波速度為3 041m／s。工件的力學(xué)性能參數(shù)分別為：彈性模量為2.16×1011Pa；泊松比為0.28；密度為7 800kg·m－3。

2.2 有限元模型

工件的整體尺寸為50mm×20mm（長×高），加載區(qū)域的中心處于上表面距左側(cè)10mm 的位置，加載區(qū)域的寬度為10 mm。缺陷中心距工件左側(cè)25mm，缺陷尺寸為1 mm×5mm（寬×深）。有兩種建立缺陷工件模型的方法：一種是直接創(chuàng)建帶有缺陷的工件，然后進行網(wǎng)格劃分；另一種是先建立完整的工件，進行網(wǎng)格劃分后，利用單元生死技術(shù)將缺陷處的單元殺死。為了便于線性計算和控制計算時間步長，選用第一種方法生成工件模型。根據(jù)超聲換能器產(chǎn)生激勵信號的形式，在加載區(qū)域內(nèi)施加的力［8］為式（3），力的單位為N，時域函數(shù)圖形與頻譜圖如圖2所示。

式中：f為激勵信號 的中心 頻率1.5 MHz，；n為周期數(shù)。取n＝3。

圖2 激勵信號的時域波形與頻譜圖

選取4個節(jié)點的2維結(jié)構(gòu)單元劃分工件，每個節(jié)點各有兩個自由度，在缺陷周圍的區(qū)域，將4節(jié)點結(jié)構(gòu)單元退化為3節(jié)點結(jié)構(gòu)單元，可以將模型劃分的更加有規(guī)律些。為滿足描述工件中彈性表面波傳播的精度要求，網(wǎng)格的大小按式（4）進行選?。?－10］。

式中：λmin為計算時所關(guān)注的最小波長。

根據(jù)表面波速度及中心頻率，可計算出：

則Δl可取為2×10－4m。

吸收邊界采用粘彈性人工邊界的形式，可等效為連續(xù)分布的并聯(lián)彈簧－阻尼器系統(tǒng)，人工邊界上法向和切向的單元同時具有彈簧剛度和阻尼系數(shù)，在此選取2維彈簧—阻尼器單元。根據(jù)以上數(shù)據(jù)建立的有限元模型如圖3所示。

圖3 有限元模型

2.3 載荷加載激勵與求解

為了確定系統(tǒng)在瞬態(tài)載荷作用下位移、應(yīng)力、應(yīng)變、力隨時間變化的規(guī)律，即系統(tǒng)的結(jié)構(gòu)動力學(xué)響應(yīng)，應(yīng)當(dāng)選用瞬態(tài)動力學(xué)分析（時間歷程分析）。結(jié)構(gòu)的動力響應(yīng)可以看作是各階模態(tài)響應(yīng)的組合，積分步長Δt應(yīng)小到能夠解出對結(jié)構(gòu)整體響應(yīng)有貢獻的最高階模態(tài)，設(shè)fmax為結(jié)構(gòu)響應(yīng)的最高階頻率（Hz），則積分步長應(yīng)為：

為了產(chǎn)生表面波，對于式（3）所描述的激勵信號，加載區(qū)域內(nèi)各節(jié)點的起振時刻需要有一定的時間差，在程序中首先根據(jù)式（2）計算好每個節(jié)點各自的延遲時間，然后按照時序與加載節(jié)點的先后順序?qū)虞d區(qū)域中每個節(jié)點進行力的加載。節(jié)點力加載完畢后，進行求解和數(shù)據(jù)處理。

2.4 結(jié)果分析與試驗

為了便于描述，將帶有缺陷的工件標(biāo)記為工件A，對A 進行有限元分析。首先說明表面波未到達缺陷前的速度，對于工件A，在上表面處距離左側(cè)16，21，36，44 mm 處分別取四個不同的節(jié)點AA、BA、CA、DA，下標(biāo)表示節(jié)點屬于工件A，顯然AA、BA處于缺陷左側(cè)而CA、DA處于缺陷右側(cè)。提取AA、BA兩個節(jié)點位移的時間歷程數(shù)據(jù)，繪制曲線，如圖4（a）所示，實線矩形框中為直通表面波信號，虛線為反射波信號。根據(jù)直通波兩波峰之間的時間差Δt1＝1.62μs以及兩節(jié)點間距Δs＝5mm，計算出表面波的速度為3 082.8m／s，與理論速度3 041m／s之間的相對誤差為1.37%，可見有限元方法計算出來的表面波速度與理論速度相符。同理，求出反射波的速度同樣為3 082.8m／s，該反射波仍為表面波，可以說明：一部分表面波信號在裂紋開口處仍以表面波的形式被反射，并沿物體表面的原路徑返回。

其次，說明表面波遇到缺陷時，會沿裂紋表面?zhèn)鞑?，并在?jīng)過裂紋后，仍能夠以表面波的形式沿工件表面繼續(xù)傳播。圖4（b）為CA、DA兩點位移的時間歷程曲線，按照上述求解速度的方法，計算出波的傳播速度為3 098.4m／s，與表面波理論速度之間的相對誤差為1.89%，說明該信號仍為表面波，因此所述的觀點能夠得以驗證。另外可利用CA、DA兩節(jié)點的起振時刻與無缺陷時的起振時刻之間的差值，估算出缺陷尺寸。

圖4 工件A 上表面節(jié)點位移的時間歷程曲線

圖5中為工件A 不同時刻的波場快照圖片，圖5（a）為5.28μs時刻的波場快照，此時R1表面波未到達缺陷；圖5（b）為10.4μs時刻的波場快照，顯示了表面波遇到缺陷后聲場發(fā)生改變的狀況。根據(jù)上述聲速計算可知，R2為遇到缺陷反射回來的表面波，R3為經(jīng)過缺陷而繼續(xù)傳播的表面波。對于信號波S1，S2，S3，L1，需要計算出波的傳播速度以確定波的類型。具體實現(xiàn)過程是，在不同時刻，分別取出相應(yīng)信號波傳播路徑中位移最大值對應(yīng)的節(jié)點，這些節(jié)點即是處于聲束軸線上的點，根據(jù)節(jié)點坐標(biāo)便可以求出節(jié)點間距，結(jié)合時間差，求出多個速度值然后取平均。將所獲取的信號波速度列于表1，其中S1，S2，S3的速度均與橫波理論速度相吻合，L1的速度與縱波速度相吻合。說明表面波在遇到缺陷時發(fā)生了波型轉(zhuǎn)換，轉(zhuǎn)換為變形縱波和變形橫波，在工件內(nèi)部傳播。

圖5 工件A 的聲場快照

表1 信號波速度的仿真求解值與理論值

通過試驗的方式對上述有限元分析結(jié)果進行說明，為了便于固定探頭和分析試驗結(jié)果，采用290mm×290mm×60mm 的45鋼工件。表面波發(fā)生波型轉(zhuǎn)換后的各類型波具有一定的指向性，即聲束軸線的角度一定，實際工件和有限元模型尺寸上有差別，厚度尺寸比例系數(shù)為3，根據(jù)比例系數(shù)求出表面波探頭與缺陷的距離為45 mm。采用漢威HS616e型便攜式探傷儀和1.5 MHz的表面波探頭激發(fā)表面波。使用型號TDS 2024C 的泰克示波器接收壓電探頭輸出的信號，主要參數(shù)為帶寬200 MHz，采樣頻率2GS／s，為了去除噪聲干擾，對采集到的信號進行8次平均。另外為了接受縱波和橫波，還使用了中心頻率均為1.5MHz的縱波探頭和K1.5橫波探頭，圖6為試驗方案，圖中矩形為表面波激勵探頭，圓形A、B、C、D 四點表示為信號接收探頭的位置，PR 為表面波探頭，PL 為縱波探頭，PS為橫波探頭。根據(jù)探頭距缺陷的距離和表面波速度，計算出表面波到達缺陷所用時間為TR＝14.8μs。

圖7（a）所示表面波到達D 處的時刻為TD＝29.6μs，與TR的差值為14.8μs，乘以表面波的速度，得到相應(yīng)波程為45mm，又因為D 點距表面裂紋42mm，因此可估算出工件表面裂紋的縱向尺寸為1.5mm。圖7（b）為B、C 處橫波探頭接收的信號，讀取試驗數(shù)據(jù)知TB＝36.6μs，TC＝50.6μs，與TR的差值即為橫波傳播所需的時間，結(jié)合B、C兩點缺陷的距離，可求出波速分別為3 303 m／s、3 296m／s，均與橫波速度吻合。圖7（c）為A 處探頭接收到的信號，從圖中得出TA＝38.2μs，按照上述方法計算出從缺陷處傳播至A 處的波速度，其數(shù)值為5 983m／s，與縱波速度值相符。

圖6 試驗方案

圖7 各探頭接收的信號

通過上述對比試驗和分析可知，根據(jù)得出的速度數(shù)據(jù)，確定信號波的類型，進而推斷出表面波遇到缺陷后整個聲場的狀況：聲場中除了含有沿著缺陷并繼續(xù)傳播的表面波外，還會生成縱波和橫波，這與有限元分析結(jié)果相一致。

3 結(jié)論

提出了等效表面波換能器有限元模型，實現(xiàn)了工件表面波檢測的有限元分析，并利用試驗進行了分析。數(shù)值計算獲得的表面波、橫波及縱波的速度均與理論值相符合。表面波遇到缺陷時一部分仍繼續(xù)傳播，另外發(fā)生的反射、波型轉(zhuǎn)換等，均能夠通過有限元模擬計算出來。利用有限元分析可以對工件檢測進行結(jié)果預(yù)測和評估，能夠為超聲波檢測提供參考和對比。并且相對于理論分析，有限元分析更加直觀方便。

［1］宗侶，王悅民，朱龍翔.超聲導(dǎo)波在彎管中傳播特性的有限元分析［J］.無損檢測，2013，35（11）：43－46.

［2］剛鐵，劉強.薄板結(jié)構(gòu)中Lamb波的檢測與仿真［J］.無損檢測，2013，35（7）：24－28.

［3］孫繼華，李書光，倪云鹿，等.固體中缺陷對聲波頻譜特性影響的有限元仿真［J］.計算機仿真，2008，25（8）：308－311.

［4］劉長福，牛曉光，李中偉，等.基于ANSYS 的超聲縱／橫波傳播仿真計算［J］.無損檢測，2011，33（6）：5－18.

［5］鄭輝，林樹青.第2版：超聲檢測［M］.北京：中國勞動社會保障出版社，2008.

［6］JASIUNIENE E，ZUKAUSKAS E.The ultrasonic wave interaction with porosity defects in welded rail head［J］.Ultragarsas，2010，65（1）：12－18.

［7］李書光，孫繼華，李樹榜，等.孔洞對透射波頻譜影響的有限元數(shù)值模擬［J］.無損檢測，2009，31（1）：72－75.

［8］BASKARAN G，RAO C L，BALASUBRAMANIAM K.Simulation of the TOFD technique using the finite element method［J］.Insight－Non－Destructive Testing and Condition Monitoring，2007，49（11）：641－646.

［9］NANDY A，MULLICK S，De S，et al.Numerical simulation of ultrasonic wave propagation in flawed domain［C］／／［s.l］：NDE，2009：160－163.

［10］GHOSE B，BALASUBRAMANIAM K，KRISHNAMURTHY C V，et al.Two dimensional FEM simulation of ultrasonic wave propagation in isotropic Solid Media using COMSOL［C］／／India：COMSOL Conference，2010.