基于廣義混合元加筋板的振動特性分析1)

2021-04-25 08:50:20陳秀濤卿光輝

力學與實踐 2021年2期

陳秀濤卿光輝

(中國民航大學航空工程學院，天津300300)

通常情況下，帶有加強筋的板殼結構具有剛度大、耗材少和重量輕等特點，因此被廣泛應用于各類實際工程中。但由于加筋結構的振動行為對整體結構的穩定性影響十分顯著，所以有關加筋結構的振動特性分析一直是國內外學者研究的熱點問題。

Mukherjee 等[1]和 Mukhopadhyay 等[2]對加筋結構振動分析的早期文獻進行了綜述：由于加筋結構的復雜性，基于位移有限元模型的數值分析方法最為普遍。Trkmen 等[3]提出了暴露于沖擊波的加筋板的位移有限元分析模型，并進行了實驗研究。Zhao等[4]使用能量方法，研究了簡單支撐的交叉加筋多層圓柱殼的自由振動。Rikards 等[5]發展了一種精度較高的三角形單元，研究了加筋復合材料層合殼的自由振動。Guo 等[6]對加筋層合板結構提出了分層有限元模型，并對這類結構進行了屈曲分析。游翔宇等[7]基于光滑有限元法對加筋板的自由振動進行了研究。

值得說明的是，一般位移有限元法是采用有限的節點代替了真實結構中無窮多的點，結果導致有限元模型比真實模型的剛度大，故而固有頻率結果大于真實解。從理論上講，由于合理地引入非協調位移項，使其近似多項式的階次完備，致使非協調元的柔性增加。系統控制方程的柔性也相應地增加，因而采用非協調元的模型可以得到精度更高的數值結果。然而，即使采用收斂性較好的非協調位移元分析加筋板的振動特性問題，其相應的有限元模型的剛度也是偏硬的，因而在有限元網格較稀疏的情況下，固有頻率結果通常是大于真實解的。

近些年來，一些研究人員應用Hamilton 正則方程半解析法分析加筋層合結構[8]或厚度不連續結構的振動特性問題[9]。Hamilton 正則方程半解析法分析板殼結構的顯著優勢是允許材料和幾何特性沿厚度方向改變，并可以方便地分析層合板殼結構類的振動問題。通常情況下，基于Hamilton 正則方程半解析法的模型，數值結果相對比較準確。但是，該方法有一個明顯的缺點，即由于層合板殼類結構的每層在面域內對位移變量和面外應力變量進行了有限元離散，因而每一層的單元組裝后導致規模較大的數值矩陣。在形成結構整體的控制方程前需要對每一層的數值矩陣分別進行關于厚度變量的指數矩陣運算，這一過程的數值運算成本高。另一方面，完成了每一層矩陣指數運算后，為保證不同材料層的層間位移變量和應力變量的連續性，又要對各層的控制方程進行乘積運算，這一計算步驟也需要一定的時間?？傊琀amilton 正則方程半解析法面對工程中龐大的結構問題，單元離散的數量大，數值運算相對于常規的有限元位移法來說，矩陣指數運算量大，時間消耗大。

傳統的混合有限元法有很多優點，例如，通常只需C0連續的插值多項式來近似表達待求的場量，位移和應力計算結果精度高。最近，Qing 等[10]結合最小勢能原理和H–R 變分原理建立了分析靜力學問題的非協調廣義混合元。文獻[11] 擴展了這種非協調廣義混合元的應用范疇。

本文在文獻 [10-11] 的基礎上，根據 Hamilton原理，首先建立考慮結構振動特性的廣義混合變分原理，然后建立了非協調廣義混合元的動力學模型，用于加筋板的振動特性分析。

1 位移元的振動特征方程

通常情況下，以位移元為基礎的無阻尼自由振動方程[12]為

式中，M表示結構的質量矩陣，為單元的一致質量矩陣，N表示形函數矩陣,ρ表示材料的密度，V表示連續體的體積，K表示結構的剛度矩陣，是單元的剛度矩陣，B為應變矩陣，D是材料的彈性矩陣。

設q(t) 為待求的解，考慮簡諧運動問題時，位移解的表達形式為

式中，Φ表示固有振型，ω表示對應固有振型的頻率。

將式 (2) 代入式 (1) 可得到無阻尼自由振動的特征方程為

對式(3) 的求解方法較多，例如，逆迭代法、子空間迭代法和行列式搜索法等。

2 廣義混合元的振動特征方程

根據廣義混合變分原理[13]，可設結構的混合能(包括應變能和余能) 為

式中，σ表示應力向量，C表示材料的剛度系數矩陣，?表示微分算子，u表示位移向量。

處于自由振動狀態下的結構的動能為

式中，ω表示固有頻率，ρ表示材料密度。

考慮動力學問題的Hamilton 變分原理為

因此含有固有頻率ω參數的廣義混合變分原理可表示為

式中，參數α的取值范圍為 0α1。

參考文獻[10-11] 的推導過程，由式(7) 可得到六面體非協調廣義混合單元的自由振動特征方程為

式中，h為單元柔度矩陣，g為單元杠桿矩陣，pe表示單元節點的應力參數，Φe表示單元固有振型參數，

各式中的形函數矩陣N和非協調項形函數矩陣Nr表達形式可參見文獻[10-11]。

由式(8) 可導出

將式 (9) 代入式 (10) 中有

其中，剛度矩陣κ=αk+(1?α)gTh?1g。

