鄧思遠， 楊其新，*， 蔣雅君， 陳　浩

(1． 西南交通大學土木工程學院交通隧道工程教育部重點實驗室， 四川 成都　610031；

2． 江蘇省交通規(guī)劃設計院股份有限公司軌道所， 江蘇 南京　210014)

?

FLAC3D流固耦合滲流模型探討

鄧思遠1， 楊其新1，*， 蔣雅君1， 陳浩2

(1． 西南交通大學土木工程學院交通隧道工程教育部重點實驗室， 四川 成都610031；

2． 江蘇省交通規(guī)劃設計院股份有限公司軌道所， 江蘇 南京210014)

摘要：FLAC3D有限元軟件內(nèi)置4種滲流模型，直接決定計算結(jié)果正確與否。滲流場在三維空間中分布復雜，難以根據(jù)計算公式定量選用何種模型，不同文獻滲流模型選擇標準也不相同。為了找到簡單、快速、合理的選擇方法，有必要對4種滲流模型進行受力分析。以一100 m×100 m×10 m各項同性彈性立方體為例，基于4種滲流模型設計4種工況。通過對比和分析各工況下總應力和孔隙水壓分布情況，主要結(jié)論如下： 1)模型A、模型C外荷載均由土顆粒骨架承擔，其余2種模型流體也參與受力； 2)流體分擔外荷載的比例與剛度系數(shù)和時間有關(guān)； 3)FLAC3D流固耦合過程正是通過調(diào)整剛度比、打開和關(guān)閉力學-流體進程來實現(xiàn)的。目前我國流固耦合計算基本基于FLAC3D有限元軟件，文中結(jié)論適用于基礎、隧道、基坑等工程，可供相關(guān)從業(yè)人士借鑒參考。

關(guān)鍵詞：流固耦合； 滲流模式； FLAC3D； 巖土工程

0引言

巖土工程施工階段引排水改變地下水分布，令土體有效應力重新分布，有效應力場反過來又影響地下水運動，這種應力場和滲流場相互作用的現(xiàn)象被稱為流固耦合[1]。20世紀中期，比奧按彈性力學理論推導出方程組，由于求解困難，當時未引起重視。近些年來，相關(guān)研究成果顯著。一方面，專家學者研制出能模擬流固耦合的模型試驗系統(tǒng)。比如： 文獻[2]通過模型試驗得到海底隧道施工中圍巖力學行為和滲流場變化規(guī)律；文獻[3]對比了模型試驗和數(shù)值計算之間的差異。另一方面，基于流固耦合的數(shù)值計算如火如荼。比如： 文獻[4]將流固耦合用于雙圓異性斷面盾構(gòu)隧道施工分析；文獻[5]研究了基坑開挖時支護變形和孔隙水壓的分布。目前，研究正向著裂隙介質(zhì)非均勻流、非飽和滲流以及多場耦合方向發(fā)展[6-7]。

雖然成果豐碩，但在流固耦合基礎應用方面還存在一些問題。我國流固耦合數(shù)值模擬主要基于FLAC3D平臺，軟件內(nèi)置4種滲流模型?；?、隧道開挖是一個復雜、循環(huán)往復的過程，既存在土體開挖卸荷，也存在引排水改變滲流場，難以根據(jù)用戶手冊選擇適合的模型，很多文獻選擇模型的標準也不盡相同。文獻[8]模擬地鐵開挖時，先計算應力場，再進行流固耦合或滲流計算；文獻[9]模擬基坑開挖時，先等滲流場穩(wěn)定后，再進行應力計算。上述處理方法均有一定道理，到目前尚未有文獻針對4種模型受力特點進行詳實論證。為了解4種模型的受力本質(zhì)，確保計算結(jié)果能指導工程，本文將以一各項同性彈性立方體為例，對各模型的受力特點進行分析。

1流固耦合概述

1.1相關(guān)理論

如圖1所示，取一微元體，建立基本微分平衡方程：

(1)

圖1　微元體受力示意圖

1.2計算原理

FLAC3D計算過程簡化如圖2所示，左邊表示應力場計算的過程，右邊表示考慮滲流場后的新增步驟。表1為FLAC3D用戶手冊滲流模型選擇標準，主要是根據(jù)特征時間和擾動類型確定的。

圖2　FLAC3D建模流程

軟件根據(jù)運動方程、平衡方程、本構(gòu)方程、相容方程以及邊界條件進行求解[1]。

1)運動方程是由達西定律v=k·i推導而來，滲流速度

(2)

