4.6　基于流固耦合分析的灘地埋管參數(shù)化模型

4.6.1　PYTHON語言對(duì)ABAQUS有限元模型參數(shù)化前處理

4.6.1.1　利用PYTHON實(shí)現(xiàn)ABAQUS參數(shù)化建模的方法

參數(shù)化設(shè)計(jì)指的用一組參數(shù)約束某幾何圖形的一組結(jié)構(gòu)尺寸序列，當(dāng)給參數(shù)序列賦予不同的數(shù)值，通過驅(qū)動(dòng)程序就可以生成新的目標(biāo)幾何圖形。參數(shù)化的優(yōu)點(diǎn)是設(shè)計(jì)的參數(shù)除描述幾何拓?fù)潢P(guān)系信息外，還能表達(dá)與處理幾何元素聞的各種設(shè)計(jì)關(guān)系和約束關(guān)系。工程設(shè)計(jì)人員在進(jìn)行產(chǎn)品設(shè)計(jì)時(shí)，就可以不需要再考慮其具體的設(shè)計(jì)細(xì)節(jié)，復(fù)雜的程序計(jì)算過程也由計(jì)算機(jī)強(qiáng)大的計(jì)算能力來完成，同時(shí)還可以通過改動(dòng)其中的一些參數(shù)和約束關(guān)系就可以達(dá)到更新設(shè)計(jì)的工作。參數(shù)化設(shè)計(jì)的方法一般有編程參數(shù)化法，代數(shù)法，人工智能方法，特征參數(shù)化法以及語言描述法等。

通常建立有限元模型的時(shí)候，往往先借助其他三維制圖軟件(Pro/E、Solid Works、CATIA，UG等)生成幾何模型，然后將幾何模型導(dǎo)入到ABAQUS CAE環(huán)境中，根據(jù)需要進(jìn)行模型修復(fù)、處理，然后進(jìn)行網(wǎng)格劃分建立有限元模型。由于軟件的兼容性以及有限元模型與幾何模型的差異性，在ABAQUS中進(jìn)行幾何模型修復(fù)也可能比較困難。于是有人連模型修復(fù)、網(wǎng)格劃分都在第三方軟件中進(jìn)行(如上述制圖軟件和Hyper Mesh等網(wǎng)格劃分軟件)。在ABAQUS中建立分析步、建立接觸對(duì)和施加約束、載荷條件后，在提交運(yùn)算之前生成INP文件以便保存模型。在借助第三方軟件建立幾何模型和網(wǎng)格劃分的情況下，過程離散，手工界面操作過多，難以進(jìn)行自動(dòng)化的參數(shù)化建模和分析，尤其是針對(duì)幾何模型改變的參數(shù)化分析。而且過程的流程性、可追溯性不強(qiáng)，對(duì)模型的修改和處理、設(shè)置往往會(huì)遺漏或混亂。直接導(dǎo)入的原有幾何模型也往往含有過多的、不必要的幾何細(xì)節(jié)，給分析過程帶來不必要的麻煩。

ABAQUS腳本接口是一個(gè)基于對(duì)象的程序庫(kù)，內(nèi)嵌腳本語言PYTHON，提供了一套應(yīng)用程序編程接口(API)來操作ABAQUS/CAE實(shí)現(xiàn)建模/后處理等功能。接口編程采用PYTHON的語法編寫腳本，但擴(kuò)展了PYTHON腳本語言，額外提供了大約500個(gè)對(duì)象模型。對(duì)象模型之間關(guān)系復(fù)雜。整個(gè)對(duì)象模型被分成三類，其中session對(duì)象用來定義視圖、遠(yuǎn)程隊(duì)列，用戶定義的視圖等;mdb對(duì)象包括計(jì)算模型對(duì)象和作業(yè)對(duì)象;odb對(duì)象包括計(jì)算模型對(duì)象和計(jì)算結(jié)果數(shù)據(jù)。每一類對(duì)象下面又包括各類子對(duì)象，比如mdb對(duì)象下面的計(jì)算模型models對(duì)象又包括很多子對(duì)象(part，sketches，section，material，loads等)。

