水平非均勻場綜合管廊地震響應(yīng)及離散元模擬研究*

黃德龍 宗鐘凌 劉強(qiáng) 岑航 湯愛平

1.江蘇海洋大學(xué)土木與港海工程學(xué)院 連云港222005

2.哈爾濱工業(yè)大學(xué)土木工程學(xué)院 150090

引言

城市地下綜合管廊的建設(shè)可以有效解決直埋管線所帶來的“馬路拉鏈”、“空中蜘蛛網(wǎng)”等“城市病”［1］。截止到2020 年，中國已有25 個(gè)試點(diǎn)城市開始推進(jìn)管廊的建設(shè)并逐步投入使用。這些城市，有50%位于Ⅶ度設(shè)防區(qū)，有20%位于Ⅷ度設(shè)防區(qū)，因此管廊的抗震研究勢在必行。

在1995 年日本阪神大地震時(shí)，地鐵車站等地下結(jié)構(gòu)遭遇嚴(yán)重毀壞，學(xué)者們才開始關(guān)注地下結(jié)構(gòu)的抗震問題［2］。管廊等地下結(jié)構(gòu)體系的破壞主要以外壁混凝土開裂以及內(nèi)容物破壞為主。其破壞的直接原因主要是土-結(jié)相互作用（SSI），與地面建筑物所產(chǎn)生的慣性效應(yīng)不同［3，4］。因此，場地振動(dòng)、地震引起土體失效以及土-結(jié)相互作用是現(xiàn)階段研究地下結(jié)構(gòu)抗震的熱點(diǎn)問題。近些年，國內(nèi)外學(xué)者們對管廊等地下結(jié)構(gòu)的抗震性能進(jìn)行了大量動(dòng)力試驗(yàn)和數(shù)值模擬，部分學(xué)者還進(jìn)行了現(xiàn)場震害調(diào)查，提出了一系列理論和技術(shù)方法。

對于試驗(yàn)研究，振動(dòng)臺試驗(yàn)是主要采取的測試手段［5-11］，得出了均勻土體中管廊的地震響應(yīng)規(guī)律以及管廊-土-地面結(jié)構(gòu)（TSF）體系的響應(yīng)規(guī)律，并對一致和非一致激勵(lì)下土-結(jié)相互作用機(jī)理進(jìn)行了深入研究。也有一些學(xué)者利用離心機(jī)試驗(yàn)進(jìn)行地下結(jié)構(gòu)抗震的研究［12-16］，探究了管廊和地震波之間的相互影響機(jī)制、土體致密化響應(yīng)特點(diǎn)，并同樣對土-結(jié)相互作用進(jìn)行了深入探索。對于數(shù)值模擬研究，學(xué)者們同樣主要關(guān)注均勻土或豎向成層土內(nèi)的土-結(jié)相互作用、管廊及節(jié)點(diǎn)的地震響應(yīng)、管廊抗震設(shè)計(jì)簡化計(jì)算方法以及管廊-土-地面建筑之間的動(dòng)力相互影響［4，17-25］。對于震害現(xiàn)場調(diào)查研究［1］，主要是針對地震作用下管廊的破壞形式進(jìn)行總結(jié)和分析。

綜合管廊可跨河床或不同地質(zhì)單元，并且只要不是膨脹土，就可選用挖出土回填，造成管廊周圍土體的水平非一致性。上述試驗(yàn)、數(shù)值模擬及現(xiàn)場震害調(diào)查中，鮮有學(xué)者考慮到穿越水平不均勻土體的管廊地震響應(yīng)機(jī)理問題，并且多是運(yùn)用有限元理論，無法精確模擬土顆粒與結(jié)構(gòu)的接觸，更無法模擬土顆粒遷移規(guī)律。本文基于振動(dòng)臺試驗(yàn)以及改進(jìn)離散元對上述不足進(jìn)行了探究。

1 試驗(yàn)與模擬方法

