吳怡凡, 沈建, 金世杰, 夏晉
(1.浙江大學結構工程研究所,杭州 310058;2.浙江大學建筑設計研究院有限公司,杭州 310028)
隨著有限元技術的不斷發展,在大型土木工程結構的模擬分析中,為了滿足整體結構的力學響應,并反應關鍵節點的塑性發展,多尺度建模成為了一種常用策略,即在整體結構的非關鍵部位使用計算效率較高的宏觀單元,而在需要精細分析的關鍵節點部位使用精度較高的微觀單元,并通過一定的連接方式將兩種類型的單元耦合起來,以達到既兼顧計算效率又保證精度的目的[1-3]。多點約束法是目前土木工程領域進行多尺度界面連接最常用的一種方法,它的核心思想是通過建立不同尺度界面節點自由度間的協調方程使兩界面達到有效的耦合[4]。在過去的研究中,學者們往往采用三種理論來構建多點約束方程,分別為:位移協調、能量協調以及位移協調與能量協調相結合。林旭川等運用平截面假定列出位移協調方程,采用MSC.MARC模擬了鋼筋混凝土框架多尺度模型在地震下的受力情況,并通過試驗驗證[5]。Yue和Wu等[6,7]利用做功相等能量協調原理,分別通過ABAQUS和ANSYS建立了鋼筋混凝土框架多尺度模型,并與實驗比較,結果表明該方法可以準確模擬結構局部破壞。但在鋼筋混凝土結構的數值模擬分析中,由于材料非線性、幾何非線性等因素,基于單一理論的多點約束法存在一定的局限性[8-10]。基于位移協調與能量協調相結合的多點約束法可以避開兩種單一理論的缺點得到更好的模擬結果,但需要插入數量龐大的約束方程,建模前處理繁瑣復雜[11]。
文中提出了一種基于連續分布耦合理論的多點約束法以實現不同尺度單元間的連接,針對鋼筋混凝土節點,使用ABAQUS比較分析了連續分布耦合與位移協調、能量協調、位移協調與能量協調相結合4種多尺度界面連接方法的優劣性。
文中以梁-實體單元的連接為例介紹不同理論實現多尺度界面連接的基本原理和方法,梁-實體單元連接如圖1所示。

圖1 梁-實體單元連接
連接界面上實體單元共有n個節點,每個節點有x、y、z3個自由度;梁單元與實體單元通過節點B連接,節點B有xB、yB、zB、θxB、θyB、θzB6個自由度。根據梁、實體單元可能產生的4種位移模式,可知界面節點間產生的變形協調如圖2所示。將精細實體單元節點的3個自由度分別用宏觀梁單元節點的6個自由度表示,一共形成3n個位移協調方程如式(1)。


圖2 梁-實體單元多尺度界面變形協調

式中,bxi、byi為實體單元節點i在X、Y坐標軸上的坐標。
在軸力、彎矩、剪力、扭矩的作用下界面節點變形示意如圖3所示,根據虛功原理可列出做功相等方程如式(2)。

圖3 梁-實體單元多尺度界面受力變形

式中,Fj為作用在宏觀單元B節點方向的力;σi1,Fj為精細單元i節點處由Fj引起的軸向應力;τi2,Fj、τi3,Fj為精細單元i節點處由Fj引起的切向應力;u1S、u2S、u3S為精細單元i節點處軸向及切向的位移。
結合材料力學公式,并用形函數來表示節點位移X={N}{X},Y={N}{Y},Z={N}{Z},代入四節點四邊形單元形函數對積分方程求解,可得到約束方程如式(3)。

式中,Cxi、Cyi、Czi、Cθxi、Cθyi為與實體單元截面尺寸及連接界面點數個數相關的系數。
基于位移協調的多點約束法在切向存在過約束的缺陷,而基于能量協調的多點約束法在轉角方向存在材料非線性和幾何非線性的缺陷。
因此將兩種方法在不同自由度上進行方程聯立,在切向上使用基于能量協調的多點約束方程,而在軸向和轉角方向上使用基于位移協調的多點約束方程,聯立得到的方程組如式(4)。

連續分布耦合示意如圖4所示,連續分布耦合允許約束面上各節點發生相對移動,通過對約束面上各節點的運動進行加權處理,使在此區域上受到的合力和合力距與施加在參考點上的力和力矩相等效。因此,使用連續分布耦合理論建立的多尺度模型在理論上可以避免位移協調理論產生的切向過約束。ABAQUS連續分布耦合有4種加權系數計算方式,包括一致、線性、二次方和三次方,計算公式如式(5)~式(8)。經過反復試算文中采用一致加權系數計算方式對界面節點進行加權處理。

圖4 連續分布耦合

