余旭東，秦 榮，沈 海，左 鵬，邵照宇

（1.北京航空航天大學(xué)，宇航學(xué)院，北京 100191；2.新加坡科技研究局，先進(jìn)再制造和技術(shù)中心，新加坡 637143；3.中國(guó)航發(fā)商用航空發(fā)動(dòng)機(jī)有限責(zé)任公司，先進(jìn)技術(shù)研究部，上海 200241）

多層復(fù)合圓管結(jié)構(gòu)，相較于單層圓管，在比強(qiáng)度、比模量、線膨脹系數(shù)、耐腐蝕性、抗疲勞性和阻尼系數(shù)等方面具有顯著優(yōu)勢(shì)，已發(fā)展成為航空航天、能源化工、船舶海洋等工業(yè)領(lǐng)域廣泛應(yīng)用的重要承力構(gòu)件［1］。因其在高溫、超壓、腐蝕等復(fù)雜苛刻的環(huán)境/條件下長(zhǎng)時(shí)服役，復(fù)合圓管結(jié)構(gòu)受疲勞、老化、應(yīng)力集中、沖擊載荷等因素影響，在結(jié)構(gòu)內(nèi)部及層間界面處易產(chǎn)生多尺度、多形態(tài)的缺陷損傷，如分層脫粘、變形、裂紋、腐蝕、凹坑等，進(jìn)而嚴(yán)重削弱系統(tǒng)整體的承載能力和服役壽命。因此，針對(duì)多層復(fù)合圓管結(jié)構(gòu)開展快速有效的無損檢測(cè)尤為重要。超聲導(dǎo)波是一種聲場(chǎng)能量能夠沿波導(dǎo)結(jié)構(gòu)傳播的彈性波，并在圓管類結(jié)構(gòu)中可遍及整個(gè)壁厚；以超聲導(dǎo)波為基礎(chǔ)的無損檢測(cè)與結(jié)構(gòu)健康監(jiān)測(cè)方法，因其具有傳播距離長(zhǎng)、檢測(cè)效率高、激發(fā)/接收方式靈活、缺陷辨識(shí)能力強(qiáng)等優(yōu)點(diǎn)，已被應(yīng)用于圓管類結(jié)構(gòu)的缺陷檢測(cè)與表征。

超聲導(dǎo)波檢測(cè)技術(shù)的發(fā)展起始于對(duì)平板中傳播導(dǎo)波的詳細(xì)研究。Lamb最早在研究無限大平板中傳播的正弦波時(shí)發(fā)現(xiàn)蘭姆波現(xiàn)象［2］，包括振動(dòng)模態(tài)沿平板中性面對(duì)稱分布的對(duì)稱S模式波（symmetric mode）與振動(dòng)模態(tài)沿中性面反對(duì)稱分布的反對(duì)稱A模式波（antisymmetric mode），并給出其嚴(yán)格的解析描述。Love［3］在研究各向同性彈性半空間上覆蓋另一種層狀介質(zhì)的結(jié)構(gòu)中的波傳播過程中發(fā)現(xiàn)水平剪切（SH）波現(xiàn)象，平板中SH波位移與蘭姆波波場(chǎng)正交，其振動(dòng)模態(tài)沿板厚呈現(xiàn)水平剪切特征。類似于板波（即Lamb波與SH波）［4］，沿著圓管及圓柱狀結(jié)構(gòu)軸向傳播的柱面導(dǎo)波亦被發(fā)現(xiàn)［5］。Cawley等［6］對(duì)管道中軸向?qū)Рǖ念l散特性和缺陷散射規(guī)律開展深入研究，并將軸向?qū)Р夹g(shù)成功應(yīng)用于油氣工業(yè)管道的長(zhǎng)距離、快速檢測(cè)。通過沿圓管軸向傳播的柱面導(dǎo)波，可從單個(gè)超聲換能器位置檢測(cè)獲取沿管道軸向約50 m（即正負(fù)方向各25 m）范圍內(nèi)的裂紋、腐蝕等缺陷，并能將其與管道固有的焊接、法蘭結(jié)構(gòu)等有效區(qū)分［7］。隨著工業(yè)用復(fù)合圓管的尺寸、層數(shù)的不斷增加以及管道所處周圍介質(zhì)（如液體、包覆層等）的影響，軸向?qū)Рㄔ诠苤袀鞑r(shí)會(huì)呈現(xiàn)較大的衰減系數(shù)，其檢測(cè)范圍和檢測(cè)能力愈受限制［8］，難以對(duì)圓管結(jié)構(gòu)中的危險(xiǎn)局部位置開展更為準(zhǔn)確的缺陷檢測(cè)與環(huán)面材料表征。而沿著圓管環(huán)面?zhèn)鞑サ闹芟驅(qū)Рㄒ嗑哂袙卟檎麄€(gè)橫截面的能力并對(duì)軸向裂紋、徑向疲勞裂紋和分層缺陷等具有較高的檢測(cè)靈敏度，可有效表征復(fù)合圓管局域內(nèi)的多種類型缺陷，相關(guān)檢測(cè)技術(shù)因此也獲得越來越多的關(guān)注與研究。

