留言板

尊敬的读者、作者、审稿人, 关于本刊的投稿、审稿、编辑和出版的任何问题, 您可以本页添加留言。我们将尽快给您答复。谢谢您的支持!

姓名
邮箱
手机号码
标题
留言内容
验证码

基于粘结单元的三维随机细观混凝土离散断裂模拟

杨贞军 黄宇劼 尧锋 刘国华

杨贞军, 黄宇劼, 尧锋, 刘国华. 基于粘结单元的三维随机细观混凝土离散断裂模拟[J]. 工程力学. doi: 10.6052/j.issn.1000-4750.2019.09.0559
引用本文: 杨贞军, 黄宇劼, 尧锋, 刘国华. 基于粘结单元的三维随机细观混凝土离散断裂模拟[J]. 工程力学. doi: 10.6052/j.issn.1000-4750.2019.09.0559
Zhen-jun YANG, Yu-jie HUANG, Feng YAO, Guo-hua LIU. THREE-DIMENSIONAL MESO-SCALE COHESIVE FRACTURE MODELING OF CONCRETE USING A PYTHON SCRIPT IN ABAQUS[J]. Engineering Mechanics. doi: 10.6052/j.issn.1000-4750.2019.09.0559
Citation: Zhen-jun YANG, Yu-jie HUANG, Feng YAO, Guo-hua LIU. THREE-DIMENSIONAL MESO-SCALE COHESIVE FRACTURE MODELING OF CONCRETE USING A PYTHON SCRIPT IN ABAQUS[J]. Engineering Mechanics. doi: 10.6052/j.issn.1000-4750.2019.09.0559

基于粘结单元的三维随机细观混凝土离散断裂模拟

doi: 10.6052/j.issn.1000-4750.2019.09.0559
基金项目: 国家自然科学基金项目(51779222,51974202)
详细信息
    作者简介:

    黄宇劼(1988−),男,福建人,工程师,博士,主要从事混凝土多尺度模拟研究(E-mail: mpcyujie@126.com)

    尧 锋(1991−),男,江西人,博士生,主要从事混凝土数值模拟研究(E-mail: yf1991x@zju.edu.cn)

    刘国华(1963−),男,福建人,教授,硕士,主要从事水工结构与优化分析研究(E-mail: zjuliu@163.com)

    通讯作者: 杨贞军(1974−),男,湖北人,教授,博士,主要从事结构工程和数值模拟研究(E-mail: zhjyang@whu.edu.cn)
  • 中图分类号: TU528