式(11)表明六面體非協調廣義混合元的剛度矩陣κ與柔度矩陣g相關，參數α的取值范圍為0～1, 并且不同的取值可得到不同的剛度矩陣。若式 (7) 中α= 1，則式 (11) 中的κ與位移元下的剛度矩陣一致，同為基于最小勢能原理的非協調位移元的剛度矩陣。α1 時，其不同的取值可以調節剛度矩陣柔性。文獻[14] 推導含參數α的廣義混合元模型時，小變形線性彈性問題是最基本的假設。另一方面，假設含參數α的廣義混合元模型的應變能與精確的應變能相等?；谶@一假設，文獻 [14] 導出了參數α= 0.75 的最優值。因此，就小變形線性彈性問題而言，廣義混合元模型中的α= 0.75 對計算結果的最主要影響是數值解更接近問題的真實解。

對式(11)求和可得整體有限元模型的振動特性方程為

關于非協調廣義混合元模型施加邊界條件的方法：最終的振動特性方程 (12) 中只含有位移變量。因此方程(12) 只需考慮位移邊界條件的引入，其引入方法與位移元法中的方法相同。

3 實例分析

以下數值分析均采用式(12)對算例進行數值分析，其中參數α=0.75。此外，在關于加筋板的算例中，為了提高計算效率，對式(12) 中的單元質量矩陣采用集中質量矩陣的形式。

3.1 層合板的自由振動特性分析

如圖 1 所示的三層板，每層均為正交異性材料，坐標軸沿彈性主方向，上、下兩層材料相同。令分別表示第一和第二層的C11值，ρ1和ρ2分別是第一層和第二層材料的密度。當δ= 1，γ= 1 時，該板即為一單層正交異性板，板長度a和寬度b，a=b，厚度h，h/a=0.1，h1/h=h3/h=0.1，h2/h=0.8。邊界條件如圖2 所示，為四邊簡支，即當x=0 和x=a時，u=w=0，當y=0 和y=b時，v=w=0。剛度系數比為

圖1 層合板尺寸圖

圖2 四邊簡支約束

以下通過改變網格密度來討論非協調廣義混合元的收斂性。

在γ= 1 和δ= 1 的模型下，在x方向和y方向劃分同等單元(從2 個單元遞增到12 個單元)，在厚度方向z劃分10 個單元，保持不變。

第一階固有頻率的收斂曲線如圖3 所示。從圖中可以明顯看出，隨著網格密度的增加，通過非協調位移元和非協調廣義混合元計算的頻率結果都逐漸趨于平穩。當網格密度為8×8×10 時，兩種方法的結果都逼近精確解，但非協調廣義混合元的結果更加精確。

圖3 不同網格密度下的1 階固有頻率

以下根據δ，γ的不同取值，分別求其一階固有頻率。網格劃分情況如圖 4 所示。兩種有限元模型的固有頻率結果見表 1，且兩種模型的 1 階振型圖一致，如圖5 所示。

比較表1 中的數據，不難看出，非協調位移元的結果高于文獻解。這是因為位移元模型的剛度相對于實際的模型偏硬。在網格模型相同的情況下，因為非協調廣義混合元模型的剛度更接近真實結構，所以本文的非協調廣義混合元的結果優于非協調位移元。這一結論與文獻[14] 中的相關論述相同。

圖4 層合板網格圖

√/C(2)11 )本文方法0.048 06 1.33 0.057 65 1.06 150.077 41相對誤差/% — 1.51 0.35

圖5 1 階振型圖

3.2 單加強筋偏心加筋板的自由振動特性分析

帶有單加強筋的偏心加筋板，尺寸圖如圖6 所示，材料參數為：彈性模量E= 68.7 GPa，泊松比ν= 0.29，密度ρ= 2823 kg/m3。幾何參數：a=b= 0.203 2 m，t= 0.001 371 6 m，w=0.006 35 m，h=0.012 7 m，邊界條件為四邊固支，如圖 7 所示，當x=0 和x=a時，u=v=w=0，當y=0 和y=a時，u=v=w=0。對該結構的振動特性進行分析。

圖6 單加強筋偏心加筋板尺寸圖

將板部分劃分成 27×27×1 個六面體單元，加強筋部分劃分成1×27×8 個六面體單元，網格劃分情況如圖8 所示。兩種有限元法的前4 階固有頻率結果見表2。前4 階固有振型見圖9，其中非協調位移元的結果由有限元商業軟件Abaqus 得到。

圖7 四邊固支約束

圖8 單加強筋偏心加筋板網格劃分圖

表2 單加強筋偏心加筋板前 4 階固有頻率(單位：Hz)

由表2 中數據可以看出，非協調廣義混合元模型的固有頻率值低于非協調位移元模型的結果數值，距離文獻解更加接近，結果精度更高。

比較圖9 中兩種方法的前4 階固有振型不難看出，非協調廣義混合元模型和非協調位移元模型的固有振型幾乎是一致的，這進一步驗證了本文方法的正確性。

3.3 雙加強筋偏心加筋板的自由振動特性分析

雙加強筋偏心加筋層合板[8]，尺寸如圖10 所示，尺寸參數：a= 1.2 m，l= 0.4 m，t= 0.01 m，t1=0.02 m，w= 0.04 m，h= 0.08 m。板部分總共三層，每層均為正交異性材料，上下兩個面層材料相同，假設面層的密度ρ1=1600 kg/m3，芯層的密度ρ1/ρ2=2，加強筋的材料參數與面層的一致。邊界條件為四邊固支，如圖11 所示，當x=0 和x=a時，u=v=w= 0，當y=0和y=a時，u=v=w= 0。對該結構進行振動分析。