基于Duffing-Holmes系統的材料非線性的量化檢測方法

2018-12-21 07:12:00劉小峰韋代平

振動與沖擊 2018年24期

柏林，唐滔，劉小峰，韋代平

(重慶大學機械傳動國家重點實驗室，重慶 400044)

針對材料非線性檢測，有研究表明，非線性Lamb波的傳播會由于材料自身的非線性以及其結構微裂紋的局部接觸而激發非線性Lamb波的非線性變化[1-2]，因此可利用高次諧波信號對材料性能進行實時程度評價。大量試驗證明，通過測量Lamb波非線性效應可以有效反映材料結構的微觀變化，以此實現早期損傷或性能退化的檢測[3-4]。但目前的研究成果中，針對非線性Lamb波的信號的特征增強以及抗噪能力的研究很少。并且在試驗中，由于Lamb波自身的擴散角以及頻散特性等，致使接收信號極易受噪聲污染從而失真。為了提高超聲波檢測技術的精度，采用非線性信號傳統的降噪方式，如Kalman濾波[5]，LMS自適應濾波[6]，以及小波閾值降噪[7]等，都已經難以滿足對微弱信號的檢測要求。鑒于非線性Lamb波對于材料特性檢測的優勢性[8]但易受外界干擾，實現非線性Lamb波抗噪量化分析對對板構件的早期性能退化識別具有重要意義。

由于杜芬混沌振子對外加擾動具有初值敏感性及噪聲免疫性，在微弱信號的檢測中得到了廣泛的運用。但目前研究都主要停留在混沌系統的對微弱信號的定性分析上，而針對非線性Lamb波的抗噪定量分析并未涉及。運用杜芬振子對微弱非線性Lamb的定量研究也尚屬于發展階段。論文利用在不同幅值激勵下杜芬系統混沌程度的差異性，并通過Lyapunov指數對混沌特性進行量化，建立激勵幅值與系統輸出的線性模型用于實現對噪聲干擾下的微弱信號量化分析。

1 非線性Lamb波傳播原理

當Lamb波經過發生換能器輸入材料中，材料在超聲波作用下內部空間受到擾動，在內部微小裂紋以及材料的非完全彈性的共同作用下，激發出高階諧波信號。對于縱向傳遞的Lamb波，即u0=A1sin(ωt)，聯立一維波動方程

(1)

式中：x為縱向傳播距離;c為縱向傳播波速;β′即為Lamb波非線性系數。當只考慮波在傳遞過程中存在的頻散現象時，式(1)有唯一解為

u=A1sin(kx-ωt)+0.125β′(kA1)2xsin 2(kx-ωt)

(2)

式中:k為基波波數;ω為基波圓頻率;A1為基波幅值。可以看出對于二次諧波幅值A2，在波的傳遞過程中材料的非線性具有累積效應，并與非線性系數β′存在線性關系

β′=(8A2)/(kA1)2x

(3)

在檢測距離x恒定的情況下，重構非線性系數β

β=A2/(A1)2

(4)

在實際非線性Lamb波檢測過程中，二次諧波幅值與基波相比非常小，常常受到噪聲信號的干擾。因此，β指數抗噪能力不佳。

2 Duffing-Holmes系統與Lyapunov指數譜

為了實現杜芬系統在不改變系統參數的情況下，對任意頻率諧波檢測的普適性。本文采用改良后的Duffing-Holmes系統，利用待測頻率ω帶系統的時間尺度進行放縮處理，其三維自治狀態方程為

(5)

式中:k為阻尼比;x3-x5為系統非線性恢復力;Fcos(ωt)為系統內周期策動力;ω即為待測頻率;F為內策動力幅值。(s(t)=Acos(ωt+φ)+n(t))外部系統擾動，其中:Acos(ωt)為待檢信號;n(t)為待檢信號中的噪聲成分。

Duffing-Holmes系統相軌跡具有如下特性：當系統阻尼k固定，內策動力幅值F從0逐漸增加臨界點F0時，相軌跡狀態由周期內軌運動過渡到倍周期分叉最后進入混沌狀態。一旦F超過F0，系統動力學特性發生突變，整體進入大周期狀態，即發生正相變。針對式(5)，待測信號s(t)中存在與內策動力同頻諧波信號，雖然同頻諧波的合成并不影響頻率成分，但改變了系統的策動力幅值與相位，觸發臨界系統正相變，即系統由混沌狀態向大周期狀態跳變。另一方面，當輸入信號為噪聲信號或不含待測成分時，將不會對系統的周期策動力F造成穩態影響，Duffing系統基本免疫。

在混沌理論中，Lyapunov指數可以在一定范圍內定量描述非線性系統的動力學行為特征，并能在一定時間范圍內對不便于觀察的微弱信號進行幅值估計。其中，Lyapunov指數沿某一方向的正負與大小，分別表征系統在該方向上的相鄰軌跡的平均發散與收斂程度。當系統的最大Lyapunov指數為正時，Duffing系統在相空間中相鄰軌跡成指數發散，系統呈現混沌狀態。

