石穎穎，趙金玲，楊樂輝,趙建平

(1.南京工業(yè)大學(xué) 機械與動力工程學(xué)院，南京 211816；2.江蘇省極端承壓裝備設(shè)計與制造重點實驗室，南京 211816)

0 引言

復(fù)合材料管道結(jié)構(gòu)廣泛應(yīng)用于軍工領(lǐng)域(航空、航天、導(dǎo)彈等)和民用領(lǐng)域(石油、化工、建筑、環(huán)保、貯罐等)。復(fù)合材料壓力容器在傳熱、傳質(zhì)、反應(yīng)等工藝過程中，不可避免地會發(fā)生損傷，其中分層損傷是層合復(fù)合材料中最常見的損傷之一[1-2]。在循環(huán)載荷作用下，分層會沿著一定的方向擴展，最終導(dǎo)致材料失效。超聲檢測可有效地對結(jié)構(gòu)的安全現(xiàn)狀進行評估，為及時進行損傷修復(fù)以及壽命預(yù)估提供依據(jù)[3-7]。

超聲檢測的開展依賴于結(jié)構(gòu)的頻散特性研究。20世紀20年代，GHOSH[8]首次根據(jù)波在空心圓柱殼中傳播的理論基礎(chǔ)，推算出其對應(yīng)的線彈性解，但沒有對應(yīng)的數(shù)值分析和試驗驗證。LIN等[9]以板殼理論為基礎(chǔ)，得到了空心圓柱殼中軸對稱模態(tài)的頻散關(guān)系。對于多層管道結(jié)構(gòu)，BARSHINGER[10]使用全局矩陣法，推導(dǎo)出導(dǎo)波在多層圓柱殼中的頻散方程。半解析有限元方法(Semi-analytical Finite Element Method,SAFE)是分析多層管道頻散特性的另一重要方法。NELSON等[11-12]利用半解析有限元法，求解彈性層狀正交各項異性板以及圓管結(jié)構(gòu)的彈性波傳播問題。HAYASHI等[13-15]基于半解析有限元法，開展了多種波導(dǎo)結(jié)構(gòu)中導(dǎo)波的傳播特性研究，包括板、圓管以及鐵軌等。最近，HAKODA等[16]提出基于Floquet邊界條件的特征頻率方法，求解了單層板、單層管的頻散曲線。

在復(fù)合材料分層損傷的檢測方面，SHOJA等[17]使用低頻導(dǎo)波對復(fù)合材料層壓板中的分層檢測進行了數(shù)值模擬；ALLEYNE等[18]研究了蘭姆波與各種缺陷的相互作用，分析了不同蘭姆波模態(tài)對缺陷的敏感性；GUO等[19]證實了蘭姆波可以用于檢測層合板中的分層，且分層的深度是影響波傳播特性的重要因素；RAMADAS等[20]觀察到分層引起的多模態(tài)轉(zhuǎn)換現(xiàn)象;SOHN等[21]設(shè)計了先進的導(dǎo)波場信號處理方法，檢測了復(fù)合材料中的分層;GAO等[22]使用電磁超聲換能器，對分層結(jié)構(gòu)中的分層進行了有效檢測;SOLEIMANPOUR等[23]研究了復(fù)合材料層合梁中分層檢測的非線性導(dǎo)波方法。

上述文獻表明，目前多層復(fù)合材料結(jié)構(gòu)頻散特性的理論求解方法主要依賴于三維彈性理論的傳遞矩陣法和半解析有限元法。傳遞矩陣法可以得到精確的頻散控制方程，但當(dāng)層數(shù)增加時，求解超越方程異常困難，且易出現(xiàn)漏根等問題。半解析有限元方法對不同的波導(dǎo)截面具有較大的適用性，但需要編寫代碼。因此，有必要建立一種高效、易實現(xiàn)的求解方法，對復(fù)合材料管道中的導(dǎo)波頻散特性進行分析。

雖然目前基于超聲導(dǎo)波的復(fù)合材料分層損傷檢測技術(shù)取得了一定的進展，但大多研究集中在復(fù)合材料板類結(jié)構(gòu)，對復(fù)合管道的分層損傷檢測較少。究其原因，發(fā)現(xiàn)復(fù)合材料管道中的導(dǎo)波傳播特性尚不明確，分層損傷對導(dǎo)波的影響機制以及檢測原理亟待研究。