在對(duì)前處理過程進(jìn)行二次開發(fā)的時(shí)候基本上都是和models對(duì)象打交道。它里面幾乎包含了建模編程需要的所有對(duì)象類型，是主要考慮的對(duì)象類型。編寫好的腳本文件提交給ABAQUS的執(zhí)行過程為:ABAQUS自帶的PYTHON解釋器解釋腳本語言;調(diào)用ABAQUS內(nèi)核執(zhí)行腳本命令;生成inp輸入文件;提交給分析器進(jìn)行有限元計(jì)算;生成結(jié)果文件。

運(yùn)行腳本的途徑有:

(1)在ABAQUS/CAE主窗口中的最下面有一個(gè)命令行接口(CLI)，可以在里面輸入腳本命令執(zhí)行。既可以輸入單行的命令執(zhí)行;也可以把編寫好的命令一整塊復(fù)制進(jìn)去執(zhí)行;還可以利用命令execfile(mainguan.pyc)運(yùn)行腳本文件。這種執(zhí)行腳本的方式有一個(gè)優(yōu)點(diǎn)，就是可以交互執(zhí)行命令，在某一條命令出錯(cuò)的情況下重新輸入正確的命令，而不需要跑到文本編輯器里面修改腳本文件，再回到ABAQUS重新運(yùn)行。

(2)啟動(dòng)ABAQUS的時(shí)候會(huì)顯示一個(gè)啟動(dòng)對(duì)話框，選擇Run Script運(yùn)行腳本文件。

(3)從ABAQUS/CAE主窗口選擇菜單File->Run Script運(yùn)行腳本文件。

(4)以命令行方式運(yùn)行ABAQUS時(shí)可以在命令行中指定要運(yùn)行的腳本文件，例如:

ABAQUScae script=mainguan.py

ABAQUS cae noGUI=mainguan.py

ABAQUS viewer noGUI=mainguan.py

execfile('mainguan.py')

利用ABAQUS程序的PYTHON腳本語言來操作ABAQUS前處理過程，實(shí)現(xiàn)自動(dòng)建模、劃分網(wǎng)格、指定材料屬性，提交作業(yè)，后處理分析結(jié)果功能。避免了在使用ABAQUS/CAE前處理過程中的完全手工操作，節(jié)省了大量時(shí)間和精力，提高了工作效率。

到目前為止，使用大型有限元分析軟件ABAQUS對(duì)地下結(jié)構(gòu)進(jìn)行流固耦合分析見諸很多文獻(xiàn)，但基于ABAQU參數(shù)化有限元分析的地下結(jié)構(gòu)還并不是很廣泛。

采用PYTHON參數(shù)化語言對(duì)灘地埋管流固耦合結(jié)構(gòu)進(jìn)行建模分析是十分有效的，由于灘地埋管長(zhǎng)達(dá)3943m，所處的地層、埋深、地下水位不同，不同位置荷載相差很大。灘地埋管流固耦合分析模型由多種材料、多種幾何形狀的部件組成，存在著管土接觸面，在改變管徑的同時(shí)需要改變接觸面，要求保持二者始終接觸。其幾何結(jié)構(gòu)和各部件之間連接關(guān)系非常復(fù)雜，有限元模型中包括成千個(gè)節(jié)點(diǎn)和單元，并且在每次建模操作過程中，面和關(guān)鍵點(diǎn)的序號(hào)也在不斷發(fā)生變化，如若只是采用ABAQUS進(jìn)行分析，需要大量重復(fù)建模，那將會(huì)花費(fèi)設(shè)計(jì)人員的大量的時(shí)間和精力。對(duì)灘地埋管流固耦合結(jié)構(gòu)來說，利用PYTHON語言建立參數(shù)化模型，只需要對(duì)灘地埋管流固耦合結(jié)構(gòu)模型的參數(shù)化文件做簡(jiǎn)單的修改，便可分別不同斷面埋管進(jìn)行分析和截面優(yōu)化設(shè)計(jì)，為設(shè)計(jì)人員減少了大量的重復(fù)性工作。這種參數(shù)化文件可以文本文件的格式保存，不會(huì)受ABAQUS軟件版本和操作系統(tǒng)平臺(tái)的限制。