對于式(5)的三維Duffing自治系統方程的Lyapunov指數，如式(2)，此時系統的Jacobian矩陣為

(6)

式中:s′(t)為系統外部擾動關于時間的一階導數;ω為待測頻率。在系統初始階段，令初始值x0=y0=t0=0，并將式(5)表示為3維連續動力學系統

(7)

此時，令Φ(t)為式(5)的基礎解空間，則該系統的線性變分方程可表示為

(8)

式中：J(t)為Duffing系統t時刻的3維Jacobian矩陣;I3為3×3的單位陣。利用4階Runge-Kutta法可分別計算變分方程在不同時刻的Jacobian矩陣。并對基礎解空間Y(t)進行QR分解，即Φ(t)=Q(t)R(t)。在文獻[9]中指出，系統的變分方程式(7)對應的Lyapunov指數滿足

(9)

式中：Rii(t)為上三角矩陣R(t)的正對角項。

3 基于諧波逼近的信號幅值估計

臨界Duffing-Holmes系統具有參數敏感性以及噪聲免疫性[10]。對于混沌振子系統，系統擾動s(t)可分為由系統頻率決定的慢變周期擾動，以及由高頻信號決定的快變周期擾動，已有研究表明慢變周期擾動能使不穩定的相軌跡進入穩定周期[11]，而快變信號僅能驅動系統的局部振蕩。對于線性放縮后的杜芬系統，Acos(ωt+φ)為慢變周期成分，噪聲成分n(t)相對于系統頻率ω仍屬于快變周期成分。鑒于Duffing系統動力學特征——Lyapunov指數——與待測慢變周期成分幅值A具有一定的相關性，利用此間關聯對噪聲下的信號成分進行幅值定量分析，具體步驟如下：

步驟1 根據Lamb波的采樣頻率fs及待測頻率ω，設置杜芬振子系統的縮放系數ω及分析步長dt=1/fs；

步驟2 根據系統分岔圖，確定混沌臨界閾值F0；

步驟3 由于接收到的Lamb波屬于瞬態信號，而杜芬系統一般用于穩定諧波檢測，因此需要對采集到的Lamb波進行周期延拓處理；

步驟4 為降低一次諧波對二次諧波檢測的影響，采用高通濾波方法濾除一次諧波，通過幅頻特性確定濾波器對二次諧波幅值影響近似忽略；

步驟5 將濾波后信號進行γ=0.4-1.6的等比放縮，并作為系統擾動信號輸入混沌系統。利用最小二乘擬合建立放縮系數γ與Lyapunov指數之間的線性方程模型

LE=a×γ+b

(10)

式中：γ即為參考信號幅值A0與待測幅值Ax之比。

步驟6 為減小噪聲帶來的局部振動的影響，通過小波軟閾值法提取噪聲信號，并與相同的采樣頻率下的參考諧波信號融合，其中，仿真諧波的幅值已知，即為式(10)中的參考幅值A0；

步驟7 測試參考信號在經過相同演化時間后Lyapunov指數，并結合式(10)，求解待測幅值Ax。

4 試驗測試

本章通過試驗驗證前文提出的方法的可行性。在諧波幅值定量分析試驗中，通過收集非線性Lamb波在非線性材料中傳播一定距離后的信號，測量其基波幅值與二次諧波幅值，以此獲得近似無噪的初始信號。

圖1 試驗設備Fig.1 Experimental equipments

試驗使用的非線性超聲測試系統為RAM-500 SNAP (RITEC Inc., Warwick, RI), 如圖1所示，具有檢測超聲衰減和波速的能力。相較于傳統線性檢測法，RAM-500具有更好的靈敏度。圖1中為本次試驗對象1.5×620×2 500 mm的無損鋁板，鋁板平行放置在海綿塊上，無固定約束。在基準線兩端分別布置壓電式發射換能器與PZT接收換能器。本次試驗距離700 mm，激發信號200 kHz，試驗周期數30，采樣頻率為10 MHz。由于換能器與鋁板之間存在空隙，試驗通過甘油作為超聲耦合劑排除間隙空氣，以便Lamb波的傳遞。

調整Duffing-Holmes系統內策動力幅值F進入臨界狀態，并設置式(5)中的系統參數k=0.5[12]，此時系統的臨界閾值F0=0.725[13]。分析試驗采集信號，提取出采集到的Lamb波包作為初始信號，并對原始信號做零相位帶通濾波后得到如圖2，分析二次諧波與基波的時域關系。

圖2 諧波時域對比圖Fig.2 Comparisons of original wavepacket and extracted harmonic

根據圖2，在有限時域采樣長度N內，由于頻散特性以及波在傳遞過程中的擴散，二次諧波在時域空間中嚴重發散。根據信號的單邊幅值譜定義

(11)

