FLAC 3D軟件在地鐵隧道地震動(dòng)力響應(yīng)的應(yīng)用

張蒙蒙 雷震宇

（同濟(jì)大學(xué)鐵道與城市軌道交通研究院，201804，上?！蔚谝蛔髡撸T士研究生）

1 FLAC 3D軟件

FLAC（Fast Lagrangian Analysisof Continua）是由美國(guó)ltasca consulting Group Inc.開(kāi)發(fā)的二維顯式有限差分法程序，采用適合于模擬大變形問(wèn)題的拉格朗日方法，即模型單元網(wǎng)格隨著變形而不斷更新。拉格朗日法按照時(shí)間步采用動(dòng)力松弛的方法來(lái)求解運(yùn)動(dòng)方程，采用混合離散化方法有效模擬模型材料的塑性破壞和流動(dòng)。

區(qū)別于其他有限元數(shù)值模擬方法，F(xiàn)LAC 3D軟件使用顯式有限差分方法求解。對(duì)于顯式法，非線性本構(gòu)關(guān)系與線性本構(gòu)關(guān)系沒(méi)有算法差別，可以方便地求出應(yīng)力增量、不平衡力，并跟蹤系統(tǒng)的演化過(guò)程。另外，在有限差分法中，基本方程組和邊界條件近似地改為差分方程表示。由于使用顯式方法按時(shí)間遞步解算代數(shù)方程避免了有限元方法中隱式矩陣解算微分方程，不需要建立剛度矩陣，省去了求解大型聯(lián)立方程組的步驟，所以占用計(jì)算機(jī)內(nèi)存少，節(jié)約了求解時(shí)間，提高了計(jì)算效率。

FLAC 3D軟件采用全動(dòng)力學(xué)運(yùn)動(dòng)方程，可以分析和計(jì)算非穩(wěn)定過(guò)程，即使在靜態(tài)模擬中，也利用動(dòng)態(tài)運(yùn)動(dòng)方程計(jì)算，使得數(shù)值模擬不存在力的不穩(wěn)定過(guò)程，消除了數(shù)值上的障礙，這是有限元方法不能解決的。

FLAC 3D軟件中單元材料可采用線性或非線性本構(gòu)模型。在外力作用下，當(dāng)材料發(fā)生屈服流動(dòng)后，網(wǎng)格能相應(yīng)發(fā)生形變及移動(dòng)，因此，在模擬大變形問(wèn)題及其他塑性問(wèn)題時(shí)，F(xiàn)LAC 3D軟件所采用的混合離散方法比有限元法中通常采用的離散集成法更為準(zhǔn)確合理。在動(dòng)力模塊中，F(xiàn)LAC 3D軟件采用完全非線性分析方法，適用于任何指定的非線性本構(gòu)模型；非線性材料定律使不同頻率的波之間可以自然地干涉和混合；可以自動(dòng)計(jì)算永久變形；采用合理塑性方程，使得塑性應(yīng)變?cè)隽颗c應(yīng)力相聯(lián)系；可以方便進(jìn)行不同本構(gòu)模型的比較；可以同時(shí)模擬壓縮波和剪切波的傳播及兩者耦合作用時(shí)對(duì)材料的影響。

FLAC 3D軟件中內(nèi)置了11種材料本構(gòu)模型，包括空單元模型、3種彈性模型、7種塑性模型。每種模型對(duì)應(yīng)1種特殊類(lèi)型巖土材料的本構(gòu)特性。模型的任何一個(gè)個(gè)體單元都可以采取不同的本構(gòu)模型或賦予不同的材料參數(shù)。

2 建模

2.1 模型建立

本文依照上海軌道交通8號(hào)線青云路段地質(zhì)狀況建模，該段主要的土層分布和土體物理參數(shù)如表1所示。地鐵隧道外徑為6m，隧道襯砌厚度為0.3m；C55混凝土，密度為2 600kg／m3，彈性模量E＝35 500 MPa；埋深分別采用9m、15m、30m。由于FLAC 3D軟件采用顯式有限差分方法進(jìn)行數(shù)值模擬，相比較于ANSYS等有限元軟件，其對(duì)于網(wǎng)格尺寸精細(xì)度的要求更低，因此可以采用較大的單元尺寸。

