一种钢混组合结构粘结滑移的有限元模拟方法与流程

本发明属于数值模拟有限元分析技术领域，尤其是涉及一种基于python编程的钢混组合结构粘结滑移的有限元模拟方法。

背景技术：

型钢与混凝土之间的粘结强度是保证两者协同工作的重要前提，但目前钢混组合结构的有限元数值模拟，大多忽略了粘结滑移的影响，所得模拟结果与实际情况相比有所差距。因此，准确模拟钢混组合结构的粘结滑移可以对研究钢与混凝土组合结构的力学性能提供技术支撑，并对钢与混凝土组合结构的设计工作具有实际指导意义。

abaqus作为国际上先进的大型通用有限元分析软件之一，具有广泛模拟性能和非线性力学分析能力，但由于研究人员专业背景和发展方向不同，研究的专业问题纷纷不一，通用软件在模拟具体专业方面不能做到面面俱到。针对这些不足，大部分的通用软件都提供了二次开发功能，以弥补软件在不同研究领域专业性问题上的不足。

目前通过限元模拟粘结滑移的方法主要有cohesive单元、弹簧单元或者umat。其中cohesive单元适用于对两个部件之间的粘性连接进行模拟，要求粘结材料的尺寸和强度均小于其粘结界面，一般用于模拟复合材料层间的开裂性能；而umat则需要编辑子程序，采用fortran语法编写，过程繁杂，容易出错，对用户编程能力要求较高；采用弹簧单元，模拟过程相对简单，计算效率高，可减少计算运行内存和时间，但abaqus自带弹簧单元为线性弹簧，模拟钢与混凝土间的粘结滑移性能精确度不高。

经检索，中国专利申请号为201910078125.9的申请案公开了一种预应力装配式混凝土结构节点半精细数值模拟方法，该申请案的其步骤如下：步骤一、构建梁柱单元网格；步骤二、预制混凝土构件通过平面应力单元模拟，预制混凝土构件间的接触面通过接触单元模拟；步骤三、预制混凝土构件中普通钢筋及预应力筋的力学变化均通过桁架单元模拟；步骤四、对于模型中在受力作用下的钢筋耗能件以及其他耗能件，用桁架单元或者特殊耗能器单元模拟；步骤五、在模型中钢筋和混凝土之间的摩擦作用，通过设定非线性弹簧方法进行模拟。

又如，中国专利申请号为201710824945.9的申请案公开了一种预应力束与混凝土间粘结性能的有限元模拟方法，该申请案首先，给出了预应力束整体滑移面及有效工作面确定原则，建立了预应力束界面微单元受力平衡关系，提出了预应力束极限粘结力计算方法；其次，给出了预应力束极限粘结力的等效原则，并据此推导了单根预应力筋等效最大粘结应力计算公式，更新了bpe粘结滑移模型，明确了模拟粘结滑移的弹簧单元的属性，实现了单根预应力筋等效粘结滑移有限元模拟；在此基础上，基于单根预应力筋等效粘结力叠加获得预应力束粘结力、单根预应力筋平均滑移作为预应力束滑移的原则，实现了预应力束粘结滑移的有限元模拟。

但上述申请案分别是针对钢筋和钢筋束(桁架单元)与混凝土(实体单元)之间的粘结滑移进行研究，而对于型钢(壳单元)和混凝土(实体单元)之间的粘结滑移，弹簧单元两端节点数量更多更复杂，因此上述方法并不能直接用于型钢和混凝土之间的粘结滑移模拟。

技术实现要素：

1.要解决的问题

本发明的目的在于克服采用现有方法对钢与混凝土之间的粘结滑移性能进行模拟时，模拟精确性相对较低的不足，提供了一种钢混组合结构粘结滑移的有限元模拟方法。本发明采用非线性弹簧对钢混组合结构之间的粘结滑移进行模拟，从而有效提高了模拟精度，进而可以为钢混组合结构的精确模拟提供了基础。

2.技术方案

为了解决上述问题，本发明所采用的技术方案如下：

本发明的一种钢混组合结构粘结滑移的有限元模拟方法，包括以下步骤：

步骤一、通过有限元分析软件创建型钢高强混凝土模型，并在装配过程中对型钢与混凝土进行分区操作；

步骤二、在型钢与混凝土接触面间添加弹簧单元，通过前处理，提交作业，生成包含有限元模型信息的inp文件；

