彈性壓應力波下軸向功能梯度變截面梁動力壓曲穩定分析

陳得良， 汪亞運， 彭旭龍， 周 露

(長沙理工大學 土木與建筑學院， 長沙 410004)

彈性壓應力波下軸向功能梯度變截面梁動力壓曲穩定分析

陳得良， 汪亞運， 彭旭龍， 周 露

(長沙理工大學 土木與建筑學院， 長沙 410004)

基于微元法以及能量守恒原理，導出了軸向功能梯度變截面梁屈曲微分控制方程及應力波波前附加邊界條件，研究了軸向功能梯度變截面梁屈曲與壓應力波耦合動力屈曲問題。采用較為簡單的數值方法，即將位移函數按Taylor級數或是Chebyshev多項式展開，從而將軸向功能梯度變截面梁屈曲問題的變系數微分控制方程轉化為含參量的線性代數方程組，進而得到了含時間參量的動力屈曲問題特征方程，隨后對軸向功能梯度變截面梁動力屈曲問題進行了數值研究，探討了變截面和材料不均勻性對系統屈曲臨界力參數的影響。研究表明，該數值方法具有很好的精度和收斂性。

功能梯度梁； 變截面； 壓應力波； Taylor級數/Chebyshev多項式； 動力屈曲

功能梯度材料(Functionally Graded Materials, FGMs)是一種廣泛應用于航天工業、土木工程以及能源等領域的新型材料。由于應用領域的廣泛性和應用環境的復雜性，其工程性能要求也變得更為苛刻，因而其結構的穩定性研究也已成為關注的焦點。但是，由于結構截面和材料特性的連續性變化，給相關研究帶來一定的困難。

對于非均勻梁的屈曲穩定性問題，學者們首先研究了靜力穩定性問題。Huang等[1]基于Euler-Bernoulli梁理論，并通過將微分控制方程轉換成含有變系數的Fredholm積分方程，研究了不同邊界條件下具變剛度的軸向非均勻梯度柱的靜力屈曲問題。Huang等[2]則通過將模態函數展開為級數的形式，研究了彈性約束下軸向非均勻變截面梁臨界屈曲荷載數值解。文獻[3-7]同樣采用不同的數值分析方法研究了軸向功能梯度變截面梁/柱的靜力屈曲問題。

上述研究工作主要針對非均勻梁的靜力屈曲問題，然而非均勻梁不僅承受靜力屈曲，很多情況下還要承受動力屈曲。因此對于梁的動力穩定問題，許多學者展開了有益的探索。Teter[8]基于不同穩定準則，利用(Analytical-Numerical Method,ANM)研究了加筋開口截面柱在矩形脈沖荷載作用下的動力耦合臨界荷載。王安穩等[9-10]基于應力波理論和失穩瞬間能量轉換及守恒原理，導出了直桿動力分岔的能量準則，為探索動力失穩問題提供了新的思路。鄭波等[11-12]以王安穩等的能量準則為基礎，運用有限元法對直桿碰撞剛性壁的過程進行模擬，研究了直桿在彈性應力波下的動力失穩問題。鐘煒輝等[13]以合理的波前連續條件代替較為常用的波前固定假設，就階躍荷載作用下不同邊界條件軸心壓桿的動力分岔屈曲進行了分析，在文中考慮了無應力桿段及遠端邊界條件對軸心壓桿沖擊屈曲的影響，得到了更為合理的結果。毛柳偉等[14]指出應力波與屈曲耦合作用下結構屈曲的真實運動與鄰近運動是不同時刻、不同區域的比較，并采用能量守恒原理，根據屈曲時刻結構能量的轉換關系，建立了彈性直桿的屈曲控制方程，并推導得到了波前邊界條件。

