基于微平面彈塑性應力-應變關系的混凝土本構模型

賈明曉

（華北水利水電學院 土木與交通學院，鄭州 450011）

微平面模型是一種跨尺度的本構模型，它通過在微平面上建立本構關系來描述材料宏觀的力學性質。微平面模型的基本思想可以追溯到1938年Taylor在研究金屬材料塑性性質時提出的金屬滑移理論。Ba?ant等將塑性滑移理論進行了改進和推廣，于1983年首次提出了適用于混凝土的微平面模型［1］，之后經過不斷的改進和發展先后歷經了5代微平面模型（簡稱M1～M5）［2～9］。該模型直接在微平面上建立矢量形式的本構關系，概念簡單并能反映混凝土的各種復雜受力行為，相對于張量形式的本構模型具有很大的優勢。

目前為止微平面模型得以實際應用的主要是M2模型，Ba?ant等［3］提出了 M2 模型，Carol等［9］建立了基于顯式算法的 M2模型，Cofer和 Kohut［10］在動力顯式有限元程序中實現了 M2模型，Ba?ant等［11］使用這一模型成功模擬了雙邊切口和開孔梁的試驗。M3及以后幾個版本的微平面模型提出了應力邊界的概念，這些改進的模型能夠擬合更多的材料受力特征，但是模型復雜程度增加，物理概念不夠明確，且存在著模型參數過多，難于確定等困難，實際應用較少。本文在M2的基礎上，通過引入微平面上體、偏、剪三個分量簡單的本構關系，建立了一個新的微平面本構模型。所建立的模型具有物理概念清楚，參數少且易于確定的特點。

1 M2 模型簡介［3］

在連續介質中任意一點，存在各個不同的方向也就是微平面，如圖1（a）所示，每個微平面上的應力、應變與宏觀意義上的應力、應變張量相聯系。在微平面本構模型中，宏觀張量間的本構關系通過微平面上的應力和應變分量之間的本構關系表達，二者通過張量與矢量間的動態約束（即投影關系）和能量原理建立聯系。圖1（a）中微平面法向和其兩個正交方向的坐標分量可表示為ni、Ii和mi，圖1（b）中的原點到正二十面體各頂點及各邊中點所對應的矢量代表微平面法矢量，共有21個微平面。

圖1 （a） 微平面及其矢量Fig.1（a）Vectors on general microplane

圖1 （b） 每半球表示21個微平面的法線方向Fig.1（b）Directions of microplane normals（circle）for system of 21 microplane per hemisphere

1.1 動態約束

任意一個微平面的應變分量和應變張量滿足動態約束，即投影關系：

其中：εN為微平面法向應變，并可分解為體積應變εV和偏應變εD，這里體積應變指平均體積應變，即εV=（εx+εy+εz）/3。εL和εM是兩個微平面切向應變，Nij、Lij和Mij為各微平面所確定投影張量，并可用微平面矢量表述，即Nij=ninj，Lij=（linj+ljni）/2，Mij=（minj+mjni）/2，其中ni、li和mi為微平面法向矢量和切向矢量，如圖1（a）所示。

1.2 虛功原理

基于虛功原理得到動態約束下的宏觀應力張量σij與微平面應力向量的關系，σij由微平面應力按照如下積分形式得到：

1.3 微平面本構關系

M2模型將微平面上各分量應力應變關系定義為：

其中：CV、CD、CT是割線模量。除壓縮體積分量外均表示為一段由上升段和軟化段表示的曲線，如圖2（a）所示，卸載模量介于初始切線模量和割線模量之間。壓縮體積分量則根據混凝土三軸壓縮力學特點定義為如圖2（b）所示的形式。

2 改進微平面本構模型

2.1 M2模型特點

限于目前的技術水平和試驗設備條件，很難從試驗室里得到混凝土的細觀力學反應數據。在M2模型中，微平面上的應力應變關系曲線是通過分析混凝土受力特征定性給出的，考慮到混凝土在宏觀上具有軟化特征，并根據直覺［4］和試驗驗證將微平面上的本構關系定義為包括上升段和下降段的光滑曲線形式，如圖2所示。實際上，M2模型正是通過這種在微平面上定義軟化段的方式來獲得混凝土宏觀的軟化現象，但這種方法在處理應力應變比值為負的象限的加卸載轉換時遇到了困難，需要引進更復雜的規則［5］，在模擬循環加載時也有很大的困難。