2)平衡方程是在小變形情況下，單位時間內(nèi)微元體含水量變化值等于流入量與流出量之差：

(3)

式中:qυ為微元體流體流入量；ζ為單元流體體積變化值。

3)本文只針對飽和土體滲流模型，其本構(gòu)方程為

(4)

式中：M為比奧模量；α為比奧系數(shù)；ε為應力場引起的體積應變。

4)應變率和速度梯度應滿足相容方程：

?！?5)

1.3工程應用難點

FLAC3D英文手冊給出滲流特征時間的計算公式，特征時間可通過式(6)—(8)[1]計算。

Lc=Vf/Af。

(6)

(7)

(8)

式中： Lc為特征長度； Vf為滲流區(qū)域體積； Af為滲流區(qū)域表面積； tc為特征時間； k為FLAC3D滲透系數(shù)； n為孔隙率； Kf為流體體積模量； K為圍巖排水條件體積模量； G為圍巖排水條件剪切模量； α為比奧系數(shù)。

對于簡單的一維滲流，容易得到特征長度和特征時間。如圖3所示，現(xiàn)有一厚度為20m的飽和土層，底部為一不透水剛體，頂部作用均布荷載，水通過頂部排除，所以特征長度取土體厚度，即20m[12]；但在實際工程中，地下水的運動是三維的，難以確定特征時間和滲流模型。

太沙基把土的滲透固結(jié)簡化為一個裝滿水的容器，容器上方有一開小孔的活塞，活塞下面用一根彈簧支撐。當荷載剛作用上去時，小孔中的水來不及排出，由水和彈簧共同承擔荷載。當荷載作用長時間后，多余的水分排出，彈簧承擔全部的荷載[13]。本文的研究思路是將土的滲流固結(jié)過程分段，建立各段與滲流模型之間的對應關(guān)系。例如對應排水固結(jié)階段的滲流模型，理論上不會出現(xiàn)超孔隙水壓。

圖3　一維滲流算例

2模型設計

2.1模型試算

為了保證建模正確性，現(xiàn)進行模型試算。根據(jù)相關(guān)文獻[3、14]，建立青島膠州灣海底隧道流固耦合模型，如圖4所示。由于部分參數(shù)(開挖時間、加固區(qū)尺寸)未在文獻中找到，所以未知參數(shù)按經(jīng)驗取值。

圖4　青島膠州灣海底隧道模型

圖5與圖6為試算結(jié)果同原文的對比，結(jié)果存在一定的差異，主要是由于部分參數(shù)取值不同而引起的，但圍巖豎向位移、內(nèi)力均在同一數(shù)量級，變形收斂趨勢是相同的。文獻[3]還就數(shù)值計算結(jié)果同模型試驗進行了對比，若不考慮裂隙水、圍巖產(chǎn)狀等因素，數(shù)值解和試驗結(jié)果相當接近。

圖5　青島膠州灣海底隧道試算結(jié)果

圖6　青島膠州灣海底隧道試算結(jié)果

流固耦合計算費時，可以沿隧道中心面取一半建模，以加快計算速度，如圖7所示。需要注意隧道中心面處孔隙水壓力邊界條件，通常設置為不透水邊界，即孔隙水壓是可變化的。

圖7　圓形隧道模型

2.2參數(shù)設計

本文的目的是為了分析4種滲流模式的受力本質(zhì)，從而與實際工程建立聯(lián)系，就如同鋼桁梁橋總體受力雖然復雜，但每根鋼桿依然準守基本的力學公式。本文以簡單模型為例，更便于分析說明，所得結(jié)論，亦適用于復雜模型。

現(xiàn)設一100m×100m×10m立方體，水位高出地面10m。初始平衡以后，在頂部施加一豎向均布荷載q，大小為500kPa。在模型中心從上到下依次取點A—J10個監(jiān)測點，監(jiān)測其初始豎向總應力σ0、初始孔隙水壓p0、加載后豎向總應力σ1、加載后孔隙水壓p1。模型示意圖見圖8。

這是卒中后吞咽障礙患者營養(yǎng)管理模式的創(chuàng)新與醫(yī)院發(fā)展戰(zhàn)略一致的大背景之一。在項目中，多學科團隊的協(xié)作，在鄭大一附院也有著深厚的基礎和濃郁的氛圍。據(jù)劉章鎖介紹，多學科管理團隊在醫(yī)院落地已有7年之久?！澳壳埃T診有28個多學科會診團隊?！?/p>

圖8　模型示意圖

