王永亮

(1. 中國礦業(yè)大學(北京)力學與建筑工程學院，北京 100083；2. 中國礦業(yè)大學(北京)煤炭資源與安全開采國家重點實驗室，北京 100083)

桿系結(jié)構(gòu)的損傷問題廣泛存在于工程實際中，工程桿系結(jié)構(gòu)大多帶裂紋工作[1-2]，裂紋損傷的存在會改變整個結(jié)構(gòu)的力學性能，影響結(jié)構(gòu)的安全性和適用性[3]。研究含多裂紋損傷梁構(gòu)件的動力特性、準確預(yù)測屈曲承載力可以有效地保障結(jié)構(gòu)在全生命周期內(nèi)安全使用[4]。曲線形梁構(gòu)件由于幾何形態(tài)復(fù)雜，容易誘發(fā)彈性屈曲失穩(wěn)[5]，精確評估各類曲線梁線型、不同曲梁夾角下深梁、淺梁的屈曲荷載成為結(jié)構(gòu)災(zāi)害分析的重要依據(jù)。曲梁中裂紋損傷的存在增加準確預(yù)測屈曲失穩(wěn)承載能力的難度，理論模型、解析方法等往往難以有效分析[6-7]。準確預(yù)測不同裂紋損傷位置、大小、數(shù)目工況下屈曲荷載承載力以及分析裂紋損傷對屈曲失穩(wěn)的影響機理[8]，成為理論研究和工程實踐的需求。

有限元法被發(fā)展和應(yīng)用于求解含裂紋損傷曲梁的彈性屈曲荷載和屈曲模態(tài)[9-11]，但解答精度依賴于網(wǎng)格劃分質(zhì)量，解答因網(wǎng)格劃分難免引入誤差[12]。有限元網(wǎng)格自適應(yīng)分析方法可有效地優(yōu)化網(wǎng)格分布，在直線梁彈性屈曲[13]、板殼振動[14]、含損傷梁振動[15]、巖體變形和斷裂[16]等問題求解中展示出很好的求解效力。本文將建立圓弧形曲梁裂紋的截面損傷缺陷比擬方案，進行裂紋大小(深度)、位置、數(shù)目的模擬，引入變截面Euler-Bernoulli梁的h型有限元網(wǎng)格自適應(yīng)分析方法[15]，求解含裂紋損傷圓弧曲梁彈性屈曲問題，得到優(yōu)化的網(wǎng)格和滿足預(yù)設(shè)誤差限的高精度屈曲荷載和屈曲模態(tài)。文中給出求解多種圓弧曲梁彈性屈曲數(shù)值算例，對網(wǎng)格自適應(yīng)劃分以及彈性屈曲解答的收斂性進行了討論，對曲梁夾角、損傷位置、數(shù)目、大小等因素影響彈性屈曲荷載和屈曲模態(tài)進行了分析，檢驗了網(wǎng)格自適應(yīng)劃分的有效性。

1 圓弧曲梁彈性屈曲

圖 1 含裂紋損傷曲梁坐標系和符號Fig. 1 Coordinate systems and symbols of cracked curved beam

圖 2 含裂紋損傷圓弧曲梁截面損傷和加載示意圖Fig. 2 Diagram of cross-section damage defect and loading for circularly curved beam with crack damage

表 1 圓弧曲梁根據(jù)徑厚比和夾角分類Table 1 Categories of circularly curved beams according to ratio of radius and thickness and subtended angle

本文研究圓弧曲梁彈性屈曲的微分控制方程為[13]：

2 裂紋損傷表征方法

本文采用Euler-Bernoulli梁理論模型研究梁高相對較小薄曲梁，若裂紋損傷發(fā)生在梁單一側(cè)，則損傷截面中性軸與無損傷截面偏移量較??；若裂紋損傷均勻發(fā)生在梁上下兩側(cè)，則損傷截面中性軸與無損傷截面完全重合。因此，本文研究的曲梁無裂紋和有裂紋橫截面僅高度不同，二者中性軸考慮為處于相同位置。曲梁中的裂紋損傷，使得梁截面產(chǎn)生弱化、梁的抗彎剛度衰減。本研究采用裂紋的截面損傷缺陷比擬方法[15]，裂紋處的截面剛度為：