本文首先利用振動(dòng)臺試驗(yàn)得到管廊的動(dòng)力響應(yīng)，并將試驗(yàn)與改進(jìn)離散元方法的結(jié)果進(jìn)行對比，以驗(yàn)證改進(jìn)離散元方法的適用性。

1.1 振動(dòng)臺試驗(yàn)

本次振動(dòng)臺試驗(yàn)在哈爾濱工業(yè)大學(xué)結(jié)構(gòu)與抗震實(shí)驗(yàn)中心進(jìn)行。試驗(yàn)首先預(yù)制縮尺管廊，見圖1，然后采用單向振動(dòng)臺，并利用疊層剪切箱，將砂土與黏土分兩側(cè)放入剪切箱，管廊埋置深度為0.4m。針對El-Centro 波、Kobe 波和人工波，分別經(jīng)過調(diào)幅（0.1g／0.2g／0.4g／0.8g／1.0g）進(jìn)行加載。試驗(yàn)中分別監(jiān)測上覆土體的豎向位移、管廊外壁土壓力以及管廊的應(yīng)變，其中在砂土和黏土內(nèi)各敷設(shè)7 個(gè)土壓力計(jì)，為E-1～E-14，用來測量管廊側(cè)壁（E-1～E-3 處于砂土，E-4～E-6 處于黏土）受到的土壓力增量以及管廊上（E-10～E-14）、下（E-8～E-10）的土壓力增量。具體的試驗(yàn)設(shè)置方案參照黃德龍等［26］關(guān)于穿越水平非均勻場地管廊振動(dòng)臺試驗(yàn)。

圖1 縮尺管廊試驗(yàn)前準(zhǔn)備Fig.1 Preparation work before scale utility tunnel test

1.2 離散元模擬

1.改進(jìn)離散元方法基本原理

離散元（PFC）的計(jì)算思想是使每個(gè)元素均滿足牛頓第二定律，用中心差分的方法求解各元素的運(yùn)動(dòng)方程，得到研究對象的整體運(yùn)動(dòng)形態(tài)，而地震波對結(jié)構(gòu)的作用是一個(gè)能量不斷積聚的過程，和PFC的思想具有一定的差異性。因此本研究中改進(jìn)的離散元方法是基于能量的思想，利用微分正交法（DQM）和Newmark-β 方法來改進(jìn)離散元，提高計(jì)算精度和計(jì)算效率。

綜合管廊外壁的總勢能V是管廊壁的應(yīng)變能U、動(dòng)能K和外部土體做的功W之和，應(yīng)變能可以寫成：

動(dòng)能可表示為：

由于土體的壓力而產(chǎn)生的外部功可以通過下式獲得：

使用Hamilton原理推導(dǎo)出結(jié)構(gòu)的控制方程，原理如下：

為了解決工程領(lǐng)域中出現(xiàn)的初始值或者邊界值問題，有很多數(shù)值方法，其中一個(gè)較好的數(shù)值方法是DQM法。與其他數(shù)值方法相比，DQM 法有以下幾個(gè)優(yōu)點(diǎn)：（1）是一種解決非線性微分方程導(dǎo)數(shù)近似的精確方法；（2）可以滿足各種邊界條件，并且需要更少的方程和編程工作；（3）準(zhǔn)確性和收斂性都很高?；谏鲜鰞?yōu)點(diǎn)，近年來，DQM法在結(jié)構(gòu)和動(dòng)力學(xué)問題分析的數(shù)值求解中越來越受歡迎。基于此方法，結(jié)構(gòu)控制方程的矩陣形式［27］可以寫成式（5）：

式中：KL、KNL、C、M、db和dd分別代表線剛度矩陣、非線性剛度矩陣、阻尼矩陣、質(zhì)量矩陣、邊界點(diǎn)和域點(diǎn)。

利用Newmark-β方法來獲得結(jié)構(gòu)的時(shí)程響應(yīng)序列，公式（5）可以改寫如下［28］：