厘清多層復(fù)合圓管中周向?qū)РＪ降念l散特性（即頻散曲線，描述導(dǎo)波相速度或能量速度隨著頻率改變的變化規(guī)律）是開展相關(guān)環(huán)面檢測(cè)的首要前提。頻散曲線可為選擇最合適的導(dǎo)波模式和激勵(lì)頻率范圍提供重要依據(jù)，以實(shí)現(xiàn)特定導(dǎo)波模式（對(duì)）的最佳激勵(lì)和波傳播條件；同時(shí)其亦能指導(dǎo)超聲導(dǎo)波檢測(cè)信號(hào)的分析與處理，實(shí)現(xiàn)缺陷的有效檢測(cè)定位以及對(duì)其他結(jié)構(gòu)特征的辨識(shí)。然而，工業(yè)用多層復(fù)合圓管呈現(xiàn)顯著的各向異性，尤其是復(fù)合材料圓管中各層不同的纖維取向亦帶來不可忽略的非均勻性。復(fù)雜的材料/結(jié)構(gòu)特征將使得多層復(fù)合環(huán)面中周向?qū)Р▎栴}的解析求解尤為艱巨。

文獻(xiàn)中關(guān)于管中周向?qū)Р▎栴}研究及頻散曲線提取工作，相較于軸向?qū)Рㄉ?，起始于Viktorov明晰圓柱狀結(jié)構(gòu)中周向?qū)Рㄅc平板結(jié)構(gòu)中導(dǎo)波傳播特性的主要差別在于圓柱狀結(jié)構(gòu)凸面可支撐表面波的傳播，其凹面則不可支撐，而平板的兩個(gè)表面均可支撐表面波的傳播，并可由此推導(dǎo)出單層環(huán)面中類蘭姆波的特征方程［9］。與平板中的導(dǎo)波傳播類似，在均勻各向同性介質(zhì)環(huán)面中，周向?qū)Рㄒ嗫煞譃橹芟蛩郊羟校–SH）波和周向類蘭姆（CLT）波。Liu等［10］基于頻散方程的解析表達(dá)式給出不同厚度/半徑比下的各向同性介質(zhì)環(huán)面中CLT模式的頻散曲線及位移場(chǎng)的數(shù)值結(jié)果。此外，Zhao等［11］利用貝塞爾函數(shù)擴(kuò)展位移場(chǎng)解析推導(dǎo)出CSH模式的頻散方程，計(jì)算了CSH模式的頻散曲線以及沿管壁厚度方向的位移和應(yīng)力場(chǎng)。然而，該方法在較大環(huán)半徑（即貝塞爾函數(shù)階數(shù)較高）情形下較難適用，高階數(shù)、小變量的貝塞爾函數(shù)求解時(shí)會(huì)出現(xiàn)計(jì)算效率低、數(shù)值上不穩(wěn)定、可計(jì)算頻散空間有限等問題［12］。為了克服上述尋根（root-finding）方法的局限性，Gridin等［13］在研究圓形鋼環(huán)中周向?qū)Рǖ膫鞑栴}時(shí)，引入多種漸近方法替代控制特征方程中的貝塞爾函數(shù)以推導(dǎo)出周向?qū)РＪ降念l散關(guān)系，但該類高頻漸進(jìn)分析方法僅僅適用于單層各向同性介質(zhì)環(huán)面。

為進(jìn)一步研究各向異性材料多層環(huán)面中周向?qū)Рǖ念l散曲線，Lowe［14］總結(jié)傳遞矩陣法（transfer matrix method）和全局矩陣法（global matrix method），構(gòu)建并求解多層結(jié)構(gòu)導(dǎo)波問題的頻散方程。然而這兩種方法均有局限性，即當(dāng)結(jié)構(gòu)層數(shù)/厚度逐漸增大時(shí)，利用傳遞矩陣法求解超越方程變得既困難又費(fèi)時(shí)，將會(huì)出現(xiàn)‘largef·d’問題（f為頻率，d為結(jié)構(gòu)厚度）。該類問題是指在導(dǎo)波求解過程中出現(xiàn)頻厚積（f·d）較大的情形時(shí)，會(huì)出現(xiàn)矩陣的病態(tài)現(xiàn)象，導(dǎo)致計(jì)算不穩(wěn)定且誤差增大。全局矩陣法對(duì)傳遞矩陣法進(jìn)行了修正，可消除傳遞矩陣法中出現(xiàn)大頻厚積情形計(jì)算精度損失的問題，但對(duì)于層數(shù)較多的情況，全局矩陣方法仍需處理較大的矩陣，計(jì)算量顯著增大。

采用上述兩種矩陣方法求解各向異性多層復(fù)合圓管中的周向?qū)Рl散問題已有相關(guān)研究。Towfighi等［15］引入傅里葉級(jí)數(shù)擴(kuò)展未知位移函數(shù)，以求解橫觀各向同性圓管中周向?qū)Рǖ念l散曲線。Yu等［16］基于線性三維彈性，利用Legendre正交多項(xiàng)式級(jí)數(shù)展開法確定均質(zhì)正交各向異性圓管中周向?qū)Рǖ念l散曲線。上述方法求解各向異性多層圓管結(jié)構(gòu)中周向?qū)Р▎栴}通常計(jì)算量較大。為進(jìn)一步提升計(jì)算效率，Lugovtsova等［17］基于比例邊界有限元法（scaled boundary finite element method，SBFEM）計(jì)算得到鋁-碳纖維多層復(fù)合圓管中周向?qū)Рǖ念l散曲線及模態(tài)。Lin等［18］將Floquet邊界條件方法與掃頻有限元方法（sweeping frequency finite element modeling，SFFEM）相結(jié)合，高效求解了正交各向異性多層復(fù)合材料圓管中的周向?qū)Рl散曲線并加以實(shí)驗(yàn)驗(yàn)證。