步骤三、通过abaqus自带的python语言接口，编写脚本，创建单元节点集合；

步骤四、在inp文件中添加弹簧刚度参数，完成对非线性弹簧的修改；

步骤五、重新在abaqus中创建作业，并输入修改后的inp文件，通过abaquscommand提交作业分析，即得到钢混组合结构粘结滑移的模拟结果。

更进一步的，步骤一中创建模型时型钢选取s4r曲壳单元，纵筋以及箍筋选取t3d2桁架单元，混凝土选取c3d8r实体单元。

更进一步的，步骤一中分区操作的具体过程为：首先创建水平与竖直方向基准轴，作为主轴；再通过原始基准轴平移得到水平、竖直两个方向初始基准面；对两个初始基准面进行偏移，即可得到沿构件全高的水平和垂直基准面，随后利用沿构件全高的水平和垂直基准面进行分区。

更进一步的，初始基准面的偏移量取值为0.04(40mm)。

更进一步的，步骤二中弹簧单元采用spring2类型弹簧，前处理参数设置如下：

(1)模拟钢混组合结构粘结滑移时设置多个分析步，包括初始分析步以及后续分析步，程序类型选取静力通用，采用增量迭代法进行计算；

(2)模拟采用单调静力轴心加载制度，创建多个荷载步，沿型钢滑移方向选取位移加载，各个分析步之间相互依次叠加；

(3)约束混凝土底部x、y、z三个方向的位移，构件加载端沿z轴方向施加50mm的位移量；

(4)在mesh模块中设置单元类型、形状及密度，取近似全局尺寸0.02，最大偏离因子0.1来为实例布置种子，然后再对部件进行网格划分。

更进一步的，步骤二中模拟钢混组合结构粘结滑移时共设置10个分析步，采用增量迭代法进行计算的最大增量步数为1000，增量步中初始设为0.1，最小设为10^-5。

更进一步的，步骤三的过程为：抽取模型中混凝土与型钢节点坐标数据，寻找重合单元节点编号，逐行读取，依次循环，按x、y、z三个方向进行排序，最终得到单元节点集合。

更进一步的，步骤二中对弹簧的方向进行定义，x轴方向为1,1；y轴方向为2,2；z轴方向为3,3；步骤四中定义弹簧3,3方向的刚度参数，该刚度参数由试验实测曲线测得。

更进一步的，修改后的inp文件包含的基本参数有：构件类型参数、几何参数、有限元参数、荷载参数、单元集合和弹簧刚度参数。

更进一步的，进行有限元模拟时，粘结滑移本构关系中的胶结粘结强度τ0、极限强度τu及残余粘结强度τr分别按下式计算：

τ0＝0.0673ft+0.153ρsv+0.109ca/ha-0.004la/ha+0.0012

τu＝0.0348ft+0.1733ρsv+0.0306ca/ha-0.0662la/ha+1.164

τr＝0.0197ft+0.0733ρsv+0.01345ca/ha-0.0142la/ha+0.605

式中，ft为混凝土抗拉强度值；ρsv为配箍率，ca/ha为型钢保护层厚度与型钢翼缘宽度的比值，la/ha为型钢埋深与翼缘宽度的比值；

当考虑型钢与混凝土之间沿锚固长度分布的粘结应力时：

τ(x)＝0.3837kεmaxe^kx

式中，τ(x)为沿型钢与混凝土之间沿锚固长度分布的粘结应力，x为计算粘结应力处与加载端处的垂直距离，k为粘结应变分布指数特征值，可按式k1＝(-0.4307la/ha+5.5634)×10^-3计算，εmax为距加载端为零的型钢等效最大应变值。

3.有益效果

相比于现有技术，本发明的有益效果为：

(1)本发明的一种钢混组合结构粘结滑移的有限元模拟方法，在abaqus中导入脚本，输出型钢与混凝土弹簧单元节点对应编号，通过修改inp文件中线性弹簧单元，并在型钢与混凝土对应弹簧单元节点编号之间输入弹簧刚度参数，成功添加非线性弹簧单元，最后在作业中提交新生成的inp文件，以此来模拟钢混组合结构粘结滑移性能，从而可以有效提高模拟精度。