THREE-DIMENSIONAL MESO-SCALE COHESIVE FRACTURE MODELING OF CONCRETE USING A PYTHON SCRIPT IN ABAQUS

  • 摘要: 开发MATLAB 和ABAQUS Python脚本程序,建立含随机三维多面体骨料的混凝土细观有限元模型。通过自编高效C++程序插设三维零厚度的粘结界面单元,以模拟复杂三维离散多裂缝的起裂与扩展。对典型混凝土试件单轴拉伸断裂的模拟,分析粘性界面单元的主要材料参数(抗拉强度与断裂能)对应力-位移曲线、断裂过程、裂缝面特征的影响。结果表明:宏观应力-位移曲线主要受砂浆、界面粘性裂缝单元的抗拉强度、断裂能绝对数值的影响;裂缝面的位置与形态决定于砂浆、界面粘结裂缝单元的抗拉强度相对比值以及断裂能相对比值;混凝土的力学响应反映其裂缝发展特征,二者既决定于断裂材料参数,也受到骨料大小、形状等细观结构因素的影响;建立的有限元模型能有效描述混凝土复杂三维断裂过程。
  • 图  1  不同形态的多面体骨料

    Figure  1.  Polyhedral aggregates of different shapes

    图  2  两个骨料之间的三种位置关系

    Figure  2.  Spatial relations of two aggregates

    图  3  在ABAQUS中通过Python脚本程序生成多面体

    Figure  3.  Polyhedra generation using Python scripts in ABAQUS

    图  4  生成的50 mm试件多面体骨料有限元模型

    Figure  4.  FEM model of a specimen with polyhedral aggregates

    图  5  细观混凝土三维粘结单元的插设算法

    Figure  5.  3D inserting algorithm of cohesive elements

    图  6  骨料表面的粘结单元(CIE_INT)

    Figure  6.  Cohesive element on the surface of an aggregate

    图  7  单轴受拉混凝土试件的边界条件

    Figure  7.  Boundary conditions under uniaxial tension

    图  8  单轴受拉混凝土试件的宏观应力-位移曲线

    Figure  8.  Average stress - displacement curves of concrete specimen under uniaxial tension

    图  9  单轴受拉开裂过程:粘结裂缝单元的受拉变形

    Figure  9.  Fracture process under uniaxial tension: the red cohesive elements with D≥0.9

    图  10  单轴受拉开裂过程

    Figure  10.  Crack propagation process under uniaxial tension

    图  11  粘结裂缝单元的抗拉强度对宏观应力-位移曲线的影响

    Figure  11.  Average stress - displacement curves: effects of fracture energy of cohesive elements

    图  12  粘结裂缝单元的断裂能对宏观应力-位移曲线的影响

    Figure  12.  Average stress - displacement curves: effects of fracture energy of cohesive elements

    图  13  粘结裂缝单元的抗拉强度对试件破坏形式的影响

    注:(a)、(b)只改变砂浆粘结裂缝单元CIE_CEM的抗拉强度;(c)、(d)只改变界面粘结裂缝单元CIE_INT的抗拉强度

    Figure  13.  Effects of tensile strength of cohesive elements on the fracture surfaces

    图  14  粘结裂缝单元的断裂能对试件破坏形式的影响

    注:(a)、(b)只改变砂浆粘结裂缝单元CIE_CEM的断裂能;(c)、(d)只改变界面粘结裂缝单元CIE_INT的断裂能

    Figure  14.  Effects of fracture energy on crack surfaces

    图  15  断裂能对开裂过程的影响

    Figure  15.  Effects of fracture energy on crack propagation processes

    图  16  断裂能对应力-位移曲线的影响

    Figure  16.  Effects of fracture energy on stress-displacement curves

    表  1  三级配骨料尺寸分布[29]

    Table  1.   Three-graded aggregate size distribution[29]

    筛孔尺寸/mm累计通过率/(%)
    12.70100
    9.5077
    4.7526
    2.360
    下载: 导出CSV

    表  2  材料参数[22]

    Table  2.   Material parameters[22]

    弹性模量E/GPa泊松比密度ρ/(×10−6 kg/mm3)初始刚度kn /(MPa/mm)抗拉强度tn /MPa断裂能Gf /(N/mm)
    骨料700.22.5
    砂浆250.22.2
    CIE_CEM2.26×10640.06
    CIE_INT2.26×10620.03
    下载: 导出CSV
  • [1] 杜修力, 金浏. 混凝土静态力学性能的细观力学方法述评[J]. 力学进展, 2011, 41(4): 411 − 426. doi:  10.6052/1000-0992-2011-4-lxjzJ2011-010

    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. (in Chinese) doi:  10.6052/1000-0992-2011-4-lxjzJ2011-010
    [2] 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): 957 − 968.
    [3] 梁诗雪, 李杰. 基于两相随机介质的混凝土破坏全过程模拟[J]. 工程力学, 2018, 35(2): 116 − 123.

    Liang Shixue, Li Jie. A two-phase random medium model for simulating the failure process of concrete [J]. Engineering Mechanics, 2018, 35(2): 116 − 123. (in Chinese)
    [4] 金浏, 苏晓, 杜修力. 钢筋混凝土梁受弯破坏及尺寸效应的细观模拟分析[J]. 工程力学, 2018, 35(10): 27 − 36. doi:  10.6052/j.issn.1000-4750.2017.06.0436

    Jin Liu, Su Xiao, Du Xiuli. Meso-scale simulations on flexural failure and size effect of reinforced concrete beams [J]. Engineering Mechanics, 2018, 35(10): 27 − 36. (in Chinese) doi:  10.6052/j.issn.1000-4750.2017.06.0436
    [5] Yang Z, Li B, Wu J. X-ray computed tomography images based phase-field modeling of mesoscopic failure in concrete [J]. Engineering Fracture Mechanics, 2019, 208: 151 − 170. doi:  10.1016/j.engfracmech.2019.01.005
    [6] 李冬, 金浏, 杜修力, 刘晶波, 张帅, 余文轩. 考虑细观组分影响的混凝土宏观力学性能理论预测模型[J]. 工程力学, 2019, 36(05): 67 − 75.

    Li Dong, Jin Liu, Du Xiuli, Liu Jinbo, Zhang Shuai, Yu Wenxuan. A theoretical prediction model of concrete macroscopic mechamoca properties considering the influence of mesoscopic composition [J]. Engineering Mechanics, 2019, 36(05): 67 − 75. (in Chinese)
    [7] Man H K, van Mier J G. Influence of particle density on 3D size effects in the fracture of (numerical) concrete [J]. Mechanics of Materials, 2008, 40(6): 470 − 486. doi:  10.1016/j.mechmat.2007.11.003
    [8] Xu W X., Chen H S. Numerical investigation of effect of particle shape and particle size distribution on fresh cement paste microstructure via random sequential packing of dodecahedral cement particles [J]. Computers & Structures, 2013, 114: 35 − 45.
    [9] Zhang Z, Song X, Liu Y, et al. Three-dimensional mesoscale modelling of concrete composites by using random walking algorithm [J]. Composites Science and Technology, 2017, 149: 235 − 245. doi:  10.1016/j.compscitech.2017.06.015
    [10] Yan P, Zhang J, Fang Q, Zhang Y, Fan J. 3D numerical modelling of solid particles with randomness in shape considering convexity and concavity [J]. Powder Technology, 2016, 301: 131 − 140. doi:  10.1016/j.powtec.2016.06.007
    [11] Xu W X, Lv Z, Chen H S. Effects of particle size distribution, shape and volume fraction of aggregates on the wall effect of concrete via random sequential packing of polydispersed ellipsoidal particles [J]. Physica A: Statistical Mechanics and its Applications, 2013, 392(3): 416 − 426. doi:  10.1016/j.physa.2012.09.014
    [12] Li X, Xu Y, Chen S. Computational homogenization of effective permeability in three-phase mesoscale concrete [J]. Construction and Building Materials, 2016, 121: 100 − 111. doi:  10.1016/j.conbuildmat.2016.05.141
    [13] Abyaneh S D, Wong H S, Buenfeld N R. Modelling the diffusivity of mortar and concrete using a three-dimensional mesostructure with several aggregate shapes [J]. Computational Materials Science, 2013, 78: 63 − 73. doi:  10.1016/j.commatsci.2013.05.024
    [14] Liu L, Shen D, Chen H, Xu W. Aggregate shape effect on the diffusivity of mortar: A 3D numerical investigation by random packing models of ellipsoidal particles and of convex polyhedral particles [J]. Computers & Structures, 2014, 144: 40 − 51.
    [15] 刘光廷, 高政国. 三维凸型混凝土骨料随机投放算法[J]. 清华大学学报(自然科学版), 2003, 43(8): 1120 − 1123.

    Liu Guanting, Gao Zhengguo. Random 3-D aggregate structure for concrete [J]. Journal of Tsinghua University (Science and Technology), 2003, 43(8): 1120 − 1123. (in Chinese)
    [16] Häfner S, Eckardt S, Luther T, Konke C. Mesoscale modeling of concrete: Geometry and numerics [J]. Computers & Structures, 2006, 84(7): 450 − 461.
    [17] Sheng P, Zhang J, Ji Z. An advanced 3D modeling method for concrete-like particle-reinforced composites with high volume fraction of randomly distributed particles [J]. Composites Science and Technology, 2016, 134: 26 − 35. doi:  10.1016/j.compscitech.2016.08.009
    [18] 方秦, 张锦华, 还毅, 张亚栋. 全级配混凝土三维细观模型的建模方法研究[J]. 工程力学, 2013, 30(1): 14 − 21.

    Fang Qin, Zhang Jinhua, Huan Yi, Zhang Yadong. The investigation into three-dimensional mesoscale modelling of fully-graded concrete [J]. Engineering Mechanics, 2013, 30(1): 14 − 21. (in Chinese)
    [19] Du C, Sun L, Jiang S, Ying Z. Numerical simulation of aggregate shapes of three-dimensional concrete and its applications [J]. Journal of Aerospace Engineering, 2011, 26(3): 515 − 527.
    [20] Huang Y, Yang Z, Ren W, Liu G, Zhang C. 3D meso-scale fracture modelling and validation of concrete based on in-situ X-ray computed tomography images using damage plasticity model [J]. International Journal of Solids and Structures, 2015, 67: 340 − 352.
    [21] Xu Y, Chen S. A method for modeling the damage behavior of concrete with a three-phase mesostructure [J]. Construction and Building Materials, 2016, 102: 26 − 38. doi:  10.1016/j.conbuildmat.2015.10.151
    [22] Caballero A, López C M, Carol I. 3D meso-structural analysis of concrete specimens under uniaxial tension [J]. Computer Methods in Applied Mechanics and Engineering, 2006, 195(52): 7182 − 7195. doi:  10.1016/j.cma.2005.05.052
    [23] Yin A, Yang X, Yang Z. 2D and 3D fracture modeling of asphalt mixture with randomly distributed aggregates and embedded cohesive cracks [J]. Procedia IUTAM, 2013, 6: 114 − 12. doi:  10.1016/j.piutam.2013.01.013
    [24] Fang J, Pan Y, Dang F, et al. Numerical reconstruction model and simulation study of concrete based on damaged partition theory and CT number [J]. Materials, 2019, 12(24): 4070. doi:  10.3390/ma12244070
    [25] Simulia D S. Abaqus 6.11 theory manual [S]. Providence, RI, USA: DS SIMULIA Corp, 2011.
    [26] Yang Z, Su X, Chen J, Liu G. Monte Carlo simulation of complex cohesive fracture in random heterogeneous quasi-brittle materials [J]. International Journal of Solids and Structures, 2009, 46(17): 3222 − 3234. doi:  10.1016/j.ijsolstr.2009.04.013
    [27] Su X, Yang Z, Liu G. Monte Carlo simulation of complex cohesive fracture in random heterogeneous quasi-brittle materials: A 3D study [J]. International Journal of Solids and Structures, 2010, 47: 2336 − 2345. doi:  10.1016/j.ijsolstr.2010.04.031
    [28] Su X, Yang Z, Liu G. Finite element modelling of complex 3D static and dynamic crack propagation by embedding cohesive elements in Abaqus [J]. Acta Mechanica Solida Sinica, 2010, 23(3): 271 − 282. doi:  10.1016/S0894-9166(10)60030-4
    [29] Hirsch T J. Modulus of elasticity of concrete affected by elastic moduli of cement paste matrix and aggregate [C]//Journal Proceedings, 1962, 59(3): 427 − 452.
    [30] Hordijk D A. Tensile and tensile fatigue behaviour of concrete: experiments, modelling and analyses [J]. Heron, 1992, 37: 1 − 79.
  • [1] 王磊, 班慧勇, 石永久, 王元清.  基于微观断裂机理的高强钢框架梁柱节点抗震性能有限元分析 . 工程力学, doi: 10.6052/j.issn.1000-4750.2017.08.0608
    [2] 段进涛, 史旦达, 汪金辉, 焦宇, 何佩珊.  火灾环境下钢结构响应行为的FDS-ABAQUS热力耦合方法研究 . 工程力学, doi: 10.6052/j.issn.1000-4750.2016.01.0061
    [3] 池寅, 黄乐, 余敏.  基于ABAQUS的钢-聚丙烯混杂纤维混凝土损伤塑性本构模型取值方法研究 . 工程力学, doi: 10.6052/j.issn.1000-4750.2016.08.0630
    [4] 赵作周, 周剑, 侯建群, 任宝双.  上下层插筋连接预制混凝土空心模剪力墙有限元分析 . 工程力学, doi: 10.6052/j.issn.1000-4750.2015.05.0411
    [5] 党像梁, 吕西林, 钱江, 周颖.  底部开水平缝预应力自复位剪力墙有限元模拟 . 工程力学, doi: 10.6052/j.issn.1000-4750.2015.12.1021
    [6] 黄维, 钱江, 周知.  考虑混凝土不同约束效应的型钢混凝土柱抗震性能模拟研究 . 工程力学, doi: 10.6052/j.issn.1000-4750.2014.10.0831
    [7] 王明强, 杨军.  HISS模型在ABAQUS中的三维实现、验证和应用 . 工程力学, doi: 10.6052/j.issn.1000-4750.2012.09.0666
    [8] 方秦, 还毅, 陈力, 柳锦春.  应变速率型RC梁柱显式分析单元及其在ABAQUS软件中的实现 . 工程力学, doi: 10.6052/j.issn.1000-4750.2011.08.0545
    [9] 陈 成, 邵永波, 杨 杰.  T型圆钢管节点抗火性能的有限元研究 . 工程力学, doi: 10.6052/j.issn.1000-4750.2011.05.0330
    [10] 聂建国, 王宇航.  ABAQUS中混凝土本构模型用于模拟结构静力行为的比较研究 . 工程力学, doi: 10.6052/j.issn.1000-4750.2011.07.0420
    [11] 王国林, 杜小伟, 朱美林, 应世洲.  轮胎成型有限元仿真研究 . 工程力学, doi: 10.6052/j.issn.1000-4750.2010.08.0599
    [12] 聂建国, 王宇航.  基于ABAQUS的钢-混凝土组合结构纤维梁模型的开发及应用 . 工程力学, doi: 10.6052/j.issn.1000-4750.2010.04.0278
    [13] 段红霞, 李守巨, 刘迎曦.  地震作用下钢结构损伤过程数值模拟 . 工程力学,
    [14] 曹 鹏, 冯德成, 田 林, 荆儒鑫.  基于弹塑性损伤理论的水泥稳定基层养生期裂缝形成机理分析 . 工程力学,
    [15] 吴 熙, 付腾飞, 吴智敏.  自密实轻骨料混凝土的双K断裂参数和断裂能试验研究 . 工程力学,
    [16] 李 强, 杨 庆, 栾茂田, 贾景超.  张开裂纹闭合的数值模拟及简化判据研究 . 工程力学,
    [17] 毛江徽, 杨显杰, 高 庆.  微电子封装钎料时相关本构模型的有限元实现 . 工程力学,
    [18] 陈 瑛, 乔丕忠, 姜弘道, 任青文.  FRP-混凝土三点受弯梁损伤粘结模型有限元分析 . 工程力学,
    [19] 沈新普, 王琛元, 周 琳.  一个钢筋混凝土损伤塑性本构模型及工程应用 . 工程力学,
    [20] 韩菊红, 赵国藩, 张雷顺.  新老混凝土粘结面断裂损伤过程区研究 . 工程力学,
  • 加载中
