考虑不同破坏模式的二维混凝土本构模型

土木工程学报,China Civil Engineering Journal,36(11), 2003,70~74

摘要：本文提出了一种可以考虑不同破坏模式的二维混凝土本构模型。该模型将混凝土按不同应力组合进行分区，按分区设定加卸载力学行为和非线性指标，可以较好地模拟混凝土在非比例加载情况下的力学行为。模型中除常见的固定裂缝模型（FCM）外，还可以选择使用转动裂缝模型（RCM），从而可以分析较复杂的混凝土开裂行为。本文还介绍了利用有限元软件MSC.MARC的二次开发功能将该模型嵌入MARC的方法。并用算例说明了本模型在解决一些复杂问题上的优势。

关键词：混凝土；本构；非线性有限元；二次开发；

A Concrete Constitutive Relationship With Various Damage Models

LU Xinzheng JIANG Jianjing

Department of Civil Engineering, Tsinghua University, Beijing, 100084

Abstract: A concrete constitutive relationship with various damage models is developed in this paper. The failure models of concrete are divided into several zones according to the principle stress ratio. Different loading/unloading behaviors and nonlinear indexes are applied to various zones. So the behavior of concrete under disproportionate loading conditions can be simulated better. Both fixed crack model and rotating crack model can be considered in this model for the complex fracture of concrete. This paper also introduces the method to embed this concrete model into the FEA software of MSC.MARC with the secondary development function of MARC. Numerical examples are shown to illustrate the advantage of this model on some complex problems.

Keywords: nonlinear finite element analysis; concrete; constitutive relationship; secondary development

1 引言

混凝土是一种拉压性能差异显著的材料，且在加载、卸载等不同应力路径下，同一应力水平也会有不同的应变量。针对混凝土这些特点，本文提出了一个综合的二维混凝土本构模型。该模型将混凝土在不同应力组合下的破坏模式进行分区，按分区设定其加卸载力学行为及非线性指标，从而可以较好地模拟混凝土在非比例加载情况下的力学行为。作者根据该本构模型编制了计算程序RCFER2002，并将其嵌入大型有限元程序MSC.MARC中，解决了一些现有有限元程序难以分析的混凝土复杂力学行为。

2 本构模型

2.1强度准则

本模型使用的破坏准则是简化的李-过和Tasuji-Slate-Nilson二轴破坏准则^[1]，其形状如图1所示，。根据主应力和的大小，将破坏曲面分成5个区，如表1所示。这样，将混凝土的二轴破坏准则分为3种破坏形式，5种应力组合分区，使程序物理意义明确，计算结果清晰。

2.2 等效单轴应力应变关系

混凝土的等效受压应力应变关系的骨架线使用Saenz公式^[2]：

文本框: 表1、破坏曲面分区
Table 1. Division of failure surface
分区破坏形态
1 <0 <0 >0.2 双压片状压坏
2 <0 <0 <0.2 单压柱状压坏
3 >0 <0 >-0.05 单压柱状压坏
4 >0 <0 <-0.05 拉坏
5 >0 >0 拉坏 (1)

这里为切线模量，为等效单向应变，为等效单向峰值压应变。为初始弹性模量，为峰值点割线弹性模量。

卸载刚度，，为此时的等效单向压应力，等效单向压应变。为系数，我们通过数值试验建议取2.0～3.5。

受拉应力应变关系采用江见鲸^[2]建议公式

（2）

为混凝土抗压强度，为混凝土开裂应变，为系数，与混凝土的断裂能相关。江见鲸^[3]建议，，这里为混凝土断裂能，为试件特征尺寸。在不知道的情况下，根据我们数值试验的结果，可以取3000～50000之间，一般建议用试算。

2. 3 裂缝模型和裂面剪力传递系数

在本文的混凝土模型中，有两种弥散裂缝模型可供选择，即固定裂缝模型(Fixed Crack Model (FCM))和转动裂缝模型(Rotating Crack Model (RCM))^[4]。在FCM中，假设裂缝一旦出现，其方向不会随主应力变化而变化。由于FCM可以较好的克服网格倾斜(Mesh bias)带来的影响^[5]，因此我们推荐在分析受拉开裂为主的问题中，FCM比较合适。但是，由于FCM可能会过高估计裂面剪力传递能力^[6]，因此，当主应力方向变化较大，裂面剪力很大的时候会带来较大误差。因此，对于这些问题，我们推荐使用RCM模型。在RCM模型中，裂缝方向始终和主拉应变方向垂直，可以降低裂面的剪力锁死问题(Stress Locking)^[6]，比较适合于主应力方向变化较大的问题。