基于DQM 及Newmark-β 的計(jì)算原理，對PFC2D程序進(jìn)行二次開發(fā)，以求解式（6）非線性離散元運(yùn)動(dòng)方程。PFC2D程序中開發(fā)自定義求解原理，需先對模型設(shè)置特有的編號，避免與標(biāo)準(zhǔn)模型沖突。在C++程序中，依次修改源文件中CM_vep：：CM_vep（）、CM_vep：：Name（）、CM_vep：：PropNames（）、CM_vep：：ReturnProp（）、CM_vep：：AcceptProp（）、CM_vep：：SaveRestore（）等成員函數(shù)，完成源文件的修改。最后，將源文件編譯成動(dòng)態(tài)鏈接庫文件（DLL 文件），可供PFC2D程序調(diào)用。

2.離散元方法的設(shè)置

圖2 為管廊的PFC2D離散元模型，此模型融入管廊及其土體的Newmark動(dòng)力求解方法，基于黃景崎等［29］的研究，在模型外側(cè)施加等效荷載來實(shí)現(xiàn)試驗(yàn)中疊層剪切箱的作用效果。對于模型的構(gòu)建，Potyondy等［30］指出：由于離散元對土體粒子數(shù)目無法按照真實(shí)情況來模擬，雖然土體顆粒的尺度因子對離散元是有影響的，但是當(dāng)L'／d＞32，即剪切箱的最大邊長與顆粒粒徑之比大于32時(shí)，此時(shí)宏觀響應(yīng)與顆粒尺寸幾乎無關(guān)，即此時(shí)顆粒的大小對管廊的響應(yīng)幾乎沒有影響。所以模擬時(shí)為了減少粒子數(shù)目和計(jì)算所需要的時(shí)間，采用比較大的球形顆粒對土粒進(jìn)行建模。模型分為兩個(gè)區(qū)，接近管廊的區(qū)域顆粒小，其最大的L'／d的值為180；離管廊較遠(yuǎn)的區(qū)域?yàn)轭w粒較大區(qū)，其最大的L'／d為125，兩者都大于32，所以此模型可用來模擬管廊及周圍土體的地震響應(yīng)。本研究粒子數(shù)為14733。

圖2 管廊的二維離散元模型Fig.2 2D discrete element simulation model of utility tunnel

整個(gè)布置區(qū)域與試驗(yàn)振動(dòng)臺的豎向截面尺寸一致，其中管廊和土箱分別采用wall來模擬。管廊周圍相近的土體孔隙比都采用0.15，砂土采用顆分的形式定義，阻尼比為0.3，黏性土的阻尼比為0.45。對于構(gòu)建砂土模型來說，顆粒與顆粒之間采用線彈性模型，如圖3a。此模型需要土體的剛度系數(shù)，以土體的基床系數(shù)為依據(jù)［31］，基于本試驗(yàn)的土體物理性質(zhì)，法向剛度和切向剛度取為3.3 ×104N／m，滑動(dòng)摩擦系數(shù)取為0.4；構(gòu)建黏土模型顆粒與顆粒之間同樣采用線彈性模型，其中法向剛度取為2.8 ×

圖3 離散元粒子與管廊接觸模型Fig.3 The contact models between discrete element particles and utility tunnel

104N／m，切向剛度取為2.0 ×104N／m，滑動(dòng)摩擦系數(shù)取為0.4；對于砂土與管廊之間的接觸（即ball 與wall 接觸）采用線性平行粘結(jié)（linearpbond）模型，如圖3b?；谑剑?）中的土壓力，可以計(jì)算接觸法向剛度和切向剛度都取為2 × 106N／m，摩擦系數(shù)取0.4，阻尼比取為0.2；對于黏土與管廊之間的接觸，也同樣采用線性平行粘結(jié)模型，接觸法向剛度和切向剛度都取為1.0 ×106N／m，摩擦系數(shù)取0.4，阻尼比均取為0.4。

3.離散元與試驗(yàn)結(jié)果對比

