基于各向異性損傷本構(gòu)模型的砌體結(jié)構(gòu)抗震分析

富秋實， 錢 江

（同濟大學(xué) 土木工程與防災(zāi)國家重點實驗室，上海200092）

1 引 言

砌體結(jié)構(gòu)是土木工程中最古老的結(jié)構(gòu)形式之一，為了有效地保護歷史建筑遺產(chǎn)和人員生命安全，針對既有砌體結(jié)構(gòu)的結(jié)構(gòu)分析手段顯得十分重要。一般，砌體結(jié)構(gòu)分析主要是使用基于經(jīng)驗的分析方法。過去的幾十年研究者對砌體結(jié)構(gòu)的數(shù)值模擬方法進行了大量研究。另一方面，由于砌體材料特殊的各向異性性質(zhì)，采用其他領(lǐng)域（如混凝土、巖石和復(fù)合材料）的數(shù)值模擬方法存在一定困難。所以，針對砌體的各向異性本構(gòu)模型的研究開發(fā)具有重要意義。

在各向異性材料本構(gòu)模型的研究領(lǐng)域，除了直接提出各種各向異性材料的屈服準(zhǔn)則函數(shù)以及損傷準(zhǔn)則函數(shù)［1］外，還有一類策略是采用基于應(yīng)力張量分量線性變換（linear transformation）的方法。使用線性變換的方法便于保持屈服準(zhǔn)則函數(shù)和損傷準(zhǔn)則函數(shù)等的外凸性，這對于在實踐中保持本構(gòu)算法的數(shù)值穩(wěn)定性十分重要。因此，這類方法得到了較為廣泛的關(guān)注。為了有效地利用已有損傷理論，開發(fā)能夠考慮砌體材料各向異性特性的損傷模型，本文采用基于應(yīng)力張量線性變換的方法。

應(yīng)力張量的線性變換方法最早由Sobotka［2］和Boehler等［3］提出。Barlat等［4，5］將各向同性屈服函數(shù)與線性變換后應(yīng)力張量的主應(yīng)力相結(jié)合，用于解決正交各向異性材料的平面應(yīng)力問題，并將這一方法應(yīng)用于三維應(yīng)力狀態(tài)。Karafillis等［6］將該方法推廣為各向同性塑性等效IPE（isotropic plasticity equivalent）理論，將線性變換方法與更加一般化的屈服函數(shù)相結(jié)合，可用于其他類型的材料。Barlat等［7］使用該方法提出了基于應(yīng)力偏量線性變換的各向異性屈服函數(shù)，用于模擬金屬板在軋制方向和橫向的各向異性塑性特性。

砌體材料各向異性本構(gòu)模型的已有研究，除了直接提出各項異性屈服函數(shù)，也有應(yīng)用張量線性變換方法的研究［8，9］。這些研究使用的方法稱為映射張量法（mapped tensor method）。

映射張量法最早由Betten［10］提出，用于解決蠕變相關(guān)問題。其假設(shè)存在一個真實的應(yīng)力空間和一個虛擬的應(yīng)力空間，兩個應(yīng)力空間之間存在映射關(guān)系，由一個四階映射張量定義。存在于真實應(yīng)力空間的各向異性材料的受力狀態(tài)可以映射到存在于虛擬應(yīng)力空間的虛擬各向同性受力狀態(tài)?？梢钥吹?，Betten的映射張量法的思想本質(zhì)上與前文介紹的線性變換方法基本一致。Oller等［11，12］對這一映射張量法進行了進一步的發(fā)展。

本文進一步完善了應(yīng)力張量的線性變換方法，提出了適應(yīng)塑性以及損傷演化程度的動態(tài)線性變換張量，并將其應(yīng)用于砌體材料本構(gòu)模型的開發(fā)。利用該本構(gòu)模型對砌體結(jié)構(gòu)的振動臺試驗進行了模擬。模擬結(jié)果驗證了本文提出的彈塑性損傷本構(gòu)模型的有效性，以及在進行結(jié)構(gòu)非線性分析方面的優(yōu)越性。

2 動態(tài)變換張量