式中,ωi為界面節點i的加權系數;ri為界面節點i到參考點的徑向距離;ω0為最遠點的徑向距離。
分析模型如圖5所示,構件高2m,截面為邊長0.4m的正方形,多尺度界面距離底部1m。圖5(a)為精細化模型,圖5(b)~圖5(e)分別為基于位移協調、能量協調、位移協調與能量協調相結合、連續分布耦合理論建立的多尺度模型。分別命名為位移法模型、能量法模型、結合法模型和DC法模型。

圖5 鋼筋混凝土節點梁-實體單元多尺度連接模型
對于實體單元部分鋼筋混凝土采用分離式建模處理,其建模原理是將混凝土和鋼筋分別建模,混凝土采用C3D8R單元,鋼筋采用T3D2單元,并利用嵌固作用將鋼筋嵌入到混凝土中,實現鋼筋與混凝土的良好結合。對于梁單元部分的鋼筋混凝土采用整體式模型,其建模原理是在同一位置建立兩個共節點B31梁單元,梁單元1為混凝土,梁單元2用箱形截面等效模擬鋼筋,其建模方式如圖6所示。混凝土采用塑性損傷本構模型,鋼筋采用無屈服點的彈塑性加工硬化兩折線模型。模型網格尺寸為50mm。

圖6 箱形截面等效配筋
在分析模型的底部設置固定約束,在頂部設置加載點通過施加位移控制加載,豎向位移從0逐漸增長到10mm。各模型在塑性階段混凝土單元的應力云圖如圖7所示,各模型在破壞階段混凝土受壓損傷云圖如圖8所示。通過觀察加載過程中應力云圖的變化可以發現,在軸力工況下,隨著豎向位移的增大,各模型的單元應力狀態逐漸由彈性進入塑性,其中構件的兩端和中部最先進入塑性,混凝土在達到極限抗壓強度后發生破壞,應力下降。在加載過程中,構件的中部逐漸形成損傷累積區并最終發生破壞。對比發現,除了位移法模型,其余各多尺度模型均能較好地反應構件在彈塑性階段的應力和損傷分布,印證了基于位移協調理論構建的多尺度模型在僅受軸力作用時,多尺度界面存在側向過約束。

圖7 軸力工況下塑性階段各模型混凝土壓應力云圖

圖8 軸力工況下各模型混凝土受壓損傷累積云圖
提取各模型在加載過程中加載點的位移和對應的軸力,作出頂點位移-軸力曲線如圖9所示。5組模型加載結束時頂點軸力分別為1510、1550、1610、1500、1530kN。相比精細化模型,各多尺度模型相對誤差分別為2.65%、5.88%、0.2%、1.12%。

圖9 頂點位移-軸力曲線
按照象限編號法則,由內到外,在各模型多尺度界面上選取16個特征單元,如圖10所示。當頂點豎向位移加載到2.15mm時構件處于塑性階段,分別提取各模型16個特征單元的Min.principle應力,統計精細化模型與多尺度模型在特征單元處的應力值并計算相對誤差。觀察發現,平均相對誤差從大到小分別為位移法模型10.1%、能量法模型1.8%、結合法模型1.6%、DC法模型4.41%。位移法模型的13、14號特征單元相對誤差最大,達到32.1%。能量法模型的1、2、3、4號特征單元相對誤差最大,達到2.5%。結合法模型的1、2、3、4號特征單元相對誤差最大,達到2.1%。DC法模型的1、2、3、4號特征單元相對誤差最大,達到5.1%。在軸力工況下,位移法模型的相對誤差率高于其他多尺度模型,印證了位移協調法在僅受軸力作用下存在側向過約束局限性的分析。DC法模型的平均相對誤差為4.4%,控制在5%之內,印證了使用連續分布耦合約束可以有效改善側向過約束降低誤差率。

圖10 多尺度界面特征單元編號
在分析模型的底部設置固定約束,在頂部設置加載點通過施加繞X軸的彎矩荷載控制加載,彎矩從0逐漸增長到40kN·m。
各模型在塑性階段混凝土單元的應力云圖如圖11所示,各模型在破壞階段混凝土受拉損傷云圖如圖12所示。通過觀察加載過程中應力云圖的變化可以發現,混凝土最大拉、壓應力分別沿受拉、受壓側邊緣分布。隨著不斷加載,構件在受拉側,沿柱身不斷形成損傷累積區并最終發生破壞。對比各模型可以發現,能量法模型相比精細化模型應力分布和損傷累積差別都較大。位移法模型和DC發模型與精細化模型在應力分布上吻合度比較高,但在多尺度界面附近的損傷累積有所差別。結合法模型與精細化模型吻合度極高,能很好地反應構件的應力分布及損傷累積。

圖11 彎矩工況下塑性階段各模型混凝土軸向應力云圖