半解析有限元（semi-analytical finite element，SAFE）方法被廣泛用于求解描述導(dǎo)波傳播的特征方程，通過采用有限元網(wǎng)格對(duì)波導(dǎo)結(jié)構(gòu)的橫截面進(jìn)行離散，而在波傳播方向上（即與波導(dǎo)截面垂直的方向）采用解析描述［19］，將三維波導(dǎo)問題轉(zhuǎn)化為二維模型進(jìn)行求解（不改變?nèi)S線彈性假設(shè)），進(jìn)而大幅減小計(jì)算成本。SAFE方法最初被用于替代一些較為傳統(tǒng)的解析求解過程（全局矩陣法、傳遞矩陣法等）以解決經(jīng)典的‘largef·d’導(dǎo)波問題。Van Velsor［20］提出SAFE方法數(shù)學(xué)框架以求解含粘彈性材料的多層各向同性圓管中周向?qū)Рǖ念l散曲線。Matuszyk等［21］提出一種將SAFE方法與完全匹配層（PML）技術(shù)相結(jié)合的手段，求解了埋地圓管中的周向?qū)Рl散曲線。Lin等［22］基于所提出的場(chǎng)變分方法，改進(jìn)了已應(yīng)用于直波導(dǎo)的半解析技術(shù)以研究各向異性圓管中周向?qū)Р▊鞑栴}，并計(jì)算了復(fù)雜周向截面結(jié)構(gòu)中周向?qū)Рǖ念l散曲線。文獻(xiàn)報(bào)道中，通過獨(dú)立編程手段執(zhí)行SAFE方法，已成功求解各向異性［23］、埋地［24］或水浸［25］的圓管結(jié)構(gòu)中周向?qū)Рǖ念l散曲線。而本課題組與合作者以及多個(gè)國(guó)內(nèi)外研究組則借助于商用有限元分析平臺(tái)COMSOL Multiphysics在笛卡爾坐標(biāo)系下執(zhí)行SAFE方法，已成功求解具有任意材料屬性［26］、任意截面形狀［27］以及任意預(yù)應(yīng)力邊界條件［28］的波導(dǎo)問題，并能求解固埋［29］、液浸［30］等情形下的頻散曲線。其可望精確計(jì)算復(fù)雜圓管結(jié)構(gòu)中周向?qū)Рǖ念l散曲線，對(duì)多層復(fù)合環(huán)面更具優(yōu)勢(shì)。

本文進(jìn)一步借助COMSOL Multiphysics商用平臺(tái)在柱坐標(biāo)系下執(zhí)行SAFE方法，通過控制特征方程數(shù)學(xué)推導(dǎo)與弱形式偏微分方程求解模塊，提出多層環(huán)面的SAFE計(jì)算模型，無需編寫冗長(zhǎng)程序代碼，即可高效、準(zhǔn)確計(jì)算具有任意材料屬性、任意曲率、任意層數(shù)的多層復(fù)合圓管結(jié)構(gòu)中周向?qū)Рǖ念l散曲線。通過與文獻(xiàn)中已有結(jié)果比較，交互驗(yàn)證該數(shù)值計(jì)算方法的準(zhǔn)確性。本文選取工業(yè)中兩種典型復(fù)合圓管結(jié)構(gòu)，即金屬-碳纖維雙層環(huán)面和多層碳纖維樹脂基復(fù)合材料環(huán)面，計(jì)算獲取其中傳播的所有周向?qū)РＪ降念l散曲線，為后續(xù)對(duì)該類結(jié)構(gòu)開展周向檢測(cè)與表征提供有效基礎(chǔ)。

1 理論框架及半解析有限元法

1.1 周向?qū)Рǖ目刂品匠?/h3>

對(duì)具有任意截面形狀波導(dǎo)結(jié)構(gòu)中沿軸向傳播的導(dǎo)波問題通常采用笛卡爾坐標(biāo)系進(jìn)行描述，而針對(duì)具有任意曲率板類結(jié)構(gòu)中周向?qū)Р▊鞑栴}，需采用柱坐標(biāo)系進(jìn)行描述。在柱坐標(biāo)系（r，θ，z）中，彈性體應(yīng)力矩陣的正應(yīng)力分量為σrr，σθθ，σzz，切應(yīng)力分量為σθz，σrz，σrθ；而應(yīng)變包括法向應(yīng)力長(zhǎng)度方向的正應(yīng)變?chǔ)舝r，εθθ，εzz以及切應(yīng)力長(zhǎng)度方向的切應(yīng)變?chǔ)纽葄，εrz，εrθ。在柱坐標(biāo)系中，ur，uθ和uz分別表示r，θ和z方向的變形，由于z軸方向的分量不隨z變化，其應(yīng)變張量分量為

取彈性波導(dǎo)中一個(gè)無限小的體積微元，由力和位移表示的平衡方程以及本構(gòu)方程分別為

式中：ρ為介質(zhì)密度；f為體力；t為時(shí)間變量；Cikjl為四階彈性系數(shù)張量。如不計(jì)體力，即fi=0，式（2）可表示為r，θ和z三個(gè)方向上應(yīng)力分量的運(yùn)動(dòng)方程。將式（1）和式（3）代入r，θ和z三個(gè)方向上應(yīng)力分量的運(yùn)動(dòng)方程中可分別得到周向?qū)Рㄑ豶，θ和z方向上位移分量的運(yùn)動(dòng)方程。

