覆冰CFRP板中水平剪切波的傳播特性

代重陽，陳智軍，王春濤，徐 君，鐘雪燕

(1.南京航空航天大學 自動化學院，南京 211106；2.中國空氣動力研究與發展中心 結冰與防除冰重點實驗室，綿陽 621000；3.南京鐵道職業技術學院，南京 210031)

飛機結冰會對飛機的氣動特性、操穩特性和起降特性產生嚴重影響[1]，因此結冰探測對于保障飛機安全飛行極為重要[2]。超聲導波憑借其可異面檢測、探測距離遠和裝置體積小等優勢，成為飛機結冰無損檢測技術中的佼佼者。目前采用超聲導波進行結冰探測的研究主要針對的仍然是鋁板結構[3-5]，而“一代材料，一代飛機”，飛機機體材料已經從以鋁合金為主跨入了以復合材料為主的新階段[6]，飛機結構復合材料化是世界航空裝備發展的趨勢[7]。碳纖維增強復合材料(CFRP)具有比強度與比模量高、抗疲勞與耐腐蝕性能好和可設計性強等優點，在航空領域得到廣泛應用[8]。飛機機體表面傳播的超聲導波屬于板波，根據偏振方向的不同可分為蘭姆波和水平剪切波。水平剪切(SH)波不存在沿板厚方向的位移，所以在飛機結冰探測中采用SH波可以消除水和防冰液等液體造成的干擾[9]。

超聲導波在波導結構中傳播時具有多模態和頻散的特點，文章利用數值計算技術，著重研究了SH波在覆冰CFRP板中的傳播特性，為在此基礎上進行的飛機結冰探測提供理論依據。采用計算速度更快的轉移矩陣法建立了波導結構模型，繪制了相速度與群速度頻散曲線，分析了結冰對CFRP板中SH波傳播特性的影響。

1 材料剛度矩陣獲取

CFRP板是一種典型的各向異性材料，其對載荷的響應與各向同性的鋁板不同，導波傳播特性也更為復雜[10]。單層CFRP板由基體中加入單向碳纖維組成(見圖1，圖中坐標軸x1為超聲導波的傳播方向，x2為水平剪切方向，x3為板厚方向，xf為碳纖維方向)。由于單向碳纖維填充排列的隨機性，單層CFRP板可視作正交各向異性材料，即垂直于xf的任何方向上的材料性質幾乎相同。

圖1 單層CFRP板結構示意

描述單層CFRP板的性能時，常采用的工程常數有拉伸彈性模量E1、E2、E3，剪切彈性模量G12、G23、G31，縱向泊松比υ12、υ21、υ13、υ31，橫向泊松比υ23、υ32，其中，E1、E2、G12、υ12和υ23是可由試驗測得的5個獨立量，在此基礎上可以分別計算出υ21和G23，又由于正交各向異性材料具有對稱性，可得E3=E2、G31=G12、υ13=υ12、υ32=υ23、υ31=υ21。

采用Einstein求和約定的應變-應力關系

εij=Sijklσkl

(1)

式中：εij與σkl分別為應變與應力矩陣：Sijkl為柔度矩陣，下標變量i、j、k、l的取值范圍均為1，2，3。

結合工程常數可以得到

(2)

式中：ε1，ε2，ε3為正應變；ε4，ε5，ε6為剪應變；σ1，σ2，σ3為正應力；σ4，σ5，σ6為剪應力。

由剛度矩陣與柔度矩陣的關系，可知單層CFRP板的材料剛度矩陣C為

C=S-1=

(3)

式(3)所示的形式僅為碳纖維方向xf與坐標系x1方向重合時的情況。但在實際應用中，總是將多個單層CFRP板通過有規律的鋪層形成層合板結構[11]，此時材料剛度矩陣將不滿足式(3)所示0°鋪層情況下的形式。在應力和應變分析中加入偏軸角度θ，可以推導出由0°鋪層到θ鋪層的材料剛度矩陣變換方法，即

Copqr=βmiβnjβokβplCijkl(o,p,q,r,i,j,k,l=1,2,3)

(4)

式中：Copqr為θ鋪層的材料剛度矩陣；Cijkl為0°鋪層的材料剛度矩陣；β為轉換矩陣。

(5)

CFRP板引入的材料各向異性問題大大增加了波導建模的復雜度，文章將材料的各向異性體現在剛度矩陣中，并提出了由CFRP板的獨立工程常數獲取材料剛度矩陣(任意角度鋪層)的方法。

2 層狀波導模型建立

從整體分析的角度出發，飛機機體表面結冰后形成了“復合材料板+冰層”的雙層結構波導。由第1部分可知，若對復合材料板進行逐層分析，則建立起的波導模型更加準確。因此，文章建立的波導模型基于層狀結構，組成部分中既存在正交各向異性的單層CFRP板，又包含各向同性的冰層。目前，針對層狀結構的建模方法主要有有限元分析法、全局矩陣法和轉移矩陣法等，其中有限元分析法在傳播距離較長時網格劃分的數量和全局矩陣法在層數較多時特征方程組系數矩陣的階數都極為龐大，所以文章提出采用轉移矩陣法建立層狀波導模型。鑒于0°和90°鋪層方式的廣泛應用，著重研究了SH波在該波導中的建模方法。