設(shè)四階動態(tài)線性變換張量為Aσ，（Vk）， ∈｛p，d｝。內(nèi)變量Vk是描述材料內(nèi)部不可逆變化進程的變量，如塑性強化參數(shù)等。上標(biāo)σ和p的含義為用于塑性本構(gòu)的應(yīng)力轉(zhuǎn)換張量，σ和d的含義為用于損傷本構(gòu)的應(yīng)力轉(zhuǎn)換張量。不同于已有的Oller和 Pelà模型［9，11，12］，本文并沒有應(yīng)變轉(zhuǎn)換張量。但為了與文獻的應(yīng)變轉(zhuǎn)換張量相區(qū)別，下文將依然保留上標(biāo)中的σ。

根據(jù)應(yīng)力轉(zhuǎn)換張量的定義，有

式中σ＊為變換后的應(yīng)力張量。設(shè)未知其具體表達的各向異性屈服／損傷面函數(shù)F為

設(shè)已知形式的各向同性本構(gòu)的屈服函數(shù)和塑性勢函數(shù)表達式為

此時根據(jù)定義顯然有

即得到了未知函數(shù)F的隱式表達。另一方面，由于不同方向各向異性強化的特性記入變換張量Aσ，（Vk）中，因此隨著強化參數(shù)Vk的演化，在F＊呈現(xiàn)各向同性強化的同時，F(xiàn)可以呈現(xiàn)各向異性強化。

3 各向異性彈塑性損傷模型建立

3.1 塑性演化

在以往研究中，對于混凝土等摩擦性材料的屈服函數(shù)一般采用 Mohr-Coulomb或Drucker-Prager屈服函數(shù)，但與試驗結(jié)果并不是很接近。本節(jié)同時采用Rankine和Drucker-Prager形式的雙屈服準(zhǔn)則，模擬砌體材料受拉及受壓屈服特性，基于動態(tài)線性變換張量，建立針對砌體材料的正交各向異性塑性修正方法。

（1）受拉Rankine形式屈服面

各向同性受拉屈服函數(shù)采用Rankine形式為

引入四階受拉塑性轉(zhuǎn)換張量Aσ，pt，則可由式（1～3）得到受拉Rankine形式屈服面的各向異性隱式表達。

（2）受壓Drucker-Prager形式屈服面

各向同性受壓屈服函數(shù)采用Drucker-Prager形式為

式中α為控制雙軸受壓強度特性的參數(shù)，單軸受壓狀態(tài)下有效應(yīng)力閾值fc＊簡化為硬化參數(shù)的線性函數(shù)。

引入四階受壓塑性轉(zhuǎn)換張量Aσ，pc，可由式（1～3）得到受壓Drucker-Prager形式屈服面的各向異性隱式表達。

3.2 損傷演化

為了反映在受拉和受壓狀態(tài)下材料的力學(xué)性能差異，損傷本構(gòu)模型需要將應(yīng)力張量σ［13］或者有效應(yīng)力張量σ［14］分解為正負(fù)張量兩部分。在基于應(yīng)變的有限元分析中，應(yīng)力張量σ在迭代收斂前是未知量，基于應(yīng)力張量的分解會極大地增加計算量和數(shù)值計算難度。另一方面，由于有效應(yīng)力張量σ和彈性應(yīng)變張量εe線性相關(guān)，基于有效應(yīng)力張量的分解更利于建立應(yīng)變的顯式表達。因此，本文采用有效應(yīng)力張量的分解方法，建立損傷演化規(guī)則：

此外，損傷修正部分只需考慮彈性損傷力學(xué)模型，更新應(yīng)力張量和損傷變量，而不需要考慮塑性變形的影響。因此，為了方便書寫，下文的理論推導(dǎo)部分將省略彈性應(yīng)變張量εe的上標(biāo)e。

（1）受拉損傷演化

受拉狀態(tài)下的各向同性損傷面函數(shù)采用Rankine形式：

τ＋＊具有應(yīng)變張量的量綱，稱為受拉等效應(yīng)力。為方便計算，損傷各向同性受拉準(zhǔn)則函數(shù)寫為

r＋＊為受拉各向同性損傷面函數(shù)對應(yīng)的等效應(yīng)力閾值。

（2）受壓損傷演化

受壓狀態(tài)下的各向同性損傷面函數(shù)采用修正的八面體應(yīng)力組合：

式中 參數(shù)K控制雙軸受壓狀態(tài)有效應(yīng)力損傷閾值與各方向單軸受壓狀態(tài)下有效應(yīng)力損傷閾值的比值，其取值由試驗確定。