周向?qū)Рㄊ窃趫A管中沿著圓周θ方向傳播的導(dǎo)波，其傳播路徑閉合；周向?qū)Р刂品匠痰慕鈛r，uθ，uz分別為r軸，θ軸和z軸三個(gè)方向上的位移分量。在均勻各向同性環(huán)面中，周向?qū)Рǚ种芟蛩郊羟校–SH）波和周向類蘭姆（CLT）波兩種。圖1顯示了圓管中周向?qū)Рǖ馁|(zhì)點(diǎn)位移，其中CSH模式的質(zhì)點(diǎn)運(yùn)動(dòng)方向（沿z軸）與波傳播θ方向垂直，且在r軸和θ軸方向的位移分量均為零，即ur=uθ=0；CLT模式質(zhì)點(diǎn)位移在r-θ面內(nèi)，位移在z軸方向的位移分量為零，即uz=0。

圖1 周向?qū)РＪ降馁|(zhì)點(diǎn)位移示意圖Fig.1 Schematics of particle displacement for circumferential guided waves modes

1.2 周向?qū)Рǖ陌虢馕鲇邢拊⊿AFE）求解方法

本文基于COMSOL Multiphysics商用有限元分析平臺(tái)執(zhí)行周向?qū)Рǖ腟AFE計(jì)算方法，通過其中的PDE weak form功能模塊求解。對(duì)于周向傳播的導(dǎo)波，彈性介質(zhì)中的波場(chǎng)不依賴于z方向，其位移矢量u只是r和θ的函數(shù)。在廣義平面應(yīng)變假設(shè)下，zθ平面中的z向延伸到無窮大；波傳播沿θ方向，并假設(shè)其作簡(jiǎn)諧波形式e-iωt，由函數(shù)eipθ描述。波導(dǎo)中位移矢量則可描述為

由于符合廣義平面應(yīng)變假設(shè)，z軸方向的分量不隨z變化，即?/?z=0。將波導(dǎo)中導(dǎo)波模式質(zhì)點(diǎn)位移分別對(duì)柱坐標(biāo)系下三個(gè)方向求偏導(dǎo)可得：

將式（5）代入r、θ和z方向的運(yùn)動(dòng)方程，等式兩邊同除ei(pθ-ωt)項(xiàng)，化簡(jiǎn)可分別得到三個(gè)方向上控制微分方程。計(jì)算真空環(huán)境下環(huán)面中周向?qū)Рǖ念l散曲線，圓管的內(nèi)外表面須滿足無應(yīng)力邊界條件，即在r=a和r=b處，σrr=0，σrθ=0，σrz=0。

控制微分方程及邊界條件對(duì)應(yīng)COMSOL Multiphysics中特征值問題的一般輸入形式為

式中：λ為待求解的特征值，即周向?qū)РＪ綄?duì)應(yīng)的角波數(shù)；u為對(duì)應(yīng)的特征向量；系數(shù)矩陣ea、da、c、α、β、a與控制方程參數(shù)對(duì)應(yīng)。將控制偏微分方程轉(zhuǎn)換為特征值問題的一般表達(dá)式，以矩陣形式表示為

式中：δij為Kronecker符號(hào)。式（7）為二次特征值問題，為將其轉(zhuǎn)換為線性形式進(jìn)行簡(jiǎn)化求解，需引入一個(gè)新變量Vj=pUj；式（6）待求解的特征向量定義為u=[U1U2U3V1V2V3]T。根據(jù)式（7）可確定COMSOL Multiphysics中待求解的特征值λ=p，以及相應(yīng)的系數(shù)矩陣。剛度矩陣采用Voigt標(biāo)注形式。對(duì)該特征值問題在頻域下求解時(shí)，通過指定單個(gè)頻率ω，可計(jì)算獲得所研究波導(dǎo)結(jié)構(gòu)中所有導(dǎo)波模式的特征值波數(shù)p及其特征向量u=[U1U2U3V1V2V3]T；在多個(gè)頻率下重復(fù)上述計(jì)算即可獲得各導(dǎo)波模式在所研究頻率范圍內(nèi)的頻散曲線。

本文中研究多層復(fù)合環(huán)面波導(dǎo)問題的SAFE計(jì)算模型如圖2所示。在COMSOL Multiphysics平臺(tái)中選擇二維軸對(duì)稱組件，采用系數(shù)型偏微分方程（coefficient PDE）模塊；模型構(gòu)建采用柱坐標(biāo)系，坐標(biāo)變量為r，θ，z；圓管內(nèi)徑為a，外徑為b，管壁厚度即為b-a。設(shè)置圓管內(nèi)外表面滿足無應(yīng)力邊界條件；圓管在z軸無限延伸，則沿z向設(shè)置波導(dǎo)截面的周期性邊界條件。周向?qū)Рㄑ卅确较騻鞑?，波函?shù)用eipθ表示，其中角波數(shù)p與環(huán)向波數(shù)k的關(guān)系為p=kb。該圓管的SAFE模型采用二階拉格朗日三角形網(wǎng)格進(jìn)行截面劃分，在z=0和z=h邊界上設(shè)置周期性邊界條件以模擬z向無限大圓管，如式（8）所定義。在該模型設(shè)置中，h可選取相對(duì)壁厚較小的高度值，減少劃分網(wǎng)格數(shù)量，并有效避免導(dǎo)波振動(dòng)模態(tài)在較大h值情形下的重復(fù)分布［31］。