表1 各土層巖體物理參數(shù)

由文獻(xiàn)［1］知，模型尺寸取隧道直徑的7～8倍時(shí)邊界影響已較小，故建模為：模型寬度100m，厚度方向60m，縱向沿隧道方向延伸10m，巖土體采用摩爾－庫(kù)倫本構(gòu)模型，隧道襯砌采用彈性體本構(gòu)模型。另外，根據(jù)FLAC 3D軟件動(dòng)力計(jì)算中網(wǎng)格的尺寸取決于輸入波形的最高頻率，通常取最高頻率對(duì)應(yīng)波長(zhǎng)的1／10～1／8的理論，以及本文所采用的震動(dòng)波在各巖層波速及波長(zhǎng)推算出的模型最大單元尺寸，為方便建模采用1m×1m網(wǎng)格計(jì)算，共生成71400個(gè)單元，79 948個(gè)節(jié)點(diǎn)。根據(jù)以往的文獻(xiàn)和計(jì)算經(jīng)驗(yàn)，這個(gè)精度是滿足要求的［3］。

2.2 地震波的矯正和輸入

FLAC 3D軟件輸入的動(dòng)力荷載可以加載在模型邊界或內(nèi)部節(jié)點(diǎn)。加載荷載的方式為向模型輸入應(yīng)力時(shí)程、加速度時(shí)程、速度時(shí)程、集中力時(shí)程。荷載的輸入方式可為FLAC 3D軟件內(nèi)嵌函數(shù)Fish語(yǔ)言定義的規(guī)律波形，也可以由table表類(lèi)讀入離散波形。

在FLAC 3D軟件動(dòng)力計(jì)算中，網(wǎng)格的尺寸取決于輸入波形的最高頻率，通常取最高頻率對(duì)應(yīng)波長(zhǎng)λ的1／10～1／8。即輸入波的最大頻率越高，網(wǎng)格的尺寸要求越小。為了保證地震波在土體中的傳播精度，同時(shí)控制網(wǎng)格尺寸及數(shù)量以節(jié)約計(jì)算時(shí)間，在波形輸入之前需要進(jìn)行濾波，把波形中的高頻成分濾掉。另外，對(duì)于各個(gè)完整波形，輸入加速度時(shí)程積分得到的最終速度和最終面積均不為0，故動(dòng)力計(jì)算過(guò)程結(jié)束時(shí)會(huì)在模型底部出現(xiàn)繼續(xù)的速度和殘余位移，為此需要進(jìn)行基線校正。SeismoSignal軟件是一款利用惠更斯原理積分的波處理軟件，可以對(duì)加速度時(shí)程的1次積分和2次積分結(jié)果進(jìn)行修正。

本次數(shù)值模擬為1995年日本阪神地震波，利用SeismoSignal軟件將地震波中高于5Hz的波形成分過(guò)濾掉，并對(duì)加速度時(shí)程進(jìn)行了基線校正。

本文采用從模型底部水平（X）方向和豎直（Z）方向同時(shí)加載。經(jīng)濾波和校準(zhǔn)的水平向、豎直向地震波如圖1和圖2所示。水平向地震波峰值為0.6 g，豎直向地震波峰值為0.2 g，作用時(shí)間為20s。

圖1 經(jīng)過(guò)濾波和校準(zhǔn)的水平向地震波

2.3 確定邊界條件

在動(dòng)力計(jì)算中，模型邊界會(huì)存在波的反射，為了更接近真實(shí)的情況和減小數(shù)值計(jì)算負(fù)擔(dān)，采用合理的模型尺寸和利用合理邊界條件吸收入射波是十分重要的。FLAC 3D軟件動(dòng)力分析的邊界條件主要分為靜力邊界（粘性邊界）條件和自由場(chǎng)邊界條件。對(duì)于地震動(dòng)力荷載，在剛性地基上可采用自由場(chǎng)邊界，模型底部不必加靜態(tài)邊界；在柔性地基上需在模型周?chē)┘幼杂蓤?chǎng)邊界，同時(shí)在模型底部施加靜態(tài)邊界條件。