（3）II型破壞損傷演化

砌塊-砂漿界面的II型破壞尚難以用損傷變量d＋和d－來描述，其滑動破壞機制也與受拉及受壓開裂有著明顯的區(qū)別。根據(jù)內(nèi)變量理論的概念，不同機制的不可逆過程應(yīng)該由不同的內(nèi)變量來表示。因此，本節(jié)在前文理論的基礎(chǔ)上，直接引入一個損傷變量dII，來描述砌體材料灰縫的II型破壞。

首先，引入II型破壞等效應(yīng)力：

式中σ12為沿灰縫滑動的剪應(yīng)力，σ22為垂直于水平灰縫的正應(yīng)力，II為II型破壞對應(yīng)的摩擦角。相應(yīng)的損傷準(zhǔn)則函數(shù)為

式中r0，II為II型破壞等效應(yīng)力的初始閾值，由剪切試驗確定。假定II型破壞損傷變量的顯式表達式為指數(shù)函數(shù)形式：

4 對比驗證

同濟大學(xué)于2011年對某砌體結(jié)構(gòu)進行了振動臺試驗，該試驗原型結(jié)構(gòu)為上海市某臨近地鐵的住宅［15，16］。由于振動臺試驗的能力限制以及砌體結(jié)構(gòu)模型試驗的特殊性，試驗采用局部代表性結(jié)構(gòu)來模擬原型結(jié)構(gòu)，從實際結(jié)構(gòu)中取出三層兩個房間作為代表性結(jié)構(gòu)進行研究。由于代表性結(jié)構(gòu)的質(zhì)量（690t）依然遠(yuǎn)超出振動臺的能力范圍，所以本試驗采用代表性結(jié)構(gòu)1／4的縮尺模型進行振動臺試驗研究。

地震波加載方向為模型的橫墻方向，即X向。模型砌體部分材料采用小磚塊和砂漿，構(gòu)造柱混凝土采用砂漿，鋼筋采用鍍鋅鐵絲，原型結(jié)構(gòu)樓屋面厚度為110mm，1／4縮尺后樓屋面板厚為27.5mm，考慮到施加質(zhì)量塊時固定砂漿會對樓板有一定程度的加強，并且樓板并不是試驗?zāi)M的重點，為使樓板能夠承受附加質(zhì)量引起的合重，樓屋面板厚取40mm，混凝土部分用砂漿，鋼筋采用鍍鋅鐵絲進行模擬。砂漿澆注時用小振搗器振搗和小鐵釬插搗，保證模型澆筑質(zhì)量。模擬砌體部分的小磚塊由MU10磚塊在三個方向上一分為四，即每塊磚切成64塊模型小磚塊，砂漿用M2.5混合砂漿。

4.1 計算模型概況

通過商業(yè)有限元程序ANSYS的Usermat子程序，將前文所建立的各向異性本構(gòu)模型嵌入ANSYS程序中，用來進行砌體結(jié)構(gòu)非線性分析。程序編寫基于平面應(yīng)力狀態(tài)，可用于平面應(yīng)力單元以及殼單元（平面應(yīng)力＋板彎曲）。

利用ANSYS程序進行三維空間建模。砌體墻、樓板及模型底座采用SHELL單元模擬，混凝土梁和柱采用BEAM單元模擬。計算模型共計單元數(shù)3392，節(jié)點數(shù)3494。計算采用的主要材料參數(shù)為x向彈性模量E1＝1806MPa，抗拉強度＝0.6MPa，抗壓強度＝4.8MPa；y向彈性模量E2＝1508MPa，抗拉強度＝0.4MPa，抗壓強度f2－＝4.0MPa。計算模型如圖1所示。

4.2 8度罕遇烈度El-Centro波作用下結(jié)構(gòu)響應(yīng)結(jié)果對比

為了盡可能使計算工況與試驗工況保持一致，計算時的加速度輸入采用相應(yīng)工況的試驗臺面記錄。

（1）樓層加速度對比

表1及圖2為8度罕遇El-Centro波下，加速度試驗與數(shù)值分析包絡(luò)結(jié)果對比。

圖3～圖6為8度罕遇El-Centro波下，加速度試驗與數(shù)值分析加速度時程結(jié)果對比。