2.2 改進方法

文獻［12］在基于有限單元方法分析巖石的破壞過程時提出，從宏觀上材料可能具有明顯的非線性性質，但在細觀層面上，局部細觀單元體的破裂性質可以假定為彈-脆性行為。文獻［13］在用有限元方法對混凝土單軸壓縮破壞過程進行細觀數值模擬時也采用了這一思想。該方法將混凝土材料分為三類單元，骨料單元、砂漿單元和界面單元，每種細觀單元都采用彈-脆性的應力-應變曲線，通過數值模擬得到了與實驗現象一致的荷載-位移曲線。

通過分析前人的研究成果和對微平面概念的認識，本文提出以下兩個與現有微平面本構模型不同的觀點：

（1）認為微平面上各分量的應力-應變關系（除體積分量外）是理想彈塑性的，應變達到一個界限值微平面發生破斷，應力消失為零；

（2）宏觀的應變軟化現象是由于微平面逐個失效后退出工作引起的。舉例說明，圖3是單軸循環壓縮試驗的軸向應力-應變滯回曲線，隨著應力循環周數的增加，材料的剛度退化，相鄰兩個滯回環中相同應變對應的應力減小，這是由于隨著應變的往復增加，微平面逐漸破斷失效，隨著微平面有效數量的減少，剛度不斷減小，微平面應力對宏觀應力貢獻的總和也隨之減少。

圖3 單軸循環壓縮試驗Fig.3 Uniaxial cyclic compression test

2.3 改進微平面本構模型

根據以上觀點，本文對混凝土微平面上各個分量的本構關系作如下改進：

（1）除壓縮體積分量外，各分量本構模型為理想彈塑性模型；

（2）提出破斷應變的概念，即微應變達到破斷應變后，微平面發生斷裂，應力消失。白衛峰等［14］在對準脆性材料破壞過程的細觀物理模型研究中曾提出斷裂應變的概念；

（3）體積分量的本構關系以靜水壓縮試驗數據為基礎進行構造；

（4）構造加卸載準則，解決原有模型中在應力應變比值為負的二四象限加卸載轉換的困難，并使得模型適用于循環加載。

改進后微平面上各分量的應力-應變曲線如圖4～6所示，下面逐個進行說明，并將各分量應力-應變關系表示為適用于增量顯式算法的形式。

2.3.1 體積分量

微平面模型中體積應變和體積應力的定義分別為εV=εkk/3和σV=σkk/3，即為宏觀的平均正應變和平均正應力，并且所有微平面上的體積分量均相同。因此，在構造體積應力應變關系σV=f（εV）時參考了現有靜水壓縮試驗結果和其它宏觀本構模型對體積應力應變關系的描述。圖4是本文構造的體積應力-應變曲線，其中e點對應的應變表示拉伸體積破斷應變。這一模型在描述土和可壓縮性泡沫材料的體積應力應變關系時也有類似應用［15］。

圖4 體積應力-應變曲線Fig.4 Volumetric strain-stress curve

2.3.2 偏分量

在構造偏分量本構關系時考慮了以下兩點：① 微平面上的偏應變定義為εD=εN-εV，當體積應變εV為零時，上式退化為εD=εN。法向分量在拉壓兩種荷載下的力學性質不同，構造偏分量時參考了法向分量的特征。② 參考了現有模型對偏分量的構造方法。M2模型中將偏分量定義為與法向分量具有相同形式的函數曲線，而且在隨后的M3和M4改進版本中法向分量和偏分量的應力邊界也具有相同的曲線形式。因此，本文也將偏分量和法向構造為相同的形式，應力應變曲線如圖5所示。圖中應力應變關系骨架的不對稱性是為了反映混凝土拉壓力學性質的不同。

將應力-應變關系表示為已知上一步應力、應變和當前步應變增量，求當前步應力的形式：