圖2 經(jīng)過(guò)濾波和校準(zhǔn)的豎直向地震波

本文模型采用自由場(chǎng)邊界。自由場(chǎng)邊界利用apply ff，apply ff＿corner，apply ff＿side命令在建立的基本模型四周及4個(gè)角生成平面網(wǎng)格和柱網(wǎng)格，通過(guò)在邊界設(shè)置切向和法向的獨(dú)立阻尼器與主體網(wǎng)格進(jìn)行耦合。自由場(chǎng)邊界中的獨(dú)立阻尼器可有效吸收來(lái)自模型內(nèi)部切向和法向的入射波，模擬真實(shí)的無(wú)限邊界中波和能量的傳播效應(yīng)。切向粘滯力fs和法向粘滯力fn分別為：

式中：

υn——邊界上的法向速度分量；

υs——邊界上的切向速度分量；

ρ——介質(zhì)的密度；

CP——P波的波速；

CS——S波的波速。

2.4 確定阻尼

FLAC 3D軟件采用動(dòng)力方程求解準(zhǔn)靜力問(wèn)題和動(dòng)力問(wèn)題，在計(jì)算中設(shè)定阻尼是必要的。對(duì)于動(dòng)力問(wèn)題中的阻尼，需要在數(shù)值模擬中重現(xiàn)自然系統(tǒng)在動(dòng)荷載加載時(shí)的阻尼大小。FLAC 3D軟件提供了瑞利阻尼（rayleigh damping）、局部阻尼（local damping）和滯后阻尼（hysteretic damping）。本文采用局部阻尼。

局部阻尼在振動(dòng)循環(huán)中通過(guò)在節(jié)點(diǎn)或結(jié)構(gòu)單元節(jié)點(diǎn)上增加或減小質(zhì)量的方法達(dá)到收斂。由于增加的單元質(zhì)量和減小的單元質(zhì)量相等，就一次總體來(lái)說(shuō)，系統(tǒng)保持質(zhì)量守恒。由于局部阻尼系數(shù)不需求解系統(tǒng)的自振頻率，且相對(duì)于瑞利阻尼不會(huì)減小時(shí)間步，所以在求解小型模型的簡(jiǎn)單問(wèn)題中有很大的優(yōu)勢(shì)。當(dāng)節(jié)點(diǎn)速度的符號(hào)改變時(shí)增加節(jié)點(diǎn)質(zhì)量，當(dāng)速度達(dá)到最大值（或最小值）時(shí)減小節(jié)點(diǎn)質(zhì)量，因此損失的能量ΔW／W 是最大瞬時(shí)應(yīng)變能Wd的一定比例，ΔW／W 是與頻率無(wú)關(guān)的。因?yàn)棣／W 是臨界阻尼比D的函數(shù)：

式中：

αL——局部阻尼系數(shù)；

D——臨界阻尼比，取值可參考瑞利阻尼的ζmin。

和祥軒的老板以搶占民女罪被縣衙抓走。最后的結(jié)局是不得不使上一筆不菲的贖金把人給贖回來(lái)。和祥軒老板仿佛是死過(guò)了一回，他親眼見(jiàn)證了那沾滿血跡的刑具是如何的寒威逼人。他只在那刑房里待了不到一個(gè)時(shí)辰，竟然就把屎尿全拉在了褲襠里。所以，以后，和祥軒的人只要瞧見(jiàn)那些個(gè)穿半頭鞋子的人都敬而遠(yuǎn)之。更要命的是，那案子還久拖而不得了結(jié)，和祥軒還得定期給官府送銀票續(xù)保。

對(duì)于巖土材料，臨界阻尼比一般取0.05，則局部阻尼系數(shù)可取0.157。

3 計(jì)算分析

3.1 初始地應(yīng)力

在施加動(dòng)力荷載前，需要對(duì)巖土層的密度、體積模量、剪切模量、黏聚力、摩擦角等物理參數(shù)賦值，并在豎直方向施加重力，solve軟件會(huì)自動(dòng)求解原始的地應(yīng)力情況以還原真實(shí)的巖土力學(xué)情況。