4.6.1.2　灘地埋管模型主要參數(shù)

本章建立的灘地埋管流固耦合結(jié)構(gòu)ABAQUS有限元參數(shù)化模型由其腳本程序PYTHONL語言描述，對(duì)不同地質(zhì)、不同埋深、不同荷載的灘地埋管截面進(jìn)行優(yōu)化分析。為了避免重復(fù)性建模，必須實(shí)現(xiàn)模型的參數(shù)化。

建立灘地埋管結(jié)構(gòu)參數(shù)化模型就要在所有的參數(shù)中確定能表征灘地埋管結(jié)構(gòu)模型的主要參數(shù)作為設(shè)計(jì)變量。PYTHON語言編程需要提取的模型參數(shù)主要包括:管厚、混凝土管內(nèi)徑、接觸面、混凝土管圓心標(biāo)高、設(shè)計(jì)水位、碎石墊層厚、混凝土墊層厚、墊層底(開挖面)標(biāo)高、地面標(biāo)高、開挖坡度、各層土的底部標(biāo)高、處理地基底部標(biāo)高、碎石墊層材料性質(zhì)和不同土層的材料性質(zhì)(強(qiáng)度指標(biāo)、密度、彈性模量、泊松比、滲透系數(shù))、荷載、孔隙水壓力邊界。改變一個(gè)幾何參數(shù)，其他參數(shù)要進(jìn)行相應(yīng)變化，保證不同材料接觸面的隨動(dòng)變化。

4.6.2　灘地埋管模型幾何參數(shù)化實(shí)現(xiàn)

為了劃分出結(jié)構(gòu)化的網(wǎng)格，先分割模型，分割后模型如圖4-7所示。分割的方式采用兩個(gè)數(shù)據(jù)點(diǎn)最短連線來確定分割面。這些點(diǎn)可以是線條相交生成的頂點(diǎn)(vertices)，也可以是定義的數(shù)據(jù)點(diǎn)(datum)。

4.6.2.1　灘地埋管模型網(wǎng)格控制與單元選擇

為了能夠在幾何參數(shù)改變時(shí)重新生成的網(wǎng)格形狀大致類似，在不同材料的邊界上撒規(guī)定數(shù)目的種子。管壁厚1m，坡度1：2，地面標(biāo)高42.5m、處理地基厚3m的網(wǎng)格如圖4-8所示。

單元類型:混凝土為不透水材料，采用CPE8平面應(yīng)變單元，土層、碎石墊層采用流固耦合單元CPE8P。

平面應(yīng)變單元形狀流固耦合形狀設(shè)為四邊形Quad，網(wǎng)格生成技術(shù)除了填土和邊坡形狀復(fù)雜，采用Free處理外，其余為Sweep。

4.6.2.2　定義接觸面

定義接觸面法向模型為硬接觸(Hard)，限制可能發(fā)生的穿透現(xiàn)象。定義接觸面水平方向允許產(chǎn)生彈性滑移變形，采取罰(Penalty)剛度算法。

4.6.3　材料本構(gòu)模型

回填土、地基土采用修正Drucker-Prage蓋帽模型。修正Drucker-Prage蓋帽模型在線性Drucker-Prage模型上增加了一個(gè)帽蓋狀屈服面，從而可以描述壓縮導(dǎo)致的屈服，同時(shí)也能控制材料在剪切作用下的無限剪脹現(xiàn)象。

對(duì)管道結(jié)構(gòu)的內(nèi)力分析，均應(yīng)按彈性體系計(jì)算，不考慮由非彈性變形所引起的塑性內(nèi)力重分布。因此，混凝土埋管材料性質(zhì)為彈性材料。