為了防止模型屈服破壞，力學模式采用理想彈性模型，滲流模式采用各項同性滲流模型，參數(shù)取整，詳見表2和表3。

表2　模型物理力學參數(shù)

表3　模型水力學參數(shù)

2.3工況設計

FLAC3D內(nèi)置4種滲流模型，各模型力學進程、流體進程、流體體積模量Kf有所差別。通常流體不考慮壓縮變形，土體的變形是由于流體遷移和土顆粒骨架空隙間的變形。FLAC3D將流體考慮成可壓縮變形的，相關(guān)資料請參考FLAC3D用戶手冊?，F(xiàn)根據(jù)不同進程和流體模量設置4種工況，見表4 。

表4　計算工況

3計算結(jié)果及分析

3.1工況1計算結(jié)果

工況1是在初始平衡后，頂部施加荷載500 kPa，維持力學進程開啟和流體進程關(guān)閉，流體體積模量Kf為0，計算至平衡。工況1應力豎向分布曲線如圖9所示。由圖9可知： 曲線p0和曲線p1重合，說明加載后孔隙水壓不變；有效應力增大約為500 kPa，即外界施加的荷載由土顆粒骨架承擔，水不分擔壓力。

圖9　工況1應力豎向分布曲線

3.2工況2計算結(jié)果

工況2同工況1類似，維持力學進程開啟和流體進程關(guān)閉。與工況1不同的是，此時需要將流體體積模量Kf設置為真實值。工況2應力豎向分布曲線見圖10。由圖10可知，外界荷載一部分由土顆粒骨架傳遞，一部分由流體承擔。

在試算時還發(fā)現(xiàn)流體承擔的荷載大小，同流體和固體的剛度比Rk有關(guān)。新增1組對比模型，將流體的體積模量降低10倍，以分析剛度比Rk對結(jié)果的影響。剛度比

(9)

式中：Kf為流體體積模量;n為孔隙率;K為固體體積模量;G為固體切變模量。

圖10　工況2應力豎向分布曲線

圖10中曲線σ1(柔性)和曲線σ1(剛性)重合。通過對比發(fā)現(xiàn)，流體分擔外荷載比例同剛度比Rk成正比。因此可作出一個假設：土顆粒和流體分擔外力原理類似圖11。剛度比Rk就是2根彈簧的彈性系數(shù)k比值，當剛度比Rk等于0時，代表流體的彈簧彈性系數(shù)k等于0，外荷載全由另一根彈簧——土顆粒骨架全部承擔。工況1正是剛度Rk等于0的特殊狀態(tài)。

圖11　流體-固體剛度比示意圖

3.3工況3計算結(jié)果

由于工況3是模擬孔隙水壓變化引起的內(nèi)力重分布，所以初始平衡之后，不再施加500 kPa荷載，而是將頂部孔隙水壓提高500 kPa。另外，工況3和工況4需要設置時間，這里將時間取一個大值，使模型滲流運動場趨于穩(wěn)定。

工況3計算過程可視作2個獨立的步驟： 第1步是孔隙水壓變化引起滲流場重分布，有效應力場不發(fā)生變化；第2步是滲流場重分布后，有效應力場變化發(fā)生。工況3應力豎向分布曲線見圖12。圖12中曲線p1(M off/F on)和曲線p1(M on/F off)重合，說明有效引力場改變并未引起孔隙水壓力重分布。結(jié)果再次證明本文3.2節(jié)假設模型的正確性——Rk=0時，流體不參與受力的。

M代表力學進程；F代表滲流進程。

圖12工況3應力豎向分布曲線

Fig. 12Distribution curves of vertical stress in Case 3

在分析實際工程時，孔隙水壓的改變通常是人為控制的。比如開挖隧道時，將開挖后的單元孔隙水壓設置為0，圍巖周圍單元的孔隙水壓會根據(jù)已知孔隙水壓力邊界條件重新分布。

3.4工況4計算結(jié)果

工況4中流體參與分擔外荷載，并且分擔的數(shù)值與時間有關(guān)。現(xiàn)在設置2組對比模型，計算時間分別為6 000 s和60 000 s，其余參數(shù)相同。工況4應力豎向分布曲線見圖13。圖13中曲線σ1(6 000 s)和曲線σ1(60 000 s)重合，隨著時間的的增加，超孔隙水壓消散，有效應力增大。結(jié)合圖9和圖10，工況4可以近似看作工況1和工況2的中間狀態(tài)。當計算時間足夠長時，結(jié)果近似于工況1；當計算時間十分短時，工況4的結(jié)果近似于工況2結(jié)果。