2.3.3 剪切分量

圖5 偏應力-應變曲線Fig.5 Deviatoric strain-stress curve on microplane

構造剪切分量時，考慮到在微平面上發生剪切變形時與方向無關，也就是說在兩個相反方向發生剪切變形的應力應變過程是相同的，因此剪切分量可以按絕對值進行構造。基于上述分析和本文提出的觀點，構造微平面剪切應力-應變關系如圖6所示。計算過程與偏分量相同。

圖6 剪切應力-應變曲線Fig.6 Shear strain-stress curve on microplane

3 參數確定與分析

模型需要確定的參數主要有剪切分量、偏分量和體積分量三個分量分別對應的彈性模量、屈服應變和破斷應變以及卸載模量等，共14個參數。

基于前述體積分量的物理意義，體積分量中的其它相關參數可以靜水壓縮試驗數據為基礎確定取值范圍。偏分量和剪切分量中的屈服應變和破斷應變參考M3及以后模型中應力應變邊界的應變取值范圍。

圖7是剪切屈服應變和破斷應變取值大小對單軸應力應變關系曲線的影響情況。從圖中可以看出，屈服應變越大峰值應力越高；破斷應變越大則峰值應力越高，且相應的峰值應變越大。這是由于屈服應變增大，單個微平面對宏觀應力的貢獻值越大；在屈服應變保持不變的情況下，隨著破斷應變的增大，峰值應變增大，同時處于最大應力狀態的微平面數目增加，使得峰值應力增加。

圖7 參數取值對應力應變曲線的影響Fig.7 The stress-strain curves impact of parameter values

4 計算實例

為了驗證文中所提觀點的合理性和模型的正確性，下面采用本文模型擬合混凝土基本實驗數據。這些數據均從文獻中得到，包括單軸壓縮試驗和單軸拉伸試驗。圖8～10中的點表示試驗數據，實線表示本文模型的計算結果。

表1～3分別給出了圖8～10中計算模型所采用的參數值。另外，模型計算時，將圖4中拉伸體積應變值取為零，ab階段體積模量取oa段的一半。

圖8 采用本文模型擬合單軸壓縮曲線（Hongnestad，1955）Fig.8 Fit of stress-strain curves measured by Hongnestad（1955）in uniaxial compression

圖9 擬合單軸壓縮曲線（Van Mier，1984）Fig.9 Fit of stress-strain curves measured by Van Mier（1984）in uniaxial compression

圖10 擬合單軸拉伸曲線（Petersson，1981）Fig.10 Fit of stress-strain curves measured by Petersson（1981）in uniaxial tension

圖8是Hongnestad［18］對三種強度混凝土所做的單軸壓縮試驗，從擬合結果來看，本文模型能夠描述不同強度混凝土從彈性階段到非線性階段的強化-軟化過程。圖9是Van Mier［19］給出的單軸壓縮應力-應變全曲線，同時給出了橫向應變隨軸向應力的變化曲線，從圖中可以看出本文模型與試驗數據符合較好。

圖10是Petersson［20］給出的混凝土單軸拉伸試驗結果。從圖中可以看出，本文模型計算結果與試驗數據擬合較好，能夠準確地描述混凝土的彈性模量和抗拉強度等基本特征，并能反應混凝土的拉伸軟化過程。

表1 模型參數取值（圖8）Tab.1 Parameter values in Fig.8

表2 模型參數取值（圖9）Tab.2 Parameter values in Fig.9

表3 模型參數取值（圖10）Tab.3 Parameter values in Fig.10

5 結論

本文在Ba?ant等的混凝土微平面本構模型M2的基礎上，提出了彈塑性微平面本構模型。本文所提出的模型，將法向應力應變分量進行體偏分解，因此包含了體、偏、剪三個分量。本文的研究結論如下：

（1）混凝土宏觀的軟化可以由微平面逐步失效后應力釋放表達。

（2）對于偏、剪切分量應力-應變關系可采用理性彈塑性模型描述，并提出了破斷應變的概念，當應變達到破斷應變后，微平面失效并釋放應力。

