基于近場動力學微分算子的變截面梁動力特性分析方法

李志遠，黃 丹，閆康昊

(河海大學工程力學系，南京 211100)

相比于等截面梁，變截面梁能夠提供更優化的強度、剛度和質量分布，可滿足特殊的結構與功能要求且經濟性好，因而廣泛應用于土木、機械、航空航天、精密儀器儀表等諸多工程領域中 [1 ?2]。國內外學者近年來采用有限元法[3?4]、微分求積法[5]、攝動法[6]、超幾何函數法[7?8]、瑞利-里茲法[9]等各種不同方法開展了變截面梁的振動特性分析。最近，LI 等[10? 11]提出了一種直接基于二維彈性理論高精度求解任意變截面梁自由振動的方法，但在處理含缺陷的變截面梁時依然存在困難。

傳統的數值方法大多是基于局部微分理論思想，如有限元法、有限差分法和有限體積法。因此，在存在間斷點和奇異點以及高階導數的情況下，傳統的數值方法面臨著巨大的挑戰。近場動力學(Peridynamics, PD)[12? 13]是一種基于空間非局部積分思想分析固體力學問題的新理論。由于在求解不連續問題時表現出獨特的優勢，該理論近二十年來受到計算力學及各工程領域相關研究者的廣泛關注[14?15]。秦洪遠等[16]和顧鑫等[17]通過改進鍵型PD 模型有效模擬了混凝土結構多裂紋擴展、沖擊破壞及侵徹問題。張鈺彬等[18]在常規態型PD 模型中加入等效水壓力項實現了頁巖水力壓裂過程的仿真。王涵等[19]建立非常規態型PD 熱黏塑性本構模型并應用于金屬材料在沖擊荷載下的熱黏塑變形與破環。馬鵬飛等[20]利用最小二乘法優化鍵型PD 中的微模量，有效地降低了邊界處的誤差，并提高了模擬裂紋擴展過程的收斂速度?；赑D 非局部相互作用的思想，MADENCI 等 [21 ?22]于2016 年提出了近場動力學微分算子(Peridynamic differential operator, PDDO)。相比于PD，PDDO 可以方便地將傳統的局部微分轉化為非局部積分形式，不需要任何參數的轉化及邊界處理，并且適用范圍更廣。目前，PDDO 已在多相流[23]、復合材料[24 ?25]、熱力耦合[26? 27]等問題的求解方面取得了初步的應用。特別是近年DORDUNCU[25]將PDDO應用于層合梁的應力分析，LI 等[27]基于PDDO 分析了變截面功能梯度梁的靜力彎曲問題。

本文在已有研究基礎上，進一步將PDDO 應用于變截面梁的動力問題求解。將變截面梁的動態微分控制方程與邊界條件通過PDDO 轉化為對應的非局部積分形式，應用變分原理和拉格朗日乘數法，將非局部積分形式的控制方程轉化為標準特征值問題表達形式，從而求得結構自振頻率與模態。通過對等截面梁的自由振動分析并與已有解析解比較，驗證了本方法良好的收斂性與準確性。其次，分析了下邊界分別以一次、二次連續變化的變截面梁的自由振動，證明了本方法在任意變截面梁自由振動問題分析中的適用性與通用性。最后，通過含孔變截面梁的自由振動分析，說明本方法在含缺陷構件振動分析和損傷識別等問題方面的潛力。

1 問題描述

考慮變截面簡支梁如圖1 所示，梁長L，左端截面高度H，上、下邊界可分別被描述為h1(x)、h2(x)。梁的本構關系(平面應力假設)為：

圖1 任意變截面簡支梁Fig. 1 Simply supported beam with arbitrarily and continuously varying cross-section

式中：σx、σy分別為x、y方向正應力；τxy為切應力；u、v分別為x、y方向位移；E為彈性模量；μ為泊松比。動態微分控制方程為：

式中，ρ 為材料密度。將式(1)代入式(2)得微分控制方程：

對于梁的自由振動分析，其周期位移分量可以表示為：

式中：U、V為振型函數；ω為自振頻率；i 為虛數單位。將式(4)代入式(3)得：