綜上所述，學者們對變截面非均勻梁靜力屈曲和等截面均勻梁在應力波傳播與屈曲耦合情況下的結構動力屈曲問題進行了廣泛的研究，但很少有文獻涉及到變截面非均勻梁的動力屈曲問題。本文將利用微元法以及王安穩等提出的能量準則，導出軸向功能梯度變截面梁動力屈曲微分控制方程以及應力波波前附加邊界條件，并提出一種簡單有效的數值方法，分析不同邊界條件下梁的靜力屈曲問題和加載端簡支遠端固支非均勻梁在應力波反射前的動力屈曲問題。

1 軸向功能梯度變截面梁動力屈曲微分控制方程及邊界條件

1.1 控制方程的推導和邊界條件

圖1所示加載端簡支遠端固支的軸向非均勻梁，梁長為L(L=L1+L2)，橫截面面積為A(x)，截面慣性矩為I(x)，材料密度為ρ(x)，彈性模量為E(x)，橫向位移為w(x,t)。在t=0時刻，在梁端受幅值為P的沖擊荷載作用，假定P的作用時間足夠長。

圖1 力學分析模型

在應力波波及區域內選取微段如圖2所示。根據一維應力波理論，基于Euler-Bernoulli梁理論，在應力波波及區域L1中，利用微元法得到軸向功能梯度變截面梁的動力屈曲微分控制方程為