為了模擬真實(shí)的地鐵隧道周?chē)鷳?yīng)力狀況，分步對(duì)左右線隧道進(jìn)行開(kāi)挖和支護(hù)。考慮到地鐵隧道管片寬度為1.3m，因此采用的開(kāi)挖模擬方法為：每次沿隧道軸向利用Null模型向前掘進(jìn)1.3m，通過(guò)Solve軟件得到開(kāi)挖后應(yīng)力狀態(tài)，然后對(duì)前一環(huán)開(kāi)挖段添加初期支護(hù)并再次通過(guò)Solve軟件計(jì)算，直到計(jì)算平衡，再逐次向前推進(jìn)；先對(duì)左線開(kāi)挖、支護(hù)，之后依照相同方法對(duì)右線開(kāi)挖并支護(hù)，以此得到支護(hù)完成后的地應(yīng)力場(chǎng)（見(jiàn)圖3、圖4）。

圖3 埋深9m時(shí)左線開(kāi)挖支護(hù)后豎直應(yīng)力云圖

圖4 埋深9m時(shí)左右線開(kāi)挖支護(hù)后豎直應(yīng)力云圖

3.2 動(dòng)力計(jì)算和分析

從模型底部由X方向和Z方向在模型X邊界和Z邊界同時(shí)加載處理過(guò)的地震波，動(dòng)力作用時(shí)間10 s，記錄每1s的最大主應(yīng)力和最小主應(yīng)力、隧道及周?chē)馏w單元X向和Z向應(yīng)力、關(guān)鍵點(diǎn)位移和速度。

3.2.1 關(guān)鍵點(diǎn)位移

取左線隧道洞頂襯砌內(nèi)側(cè)和外側(cè)點(diǎn)、隧道肩部點(diǎn)、隧道左右兩側(cè)襯砌內(nèi)外側(cè)點(diǎn)作為主要監(jiān)測(cè)點(diǎn)（見(jiàn)圖5），記錄埋深為9m、12m、30m時(shí)，隧道頂、肩和兩側(cè)監(jiān)測(cè)點(diǎn)在震動(dòng)輸入的豎直方向和水平方向每秒鐘的位移情況，如圖6～11所示。其中，GP1為隧道頂部襯砌內(nèi)側(cè)一點(diǎn)；GP2為隧道頂部襯砌外側(cè)圍巖上一點(diǎn)；GP5為隧道肩部襯砌外圍巖上一點(diǎn)；GP10為隧道襯砌壁遠(yuǎn)離相鄰隧道方向的側(cè)面上一點(diǎn)；GP11為隧道襯砌壁相近相鄰隧道方向的側(cè)面上一點(diǎn)。

圖5 左線隧道各考察點(diǎn)位置

圖6 埋深9m時(shí)各考察點(diǎn)豎直位移變化

圖7 埋深9m時(shí)各考察點(diǎn)水平位移變化

圖8 埋深12m時(shí)各監(jiān)測(cè)點(diǎn)豎直位移變化

圖9 埋深12m時(shí)各監(jiān)測(cè)點(diǎn)水平位移變化

圖10 埋深30m時(shí)各監(jiān)測(cè)點(diǎn)豎直位移變化

圖11 埋深30m時(shí)各監(jiān)測(cè)點(diǎn)水平位移變化

從圖中可以看出，埋深不同，對(duì)于各點(diǎn)在豎直和水平地震波傳播方向的位移影響較大，隨著埋深增加，隧道頂、肩、兩側(cè)位移更加趨于一致，位移量減小，隧道洞壁和襯砌變形變小，更加不易破壞。這說(shuō)明：隨著上覆土層體積的增加，重力及周?chē)鷰r體性質(zhì)對(duì)隧道襯砌的約束增強(qiáng)，使得隧道襯砌對(duì)地震動(dòng)力的響應(yīng)減弱。

另外，隧道襯砌在整個(gè)地震動(dòng)力響應(yīng)過(guò)程中的最大位移響應(yīng)出現(xiàn)在第5s，這為之后的分析提供了依據(jù)。