圖1 有限元模型Fig.1 FEM model

表1 加速度試驗與數(shù)值分析包絡(luò)結(jié)果對比Tab.1 Comparison of acceleration results obtained from test and simulation

可以看出，除了模型底板處的動力響應(yīng)結(jié)果與試驗偏差略大以外，整體結(jié)構(gòu)的相應(yīng)幅值與試驗結(jié)果比較接近；相位結(jié)果有一定偏差，其原因可能為計算模型未考慮混凝土底梁的非線性行為，但在試驗過程中混凝土底梁進入了非線性狀態(tài)。

圖2 加速度試驗與數(shù)值分析包絡(luò)結(jié)果對比Fig.2 Comparison of acceleration results obtained from test and simulation

圖3 一層加速度反應(yīng)時程結(jié)果對比Fig.3 Comparison of acceleration time history curves of the first floor

圖4 二層加速度反應(yīng)時程結(jié)果對比Fig.4 Comparison of acceleration time history curves of the second floor

圖5 三層加速度反應(yīng)時程結(jié)果對比Fig.5 Comparison of acceleration time history curves of the third floor

圖6 頂板層加速度反應(yīng)時程結(jié)果對比Fig.6 Comparison of acceleration time history curves of the roof

（2）樓層相對位移對比

表2和圖7為8度罕遇El-Centro波作用下，樓層相對位移試驗與數(shù)值分析包絡(luò)結(jié)果對比。兩者結(jié)果除了2層處相差較遠(yuǎn)外，其他位置結(jié)果比較接近。考慮本模型在2層處沒有附加質(zhì)量和剛度突變，計算所得的位移包絡(luò)曲線形狀更加符合結(jié)構(gòu)概念，因此，認(rèn)為計算得到的結(jié)果具有較好的可信度。

（3）損傷分布情況對比

計算結(jié)構(gòu)損傷分布情況如圖8所示，損傷主要分布在模型的底層。試驗結(jié)果顯示，在8度罕遇El-Centro波作用下，試驗?zāi)Ｐ椭饕茐哪Ｊ綖闄M墻1／2高度處起始的斜向階梯型裂縫（試驗過程中能聽到明顯的開裂聲音），相應(yīng)的試驗?zāi)Ｐ土芽p照片和計算模型損傷分布對比情況如圖9所示，可見計算模型的損傷分布與試驗?zāi)Ｐ土芽p開展情況基本相符。同時，縱墻的墻體底部也發(fā)生了較為輕微的水平開裂，試驗?zāi)Ｐ土芽p照片和計算模型損傷分布對比情況如圖10所示，兩者結(jié)果也基本相符。

表2 樓層相對位移試驗與數(shù)值分析包絡(luò)結(jié)果對比Tab.2 Comparison of maximum story relative displacements

圖7 樓層相對位移試驗與數(shù)值分析包絡(luò)結(jié)果對比Fig.7 Comparison of maximum story relative displacements

圖8 計算模型整體損傷分布情況Fig.8 Damage distribution of FEM model

圖9 底層橫墻裂縫／損傷分布對比Fig.9 Crack／damage distribution of bottom lateral wall

5 結(jié) 論

本文進一步完善了應(yīng)力張量的線性變換方法，提出了適應(yīng)塑性以及損傷演化程度的動態(tài)線性變換張量，建立了以該變換張量為基礎(chǔ)發(fā)展各向異性彈塑性本構(gòu)模型以及損傷本構(gòu)模型的一般方法，并將其應(yīng)用于砌體材料本構(gòu)模型的開發(fā)。

基于開發(fā)的砌體各向異性彈塑性損傷本構(gòu)模型，本文對砌體結(jié)構(gòu)的振動臺試驗進行了模擬。模擬結(jié)果進一步驗證了本文提出的彈塑性損傷本構(gòu)模型的有效性，以及在進行結(jié)構(gòu)非線性分析方面的優(yōu)越性，即除了能夠獲取結(jié)構(gòu)的反應(yīng)外，還能實時提供結(jié)構(gòu)中的損傷分布狀態(tài)信息，找出結(jié)構(gòu)的薄弱部位，并據(jù)此合理地設(shè)計結(jié)構(gòu)或進行相應(yīng)的結(jié)構(gòu)修復(fù)。