基于Zhang‒Hou损伤计算模型的不同强度混凝土裂缝贯通与损伤演化规律

Crack Coalescence and Damage Evolution Laws of Concretes with Different Strengths Based on the Zhang‒Hou Damage Calculation Model

摘要

混凝土裂缝演化规律的研究一直具有挑战性。为分析不同强度混凝土裂缝贯通与损伤演化的内在规律，基于混凝土塑性损伤（concrete damaged plasticity，CDP）模型的理论框架，结合Sidoroff损伤理论，建立Zhang‒Hou损伤计算模型。该模型能够定量描述混凝土的单轴压缩和单轴拉伸损伤行为，有利于深入研究混凝土裂缝扩展及其损伤演化过程。研究发现：不同强度混凝土的单轴压缩和单轴拉伸损伤演化规律基本一致，并且在损伤演化过程中发现了单轴压缩和单轴拉伸损伤交叉点，对应的损伤分别为0.57和0.62。损伤交叉点是损伤演化曲线从快速发展转向缓慢发展的一个临界点。不同强度混凝土在损伤交叉点前后的损伤演化规律一致。在交叉点前，强度越高，对应的损伤越小；交叉点之后则呈现相反的情况。基于以上发现，根据突变性判据的内涵，将发现的损伤交叉点作为混凝土出现裂缝区域性贯通的临界判据，对应的损伤被解释为量化混凝土裂缝区域性贯通的损伤阈值。进一步地，利用已验证的钢筋混凝土（reinforced concrete，RC）有限元模型，对不同加载位移下RC构件的损伤状态和裂缝贯通区域进行系统分析，验证损伤阈值的合理性。所得结论可为基于连续介质力学的混凝土损伤分析和裂缝模拟提供理论依据。

Abstract

The study of concrete crack propagation is always challenging. A Zhang‒Hou damage calculation model is developed based on the theoretical framework of the concrete damaged plasticity (CDP) model and incorporates Sidoroff's damage theory to analyze the inherent patterns between crack coalescence and damage evolution in concrete of different strengths. This model quantitatively describes the uniaxial compressive and tensile damage behavior of concrete, enabling an in-depth investigation of crack propagation and damage evolution in concrete. The study finds that the damage evolution patterns under uniaxial compression and tension are generally consistent across concretes of different strengths. A critical damage crossover point is identified, where the damage values for compression and tension are 0.57 and 0.62, respectively. This crossover point marks a transition from rapid to slow damage evolution. The damage evolution patterns before and after this point remain consistent across different strengths: before the crossover, higher-strength concrete exhibits lower damage, while after the crossover, the trend reverses. Based on these findings and the implications of the catastrophe criterion, the identified damage crossover point is proposed as a critical criterion for the crack coalescence of concrete, and the corresponding damage values are interpreted as a damage threshold that characterizes the crack coalescence in concrete. In addition, a validated finite element model of reinforced concrete (RC) is utilized to comprehensively analyze its damage state and crack coalescence regions under various loading displacements, confirming the validity of the proposed damage threshold. The findings provide a theoretical basis for concrete damage analysis and crack simulation within the framework of continuum mechanics.

Graphical abstract

关键词

Zhang-Hou损伤计算模型 / 压缩损伤 / 拉伸损伤 / 裂缝贯通 / 损伤阈值

本刊网刊

随着中国二次能源结构逐渐优化，混凝土在中国水利工程建设中得到广泛应用^[1]。然而，混凝土的复杂材料属性使其在形成初期就存在微裂缝和微孔洞等缺陷。单轴荷载下的混凝土应力‒应变曲线表现出上升段和下降段，这是由混凝土细观裂缝发展和宏观裂缝形成引起的，也是材料连续损伤和损伤逐步累积的宏观表现。混凝土的损伤发展是微观裂缝萌生、扩展、贯通的连续损伤演化过程，最终破坏是裂缝出现贯通直至最后失效的阶段^[2]。鉴于混凝土损伤过程的复杂性，有必要将混凝土的裂缝扩展与损伤演化过程紧密联系起来，以进一步研究混凝土的破坏过程。

目前已经有许多关于混凝土破坏机理的研究，其中包括：通过原位XCT扫描提出具有真实内部特征的混凝土结构2维细观模型，并分析裂缝发展对混凝土断裂损伤行为的影响^[3]；采用离散元法建立包含骨料、砂浆和界面过渡区的非均质模型，模拟混凝土微观裂纹及损伤演化过程^[4]；建立基于自适应图像的多尺度建模和计算方法，模拟混凝土结构从细观到宏观损伤直至结构破坏的跨尺度动态破坏过程^[5‒6]等。同时，由于混凝土破坏过程具有跨尺度、非线性和不可逆等基本特点，不同的非线性理论（如分形几何、损伤和断裂力学等）被用来探索混凝土从细观到宏观破坏的演化过程。例如：殷亚娟等^[2]提出一种基于裂纹分形维的损伤变量计算方法，并利用分形维对混凝土开裂过程中的裂纹分布情况进行表征；Yin^[7]和Zhang^[8]等通过建立细观骨料模型研究了混凝土在单轴荷载下的宏观应力‒应变曲线及裂纹损伤状态；商效瑀^[9]和韩燕华^[10]等采用CT图像扫描技术研究了混凝土内部细观破裂特征，为混凝土内部微裂缝的研究提供了新思路。此外，模拟混凝土裂缝的材料模型也得到相继发展。如：刘智等^[1]建立了考虑滞回效应的混凝土损伤模型，并分析了重力坝的震害损伤分布；焦延涛等^[11]建立了新的弹塑性损伤本构模型，模拟了双边开口梁的裂纹扩展路径；虢成功等^[12]通过混凝土纳‒微‒细观随机损伤模型预测了混凝土的单轴应力‒应变曲线。