圖2 COMSOL Multiphysics平臺(tái)中建立的SAFE計(jì)算模型Fig.2 Diagrams for SAFE modelling via COMSOL Multiphysics

通過環(huán)面SAFE模型計(jì)算提取周向?qū)Рǖ念l散關(guān)系時(shí)，先根據(jù)上述求解方法在頻率-波數(shù)域求解特征方程的根p，由角波數(shù)p計(jì)算出在任意半徑R處的環(huán)向波數(shù)k（k=p/R），然后分別求出多層結(jié)構(gòu)中任意半徑處的相速度（cp=ω/k=ωR/p）和群速度（cg=?ω/?k），以獲得其頻散曲線。其中，相速度cp從環(huán)的內(nèi)徑到外徑線性增大。這種相速度隨半徑增大而增大的現(xiàn)象是在環(huán)面沿厚度方向保持相位不變的必要條件，且屬于波沿著曲面?zhèn)鞑r(shí)的獨(dú)特現(xiàn)象。

通過周向SAFE計(jì)算，對(duì)環(huán)面周向?qū)Рǖ哪B(tài)特性進(jìn)行研究，旨在識(shí)別可用于篩選的適當(dāng)模式-頻率組合，并預(yù)測(cè)其對(duì)多層復(fù)合環(huán)面不同區(qū)域缺陷的檢測(cè)靈敏度。為了選擇感興趣的解決方案，在網(wǎng)格的每個(gè)節(jié)點(diǎn)位置使用Pξi=-vjτji t計(jì)算模式的時(shí)間平均功率流。式中，Pξi為坡印亭（Poynting）矢量在ξi方向上分量為一段時(shí)間內(nèi)的時(shí)間平均值；vi=-iωUi為速度分量；τji為應(yīng)力張量；下標(biāo)i，j=1，2，3分別表示r，θ，z方向。

2 周向?qū)Р⊿AFE計(jì)算模型與驗(yàn)證

本節(jié)基于1.2節(jié)周向SAFE方法求解單層各向同性圓管中周向?qū)Рǖ念l散曲線，計(jì)算結(jié)果與公開文獻(xiàn)中所報(bào)道的頻散曲線結(jié)果比較，以驗(yàn)證周向SAFE方法及其計(jì)算執(zhí)行的準(zhǔn)確性。所用驗(yàn)證算例中各向同性圓管的幾何與材料參數(shù)見表1。用于計(jì)算環(huán)向波數(shù)k的半徑對(duì)應(yīng)空心圓管的外徑b；圓管的內(nèi)外徑比定義為η=IR/OR=a/b。

表1 公開文獻(xiàn)報(bào)道算例中空心圓管的幾何與材料參數(shù)

首先將周向SAFE方法計(jì)算得到的η=0.1鋁管中CLT模式的頻散曲線與Towfighi等論文［15］中報(bào)道的頻散曲線進(jìn)行比較，結(jié)果如圖3所示。其中，量綱一化波數(shù)為，量綱一化角頻率為式中，cT為各向同性介質(zhì)中傳播的體橫波波速；μ為介質(zhì)的拉梅常數(shù)。由圖3中兩種方法計(jì)算結(jié)果的比較可知，周向SAFE方法計(jì)算得到的頻散曲線與文獻(xiàn)中解析方法計(jì)算得到的頻散曲線結(jié)果高度吻合，在所研究的頻率范圍內(nèi)各CLT模式波數(shù)對(duì)應(yīng)的最大相對(duì)誤差均小于1.7%。

圖3 內(nèi)外徑比為0.1鋁制空心圓管中CLT導(dǎo)波模式的量綱一角頻率與量綱一波數(shù)的頻散曲線Fig.3 Dispersion curves of relation between dimensionless frequency and dimensionless wavenumber for an aluminum circular tube with a ratio of inner diameter to outer diameter of 0.1

類似地，本節(jié)采用周向SAFE方法進(jìn)一步計(jì)算了內(nèi)外徑比η=0.5的鐵制空心圓管中CLT和CSH模式導(dǎo)波的頻散曲線，并與Rose等［12］在公開文獻(xiàn)中報(bào)道的對(duì)應(yīng)kd-ω頻散曲線進(jìn)行比較，結(jié)果如圖4所示。CSH及CLT導(dǎo)波模式的數(shù)值與解析計(jì)算結(jié)果高

圖4 內(nèi)外徑比為0.5鋼制空心圓管中周向?qū)РＪ降慕穷l率與波數(shù)-厚度積的頻散曲線Fig.4 Dispersion curves of relation between angular frequency and wavenumber-thickness product for a steel circular tube with a ratio of inner diameter to outer diameter of 0.5

Tab.1 Geometric and material properties for circu

lar tubes reported in the literature度吻合，最大相對(duì)誤差為2.4 %，有效驗(yàn)證了周向SAFE計(jì)算方法的準(zhǔn)確性。相較于平板中的蘭姆波，環(huán)面中CLT波的角頻率-波數(shù)頻散曲線圖原點(diǎn)處出現(xiàn)新的導(dǎo)波模式，即CLT0模式，如圖4 b所示。該模式始于原點(diǎn)，然后迅速到達(dá)一個(gè)峰值，之后降低并與kd軸相交于點(diǎn)（kd=1-η），這與Rose公開文獻(xiàn)［12］中觀察到CLT波的特征一致。CLT模式的產(chǎn)生是由于環(huán)的內(nèi)表面不能支撐表面波，而外表面可支撐表面波傳播。因此，任何在環(huán)面內(nèi)表面產(chǎn)生的表面波都會(huì)迅速泄漏到環(huán)內(nèi)，然后反射到外表面，從而形成新的導(dǎo)波模式。在平板中，由于板的兩個(gè)表面都能支持表面波傳播，因而不存在該種模式［10］。當(dāng)環(huán)的內(nèi)外徑比值逐漸減小時(shí)，類似于S0板波的環(huán)面CLT模式將與其逐漸產(chǎn)生顯著差異，該現(xiàn)象說明在具有較小內(nèi)外徑比值的環(huán)面上不能支持對(duì)稱的波位移剖面。

