基于半解析法的圓柱殼結構自由振動特性分析

李海超， 龐福振， 田宏業， 劉江濤

(哈爾濱工程大學 船舶工程學院，哈爾濱 150001)

柱殼結構被廣泛應用于航空航天、船舶、化工、機械等領域，因此開展一般邊界條件下柱殼結構自由振動特性分析，明確其自由振動特性規律，對豐富圓柱殼結構基礎理論及指導工程應用具有重要的意義。在此方面，Rayleigh[1]對光圓柱殼的拉伸振動和彎曲振動進行了研究，給出了無限長圓柱殼在真空環境中的自由振動固有頻率計算公式。Galletly等[2]采用數值算法求解不同端部結構的圓柱殼的自由振動。Cheng等[3]以圓柱殼-圓板為模型，采用Rayleigh-Ritz法求解特征值方程從而求解其自由振動特性。駱東平等[4]采用Flügge和有限差分法，以環肋增強圓柱殼為研究對象。Yu[5]討論了簡支和固支邊界條件下有限長圓柱殼的自由振動特性。王宇等[6]基于Love殼體理論對固支-自由約束條件下受徑向載荷的薄壁圓柱殼構件開展受迫振動響應特征分析。汪志強等[7]采用Flügge經典薄殼理論和波傳播方法討論了正交各向異性圓柱殼的自由振動問題，且所提出方法可考慮正交各向異性圓柱殼在復雜和受外力的情況。向宇等[8]提出了一種半解析方法，求解飽和多孔介質圓柱殼動力學問題。Li等[9]依據Flügge薄殼理論建立振動分析模型，采用Jacobi-Ritz法分析了均厚度以及變厚度圓柱殼、球殼等組合殼體結構的自由及受迫振動特性。龐福振等[10]基于Flügge殼體振動理論，將改進精細傳遞矩陣法應用于水下加筋柱殼聲輻射問題，分析了經典邊界條件、結構損耗因子、流體介質以及殼體厚度對結構聲輻射的影響。Pang等[11]提出了一種求解雙曲率殼自由振動的半解析方法。對于區域能量分解法，最重要的是廣義變分原理理論的成熟，我國著名數學和力學家錢偉長[12]對廣義變分原理進行了系統的研究，對廣義變分原理成熟應用起到了極大的推動作用，也為區域分解法在彈性體中的應用提供了堅實的基礎。Zienkiewicz等[13]在其《有限元法》一書對變分約束條件進行了講解，確定拉格朗日乘子系數，并應用到彈性力學中。宋文煜[14]應用區域分解法，分別研究了鉆柱縱向、扭轉、橫向以及耦合振動的固有頻率及長度、鉆柱壁厚等對鉆柱固有頻率的影響。Zhou等[15]提出了一種求解圓柱和實心圓柱自由振動的通用方法。分析過程基于小應變、線性和精確彈性理論。利用Chebyshev多項式級數乘以邊界函數滿足幾何邊界條件作為可容許函數，應用Ritz法導出了圓柱的頻率方程。

由以上分析可知：一方面現有文獻對任意邊界條件下圓柱殼自由振動特性研究較少；另一方面，現有研究方法尚未形成統一的形式，且現有半解析方法仍有待進一步豐富。為此，本文基于區域能量分解法，開展圓柱殼自由振動特性分析，旨在提出統一的求解公式，為任意邊界條件圓柱殼自由振動特性分析提供數據積累和方法依據。此外，運用區域能量分解法時，可通過改變邊界控制參數可快速地分析不同組合邊界條件下殼體的振動特性，無須重新形成整個殼體的質量和剛度矩陣，從而大幅提高計算效率。

1 柱殼結構區域能量分析模型的建立

本文基于區域能量分解方法的基本原理，對其進行深度擴展，并運用到求解柱殼結構的振動特性當中。圓柱殼示意圖如圖1所示，結構的長度為L，半徑為R，殼體厚度為h，結構的坐標系統如圖1所示。

圖1 柱殼結構的理論模型Fig.1 Theoretical model of cylindrical shell structure

1.1 柱殼結構的能量泛函建立

將圓柱殼體沿圓柱的軸線方向均勻地截斷成NL段，即每一段的長度Li=L/NL。根據修正的Hamilton原理，考慮到每一段的能量和相鄰兩段之間的影響，圓柱殼體的總勢能為

(1)

式中：TL,i，UL,i，WL,i分別為圓柱殼體的第i段的動能、應變能、外力功和附加能量泛函；ΠK,L為相鄰分段i和i+1之間的附加界面勢能。

當忽略圓柱殼體旋轉慣性的條件下，圓柱殼體的第i段的動能可以表示為

(2)