上述研究在试验技术、裂缝建模和材料模型等方面为混凝土破坏研究提供了有力支持。然而，尚缺乏定量分析混凝土的裂缝扩展及其损伤演化等方面的研究，对于混凝土出现宏观断裂时的定量描述存在差异。谢和平^[13]指出应该以材料出现宏观断裂时的损伤阈值作为破坏判据，简单地将损伤变量D=\mathrm{1}作为破坏判据并不合理。殷亚娟等^[2]和吕从聪^[14]将损伤阈值设定为0.8，王高辉^[15]和李祚华^[16]假定其为0.75。但这些研究尚未区分混凝土在压缩和拉伸损伤‒破坏方面的不同。此外，目前研究未分析不同强度混凝土裂缝扩展及其损伤演化的内在规律。

本文基于新建立的Zhang‒Hou损伤计算模型，定量分析了不同强度混凝土裂缝扩展与损伤演化过程。同时，通过有限元模拟和试验结果验证了研究结果。

1 研究方法

1.1 CDP模型损伤表述

混凝土塑性损伤（concrete damage plasticity，CDP）模型^[17]根据塑性增量理论将应变张量\mathit{\epsilon }分解为弹性应变张量{\mathit{\epsilon }}^{\mathrm{e}}和塑性应变张量{\mathit{\epsilon }}^{\mathrm{p}}，并采用有效应力进行损伤建模，同时使用标量退化损伤参数D体现各向同性损伤，相应的损伤本构关系表述如下：

\overline{\mathit{\sigma }}=\mathit{E}:{\mathit{\epsilon }}^{\mathrm{e}}=(\mathrm{1}-D){\mathit{E}}_{\mathrm{0}}:(\mathit{\epsilon }-{\mathit{\epsilon }}^{\mathrm{p}})（1）

\mathit{E}=(\mathrm{1}-D){\mathit{E}}_{\mathrm{0}}（2）

式（1）~（2）中，\overline{\mathit{\sigma }}为有效应力，{\mathit{E}}_{\mathrm{0}}为初始弹性刚度张量，\mathit{{\rm E}}为4阶弹性刚度张量。

若将式（1）退化为1维问题，可得单轴混凝土损伤本构关系，具体如下：

{\sigma }_{\mathrm{c}}=(\mathrm{1}-d_{\mathrm{c}})E_{\mathrm{0}}({\epsilon }_{\mathrm{c}}-{\epsilon }_{\mathrm{c}}^{\mathrm{p}\mathrm{l}})（3）

{\sigma }_{\mathrm{t}}=(\mathrm{1}-d_{\mathrm{t}})E_{\mathrm{0}}({\epsilon }_{\mathrm{t}}-{\epsilon }_{\mathrm{t}}^{\mathrm{p}\mathrm{l}})（4）

式（3）~（4）中：{\sigma }_{\mathrm{c}}和{\sigma }_{\mathrm{t}}为单轴压缩和拉伸应力；{\epsilon }_{\mathrm{c}}和{\epsilon }_{\mathrm{t}}为单轴压缩和拉伸应变；E_{\mathrm{0}}为材料初始弹性刚度；d_{\mathrm{c}}和d_{\mathrm{t}}分别为单轴压缩和拉伸损伤演化参数；{\epsilon }_{\mathrm{c}}^{\mathrm{p}\mathrm{l}}和{\epsilon }_{\mathrm{t}}^{\mathrm{p}\mathrm{l}}分别为压缩和拉伸等效塑性应变，由式（5）和（6）计算。

{\epsilon }_{\mathrm{c}}^{\mathrm{p}\mathrm{l}}={\epsilon }_{\mathrm{c}}^{\mathrm{i}\mathrm{n}}-\frac{d_{\mathrm{c}}}{\mathrm{1}-d_{\mathrm{c}}}\cdot \frac{{\sigma }_{\mathrm{c}}}{E_{\mathrm{0}}},{\epsilon }_{\mathrm{c}}^{\mathrm{i}\mathrm{n}}={\epsilon }_{\mathrm{c}}-{\sigma }_{\mathrm{c}}/E_{\mathrm{0}}（5）

{\epsilon }_{\mathrm{t}}^{\mathrm{p}\mathrm{l}}={\epsilon }_{\mathrm{t}}^{\mathrm{c}\mathrm{k}}-\frac{d_{\mathrm{t}}}{\mathrm{1}-d_{\mathrm{t}}}\cdot \frac{{\sigma }_{\mathrm{t}}}{E_{\mathrm{0}}},{\epsilon }_{\mathrm{t}}^{\mathrm{c}\mathrm{k}}={\epsilon }_{\mathrm{t}}-{\sigma }_{\mathrm{t}}/E_{\mathrm{0}}（6）

式（5）~（6）中，{\epsilon }_{\mathrm{c}}^{\mathrm{i}\mathrm{n}}为非弹性应变，{\epsilon }_{\mathrm{t}}^{\mathrm{c}\mathrm{k}}为开裂应变。

1.2 基于Sidoroff损伤理论的弹塑性损伤模型

在连续损伤力学中，损伤及其增长对材料变形的影响可以用有效应力表示。Sidoroff损伤理论^[18‒19]的能量等效假设原理指出，无损材料的弹性余能\rho {\psi }_{\mathrm{e},\mathrm{0}}和等效损伤材料的弹性余能\rho {\psi }_{\mathrm{e},\mathrm{d}}相等，即：

