南水北调东线穿黄河工程建设理论与实践
上QQ阅读APP看本书,新人免费读10天
设备和账号都新为新人

4.6 基于流固耦合分析的滩地埋管参数化模型

4.6.1 PYTHON语言对ABAQUS有限元模型参数化前处理

4.6.1.1 利用PYTHON实现ABAQUS参数化建模的方法

参数化设计指的用一组参数约束某几何图形的一组结构尺寸序列,当给参数序列赋予不同的数值,通过驱动程序就可以生成新的目标几何图形。参数化的优点是设计的参数除描述几何拓扑关系信息外,还能表达与处理几何元素闻的各种设计关系和约束关系。工程设计人员在进行产品设计时,就可以不需要再考虑其具体的设计细节,复杂的程序计算过程也由计算机强大的计算能力来完成,同时还可以通过改动其中的一些参数和约束关系就可以达到更新设计的工作。参数化设计的方法一般有编程参数化法,代数法,人工智能方法,特征参数化法以及语言描述法等。

通常建立有限元模型的时候,往往先借助其他三维制图软件(Pro/E、Solid Works、CATIA,UG等)生成几何模型,然后将几何模型导入到ABAQUS CAE环境中,根据需要进行模型修复、处理,然后进行网格划分建立有限元模型。由于软件的兼容性以及有限元模型与几何模型的差异性,在ABAQUS中进行几何模型修复也可能比较困难。于是有人连模型修复、网格划分都在第三方软件中进行(如上述制图软件和Hyper Mesh等网格划分软件)。在ABAQUS中建立分析步、建立接触对和施加约束、载荷条件后,在提交运算之前生成INP文件以便保存模型。在借助第三方软件建立几何模型和网格划分的情况下,过程离散,手工界面操作过多,难以进行自动化的参数化建模和分析,尤其是针对几何模型改变的参数化分析。而且过程的流程性、可追溯性不强,对模型的修改和处理、设置往往会遗漏或混乱。直接导入的原有几何模型也往往含有过多的、不必要的几何细节,给分析过程带来不必要的麻烦。

ABAQUS脚本接口是一个基于对象的程序库,内嵌脚本语言PYTHON,提供了一套应用程序编程接口(API)来操作ABAQUS/CAE实现建模/后处理等功能。接口编程采用PYTHON的语法编写脚本,但扩展了PYTHON脚本语言,额外提供了大约500个对象模型。对象模型之间关系复杂。整个对象模型被分成三类,其中session对象用来定义视图、远程队列,用户定义的视图等;mdb对象包括计算模型对象和作业对象;odb对象包括计算模型对象和计算结果数据。每一类对象下面又包括各类子对象,比如mdb对象下面的计算模型models对象又包括很多子对象(part,sketches,section,material,loads等)。

在对前处理过程进行二次开发的时候基本上都是和models对象打交道。它里面几乎包含了建模编程需要的所有对象类型,是主要考虑的对象类型。编写好的脚本文件提交给ABAQUS的执行过程为:ABAQUS自带的PYTHON解释器解释脚本语言;调用ABAQUS内核执行脚本命令;生成inp输入文件;提交给分析器进行有限元计算;生成结果文件。

运行脚本的途径有:

(1)在ABAQUS/CAE主窗口中的最下面有一个命令行接口(CLI),可以在里面输入脚本命令执行。既可以输入单行的命令执行;也可以把编写好的命令一整块复制进去执行;还可以利用命令execfile(mainguan.pyc)运行脚本文件。这种执行脚本的方式有一个优点,就是可以交互执行命令,在某一条命令出错的情况下重新输入正确的命令,而不需要跑到文本编辑器里面修改脚本文件,再回到ABAQUS重新运行。

(2)启动ABAQUS的时候会显示一个启动对话框,选择Run Script运行脚本文件。

(3)从ABAQUS/CAE主窗口选择菜单File->Run Script运行脚本文件。

(4)以命令行方式运行ABAQUS时可以在命令行中指定要运行的脚本文件,例如:

ABAQUScae script=mainguan.py

ABAQUS cae noGUI=mainguan.py

ABAQUS viewer noGUI=mainguan.py