由此得到這樣一個結(jié)論——模型D適用于任何情況。分析透水性好的砂土排水固結(jié)或是透水性差的黏土不排水固結(jié)，都可以選用模型D，只要時間按真實情況取值即可；但流固耦合相當費時，所以模型D只是“理論上”適用于各種情況。根據(jù)大量的試算還發(fā)現(xiàn)，不僅合理的滲流模型能節(jié)約時間，采取合理的網(wǎng)格形狀，避免單元尺寸過小，亦能減小單元數(shù)量和增大時間步time-step，從而縮短時間。

圖13　工況4應力豎向分布曲線

4結(jié)論與建議

1)模型A、模型C外荷載全部由土顆粒骨架承擔，模型B、模型D流體也分擔外荷載。

2)流體分擔外荷載比例與剛度比Rk成正比，模型D流體分擔比例還與作用時間有關(guān)。

如透水砂層在外荷載作用下長期變形，宜選用模型A；透水不良黏土層突然受外力作用，宜選用模型B；抽水引起基礎變形，宜選用模型C；黏土層排水固結(jié)過程可采用模型D。如果將土的滲透固結(jié)過程看作是用手在壓一個裝滿水的容器的話[13]，那么FLAC3D中4個模式對應關(guān)系參考圖14(筆者注： 圖11與圖14的模型是在大量試算基礎上，為了便于理解而假設的定性模型)。

圖14　對應關(guān)系示意圖

3)FLAC3D流固耦合正是通過調(diào)整剛度比Rk，打開和關(guān)閉力學-流體進程來實現(xiàn)的；但是由于流固耦合十分費時，一般通過調(diào)整流體體積模量Kf或比奧模量M減少計算時間。此外，手冊還指出計算滲流穩(wěn)定狀態(tài)時，流體比奧模量要滿足M>αLzρωg/n[1]，以確保流體自由面穩(wěn)定。相關(guān)理論還待作進一步研究。

參考文獻(References)：

[1]Itasca Consulting Group Minneapolis. Fast language analysis of continua in 3 Dimensions[EB/OL].Version 3.0, User’s manual.Minneapolis: Itasca Consulting Group,2005.

[2]李術(shù)才,宋曙光,李利平,等. 海底隧道流固耦合模型試驗系統(tǒng)的研制及應用[J]. 巖石力學與工程學報,2013,32(5): 883-890.(LI Shucai, SONG Shuguang, LI Liping, et al. Development on subsea tunnel model test system for solid-fluid coupling and its application [J]. Chinese Journal of Rock Mechanics and Engineering, 2013,32(5): 883-890. (in Chinese))

[3]蔚立元,李術(shù)才,徐幫樹,等. 水下隧道流固耦合模型試驗與數(shù)值分析[J]. 巖石力學與工程學報,2011,30(7): 1467-1474.(YU Liyuan, LI Shucai, XU Bangshu, et al. Study of solid-fluid coupling model test and numerical analysis of underwater tunnels [J]. Chinese Journal of Rock Mechanics and Engineering, 2011,30(7): 1467-1474. (in Chinese))

[4]陳宇,朱繼文. 雙圓盾構(gòu)隧道的流固耦合分析[J]. 地下空間與工程學報,2010,6(1): 21-27.(CHEN Yu, ZHU Jiwen. Coupled fluid-mechanical analysis of DOT shield tunnel [J]. Chinese Journal of Underground Space and Engineering,2010,6(1): 21-27. (in Chinese))

[5]馮懷平，張伏光，岳祖潤. 考慮流固耦合作用的深基坑有限元分析[J]. 地下空間與工程學報,2012,8(2): 286-291.(FENG Huaiping, ZHANG Fuguang, YUE Zurun. Numerical analysis of excavation considering the coupling of seepage field and stress field [J]. Chinese Journal of Underground Space and Engineering,2012,8(2): 286-291.(in Chinese))

[6]劉泉聲,劉學偉. 多場耦合作用下巖體裂隙擴展演化關(guān)鍵問題研究[J]. 巖土力學,2014,35(2): 305-320.(LIU Quansheng, LIU Xuewei. Research on critical problem for fracture network propagation and evolution with multifield coupling of fractured rock mass [J]. Rock and Soil Mechanics,2014,35(2): 305-320. (in Chinese))