基于DQM 計(jì)算方法，加速度、位移和應(yīng)變的響應(yīng)結(jié)果與試驗(yàn)對比如圖4 所示?？梢园l(fā)現(xiàn)，加速度和位移的計(jì)算結(jié)果吻合較好。

圖4 離散元模擬和試驗(yàn)結(jié)果對比（El-Centro，PGA ＝1.0g）Fig.4 Comparison of numerical simulation and experimental results（El-Centro，PGA ＝1.0g）

在模擬中，由于是二維模型并且土體顆粒尺寸與真實(shí)顆粒相差較大，導(dǎo)致顆粒之間的剛度和阻尼無法按照真實(shí)情況模擬。這種情況下，粒子的宏觀運(yùn)動(dòng)特性可以準(zhǔn)確實(shí)現(xiàn)，微應(yīng)變等微觀變形響應(yīng)的時(shí)程曲線將很難與試驗(yàn)準(zhǔn)確對應(yīng)，但是其幅值均在可接受范圍內(nèi)。

2 試驗(yàn)與模擬結(jié)果討論

2.1 豎直位移

圖5 為輸入PGA 從0.1g 調(diào)幅至1.0g（增幅0.2g）的三種地震波后，利用激光位移計(jì)所測得的砂土和黏土豎向位移。可以看出砂土數(shù)值增大明顯，砂土區(qū)發(fā)生了剪沉現(xiàn)象（沉降），而黏土區(qū)數(shù)值呈減小趨勢，發(fā)生了剪脹現(xiàn)象。隨著加速度的增大，剪沉和剪脹增加量逐漸減小，即剪沉和剪脹雖然逐漸累積，但累積速度在放緩。到最后加載PGA ＝1.0g的人工波后，發(fā)現(xiàn)砂土依然在剪沉，黏土同樣也發(fā)生了輕微剪沉現(xiàn)象，但是其量值很小。說明隨著加載的進(jìn)行，砂土的密度逐漸增大，而黏土的密度逐漸減小到一定值后還有增加的趨勢。對于三種波在PGA ＝0.8g時(shí)出現(xiàn)了與其他地震動(dòng)強(qiáng)度不同的情況，即砂土出現(xiàn)了明顯的變形上浮，主要是因?yàn)槿N波在PGA ＝0.4g均發(fā)生了明顯的剪沉。

圖5 不同土體表面的豎向位移Fig.5 Vertical displacement of different soil surfaces

土體變形將會導(dǎo)致土壓力發(fā)生變化，而管廊所受的土壓力又是管廊發(fā)生變形的直接因素。圖6 為試驗(yàn)中砂土發(fā)生了剪沉現(xiàn)象，可以發(fā)現(xiàn)砂土流向與剪切箱接觸邊的中心，并逐漸向下流動(dòng)。

圖6 試驗(yàn)觀測砂土沉降Fig.6 Sand settlement observation

圖7 為土體發(fā)生豎向位移的改進(jìn)離散元二維截面模擬情況。此程序無法實(shí)現(xiàn)內(nèi)部管廊的轉(zhuǎn)動(dòng)作用，但是從圖7a 中可以發(fā)現(xiàn)左側(cè)砂土區(qū)發(fā)生了剪沉。土體與管廊之間出現(xiàn)了明顯的空隙，之后砂土區(qū)左側(cè)土體補(bǔ)充此空隙，但是有明顯的滯后效應(yīng)。這與圖6 砂土區(qū)剪切箱邊砂土發(fā)生的流動(dòng)變化相一致。并且圖7b 土顆粒的運(yùn)動(dòng)也可以說明砂土主要是向下沉降，而黏土是向側(cè)壁擴(kuò)展，進(jìn)而導(dǎo)致剪脹現(xiàn)象的發(fā)生。

圖7 土體豎向移動(dòng)模擬Fig.7 Simulation of vertical movement of the soil

2.2 管廊上部土壓力增量