3 工業(yè)用多層復(fù)合圓管結(jié)構(gòu)中的周向?qū)Рl散特性

相較于單層各向同性圓管，工業(yè)中常用多層復(fù)合圓管結(jié)構(gòu)，其在力學(xué)性能方面具有更大的裁剪自由度，但在材料組分與結(jié)構(gòu)構(gòu)型方面亦更為復(fù)雜；因此對(duì)該類多層復(fù)合圓管的加工制造、質(zhì)量監(jiān)控及服役性能評(píng)估任務(wù)均更為艱巨。本文選取了工業(yè)中兩種典型復(fù)合圓管構(gòu)件，即鋼-碳纖維雙層復(fù)合圓管和多層碳纖維樹脂基復(fù)合材料圓管。其中，鋼-碳纖維復(fù)合圓管常作為III型復(fù)合材料壓力容器以用于儲(chǔ)存氫氣［17］；因其具有重量輕、可提高壓縮氣體儲(chǔ)存效率等優(yōu)點(diǎn)，正逐漸代替I型全金屬壓力容器而應(yīng)用。而多層碳纖維樹脂基復(fù)合材料圓管構(gòu)件具有強(qiáng)度高，壽命長(zhǎng)、耐腐蝕，質(zhì)量輕、密度低等優(yōu)點(diǎn)，已被用作無人機(jī)、機(jī)器人、衛(wèi)星支承等工業(yè)設(shè)備中的關(guān)鍵承載結(jié)構(gòu)件。沿該類多層復(fù)合圓管中傳播的周向?qū)Р▽y帶多層環(huán)面中的材料/結(jié)構(gòu)/缺陷信息；深入研究并厘清其中周向?qū)Рǖ念l散特性則可望形成高效、準(zhǔn)確的無損檢測(cè)與評(píng)估技術(shù)。在實(shí)際檢測(cè)過程中，通常易測(cè)量導(dǎo)波模式沿空心圓管表面?zhèn)鞑サ南嗨俣然蛉核俣?，因此本?jié)選取相速度與頻率的關(guān)系表征周向?qū)РＪ降念l散曲線。

基于第2節(jié)驗(yàn)證的周向SAFE計(jì)算方法，進(jìn)一步建模計(jì)算鋼-碳纖維雙層復(fù)合圓管中的周向?qū)Р?。該雙層復(fù)合圓管由各向同性鋼制環(huán)面與各向異性碳纖維復(fù)合材料環(huán)面疊壓而成。鋼-碳纖維雙層復(fù)合圓管的幾何與材料參數(shù)：鋼層的內(nèi)徑為30 mm，外徑為32 mm，密度ρ=7 850 kg·m-3，體縱波波速cL=5 960 m·s-1，體橫波波速cT=3 230 m·s-1；碳纖維層的內(nèi)徑為32 mm，外徑為36 mm，密度ρ=1 600 kg·m-3。纖維方向沿圓周θ方向，表2給出表征碳纖維復(fù)合材料層正交各向異性對(duì)稱性所需的9個(gè)獨(dú)立彈性常數(shù)，然后可依據(jù)正交各向異性彈性常數(shù)轉(zhuǎn)換得到碳纖維層的剛度矩陣。值得注意的是，復(fù)合圓管的外層碳纖維樹脂基環(huán)面介質(zhì)采用均勻（homogeneous）假設(shè)，認(rèn)為其正交各向異性剛度矩陣在該層內(nèi)各個(gè)位置是保持一致的。

表2 碳纖維正交各向異性復(fù)合材料的彈性常數(shù)Tab.2 Elastic properties for carbon fiber orthotropic composites

通過周向SAFE模型計(jì)算其特征值波束，圖5則給出該雙層復(fù)合環(huán)面中周向?qū)Рǖ南嗨俣阮l散曲線，包含所研究頻率范圍內(nèi)所有CLT模式與CSH模式。其中環(huán)向波數(shù)及相速度對(duì)應(yīng)的半徑選取該環(huán)面內(nèi)層的外徑值，圖中虛線表示周向水平剪切（CSH）波的頻散曲線，實(shí)線表示周向類蘭姆（CLT）波的頻散曲線。不同于平板中SH0模式，環(huán)面中CSH0模式出現(xiàn)頻散現(xiàn)象。其他高階CSH模式亦發(fā)生頻散，隨著頻率的增加頻散降低，均存在截止頻率。CLT1模式隨著頻率增加頻散逐漸降低；其他高階CLT模式亦均存在截止頻率。針對(duì)頻散較低的周向?qū)РＪ?頻率組合，可激勵(lì)對(duì)應(yīng)導(dǎo)波模式，通過缺陷反射信號(hào)到達(dá)時(shí)間以定位環(huán)面缺陷，并通過導(dǎo)波透射或反射系數(shù)推測(cè)缺陷尺寸。針對(duì)頻散較高的周向?qū)РＪ?頻率組合，可利用頻厚積微小變化對(duì)應(yīng)顯著速度差的關(guān)系，結(jié)合導(dǎo)波層析成像反演方法［32］，可定量表征環(huán)面上的腐蝕減薄特征或應(yīng)用于其環(huán)面材料表征［33］。