針對不同截面梁，即不同的上、下邊界表述函數hk(x) (k=1, 2)，上、下邊界條件為：

2 變截面梁振動分析的PDDO 求解

2.1 PDDO 簡述

將二維標量函數f(x)=f(x,y)，二階泰勒展開并忽略余項可得[21? 22]：

2.2 非局部化

將式(9)代入式(5)，可得非局部形式的控制方程：

可采用類似的方法得到本構方程和邊界條件非局部形式，此處不作贅述。

2.3 離散格式

求解域D可被均勻離散為N個物質點，如圖3所示。Δx為物質點離散間距，A=(Δx)2為物質點代表面積。物質點x(i)附近δ 距離內的物質點組成近場范圍H(i)，可表達為：

圖3 均勻PD 離散Fig. 3 Illustrations of uniform PD discretization

從圖2 可以看出，近場范圍可以為任意形狀。傳統PD 中一般采用圓形近場范圍，并需要額外的體積修正步驟來提高精度。為避免考慮額外的步驟，本文采用正方形近場范圍，并選取近場范圍尺寸δ=3Δx，如圖3 中物質點x(i)的近場范圍H(i)。對于含缺陷的結構，則需將相互作用跨過缺陷的物質點從原近場范圍物質點集合中去除。靠近邊界的物質點x(j)、x(k)自然會有非對稱的近場范圍H(j)、H(k)。邊界條件可直接施加在邊界物質點上。

圖2 物質點間的相互作用Fig. 2 Interaction of material points with arbitrary family

同樣也可以得到本構方程和邊界條件離散形式。

2.4 求解體系

將控制方程和邊界條件表達為代數方程組：

式中：L和c分別為離散控制方程和邊界條件中PD 函數所產生的系數矩陣；u為待求位移向量，可表示為：

于是，式(24)可轉化為標準特征值問題：

其中：

3 數值算例

本節應用所提出的方法開展多種變截面梁的自由振動分析。首先，通過等截面梁的自由振動分析，驗證本文方法的準確性與收斂性。其次，分別對下邊界一次、二次連續變化梁的自由振動分析，檢驗本方法對于任意變截面梁自由振動分析的適用性。最后，考慮含孔洞的下邊界二次變化梁，分析孔徑對變截面梁自由振動的影響。本節所有算例，如無特殊說明，參數均取為：梁長L=10 m，彈性模量E= 206 GPa，泊松比μ= 0.3，密度ρ = 7800 kg/m3，上邊界h1(x) = 0，物質點離散間距Δx= 25 mm。采用無量綱自振頻率的表達式：

3.1 等截面梁

考慮等截面簡支梁，如圖4 所示，截面高度H= 1 m，下邊界h2(x) =H。

圖4 等截面簡支梁Fig. 4 The simply supported beam with constant cross-section

表1 中給出了等截面梁前8 階無量綱自振頻率，并與解析解[10]作對比。可以看出本文解與解析解具有很好的一致性，最大相對誤差為0.946%，說明了本方法的準確性與高精度。

表1 等截面梁的無量綱自振頻率Table 1 Non-dimensional frequencies of the beam with constant cross-section

為進一步驗證本方法的收斂性，分別考慮物質點離散間距Δx= 25 mm、31.25 mm、40 mm、50 mm，前3 階無量綱自振頻率的相對誤差如圖5所示。可以看出，隨著物質點離散間距的減小，相對誤差迅速減小，說明了本方法良好的收斂性。

圖5 各離散間距下前3 階無量綱自振頻率的相對誤差Fig. 5 Relative error of first 3 non-dimensional frequencies for different dense mesh

3.2 變截面梁

考慮下邊界一次變化的梁和下邊界二次變化的梁，如圖6 和圖7 所示，左端截面高度H= 0.5 m。下邊界一次變化的梁右端截面高度為Ha，下邊界可表述為ha(x) = (Ha?H)x/L+H。下邊界二次變化的梁跨中截面高度為Hb，下邊界可描述為hb(x) =(H?Hb)(2x/L?1)2+Hb。