為了解決上述問題，本文結(jié)合Floquet周期性邊界條件和有限元特征頻率求解法，高效、簡便地分析復(fù)合材料管道中超聲波的頻散特性；進一步從數(shù)值仿真角度，針對反射和模態(tài)轉(zhuǎn)換兩種損傷因子，研究分層損傷對導(dǎo)波傳播特性的調(diào)制規(guī)律，為后續(xù)試驗研究提供堅實的理論指導(dǎo)。

1 基于Floquet邊界條件的特征頻率方法

1.1 理論模型

Floquet邊界條件通常作用于平行表面，定義為：

udst=usrce-ikF(rdst-rsrc)

(1)

式中，u為位移場；dst，src分別為目標邊界、源邊界；kF為波數(shù)；r為空間位置。

假設(shè)一個無限長的管道波導(dǎo)(見圖1)，采用Floquet邊界條件，將域從圖1縮減到圖2，用元胞結(jié)構(gòu)進行等效求解。具體原理如式(2)～(10)所示。

圖1 無限長管道示意Fig.1 Infinitely large pipeline

圖2 元胞結(jié)構(gòu)與Floquet邊界條件Fig.2 Cell structure and Floquet boundary conditions

無限長的管道波導(dǎo)結(jié)構(gòu)的位移為:

u(r,t)=U(r)e-i(mθ+kzz)+iωt

(2)

圖2中的元胞結(jié)構(gòu)詳細圖解如圖3所示，幾何變量和物理變量的變化如式(3)～(8)所示。

r=r′

(3)

θ=θ′+mθθs

(4)

z=z′+mzsz

(5)

kr=0

(6)

(7)

(8)

將式(3)～(8)代入式(2)得：

u(r,t)=e-i(k′θmθR0θs+k′zmzsz)U(z′)e-i(mθ′+kzz′)+iωt

(9)

如果e-i(k′θmθR0θs+k′zmzsz)是連續(xù)的，則有：

u(r,t)=U(r)e-i(mθ′+kzz′)+iωt

(10)

將式(10)與式(2)對比可知，可以用小的元胞結(jié)構(gòu)代替整個無限大的管道結(jié)構(gòu)(如圖3所示)，使得整個管道結(jié)構(gòu)表示為元胞結(jié)構(gòu)的周期性延拓，從而簡化計算。

圖3 元胞結(jié)構(gòu)示意Fig.3 Schematic diagram of cell structure

1.2 等效的COMSOL模型

采用COMSOL MULTIPHYSIC有限元仿真軟件中的固體力學(xué)模塊，考慮R方向厚度為d的管道結(jié)構(gòu)，基于上述理論建立如圖4所示的元胞結(jié)構(gòu)(以8層管道為例)。在Z方向上加Floquet周期性邊界條件，在T方向上加連續(xù)性邊界條件，厚度方向上使用自由邊界條件。

圖4中z是源邊界和目標邊界的距離，即Floquet 周期性邊界條件中的晶格平移量，是由Z方向的本征波矢kz和邊界的相對距離來確定的。需要注意的是，引入Floquet周期性邊界條件會使得頻散曲線在布洛赫邊界發(fā)生折疊，對于管道結(jié)構(gòu)來說這個折疊是虛假的?？紤]到布里淵邊界處波矢量的最大取值為π/z，可以使用較小的z值，將求解的特征頻率折疊推到極高頻，然后截斷，即可得到所需的頻散關(guān)系。使用較小的z值還可以簡化模型的復(fù)雜性，從而減少計算時間。根據(jù)經(jīng)驗可知，z的取值為管道厚度d的1/10。θs為所取元胞結(jié)構(gòu)的弧度，取值也相對小一點，可以減少網(wǎng)格的數(shù)量，從而減少計算時間。

圖4 8層管道元胞結(jié)構(gòu)Fig.4 Cell structure of an eight-layer pipeline