圖5 鋼-碳纖維雙層復(fù)合環(huán)面中周向?qū)Рǖ南嗨俣阮l散曲線Fig.5 Phase velocity dispersion curves of circumferential guided waves in steel-CFRP double-layered composite circular tube

通過周向SAFE計(jì)算進(jìn)一步提取特征波數(shù)對(duì)應(yīng)的特征向量，可獲取對(duì)應(yīng)周向?qū)Рǖ膫鞑ツB(tài)。典型地，圖6分別顯示100 kHz頻率下鋼-碳纖維雙層復(fù)合環(huán)面中CLT1模式、CSH0模式和CLT2模式的波結(jié)構(gòu)及坡印亭矢量分布情形，其中顏色幅值表示導(dǎo)波模式能流密度的相對(duì)幅值（低至高：黑色至白色）；箭頭表示導(dǎo)波模式的粒子位移（幅值與方向）；振動(dòng)模態(tài)沿徑向壁厚提取r，θ，z方向的位移分量與各應(yīng)力分量（下同）。由結(jié)果可知，不同周向?qū)РＪ降穆晥?chǎng)能量匯聚于鋼-碳纖維復(fù)合圓管的不同區(qū)域，可依據(jù)能量分布以預(yù)測(cè)周向?qū)РＪ綄?duì)鋼-碳纖維雙層復(fù)合環(huán)面中不同區(qū)域的檢測(cè)能力。例如，其中CSH0模式的能量主要集中于內(nèi)環(huán)（鋼層），可能對(duì)界面狀態(tài)改變不靈敏；而CLT1模式和CLT2模式部分能量局域于層間界面處，對(duì)界面處的分層脫粘缺陷將呈現(xiàn)較高的檢測(cè)靈敏度。同時(shí)，沿雙層環(huán)面的徑向壁厚提取了各應(yīng)力分量的分布情形，亦可為各導(dǎo)波模式的缺陷檢測(cè)靈敏度提供參考依據(jù)。

圖6 鋼-碳纖維雙層復(fù)合圓管中周向?qū)РＪ皆?00 kHz頻率下典型波結(jié)構(gòu)及坡印亭矢量分布示意圖Fig.6 Diagrams of typical wave structure for varied circumferential guided wave modes in a steel-CFRP compound circular annuli,calculated by the SAFE method at 100 kHz

為觀察不同模式在r，θ，z方向上的位移分量占比情形，可沿圓管徑向厚度（即虛線截取線）提取聲場(chǎng)位移分量，可知CLT1模式的徑向位移ur占主導(dǎo)，內(nèi)層（鋼層）徑向位移ur振幅略高于外層（碳纖維層）徑向位移幅值，周向位移uθ和軸向位移uz幾乎為0，這與平板中反對(duì)稱式（A）導(dǎo)波模式的模態(tài)類似。而CSH0模式的軸向剪切位移uz占主導(dǎo)，內(nèi)層（鋼層）軸向位移uz幅值較低于外層（碳纖維復(fù)合層）軸向位移幅值，其徑向位移ur和周向位移uθ幾乎為0，這與平板中水平剪切型（SH）導(dǎo)波模式的振動(dòng)模態(tài)相似。CLT2模式的周向位移uθ占主導(dǎo)，其軸向位移uz為0，模態(tài)與平板中對(duì)稱式（S）導(dǎo)波模式相似。

第二類經(jīng)典環(huán)面結(jié)構(gòu)選取12層T800/924碳纖維增強(qiáng)樹脂基復(fù)合材料圓管結(jié)構(gòu)，其纖維鋪層順序?yàn)椋郏?45°/-45°）3］s，其中s表示為纖維鋪層對(duì)稱，每層厚度為1 mm，疊壓成形后的環(huán)面內(nèi)徑和外徑分別為138 mm和150 mm，密度ρ=1 500 kg·m-3。通過周向SAFE方法可計(jì)算該多層復(fù)合環(huán)面中周向?qū)Рǖ念l散曲線，這里注意模型采用非均勻（heterogeneous）假設(shè)，即各層定義隨纖維取向改變的單方向預(yù)浸料的剛度矩陣，其中纖維沿主軸0°方向定義剛度矩陣如表3所示。

表3 單方向碳纖維增強(qiáng)樹脂基T800/924復(fù)合材料的剛度系數(shù)Tab.3 Elastic constants of unidirectional T800/924 CFRP composites GPa

借助于商用有限元平臺(tái)COMSOL Multiphysics，基于所提出的周向SAFE方法獲得T800/924多層復(fù)合材料環(huán)面中周向?qū)Рǖ念l散曲線，如圖7所示，其中環(huán)向波數(shù)、相速度的半徑選取第四層環(huán)面的外徑值。圖7包含所研究頻率范圍內(nèi)所有CLT模式與CSH模式。環(huán)面中CSH波的頻散特性及振動(dòng)模態(tài)均與板中SH波相似。這里CSH0模式與SH0板波類似，均具有非頻散特性。工業(yè)中非頻散的導(dǎo)波模式能夠克服導(dǎo)波頻散特性給信號(hào)處理帶來的困難，最易直接應(yīng)用于無損檢測(cè)。高階CSH模式均發(fā)生頻散，且隨著頻率增加，其頻散程度降低，均存在對(duì)應(yīng)的截止頻率。所研究頻率范圍內(nèi)，各CLT模式頻散曲線亦給出，其頻散曲線并無相交。該現(xiàn)象在多層平板中亦然［34］。