根據平面簡諧波的假設有

ui=Aieik(x1+αx3-cpt)

(6)

式中：ui為位移矢量的分量；Ai為振幅；k為波數；x1和x3為對應方向上的位置；α為x3方向上的衰減系數；cp為相速度；t為時間。

由牛頓第二定律和應力-應變、應變-位移的關系可推導出每一層結構中的波動方程

(7)

式中：ρ為該層材料的密度。

將位移分量代入波動方程，可得到Christoffel方程

KA=0

(8)

式中：矩陣A為[A1A2A3]T；矩陣K為3×3的系數矩陣。

各分量的具體表達式為

(9)

結合第1部分可知，冰層與0°/90°鋪層CFRP板材料剛度矩陣中的非零項都分布在相同的位置處，所以式(9)中的K12、K21、K23和K32都為0。因此，可以將式(8)分解為

(10)

(11)

式(10)和(11)表明，在由特殊角度(0°和90°)鋪層的CFRP板與各向同性的冰層所組成的波導介質中，偏振方向在x1-x3平面內的蘭姆波和偏振方向在x2方向上的SH波是解耦的，這種解耦現象為單獨研究SH波的傳播特性創造了條件。同時也揭示了蘭姆波和SH波的解耦不僅會發生在各向同性的介質中，也會發生在一些特殊的各向異性材料中。

由式(11)可知，當振幅A2具有非零解時，K22的值為零。K22是關于α的二次多項式，依據部分波理論，每一層結構可得到兩個SH型部分波的解。在確定了α的值后，可將每一層結構中x2方向上的位移分量大小u2描述為兩個部分波的線性組合

exp(ik(x1-cpt))

(12)

式中：αm為K22等于零的條件下得到的x3方向的衰減系數；A2m為對應的振幅；Bm為兩個部分波對應的加權系數。

式(12)中的加權系數Bm需要結合邊界條件進行確定，邊界條件涉及位移場和應力場，將每一層結構中x3-x2平面內的應力分量σ4表示成部分波的解

(ik)exp(ik(x1-cpt))

(13)

在建立層狀結構波導模型時，需要對聲場分量進行構造，以便于后續矩陣計算時的消元。在忽略公共簡諧項exp[ik(x1-cpt)]的前提下，用Γ矩陣來表述聲場矢量，由式(12)和式(13)可以得到

(14)

將Γ矩陣進一步分解為三部分

Γ=XWB

(15)

式中：X為2×2的矩陣，是波數k和相速度cp的函數；W也為2×2的矩陣，是波數k、相速度cp和x3方向上位置的函數；B為2×1的矩陣，是加權系數Bm的函數。

通過轉移矩陣法建立N層CFRP板+冰層結構波導模型的方法如圖2所示[n為多層板中的任意一層，N為多層板中的最后一層；上標“+”表示上表面，“-”表示下表面；上標“(ice)”表示冰層]。

圖2 轉移矩陣法建立N層CFRP板+冰層結構波導模型示意

在該N+1層結構中建立N+1個局域坐標系，則每層結構底部的x3值為0，這將使得W矩陣成為單位矩陣，從而可簡化聲場分量為

Γ(n)-=X(n)B(n)

(16)

同時，在每層結構的頂部，聲場分量又可以表述為

Γ(n)+=X(n)(W(n)|h(n))B(n)

(17)

式中：h(n)為每層結構的厚度。

由式(16)和式(17)可知，每一層上表面的聲場可由下表面的聲場表示，即

Γ(n)+=Ψ(n)Γ(n)-

(18)

式中：Ψ(n)矩陣為第n層的轉移矩陣。

由于相鄰兩層結構的連接處，界面位移和應力滿足連續性條件，故可以得到

Γ(n)+=Γ(n+1)-

(19)

將式(19)代入式(18)可進一步得到

Γ(n+1)-=Ψ(n)Γ(n)-

(20)

將式(20)所示的關系應用于N+1層結構內部的N個連接處，有

Γ(ice)+=Ψ(ice)Γ(n)+=Ψ(ice)Ψ(n)…Ψ(1)Γ(1)-=

(21)

式(21)體現了轉移矩陣法的優勢：總轉移矩陣Ψ的階數與層數無關，始終為2×2。

在此基礎上，再將N+1層結構上、下表面處，界面應力滿足自由邊界條件考慮進來，得到

(22)

拆分出式(22)中的一部分，有

(23)

式(23)也稱作SH波傳播的特征方程。由于u2(1)-有非零解，當波數k和相速度cp滿足Ψ21等于零時，就能得到頻散曲線。