3 網(wǎng)格自適應(yīng)細分加密

有限元計算存在相比當前網(wǎng)格解答具有更高收斂階的超收斂點[19]，利用超收斂點結(jié)合單元拼片、高階形函數(shù)插值技術(shù)，可以提高當前有限元解的精度，得到全域的超收斂解[15,20-21]。本文對于圓弧曲梁的彈性屈曲問題，求得當前網(wǎng)格下屈曲模態(tài)(位移)的有限元解后，利用有限元后處理超收斂拼片恢復(fù)方法，得到屈曲模態(tài)的超收斂解：

利用屈曲模態(tài)誤差估計，網(wǎng)格可以進行優(yōu)化處理來降低和控制屈曲模態(tài)的誤差，達到預(yù)設(shè)的解答精度。本文方法對每個有限元單元e上的振型誤差進行判斷，如果誤差控制式(12)不滿足，則表明該單元上屈曲模態(tài)解答的誤差過大，需要通過進行網(wǎng)格優(yōu)化處理，本文采用單元均勻細分加密的h型網(wǎng)格自適應(yīng)方式來增加模型自由度、降低單元上解答的誤差[15]。當前單元細分生成的新單元長度與目前誤差和單元階次相關(guān)，即利用當前誤差可以估計新單元長度：

4 數(shù)值算例

本文方法已經(jīng)編制相應(yīng)的Fortran 90語言程序代碼，程序開發(fā)實施基于Microsoft Visual Studio和Intel Visual Fortran編程軟件平臺。本節(jié)給出求解具有代表性的多種圓弧曲梁彈性屈曲數(shù)值算例，對網(wǎng)格自適應(yīng)劃分以及彈性屈曲解答的收斂性進行了討論，對曲梁夾角、損傷位置、數(shù)目、大小等因素影響彈性屈曲荷載和屈曲模態(tài)進行了分析，檢驗了網(wǎng)格自適應(yīng)劃分的有效性。本節(jié)所有算例均采用3次元，初始網(wǎng)格采用2個單元，給定的初始誤差限為。

4.1 圓弧曲梁屈曲分析的網(wǎng)格優(yōu)化

考慮一兩端簡支的圓弧曲梁，梁的幾何和物理參數(shù)如下：

本文方法對該曲梁彈性屈曲進行求解，得到屈曲荷載和屈曲模態(tài)的解答。文獻[23]采用層狀組合梁模型及有限元方法、文獻[7]采用理論模型解析法分別對該梁彈性屈曲進行分析，得到無裂紋損傷梁屈曲荷載值。為進行對比分析，將上述各方法求得的屈曲荷載值均列于表2。同時，為討論網(wǎng)格數(shù)目對解答收斂的影響，本研究采用常規(guī)有限元法在4個、5個、6個單元(稀疏均勻網(wǎng)格)以及25個、50個、100個單元(密集均勻網(wǎng)格)上分別進行求解，得到屈曲荷載值；可知隨著單元增多，解答趨于穩(wěn)定，在接近100個單元時得到收斂解答。使用本文自適應(yīng)有限元方法進行求解，得到僅為16個單元的優(yōu)化網(wǎng)格和在此網(wǎng)格下的收斂解答，該屈曲荷載解答(64.406 kN/m)與理論模型解析解(64.966 kN/m)具有很好的吻合度。需要說明的是，組合梁模型使用4個、5個、6個單元得到的解答與解析解相差較大，為梁模型的層狀建模引起。

表 2 單元數(shù)目與屈曲荷載結(jié)果收斂性Table 2 Convergence for number of elements and buckling loads results

