基于ABAQUS二次开发的可靠度计算程序的设计与应用

2021-09-08徐志纯

人民珠江 2021年9期

邹宏，伍剑，徐志纯，杨立

(成都理工大学地质灾害防治与地质环境保护国家重点实验室，四川成都 610059)

结构的安全性、适用性与耐久性对结构的可靠性提出了较高的要求，保证结构在规定的使用期内能够承受设计考虑的各种作用，满足设计要求的各项使用功能，这是工程结构可靠性的基本内容。为了对工程结构的可靠性进行定量分析，将结构在规定的条件下和规定时间内完成规定功能的概率定义为结构的可靠度[1]。在结构可靠度的计算方法当中，Monte Carlo方法是最基本的方法，当计算数据量足够时，也是相对比较准确的方法。如果通过真实试验进行Monte Carlo模拟计算可靠度，相关试验可能会因为条件不足或代价高昂而难以实施。可使用有限元模型替代试验模型，即采用计算机模拟或预测一个结构关于结构参数或设计变量的性能或响应[2]，使用Monte Carlo方法计算其可靠度。

在有限元模型的可靠度计算方面，近些年已开展了一些研究工作。徐军等[3]将可靠度计算与响应面有限元直接耦合，提出了可靠度响应面有限元法。蔡阳等[4]建立了重力式挡土墙可靠度分析的数学模型，并叙述了其Monte Carlo模拟的计算思路。伍国军等[5]编制了Matlab-Abaqus联合计算可靠度程序，对圆形隧洞开挖锚固承载力进行了可靠度分析。任斌斌等[6]使用Python语言获取离散随机场，与ABAQUS有限元模型相结合，计算了边坡的可靠度。蔡德咏等[7]通过ABAQUS二次开发实现了可靠性分析与有限元程序相结合，并将其应用于复合材料定向管优化设计方面。上述研究在有限元模型的可靠度计算方面取得了一定的进展，但仍然存在不足，首先，上述程序设计往往针对单一结构类型进行，不具有普适性，其他人员使用时需重新编写；其次，由于有限元脚本语言与编写程序所用语言不一致，相结合后可能会降低计算效率。

本文运用ABAQUS二次开发技术，采用Python语言将Monte Carlo模拟所需的重复模拟过程程序化，根据设计参数的概率分布自动生成随机参数表，程序读取随机参数表自动完成模型建立、作业提交、结果处理等工作。对于工程设计人员，该方法简化了重复建模过程的复杂程度，缩短了计算周期，提高了可靠度计算的效率和自动化程度；对于研究人员，模块化结构使其可便捷地替换可靠度计算对象，添加其他可靠度计算方法，具有很强的可移植性和可扩展性。

1 ABAQUS二次开发与Python脚本

ABAQUS在各个仿真领域都日益得到越来越广泛的应用，特别是其开放的二次开发功能。ABAQUS二次开发可分为子程序开发和用户图形界面程序开发两类[8]。子程序开发使用Fortran语言，主要用于材料本构关系、自定义单元等子程序的编写。用户图形界面程序开发基于Python语言，主要用于对原有ABAQUS/CAE界面的继承和扩展，开发专用的前后处理模块以及GUI工具等。

作为ABAQUS二次开发的工具语言，Python具有强大的功能。Python程序不仅可以实现ABAQUS/CAE中的所有前后处理操作，还可以实现许多超出ABAQUS基本功能的操作，同时也能减少很多重复性工作，大大提高计算效率。Python具有面向对象、适应性强、可扩展性强等特点[9]，加上其丰富且强大的扩展程序库，无论是数据处理还是科学计算都具有显著的优势。

2 可靠度计算程序的设计

ABAQUS模型可靠度计算程序主要包括参数生成、模型建立、结果处理、循环执行、可靠度计算5个模块，具体的程序流程见图1。

图1 有限元法可靠度计算流程

2.1 参数生成模块