同時，在利用有限元方法進行求解的過程中，由于周期性條件的引入會導(dǎo)致計算收斂性變差。為了提高收斂性，在進行網(wǎng)格劃分時，需要注意使互為周期性的兩個邊界上的網(wǎng)格完全一致。通過波矢kz進行參數(shù)化掃描，可以求解出質(zhì)點振動的本征方程對應(yīng)于不同波矢的特征頻率，進而求解多層復(fù)合管道的頻散曲線。

2 復(fù)合材料管道中的頻散特性

以圖4的8層管道為例，管道內(nèi)徑為60 mm，長度為500 mm，采用上述方法對其頻散特性進行分析。表1列出了單層材料屬性，結(jié)構(gòu)鋪層順序為[0/45/-45/90]s,每層厚度為0.15 mm。在COMSOL中建立相應(yīng)的模型，元胞高度為管道原厚度，設(shè)置為1.2 mm，源邊界與目標邊界的距離z設(shè)置為0.12 mm，θs為0.3°。網(wǎng)格為自由四面體網(wǎng)格，網(wǎng)格尺寸采用預(yù)定義里的coarser(粗糙的)。得到的頻散曲線見圖5。

表1 材料參數(shù)Tab.1 Material parameters

圖5 [0/45/-45/90]s復(fù)合材料管道中的頻散曲線

進一步對模態(tài)的波結(jié)構(gòu)進行分析。取如下3個模態(tài)：k1=942.48 m-1，f1=181 kHz；k2=314.16 m-1，f2=166 kHz；k3=314.16 m-1，f3=279 kHz，在COMSOL中提取對應(yīng)的特征向量，繪制波結(jié)構(gòu)如圖6所示?？梢钥闯?，圖6(a)是以ur為主導(dǎo)，判斷其為L(0,1)的波結(jié)構(gòu)；圖6(b)是以uθ為主導(dǎo)，判斷其為T(0,1)的波結(jié)構(gòu)；圖6(c)是以uz為主導(dǎo)，判斷其為L(0,2)的波結(jié)構(gòu)。

使用半解析有限元法[24]對上述方法求解得到的頻散曲線進行驗證。選用一維二次單元在厚度方向進行有限元離散，每層厚度取2個單元，可以分析多層復(fù)合管道中導(dǎo)波的頻散特性。圖5初步驗證了本文提出的特征頻率方法求解頻散曲線的正確性。

圖6 波結(jié)構(gòu)示意

采用Abaqus動力學(xué)分析模塊，對200 kHz下的T(0,1)，L(0,1)，L(0,2)模態(tài)進行模擬，提取其傳播速度，與上述方法作比較。

在Abaqus中建立對應(yīng)的三維管道，長度為0.5 m，單元尺寸為1 mm，單元類型為C3D8R。采用隱式動力學(xué)分析步，設(shè)置固定增量步為2.5×10-7s。在管道的一端施加固定約束，另一端施加壓力載荷，大小為五波峰調(diào)制的正弦信號:

(11)

式中，N為周期數(shù)，N=5；f0為中心頻率，kHz，f0=200 kHz。

依據(jù)上述波結(jié)構(gòu)分析結(jié)果，分別用3個方向上的力來加載，激勵出L(0,1),L(0,2),T(0,1)模態(tài)，不同模態(tài)傳播的云圖如圖7所示。

圖7 不同模態(tài)傳播的位移云圖

管道軸向一條路徑上所有節(jié)點的位移時域信號，通過一維傅里葉變換，將信號從時域轉(zhuǎn)換至波數(shù)域，結(jié)果見圖8。計算可得到相應(yīng)的相速度。由表2可以看出，L(0,2),T(0,1)模態(tài)的COMSOL解(基于Floquet邊界條件的有限元特征頻率法)和Abaqus解(動力學(xué)方法)相差很小，但L(0,1)模態(tài)的相對誤差較大。究其原因，L(0,1)模態(tài)的波長約為6 mm，仿真中單元尺寸設(shè)為1 mm，不滿足一個波長內(nèi)至少8～10個單元的要求，可能是導(dǎo)波L(0,1)模態(tài)的速度誤差較大的原因。

圖8 不同模態(tài)的位移空間信號與波數(shù)域信號

表2 COMSOL求解值與Abaqus求解值對比Tab.2 Comparison of COMSOL solved values withAbaqus solved values

注：COMSOL解為基于Floquet邊界條件的有限元特征頻率法;Abaqus解為動力學(xué)方法。