圖12～17表示在位移響應(yīng)最大的第5s時(shí)隧道在埋深分別為9m、12m、30m情況下的豎直應(yīng)力云圖和水平應(yīng)力云圖（單位：Pa）。

表2是隧道頂、肩、兩側(cè)考察點(diǎn)周?chē)鷨卧谖灰浦捣逯?s時(shí)，埋深分別為9m、12m、30m情況下的豎直應(yīng)力值（szz）和水平應(yīng)力值（sxx），以及埋深增加后（9m→12m、12m→30m）應(yīng)力值增長(zhǎng)百分比。

圖12 埋深9m，5s時(shí)豎直應(yīng)力云圖

圖13 埋深9m，5s時(shí)水平應(yīng)力云圖

圖14 埋深12m，5s時(shí)豎直應(yīng)力云圖

圖15 埋深12m，5s時(shí)水平應(yīng)力云圖

圖16 埋深30m，5s時(shí)豎直應(yīng)力云圖

圖17 埋深30m，5s時(shí)水平應(yīng)力云圖

由圖12～17及表2可以看出：① 在埋深增加的情況下，隧道頂部在豎直方向拉應(yīng)力增長(zhǎng)，壓應(yīng)力區(qū)向兩肩擴(kuò)散，且數(shù)值逐漸增大；② 埋深越大，隧道單元應(yīng)力值越大，但增長(zhǎng)百分比越小，增長(zhǎng)幅度明顯減小。這說(shuō)明：在埋深較深的地方圍巖對(duì)隧道的約束越大，地震波對(duì)隧道襯砌的內(nèi)力影響相對(duì)越弱。

3.2.3 變形云圖

FLAC 3D軟件的后處理功能可以繪出模型的變形云圖，用以形象地表征變形的發(fā)生方向和相對(duì)大小。圖18～20采用模型的變形云圖、變形后網(wǎng)格疊加的方式顯示了埋深分別為9m、12m、30m，震動(dòng)持 續(xù)5s時(shí)（位移變形最大時(shí)刻）的位移變形情況。

表2 不同埋深下隧道各考察點(diǎn)周?chē)鷨卧獞?yīng)力值

圖18 埋深9m變形云圖

圖19 埋深12m變形云圖

圖20 埋深30m變形云圖

從圖中可以看出，上覆土層厚度為12m時(shí)的隧道變形情況較埋深為9m、30m時(shí)更加嚴(yán)重。這主要是由于隧道所處土層造成的，此時(shí)隧道處于淤泥質(zhì)粉砂黏土層和黏土層，相對(duì)于其他情況所處的粉質(zhì)黏土層，其彈性模量E、摩擦角f及內(nèi)聚力c較小，由于周?chē)鷰r土控制地震波的傳遞和衰減特性，地震波對(duì)于軟弱地層造成的影響顯然比堅(jiān)硬地層大得多。

4 結(jié)語(yǔ)

本文通過(guò)FLAC 3D軟件模擬地鐵雙線隧道在不同埋深情況下，在豎直方向和水平方向同時(shí)加載地震波時(shí)的應(yīng)力和位移變形情況，得出以下結(jié)論：

（1）上覆土層厚度12m時(shí)，隧道的變形情況比其他情況下更加嚴(yán)重。由于隧道處于淤泥質(zhì)粉砂黏土層和黏土層，其E，f，c相對(duì)于其他情況所處粉質(zhì)黏土層較小，根據(jù)周?chē)鷰r土控制地震波的傳遞和衰減的特性，因此軟弱地層地震波造成的影響顯然要大的多。

（2）在埋深增加的情況下，隧道頂部在豎直方向拉應(yīng)力增長(zhǎng)，壓應(yīng)力區(qū)向兩肩擴(kuò)散，且數(shù)值逐漸增大。

（3）埋深越大，隧道單元應(yīng)力值越大，但增長(zhǎng)百分比越小，增長(zhǎng)幅度明顯減小，說(shuō)明在埋深較深的地方圍巖對(duì)隧道的約束越大，地震波對(duì)隧道襯砌的內(nèi)力影響越弱。