(2)本发明的一种钢混组合结构粘结滑移的有限元模拟方法，实现了快速寻找钢与混凝土弹簧单元节点对应编号，摒弃了通过切片观察混凝土内部节点编号的方法，从而有效缩短了模拟过程后期添加非线性弹簧的工作量，为科研人员模拟钢混组合结构粘结滑移带来新思路和方法。

(3)本发明的一种钢混组合结构粘结滑移的有限元模拟方法，通过对型钢与混凝土间的粘结滑移本构关系中的参数计算公式进行优化，从而有利于进一步提高模拟的精度和准确性。

附图说明

图1为本发明的基于python编程模拟钢混组合结构粘结滑移的流程图。

图2为有限元模型分区示意图。

图3为添加及生成spring2类型弹簧的示意图。

图4为python编程代码和节点单元集合的示意图。

图5为混凝土内部视图切面和型钢应力云图。

具体实施方式

abaqus有限元软件中自带的弹簧单元为线性弹簧，模拟钢混组合结构中的粘结滑移与实际粘结滑移行为差距较大，目前绝大多数组合结构数值模拟中均忽略了粘结滑移的影响，仅采用嵌入方式传递应力，模拟结果较为粗略。而本发明在abaqus中导入脚本，输出型钢与混凝土弹簧单元节点对应编号，通过修改inp文件中线性弹簧单元，并在型钢与混凝土对应弹簧单元节点编号之间输入弹簧刚度参数，成功添加非线性弹簧单元，最后在作业中提交新生成的inp文件，以此模拟钢混组合结构粘结滑移性能，从而有效提高了模拟精度。

传统添加非线性弹簧单元的方法是对照实体模型，记录型钢与混凝土每个单元对应的重合节点编号，混凝土内部节点编号还需要视图剖切进行寻找，过程繁琐，工作量大，因此迫切需要寻求一种更加高效快捷的添加弹簧单元的方法。本发明实现了快速寻找钢与混凝土弹簧单元节点对应编号，摒弃通过切片观察混凝土内部节点编号的方法，有效缩短了模拟过程后期添加非线性弹簧的工作量，为科研人员模拟钢混组合结构粘结滑移带来新思路和方法。

如图1所示，本发明的钢混组合结构粘结滑移的有限元模拟方法，具体包括以下步骤：

步骤一、通过有限元分析软件创建型钢高强混凝土模型，并在装配过程中对构件进行分区操作；

试件尺寸：混凝土柱：300mm×300mm×600mm，内置型钢：150mm×75mm×7mm×4mm，纵筋为4根公称直径为12mm的三级钢，箍筋采用间距为100公称直径为8mm的一级钢。

混凝土采用塑性损伤模型，本构关系根据混凝土结构设计规范定义的单轴本构关系确定，钢筋和型钢采用理想弹塑性本构模型。

型钢选取s4r曲壳单元，适应性较好；纵筋以及箍筋选取t3d2桁架单元，拥有良好塑性，能有效模拟钢筋在混凝土中的受力和应变情况；混凝土选取c3d8r实体单元，适用于大变形分析和划分细网格。

在装配过程中对模型进行分区处理，方便后期输出对应节点编号，分区情况如图2所示。分区处理的具体操作为：首先创建水平与竖直方向基准轴，作为主轴，再通过原始基准轴平移得到水平、竖直两个方向基准面，对两个基准面进行偏移，偏移量取值0.04，即可得到沿构件全高的水平和垂直基准面，随后利用基准面进行分区。

步骤二、在相互作用中对型钢分区节点添加弹簧单元，如图3所示。

弹簧单元类型选用spring2弹簧，型钢与混凝土节点采用三个spring2弹簧进行连接，表示模型在三维空间上三个方向的粘结性能。对弹簧的方向进行定义，x轴方向为1,1；y轴方向为2,2；z轴方向为3,3，一次只能添加一个方向的弹簧，定义另一个弹簧方向时，单元编号需累积增加，防止后加弹簧覆盖前面已添加的弹簧。垂直于型钢推出方向的两个弹簧(1,1方向和2,2方向)用来模拟推出过程中四周混凝土对型钢的握裹力，设为大刚度弹簧，定义型钢和混凝土接触面的法向接触。

步骤三、创建相互作用、荷载和边界条件、进行网格划分、提交作业分析，生成包含有限元模型信息的inp文件，前处理参数设置如下：