［1］ Ba?ant Z P，Oh B H.Microplane model for fracture analysis of concrete structures.Proc.，Symp.On the Interactions of Non-Nuclear Munitions with Structure［C］//U.S.Air Force Academy，Colorado Springs，Colo.，1983：49-53.

［2］ Ba?ant Z P，Oh B H.Microplane model for progressive fracture of concrete and rock［J］.Journal of Engineering Mechanics，1985，114（4）：559-582.

［3］ Ba?ant Z P，Prat P C.Microplane model for brittle plastic material：I.Theory［J］.Journal of Engineering Mechanics，1988，114：1672-1688.

［4］ Ba?ant Z P，Prat P C.Microplane model for brittle plastic material： II. Verification［J］. JournalofEngineering Mechanics，1988，114：1689 –1702.

［5］Ba?ant Z P，Xiang Y，Prat P C.Microplane model for concrete.I.Stress-strain boundariesand finite strain［J］.Journalof Engineering Mechanics，1996，122（3）：245-254.

［6］ Ba?ant Z P，Caner F C，Carol I，et al.Microplane model M4 for concrete：I.Formulation with work-conjugate deviatoric stress［J］. JournalofEngineering Mechanics，2000，126（9）：944-953.

［7］ Ba?ant Z P，Caner F C.Microplane model M5 with kinematic and static constraints for concrete fracture and anelasticity.I：Theory［J］.Journal of Engineering Mechanics，2005，131（1）：31-40.

［8］ Ba?ant Z P，Caner F C.Microplane model M5 with kinematic and static constraints for concrete fracture and anelasticity.II：Computation［J］.Journal of Engineering Mechanics，2005，131（1）：41-47.

［9］Carol I，Prat P C，Ba?ant Z P.New explicit microplane model for concrete： Theoretical aspects and numerical implementation［J］.International Journal of Solids and Structures，1992，29（9）：1173-1191.

［10］ Cofer W F，Kohut S W.A general nonlocal microplane concrete material model for dynamic finite element analysis［J］.Computers ＆ Structures，1994，53（1）：189-199.

［11］ Pivonka P，Ozbolt J，Lackner R，et al.Comparative studies of 3D-constitutive models for concrete：Application to mixedmode fracture［J］. International Journal for Numerical Methods in Engineering，2004，60（2）：549-570.

［12］ Tang C A.Numerical simulation on progressive failure leading to collapse and associated seismicity［J］.International Journal of Rock Mechanics and Mining Science，1997，34（2）：249-261.

［13］ Dang F N，Han W T，Tian W，et al.Three-dimensional numerical simulation on meso-level of uniaxial failure process in concrete［J］.Journal of Xi＇an University of Technology，2006，22（2）：113-118.

［14］白衛峰，陳健云，胡志強，等.準脆性材料單軸拉伸破壞全過程物理模型研究［J］.巖石力學與工程學報，2007，26（4）：670-681.

［15］ LSTC.LS-DYNA USER＇S MANUAL（Version 970）.2003.

［16］ Ba?ant Z P，Caner F C，Adley M D，et al.Fracturing rate effect and creep in microplane for dynamics［J］.Journal of Engineering Mechanics，2000，126（9）：962-970.

［17］Yan D M，Lin G.Dynamic properties of concrete in direct tension［J］.Cement and Concrete Research，2006，36：1371-1378.

［18］ Hognestad E，Hanson N W，McHenry.Concrete stress distribution in ultimate strength design［J］.D.J.Am.Concr.Inst，1955，52（4）：455-477.

［19］ Van Mier J G M.Strain-softening of concrete under multiaxial loading conditions［D］.Thesis presented to University Eindhoven，the Netherlands，in partial fulfillment of the requirements for the degree of Doctor of Philosophy，1984.

［20］ Petersson P E.Crack growth and development for fracture zones in plain concrete and similar materials［R］.Report No.TVBM 1006，Lund Inst.of Tech，.Sweden，1981.

［21］ Dilger W H，Koch R，Kowalczyk R.Ductility of plain and confined concrete under different strain rates［J］.Am.Concrete Inst.Spec.Publ.，Detroit，Mich.R，1984，81（1）：73-81.