圖3給出本文方法利用自適應(yīng)網(wǎng)格求解得到的屈曲模態(tài)解答，為方便直觀顯示和分析，屈曲模態(tài)結(jié)果均進行歸一化處理(令最大模態(tài)值為1)。由于該梁的物理性質(zhì)和幾何形態(tài)具有左右對稱性，得到圖3(a)所示的屈曲模態(tài)亦為左右對稱形式；同時，水平坐標軸上給出了本文的自適應(yīng)網(wǎng)格(單元端節(jié)點)分布，可知單元分布適應(yīng)模態(tài)變化、同樣具有左右對稱性，且在兩端邊界部分使用了相對細密的網(wǎng)格。算法自動優(yōu)化出非均勻網(wǎng)格，在屈曲模態(tài)變化平緩區(qū)域使用稀疏網(wǎng)格、在屈曲模態(tài)變化劇烈處采用了相對細密的網(wǎng)格，避免了全域使用一致細密網(wǎng)格的冗余性。為了檢驗屈曲模態(tài)解答的精確性，本研究使用常規(guī)有限元法2500個單元(高密集均勻網(wǎng)格)求解得到高精度屈曲模態(tài)作為解析解；該高精度屈曲模態(tài)形態(tài)與圖3(a)相同，這里不再給出。圖3(b)所示為使用本文自適應(yīng)方法的屈曲模態(tài)與高精度屈曲模態(tài)的差值曲線分布，該差值曲線最大值2.17×10-5小于預(yù)設(shè)誤差限10-4，驗證了本文方法求解解答的精確性。

常規(guī)有限元采用4個、5個、6個、25個、50個、100個單元求解得到的屈曲模態(tài)解答與高精度屈曲模態(tài)解答的差值分布分別如圖4所示，并在圖中水平坐標軸上給出各網(wǎng)格分布。各單元下差值在全域上的最大值分別為1.23×10-1、1.23×10-1、4.12×10-3、7.70×10-3、1.62×10-3、2.63×10-6，可見隨著網(wǎng)格的均勻加密，屈曲模態(tài)解答的誤差呈現(xiàn)逐漸降低的趨勢；直到提供足夠多的單元(100個單元)時，才接近獲得滿足預(yù)設(shè)誤差限10-4的解答。本文自適應(yīng)方法僅采用16個單元的非均勻分布網(wǎng)格，即可避免高密度均勻分布網(wǎng)格的單元冗余性，提高了計算效率。

圖 3 自適應(yīng)網(wǎng)格下屈曲模態(tài)Fig. 3 Buckling mode on adaptive refinement mesh

4.2 變曲梁夾角下的屈曲荷載

為檢驗本文方法求解不同幾何形式圓弧曲梁彈性屈曲的適用性，本研究對不同圓弧夾角曲梁(淺梁、深度梁)進行分析?？紤]一兩端簡支的圓弧曲梁，梁的幾何和物理參數(shù)如下：

4.3 彈性屈曲的裂紋損傷位置影響

圖 4 加密網(wǎng)格下屈曲模態(tài)收斂情況Fig. 4 Convergence of buckling modes on refined meshes

表 3 兩端簡支曲梁不同夾角下屈曲荷載值Table 3 Buckling loads of of curved beam with hinged-hinged supports under different subtended angles

表 4 不同裂紋損傷位置下曲梁屈曲荷載值Table 4 Buckling loads of curved beam with crack damage at different locations

圖 5 不同裂紋損傷位置下圓弧曲梁彈性屈曲Fig. 5 Elastic buckling of circularly curved beam with different locations of crack damage

4.4 彈性屈曲的裂紋損傷大小影響

本節(jié)分析圓弧曲梁裂紋損傷大小(深度)對彈性屈曲的影響，固定裂紋損傷位于跨中(θc=π/2)，分析裂紋損傷大小hc/h值為0.1、0.2、0.3、0.4、0.5時的彈性屈曲。使用本文方法分別進行求解，得到彈性屈曲荷載解答，結(jié)果列于表5?？芍S著裂紋損傷hc增大(hc/h=0.1 →0.5 )，qˉc逐漸降低。

表5 不同裂紋損傷大小下曲梁屈曲荷載值Table5 Buckling loadsof curved beam with crack under different magnitudes