圖7 T800/924多層復(fù)合環(huán)面中周向?qū)Рǖ南嗨俣阮l散曲線（虛線表示周向水平剪切（CSH）波，實(shí)線表示周向類蘭姆（CLT）波）Fig.7 Phase velocity dispersion curves of circumferential guided waves in T800/924 multilayered composite circular tube(the dashed line representing the circumferential horizontal shear(CSH)wave,and the solid line representing the circumferential lamb-like(CLT)wave)

平板中水平剪切（SH）模式位移與蘭姆波波場(chǎng)正交，其振動(dòng)模態(tài)沿板厚呈現(xiàn)水平剪切特征。類似于平板中的SH模式，周向水平剪切（CSH）模式位移與周向類蘭姆（CLT）模式波場(chǎng)正交，其振動(dòng)模態(tài)沿環(huán)面厚度呈現(xiàn)水平剪切特征，其粒子位移振動(dòng)方向（沿z軸）與波傳播θ方向垂直，且在r軸和θ軸方向的位移分量均為零。圖8顯示了在300 kHz下由周向SAFE方法計(jì)算得到的一系列高階CSH模式的波結(jié)構(gòu)及坡印亭矢量分布情況。CSH模式的波能量集中于多層復(fù)合材料環(huán)面中的不同區(qū)域及層間界面且分布連續(xù)，其沿環(huán)面?zhèn)鞑ミ^程沿厚度方向傳遞剪切應(yīng)力，可望對(duì)各界面處的分層脫粘缺陷具有較高的檢測(cè)靈敏度。此外，CSHn導(dǎo)波模式的順序根據(jù)其從低頻到高頻出現(xiàn)在頻散曲線中的順序進(jìn)行索引，亦與圓管厚度方向上水平剪切振動(dòng)的節(jié)點(diǎn)數(shù)量對(duì)應(yīng)。

圖8 T800/924多層復(fù)合圓管中周向?qū)РＪ皆?00kHz頻率下典型的高階CSH模式的波結(jié)構(gòu)及坡印亭矢量分布Fig.8 Diagrams of typical wave structure for varied circumferential guided wave modes in a T800/924 multilayered composite annuli,calculated by the SAFE method at 300 kHz

4 結(jié)論

本文提出了一種基于商用有限元平臺(tái)COMSOL Multiphysics的半解析有限元（SAFE）執(zhí)行方法，可高效準(zhǔn)確地計(jì)算多層復(fù)合環(huán)面中周向?qū)Рǖ念l散特性。文中首先從各向異性彈性波導(dǎo)的平衡方程和幾何關(guān)系出發(fā)，結(jié)合本構(gòu)方程，推導(dǎo)了柱坐標(biāo)系下圓管結(jié)構(gòu)中周向?qū)Рǖ目刂莆⒎址匠?；基于周向SAFE法及特征值問題數(shù)學(xué)框架，借助COMSOL Multiphysics中弱形式系數(shù)偏微分方程模塊，實(shí)現(xiàn)了對(duì)任意圓管中周向?qū)РＪ筋l散關(guān)系的精確求解，并交互驗(yàn)證該數(shù)值計(jì)算方法的準(zhǔn)確性。在此基礎(chǔ)上，本文將所提出的周向SAFE方法應(yīng)用于工業(yè)中兩種典型復(fù)合圓管結(jié)構(gòu)，即金屬-碳纖維雙層環(huán)面和多層碳纖維樹脂基復(fù)合材料環(huán)面，有效計(jì)算獲取其中傳播的所有周向?qū)РＪ降念l散曲線，并評(píng)估了不同周向?qū)РＪ綉?yīng)用于環(huán)面缺陷檢測(cè)與材料表征的潛力。

本文所提出的周向SAFE方法可計(jì)算任意材料屬性、任意曲率及任意層數(shù)的復(fù)合圓管中周向?qū)Рǖ念l散特性，并可推廣至具有任意環(huán)向截面形狀的復(fù)雜結(jié)構(gòu)中的導(dǎo)波問題求解。該方法的提出為獲取復(fù)雜材料與結(jié)構(gòu)中周向?qū)Рǖ念l散曲線提供了高效且便捷的手段。在此工作的基礎(chǔ)上，可以進(jìn)一步開展周向?qū)Рz測(cè)技術(shù)研究，并結(jié)合超聲換能器/傳感器布置策略構(gòu)建自診斷的智能材料與結(jié)構(gòu)，實(shí)現(xiàn)對(duì)該類工程多層圓管結(jié)構(gòu)的質(zhì)量控制與健康狀態(tài)監(jiān)測(cè)。

作者貢獻(xiàn)聲明：

余旭東：提出研究選題、設(shè)計(jì)研究方案、修訂論文、獲取研究經(jīng)費(fèi)、指導(dǎo)性支持。

秦榮：實(shí)施研究過程、理論推導(dǎo)與計(jì)算、統(tǒng)計(jì)分析、起草論文。

沈海：實(shí)施研究過程、數(shù)值計(jì)算、調(diào)研整理文獻(xiàn)。

左鵬：設(shè)計(jì)研究方案、修訂論文、指導(dǎo)性支持。

邵照宇：技術(shù)與材料支持。