当裂缝方向的时，我们认为此时裂缝受压闭合。

开裂混凝土的折算剪切模量 ^[7]，在本文的混凝土模型中，裂面剪力传递系数取值为：

(3)

这里，分别为混凝土刚刚开裂和裂缝充分发展时的剪力传递系数。为常数，反映降低的速度。该关系具有比较好的通用性，如取，则为一定值，即为常用的固定剪力传递系数模型。根据我们数值试验的结果：这种模型，对于混凝土剪力墙，推荐取固定剪力传递系数为0.125；对于混凝土深梁，取0.25。一般梁板结构取为0.5。如果取，，，则与AL-Manaidi^[7]建议的裂面剪力传递系数变化关系相接近，裂面剪力传递系数将随裂缝开展宽度增加而降低。

2. 4 加卸载关系

当材料处于破坏曲面以内时，程序使用应力空间来定义加卸载关系。在破坏曲面以外，则用应变空间来定义加载和卸载。在本混凝土模型中，共使用了6个非线性指标 ~ ，其各自含义见表2。定义表示各非线性指标曾经达到的最大值，则如则处于卸载或卸载后再加载阶段，，则处于沿骨架线加载阶段。

表2 非线性指标

非线性指标
含义	双压破坏非线性指标	单压破坏非线性指标	未出现破坏时的非线性指标	有一条裂缝时，对应该裂缝的非线性指标	有两条裂缝时，对应第一条裂缝的非线性指标	有两条裂缝时，对应第二条裂缝的非线性指标

各非线性指标的具体定义为：

（4）

这里为消除泊松比影响后的主压应变。为混凝土峰值压应变。

（5）

这里为混凝土主压应变。

根据应力组合分区不同而有所不同，如图2所示。

(6)

(7)

(8)

这里为消除泊松比影响后的主拉应变，为混凝土开裂应变。

在应力空间计算时，加卸载的判别是对屈服面求偏导，如，则为加载，，则为卸载。而在有限元位移法中每次迭代只能得到，，为本构矩阵。矩阵的取值又依赖于加卸载的判断，因此矩阵和相互依赖，需要迭代判别。在文献[7]中建议反复迭代直至预测加卸载状态和实际加卸载状态完全一致。但是我们通过数值试验发现这种方法很可能因为加卸载状态的反复摆动而使得计算效率很低。同时我们发现经过一次迭代修正后其加卸载状况和实际情况一般已经很接近了。因此在程序中，对加卸载状态的迭代修正都取为1次。

3 RCFER2002和MARC的衔接

由MSC公司开发的有限元软件MARC具有强大的非线性分析能力，完整的单元库和便捷的前后处理系统，在科研和生产领域内都有很多成功应用。而其提供的用户子程序^[9]（user subroutine）功能允许用户根据自己需要编制相应的子程序模块，编译后嵌入MARC中，实现某些特殊功能，更是给用户的科研和二次开发带来了很大的方便。通过将该混凝土本构模型嵌入MARC，取代MARC中的混凝土材料模型。同时利用MARC强大的前后处理能力和出色的非线性求解模块，可以用于解决一些复杂的混凝土结构问题。

3.1 本构子程序

在MARC的user subroutine中，提供了多种添加材料本构模型的子程序。经过比较，我们认为，其中的超弹性（hyperelastic）本构子程序HYPLEA稳定性和接口都较好，比较适合用来和RCFER2002相衔接，该子程序的接口为：

SUBROUTINE HYPELA (D, G, E, DE, S, TEMP, DTEMP, NGENS, N, NN, KC, MATS, NDI, NSHEAR)