為更直觀顯示和分析變化趨勢，利用上述不同裂紋損傷大小下屈曲荷載結(jié)果繪制圖6(a)所示變化曲線。可見屈曲荷載隨裂紋損傷加深，出現(xiàn)準線性的降低趨勢；在裂紋擴展到梁高的一半(hc/h=0.5)時，屈曲荷載承載力降低約10%。圖6(b)所示為不同裂紋損傷大小時各屈曲模態(tài)與無損傷屈曲模態(tài)的差值曲線，可以看出裂紋損傷大小對屈曲模態(tài)變化幅度具有重要影響，裂紋損傷程度越大，越容易誘發(fā)屈曲模態(tài)的大幅度變化。同時，水平坐標軸上給出本文方法求解hc/h=0.5時的自適應(yīng)最終非均勻網(wǎng)格。

4.5 多裂紋損傷下屈曲模態(tài)和網(wǎng)格劃分

為分析多裂紋損傷對曲梁彈性屈曲的影響，采用圖7 所示的兩端簡支圓弧曲梁，該曲梁夾角為θ=π，其余基本幾何和物理參數(shù)同式(12)。該曲梁包含3條裂紋損傷，考慮I(裂紋損傷沿曲梁均勻分布)、II(裂紋損傷集中于曲梁端部)、III(裂紋損傷集中于曲梁中部)等3種工況，各工況中裂紋分布角坐標值如表6所示，裂紋分布如圖7所示。

使用本文方法計算該曲梁多裂紋工況的結(jié)果如表6所示，工況III的屈曲荷載值最小，即多裂紋損傷越集中于跨中，越容易誘發(fā)彈性失穩(wěn)，這與上文單一裂紋損傷位于跨中易于導致彈性失穩(wěn)的結(jié)論一致。

圖8所示為多裂紋工況各屈曲模態(tài)與無損傷屈曲模態(tài)的差值曲線，可以看出各裂紋損傷所在局部區(qū)域?qū)ηB(tài)變化有重要影響，多裂紋損傷缺陷導致屈曲模態(tài)劇烈變化。同時，本文算法在各裂紋損傷附近區(qū)域使用了相對密集的網(wǎng)格，用于適應(yīng)裂紋損傷引起屈曲模態(tài)的變化，形成優(yōu)化的非均勻網(wǎng)格、確保解答的可靠性。

表6 多裂紋損傷不同位置下曲梁屈曲荷載值Table 6 Buckling loadsof curved beam with multiple cracks at different locations

圖6 不同裂紋損傷大小下圓弧曲梁彈性屈曲Fig.6 Elastic buckling of circularly curved beam with different magnitudes of crack damage

圖7 各工況多裂紋分布示意圖Fig.7 Diagram of casesfor multiple cracksdistributions

圖8 多裂紋損傷下網(wǎng)格分布和屈曲模態(tài)差Fig.8 Mesh distribution and buckling modes of circularly curved beam with multiple cracks

5 結(jié)論

本文建立了圓弧形曲梁裂紋的截面損傷缺陷比擬方案，實現(xiàn)裂紋的大小(深度)、位置、數(shù)目的模擬；針對含裂紋損傷圓弧曲梁彈性屈曲，引入有限元網(wǎng)格自適應(yīng)分析方法，得到了優(yōu)化的網(wǎng)格和滿足預(yù)設(shè)誤差限的高精度屈曲荷載和屈曲模態(tài)解答。經(jīng)數(shù)值算例檢驗，本文方法對淺梁、深度梁等各類變化幾何形式圓弧曲梁的彈性屈曲求解具有良好適用性，解答與解析解具有較高吻合度。研究發(fā)現(xiàn)：隨著裂紋損傷增大、接近跨中，均會不同程度降低屈曲荷載，越容易誘發(fā)彈性失穩(wěn)。裂紋損傷將誘發(fā)屈曲模態(tài)變化，本文自適應(yīng)算法可劃分出非均勻網(wǎng)格，在裂紋附近區(qū)域使用了相對密集的網(wǎng)格來適應(yīng)裂紋損傷引起屈曲模態(tài)的變化。