式中：ui，vi，wi分別為不同方向的位移矢量；ρ為結構的質量密度；hi為第i段的結構厚度；Si為結構中面面積。

根據Reissner-Naghdi’s線性薄殼理論，i分段的最大結構應變能可以表示為

(3)

假設外部載荷全部作用在中面位置處，圓柱殼結構的第i分段分布著沿x方向、θ方向、z方向外力fu,i,fv,i,fw,i，此時結構第i分段的外力作功為

(4)

圓柱殼體的第i段、第i+1段的界面上的附加能量泛函為

(5)

式中：λ,β,?,ψ分別為柱殼結構分段i與i+1交界面處的未知拉格朗日乘子；Θu,Θv,Θw和Θr分別為柱殼結構分段i與i+1交界面處位移協調方程，它們可以表示為

(6)

將式(2)～式(5)代入式(1)，并且根據廣義變分原理，對ui，vi，wi，ui+1，vi+1，wi+1，?wi/?x，?wi+1/?x，λ，β，?和ψ做變分運算，可得到

(7)

(8)

(9)

(10)

將式(7)～式(10)代入式(1)后得到新的能量泛函

(11)

式中：Nx=λ；Nθ=β;Qx=?;-Mx=ψ。

為了保證數值算法的計算穩定性，在式(11)的基礎上添加一項子結構交界面位移連續方程的最小二乘加權參數殘值Πκ,L，此時結構的完整能量泛函表示為

(12)

(13)

式中：κu,κv,κw和κr分別為柱殼分區后第i段和第i+1段分區界面加權參數。

在式(13)中引入邊界條件控制參數?u，?v,?w,?r, 通過控制參數的選取來控制邊界條件。此時的能量泛函可以表示為

(14)

表1給出了邊界條件控制權參數，通過選擇不同的權參數來得到不同的邊界條件。

表1 不同邊界條件對應的控制參數 ?t(t=u,v,w,r)Tab.1 Control parameters ?t(t=u,v,w,r)corresponding to different boundary conditions

1.2 位移函數的表達式

對于柱殼結構，柱殼結構的位移函數采用Chebyshev行列式[16]和傅里葉級數進行展開，系統的位移可以寫成

(15)

(16)

(17)

式中：P為圓柱殼母線方向上位移分量的Chebyshev多項式截斷值;N為殼體周向位移分量的Fourier級數截取階數；Φp(x)為Chebyshev多項式，表示圓柱殼體軸向階數;fni(θ)為傅里葉級數，表示圓柱殼體周向波數，它們分別表示為

Φ0(x)=1,Φ1(x)=x,Φi+2(x)=2xΦi+1(x)-Φi(x)

(18)

fn1(θ)=cos(nθ)+sin(nθ)

(19)

fn2(θ)=sin(nθ)+cos(nθ)

(20)

fn3(θ)=cos(nθ)+sin(nθ)

(21)

2 柱殼結構的振動特性分析

2.1 柱殼結構自由振動特性方程的求解

由“1.1”節可知，結構的整體能量泛函可表示為

(22)

對自由振動而言

(23)

對未知位移系數進行變分運算，此時可得到結構的自由振動特性方程

(24)

(25)

其中，

(26)

(27)

(28)

(29)

對式(29)求解，即可得結構的模態陣型。為統一表述，下文無量綱公式定義為

Ω=ωR(ρ(1-μ2)/E)1/2

(30)

2.2 最小二乘加權參數與分段數的收斂性研究

本節選取結構參數為：L=6 m，R=1 m，h/R=0.01，E=210 GPa，ρ=7 800 kg/m3，μ=0.3。彈性參數取值設為ku=kv=2×108N/M[17]，位移容許函數在母線方向上的切比雪夫正交多項式的截斷值為8，結構的分段數NL=2。F-E邊界條件下最小二乘加權參數收斂性分析如表2所示。

由表2可知，當加權參數增大時，結構的固有頻率也隨之增大，然而當加權參數達到某一個臨界取值時，結構的計算結果已趨于穩定，在實際數值計算中，加權參數的取值不可能無窮大，因此，加權參數取臨界值即可。從表2可以看出，當κ取值1×1012時，已達到臨界值，但為保證計算的絕對數值穩定性，在后續的計算中，最小二乘殘差加權系數取值為κ=1×1014。

表2 最小二乘加權參數收斂性分析Tab.2 Convergence analysis of least squares weighted parameters

為兼顧求解精度和計算效率，需開展區域分段收斂性分析。選取結構參數與表2保持一致，并且加權參數取值為κ=1×1014。同時，為驗證本文方法有效性，將本文計算方法與有限元法進行了對比，不同分段數下結構無量綱化固有頻率Ω如表3所示。