(1)模拟钢混组合结构粘结滑移时需要设置多个分析步，包括初始分析步以及后续分析步，为了计算结果易于收敛，共设立10个分析步，程序类型选取静力通用，采用增量迭代法进行计算，最大增量步数为1000，增量步中初始设为0.1，最小设为10^-5。

(2)模拟采用单调静力轴心加载制度，创建多个荷载步，沿型钢滑移方向选取位移加载，各个分析步之间相互依次叠加。

(3)约束混凝土底部x、y、z三个方向位移，构件加载端沿z轴方向施加50mm的位移量。

(4)为了更加直观的反映试验现象，尽量使型钢与混凝土划分网格尺寸一致，网格数目越多，模拟结果越精确，此时需要在mesh模块中设置合理的单元类型、形状及密度，型钢选取s4r曲壳单元，纵筋以及箍筋选取t3d2桁架单元，混凝土选取c3d8r实体单元。单元形状皆为四边形，取近似全局尺寸0.02，最大偏离因子0.1来为实例布置种子，然后再对部件进行网格划分。

inp文件中包含的基本参数包括单元类型参数、几何参数、有限元参数、荷载参数，其中单元类型参数包括型钢和混凝土的单元类型；几何参数包括型钢和混凝土的截面尺寸和模型高度，钢筋骨架的外形尺寸；有限元参数包括材料的编号、节点坐标；弹簧单元的方向、个数、刚度；荷载参数包括荷载位置、大小及方向。

步骤四、通过abaqus自带的python语言接口，编写脚本，创建单元节点集合。

抽取模型中混凝土与型钢节点坐标数据，通过以下代码if(nx1＝＝nx2)and(ny1＝＝ny2)and(nz1＝＝nz2)实现，若两者坐标一致，则输出节点编号，依次循环，直至结束。

单元集合中实例名由部件名加编号组成，通过'*instance,name＝'+partname+'-1'实现。

节点集合按以下处理，allnodespair+＝[[node1,node2],]，最终依照allnodespair＝buildnodespair(dic1,dic2)输出所有单元集合，按x、y、z三个方向依次排序，如图4所示。

python编程将模型中实例部件输入数据、节点输入数据、弹簧单元数据按照预设的数据输入语法格式写入文件并按照预定路径保存。在cae菜单栏file下点击运行脚本(runscript)，弹出对话框，选择编写的脚本文件，导入脚本，自动输出型钢与混凝土弹簧单元节点对应编号。

步骤五、将脚本输出编号写入inp文件中，添加非线性弹簧单元，定义非线性弹簧刚度参数，重新生成inp文件；

型钢与混凝土间的粘结力分布由上升段、曲线上升段、曲线下降段和平稳段四部分组成，其粘结滑移本构关系中的胶结粘结强度τ0、极限强度τu及残余粘结强度τr可按下式计算：

τ0＝0.0673ft+0.153ρsv+0.109ca/ha-0.004la/ha+0.0012

τu＝0.0348ft+0.1733ρsv+0.0306ca/ha-0.0662la/ha+1.164

τr＝0.0197ft+0.0733ρsv+0.01345ca/ha-0.0142la/ha+0.605

式中，ft为混凝土抗拉强度值，单位为kn/mm²；ρsv为配箍率，ca/ha为型钢保护层厚度与型钢翼缘宽度的比值，la/ha为型钢埋深与翼缘宽度的比值。

当考虑型钢与混凝土之间沿锚固长度分布的粘结应力时：

τ(x)＝0.3837kεmaxe^kx

式中，τ(x)为沿型钢与混凝土之间沿锚固长度分布的粘结应力，单位为kn/mm²，x为计算粘结应力处与加载端处的垂直距离，单位为mm；k为粘结应变分布指数特征值，可按式k1＝(-0.4307la/ha+5.5634)×10^-3计算，εmax为距加载端为零的型钢等效最大应变值。

不同类型钢混组合结构中添加的弹簧刚度参数各不相同，具体参数亦可由试验测得的数据确定。最终生成inp文件中包含的基本参数有：单元类型参数、几何参数、有限元参数、荷载参数、单元集合和非线性弹簧刚度参数(弹簧3,3方向的刚度参数)。

步骤六、上述步骤操作结束，成功添加非线性弹簧，保存文件，重新在abaqus中创建作业，输入修改后的inp文件，通过abaquscommand提交作业分析，输出结果，即得到钢混组合结构粘结滑移变形云图，如图5所示。