（4）埋深不同，對(duì)于各點(diǎn)在豎直和水平地震波傳播方向的位移影響較大，隨著埋深增加，隧道頂、肩、兩側(cè)位移更加趨于一致，位移量減小，隧道洞壁和襯砌變形變小，更加不易破壞。

綜合各結(jié)果可知，隨著隧道上覆土層厚度的增加，隧道的破壞程度逐漸減小。這是由于土體重度增加和對(duì)隧道襯砌產(chǎn)生影響的巖體量增加，約束增強(qiáng)，使得地震波的動(dòng)力作用對(duì)隧道襯砌的位移和內(nèi)力影響減弱。

利用FLAC 3D軟件對(duì)地鐵隧道進(jìn)行初始應(yīng)力場(chǎng)、開(kāi)挖和簡(jiǎn)單支護(hù)、地震波作用模擬，可直觀地體現(xiàn)地鐵隧道在動(dòng)力作用下圍巖位移和應(yīng)力分布情況，并據(jù)此進(jìn)行一系列探討。相比其它有限元分析軟件，F(xiàn)LAC 3D軟件更加節(jié)約計(jì)算機(jī)內(nèi)存和計(jì)算時(shí)間，其強(qiáng)大的后處理功能使分析更加直觀和易得。因此，F(xiàn)LAC 3D軟件在隧道地震動(dòng)力響應(yīng)方面的應(yīng)用具有顯而易見(jiàn)的適用性和優(yōu)越性。由于時(shí)間等因素的限制，本文所采取的模型仍有許多值得改進(jìn)之處：① 模型底部邊界若采用靜態(tài)邊界限制或可得到更精確的分析結(jié)果；② 隧道的襯砌只采用簡(jiǎn)單的初期支護(hù)，或可采用更加精細(xì)的支護(hù)模擬以更加接近真實(shí)工況；③ 本文采取震動(dòng)每秒結(jié)束時(shí)的位移和應(yīng)力情況進(jìn)行分析，或可采用峰值進(jìn)行分析對(duì)比。這些問(wèn)題將在以后的研究工作中逐步完善。

［1］王鐵男.FLAC在地鐵隧道數(shù)值模擬中的應(yīng)用［J］.沈陽(yáng)大學(xué)學(xué)報(bào)，2010，22（1）：11.

［2］唐益群，欒長(zhǎng)青.地鐵振動(dòng)荷載作用下隧道土體變形數(shù)值模擬［J］.地下空間與工程學(xué)報(bào)，2008，4（1）：105.

［3］陳育民，徐鼎平.FLAC／FLAC 3D基礎(chǔ)與工程實(shí)例［M］.北京：中國(guó)水利水電出版社，2011：75－88，297－312.

［4］彭文斌.FLAC 3D實(shí)用教程［M］.北京：機(jī)械工業(yè)出版社，2011：207－249.

［5］李圍.隧道及地下工程FLAC解析方法［M］.北京：中國(guó)水利水電出版社，2009：238－245.

［6］陳貴紅.沉管隧道抗震數(shù)值分析［D］.成都：西南交通大學(xué)，2002.

［7］鄔玉斌，劉如山，楊學(xué)山.地鐵交通隊(duì)某實(shí)驗(yàn)樓的環(huán)境振動(dòng)的模擬分析［J］.世界地震工程，2008，24（4）：70.

［8］馬樂(lè)元，丁毅，郭傳創(chuàng).基于FLAC 3D的移動(dòng)荷載作用下基坑結(jié)構(gòu)響應(yīng)的分析［J］.工程與建設(shè)，2011，25（1）：90.

［9］曹恒將.考慮動(dòng)載的巷道頂板離層及控制的FLAC 3D模擬研究［D］.青島：山東科技大學(xué)，2010.

［10］劉波，葉圣國(guó)，陶龍光.地鐵盾構(gòu)施工引起鄰近基礎(chǔ)沉降的FLAC元數(shù)值模擬［J］.煤炭科學(xué)技術(shù)，2002，30（10）：9.

［11］蘇麗娟.高速鐵路隧道圍巖支護(hù)參數(shù)優(yōu)化設(shè)計(jì)［D］.北京：北京交通大學(xué)，2011.