[7]謝東海,馮濤,趙延林,等. 裂隙煤巖體的流固耦合精細模型[J]. 中南大學學報(自然科學版),2013,44(5): 2014-2021.(XIE Donghai, FENG Tao, ZHAO Yanlin,et al. Fluid-solid coupling detailed model of rock-coal mass [J]. Journal of Central South University(Science and Technology),2013,44(5): 2014-2021.(in Chinese))

[8]陳育民,徐鼎平. FLAC/FLAC3D基礎與工程實例[M]. 北京: 中國水利水電出版社,2013. (CHEN Yumin, XU Dingping. FLAC/FLAC3D basis and engineering instance[M]. Beijing: China Water & Power Press, 2013. (in Chinese))

[9]彭文斌. FLAC3D實用教程[M]. 北京: 機械工業(yè)出版社,2008. (PENG Wenbin. FLAC3D practical tutorial [M]. Beijing: China Machine Press, 2008.(in Chinese))

[10]李廣信. 高等土力學[M]. 北京: 清華大學出版社,2004. (LI Guangxin. Advanced soil mechanics [M]. Beijing: Tsinghua University Press, 2004.(in Chinese))

[11]王光欽. 彈性力學[M]. 北京: 中國鐵道出版社,2008. (WANG Guangqin. Elastic mechanics [M]. Beijing: China Railway Publishing House, 2008. (in Chinese))

[12]劉波,韓彥輝. FLAC3D原理實例與運用指南[M]. 北京: 人民交通出版社,2005.(LIU Bo, HAN Yanhui. FLAC3D principle examples and application guide[M]. Beijing: China Communications Press, 2005. (in Chinese))

[13]劉成宇. 土力學[M]. 北京: 中國鐵道出版社,2002. (LIU Chengyu. Soil mechanics [M]. Beijing: China Railway Publishing House, 2002. (in Chinese))

[14]馬棟,譚忠盛. 青島膠州灣海底隧道施工及機械配套技術(shù)[J]. 中國工程科學,2009,11(7): 45-52.(MA Dong, TAN Zhongsheng. Key tunneling technology and configuring big machines in Qingdao Jiaozhouwan subsea tunnel [J].Engineering Science,2009,11(7): 45-52. (in Chinese))

Discussion on Seepage Model Based on FLAC3D Solid-fluid Coupling

DENG Siyuan1, YANG Qixin1，*, JIANG Yajun1, CHEN Hao2

(1.KeyLaboratoryofTransportationTunnelEngineeringofMinistryofEducation,SchoolofCivilEngineering,SouthwestJiaotongUniversity,Chengdu610031,Sichuan,China; 2.TrackwayInstitute,JiangsuProvincialCommunicationsPlanningandDesignInstituteCo.,Ltd.,Nanjing210014,Jiangsu，China)

Abstract:There are 4 seepage models in the FLAC3D finite element software; it is decided by selecting proper seepage models whether the result of numerical calculation is right or not. The seepage field distribution in three-dimensional space is complicate; as a result, it is difficult to decide a quantitative yardstick. There are different selection methods in different references. It is necessary to make stress analysis on the 4 seepage models, so as to choose a simple, rapid and rational model selection method. Four construction cases based on the 4 seepage models are designed, with a 100 m×100 m×10 m isotropic elastic cube as an example. The total stresses and seepage field distributions of the 4 construction cases are compared. The conclusions drawn are as follows: 1) For model A and model C, the soil skeleton carries all loads; for other 2 models， the fluid carries partial load too. 2) The load sharing ratio of fluid is based on the stiffness ratio and time. The solid-fluid coupling process in FLAC3D is realized by adjusting the stiffness ratio and switching mechanics-fluid process. So far, the calculation of solid-fluid coupling is almost based on FLAC3D in China. The conclusions can be applied to foundation engineering and underground engineering, and the paper can provide reference for similar projects.

Keywords:solid-fluid coupling; seepage model; FLAC3D; geotechnical engineering

中圖分類號：U 45

文獻標志碼：A

文章編號：1672-741X(2016)02-0179-07

DOI:10.3973/j.issn.1672-741X.2016.02.009

作者簡介：第一 鄧思遠(1990—)，男，四川宜賓人，西南交通大學橋梁與隧道專業(yè)在讀碩士，主要研究方向為隧道及地下工程設計和施工技術(shù)。E-mail: dengsiyuan_2715@126.com。*通訊作者： 楊其新，E-mail: yangqixin113@163.com。

基金項目：四川省交通科技項目計劃任務(2012C14-1)； 中央高?；究蒲袠I(yè)務費專項資金資助(2682014CX065)

收稿日期：2015-11-15; 修回日期: 2015-11-28