圖8 為El-Centro 波，PGA ＝0.4g 和1.0g，分別對應(yīng)砂土和黏土區(qū)20cm和40cm深度處（管廊上部）的土壓力增量ΔσE時(shí)程?？梢园l(fā)現(xiàn)深度越深，土壓力增量越大；小震時(shí)砂土區(qū)的土壓力較大，而大震時(shí)，黏土區(qū)的土壓力較大；其中部分土壓力增量迅速達(dá)到一個(gè)平衡階段，然后在該位置震蕩直到結(jié)束，見圖8a PGA ＝0.4g 的曲線，這與Cilingir 和Madabhushi［32］發(fā)現(xiàn)的規(guī)律一致。Hushmand等［14］給出的解釋是土體致密化會導(dǎo)致土壓力增量增大，土壓力最后可能有殘余值。結(jié)合圖5也可以發(fā)現(xiàn)，隨著加載的進(jìn)行，在PGA 增大的同時(shí)，砂土和黏土各自發(fā)生剪沉和剪脹現(xiàn)象，其殘余內(nèi)力在增大，表現(xiàn)為土壓力增量在增大。圖9 為利用改進(jìn)的離散元程序?qū)虞d1.0g El-Centro波的土體內(nèi)部不平衡力變化過程進(jìn)行模擬，此不平衡力可以反映土壓力的大小。由于El-Centro波在3s左右加速度達(dá)到最大峰值，因此對模擬的前3s中6 個(gè)時(shí)刻的土壓力進(jìn)行研究（限于篇幅，圖9 中僅給出3 個(gè)代表時(shí)刻）?？梢园l(fā)現(xiàn)隨著加載時(shí)刻的進(jìn)行，加速度逐漸增大，砂土區(qū)管廊與土體之間的空隙逐漸增大。由于上部土體的補(bǔ)充，此空隙逐漸減小。可以看出土壓力整體在逐漸增大。開始時(shí)，由于砂土區(qū)下部土體沉降，此區(qū)域土體內(nèi)部不平衡力逐漸增大，即土壓力增大。隨著加速度的增大，管廊上覆土體的土壓力也在增大，達(dá)到加速度峰值1.0g時(shí)，管廊上部土壓力超過下部，在邊界以及兩種土體分界區(qū)土壓力最大。

圖8 不同土體的土壓力增量（El-Centro 波）Fig.8 Earth pressure increments for different soils（El-Centro）

圖10 表示利用改進(jìn)的離散元程序?qū)GA ＝1.0g時(shí)三種波作用下的土壓力峰值，可以發(fā)現(xiàn)Kobe波的土壓力峰值最小。

圖10 三種地震波PGA ＝1.0g 時(shí)土壓力峰值模擬結(jié)果（單位：kPa）Fig.10 Simulation of peak earth pressure at 1.0g for three waves（unit：kPa）

2.3 管廊側(cè)壁土壓力增量

圖11 分別為砂土區(qū)和黏土區(qū)管廊側(cè)壁沿縱向受到的土壓力增量ΔσE時(shí)程，可以發(fā)現(xiàn)砂土區(qū)土體對側(cè)壁土壓力增量略大。并且側(cè)壁中間（砂土E-2 和黏土E-5）土壓力增量較大，其原因是處于側(cè)壁中間處的土體受到側(cè)壁的約束，殘余內(nèi)力無法自由釋放，而導(dǎo)致對側(cè)壁土壓力的增大。同樣可以發(fā)現(xiàn)砂土區(qū)E-2 發(fā)生了土壓力增量的負(fù)向殘余值，其原因是縱向運(yùn)動(dòng)的過程中，砂土與管廊有垂直側(cè)壁方向的脫離趨勢。

圖11 沿豎直方向分布土壓力增量響應(yīng)（El-Centro，PGA ＝1.0g）Fig.11 Response of earth pressure increments along vertically distributed（El-Centro，PGA ＝1.0g）