相比于計算繁瑣的全局矩陣法[12]，文章所采用的轉移矩陣法可通過矩陣相乘的方式消除內部的連續性條件，只保留上、下表面的邊界條件，特征方程的階數與層狀結構的層數無關，顯著提高了計算速度。

3 導波傳播特性分析

超聲導波多模態和頻散的特點可以很直觀地體現在頻散曲線上，因此通過建立“CFRP板+冰層”波導模型，結合繪制的頻散曲線對SH波的傳播特性進行分析。為了提高模型的實用價值，以IM7/977-3(被廣泛應用于飛機機身蒙皮)CFRP材料為例，其ρ=1.608 g/cm3，E1=172 GPa，E2=9.8 GPa，G12=6.1 GPa，υ12=0.37，υ23=0.55[13]。

模型中的CFRP層合板為8層結構，每個單層板厚度為0.2 mm，采用一種典型的鋪層方式，鋪層代碼為[(0/90)s]2。由兩種不同鋪層角度的8層CFRP板和冰層組成的層狀結構如圖3所示。

圖3 8層CFRP板+冰層結構

飛行過程中，飛機表面形成的冰層主要有明冰、霜冰和混合冰，其力學參數如表1所示[14]。

表1 3類冰層的力學參數

3類冰層中，以明冰對飛行的危害最為嚴重，所以首先分析明冰對CFRP板中SH波傳播特性的影響。CFRP板表面形成明冰前后對應的頻散曲線如圖4所示(f為頻率,cg為群速度)，其中表面明冰層厚度為1 mm。

圖4 CFRP板表面有無冰層時SH波的頻散曲線

由圖4可知：① 在相同的頻率和速度范圍內，CFRP板表面結冰與否所對應的頻散曲線數目是不同的；② 在相速度頻散曲線中，高階SH模態導波的截止頻率會受冰層的影響而顯著減小；③ 當CFRP板表面無冰時，其SH0模態導波的相速度與群速度都為常數，但結冰會改變這一現象，并且群速度的變化更為明顯；④ 結冰所引起的SH0模態導波變化主要發生在頻率較高處，而在低頻處仍表現為常數，不過該常數值的大小會發生微小改變。

對于冰層對CFRP板中SH波傳播特性的影響，不僅可進行定性分析，也可進行定量研究。圖5為CFRP板表面明冰層厚度分別為0.25，0.50，0.75 mm時的頻散曲線。

圖5 CFRP板表面冰層厚度不同時SH波的頻散曲線

由圖5可知：① 在明冰層厚度由0.25 mm開始均勻遞增至0.75 mm的過程中，相同區間內頻散曲線數目增加了1條；② 以SH1模態導波為例，其相速度截止頻率由464.1 kHz經424.7 kHz減小至384.7 kHz，即明冰層厚度每增加0.1 mm會引起約16 kHz的截止頻率變化量；③ SH0模態導波的相速度和群速度都會隨著明冰層厚度的增加而減小，例如在頻率為1.5 MHz處，3種明冰層厚度條件下的群速度相較于未結冰時分別減小了約17，55，92 m/s；④ 在低頻的50 kHz附近，SH0模態導波的速度都保持不變且相速度與群速度相等，但不同明冰層厚度對應的速度分別約為1.940，1.934，1.928 km/s。

除了冰層厚度，CFRP板表面冰層類型的差異也會引起傳播特性的變化。當CFRP板表面冰層厚度都為0.25 mm時，冰層類型為明冰、霜冰和混合冰所對應的頻散曲線如圖6所示。

圖6 CFRP板表面冰層類型不同時SH波的頻散曲線

由圖6可以看出，明冰與混合冰對SH波傳播特性的影響較為接近，僅在SH0模態高頻處以及高階模態處具有明顯差異，而霜冰的影響卻與二者不同，尤其反映在高階模態的相速度和SH0模態的群速度上，隨著模態的增高，霜冰對應的相速度截止頻率與明冰所對應截止頻率的差距逐漸擴大；對于SH0模態的群速度，霜冰產生的影響是一種截然不同的趨勢，在0.5～1.5 MHz頻率區間內，群速度劇烈減小。

這些反映在頻散曲線上的差異可以指導基于SH波的結冰探測研究，合理利用這些傳播特性的變化規律，能夠準確地對有無冰層、冰層厚度和冰層類型進行探測。

4 結語

針對以CFRP層合板為主的復合材料，提出了一種由工程常數獲得任意角度鋪層的單層CFRP板材料剛度矩陣的方法。將轉移矩陣法應用于建立SH波在CFRP板+冰層介質中的傳播模型，利用矩陣相乘的方式消除層狀結構內部的連續性邊界條件，使得計算時矩陣的階數與層狀結構的層數無關，顯著提高了計算速度。以材料IM7/977-3、鋪層方式為[(0/90)s]2的典型CFRP層合板作為示例，結合繪制的頻散曲線，分析了有無結冰、冰層厚度和冰層類型對SH波傳播特性的影響；歸納了這些待測量引起頻散曲線變化的規律，為基于SH波的航空復合材料結冰探測研究提供了理論基礎。