基于ABAQUS的水工混凝土细观分析周期性边界条件自动施加

黄叶飞 徐 磊 刘 杰 荆帅召 周昌巧 金永苗

(1.中冶华天南京工程技术有限公司,南京 210019;2.河海大学 水利水电工程学院,南京 210098;3.重庆市水利电力建筑勘测设计研究院,重庆 400020)

随着数值分析方法的发展和计算机技术的进步,在混凝土结构分析中开展宏细观相结合的多尺度结构分析方法已成为当前的研究热点之一.

相对于基于给定的混凝土材料宏观本构模型并在单一宏观尺度上开展混凝土结构分析,在多尺度分析中混凝土材料的宏观力学行为是直接基于细观分析结果给出的,理论上具有更高的精度.在多尺度分析中为实现宏观尺度与细观尺度间的连接,通常需要为宏观尺度上的任一材料点赋予一个细观尺度上的代表性体积单元(RVE),继而通过降尺度技术将宏观应变转化为细观RVE的边界条件,并开展细观分析得到RVE的细观力学响应,在此基础上,通过升尺度技术(计算均匀化),将细观尺度上的计算结果转化为宏观尺度计算所需的本构行为[1-2].

如前所述,将宏观应变转化为细观RVE的边界条件是实现尺度连接的关键环节之一.目前,常用的细观RVE边界条件主要有均匀应变(Dirichlet)边界条件、均匀应力(Neuman)边界条件以及周期性边界条件.研究表明[3],随着RVE尺寸的增大,Dirichlet和Neuman边界条件下的结果会分别从上限和下限逼近周期性边界条件下的结果.对于固定尺寸的RVE,相比其他两种边界条件,周期性边界条件得到的材料宏观有效性能更加精准.因而,在混凝土结构的多尺度分析中,周期性边界条件是较为理想的选择.

现阶段,周期性边界条件的施加方法多采用手动人工操作的方法.对于水工混凝土而言,由于其骨料尺寸相对较大,相应的RVE尺寸也较大,这意味着其边界结点数量众多,当采用常规的手动方法对RVE施加周期性边界条件时,需要将相对边界上的节点一一匹配,逐点施加边界约束.在混凝土细观数值分析过程中,该方法存在操作繁琐,工作量巨大的弊端.而且在当前的研究工作中,对此问题仍然缺乏有效的解决方案.鉴于此,本文以周期性边界条件理论为基础,并结合有限元计算软件ABAQUS平台,提出了周期性边界条件自动施加的方法,给出了具体的实施流程,开发了基于ABAQUS的混凝土细观分析周期性边界条件自动施加程序.最后,通过算例分析验证了该方法的可行性.

1 周期性边界条件

在多尺度分析方法中,假设混凝土结构在宏观尺度上的任一积分点X都有相应的RVE与之关联,其中RVE的区域是V,边界是∂V,如图1所示.

图1 RVE选取示意图

在均匀化理论中,RVE上所施加的边界条件应满足宏观、细观变形能量相等的条件,即Hill-Mandel条件.Hill-Mandel条件如下:

根据虚功原理,外力在容许位移上做的功等于应力在容许应变上做的功.在不考虑体力的情况下,公式为:

故可将Hill-Mandel条件改写为:在细观分析中,RVE边界上的细观位移场可分为均匀应变场和位移波动场,即

将式(4)代入到(3)式,可得Hill-Mandel的简化公式:

对于一个二维的RVE几何结构,如图2所示,V1,V2,V3,V4为RVE的4个顶点,∂V+和∂V-为边界∂V的正负两个部分.在RVE的边界上施加周期性边界条件,要求其相对边界上的节点分布相同,同时满足位移连续和应力连续条件,即相对边上的节点位移相同,应力矢量大小相等,方向相反.周期性边界条件的特征可用下式表示:

图2 RVE的几何结构

可以发现将式(6)和式(7)代入到式(5)的Hill-Mandel定理中等式仍然成立,说明在RVE上施加周期性边界是满足Hill-Mandel条件的.

2 周期性边界的自动施加

对于RVE模型而言,其边界条件是由宏观应变场转换成的细观尺度上边界节点的位移来决定的[4].在定义RVE边界条件时,顶点V1在x和y方向的位移被约束,通过在其他边界节点上施加相应的位移约束来定义边界条件.对于图2所示的RVE结构,边界节点的具体位移值如下所示:

其中,u表示细观尺度上节点的位移,上标V1,V2,V3,V4表示RVE 4个顶点的节点编号,从左下角开始逆时针方向排序,而LR和TB分别表示RVE左右边界和上下边界上的节点编号;下标1,2表示节点的自由度;ε和γ表示宏观的应变分量;L x和L y分别表示RVE水平和竖直方向的尺寸.

根据式(4)和(6),可将周期性边界条件改写为:

在ABAQUS中,周期性边界条件的施加,一般是通过在RVE相对边界上的相应节点处使用约束方程来实现的(可参考ABAQUS帮助文档),约束方程的具体形式如下:

其中,A1,A2,A N是定义节点相对运动的系数;P,Q,R是节点的编号;i,j,k是节点自由度,如1,2,3分别表示x,y,z方向;^u是指定的位移约束,作为边界条件使用.

为了方便在ABAQUS定义周期性边界条件,引入了虚拟节点这个概念[5].将式(15)改写为:

其中表示虚拟节点Z在s方向的位移值.定义虚拟节点时,应注意的是虚拟节点Z的编号要足够大,不能与模型中的节点编号重复,且虚拟节点的坐标位置不能与模型的任何部分相连.

在ABAQUS的input文件中,使用命令*E-quation去定义这个约.

周期性边界条件的自动施加是建立在混凝土细观有限元模型的基础上.首先需要将混凝土细观结构进行网格划分,获取所有单元节点的编号以及坐标等信息.有限元网格划分时,相对边界上的网格密度应相同.

在获取有限元模型的数据信息后,本文基于周期性边界条件理论以及ABAQUS的input计算文件的格式要求,编制了相应程序用以实现周期性边界条件的自动施加,图3为相应程序的流程图.

图3 程序流程图

下面对程序中的约束方程进行详细说明.为了在有限元分析中应用周期性边界条件,程序引入了3个虚拟节点Set-dummy-XX、Set-dummy-XY、Set-dummy-YY,节点编号设为100000、100001、100002,坐标设为(-10,10)、(-10,0)、(10,0).然后使用下述命令行对边界节点施加约束方程.

最后,定义一个加载步,在该加载步中指定虚拟节点所对应的宏观应变分量.具体命令行如下:

3 算例验证

为验证2中所提到的周期性边界自动施加程序的正确性,取一个300 mm×300 mm的二维水工混凝土试件进行算例分析.

在本项研究中,认为混凝土在细观尺度上是由骨料、砂浆以及界面过渡区(ITZ)3部分构成的.考虑到水工混凝土的实际结构,骨料的几何形状被定为不规则多边形,界面过渡区则是用骨料和砂浆之间具有一定厚度的界面表示,砂浆被认为是连续均匀材料.混凝土试件中骨料含量为50%,骨料级配采用常用的混凝土三级配曲线,其中小石∶中石∶大石的比例为3∶3∶4.根据蒙特卡洛原理以及take-and-place方法,生成混凝土试件的随机骨料结构.详细步骤见文献[6].混凝土试件的随机骨料结构如图4所示.

图4 混凝土试件的几何模型

在材料的本构模型方面,对于骨料采用各向同性线弹性模型,对于砂浆和ITZ则采用CDP模型,而ITZ是混凝土材料中的薄弱环节,其相应材料参数均取为砂浆参数的0.75倍.

为了使混凝土细观材料参数更加符合实际,依据《混凝土结构设计规范(GB50010-2010)》[7]中给出的C20混凝土单轴受拉应力应变曲线,反演出砂浆和ITZ的力学参数,表1给出了混凝土细观组成的力学参数.

表1 混凝土细观组成的力学参数

根据混凝土的细观结构,借助于有限元软件平台ABAQUS,可实现混凝土细观结构的网格剖分,图5为混凝土试件的有限元模型,其中骨料位置用红色突出显示.网格密度设为骨料最小粒径的0.4倍,保证了相对边界上节点一一对应.

图5 混凝土试件的有限元模型

采用2中提出的自动施加周期性边界条件的方法,对混凝土有限元模型施加3种不同的约束,包括主轴方向的拉伸和平面内的剪切,这3种工况边界条件下的宏观应变分量见表2.

表2 RVE的边界条件

图6～7为工况1边界条件下的变形图和等效应力图,图8～9为工况2边界条件下的变形图和等效应力图,图10～11为工况3边界条件下的变形图和等效应力图.

图6 工况1变形图

图7 工况1等效应力图

图8 工况2变形图

图9 工况2等效应力图

图10 工况3变形图

图11 工况3等效应力图

由图6、图8和图10可以发现,在3种不同工况下,混凝土试件相对边界都表现出变形形态相同,位移差相等的特点,满足周期性边界位移连续的要求;由图7、图9和图11可以发现,在3种不同工况下,混凝土试件的应力分布并不均匀,但是相对边界的对应区域呈现应力一致的特点,满足周期性边界应力连续的要求.综上,可以证明本文所提出的方法满足周期性边界正确施加的两个基本特征:位移连续和应力连续,验证了本文方法的正确性.

与此同时,相比手动施加周期性边界条件方法的繁琐耗时,由于本文采用的是自动施加,在操作上更加便捷,施加速度也具有明显的优势.

4 结 语

针对水工混凝土细观分析中周期性边界条件施加存在操作繁琐,耗时过长等问题,本文以周期性边界条件理论为基础,结合有限元计算软件ABAQUS,提出了周期性边界条件自动施加方法,并通过算例对本文所提方法进行了验证.在3种不同工况下,使用本文方法对混凝土细观数值模型施加周期性边界条件.结果表明,本文所提方法满足周期性边界条件位移连续和应力连续的基本要求,且操作便捷,周期性边界条件施加速度显著提升,在水工混凝土细观分析中具有较好的可行性.