execfile('mainguan.py')

利用ABAQUS程序的PYTHON脚本语言来操作ABAQUS前处理过程,实现自动建模、划分网格、指定材料属性,提交作业,后处理分析结果功能。避免了在使用ABAQUS/CAE前处理过程中的完全手工操作,节省了大量时间和精力,提高了工作效率。

到目前为止,使用大型有限元分析软件ABAQUS对地下结构进行流固耦合分析见诸很多文献,但基于ABAQU参数化有限元分析的地下结构还并不是很广泛。

采用PYTHON参数化语言对滩地埋管流固耦合结构进行建模分析是十分有效的,由于滩地埋管长达3943m,所处的地层、埋深、地下水位不同,不同位置荷载相差很大。滩地埋管流固耦合分析模型由多种材料、多种几何形状的部件组成,存在着管土接触面,在改变管径的同时需要改变接触面,要求保持二者始终接触。其几何结构和各部件之间连接关系非常复杂,有限元模型中包括成千个节点和单元,并且在每次建模操作过程中,面和关键点的序号也在不断发生变化,如若只是采用ABAQUS进行分析,需要大量重复建模,那将会花费设计人员的大量的时间和精力。对滩地埋管流固耦合结构来说,利用PYTHON语言建立参数化模型,只需要对滩地埋管流固耦合结构模型的参数化文件做简单的修改,便可分别不同断面埋管进行分析和截面优化设计,为设计人员减少了大量的重复性工作。这种参数化文件可以文本文件的格式保存,不会受ABAQUS软件版本和操作系统平台的限制。

4.6.1.2 滩地埋管模型主要参数

本章建立的滩地埋管流固耦合结构ABAQUS有限元参数化模型由其脚本程序PYTHONL语言描述,对不同地质、不同埋深、不同荷载的滩地埋管截面进行优化分析。为了避免重复性建模,必须实现模型的参数化。

建立滩地埋管结构参数化模型就要在所有的参数中确定能表征滩地埋管结构模型的主要参数作为设计变量。PYTHON语言编程需要提取的模型参数主要包括:管厚、混凝土管内径、接触面、混凝土管圆心标高、设计水位、碎石垫层厚、混凝土垫层厚、垫层底(开挖面)标高、地面标高、开挖坡度、各层土的底部标高、处理地基底部标高、碎石垫层材料性质和不同土层的材料性质(强度指标、密度、弹性模量、泊松比、渗透系数)、荷载、孔隙水压力边界。改变一个几何参数,其他参数要进行相应变化,保证不同材料接触面的随动变化。

4.6.2 滩地埋管模型几何参数化实现

为了划分出结构化的网格,先分割模型,分割后模型如图4-7所示。分割的方式采用两个数据点最短连线来确定分割面。这些点可以是线条相交生成的顶点(vertices),也可以是定义的数据点(datum)。

img

图4-7 分割后模型

4.6.2.1 滩地埋管模型网格控制与单元选择

为了能够在几何参数改变时重新生成的网格形状大致类似,在不同材料的边界上撒规定数目的种子。管壁厚1m,坡度1:2,地面标高42.5m、处理地基厚3m的网格如图4-8所示。

单元类型:混凝土为不透水材料,采用CPE8平面应变单元,土层、碎石垫层采用流固耦合单元CPE8P。

平面应变单元形状流固耦合形状设为四边形Quad,网格生成技术除了填土和边坡形状复杂,采用Free处理外,其余为Sweep。

4.6.2.2 定义接触面

定义接触面法向模型为硬接触(Hard),限制可能发生的穿透现象。定义接触面水平方向允许产生弹性滑移变形,采取罚(Penalty)刚度算法。

img

图4-8 管壁厚1m,坡度1:2,地面标高42.5m、处理地基厚3m网格

4.6.3 材料本构模型

回填土、地基土采用修正Drucker-Prage盖帽模型。修正Drucker-Prage盖帽模型在线性Drucker-Prage模型上增加了一个帽盖状屈服面,从而可以描述压缩导致的屈服,同时也能控制材料在剪切作用下的无限剪胀现象。

对管道结构的内力分析,均应按弹性体系计算,不考虑由非弹性变形所引起的塑性内力重分布。因此,混凝土埋管材料性质为弹性材料。