圖12 彎矩工況下各模型混凝土受拉損傷累積云圖
提取各模型在加載過程中加載點繞X軸的彎矩和對應繞X軸的轉角,作出頂點轉角-彎矩曲線如圖13所示。5組模型加載結束時頂點轉角分別為0.00718、0.00923、0.01768、0.00676、0.00493。相比精細化模型,各多尺度模型相對誤差分別為28.55%、91.55%、5.85%、31.94%。

圖13 頂點轉角-彎矩曲線
當頂點彎矩加載到36.8kN·m時構件處于塑性階段,分別提取各模型16個特征單元的S33應力,統計精細化模型與多尺度模型在特征單元處的應力值并計算相對誤差,觀察發現,平均相對誤差從大到小分別為能量法模型78.5%、位移法模型27.7%、DC法模型6.6%、結合法模型0.7%。能量法模型的3、4號特征單元相對誤差最大,達到239.2%。位移法模型的1、2、3、4號特征單元相對誤差最大,達到50%。DC法模型的1、2、3、4號特征單元相對誤差最大,達到7.1%。結合法模型的9、10、11、12號特征單元相對誤差最大,達到2.3%。在彎矩工況下,能量法模型的相對誤差率遠大于其他多尺度模型,印證了能量協調法在彎矩作用下局限性的分析,既當構件進入塑性階段后,彎矩產生的正應力呈非線性分布,基于能量協調的多點約束方程不再適用。
在分析模型的底部設置固定約束,在頂部設置加載點通過施加沿X軸的剪力荷載控制加載,剪力從0逐漸增長到20kN。
各模型在塑性階段混凝土單元的應力云圖如圖14所示,各模型在破壞階段混凝土受拉損傷云圖如圖15所示。通過觀察加載過程中應力云圖的變化可以發現,構件在端部剪力的作用下,底部混凝土首先進入塑性,并在拉壓兩側產生應力集中。隨著不斷加載,構件底部發生破壞,塑性變形快速增大。對比各模型可以發現,位移法模型與結合法模型的應力分布、損傷累積計算結果均較貼合精細化模型所得結果。能量法模型在多尺度界面出現了明顯應力突變,損傷積累也與精細化模型有較大差別。

圖14 剪力工況下塑性階段各模型混凝土軸向應力云圖
提取各模型在加載過程中加載點沿X軸的剪力和對應沿X軸的位移,作出頂點位移-剪力曲線如圖16所示。5組模型加載結束時頂點位移分別為1.92、1.89、1.89、1.87、1.89mm。相比精細化模型,各多尺度模型相對誤差分別為1.22%、1.22%、2.25%。

圖16 頂點位移-剪力曲線
當頂點剪力加載到20kN時構件處于塑性階段,分別提取各模型16個特征單元的S33應力,統計精細化模型與多尺度模型在特征單元處的應力值并計算相對誤差見表1,觀察發現,在剪力荷載工況下,平均相對誤差從大到小分別為能量法模型142.2%、位移法模型23.7%、DC法模型15.2%、結合法模型3.6%。能量法模型的1號特征單元相對誤差最大,達到341.6%。位移法模型的1、2、3、4號特征單元相對誤差最大,達到33.3%。DC法模型的1、2、3、4號特征單元相對誤差最大,達到45.1%。結合法模型的6、7、8號特征單元相對誤差最大,達到5.4%。在剪力工況下,能量法模型的相對誤差率遠大于其他多尺度模型,是因為在端部剪力的作用下構件底部產生彎矩,混凝土進入塑性后彎矩產生的正應力呈非線性分布,基于能量協調的多點約束方程不再適用。

表1 鋼筋混凝土節點平均應力相對誤差統計
文中針對鋼筋混凝土節點,基于位移協調、能量協調、位移協調與能量協調相結合以及基于連續分布耦合理論實現多尺度界面的連接并進行算例分析,所得具體研究成果如下:
(1) 基于位移協調實現多尺度模型界面連接時,除了在僅受軸力這種工況外,可以較為合理地反應構件的應力分布和損傷累積。
(2) 基于能量協調實現多尺度模型界面連接時,除了在軸力工況下得到的結果較為精確外,在剪力工況、彎矩工況下,由于材料進入塑性,彎矩產生的正應力分布公式不再成立,約束方程不再適用,塑性區域單元迭代運算逐漸出現誤差累積。
(3) 基于位移協調和能量協調相結合實現多尺度模型界面連接時,在各工況下均可合理地反應構件應力分布和損傷累積,得到的結果也最接近精細化模型,但需要輸入大量約束方程,建模前處理較為復雜。
(4) 基于連續分布耦合理論實現多尺度模型界面連接時,在各工況下均可合理地反應構件的應力分布和損傷累積,相比位移協調理論不僅可以有效改善側向過約束,在軸力工況下得到較為精確的結果,還能避免輸入大量的約束方程等繁瑣的建模前處理。