(0

(1)

式中,m(x)=ρ(x)A(x)，令D(x)=E(x)I(x)，w(x,t)=W(x)·T(t)，并引用無量綱坐標ξ=x/L1，則有

(2)

式中：ω為梁的固有頻率；λ*=PL12/(E0I0)、α*=ω2L14分別為臨界荷載參數和動力參數；к=m(x)/m0為材料線密度沿軸向的變化率；γ2=E0I0/m0；E0、I0和m0為常數。

相應邊界條件

簡支端(S)

W(ξ)=0，W″(ξ)=0

(3a)

固支端(C)

W(ξ)=0，W′(ξ)=0

(3b)

自由端(F)

(3c)

根據應力波準則，在壓縮波前ξ=1處有

W(ξ)=0，W′(ξ)=0

(4)

由控制方程式(2)和邊界條件式(3)、式(4)確定特征參數λ*和α*的值還缺乏一個補充條件。

圖2 梁微段受力分析

1.2 軸向非均勻變截面梁波前附加邊界條件

(5)

動能增量為

(6)

外力功增量為

ΔW=P·u1(0,t)

(7)

定義

(8)

由能量守恒定律ΔW=ΔU+ΔK，以及式(8)，有

Upre=Ubuc+ΔK

(9)

式(9)表示失穩能量轉化過程服從守恒定律，其中Upre和Ubuc分別表示失穩時釋放出的部分變形能和新的屈曲變形能。將式(9)對時間t求導，即可得能量轉換率守恒條件

(10)

在動力失穩中ΔK2相對ΔK1為小量，因而可以忽略ΔK2。將式(6)和式(8)代入臨界條件式(9)，并采用分離變量w(x,t)=W(x)·T(t)以及無量綱坐標ξ=x/L1，可以得到式(2)以及相應邊界條件式(3)和式(4)。將式(6)和式(8)代入臨界條件式(10)，即可得到非均勻梁波前附加邊界條件

W″(ξ)=0，(ξ=1)

(11)

邊界條件式(3)、式(4)和式(11)構成微分控制方程式(2)第二式的完備定解條件。

2 數值求解方法

式(2)為變系數高階微分控制方程，無法直接求解。在此，將位移函數以Taylor級數或第一類Chebyshev正交多項式展開，其統一形式為

(12)

式中:li為未知系數;N為級數展開項數;Ti(ξ)為區間[0,1]上的第一類Chebyshev正交多項式或Taylor展開式。

將式(12)代入微分控制方程式(2)第二式中，則有

(13)

在式(13)兩邊同乘ξt(t=0,1,…,N-5)，并對ξ從0～1積分，則得到N-4個獨立的線性代數方程

(14)

其中,

為確定N個未知系數li，需引入另外4個獨立的方程。將式(12)代入邊界條件式(3)和式(4)，并采用與式(14)相同的形式表達，即和式(14)形成封閉的線性代數方程組，其矩陣表達式為

(15)

式(15)為所求問題的特征方程。因方程存在非平凡解，則式(15)的系數矩陣行列式為0，即

(16)

在靜力屈曲問題中α*=0，可由式(16)直接得到靜力屈曲臨界力參數。對于屈曲和應力波耦合動力壓曲穩定，還需利用波前附加約束邊界條件式(11)，其具體分析方法為：對于給定沖擊荷載P，逐步增大動力參數α*的值進行特征值分析，當滿足波前附加條件式(11)時，即可得到相應各階屈曲臨界力參數和模態。

3 數值算例及討論

3.1 等截面均勻梁

為了說明本文方法的可靠性和收斂性，首先對等截面均勻梁的靜力和動力壓曲穩定進行分析，并與有關文獻結果進行對比分析。在動力參數α*=0時，對應靜力壓曲穩定的情況。在此處，考慮不同邊界條件均勻梁的靜力壓曲，其中屈曲臨界力參數λ*=PL2/EI，E和I為常數，相應臨界力參數數值解列于表1。

利用表1，比較本文方法所得結果與文獻[2]給出的精確解和數值解可知，該方法能快速收斂，且隨位移函數展開項數N的增加，所得結果與精確解的誤差將大大降低。說明該方法求解靜力屈曲問題是精確有效的。

對于動力穩定分析，動力參數α*>0恒成立。為了驗證所提出的方法求解加載端簡支遠端固支梁動力屈曲問題的有效性，此處采用Taylor級數展開的方法計算文獻[11]中的模型。計算模型長度L=5 m，彈性模量為E=207 GPa，橫截面面積A=4.0×10-4m2，材料密度ρ=7 778 kg/m3，截面慣性矩I=1.333×10-8m4。并假定壓應力波在傳到L1=2.5 m時，梁發生動力失穩，前3階動力屈曲臨界荷載列于表2。其動力屈曲模態和彎矩圖則分別如圖3和圖4所示。

表1 不同邊界均勻梁的屈曲臨界力參數(α*=0)

表2 應力波在固定端反射前對應的前三階動力屈曲臨界荷載

表2表明，當取N=16時，所得動力屈曲臨界荷載值和文獻相比非常接近，誤差小于0.13%。本文所得到的應力波反射前的動力失穩模態和彎矩圖也與文獻[9-11]一致。說明本文的方法不僅能有效分析靜力穩定性問題，同時也能有效用于動力穩定性問題的研究，且具有很好的收斂性和可靠性。

圖3 壓應力波未反射時的動力屈曲模態

圖4 壓應力波未反射前梁的彎矩圖

在后續工作中將采取Taylor級數展開的方法并取N=16分析計算變截面梁和軸向材料非均勻梁的靜力屈曲和動力屈曲問題。

3.2 變截面梁

對變截面梁靜力與動力屈曲的研究，取如下參數，彈性模量E和密度ρ為常量，其截面變化的情況如下[15]：

情況1A/A0=I/I0=1+aξ，截面高度一定，寬度沿軸向線性變化；

情況2A/A0=1+aξ，I/I0=(1+aξ)3，截面寬度一定，高度沿軸向線性變化。

其中a為幾何參數，即截面延梁長度L的變化率。在長度L1上，幾何參數則變為a*=(L1/L)a。兩種情況下，梁在應力波反射前的前三階動力屈曲臨界力參數分別列于表3和表4。圖5則給出情況1中a=0.1時應力波傳播長度對屈曲臨界荷載Pcr=λ*E0I0/L2的影響圖，圖5中Pcr max取幾種情況臨界荷載中的最大值。

表3和表4分別研究了矩形截面梁截面沿軸向變化的兩種情況下的靜力屈曲臨界力參數和動力屈曲前三階臨界力參數。表3和表4表明，隨著截面幾何參數的增加，結構所能承受的屈曲臨界力參數也將隨之增加。等效長度梁的動力屈曲臨界力參數大于靜力屈曲臨界力參數。比較表3和表4，可知兩種不同的截面變化方式對結構屈曲臨界力參數的影響不同，其中以情況2更能提高結構的抗屈曲能力，這與工程實際相符合。此外，從表3、表4及圖5可知梁的等效長度對其抗屈曲能力的影響，長度越長梁屈曲臨界荷載越小。

表3 變截面梁A/A0=I/I0=1+aξ時靜力和動力屈曲臨界力參數

表4 變截面梁A0=1+aξ，I/I0=(1+aξ)3時靜力和動力屈曲臨界力參數

3.3 軸向材料非均勻梁

本節將討論材料性能沿軸向變化對等截面梁動力穩定的影響，設材料密度ρ(ξ)和彈性模量E(ξ)有如下幾種形式：

方案1 材料密度為常量，彈性模量取特定的形式

方案2 材料密度為線性變化，彈性模量取特定的形式

圖5 壓應力波傳播長度對屈曲臨界荷載的影響

Fig.5 Effect of compression wave propagation length on the critical buckling load

在上述所取假定中，ρ0、E0和I0為常數；aj(j=0,1)、bjh(j=0,1,h=0,…,j+4) 為任意常數，具體取值如下

對于加載端簡支遠端固支梁，兩種方案得到的靜力屈曲臨界力參數和第一階動力屈曲臨界力參數列于表5。

表5 壓應力波反射前加載端簡支遠端固支梁的動力屈曲臨界力參數

從表5中數據可以看出，材料特性以不同形式的多項式表示時對結構抗屈曲能力有較大的影響，并隨多項式項數變化而發生波動。因此，在工程應用中結構外形限定時，可以通過改變材料特性的分布規律來改善結構性能以滿足工程需要。

4 結 論

(1) 本文主要就軸向功能梯度變截面梁靜力屈曲問題和加載端簡支遠端固支非均勻梁在應力波反射前的動力屈曲問題進行了分析，通過建立軸向非均勻梁屈曲條件，求得了應力波在固定端反射前的屈曲臨界力參數和相應屈曲模態。

(2) 基于能量準則，推導了軸向功能梯度變截面梁應力波波前附加邊界條件。并提出一種簡單有效的數值方法，從而避免了直接求解含變系數的高階微分方程。

(3) 通過和已有文獻數據的比較，驗證了方法的精度和有效性，研究了變截面和軸向材料非均勻性對屈曲臨界力參數的影響。文章研究內容對工程實際應用具有重要意義。

[1] HUANG Y, LI X F. Buckling analysis of non-uniform and axially graded columns with varying flexural rigidity[J]. Journal of Engineering Mechanics, 2011, 137(1): 73-81.

[2] HUANG Y, LUO Q Z. A simple method to determine the critical buckling loads for axially inhomogeneous beams with elastic restraint[J]. Computers and Mathematics with Applications, 2011, 61(9): 2510-2517.

[3] BELLMAN R, CASTI J. Differential quadrature and long-term integration[J]. Journal of Mathematical Analysis and Applications, 1971, 34(2): 235-238.

[4] SHAHBA A, RAJASEKARAN S. Free vibration and stability of tapered Euler-Bernoulli beams made of axially functionally graded materials[J]. Applied Mathematical Modelling, 2012, 36(7): 3094-3111.

[5] YILMAZ Y, GIRGIN Z, EVRAN S. Buckling analyses of axially functionally graded nonuniform columns with elastic restraint using a localized differential quadrature method[J]. Mathematical Problems in Engineering, 2013(2): 1-12.

[6] SINGH K V, LI G. Buckling of functionally graded and elastically restrained non-uniform columns[J]. Composites Part B: Engineering, 2009, 40(5): 393-403.

[7] SHVARTSMAN B, MAJAK J. Numerical method for stability analysis of functionally graded beams on elastic foundation[J]. Applied Mathematical Modelling, 2016, 40(5/6): 3713-3719.

[8] TETER A. Dynamic critical load based on different stability criteria for coupled buckling of columns with stiffened open cross-sections[J]. Thin-Walled Structures, 2011, 49(5): 589-595.

[9] 王安穩. 彈性壓應力波下直桿動力失穩的機理和判據[J]. 力學學報, 2001, 33(6): 812-820.

WANG Anwen. Mechanism and criterion for dynamic instability of bars under elastic compression wave[J]. Acta Mechanica Sinica, 2001, 33(6): 812-820.

[10] WANG Anwen, TIAN Wenying. Characteristic-value analysis for plastic dynamic buckling of columns under elastoplastic compression waves[J]. International Journal of Non-linear Mechanics, 2003, 38(5): 615-628.

[11] 鄭波, 王安穩. 彈性應力波下直桿分叉動力失穩特征值有限元分析[J]. 工程力學, 2006, 23(12): 36-40.

ZHENG Bo, WANG Anwen. Finite element character-value analysis for dynamic buckling of bars under elastic compression wave[J]. Engineering Mechanics, 2006, 23(12):36-40.

[12] 鄭波, 王安穩. 直桿碰撞剛性壁彈塑性動力后屈曲有限元分析[J]. 爆炸與沖擊, 2007, 27(2): 126-130.

ZHENG Bo, WANG Anwen. Finite element analysis for elastic-plastic dynamic postbuckling of bars subjected to axial impact[J]. Explosion and Shock Waves, 2007, 27(2): 126-130.

[13] 鐘煒輝, 郝際平, 雷蕾, 等. 考慮應力波的軸心壓桿沖擊分岔屈曲研究[J]. 振動與沖擊, 2010, 29(10): 201-205.

ZHONG Weihui, HAO Jiping, LEI Lei, et al. Impact bifurcation buckling of an axial compression bar with stress wave[J]. Journal of Vibration and Shock, 2010, 29(10): 201-205.

[14] 毛柳偉, 王安穩. 應力波傳播與屈曲耦合情況下結構動力屈曲控制方程及波前邊界條件的探討[J]. 固體力學學報, 2013, 34(2): 152-157.

MAO Liuwei, WANG Anwen. Derivation of governing equation for dynamic buckling of structures under coupled stress wave propagating and dynamic buckling[J]. Acta Mechanica Solida Sinica, 2013, 34(2): 152-157.

[15] HUANG Y, LI X F. A new approach for free vibration analysis of axially functionally graded beams with non-uniform cross-section[J]. Journal of Sound and Vibration, 2010, 329(11): 2291-2303.

Dynamic buckling of axially functionally-graded beams with non-uniform cross-section under elastic compression stress wave

CHEN Deliang, WANG Yayun, PENG Xulong, ZHOU Lu

(School of Civil Engineering and Architecture, Changsha University of Science and Technology, Changsha 410004, China)

Here, the dynamic buckling problem during the buckling of axially functionally graded beams with variant cross-section being coupled with compression stress wave was investigated. The buckling differential governing equation and the wavefront additional boundary conditions of compression stress wave for an axially functionally graded beam with variant cross-section were established based on the differential-element method and the principle of energy conservation. A simpler numerical method was introduced to transfer this varying-coefficient governing differential equation into a set of linear algebraic equations with the displacement function being expanded using Taylor series or Chebyshev polynomials. Then the eigen-equation for the axially functionally graded beam with non-uniform cross-section was obtained. Moreover, a numerical investigation for the dynamic buckling of the axially functionally graded beams was conducted to discuss the effects of variant cross-section and material inhomogeneity on the system’s critical buckling force parameters. The results showed that the proposed method has good accuracy and convergence.

axially functionally graded beam; non-uniform cross-section; compression stress wave; Taylor series or Chebyshev polynomials; dynamic buckling

國家自然科學基金資助(11172051; 11202038); 湖南省自然科學基金(2015JJ4006);長沙理工大學研究生科研創新(CX2016SS02)

2016-02-29 修改稿收到日期：2016-05-19

陳得良 男，博士，教授， 1971年生

O347

A

10.13465/j.cnki.jvs.2017.13.004