顆粒離散元法在框架結構中的應用

齊 念

(中冶華天南京工程技術有限公司，江蘇 南京 210019)

顆粒離散元法(Discrete Element Method，DEM)是20世紀70年代由美國學者Cundall教授[1]首先提出的，該方法最初應用于土體、巖石等材料的力學性能分析。其基本思想是把研究對象離散為剛性單元的集合，單元與單元之間通過彈簧連接，接觸力與接觸位移之間的關系構成了DEM的接觸本構模型。各個剛性單元滿足牛頓運動方程，DEM法分析時允許單元發(fā)生相對運動，不要求必須滿足位移連續(xù)和變形協(xié)調條件，尤其適合非線性問題的研究。目前DEM法已廣泛用于計算散體力學領域[2,3]，如粉末的加工、研磨，藥品的運輸?shù)取１疚耐ㄟ^理論推導和數(shù)值分析，將DEM法推廣用于框架結構幾何非線性大變形問題。然后通過一個經典算例的大變形行為進行模擬和分析并與其他方法的計算結果進行比較。結果表明，DEM方法適用于框架結構大變形問題的分析。

1 顆粒DEM法基本原理

1.1 單元之間的接觸力分析

顆粒離散元，以顆粒為基本研究對象，將物體離散作為具有代表性的數(shù)個單元，利用顆粒流模型構建物體的力學性質。以平面問題為例，假定單元形狀為圓盤，單元之間的接觸采用柔性接觸，即允許接觸處產生重疊，圖1為任意兩個相鄰圓A和圓B間的接觸模型，C點為接觸中心。

當圓A與圓B之間產生接觸時，接觸重疊量Un可表示為：

Un=R[A]+R[B]-D

(1)

其中，R[A]，R[B]，D分別為圓A、圓B的半徑及兩圓之間的中心距。

在離散元方法中，相鄰圓之間的相互接觸作用力是通過接觸彈簧來實現(xiàn)[4]。如圖1所示，在接觸平面內，接觸力F可分解成垂直于接觸面的法向力Fn和平行于接觸面的切向力Fs，即：

F=Fn+Fs=Fn×n+Fs×s

(2)

其中，F(xiàn)n，F(xiàn)s分別為法向力和切向力的大小；n，s分別為接觸平面的法向和切向單位矢量。

此外，接觸點C處因兩圓盤發(fā)生相對轉動還會引起大小為Mz的接觸力矩，這與散粒體離散元中的接觸模型不同。

將接觸力分別向坐標軸進行投影，得：

Fi=Fn×ni+Fs×si(i=1,2)

(3)

根據(jù)力系平移定理，將Fi移至各單元的中心處，有：

(4)

(5)

1.2 接觸力和接觸力矩的計算

接觸力和接觸力矩的計算一般是通過力與位移之間的關系來表達。對于法向接觸力Fn：

ΔFs=-KnΔUn

(6)

切向接觸力Fs，用增量的形式進行計算：

ΔFs=-KsΔUs

(7)

ΔUs=VsΔt

(8)

其中，Kn,Ks分別為彈簧法向剛度和切向剛度系數(shù);ΔUs為相對切向位移增量;Δt為計算時步;Vs為接觸點C處的切向速度。

對于接觸力矩M3的計算，也采用增量形式來表達：

(9)

其中，Kθ為彈簧轉動剛度;Δθ為相對轉角增量;ω[B],ω[A]分別為圓A和圓B的轉動角速度。

然后，按以下方式進行接觸作用力的更迭計算：接觸形成時，對Fn,Fs,M3初始化為零，然后進行下一時步的計算，將由相對位移引起的接觸力增量、相對轉角引起的接觸力矩增量累加到現(xiàn)值上，即：

(10)

1.3 運動控制方程

(11)

對于運動方程的求解，DEM通常采用顯示積分算法，本文采取中心差分格式。

2 算例驗證

平面正方形剛架大變形分析。

如圖2所示，一平面正方形剛架在對邊受拉力作用下發(fā)生大位移大轉動行為，采用DEM法模擬并與文獻[5]中采用橢圓積分計算的結果進行對比。剛架各構件的材料和尺寸均相同，相關幾何及物理參數(shù)為：L=10 cm，A=0.5 cm2，I=0.010 4 cm4，E=1.6×107N/cm2，材料為線彈性。DEM建模時，將剛架各邊構件均離散成11個半徑為1.0 cm的圓盤。計算時步Δt=1.0×10-3s，阻尼系數(shù)ζ取0.6。

表1為用DEM法計算得到的剛架在對邊受拉情形下的荷載—變形結果，定義荷載無量綱因子：λ=PL2/EI。與文獻中KJELL采用橢圓積分分析進行對比發(fā)現(xiàn)，兩種方法計算的結果相當接近，最大誤差不超過0.9%，表明本文方法具有較高的數(shù)值精度。

表1 平面正方形剛架對邊受拉下的荷載—變形結果

圖3為不同荷載大小作用下平面正方形剛架的最終變形圖。從圖3中看出，荷載值越大，結構變形越顯著，用DEM分析時很容易獲得結構的最終變形，因為它計算時會自動考慮幾何非線性的影響，而不去區(qū)分是大變形還是小變形，因此求得的解就是它真實的狀態(tài)。

3 結語

本文將散體力學中應用廣泛的顆粒DEM法推廣用于求解平面框架結構幾何非線性大變形問題。算例分析結果表明：本文方法計算精度較高，對框架結構幾何非線性大變形(大位移、大轉動)問題能很好的處理。

DEM法是以牛頓第二定律為力學基礎，在求解結構大變形問題時，不用組集單元的剛度矩陣，不需要迭代，流程相對簡單。另外，該方法不要求滿足位移連續(xù)和變形協(xié)調條件，突破了傳統(tǒng)非線性有限元在結構非線性問題中遇到的困難，為今后求解結構大變形這類問題提供了一種新的研究思路和方法。

[1] Cundall P A,Strack ODL.A discrete numerical model for granular assemblies[J].Geotechnique,1979,29(1):47-65.

[2] 徐 泳，孫其誠，張 凌.顆粒離散元法研究進展[J].力學進展，2003，33(2)：251-260.

[3] 唐志平.激光輻照下充壓柱殼失效的三維離散元模擬[J].爆炸與沖擊，2001，21(1)：1-7.

[4] Itasca Consulting Group Inc.,Particle Flow Code in 2-Dimensions(PFC2D),Version 3.10,Minneapolis,Minnesota,2004.

[5] KJELL MATTIASSON.Numerical results from large deflection beam and frame problems analysed by means of elliptic integrals[J].International journal for numerical methods in engineering,1981,17(1):145-153.