计量
  • 文章访问数:  104
  • HTML全文浏览量:  36
  • PDF下载量:  1
  • 被引次数: 0
出版历程
  • 收稿日期:  2019-09-24
  • 修回日期:  2020-02-12
  • 网络出版日期:  2020-06-02

基于粘结单元的三维随机细观混凝土离散断裂模拟

doi: 10.6052/j.issn.1000-4750.2019.09.0559
    基金项目:  国家自然科学基金项目(51779222,51974202)
    作者简介:

    黄宇劼(1988−),男,福建人,工程师,博士,主要从事混凝土多尺度模拟研究(E-mail: mpcyujie@126.com)

    尧 锋(1991−),男,江西人,博士生,主要从事混凝土数值模拟研究(E-mail: yf1991x@zju.edu.cn)

    刘国华(1963−),男,福建人,教授,硕士,主要从事水工结构与优化分析研究(E-mail: zjuliu@163.com)

    通讯作者: 杨贞军(1974−),男,湖北人,教授,博士,主要从事结构工程和数值模拟研究(E-mail: zhjyang@whu.edu.cn)
  • 中图分类号: TU528

摘要: 开发MATLAB 和ABAQUS Python脚本程序,建立含随机三维多面体骨料的混凝土细观有限元模型。通过自编高效C++程序插设三维零厚度的粘结界面单元,以模拟复杂三维离散多裂缝的起裂与扩展。对典型混凝土试件单轴拉伸断裂的模拟,分析粘性界面单元的主要材料参数(抗拉强度与断裂能)对应力-位移曲线、断裂过程、裂缝面特征的影响。结果表明:宏观应力-位移曲线主要受砂浆、界面粘性裂缝单元的抗拉强度、断裂能绝对数值的影响;裂缝面的位置与形态决定于砂浆、界面粘结裂缝单元的抗拉强度相对比值以及断裂能相对比值;混凝土的力学响应反映其裂缝发展特征,二者既决定于断裂材料参数,也受到骨料大小、形状等细观结构因素的影响;建立的有限元模型能有效描述混凝土复杂三维断裂过程。