表3 區域分段數收斂性分析Tab.3 The convergence analysis of the number of regional segments

從表3可知，當NL=8時，其計算結果已收斂，且與有限元法相比具有較高的計算精度。因此，在后續的計算中，如不做特別說明，其分段數取值為NL=8。

3 柱殼結構自由振動特性分析

3.1 經典邊界條件

在前文研究基礎上，本節主要對經典邊界條件(C-C)下的結構自由振動特性進行計算分析。結構的材料參數為：E=210 GPa，ρ=7 800 kg/m3，μ=0.3。結構參數與表2相同。兩端固支邊界條件下結構頻率參數計算結果如表4所示。

從表4可知，在C-C邊界條件下的計算結果與文獻的計算結果和有限元仿真結果吻合良好，為進一步驗證本文方法的有效性，計算S-S,F-F,F-S,C-S邊界條件下的頻率參數Ω，并與其他文獻結果對比，其對比結果如表5所示(圓柱殼幾何參數及材料參數與文獻中圓柱殼相同，h/R=0.05，R/L=0.05，m=1)。

表4 圓柱殼在C-C邊界條件下結構頻率參數ΩTab.4 Structural frequency parameters Ω of a cylindrical shell under C-C boundary condition

表5 圓柱殼在不同邊界條件下文獻對比結果Tab.5 Comparison of cylindrical shells under different boundary conditions

從表5中對比結果可以看出，在不同邊界條件下，本文方法計算結果與文獻對比結果吻合良好，更進一步說明本文方法的有效性。

3.2 一般邊界條件

在前文研究基礎上，本節主要對圓柱殼結構一般彈性邊界條件下的自由振動特性進行研究。材料參數及結構參數與“3.1”節相同。使用區域能量分解法，可通過改變邊界控制參數快速形成不同組合邊界條件，這也是本方法的優勢之一。首先定義三種彈性邊界：E1彈性邊界，其中只有u向為彈性約束，其余為固定邊界條件(u≠0,v=w=?w/?x=0)；E2彈性邊界，其中只有v向為彈性約束，其余為固定邊界條件(v≠0,u=w=?w/?x=0)；E3彈性邊界，其中只有v向與u向為彈性約束，其余為固定邊界條件(v≠0,u≠0,w=?w/?x=0)。在u向的彈性約束剛度值大小為κu=4×106N/m，在v向的彈性約束剛度值大小為κv=8×106N/m。表6～表8給出了E1-E2,E1-E3,E2-E3邊界條件下結構的無量綱頻率參數。

表6 圓柱殼在E1-E2邊界條件下頻率參數ΩTab.6 The frequency parameters Ω of a cylindrical shell under the E1-E2 boundary condition

表7 圓柱殼在E1-E3邊界條件下頻率參數ΩTab.7 The frequency parameters Ω of a cylindrical shell under the E1-E3 boundary condition

表8 圓柱殼在E2-E3邊界條件下頻率參數ΩTab.8 The frequency parameters Ω of a cylindrical shell under the E2-E3 boundary condition

由表6～表8可知，在彈性邊界條件下，本文的計算結果與有限元計算結果吻合良好，說明本文方法可應用于一般邊界條件下圓柱殼結構自由振動特性分析。

在此基礎上，計算另外三組邊界條件(SS-F，SS-SD，SS-C)圓柱殼結構自由振動特性，為一般邊界條件下圓柱殼結構自由振動特性分析提供數據積累(見表9)。

表9 SS-F，SS-SD，SS-C邊界條件下圓柱殼結構頻率參數ΩTab.9 The frequency parameters Ω under the boundary conditions of SS-F, SS-SD and SS-C

從表9可知，母線與半徑長度比對計算結果影響較大，在相同半徑長度與厚度比情況下，母線長度增大，結構的固有頻率會隨之增大。同理，在相同的母線長度與半徑的長度比下，結構厚度增大，結構的固有頻率也隨之增大，但在不同軸向波數情況下，其增長率并不一致。

4 結 論

本文基于區域能量分解法，位移容許函數采用第一類正交切比雪夫多項式與三角級數表示，基于廣義變分原理和最小二乘殘差原理，建立了一般邊界條件下圓柱殼結構自由振動分析模型，并對最小二乘殘差加權系數與區域分段數對于數值計算的穩定性進行研究分析，在此基礎上，從固有頻率參數出發，對經典邊界條件及一般彈性邊界下圓柱殼結構自由振動特性進行了分析，驗證了本文方法的有效性。為柱殼結構的前期快速設計提供一種可靠的手段。研究成果可為一般邊界條件下圓柱殼結構自由振動特性分析提供數據積累和方法依據。