钻孔ZKM204(桩号3+400)所在剖面的几何和材料参数见表4-3。回填土修正Drucker-Prage盖帽模型参数见表4-4。

表4-3 钻孔ZKM204所在剖面几何和材料参数

img

续表

img

表4-4 回填土修正Drucker-Prage盖帽模型参数

img

4.6.4 荷载、荷载组合和分析工况

4.6.4.1 地表淤积土压力

滩地淤积按黄河河床每年平均淤积10cm计,运行期考虑黄河滩地30年、50年淤积两种情况。30年淤积3m。取壤土湿密度1.73t/m3

4.6.4.2 正常土压力

地面标高按42.5m、40.5m,静止土压力系数可按式(4-19)计算。

埋管上垂直土压力侧向土压力的作用分项系数,当其作用效应对管体结构不利时应采用1.1;有利时应采用0.9。

img

若墙厚填土为正常固结黏土,K0也可由式(4-20)计算:

img

4.6.4.3 内水压力

正常运行取进口前水位38.91m。

4.6.4.4 外水压力

正常地下水位35.40m。

汛期水位:黄河位山段大堤设防标准为11000m3/s,初步考虑黄河的淤积情况,2050年的相应洪水位为50.40m,滩地埋管按此洪水标准设计。

4.6.4.5 施工荷载

施工荷载为4kN/m2

按照《水工混凝土结构设计规范》(SL191—2008),有限元模型的荷载和分析步考虑了不同运行工况、荷载组合、荷载分项系数,见表4-5。

表4-5 有限元分析荷载、荷载组合及作用分项系数

img

注 表内括号内为荷载分项系数,在程序中定义埋管自重、填土自重的重力加速度为-9.81m/s2乘以分项系数后的组合值。斜坡荷载和地基自重的重力加速度为-9.81m/s2。淤积土压力3m为30年淤积土压力组合值,淤积土压力2m为20年内淤积土压力组合值。

4.6.5 边界条件及参数化实现

4.6.5.1 位移边界条件

模型左右边界施加水平方向的约束,底面边界施加竖直方向的约束。

4.6.5.2 孔隙水压力边界

(1)施工阶段降水孔隙水压力。施工阶段降水,地下水位应控制在建基面以下0.5m。

(2)设计外水孔隙水压力。标高35.40m,孔隙水压力边界为0;模型下边界孔隙水压力为设计外(水位35.40m)产生的孔隙水压力,为354 kPa。

(3)汛期孔隙水压力。

模型上表面标高为40.5m时,上表面孔隙水压力50.4-40.5=9.9m,即99kPa。

模型上表面标高为42.5m时,上表面孔隙水压力50.4-42.5=7.9m,即79kPa。

模型下边界孔隙水压力为洪水(水位50.40m)产生的孔隙水压力,为504kPa。

4.6.6 分析步设置

流固耦合有限元分析,首先应平衡地应力,建立初始的平衡状态,然后进行耦合分析。

4.6.6.1 地应力平衡

在运用ABAQUS软件分析岩土工程问题中,第一步是首先应建立初始平衡状态,包括水平方向和竖直方向的平衡,给岩土体设置初始有效应力、初始孔隙比、初始孔隙水压力和重力载荷等参数。

4.6.6.2 渗流-应力耦合分析

ABAQUS中流固耦合的分析有稳态分析和瞬态分析两种形式,稳态分析认为流体的流动速度,体积不随时间变化,瞬态分析可求解孔隙水压力、沉降随时间的变化过程。本次计算采用稳态分析。

在进行渗流-应力耦合分析时,解的形式有两种:总孔隙压力(Total pore pressure)和超孔隙压力(Excess pore pressure)。如定义重力载荷时用Gravity分布载荷类型(distribution Type),ABAQUS提供的是总孔隙压力解;如定义重力载荷时采用其他类型,ABAQUS提供的是超孔隙压力解。

采用总孔隙压力解的求解方式,用Material→Density定义材料的干密度,用Permeability→Specific重量参数定义流体的重度,用Gravity定义重力载荷。

4.6.7 PYTHON语言对ABAQUS有限元分析结果后处理