觀察圖5，基于3種方法求解的導(dǎo)波頻散特性吻合度很高，驗證了本文提出的有限元特征頻率求解法求解頻散曲線的正確性，同時也驗證了本文建立的Abaqus模型的正確性。

3 分層損傷對導(dǎo)波傳播特性的影響

3.1 損傷模型的建立

為了進一步研究分層損傷對導(dǎo)波傳播特性的調(diào)制規(guī)律，在圖7模型的基礎(chǔ)上建立分層損傷。分層損傷通常認為是零體積的，且無法傳遞面內(nèi)剪力。有限元軟件中具體的設(shè)置方法是將有損傷的地方解除tie約束。

采用激勵頻率200 kHz下L(0,2)模態(tài)檢測損傷，以圖9為例，損傷大小為2 mm×2 mm，距離左端部250 mm處，位于第4層和第5層中間位置。傳播的云圖見圖9，觀察到經(jīng)過損傷后有4股波出現(xiàn)，從左到右依次為經(jīng)損傷反射的L(0,2)R、經(jīng)損傷模態(tài)轉(zhuǎn)換反射的L(0,1)R、經(jīng)損傷模態(tài)轉(zhuǎn)換透射的L(0,1)T以及經(jīng)損傷透射的L(0,2)T。

圖9 含損傷管道反射波、透射波云圖Fig.9 Nephograms of reflected and transmitted wavesin a pipe with delamination damage

(a)空間域

(b)波數(shù)域圖10 含損傷多層復(fù)合材料管道結(jié)構(gòu)中的導(dǎo)波信號Fig.10 Wave guided signals in a multilayer composite pipewith delamination damage

選取管道上的一條路徑(如圖9所示)，提取t=83.75 μs時刻下的空間域信號(見圖10(a))，作一維傅里葉變換(見圖10(b))。

由圖10(b)可以看出兩個波數(shù)峰值：k1=36 m-1，k2=177 m-1，分別對應(yīng)模態(tài)為L(0,2)和L(0,1)，因此，可判斷發(fā)生了從L(0,2)向L(0,1)的模態(tài)轉(zhuǎn)換現(xiàn)象。

3.2 損傷深度對導(dǎo)波傳播特性的影響

固定復(fù)合管道長度為500 mm，損傷尺寸為20 mm×20 mm，損傷發(fā)生在距離管道端部250 mm 處，通過改變分層發(fā)生的界面位置，分析其對導(dǎo)波傳播特性的影響規(guī)律。不同損傷位置位移的空間域信號如圖11所示。

圖11 不同損傷深度的導(dǎo)波空間域信號Fig.11 Spatial domain signals of guided waves withdelamination damage in different depths

定義4種歸一化的損傷因子IL(0,2)T，IL(0,1)T，IL(0,1)R，IL(0,2)R：

(12)

(13)

(14)

(15)

式中，IL(0,2)T，IL(0,1)T，IL(0,1)R，IL(0,2)R分別為L(0,2)模態(tài)透射因子、L(0,1)模態(tài)透射因子、L(0,1)模態(tài)反射因子以及L(0,2)模態(tài)反射因子；AL(0,2)T，AL(0,1)T，AL(0,1)R，AL(0,2)R分別為L(0,2)模態(tài)透射波幅值最大值、L(0,1)模態(tài)透射波幅值最大值、L(0,1)模態(tài)反射波幅值最大值以及L(0,2)模態(tài)反射波幅值最大值；A0為入射波幅值最大值。

從圖12可以觀察到，L(0,2)T隨損傷深度變換呈波動變化，L(0,2)R和發(fā)生模態(tài)轉(zhuǎn)換的L(0,1)T都隨深度增大,呈現(xiàn)先減小、后增大的趨勢，且都是損傷在第4層和第5層之間時達到最小值。隨深度增大，L(0,1)R呈現(xiàn)先增大、后減小的趨勢。

圖12 不同損傷位置下的透射、反射、模態(tài)轉(zhuǎn)換因子Fig.12 Transmission,reflection and mode conversionindicators with delamination damage in different depths