根据随机变量的概率分布，产生足够多的样本值即随机数，这一过程称为对该随机变量的随机抽样[2]。参数生成的过程就是利用随机抽样原理，将重要的设计参数按其各自的分布类型生成许多组随机数，参数生成的主要步骤如下。

a)筛选并列举出重要的设计参数。设计参数主要分为两大类：一类是结构的基本属性，如几何尺寸、材料属性、接触条件等；另一类是施加在结构上的直接作用或间接作用，如重力荷载、温度作用、地震作用等。

b)定义随机数生成函数。Python语言的扩展程序库NumPy提供了常见的随机数生成函数，包括但不限于均匀分布、高斯分布、对数正态分布等常见分布类型，如果提供的类型不能满足要求，可根据相关公式自行编写。

c)通过提前对实验数据或统计资料的分析，确定设计参数的分布类型、均值、标准差等统计数据，根据参数的分布类型调用对应的随机数生成函数。

d)保存抽样结果为文本文件，方便后续调用。新建一个名称为input.txt的文本文件，将第三步生成的随机数按列写入此文件，每列之间用空格间隔。生成的文本文件的列数等于设计参数的个数，行数等于预先设定的随机数抽样次数。参数生成模块GrowthParameter.py的部分程序如下。

L =100

#定义参数

def GAUS(Meanvalue,Standarddeviation,size):

#定义分布函数

L = np.random.normal(loc=Meanvalue,

scale=Standarddeviation,size=size)

return(L)

LR=UNIF1(L,0.1*L,10)

#调用分布函数

file1 = open(′input.txt′,′w′)

#保存为文本文件

for i in range(len(LR)):

file1.write(str(LR[i])+′ ′)

file1.close()

2.2 模型建立模块

在结构的可靠度计算中，对结果影响最大的是模型的准确性，准确的有限元模型将会得到优良的可靠度计算结果。同时，模型的建立也是可靠度计算过程中最复杂的步骤。在有限元模型的建立和分析中，使用参数化有限元分析方法，利用参数来描述结构特征，通过参数来表征分析过程，从而实现可变结构参数的有限元分析[10]。

模型建立模块的主要工作是按照ABAQUS分析流程编写可变参数的有限元分析的命令流文件，提交命令流文件并生成结果数据库文件。可直接编辑生成命令流文件，但对使用者的能力要求较高，也可通过界面输入方式完整地完成一次有限元分析流程，在此流程中获取命令流文件。提交命令流文件的方法有2种：①生成并修改后缀名为inp的命令流文件，将其提交到ABAQUS Command中完成计算并得到结果数据库文件；②生成并修改后缀名为py的命令流文件，将其提交到ABAQUS GUI中完成计算并得到结果数据库文件。

上面2种方法都能获得相同的结果数据库文件，但计算速度存在差别。通过第2种方式提交计算时，是先将Python命令流运行后生成后缀名inp的输入文件，再将此文件提交计算，因此第2种方式的计算速度较慢。此外，第1种方法可同时批量计算多个有限元模型，相较于第2种更具有优势。

2.3 结果处理模块

可靠度计算模块需要从有限元模拟结果中获得结构的响应，这个响应是一个或多个具体的值，例如最大应力、最大应变等。ABAQUS输出的结果数据来自于整个模型或者模型的大部分区域，不但包含的数据信息量非常大，而且分为2种数据保存类型，包括以分析步划分的场输出数据和以帧划分的历史输出数据[11]。因此，为了提高结果提取效率，需编写专门的模块读取结果数据库文件，根据结构响应的数据类型及相关要求从数据库当中输出特定的数据。

ABAQUS的结果数据保存在工作目录下后缀名为odb的数据库文件当中，可通过Python脚本读取结果数据。例如，获取特定区域中的最大应力(应变)值，首先获取特定区域的节点(单元)号并创建节点(单元)集，按照节点(单元)集中的每一个节点(单元)的编号到数据库中读取其对应的应变(应力)，通过比较大小获得最大应力(应变)值，最后删除ABAQUS工作目录的所有文件，避免后续新的模型文件因权限问题而无法创建。

结果处理模块ExtractionData.py的部分程序如下。