鉆孔ZKM204(樁號(hào)3+400)所在剖面的幾何和材料參數(shù)見表4-3。回填土修正Drucker-Prage蓋帽模型參數(shù)見表4-4。

4.6.4　荷載、荷載組合和分析工況

4.6.4.1　地表淤積土壓力

灘地淤積按黃河河床每年平均淤積10cm計(jì)，運(yùn)行期考慮黃河灘地30年、50年淤積兩種情況。30年淤積3m。取壤土濕密度1.73t/m3。

4.6.4.2　正常土壓力

地面標(biāo)高按42.5m、40.5m，靜止土壓力系數(shù)可按式(4-19)計(jì)算。

埋管上垂直土壓力側(cè)向土壓力的作用分項(xiàng)系數(shù)，當(dāng)其作用效應(yīng)對(duì)管體結(jié)構(gòu)不利時(shí)應(yīng)采用1.1;有利時(shí)應(yīng)采用0.9。

若墻厚填土為正常固結(jié)黏土，K0也可由式(4-20)計(jì)算:

4.6.4.3　內(nèi)水壓力

正常運(yùn)行取進(jìn)口前水位38.91m。

4.6.4.4　外水壓力

正常地下水位35.40m。

汛期水位:黃河位山段大堤設(shè)防標(biāo)準(zhǔn)為11000m3/s，初步考慮黃河的淤積情況，2050年的相應(yīng)洪水位為50.40m，灘地埋管按此洪水標(biāo)準(zhǔn)設(shè)計(jì)。

4.6.4.5　施工荷載

施工荷載為4kN/m2。

按照《水工混凝土結(jié)構(gòu)設(shè)計(jì)規(guī)范》(SL191—2008)，有限元模型的荷載和分析步考慮了不同運(yùn)行工況、荷載組合、荷載分項(xiàng)系數(shù)，見表4-5。

4.6.5　邊界條件及參數(shù)化實(shí)現(xiàn)

4.6.5.1　位移邊界條件

模型左右邊界施加水平方向的約束，底面邊界施加豎直方向的約束。

4.6.5.2　孔隙水壓力邊界

(1)施工階段降水孔隙水壓力。施工階段降水，地下水位應(yīng)控制在建基面以下0.5m。

(2)設(shè)計(jì)外水孔隙水壓力。標(biāo)高35.40m，孔隙水壓力邊界為0;模型下邊界孔隙水壓力為設(shè)計(jì)外(水位35.40m)產(chǎn)生的孔隙水壓力，為354 kPa。

(3)汛期孔隙水壓力。

模型上表面標(biāo)高為40.5m時(shí)，上表面孔隙水壓力50.4-40.5=9.9m，即99kPa。

模型上表面標(biāo)高為42.5m時(shí)，上表面孔隙水壓力50.4-42.5=7.9m，即79kPa。

模型下邊界孔隙水壓力為洪水(水位50.40m)產(chǎn)生的孔隙水壓力，為504kPa。

4.6.6　分析步設(shè)置

流固耦合有限元分析，首先應(yīng)平衡地應(yīng)力，建立初始的平衡狀態(tài)，然后進(jìn)行耦合分析。

4.6.6.1　地應(yīng)力平衡

在運(yùn)用ABAQUS軟件分析巖土工程問題中，第一步是首先應(yīng)建立初始平衡狀態(tài)，包括水平方向和豎直方向的平衡，給巖土體設(shè)置初始有效應(yīng)力、初始孔隙比、初始孔隙水壓力和重力載荷等參數(shù)。

4.6.6.2　滲流-應(yīng)力耦合分析

ABAQUS中流固耦合的分析有穩(wěn)態(tài)分析和瞬態(tài)分析兩種形式，穩(wěn)態(tài)分析認(rèn)為流體的流動(dòng)速度，體積不隨時(shí)間變化，瞬態(tài)分析可求解孔隙水壓力、沉降隨時(shí)間的變化過程。本次計(jì)算采用穩(wěn)態(tài)分析。