English Abstract

杨贞军, 黄宇劼, 尧锋, 刘国华. 基于粘结单元的三维随机细观混凝土离散断裂模拟[J]. 工程力学. doi: 10.6052/j.issn.1000-4750.2019.09.0559
引用本文: 杨贞军, 黄宇劼, 尧锋, 刘国华. 基于粘结单元的三维随机细观混凝土离散断裂模拟[J]. 工程力学. doi: 10.6052/j.issn.1000-4750.2019.09.0559
Zhen-jun YANG, Yu-jie HUANG, Feng YAO, Guo-hua LIU. THREE-DIMENSIONAL MESO-SCALE COHESIVE FRACTURE MODELING OF CONCRETE USING A PYTHON SCRIPT IN ABAQUS[J]. Engineering Mechanics. doi: 10.6052/j.issn.1000-4750.2019.09.0559
Citation: Zhen-jun YANG, Yu-jie HUANG, Feng YAO, Guo-hua LIU. THREE-DIMENSIONAL MESO-SCALE COHESIVE FRACTURE MODELING OF CONCRETE USING A PYTHON SCRIPT IN ABAQUS[J]. Engineering Mechanics. doi: 10.6052/j.issn.1000-4750.2019.09.0559
  • 传统的宏观混凝土模型假设均质材料,不能反映细观组分空间分布的随机性和非均质性引起的损伤与断裂局部化,难以可靠预测材料的力学响应与裂缝扩展过程。因此进行混凝土细观断裂模拟研究十分必要[1-5]。在三维细观混凝土模拟中,假设骨料为球体是一种常用的简化建模方法[4, 6-7],虽然可以高效地生成细观结构,但无法反映骨料形态的影响。近年来,采用更为复杂的椭球和凸多面体来模拟骨料日趋常见[8-9]。在应用上,复杂骨料模型大多用于几何特征研究[8, 10]、边界效应分析[11]、渗透性能研究[12-14]以及线弹性应力分析[15-17];在三维复杂非线性损伤和断裂的模拟方面仍非常有限,且大多采用连续介质损伤塑性模型[18-21],难以定量预测裂缝宽度,因此有必要采用离散裂缝模型进行更精细的断裂模拟。另一方面,在骨料几何建模方法上,一般通过自编程序进行参数化建模[18, 21],或基于Voronoi图生成细观结构[22],或基于ANSYS的ADPL命令流建模[23]。但是,对三维复杂细观结构进行网格划分并不容易,较为直接的方法是将细观结构映射到规则的立方体背景网格上[4, 18-20, 24],但要保证精度则会导致巨大的单元数量。ABAQUS具有强大的前处理模块,可以灵活有效地对复杂结构进行网格自动划分,不过尚未见到将其用于三维混凝土细观结构建模与网格划分的研究。

    本文首先采用Python编写程序,利用ABAQUS脚本接口[25]进行前处理二次开发,考虑骨料实际级配,有效地建立多面体随机骨料几何模型并进行有限元网格划分。然后通过自编C++程序在骨料-砂浆界面上以及砂浆中高效插设零厚度的离散粘结裂缝单元[26-28],模拟受拉试件复杂细观三维多裂缝的起裂与扩展直至破坏的全过程,并对粘性裂缝单元的主要材料参数(抗拉强度与断裂能)对荷载-位移曲线、断裂过程、裂缝面特征的影响进行详细分析,以加深对混凝土细观断裂机理的认识。应该指出,笔者论文[26-28]也采用了插设界面单元的方法,但均非基于本文提出的随机骨料模型;虽然文献[26-27]采用随机场来间接表征混凝土的空间非均质性,但和文献[28]一样,仅提出了均质材料中插设二维、三维界面的方法。本文是前述工作的扩展,并将Matlab程序改进为C++,大大提高了效率,更适合大规模三维有限元网格的处理。

    • 算法分两步,第一步建立几何模型,第二步建立有限元模型。

      首先编写MATLAB程序进行多面体随机骨料的生成与投放,建立几何模型。采用基于球体的方式产生凸多面体,即凸多面体的各顶点位于球面上。该球体称为基球,其球面上随机产生8个~25个点作为多面体的顶点,顶点i的笛卡尔坐标用Euler空间的极角αi和方位角βi表示为:

      $$\tag{1a}x=\cos ({\alpha _i})\sin ({\beta _i}){r_0} + {x_0}$$
      $$\tag{1b}y=\sin ({\alpha _i})\cos ({\beta _i}){r_0} + {y_0}$$
      $$\tag{1c}{\textit {z}}=\cos ({\beta _i}){r_0} + {{\textit {z}}_0}\quad\quad\quad$$

      式中:(x0, y0, z0)是球心的坐标;r0是球半径;αi在[0, 2π]、βi在[0, π]范围内按均匀分布随机产生。然后采用MATLAB算法“convhulln”基于顶点建立三角形面,这些面构成骨料多面体。为避免顶点距离太近而造成三角形面畸形,需保证球面上的顶点间距不小于ξr0,本文取ξ=0.5。多面体的生成较为耗时;为提高效率,本文不采用边投放、边生成骨料的方式,而是通过生成足够数量(本文为20000)的直径为单位一的多面体,预先建立骨料形态数据库,图1为一些骨料例子。

      图  1  不同形态的多面体骨料

      Figure 1.  Polyhedral aggregates of different shapes

      针对给定模拟区域,具体算法如下:

      1) 根据所需骨料含量,采用实际骨料级配表1求解各级配的骨料体积,从骨料直径最大的级配开始分级配生成并投放骨料;

      表 1  三级配骨料尺寸分布[29]

      Table 1.  Three-graded aggregate size distribution[29]

      筛孔尺寸/mm累计通过率/(%)
      12.70100
      9.5077
      4.7526
      2.360

      2) 对每个级配的操作步骤:随机产生一个球心坐标作为投放位置,从骨料数据库中随机选取一个骨料,并在此级配区间内按均匀分布随机赋予其直径;

      3) 若该骨料不与已有骨料相交或重叠,则该骨料投放成功,进行下一个投放,否则重新进行本次投放;

      4) 本级配骨料体积达到后,进行下一级配投放;

      5) 循环2)步~4)步,直至达到总骨料含量。

      在上述算法第3)步中,新生成骨料与已有骨料之间可能存在如图2所示的3种位置关系:不相交也不重叠、相交、重叠(包含)。

      图  2  两个骨料之间的三种位置关系

      Figure 2.  Spatial relations of two aggregates

      通过下述算法判断两个多面体相交或重叠:如图2(a)所示,对于新生成的多面体A和已有多面体B,要使二者不相交且不重叠,则B的各面需在A各面的同一侧。A的顶点A1~A3所构成平面A123的方程是:

      $$F=\left| {\begin{array}{*{20}{c}} x&y&{\textit{z}}&1 \\ {{x_{{\rm{A}}1}}}&{{y_{{\rm{A}}1}}}&{{{\textit{z}}_{{\rm{A}}1}}}&1 \\ {{x_{{\rm{A}}2}}}&{{y_{{\rm{A}}2}}}&{{{\textit{z}}_{{\rm{A}}2}}}&1 \\ {{x_{{\rm{A}}3}}}&{{y_{{\rm{A}}3}}}&{{{\textit{z}}_{{\rm{A}}3}}}&1 \end{array}} \right|=0$$ (2)

      式中:(xAi, yAi, zAi)是顶点Ai(i=1~3)的坐标,则对于平面A123外的点(x, y, z),有F(x, y, z) > 0或F(x, y, z)<0。则上述位置判断问题可以等效为:遍历A的各面,如果存在一个面(例如A123),使得A的内部各点(凸多面体可只选取其形心)与B的任意顶点(xBj, yBj, zBj)位于面A123的不同侧,则可判定A与B既不相交也不重叠。相比于遍历A、B的顶点以保证各自顶点不在另一个多面体内部的算法,上述位置判断算法更为直观和有效。若A的形心坐标为(xA0, yA0, zA0),A与B既不相交也不重叠时只需满足下式:

      $$F({x_{{\rm{A}}0}},{\rm{ }}{y_{{\rm{A}}0}},{\rm{ }}{{\textit{z}}_{{\rm{A}}0}}) \cdot F({x_{{\rm{B}}j}},{\rm{ }}{y_{{\rm{B}}j}},{\rm{ }}{{\textit{z}}_{{\rm{B}}j}})<0$$ (3)

      骨料与骨料、骨料与边界之间预设最小间距(本文设置为较大骨料直径的0.05倍)以方便网格生成。

      第二步,编写Python脚本程序,进行ABAQUS前处理二次开发,建立多面体骨料并完成模型网格划分。具体算法如下:

      1) 建立模型数据库;

      2) 生成混凝土试件立方体部件(part);

      3) 读取随机骨料几何模型中所有骨料的信息:基球的球心坐标、顶点坐标、各面的顶点构成;

      4) 对于每个骨料,生成三维球体部件;

      5) 延拓骨料的各面,切割该球体,过程见图3

      6) 在assembly模块中创建部件实例(instance),采用布尔操作合并生成的多面体与试件立方体;

      7) 重复上述步骤,直至创建所有骨料,见图4

      图  3  在ABAQUS中通过Python脚本程序生成多面体

      Figure 3.  Polyhedra generation using Python scripts in ABAQUS

      图  4  生成的50 mm试件多面体骨料有限元模型

      Figure 4.  FEM model of a specimen with polyhedral aggregates

      8) 布置全局网格种子,采用“自由网格”技术用四面体单元进行网格划分。

      以边长为50 mm立方体试件和骨料含量30%为例,上述算法建立的骨料有限元模型见图4。该模型含有423个骨料,网格单元平均尺寸为2 mm。

    • 本文将ABAQUS 零厚度三维六节点粘结裂缝单元COH3D6插设到砂浆的有限元实体单元之间,以及骨料-砂浆界面上,用于模拟混凝土复杂非线性断裂过程。采用损伤系数D [0, 1]描述裂缝单元达到强度后的损伤程度,反映刚度退化,该系数是裂缝单元有效相对位移δm的函数:

      $${\delta _{\rm{m}}}=\sqrt {<{\delta _{\rm{n}}}{ > ^2} + \delta _{\rm{s}}^2 + \delta _{\rm{t}}^2} $$ (4)

      式中:δnδtδs分别是法向、切向的相对位移;< >为Macaulay括号,表示为:

      $$ < {\delta _{\rm{n}}} > = \left\{ {\begin{array}{*{20}{c}} {{\delta _{\rm{n}}}{\rm{ , }}\;\;\;{\delta _{\rm{n}}} \geqslant 0}\\ {0{\rm{ , }}\;\;\;{\delta _{\rm{n}}} < 0} \end{array}} \right.{\rm{ }}\begin{array}{*{20}{c}} {{\text{受拉}}}\\ {{\text{受压}}} \end{array} $$ (5)

      以线性软化准则为例,损伤系数可以写成:

      $$D=\frac{{\delta _{{\rm{mf}}}^{}(\delta {{_{\rm{m}}^{}}_{{\rm{,max}}}} - \delta _{{\rm{m}}0}^{})}}{{\delta {{_{\rm{m}}^{}}_{,\max }}(\delta _{{\rm{mf}}}^{} - \delta _{{\rm{m}}0}^{})}}$$ (6)

      式中:δm,max是加载历史中的最大有效相对位移;δm0δmf 分别是裂缝起裂和完全破坏时的有效相对位移。

      法向刚度kn与切向刚度kskt分别表示为:

      $$\tag{7a}{k_{\rm{n}}}=(1 - D){k_{{\rm{n}}0}}$$
      $$\tag{7b}{k_{\rm{s}}}=(1 - D){k_{{\rm{s}}0}}$$
      $$\tag{7c}{k_{\rm{t}}}=(1 - D){k_{{\rm{t}}0}}$$

      相应的应力则可以表示为:

      $$\tag{8a}{t_{\rm{n}}}=\left\{ {\begin{array}{*{20}{c}} {(1 - D){k_{{\rm{n}}0}}{\delta _{\rm{n}}},\;\;\;\;\;\;\;\;\;{\delta _{\rm{n}}} \geqslant {\rm{0 }}} \\ {{k_{{\rm{n}}0}}{\delta _{\rm{n}}},\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;{\rm{ }}{\delta _{\rm{n}}}<{\rm{0 }}} \end{array}} \right.$$
      $$\tag{8b}{t_{\rm{s}}}=(1 - D){k_{{\rm{s}}0}}{\delta _{\rm{s}}}\;\;\qquad\qquad\qquad$$
      $$\tag{8c}{t_{\rm{t}}}=(1 - D){k_{{\rm{t}}0}}{\delta _{\rm{t}}}\;\;\qquad\qquad\qquad$$

      采用名义应力平方准则作为起裂准则:

      $${\left\{ {\frac{{\langle {t_{\rm{n}}}\rangle }}{{{t_{{\rm{n0}}}}}}} \right\}^2} + {\left\{ {\frac{{{t_{\rm{s}}}}}{{{t_{{\rm{s0}}}}}}} \right\}^2} + {\left\{ {\frac{{{t_{\rm{t}}}}}{{{t_{{\rm{t0}}}}}}} \right\}^2}=1$$ (9)

      式中,tntstt分别为法向、切向的应力。

      文献[26-27]报道了在二维、三维实体有限元单元网格中插设粘结单元的算法和MATLAB程序实施,用于均质材料。本文将该算法加以拓展,通过编写C++程序,在混凝土细观结构中插设三维粘结单元,图5给出了算法示意图。假设骨料不会开裂,骨料内无需插设粘结单元,因此只在骨料-砂浆界面以及砂浆中采用粘结单元(图5(c)),分别标记为CIE_INT与CIE_CEM。

      图  5  细观混凝土三维粘结单元的插设算法

      Figure 5.  3D inserting algorithm of cohesive elements

      具体插设算法如下:

      1) 读取初始网格,图5(a)显示了网格中与骨料表面相连的2个砂浆四面体单元A与B,以及1个与砂浆单元B相连的砂浆单元C。

      2) 如图5(b)所示,搜寻单元的面,骨料的面存入结构体FACE_INT中,砂浆的面存入FACE_CEM中。

      3) 如图5(b)所示,将节点按其所在单元的数量拆开(即不改变坐标地进行节点的原位复制)并重新编号,再分配给各单元。

      4) 如图5(c)所示,根据FACE_INT、FACE_CEM分别建立粘结单元CIE_INT、CIE_CEM,从而完成粘结单元的插设,图6给出一个骨料表面插设的零厚度粘结单元CIE_INT作为示例。

      图  6  骨料表面的粘结单元(CIE_INT)

      Figure 6.  Cohesive element on the surface of an aggregate

      5) 输出插设粘结单元后的整体网格信息,作为ABAQUS可读取的*.inp输入文件供计算使用。

      对于图4中约有3.2万个节点、18.1万个四面体单元(C3D4)以及36.7万个单元面的原始模型,采用原有MATLAB程序[27-28]进行插设耗时5 h,而本文C++程序不到3 min。插设粘结单元之后,模型大约有53.1万个节点数、27.0万个粘结单元,其中骨料-砂浆界面上大约有2.6万个粘结单元。

    • 骨料和砂浆假设为线弹性材料,粘结单元采用线性受拉/受剪软化关系和名义应力平方起裂准则进行模拟。采用文献[22]中的材料参数(表2),并假设受剪与受拉断裂参数相同(偏于保守)。粘结单元的初始弹性刚度应足够高以模拟未开裂初始状态,但刚度过高会导致算法不收敛。经过多次尝试,将其取为6×106 MPa/mm。

      表 2  材料参数[22]

      Table 2.  Material parameters[22]

      弹性模量E/GPa泊松比密度ρ/(×10−6 kg/mm3)初始刚度kn /(MPa/mm)抗拉强度tn /MPa断裂能Gf /(N/mm)
      骨料700.22.5
      砂浆250.22.2
      CIE_CEM2.26×10640.06
      CIE_INT2.26×10620.03

      对生成的50 mm的细观混凝土立方体试件(图4)进行单轴受拉模拟。边界条件如图7所示,左端面全约束,右端面受均匀分布水平(x向)位移荷载控制,最终位移0.06 mm。使用ABAQUS/Explicit显式求解器进行准静态计算分析,设定荷载步时间为0.005 s。

      图  7  单轴受拉混凝土试件的边界条件

      Figure 7.  Boundary conditions under uniaxial tension

      图8给出了该模型采用两组不同断裂能时的平均应力-位移曲线,其中一组采用表2的参数(即参照组),另一组Gc0.06Gi0.01表示界面粘结单元的断裂能为0.01 N/m,其余参数同表2。图中也给出了典型单轴受拉试验数据[30]作为对比,该试验试件的尺寸是50 mm×50 mm×60 mm,与模型比较接近。可见Gc0.06Gi0.01得到的模拟结果与试验吻合较好,表明模型能够有效预测混凝土力学响应。但由于试验试件与本文数值模型在细观组分方面并不相同,该试验数据仅作参考。

      图  8  单轴受拉混凝土试件的宏观应力-位移曲线

      Figure 8.  Average stress - displacement curves of concrete specimen under uniaxial tension

    • 图9显示了采用参考材料参数(表2)模拟的试件正面和背面视角的开裂过程(对应于图8曲线的A点~F点)。其中红色的单元为损伤因子D≥0.9的粘结单元,在本文中假设为已完全损坏的离散裂缝。为了观察裂缝表面形态,图10显示了去掉这些损坏裂缝单元后的相应变形过程。为了方便观察,变形放大系数均设为200。

      图  9  单轴受拉开裂过程:粘结裂缝单元的受拉变形

      Figure 9.  Fracture process under uniaxial tension: the red cohesive elements with D≥0.9

      图  10  单轴受拉开裂过程

      Figure 10.  Crack propagation process under uniaxial tension

      图9图10所示,在早期加载阶段(A点之前),由于界面粘结单元的断裂参数比砂浆低,骨料-砂浆界面上出现大量微裂缝,但其损伤因子D仍较低,反映在裂缝单元的张开位移较小,以及荷载-位移曲线在加载位移较小时即显示出非线性特征;随着加载位移的增加直至试件达到抗拉强度(A点),微裂缝宽度缓慢发展,但裂缝并不显著;当试件进入软化段时(B点~D点),一些微裂缝的宽度迅速增加进而局部化形成主裂缝,并与砂浆中新生成的裂缝连通,而另外一些微裂缝则逐渐卸载闭合;如图中阶段(E点~F点)所示,随着位移进一步增加,连通的裂缝迅速张开变宽,形成一条贯穿主裂缝,试件拉裂。这和混凝土单向受拉实验[30]和其他数值模拟[22, 28]的断裂扩展过程一致。

    • 现有研究表明[26-27],粘结裂缝单元的材料参数,特别是抗拉强度与断裂能,对模拟结果影响较大。本文对砂浆中与骨料-砂浆界面上的粘结单元选取不同组别的抗拉强度与断裂能进行研究。保持断裂能不变(见表2),粘结单元取不同抗拉强度时,使用TcTi的编号规则,如Tc4Ti2表示砂浆和界面粘结单元的抗拉强度分别为4 MPa和2 MPa;保持表2中抗拉强度不变,取不同断裂能,使用GcGi的编号规则,如Gc0.06Gi0.03表示砂浆和界面粘结单元的断裂能分别为0.06 N/mm和0.03 N/mm。将采用表2材料参数的模型设置为参考组(编号为Tc4Ti4_REF与Gc0.06Gi0.03_REF)。

      取8组粘结单元材料参数进行模拟,对每一组,只变动TcTiGcGi参数中的一个,其余参数与参考组相同。得到如图11图12所示试件宏观平均应力-位移曲线,以及如图13图14所示的试件破坏形式。

      图  11  粘结裂缝单元的抗拉强度对宏观应力-位移曲线的影响

      Figure 11.  Average stress - displacement curves: effects of fracture energy of cohesive elements

      图  12  粘结裂缝单元的断裂能对宏观应力-位移曲线的影响

      Figure 12.  Average stress - displacement curves: effects of fracture energy of cohesive elements

      图  13  粘结裂缝单元的抗拉强度对试件破坏形式的影响

      Figure 13.  Effects of tensile strength of cohesive elements on the fracture surfaces

      图  14  粘结裂缝单元的断裂能对试件破坏形式的影响

      Figure 14.  Effects of fracture energy on crack surfaces

      图11可见,粘结单元的抗拉强度对试件应力-位移曲线峰值以及峰值前非线性段影响显著:随着粘结单元强度的增加,试件强度增加;当骨料-砂浆界面裂缝单元强度降低而砂浆裂缝单元强度增加时,试件强度亦增加,表明砂浆粘结单元的抗拉强度对试件承载力起控制作用。图11也表明:随着试件强度的增加,其脆性也逐渐增加。由图12可见,粘结单元的断裂能对试件应力-位移曲线软化段影响显著,随着断裂能的提高,试件的延性增加,但对试件强度以及曲线峰前段影响不大。

      图13可见,粘结单元抗拉强度对裂缝面形态有较大影响,该影响主要体现于砂浆、界面粘结单元的强度相对比值γT=Tc/Ti。当γT≠1时,断裂面的位置基本相同,γT主要影响断裂面的形态;当γT=4时(图13(c)设计组Tc4Ti1),断裂面上出现的大骨料数量最多,其次是γT=3(图13(b)设计组Tc6Ti2),然后是γT=2(图13参考组Tc4Ti2)。由于本文粘结单元以名义拉应力作为起裂准则,随着砂浆、界面粘结单元强度差别的增大,界面愈发成为薄弱环节。由于大骨料表面的界面单元面积较大,较容易成为裂缝起裂与扩展通道,使试件倾向于在大骨料附近形成裂缝面。另外,当γT=1时,如图13(a)图13(d)所示,试件断裂面较平滑且以小骨料数量居多;其中设计组Tc2Ti2的粘结单元强度较低,离加载端最近一层裂缝单元因应力集中而破坏。

      图14显示砂浆、界面粘结单元的断裂能相对比值γG=Gc/Gi对裂缝形态的影响。当γG≤3时(图14(a)、图14(b)、图14(d)和参考组),试件破坏时的位移和断裂面形态基本一致;随着断裂能的提高,裂缝面略趋于曲折;但当γG=6时(图14(c)设计组Gc0.06Gi0.01),砂浆裂缝单元的断裂能显著高于界面单元,界面成为薄弱环节,由于大骨料表面的界面单元面积较大,较容易成为裂缝扩展通道因而影响了裂缝扩展路径,使最终断裂面的位置与形态与其他组差异较大。

      进一步将Gc0.06Gi0.01与参考组Gc0.06Gi0.03_REF比较,分析断裂能、宏观应力-位移曲线以及开裂过程之间的关联。如图15图16所示,在这两个模型的应力-位移曲线的软化段中选取4个应力水平接近的阶段(A点~D点)进行分析。由图可见,由于Gc0.06Gi0.01的裂缝沿着图14(c)黑色箭头所指的大骨料(大骨料的界面作为薄弱环节)扩展,导致其裂缝面的形成比Gc0.06Gi0.03_REF更为曲折复杂,因此其应力-位移曲线的软化段由于断裂耗能的需要变得较为平缓,而没有呈现出Gi减小导致的脆性特征。由此可见,混凝土的力学响应能够反映其裂缝发展特征,二者既决定于断裂材料参数,也受到骨料大小、形状等细观结构因素的影响,显示出十分复杂的破坏机理。

      图  15  断裂能对开裂过程的影响

      Figure 15.  Effects of fracture energy on crack propagation processes

      图  16  断裂能对应力-位移曲线的影响

      Figure 16.  Effects of fracture energy on stress-displacement curves

    • 本文采用Python编程,利用ABAQUS脚本接口进行前处理二次开发,建立混凝土三维随机多面体骨料的细观有限元模型;通过高效插设离散粘结裂缝单元,成功模拟了单向受拉复杂三维多裂缝的起裂与扩展。主要结论如下:

      (1) 在早期加载阶段,由于界面粘结裂缝单元的断裂参数比砂浆低,骨料-砂浆界面上出现大量微裂缝,荷载-位移曲线较早显示出非线性特征;随着位移的增加直至试件达到抗拉强度,微裂缝宽度缓慢增加;进入软化段后,一些微裂缝宽度迅速增加,并与砂浆中新生成裂缝连通形成局部化的主裂缝,而另外一些微裂缝则逐渐卸载闭合。

      (2) 宏观应力-位移曲线主要受砂浆、界面粘性裂缝单元的抗拉强度和断裂能绝对数值的影响:砂浆粘结裂缝单元的抗拉强度对试件承载力起控制作用,试件承载力随着其强度的提高而增大;而断裂能对非线性软化段影响显著,断裂能的提高增加了试件的延性。

      (3) 裂缝面的位置和形态主要受砂浆、界面粘性裂缝单元的抗拉强度、断裂能相对比值的影响。比值较大时(γT>2或者γG>3),即界面粘结单元的力学性能远弱于砂浆粘结单元时,由于大骨料表面的界面单元面积较大,较容易成为裂缝起裂与扩展通道,使试件倾向于在大骨料附近形成裂缝面。

      (4) 混凝土的力学响应反映其裂缝发展特征,二者既决定于断裂参数,也受到细观结构的影响,本文建立的模型能够有效地描述混凝土复杂三维断裂过程。

目录

    /

    返回文章
    返回