圖6 下邊界一次變化的簡支梁Fig. 6 The beam with linearly varying lower surface

圖7 下邊界二次變化的簡支梁Fig. 7 The beam with parabolic convex lower surface

表2 列出右端截面高度Ha= 0.75 m、1.0 m 情況下，下邊界一次變化梁的前8 階無量綱自振頻率。表3 所示為不同跨中截面高度Hb= 0.8 m、 1.0 m情況下，下邊界二次變化梁的前8 階無量綱自振頻率。從兩表中對比可見，本文解與解析解[10]吻合較好，驗證了本方法對于邊界一次、二次變化梁自由振動求解有效性，也間接說明了本方法對于任意變截面梁自由振動分析的適用性。

表2 下邊界一次變化梁的無量綱自振頻率Table 2 Non-dimensional frequencies of the beam with linearly varying lower surface

表3 下邊界二次變化梁的無量綱自振頻率Table 3 Non-dimensional frequencies of the beam with parabolic convex lower surface

圖8 所示為下邊界一次、二次變化變截面梁(Ha/H=Hb/H= 2)在y= 0.4H水平線處的前3 階軸向與橫向振型圖。為方便直觀顯示，圖8 中結果均進行歸一化處理，可以發現前3 階振型均以橫向振動為主。由于結構的對稱性，下邊界二次變化梁的橫向振幅也是對稱的，而下邊界一次變化梁的橫向振幅隨著梁厚度的增加在逐漸減小。

圖8 下邊界一次、二次變化梁的前3 階振型(Ha/H = Hb/H = 2, y = 0.4H)Fig. 8 First 3 mode shapes in x and y directions at y = 0.4H for the beam with linearly varying and parabolic convex lower surface (Ha/H = Hb/H = 2, y = 0.4H)

3.3 含孔變截面梁

在3.2 節下邊界二次變化梁(Hb/H= 2)的基礎上，考慮跨中位置存在一個圓孔，如圖9 所示，圓心Q(5, 0.5)，半徑r。

表4 列出不同孔徑r= 0.2 m、0.3 m、0.4 m 下，下邊界二次變化的含孔梁的前8 階無量綱自振頻率。由表4 中結果可見，隨著孔徑的增加，第2 階、4 階、5 階、7 階頻率減小，第1 階、3 階、6 階頻率變化較小，第8 階頻率增加。圖10 所示為下邊界二次變化含孔梁(r/H= 0.6)的前4 階軸向應力模態。為方便分析，軸向應力結果也同樣進行了歸一化處理，并呈現在對應歸一化處理后的同階振型上。可以看出，前4 階振型均以橫向振動為主，孔邊也存在明顯的應力集中現象。本算例初步說明了基于PDDO 的非局部方法在含缺陷構件振動分析和損傷識別問題方面的潛力，為相關問題分析提供可以進一步深入研究的新思路。

圖10 下邊界二次變化含孔梁的前4 階軸向應力模態(r/H = 0.6)Fig. 10 First 4 stress modes in x direction for the beam with a hole and parabolic convex lower surface (r/H = 0.6)

表4 下邊界二次變化含孔梁的無量綱自振頻率Table 4 Non-dimensional frequencies of the beam with a hole and parabolic convex lower surface

4 結論

本文嘗試基于非局部PDDO 思想開展變截面梁的動力特性分析，主要初步結論如下：

(1)可應用PDDO 將變截面梁的動態微分控制方程與邊界條件轉化為對應的非局部積分形式，并進一步應用拉格朗日乘數法與變分原理求解得到自振頻率與模態。

(2)通過對等截面梁和下邊界不同變化梁的自由振動分析，驗證了本方法的精度和收斂性，以及對于任意變截面梁自由振動問題的適用性。

(3)通過對含孔變截面梁的動力特性分析，說明了本方法在含缺陷結構的振動分析和損傷識別等工程問題方面的潛力，可為含缺陷變截面構件的動力分析問題提供新思路?；谶@一思路，或可進一步研究不規則、非均質、不連續、含缺陷構件的動力學問題和結構損傷識別。