對于瞬態信號，其對應頻率的幅值由信號采樣長度N以及瞬態信號長度N1共同影響。但相比與基波長度N，二次諧波長度N1難以確定。為防止由于截斷時長帶來的幅值譜差異，對待測信號統一補零到相同時窗長度后計算幅值譜。取接收到的Lamb波包信號進行分析，如圖3(a)所示；為模擬試驗干擾，在接收信號中加入標準差為1×10-3的白噪聲，即圖3(b)。

(a)原始信號

(b)加噪信號圖3 原始信號與加噪信號對比圖Fig.3 Comparison of the received signal and noised signal

從圖3(a)可知，接收到的Lamb波包中二次諧波非常微弱，與一次諧波幅值相差2個數量級。在試驗中，即使存在微弱干擾，也會使得諧波幅值劇烈波動。在圖3(b)中二次諧波幅值依然清晰可見，但根據表1可知，此時二次諧波幅值明顯失真，傳統非線性指數β在干擾情況下已失效。

采用本文方法首先對Lamb波包進行整周期延拓濾波后，并乘以縮放系數γ=1后作為系統擾動參數s(t)輸入系統。通過4～5階變步長Runge-Kutta法求解，得出系統的相軌跡圖與Lyapunov指數譜如圖4(a)所示。

表1 不同算法條件下的β指數對比Tab.1 Comparison of calculated β with different algorithms

在圖4(b)中，當系統無擾動輸入時，混沌系統存在奇怪吸引子。對于三維Duffing-Holmes自治系統，LEx，LEy和LEt分別表示系統沿x軸，y軸，時間軸t方向的Lyapunov指數。

(a)γ=1時杜芬系統輸出

(b) 無輸入時杜芬系統輸出圖4 系統輸出響應對比Fig.4 Comparison of chaotic system outputs

此時沿某一方向相軌跡必發散，即圖中沿x軸方向LEx>0。在圖4(a)中，由于系統擾動項中存在穩定待測頻率成分，混沌系統出現普通吸引子——極限環，在Lyapunov指數譜中，有且僅存沿時間軸方向動力學指數為0，其余指數皆小于0。經過試驗分析，隨著待測成分幅值Ax增大，杜芬系統整體趨向穩定速度加快，在相同演化時間內與沿x軸方向Lyapunov指數呈現負相關趨勢。本文利用LEx分析系統動力學特性。

調整擾動信號的縮放系數γ∈[0.4,1.6]，在經歷相同的演化時間后，記錄對應的沿x軸方向的Lyapunov指數LEx，如表2所示。

表2 放縮系數γ與Lyapunov指數對應關系Tab.2 The relation between the coefficient of contraction γ and the Lyapunov exponent

圖5中，利用最小二乘法對得到的離散數據進行一次線性擬合，得到放縮系數γ與LE之間的線性模型(擬合確定度R-square: 0.968 1)

LE=-0.004 52γ-0.078 1

(12)

圖5 γ與LE的一階線性擬合Fig.5 Linear fitting between γ and LE

在整個時域空間內，混沌系統的整體趨勢主要受到待測頻率成分幅值影響，而噪聲成分僅造成局部的波動。因此，對于擾動信號的放縮，可以等價視為不同幅值的無噪諧波信號直接作用于混沌系統，此時系統狀態將表現出于待測信號縮放系數之間的直接相關性，即式(12)，系統狀態差異通過Lyapunov指數進行量化。由式(10)可知，當系統線性模型確定，為定量分析噪聲下的準確二次諧波幅值Ax，先引入參考信號(波數30，頻率0.4 MHz，采樣頻率10 MHz的正弦諧波信號)，此時仿真信號并沒有考慮波的擴散，將參考信號通過補零的方式延拓至與原始信號相同的采樣長度，分析此時幅值譜，如圖6所示。

圖6 同頻率基準信號Fig.6 Reference signal of the same frequency

圖7中，在相同的采樣長度下，參考信號幅值為5.508×10-5V。考慮噪聲帶來的局部振動，對原始信號進行小波軟閾值信噪分離提取噪聲成分，如圖7所示。

圖7 提取噪聲信號Fig.7 Extracted noise

將參考信號與噪聲成分結合延拓后作為系統擾動輸入杜芬系統，經過相同的演化時間后，沿X軸方向動態指數LEx=-0.083 8，如圖8所示。聯合線性方程式(10)與式(12)，代入參考幅值A0=5.508×10-5V。代入求解參考信號與真實信號之間的放縮系數γ=A0/Ax=1.26，即預測二次諧波幅值為Ax=4.371 4×10-5V。

圖8 參考信號的Lyapunov指數Fig.8 Lyapunov exponent of the reference signal

通過表1，當加入微量的噪聲，此時基波幅值基本維持不變，但二次諧波參數劇烈波動。通過小波閾值法可以看出，當噪聲信號與待測信號出現在同一頻率上時，小波并不能實現很好的分離。Kalman濾波器在信噪分離階段，分離出的信號并不是最佳估計，信號失真嚴重。LMS自適應濾波收斂速度較慢，信噪分離效果不佳。而此時混沌預測達到了較好的估計效果。

5 結論