在進(jìn)行滲流-應(yīng)力耦合分析時(shí)，解的形式有兩種:總孔隙壓力(Total pore pressure)和超孔隙壓力(Excess pore pressure)。如定義重力載荷時(shí)用Gravity分布載荷類型(distribution Type)，ABAQUS提供的是總孔隙壓力解;如定義重力載荷時(shí)采用其他類型，ABAQUS提供的是超孔隙壓力解。

采用總孔隙壓力解的求解方式，用Material→Density定義材料的干密度，用Permeability→Specific重量參數(shù)定義流體的重度，用Gravity定義重力載荷。

4.6.7　PYTHON語言對(duì)ABAQUS有限元分析結(jié)果后處理

采用PYTHON開發(fā)的程序代碼，訪問ABAQUS輸出數(shù)據(jù)庫(kù)ODB文件，提取灘地埋管拱頂、設(shè)計(jì)水位線與埋管交接處截面、拱腰、側(cè)拱與管底交界處、管底截面的每一分析步最后一幀的各剖面局部坐標(biāo)下內(nèi)力，并輸出到文件maiguanforcea.CSV。

當(dāng)埋管厚為0.95m時(shí)，樁號(hào)4+470斷面各工況各截面內(nèi)力及主應(yīng)力分布如圖4-9～圖4-16所示。圖中帶下劃線的數(shù)字為彎矩，不帶下劃線的數(shù)字為軸力。

4.6.8　MATLAB程序配筋計(jì)算

采用MATLAB開發(fā)的程序代碼，根據(jù)PYTHON后處理程序提取的灘地埋管拱頂、設(shè)計(jì)水位線與埋管交接處截面、拱腰、側(cè)拱與管底交界處、管底截面的內(nèi)力，按照《水工混凝土結(jié)構(gòu)設(shè)計(jì)規(guī)范》(SL 191-2008)進(jìn)行配筋計(jì)算和裂縫寬度驗(yàn)算，滿足承載能力極限狀態(tài)、正常使用極限狀態(tài)和構(gòu)造要求。程序自動(dòng)判斷截面大偏壓、小偏壓、大偏拉、小偏拉以及受拉側(cè)。

MATLAB配筋程序包含1個(gè)主程序Design Main.m和6各子程序。Design Main.m讀取內(nèi)力文件maiguanforce.CSV存放的各截面的每個(gè)分析步最后一幀的內(nèi)力，子程序peijin.m判斷構(gòu)件受力狀態(tài)(大小偏壓、大小偏拉)，并分別調(diào)用Design Main.m(偏拉構(gòu)件配筋計(jì)算子程序)，XiaoPian Ya.m(小偏壓構(gòu)件配筋計(jì)算子程序)、DaPian Ya.m(大偏壓構(gòu)件配筋計(jì)算子程序)。Crackwidth.m用來計(jì)算正常使用極限狀態(tài)荷載組合情況埋管的裂縫寬度。Crackcheck.m根據(jù)peijin.m配筋結(jié)果，調(diào)用Crackwidth.m計(jì)算得到的裂縫寬度值，判斷裂縫寬度是否滿足要求(0.25mm)，若不滿足要求，則增加受拉側(cè)配筋面積。

樁號(hào)4+470斷面各工況各截面內(nèi)力配筋見表4-6～表4-12，為管厚為0.95m時(shí)各工況內(nèi)力。完建期、檢修期、汛期均為偏壓狀態(tài)，各運(yùn)行期均為偏拉狀態(tài)。淤積土壓力的增加，使得運(yùn)行期管頂(截面1)、管側(cè)(截面3)和管底(截面5)的軸向拉力減小，但彎矩增加。

可以從表4-6～表4-12看出，樁號(hào)4+470斷面運(yùn)行期3(Step 7)各截面內(nèi)力為最不利。配筋按此內(nèi)力組合進(jìn)行。運(yùn)行期3截面的主應(yīng)力分布如圖4-17所示?？梢钥闯?，管頂內(nèi)側(cè)、管底內(nèi)側(cè)都產(chǎn)生很大主拉應(yīng)力，是受力的最薄弱截面。