D(NGENS, NGENS)为材料的切线刚度阵。G(NGENS)为温度引起的应力增量。E(NGENS)为此时的应变。DE(NGENS)为应变增量。S(NGENS)为此时的应力。TEMP为此时的温度。DTEMP为温度增量。NGENS为应力分量的数量，对于平面应力问题NGENS=3，对于空间问题NGENS=6。N为单元号。NN为积分点号。KC为材料层号（对于多层材料）。MATS为材料编号。NDI为正应力分量数量。NSHEAR为剪应力分量数量。

由于该子程序是针对超弹性材料的，即材料的状态和应力历史无关，这和混凝土的实际情况有较大差距。为了解决这个问题，我们利用Fortran 90语言提供的模块module功能，在内存中申请一个模块，用来存放混凝土材料的应力历史，从而可以有效解决这个问题。

具体做法是：

在用户程序中添加以下语句：

module Varies

type(typ_Concrete) :: Concrete (NElem, NIntPoint)

end module

其中typ_Concrete为混凝土材料和应力历史参数的数据结构。NElem为混凝土单元数，NIntPoint为每个单元的高斯积分点数。这样就实现了在内存中开辟数据存储空间的功能。

只要在MARC的本构子程序HYPLEA开头添加以下语句：

use Varies

就可以实现对Varies模块的访问。

3.2 后处理子程序

利用MARC user subroutine中的POSTV子程序，可以很方便的进行绘制开裂应变、应力状态等后处理操作。由于混凝土的各项应力历史参数都已经存放在Varies模块中，因此只要在POSTV子程序中引用该模块，并把相应的变量赋予相应的后处理参数，就可以十分方便的使用MARC强大的后处理功能。

4 算例

图3、图4为RCFER2002嵌入MARC后得到的在单向往复受压和往复拉压作用下的应力应变关系。与文献[3]中给出的试验结果吻合良好。可见本模型可以较好反映混凝土在往复荷载作用下的力学行为。

纤维增强材料（FRP）布与混凝土之间的粘结锚固性能一直是混凝土有限元分析的一个难点。如果使用MARC自带的混凝土模型，由于其使用固定裂缝（FCM）模型，且剪力传递系数为一常数，很难模拟出纤维布剥离的效果。其他一些有限元软件如ANSYS等也难以实现对剥离的模拟^[10]。而在RCFER2002中，由于使用了转动裂缝模型（RCM），加上程序中比较复杂的裂面行为模型，使得计算结果和试验结果^[9]（试件PG1-11和PG1-12）吻合良好。并且利用MARC强大的非线性分析功能，得到了整个纤维布剥离的全过程。有限元模型如图5所示，计算结果对比如图6所示，最终了裂缝和变形分布如图7所示。

5 结论

本文提出了一种可以考虑不同破坏模式的混凝土二维本构模型，并将该模型嵌入MSC.MARC，可以分析某些特殊应力条件下的混凝土构件，尤其是那些主应力方向变化显著或加卸载行为显著的混凝土问题，计算结果要优于一些现有有限元软件中混凝土本构模型的计算结果。

参考文献

1 过镇海，钢筋混凝土原理[M]，北京：清华大学出版社，1999

2 江见鲸，钢筋混凝土结构非线性有限元分析[M]，西安：陕西科学技术出版社，12-13，1984

3 沈聚敏，王传志，江见鲸，钢筋混凝土有限元与板壳极限分析[M]，北京：清华大学出版社，234，1993

4 Z. P. Bazant, Concrete fracture models: testing and practice[J], Engineering Fracture Mechanics, Vol.69, 165-205, 2002

5 Y. J. Li, Th. Zimmerman, Numerical evaluation of rotating crack model[J], Computers and Structures, Vol.69, 487-497,1998

6 Classification of smeared crack models based on material and structural properties[J], Int. J. Solid Structures. Vol. 35 No.12, 1289-1308, 1998

7 吕西林，金国芳，吴晓涵，钢筋混凝土结构非线性有限元理论与应用[M]，上海：同济大学出版社，173，1997

8 MSC.Marc User’s Manual, Volume D (User Subroutines and Special Routines), MSC. Software Corporation, 2001

9 谭壮，GFRP布加固混凝土梁受力性能的试验研究[D]，清华大学硕士学位论文，2002

10 陆新征，严吉洁，韦韩等，FRP布-混凝土界面有限元分析技术难点讨论，全国土木工程研究生学术论坛论文集（已接受）