\rho {\psi }_{\mathrm{e},\mathrm{0}}=\frac{{\sigma }^{\mathrm{2}}}{\mathrm{2}{E}_{\mathrm{0}}}=\rho {\psi }_{\mathrm{e},\mathrm{d}}=\frac{{\overline{\sigma }}^{\mathrm{2}}}{\mathrm{2}{E}_{\mathrm{d}}}=\frac{{\sigma }^{\mathrm{2}}}{\mathrm{2}{E}_{\mathrm{0}}{(\mathrm{1}-D)}^{\mathrm{2}}}（7）

式中，\sigma为应力，E_{\mathrm{d}}为等效损伤材料的弹性模量。

由式（7）可得损伤材料的应力‒应变关系：

\sigma =E_{\mathrm{d}}\epsilon =E_{\mathrm{0}}{(\mathrm{1}-D)}^{\mathrm{2}}\epsilon（8）

式中，\epsilon为弹性应变。

因此，标量退化损伤参数D可表示为：

D=\left\{\begin{array}{l}{d}_{\mathrm{c}}=\mathrm{1}-\sqrt[]{\frac{{\sigma }_{\mathrm{c}}}{E_{\mathrm{0}}{\epsilon }_{\mathrm{c}}}},\\ d_{\mathrm{t}}=\mathrm{1}-\sqrt[]{\frac{{\sigma }_{\mathrm{t}}}{E_{\mathrm{0}}{\epsilon }_{\mathrm{t}}}}\end{array}\right.（9）

显然，式（9）与总应变（{\epsilon }_{\mathrm{c}}和{\epsilon }_{\mathrm{t}}）相关，没有考虑塑性应变对损伤演化的影响，属于弹性损伤计算模型。为反映塑性应变对损伤演化的影响，根据CDP模型的应变分解方法，分别采用{\epsilon }_{\mathrm{c}}-{\epsilon }_{\mathrm{c}}^{\mathrm{p}\mathrm{l}}和{\epsilon }_{\mathrm{t}}-{\epsilon }_{\mathrm{t}}^{\mathrm{p}\mathrm{l}}代替方程（9）中的{\epsilon }_{\mathrm{c}}和{\epsilon }_{\mathrm{t}}，则可得到新的弹塑性损伤计算模型：

D=\left\{\begin{array}{l}{d}_{\mathrm{c}}=\mathrm{1}-\sqrt[]{\frac{{\sigma }_{\mathrm{c}}}{E_{\mathrm{0}}({\epsilon }_{\mathrm{c}}-{\epsilon }_{\mathrm{c}}^{\mathrm{p}\mathrm{l}})}},\\ d_{\mathrm{t}}=\mathrm{1}-\sqrt[]{\frac{{\sigma }_{\mathrm{t}}}{E_{\mathrm{0}}({\epsilon }_{\mathrm{t}}-{\epsilon }_{\mathrm{t}}^{\mathrm{p}\mathrm{l}})}}\end{array}\right.（10）

由式（5）和（6）可知，变量{\epsilon }_{\mathrm{c}}^{\mathrm{p}\mathrm{l}}和{\epsilon }_{\mathrm{t}}^{\mathrm{p}\mathrm{l}}与损伤变量d_{\mathrm{t}}和d_{\mathrm{c}}相耦合。为求解变量{\epsilon }_{\mathrm{t}}^{\mathrm{p}\mathrm{l}}和{\epsilon }_{\mathrm{c}}^{\mathrm{p}\mathrm{l}}，假定等效塑性应变与非弹性应变/开裂应变之间存在线性关系：

{\epsilon }_{\mathrm{c}}^{\mathrm{p}\mathrm{l}}={\beta }_{\mathrm{c}}{\epsilon }_{\mathrm{c}}^{\mathrm{i}\mathrm{n}},{\epsilon }_{\mathrm{t}}^{\mathrm{p}\mathrm{l}}={\beta }_{\mathrm{t}}{\epsilon }_{\mathrm{t}}^{\mathrm{c}\mathrm{k}}（11）

式中，{\beta }_{\mathrm{c}}和{\beta }_{\mathrm{t}}分别为等效压缩应变因子和等效拉伸应变因子，，。

为简便，将上述新的弹塑性损伤模型命名为Zhang‒Hou损伤计算模型，即：

\left\{\begin{array}{l}{d}_{\mathrm{c}}=\mathrm{1}-\sqrt[]{\frac{{\sigma }_{\mathrm{c}}}{E_{\mathrm{0}}({\epsilon }_{\mathrm{c}}-{\beta }_{\mathrm{c}}{\epsilon }_{\mathrm{c}}^{\mathrm{i}\mathrm{n}})}},\\ d_{\mathrm{t}}=\mathrm{1}-\sqrt[]{\frac{{\sigma }_{\mathrm{t}}}{E_{\mathrm{0}}({\epsilon }_{\mathrm{t}}-{\beta }_{\mathrm{t}}{\epsilon }_{\mathrm{t}}^{\mathrm{c}\mathrm{k}})}}\end{array}\right.（12）

由式（12）可知，d_{\mathrm{t}}和d_{\mathrm{c}}的计算需要利用混凝土的单轴本构关系。其中，单轴受压本构关系采用3段式描述，分别为线弹性段（式（13））、强化段（式（14））^[20]和软化段（式（15））^[21]，具体可表示为：

{\sigma }_{\mathrm{c}\mathrm{1}}=E_{\mathrm{0}}{\epsilon }_{\mathrm{c}}（13）

{\sigma }_{\mathrm{c}\mathrm{2}}=\frac{E_{\mathrm{c}\mathrm{i}}\frac{{\epsilon }_{\mathrm{c}}}{f_{\mathrm{c}\mathrm{m}}}-{\left(\frac{{\epsilon }_{\mathrm{c}}}{{\epsilon }_{\mathrm{c}\mathrm{m}}}\right)}^{\mathrm{2}}}{\mathrm{1}+\left(E_{\mathrm{c}\mathrm{i}}\frac{{\epsilon }_{\mathrm{c}\mathrm{m}}}{f_{\mathrm{c}\mathrm{m}}}-\mathrm{2}\right)\frac{{\epsilon }_{\mathrm{c}}}{{\epsilon }_{\mathrm{c}\mathrm{m}}}}{f}_{\mathrm{c}\mathrm{m}}（14）

{\sigma }_{\mathrm{c}\mathrm{3}}={\left(\frac{\mathrm{2}+{\gamma }_{\mathrm{c}}{f}_{\mathrm{c}\mathrm{m}}{\epsilon }_{\mathrm{c}\mathrm{m}}}{\mathrm{2}{f}_{\mathrm{c}\mathrm{m}}}-{\gamma }_{\mathrm{c}}{\epsilon }_{\mathrm{c}}+\frac{{\epsilon }_{\mathrm{c}}^{\mathrm{2}}{\gamma }_{\mathrm{c}}}{\mathrm{2}{\epsilon }_{\mathrm{c}\mathrm{m}}}\right)}^{-\mathrm{1}}（15）

其中，

{\gamma }_{\mathrm{c}}=\frac{{\mathrm{\pi }}^{\mathrm{2}}{f}_{\mathrm{c}\mathrm{m}}{\epsilon }_{\mathrm{c}\mathrm{m}}}{\mathrm{2}{\left[\frac{G_{\mathrm{c}}}{l_{\mathrm{e}\mathrm{q}}}-\mathrm{0.5}{f}_{\mathrm{c}\mathrm{m}}\left({\epsilon }_{\mathrm{c}\mathrm{m}}\left(\mathrm{1}-b\right)+\mathrm{0.9}\frac{f_{\mathrm{c}\mathrm{m}}}{E_{\mathrm{0}}}\right)\right]}^{\mathrm{2}}}（16）

E_{\mathrm{0}}=(\mathrm{0.8}+\mathrm{0.2}\frac{f_{\mathrm{c}\mathrm{m}}}{\mathrm{88}})E_{\mathrm{c}\mathrm{i}}（17）

E_{\mathrm{c}\mathrm{i}}=\mathrm{21}\mathrm{500}{\left(\frac{f_{\mathrm{c}\mathrm{m}}}{\mathrm{10}}\right)}^{\mathrm{1}/\mathrm{3}}（18）

式（13）~（18）中：{\sigma }_{\mathrm{c}\mathrm{1}}、{\sigma }_{\mathrm{c}\mathrm{2}}和{\sigma }_{\mathrm{c}\mathrm{3}}分别为线弹性段、强化段和软化段的应力；E_{\mathrm{c}\mathrm{i}}为初始切割模量；{\epsilon }_{\mathrm{c}\mathrm{m}}为峰值应力对应的应变；f_{\mathrm{c}\mathrm{m}}为圆柱体抗压平均值；{\gamma }_{\mathrm{c}}为软化参数；l_{\mathrm{e}\mathrm{q}}为单元特征长度；b为等效塑性应变与非弹性应变之比，初始取值为0.9；G_{\mathrm{c}}为压溃断裂能^[22]。

另有：

G_{\mathrm{F}}=\mathrm{0.073}(f_{\mathrm{c}\mathrm{m}}{)}^{\mathrm{0.18}},G_{\mathrm{c}}={\left(\frac{f_{\mathrm{c}\mathrm{m}}}{f_{\mathrm{t}\mathrm{m}}}\right)}^{\mathrm{2}}{G}_{\mathrm{F}}（19）

f_{\mathrm{c}\mathrm{m}}=f_{\mathrm{c}\mathrm{k}}+\mathrm{8},f_{\mathrm{t}\mathrm{m}}=\mathrm{0.301}\mathrm{6}{f}_{\mathrm{c}\mathrm{k}}^{\mathrm{2}/\mathrm{3}}（20）

式（19）~（20）中，G_{\mathrm{F}}为开裂断裂能^[20]，f_{\mathrm{t}\mathrm{m}}为圆柱体抗拉强度平均值，f_{\mathrm{c}\mathrm{k}}为圆柱体抗压强度标准值。

单轴受拉本构关系采用李晓琴等^[22]提出的应力‒开裂位移关系，具体如下：

{\sigma }_{\mathrm{t}}(w)=f_{\mathrm{t}\mathrm{m}}{\mathrm{e}}^{-cw}（21）

c=\left\{\begin{array}{l}\mathrm{2.135}\mathrm{7}{f}_{\mathrm{c}\mathrm{u}}^{\mathrm{0.527}\mathrm{2}},f_{\mathrm{c}\mathrm{u}}\le \mathrm{50}\mathrm{M}\mathrm{P}\mathrm{a};\\ \mathrm{2.501}\mathrm{5}{f}_{\mathrm{c}\mathrm{u}}^{\mathrm{0.506}\mathrm{8}},f_{\mathrm{c}\mathrm{u}}>\mathrm{50}\mathrm{M}\mathrm{P}\mathrm{a}\end{array}\right.（22）

式（21）~（22）中，{\sigma }_{\mathrm{t}}(w)为拉伸应力，w和w_{\mathrm{0}}分别为开裂位移和最大开裂位移，w_{\mathrm{0}}=\mathrm{5}{G}_{\mathrm{F}}/f_{\mathrm{t}\mathrm{m}}，c为软化参数，f_{\mathrm{c}\mathrm{u}}为立方体抗压强度。

根据裂缝带理论^[23]，可利用式（23）将应力‒开裂位移曲线转换为应力‒应变曲线：

{\epsilon }_{\mathrm{t}}=f_{\mathrm{t}\mathrm{m}}/E_{\mathrm{0}}+w/l_{\mathrm{e}\mathrm{q}}（23）

图1为弹性损伤计算模型和Zhang‒Hou损伤计算模型的损伤演化曲线对比，使用到的材料参数为：f_{\mathrm{c}\mathrm{k}}=\mathrm{24.5}\mathrm{M}\mathrm{P}\mathrm{a}，f_{\mathrm{t}\mathrm{m}}=\mathrm{1.95}\mathrm{M}\mathrm{P}\mathrm{a}，E_{\mathrm{0}}=\mathrm{24}\mathrm{789}\mathrm{M}\mathrm{P}\mathrm{a}，E_{\mathrm{c}\mathrm{i}}=\mathrm{28}\mathrm{972}\mathrm{M}\mathrm{P}\mathrm{a}，G_{\mathrm{F}}=\mathrm{0.13}\mathrm{N}/\mathrm{m}\mathrm{m}，G_{\mathrm{c}}=\mathrm{20.39}\mathrm{N}/\mathrm{m}\mathrm{m}，l_{\mathrm{e}\mathrm{q}}=\mathrm{150}\mathrm{m}\mathrm{m}，w_{\mathrm{0}}=\mathrm{0.342}\mathrm{m}\mathrm{m}，c=\mathrm{13.90}，{\beta }_{\mathrm{t}}=\mathrm{0.5}，{\beta }_{\mathrm{c}}=\mathrm{0.5}。由图1可知，两种模型的单轴压缩/拉伸损伤演化曲线变化规律基本一致，开始的线弹性阶段损伤为0，随后损伤进入快速发展阶段，最后损伤增长变慢并趋于1。此外，考虑塑性应变影响的Zhang‒Hou损伤计算模型的损伤演化曲线发展速度比弹性损伤模型更缓慢。

采用单轴循环试验结果^[24]验证Zhang‒Hou损伤计算模型，并与弹性损伤计算模型对比，结果见图2。由图2可知，Zhang‒Hou损伤计算模型能准确反映试验的加载—卸载—再加载特征，而弹性损伤计算模型高估了试验卸载—再加载斜率。结果证实了Zhang‒Hou损伤计算模型的合理性，其具有较强的普适性。

1.3 突变性判据

系统处于极限平衡状态意味着它是一种平衡状态向另一种平衡状态或非平衡状态的转变，即系统的状态发生了突变^[2]。突变性判据认为任何反映系统状态突变的现象都可以被视为失稳判据，如：位移突变判据、能量突变判据等。本文将混凝土裂缝出现区域性贯通理解为宏观断裂的起点并利用不同强度混凝土在单轴荷载下裂缝扩展过程中损伤演化的共性规律来探究混凝土裂缝出现区域性贯通时的临界判据。

2 混凝土裂缝贯通与损伤演化

2.1 混凝土的单轴力学性质及其裂缝状态

混凝土单轴应力‒应变曲线描述了混凝土在加载过程中的力学性质，其发展过程对应混凝土裂缝的萌生、扩展、贯通直至最后失效的过程。如图3所示，单轴压缩应力‒应变曲线有4个阶段：初始线性段、非线性段、峰值后的下降段和平稳段；单轴拉伸应力‒应变曲线主要包括初始线性段、峰值后下降段和平稳段。

混凝土的破坏是由材料内部裂缝扩展引起的，单一的单轴试验难以研究材料内部破坏机理。细观破坏模拟是研究材料内部破坏机理的重要手段^[26]，殷亚娟^[2]、Sun^[5]等建立的非均质细观模型为混凝土在单轴拉伸和压缩荷载下的全过程裂缝扩展研究提供了详实结果，这些结果（如应力‒应变曲线和裂缝破坏形态）与其他数值模拟结果^[27-29]和试验结果^[25,30]相符，这种细观数值模拟方法的准确性在工程结构模拟^[31]中也得到了验证。本文总结了混凝土在单轴压缩和单轴拉伸荷载下的裂缝扩展结果，如图4~5所示。

这些结果显示了混凝土在单轴应力下的裂缝扩展特征：在B点处，裂缝趋于局部化发展；随后，裂缝逐渐扩展，大约在C点处出现裂缝区域性贯通，根据文献^[2,5,25,32]的研究结果，该点对应的单轴压缩应力为峰值应力的30%~60%，单轴拉伸应力为峰值应力的20%~51%。之后，贯通裂缝继续扩展并伴随新的宏观裂缝的形成，直至最后结构失效（D点）。

2.2 不同强度混凝土的单轴应力‒应变曲线及其损伤演化结果

为研究不同强度混凝土裂缝发展与损伤演化的内在规律，对强度等级在C30至C80范围内的混凝土进行分析。不同强度混凝土的材料参数见表1，相应的分析结果如图6所示。由图6（a）、（b）可以看出，混凝土强度等级越高，下降段曲线坡度越陡峭，并且最大压缩和拉伸应力对应的应变有轻微提高。由图6（c）、（d）可以看出：1）不同强度混凝土的单轴压缩和单轴拉伸损伤演化规律基本一致，并且在损伤演化过程中发现了单轴压缩和单轴拉伸损伤交叉点，对应的损伤分别为0.57和0.62。2）不同强度混凝土在损伤交叉点前后的损伤演化规律一致。在交叉点前，强度越高，对应的损伤越小；交叉点之后则相反。3）观察损伤交叉点前后损伤演化曲线的变化形态，该交叉点是不同强度混凝土损伤演化曲线从快速发展转向缓慢发展的一个临界点。因此，根据突变性判据的内涵，该损伤交叉点可作为量化混凝土裂缝网格状态的一个有力指标。

表2为不同强度混凝土在损伤交叉点处的单轴应力分析结果。

表2中， f_{\mathrm{c}}^{\mathrm{C}}和f_{\mathrm{t}}^{\mathrm{C}}分别表示损伤交叉点（C点）处的单轴压缩和拉伸应力，f_{\mathrm{c}}^{\mathrm{C}}/f_{\mathrm{c}\mathrm{k}}和 f_{\mathrm{t}}^{\mathrm{C}}/f_{\mathrm{t}\mathrm{m}}分别表示C点处单轴压缩和拉伸应力相对峰值应力的百分比。由表2可知，f_{\mathrm{c}}^{\mathrm{C}}/f_{\mathrm{c}\mathrm{k}}和f_{\mathrm{t}}^{\mathrm{C}}/f_{\mathrm{t}\mathrm{m}}分别在30.40%~52.97%和32.95%~51.41%之间，对应的平均值分别为39.16%和41.07%，与图3中根据文献^[2,5,25,32]综述的C点应力结果基本吻合。

综上可知，损伤交叉点可作为混凝土出现裂缝区域性贯通的临界判据，相应的损伤可被解释为量化混凝土裂缝区域性贯通的损伤阈值。

3 RC构件模拟与损伤分析

3.1 有限元模拟与验证

为验证损伤阈值的合理性，采用ABAQUS软件对王文炜等^[33]开展的RC简支梁试验（CL30）进行数值模拟。构件的几何尺寸和配筋信息如图7所示，试验材料性质参数见表3。

钢筋采用T3D2单元模拟，混凝土采用C3D8R单元模拟。钢筋和混凝土分别采用弹塑性模型和CDP模型模拟。根据已有认识^[34]，对于具有对称性的单调加载情况，当结构破坏主要由对称裂缝主导时，黏结滑移对模拟结果影响不大。因此，采用Embedded模拟钢筋与混凝土的约束关系，即不考虑黏结滑移。RC简支梁有限元模型的两端采用铰约束，在模型的3分点位置处施加位移荷载。模型同时考虑几何和材料非线性，采用自动增量Newton‒Raphson迭代法求解。

为确保有限元模型的可靠性，进行有限元网格敏感性分析，结果如图8所示。由图8可知：75 mm网格模型的精度不够，使得相邻单元之间的云图变化不连续，同时承载力预测偏高。然而，随着网格尺寸的细化，25 mm和50 mm网格模型的模拟结果趋于一致，满足了网格收敛性的要求。考虑到25 mm网格模型的单元数量较多且计算成本较高，本文选择50 mm网格模型作为基准分析模型。

50 mm网格模型模拟与试验的荷载‒位移曲线对比结果如图9所示，可以看出两曲线的整体形态吻合较好。构件的屈服荷载和极限荷载分别为85 kN和102 kN，对应的位移分别为10.02 mm和33.40 mm。在相同位移处，模拟的承载力分别为85.2 kN和98.6 kN，与试验的误差在0.24%~3.33%之间。试验和模拟能量耗散分别为2 727.9和2 756.4 kN·mm，误差为1.04%，验证了模型的可靠性。

3.2 构件损伤结果分析

图10为RC构件有限元模型在不同变形阶段（图9中的P1至P3点，下同）的压缩损伤和拉伸损伤云图及试验结果^[33]。图10中：灰色区域为超过损伤阈值的区域，表示裂缝区域性贯通和贯通裂缝持续扩展直至失效；双箭头表示方向垂直于单元的开裂面，双箭头指向外侧表示单元受拉开裂，双箭头指向内侧表示单元受压开裂。由图10可知，在加载位移为3.21 mm时，梁底跨中的部分单元显示灰色，并且梁底的双箭头指向外侧，基本平行于梁底。这表明在梁底跨中出现了竖向裂缝且发生了区域性的贯通。然而，梁顶跨中的最大压缩损伤为0.19，从图6（a）和（c）可知，当前阶段混凝土处于峰值应力左右，单元尚未出现裂缝贯通。当加载位移增加至10.02 mm时，梁底跨中的灰色云图沿着梁的两端扩展，同时梁跨中的灰色云图也沿着截面高度方向扩展。此时，梁顶跨中的最大压缩损伤为0.47，表明裂缝处于继续发展状态但尚未贯通。当加载位移增加至33.40 mm时，梁底部的灰色云图区域扩大明显，梁底跨中高度方向有3层单元进入裂缝区域性贯通状态；同时，梁顶跨中有2层单元出现局部压溃现象。此外，观察云图中的箭头指向变化可以发现，梁底跨中的箭头指向外侧并与梁底面平行，这表明在梁底跨中部位以竖向裂缝为主。随着加载位移的增加，靠近梁端的箭头指向外侧并与梁轴线的夹角逐步减小，这表明沿着梁端轴线方向的裂缝倾角逐渐减小并以斜向剪切裂缝为主。这些模拟结果与试验结果（图10（d））^[33]相当接近，验证了损伤阈值的合理性。

图11为不同加载位移下模型底层单元从梁左端到梁右端的拉伸损伤与最大主应力的变化过程。由图11可知，受边界条件的影响，梁端附近的单元出现负应力，因此可忽略此处的应力结果。当加载位移为3.21 mm时，梁底跨中有14个单元出现裂缝贯通，相应的最大主应力衰减至1.3 MPa左右；当加载位移增加至10.02 mm时，梁底裂缝贯通单元数量扩展至36个，此时最大主应力衰减至0.6 MPa左右；当加载位移增加至33.40 mm时，梁底裂缝贯通单元数量扩展至42个，此时最大主应力衰减至零应力状态。

图12为不同加载位移下模型顶层单元从梁左端到梁右端的压缩损伤与最大主应力的变化过程。由图12可知，加载施加在模型单元上，导致加载附近的单元出现损伤和应力跳跃，因此可忽略此处的结果。当加载位移为3.21 mm时，梁顶跨中没有单元出现压溃，相应的最大主应力为12.5 MPa；当加载位移为10.02 mm时，梁顶跨中仍没有单元出现压溃，相应的最大主应力从峰值应力衰减至20.0 MPa；当加载位移为33.40 mm时，梁顶跨中14个单元出现压溃，相应的最大主应力从峰值应力衰减至15.0 MPa。

综合图11和12可知：随着加载位移的增加，模型底层单元的拉伸损伤逐渐增加，最大主应力逐渐降低；模型顶层单元中的压缩损伤也随之增加，最大主应力在一定加载范围内先增加后减小。这些结果合理地反映了混凝土在加载过程中的损伤发展和应力变化特征，说明Zhang‒Hou损伤计算模型及损伤阈值能够有效模拟混凝土裂缝区域性贯通的起点及裂缝贯通后的扩展特征。

4 结 论

本文利用Zhang‒Hou损伤计算模型定量分析了不同强度混凝土的裂缝扩展及其损伤演化过程，主要得到以下结论：

1）考虑塑性应变影响的Zhang‒Hou损伤计算模型能够较好地描述混凝土的损伤演化特征，具有较强的普适性。

2）不同强度混凝土的单轴压缩和单轴拉伸损伤演化规律基本一致，在损伤演化过程中发现了单轴压缩和单轴拉伸损伤交叉点，对应的损伤分别为0.57和0.62。该交叉点可作为混凝土出现裂缝区域性贯通的临界判据，相应的损伤可被解释为量化混凝土裂缝区域性贯通的损伤阈值。

3）本文建立的RC简支梁模型的模拟结果满足网格收敛性要求。与试验结果相比，模拟的承载力和能量耗散误差分别小于3.33%和1.04%，验证了有限元模型的可靠性。

4）Zhang‒Hou损伤计算模型和损伤阈值能够有效模拟混凝土裂缝区域性贯通的起点及裂缝贯通后的扩展特征。

参考文献

[1] Liu Zhi, Zhao Lanhao, Wu Xiaobin,et al.Damage model of concrete considering hysteretic effect under cyclic loading[J].Advanced Engineering Sciences,2020,52(4):117‒123.

[2] 刘智,赵兰浩,吴晓彬,等.循环荷载下考虑滞回效应的混凝土损伤模型[J].工程科学与技术,2020,52(4):117‒123.

[3] Yin Yajuan, Ren Qingwen, Shen Lei,et al.Study on crack propagation and damage evolution process of concrete based on fractal dimension[J].Journal of Hydraulic Engineering,2021,52(11):1270‒1280.

[4] 殷亚娟,任青文,沈雷,等.基于分形维的混凝土裂纹扩展及损伤演化过程研究[J].水利学报,2021,52(11):1270‒1280.

[5] Li Junsheng, Ren Wenyuan, Zhang Aijun,et al.Two dimensional meso-compression fracture simulation of concrete based on XCT images[J].Journal of Water Resources and Water Engineering,2021,32(6):170‒177.

[6] 李俊生,任文渊,张爱军,等.基于XCT扫描图像的混凝土二维细观受压断裂模拟[J].水资源与水工程学报,2021,32(6):170‒177.

[7] Ren Huilan, Song Shuizhou, Ning Jianguo.Damage evolution of concrete under tensile load using discrete element modeling[J].Theoretical and Applied Fracture Mechanics,2022,122:103622. doi:10.1016/j.tafmec.2022.103622

[8] Sun Bin, Li Zhaoxia.Adaptive image-based method for integrated multi-scale modeling of damage evolution in heterogeneous concrete[J].Computers & Structures,2015,152:66‒81. doi:10.1016/j.compstruc.2015.02.015

[9] Sun Bin.A numerical method for improving seismic resistance of concrete based on its meso-scale damage mechanisms[J].International Journal of Computational Materials Science and Engineering,2019,8(1):1850031.

[10] Yin Yajuan, Ren Qingwen, Shen Lei.Study on the effect of aggregate distribution on mechanical properties and dama-ge cracks of concrete based on multifractal theory[J].Construction and Building Materials,2020,262:120086. doi:10.1016/j.conbuildmat.2020.120086

[11] Zhang Yuhang, Chen Qingqing, Wang Zhiyong,et al.3D me-soscale fracture analysis of concrete under complex loading[J].Engineering Fracture Mechanics,2019,220:106646. doi:10.1016/j.engfracmech.2019.106646

[12] Shang Xiaoyu, Yang Jingwei, Li Jiangshan.Fractal characteristics of meso-failure crack in recycled coarse aggregate concrete based on CT image[J].Acta Materiae Compositae Sinica,2020,37(7):1774‒1784.

[13] 商效瑀,杨经纬,李江山.基于CT图像的再生混凝土细观破坏裂纹分形特征[J].复合材料学报,2020,37(7):1774‒1784.

[14] Han Yanhua, Xie Zuowei, Liu Fang,et al.Internal micro-deformation analysis of recycled aggregate concrete based on CT images[J].Journal of South China University of Tec-hnology(Natural Science Edition),2019,47(5):73‒80.

[15] 韩燕华,谢作为,刘方,等.基于CT图像的再生混凝土内部微变形演化分析[J].华南理工大学学报(自然科学版),2019,47(5):73‒80.

[16] Jiao Yantao, Cheng Liping.A new anisotropic plastic damage model and its numerical implementation for plain concrete[J].Engineering Mechanics,2022,39(8):122‒137.

[17] 焦延涛,程立平.一种新的混凝土各向异性弹塑性损伤本构模型及其数值实施[J].工程力学,2022,39(8):122‒137.

[18] Guo Chenggong, Li Jie.A new stochastic damage constitutive model of concrete considering strain rate effect[J].Chinese Journal of Theoretical and Applied Mechanics,2022,54(12):3456‒3467.

[19] 虢成功,李杰.考虑应变率效应的混凝土随机损伤本构模型研究[J].力学学报,2022,54(12):3456‒3467.

[20] 谢和平.岩石混凝土损伤力学[M].徐州:中国矿业大学出版社,1990.

[21] Congcong Lü.Study on some problems of numerical model of elastic-plastic damage and permeability coupling of concrete[D].Yangling:Northwest A & F University,2017.

[22] 吕从聪.混凝土弹塑性损伤与渗透耦合数值模型若干问题研究[D].杨凌:西北农林科技大学,2017.

[23] Wang Gaohui.Dynamic response behavior and damage me- chanism of concrete gravity dam under extreme load[D].Tianjin:Tianjin University,2014.

[24] 王高辉.极端荷载作用下混凝土重力坝的动态响应行为和损伤机理[D].天津:天津大学,2014.

[25] Li Zuohua.Damage model and earthquake damage desc-ription of high-rise reinforced concrete structure[D].Harbin:Harbin Institute of Technology,2010.

[26] 李祚华.高层钢筋混凝土结构损伤模型及地震损伤描述[D].哈尔滨:哈尔滨工业大学,2010.

[27] Lee J, Fenves G L.Plastic-damage model for cyclic loading of concrete structures[J].Journal of Engineering Mechanics,1998,124(8):892‒900. doi:10.1061/(asce)0733-9399(1998)124:8(892)

[28] Sidoroff F.Description of anisotropic damage application to elasticity[C]//Physical Non-Linearities in Structural Analysis.Berlin:Springer,1981:237‒244. doi:10.1007/978-3-642-81582-9_35

[29] Cordebois J P, Sidoroff F.Damage induced elastic anisotropy[C]//Mechanical Behavior of Anisotropic Solids/Comportment Méchanique des Solides Anisotropes.Dordrecht:Springer Netherlands,1982:761‒774. doi:10.1007/978-94-009-6827-1_44

[30] FIP C. Model code for concrete structures 2010 [S].Ernst & Sohn,2013. doi:10.1002/9783433604090

[31] Krätzig W B, Pölling R.An elasto-plastic damage model for reinforced concrete with minimum number of material parameters[J].Computers & Structures,2004,82(15/16):1201‒1215. doi:10.1016/j.compstruc.2004.03.002

[32] Li Xiaoqin, Zhang Tian.Calculation method of cracking-closing behaviour of concrete under cyclic loading[J].Jou-rnal of Vibration and Shock,2021,40(9):254‒263.

[33] 李晓琴,张田.循环荷载下混凝土开裂‒闭合行为计算方法研究[J].振动与冲击,2021,40(9):254‒263.

[34] Bažant Z P, Oh B H.Crack band theory for fracture of concrete[J].Matériaux et Construction,1983,16(3):155‒177. doi:10.1007/bf02486267

[35] Sinha B P, Gerstle K H, Tulin L G.Stress‒strain relations for concrete under cyclic loading[J].ACI Journal Proceedings,1964,61(2):195‒212. doi:10.14359/7775

[36] Ren Xiaodan.Experimental study on random damage constitutive relation of concrete[D].Shanghai:Tongji University,2006.

[37] 任晓丹.混凝土随机损伤本构关系试验研究[D].上海:同济大学,2006.

[38] Du Xiuli, Jin Liu.A review on meso-mechanical method for studying the static-mechanical properties of concrete[J].Advances in Mechanics,2011,41(4):411‒426.

[39] 杜修力,金浏.混凝土静态力学性能的细观力学方法述评[J].力学进展,2011,41(4):411‒426.

[40] Trawiński W, Bobiński J, Tejchman J.Two-dimensional si-mulations of concrete fracture at aggregate level with cohesive elements based on X-ray μCT images[J].Engineering Fracture Mechanics,2016,168:204‒226. doi:10.1016/j.engfracmech.2016.09.012

[41] Grassl P, Jirásek M.Meso-scale approach to modelling the fracture process zone of concrete subjected to uniaxial tension[J].International Journal of Solids and Structures,2010,47(7/8):957‒968. doi:10.1016/j.ijsolstr.2009.12.010

[42] Ren Wenyuan, Yang Zhenjun, Sharma R,et al.Two dimensional X-ray CT image based meso-scale fracture modelling of concrete[J].Engineering Fracture Mechanics,2015,133:24‒39. doi:10.1016/j.engfracmech.2014.10.016

[43] Skarżyński, Tejchman J.Experimental investigations of fr- acture process in concrete by means of X-ray micro-com-puted tomography[J].Strain,2016,52(1):26‒45. doi:10.1111/str.12168

[44] Sun Bin, Li Zhaoxia.Adaptive concurrent three-level mu-ltscale simulation for trans-scale process from material mesodamage to structural failure of concrete structures[J].International Journal of Damage Mechanics,2016,25(5):750 ‒769. doi:10.1177/1056789516648371

[45] 李杰,吴建营,陈建兵.混凝土随机损伤力学[M].北京:科学出版社,2014.

[46] Wang Wenwei, Zhao Guofan, Huang Chengkui,et al.An experimental study of strengthening of initially loaded reinforced concrete beams using CFRP sheets[J].Engineering Mechanics,2004,21(4):172‒178.

[47] 王文炜,赵国藩,黄承逵,等.碳纤维布加固已承受荷载的钢筋混凝土梁抗弯性能试验研究及抗弯承载力计算[J].工程力学,2004,21(4):172‒178.

[48] Li Xiaoqin, Zhang Tian, Ding Zude,et al.Numerical analysis of normal concrete lining strengthening methods under different damage levels[J].Structure and Infrastructure Engineering,2021,17(12):1597‒1611. doi:10.1080/15732479.2020.1817106