Odb=openOdb(r′Job-1.odb′)

#打开ODB文件

inX=Odb.rootAssembly.instances[′PART-1-1′]

endNode = inX.nodeSets[′SET-1′]

#获取节点集的编号

Var= Odb.steps[′Step-1′].frames[-1].fieldOutputs[′U′]

#提取位移

nset_val =Var.getSubset(region=endNode).values

stress_data = map(lambda x:[

x.nodeLabel,x.data[1]],nset_val)

Res=stress_data[0][1]

del mdb.jobs[′Job-1′]

# 删除作业

2.4 循环执行模块

在使用Monte Carlo方法时，为批量完成模型建立、任务提交、结果输出等过程，需要设置专门的循环执行模块。该模块的主要任务是读取参数生成模块所生成的随机参数表，将其中的变量作为参数输入模型建立模块中，调用结果处理模块读取模型建立模块生成的结果数据库文件，保存结构的响应。

循环执行模块的具体流程如下：读取参数生成模块输出的input.txt文件，通过按行读取以实现循环，通过while函数判断是否结束循环;将读取的单行数据处理后，分别对模型建立模块生成的命令流文件当中的参数符号进行赋值，将修改过后的命令流文件提交到ABAQUS中，计算得到后缀名为odb的结果数据库文件，调用结果处理模块处理数据库文件，获取结构的响应，并将文件保存为output.txt文件。

循环执行模块LoopExecution.py的部分程序如下。

from ExtractionData import output

#导入结果处理模块

from ModelBuilding import create

#导入模型建立模块

path = r″input.txt″

#读取参数数值表

file = open(path,″r″)

mystr = file.readline()

#一次读取一行

L = float(mystr.split()[0])

#为参数赋值

create(L)

#创建模型

value=output()

#数据输出

file2= open(′output.txt′,′w′)

#结果保存

file2.write(value+′ ′)

2.5 可靠度计算模块

使用Monte Carlo模拟法计算可靠度的基本思路是：当重要设计参数x1,x2,…,xn(n为参数数量)的概率分布类型、均值、标准差等数据已知时，利用算法产生符合相应重要参数概率分布的随机数矩阵，矩阵的列数等于参数的数量，矩阵的行数等于每个参数生成的随机数的数量，每次从中抽取一行随机数X=(xn1,xn2,…,xnm)T(m为参数的抽样数)组成随机样本输入结构的参数化有限元模型，提交计算后提取得到一组随机抽样值S=(sn1,sn2,…,snk)T(k为响应的数量)，设极限状态函数：

Znk=sn0-snk

(1)

式中sn0——失效状态值。

当snk超过sn0，即Znk<0,则认为结构失效。由大数定律中的伯努利定理可知，当抽样次数足够大时，随机事件出现的频率近似于它的概率[12]。因此，将结构失效的次数n与总模拟次数m之比n/m近似为结构的失效概率pf，查表可得可靠度指标β。

3 可靠度计算程序的应用

3.1 模型介绍

某简支钢桁架桥二维模型[13]见图2，该钢桁架桥模型由23根杆件组成，所有水平杆都具有完全相同的弹性模量和横截面积，斜杆也是如此。该钢桁架桥的长度为24 m，高度为2 m，支座A为固定铰支座，支座B为滑动铰支座。荷载为6个相互独立的集中荷载，大小为50 kN。

图2 钢桁架桥模型

3.2 参数统计

钢桁架桥在施工和运营过程中，构件的几何误差、材料的随机性、荷载的不确定性都会对桥梁结构的可靠度产生影响，在设计中必须考虑这些因素以确保结构的安全性[14]。将水平杆弹性模量E1、斜杆弹性模量E2、水平杆截面面积A1、斜杆截面面积A2以及各集中力的大小P1～P6作为随机参数。各随机参数间相互独立，除集中力的大小P1～P6的概率分布为Ⅰ型极值分布外，其他随机参数均满足对数正态分布，采用Latin超立方抽样法[15]进行抽样。采用表1所示的参数统计数据。