柔度比F表示結(jié)構(gòu)周圍土體與結(jié)構(gòu)之間的相對剛度［2］。Wang［3］和Hashash等［4］都對矩形地下結(jié)構(gòu)的柔度比進(jìn)行了深入研究，其計(jì)算公式如下：

式中：Gm為在自由場中沿管廊高度土體的平均應(yīng)變所對應(yīng)的剪切模量；K 為管廊的剛度系數(shù)，可通過結(jié)構(gòu)力學(xué)的方法確定，即在限制管廊底部水平運(yùn)動(dòng)的同時(shí)，在管廊頂部施加單位力所引起的管廊頂部橫向位移的倒數(shù)；W為管廊截面的寬度；H為管廊截面的高度。

本研究的工況只有砂土和黏土兩種，只涉及兩個(gè)柔度比，即砂土區(qū)柔度比F1和黏土區(qū)柔度比F2。由于本研究砂土的剪切模量Gm-sand近似為7.7MPa，黏土的剪切模量Gm-clay近似為3.5MPa，而兩個(gè)工況柔度比所涉及的其他參數(shù)均相同，因此F1＞F2。

圖12 為砂土和黏土區(qū)沿管廊高度分布土壓力增量峰值變化曲線?？梢园l(fā)現(xiàn)隨著PGA 的增大，土壓力增量峰值有增大的趨勢，但是由于土體和管廊有分離的趨勢，部分土壓力增量減小。小震PGA ＝0.1g時(shí)，在側(cè)壁中間高度處，土壓力增量都要比PGA ＝0.2g時(shí)大，這是由于地震波以PGA ＝0.1g開始加載，當(dāng)突然加載或波形（頻率）發(fā)生改變時(shí)，土壓力均將發(fā)生顯著變化。大震時(shí)砂土區(qū)土壓力增量要比黏土區(qū)大，即柔度比在一定范圍內(nèi)越大時(shí)，管廊相對于土體越柔，受到的土壓力越大。還可以發(fā)現(xiàn)隨著加速度的增大，土壓力增量分布近似由底部（頂部）較大向中間位置較大變化。Hushmand 等［14］指出隨著地震作用的增大，土壓力增量的分布近似由三角形變?yōu)閽佄锞€形，柔度比可以決定側(cè)壁土壓力的形狀［2］。本文得出的結(jié)果與上述結(jié)論相符，但是兩者不同點(diǎn)在于砂土區(qū)管廊下部的土壓力增量峰值要大于上部，而黏土區(qū)則相反。因此砂土區(qū)要著重關(guān)注管廊下部，而黏土區(qū)上部較危險(xiǎn)，這與2.1 節(jié)土體豎向變形有密切關(guān)系。

圖12 沿管廊側(cè)壁高度土壓力增量峰值變化Fig.12 Variations of earth pressure increment peaks along the height of utility tunnel

在地震作用下，地基土體會發(fā)生沉降變形。管廊地基土體的沉降變形是破壞管廊穩(wěn)定性的直接原因。土體一旦出現(xiàn)沉降變形就會使管廊受力發(fā)生變化。受力超過允許值，管廊將遭到破壞。由地震作用所產(chǎn)生的管廊側(cè)壁水平向動(dòng)土壓力為［33］：

這表明擋板上的土壓力水平分量py隨深度增加趨于一個(gè)極值，而傳統(tǒng)的Coulomb 土壓力理論和Rankine土壓力理論都認(rèn)為土壓力沿深度線性增加，這說明地震動(dòng)力作用下的擋板土壓力分布與Coulomb土壓力理論和Rankine 土壓力理論有較大差別。這與圖12 管廊側(cè)壁土壓力呈現(xiàn)中間大兩端小的極值曲線相接近。

分析式（8）可以發(fā)現(xiàn)，土壓力的極值與側(cè)壓力系數(shù)K、土體重度γ、土體內(nèi)摩擦角φ、擋板與土體摩擦角δ、地震加速度系數(shù)kv、kh和剪切箱-管廊側(cè)壁凈距L有關(guān)。