從圖12還可以看出，發(fā)生模態(tài)轉(zhuǎn)換得到的L(0,1)T的透射因子明顯大于L(0,2)R的反射因子，因此相對于L(0,2)R，經(jīng)由模態(tài)轉(zhuǎn)換得到的L(0,1)T對分層損傷檢測靈敏度更高。從以上的分析可以得出結(jié)論，L(0,1)T對損傷的敏感度明顯高于L(0,2)R，并呈現(xiàn)一定的規(guī)律，這也為后續(xù)開展定位損傷檢測提供基礎(chǔ)。

3.3 損傷尺寸對導(dǎo)波傳播特性的影響

固定損傷位置在第4層和第5層之間，固定軸向尺寸為20 mm，改變損傷的周向尺寸l，分析其對導(dǎo)波傳播特性的影響規(guī)律。

由圖13可以觀察到，L(0,2)T隨損傷尺寸的改變呈現(xiàn)波動變換；L(0,2)R反射因子一開始隨著尺寸的變大而變大，當(dāng)尺寸與波長的比值接近0.4時，反射因子達到最大，隨著損傷尺寸的進一步變大，L(0,2)R的反射因子隨損傷尺寸的變大而呈平緩上升的趨勢；L(0,1)R反射因子一開始隨尺寸的增大而逐漸減小，在損傷尺寸與波長的比值接近0.4時，L(0,1)R反射因子達到最小，隨著尺寸的進一步增大，L(0,1)R反射因子隨損傷尺寸的變大呈平緩上升的趨勢；L(0,1)T透射因子一開始隨尺寸的增大而逐漸減小，在損傷尺寸與波長的比值接近0.4時，L(0,1)R反射因子達到最小，隨著尺寸的進一步增大，L(0,1)R反射因子隨損傷尺寸的變大呈平緩上升的趨勢。

圖13 不同損傷尺寸下的透射、反射、模態(tài)轉(zhuǎn)換因子Fig.13 Transmission,reflection and mode conversionindicators with delamination damage of different sizes

從圖13還可以看出，損傷尺寸與波長之比小于0.5時，L(0,2)R反射因子顯然大于經(jīng)模態(tài)轉(zhuǎn)換得到的L(0,1)的透射因子與反射因子，相對于L(0,1)T和L(0,1)R，L(0,2)R對分層損傷檢測靈敏度更高，故此時適宜使用L(0,2)R進行損傷檢測。但在損傷尺寸與波長之比介于0.5與1.6之間時，經(jīng)模態(tài)轉(zhuǎn)換的L(0,1)R的反射因子顯然大于其透射因子與L(0,2)R的反射因子，此時則應(yīng)使用經(jīng)模態(tài)轉(zhuǎn)換得到的L(0,1)R進行缺陷檢測，這也為后續(xù)開展定量損傷檢測提供了依據(jù)。

4 結(jié)論

本文基于Floquet周期性邊界條件的有限元特征頻率法，高效、簡便地分析了多層復(fù)合材料管道中的超聲導(dǎo)波頻散特性；并建立了含分層損傷的復(fù)合材料管道中的超聲導(dǎo)波動力學(xué)模型，針對反射和模態(tài)轉(zhuǎn)換兩種因子，研究了損傷深度和尺寸對導(dǎo)波傳播特性的影響，得出以下結(jié)論。

(1)基于有限元特征頻率法，引入Floquet周期性邊界條件，將復(fù)合材料管道中的超聲波傳播的頻散控制方程轉(zhuǎn)化為元胞結(jié)構(gòu)的特征頻率方程，高效地分析了導(dǎo)波頻散特性。該方法不需要編程，且對不同波導(dǎo)結(jié)構(gòu)具有普適性。

(2)在分層損傷處會發(fā)生L(0,2)到L(0,1)的模態(tài)轉(zhuǎn)換。當(dāng)損傷發(fā)生在不同深度時，經(jīng)損傷發(fā)生模態(tài)轉(zhuǎn)換的L(0,1)T對損傷的敏感度最大。損傷尺寸與波長比值小于0.5時，L(0,2)R對損傷的敏感度最大；當(dāng)損傷尺寸與波長比值大于0.5時，經(jīng)損傷發(fā)生模態(tài)轉(zhuǎn)換的L(0,1)R對損傷的敏感度最大。

(3)多損傷因子分析可為分層損傷的定量檢測提供依據(jù)。