采用PYTHON开发的程序代码,访问ABAQUS输出数据库ODB文件,提取滩地埋管拱顶、设计水位线与埋管交接处截面、拱腰、侧拱与管底交界处、管底截面的每一分析步最后一帧的各剖面局部坐标下内力,并输出到文件maiguanforcea.CSV。

当埋管厚为0.95m时,桩号4+470断面各工况各截面内力及主应力分布如图4-9~图4-16所示。图中带下划线的数字为弯矩,不带下划线的数字为轴力。

img

图4-9 桩号4+470断面正常使用极限状态(Step 2)截面内力和主应力(管厚0.95m)

img

图4-10 桩号4+470断面正常使用极限状态(Step 3)截面内力和主应力(管厚0.95m)

img

图4-11 桩号4+470断面正常使用极限状态(Step 4)截面内力和主应力(管厚0.95m)

img

图4-12 桩号4+470断面运行期1(Step 5)截面内力和主应力(管厚0.95m)

img

图4-13 桩号4+470断面运行期2(Step 6)截面内力和主应力(管厚0.95m)

img

图4-14 桩号4+470断面运行期3(Step 7)截面内力和主应力(管厚0.95m)

img

图4-15 桩号4+470断面检修期(Step 8)截面内力和主应力(管厚0.95m)

img

图4-16 桩号4+470断面汛期(Step 9)截面内力和主应力(管厚0.95m)

4.6.8 MATLAB程序配筋计算

采用MATLAB开发的程序代码,根据PYTHON后处理程序提取的滩地埋管拱顶、设计水位线与埋管交接处截面、拱腰、侧拱与管底交界处、管底截面的内力,按照《水工混凝土结构设计规范》(SL 191-2008)进行配筋计算和裂缝宽度验算,满足承载能力极限状态、正常使用极限状态和构造要求。程序自动判断截面大偏压、小偏压、大偏拉、小偏拉以及受拉侧。

MATLAB配筋程序包含1个主程序Design Main.m和6各子程序。Design Main.m读取内力文件maiguanforce.CSV存放的各截面的每个分析步最后一帧的内力,子程序peijin.m判断构件受力状态(大小偏压、大小偏拉),并分别调用Design Main.m(偏拉构件配筋计算子程序),XiaoPian Ya.m(小偏压构件配筋计算子程序)、DaPian Ya.m(大偏压构件配筋计算子程序)。Crackwidth.m用来计算正常使用极限状态荷载组合情况埋管的裂缝宽度。Crackcheck.m根据peijin.m配筋结果,调用Crackwidth.m计算得到的裂缝宽度值,判断裂缝宽度是否满足要求(0.25mm),若不满足要求,则增加受拉侧配筋面积。

桩号4+470断面各工况各截面内力配筋见表4-6~表4-12,为管厚为0.95m时各工况内力。完建期、检修期、汛期均为偏压状态,各运行期均为偏拉状态。淤积土压力的增加,使得运行期管顶(截面1)、管侧(截面3)和管底(截面5)的轴向拉力减小,但弯矩增加。

表4-6 桩号4+470断面施工期(Step2)截一面内力(管厚0.95m)

img

表4-7 桩号4+470断面完建期(Step3)截面内力(管厚0.95m)

img

表4-8 桩号4+470断面运行期1(Step 5)截面内力(管厚0.95mm)

img

表4-9 桩号4+470断面运行期2(Step 6)截面内力(管厚0.95m)

img

表4-10 桩号4+470断面运行期3(Step 7)截面内力(管厚0.95m)

img

表4-11 桩号4+470断面检修期(Step 8)截面内力(管厚0.95m)

img

表4-12 桩号4+470断面汛期(Step 9)截面内力(管厚0.95m)

img

可以从表4-6~表4-12看出,桩号4+470断面运行期3(Step 7)各截面内力为最不利。配筋按此内力组合进行。运行期3截面的主应力分布如图4-17所示。可以看出,管顶内侧、管底内侧都产生很大主拉应力,是受力的最薄弱截面。

img

图4-17 桩号4+470断面运行期3(Step 7)截面主拉应力分布(管厚0.95m)