2.4 地表豎向位移與管廊側(cè)壁土壓力的關(guān)系

圖13 為試驗(yàn)中地表豎向位移增量與管廊側(cè)壁土壓力增量的關(guān)系?？梢园l(fā)現(xiàn)，砂土區(qū)隨著豎向沉降增量峰值的增大，側(cè)壁土壓力增量在開始階段較大，砂土區(qū)側(cè)壁中間E-2 的土壓力增量在沉降增量為0.4mm～0.8mm 附近出現(xiàn)極小值。說明在發(fā)生沉降初期，側(cè)壁土壓力變化明顯，側(cè)壁中間土壓力增量最大，并且存在極小值。而黏土區(qū)的剪脹與土壓力增量關(guān)系存在極大值，當(dāng)剪脹量在0.2mm～0.3mm 附近時(shí)，側(cè)壁土壓力增量出現(xiàn)極大值。

圖13 管廊側(cè)壁土壓力增量峰值與地表豎向位移之間的關(guān)系Fig.13 Relationship between the vertical displacement of the ground surface and peak value of the earth pressure on the side wall of the tunnel

圖14 為利用改進(jìn)的離散元程序?qū)GA ＝1.0g時(shí)El-Centro波作用下砂土和黏土區(qū)管廊側(cè)壁土壓力峰值的模擬?？梢园l(fā)現(xiàn)砂土區(qū)土壓力整體要比黏土區(qū)土壓力大，并且側(cè)壁中部土壓力大于側(cè)壁上下兩側(cè)的土壓力，與圖12 曲線相一致。還可以發(fā)現(xiàn)管廊側(cè)壁的土壓力要比管廊上下壁土壓力大，管廊的破壞主要受側(cè)壁土壓力控制。

圖14 El-Centro 波管廊側(cè)壁土壓力峰值模擬（PGA ＝1.0g）（單位：kPa）Fig.14 Simulation of peak earth pressure corresponding to El-Centro wave（PGA ＝1.0g）（unit：kPa）

3 結(jié)論

本文基于振動(dòng)臺試驗(yàn)以及改進(jìn)的離散元數(shù)值方法對水平非均勻場地綜合管廊地震響應(yīng)進(jìn)行探究，通過利用改進(jìn)的離散元程序?qū)芾燃捌渲車馏w進(jìn)行模擬，并與試驗(yàn)進(jìn)行對比，得出如下結(jié)論：

1.基于能量思想，并利用微分正交法以及Newmark-β方法得到模型體系的運(yùn)動(dòng)控制方程，以此改進(jìn)離散元程序，建立PFC2D管廊模型。模擬結(jié)果與試驗(yàn)結(jié)果的誤差均在20%以內(nèi)，即二次開發(fā)的離散元數(shù)值模擬方法是準(zhǔn)確的；

2.砂土區(qū)土體發(fā)生了剪沉，黏土區(qū)發(fā)生了剪脹現(xiàn)象，并且伴隨著土體密度的變化，管廊發(fā)生了輕微轉(zhuǎn)動(dòng)，砂土區(qū)剪沉現(xiàn)象的發(fā)生是為了補(bǔ)充管廊下部砂土與管廊的空隙；

3.土體致密化會導(dǎo)致土壓力增大，土壓力增量可能會出現(xiàn)殘余值，土壓力的大小與所輸入PGA有關(guān)，不均勻土體對土壓力的影響較大。管廊側(cè)壁的土壓力增量呈中間大兩端小的極值形，并且在一定范圍內(nèi)，柔度比越大，管廊所受到的土壓力越大，管廊側(cè)壁變形也越大；

4.砂土區(qū)在發(fā)生沉降初期，側(cè)壁土壓力變化明顯，且之后存在極小值；黏土區(qū)發(fā)生剪脹，隨著剪脹量的減小，土壓力增量存在極大值。