張浩 呂達 陳飛 宋家才 韓金柯 張宇 孫明源 劉鵬 范佳林

1深圳市燃氣集團股份有限公司

2中國石油大學(xué)（華東）儲運與建筑工程學(xué)院

2020年5月，某燃氣管道由于第三方施工導(dǎo)致管線埋土下沉，管線承受較大作用力，致使管道存在斷裂風(fēng)險。受影響的燃氣管道位于人行道與車行道交界處，直徑D508，材質(zhì)為L360 鋼，管道厚度7.9 mm，運行壓力1.5 MPa，于2004 年投產(chǎn)運行。第三方施工項目開挖地基與D508 次高壓燃氣管線較近，最近距離約8 m（圖1）。項目施工前期，為保證基坑及周邊安全，施工單位對燃氣管道的變形與安全狀況進行了監(jiān)測。監(jiān)測數(shù)據(jù)顯示，次高壓燃氣管道間接沉降超出預(yù)定報警值，最大差異沉降為0.63‰，沉降未趨于收斂穩(wěn)定，存在安全隱患。為了控制管道沉降，避免管道事故的發(fā)生，施工單位進行了管底注漿加固工作，并進一步開展了管道應(yīng)力檢測，發(fā)現(xiàn)管道部分位置應(yīng)力已經(jīng)超過管道材料許用應(yīng)力。

圖1 第三方施工項目及受影響管道圖Fig.1 Diagram of third party construction projects and affected pipelines

由上述案例可知，在管道運營過程中，因地下工程、降雨、基坑開挖等因素不可避免的會造成土體下沉，嚴重威脅管道安全，輕則導(dǎo)致管道變形，重則導(dǎo)致管道斷裂。鑒于沉降引起管道事故的嚴重性，有必要研究沉降作用下埋地管道力學(xué)響應(yīng)。

理論分析、試驗?zāi)M和數(shù)值模擬是研究沉降作用下管道力學(xué)響應(yīng)、評價管道安全狀態(tài)的基本方法，國內(nèi)外學(xué)者已開展大量研究。在理論分析上，WINKLER[1]提出計算管道的變形受力的方法，將埋地管道看作彈性地基梁，管道周圍土體看作具有一定抗壓剛度系數(shù)的彈簧。日本學(xué)者高田至郎[2]通過選取沉陷區(qū)和非沉陷區(qū)的管道作為分析對象，運用彈性地基梁模型分析，得出地埋管線在塌陷情況下的簡化分析公式。高惠瑛[3]等利用三次曲線模擬沉陷區(qū)管道的幾何大變形，用彈性地基梁模型模擬非沉陷區(qū)的管道變形，推導(dǎo)出沉陷區(qū)的管道的位移及內(nèi)力計算公式。鄧道明[4]將災(zāi)害區(qū)外部管道看作半無限長的彈性地基梁模型，基于災(zāi)害體內(nèi)外管道內(nèi)力和變形的連續(xù)性，推導(dǎo)出管道應(yīng)力和位移的計算方法。張土喬[5]等分析管道的受力特性，把埋地管道視為連續(xù)彈性地基梁，計算出管道上的最大撓度、轉(zhuǎn)角和彎矩。關(guān)惠平等[6]對非沉陷區(qū)和沉陷區(qū)的管道分別構(gòu)建了計算模型，得出不同塌陷區(qū)長度、不同沉降量的管道內(nèi)部最大軸向應(yīng)力。在試驗?zāi)M上，陳志磊[7]、張旭[8]在室內(nèi)土工槽內(nèi)建立試驗?zāi)Ｐ?，模擬了由于土體塌陷或沉降引起管線變形、影響管線正常使用這一過程。王樂[9]通過實驗研究了管道在不均勻沉降作用下的力學(xué)性能，得到了不均勻沉降過程中管道的應(yīng)變值和管道受力特征。在數(shù)值模擬上，RAJEEV 等[10-12]使用有限差分軟件FLAC3D分析了膨脹土中管道的應(yīng)力變化，同時在大型管箱中對塑料管進行了實驗室實驗，實驗結(jié)果和模型預(yù)測之間有相當好的一致性，驗證了軟件分析的準確性。吳韶艷[13]等結(jié)合工程實際，比較基于有限差分軟件FLAC3D的數(shù)值模擬結(jié)果和現(xiàn)場變形監(jiān)測結(jié)果，驗證了FLAC3D 模擬的可行性。趙瀟[14]、王強[15]應(yīng)用有限元模擬軟件ABAQUS 研究了煤礦開采過程中埋地管道的變形和受力狀況。

由于目前研究中，大多僅提出理論公式，并未分析公式的可靠性，也未討論公式應(yīng)用到工程實際的可行性。為此，在基于燃氣管道沿線地面沉降的實際情況基礎(chǔ)上，將力學(xué)模型計算與數(shù)值模擬方法相結(jié)合，探究土體沉降作用下管道在不同工況時的受力特征，并對已有力學(xué)模型進行修正，提出適用于不同工程實際的力學(xué)模型計算方法。

1 沉降區(qū)管道力學(xué)計算

1.1 力學(xué)模型

對于沉降區(qū)埋地管道，管道兩端受到土體支撐且外側(cè)管道受到穩(wěn)定土體的約束，沉降區(qū)管道受到由于土體沉降造成的豎向載荷，根據(jù)Winkler 假設(shè)，鄧道明[4]將沉降區(qū)埋地管道簡化為彈性地基梁模型（圖2）。

圖2 管道力學(xué)模型Fig.2 Pipe mechanics model

對于沉降區(qū)內(nèi)管道受力計算公式為

式中：E為管材的彈性模量，Pa；I為管道的截面慣性矩，m4；y為管道撓度，m；x為管道位置，m；M0為管道在x=0處截面的彎矩，N·m；υ0為管道在x=0處的撓度，m；q為作用于管道上方的均布荷載，N/m；L為管道穿越沉降區(qū)的長度，m；S0為管道的當量軸向力，N。

計算方法見式（2）。

式中：N0為管道在x=0 截面的軸力，N；p為管道輸送壓力，Pa；d為管道內(nèi)徑，m。

結(jié)合邊界條件求解方程，沉降區(qū)邊界管道的截面彎矩計算見公式（3）。

式中：M0為沉降區(qū)邊界管道的截面彎矩，N·m；k為系數(shù)；β為特征系數(shù)；D為管道外徑，m。

沉降區(qū)中心管道的撓度y計算見公式（4）。

管道在滑坡段中點的截面彎矩計算見公式（5）。

式中：Mc為滑坡段中點的截面彎矩，N·m。

其中：

式中：Cy0為管道在沉降區(qū)外穩(wěn)定土層的縱向阻力系數(shù)，N/m3。

以沉降區(qū)內(nèi)管道為研究對象，利用L管段軸向應(yīng)變的幾何非線性關(guān)系，得到管道在x=0 處的撓度，如公式（8）所示：

式中：A為管道截面積，m2；α為鋼材的線膨脹系數(shù)，m/(m·℃)；ΔT為管體正溫差，℃；μ為泊松比，取0.3；σh為內(nèi)壓引起的管道環(huán)向應(yīng)力，MPa。

以沉降區(qū)外穩(wěn)定土層中管段為研究對象，由軸向位移的連續(xù)條件，利用半無限長埋地直管段的軸向位移一般解，確定υ0。

當(N0-Nr)γ/tf≤1 時，υ0計 算如公式（9）所示：

式中：γ，Nr為系數(shù)；tf為沉降區(qū)外土層的極限抗剪強度，N/m。

當(N0-Nr)γ/tf＞1 時，υ0計算如公式（10）所示：

Nr的計算如公式（11）所示：

γ的計算如公式（12）所示：

式中：cx0為沉降區(qū)外土層的軸向位移阻力系數(shù)，N/m3。

1.2 沉降區(qū)管道承受載荷

上述公式計算結(jié)果的準確性很大程度上取決于輸入?yún)?shù)的準確性。其中沉降區(qū)管道的承受荷載起到關(guān)鍵性作用。目前計算中對管道上方荷載通常僅考慮重力，而結(jié)合實際沉降過程中的管土相互作用分析可知，土體黏聚力和摩擦力的變化對于管土相互作用具有重要影響，據(jù)此在考慮土體、管道及管內(nèi)介質(zhì)重力的基礎(chǔ)上，考慮土體摩擦力以及黏聚力的作用[16]，對載荷q的計算進行修正。

自重計算如公式（13）所示：

式中：q為重力作用下管道受到的總荷載壓力，N/m；q1為管道自重作用于管道上的荷載壓力，N/m；q2為輸送介質(zhì)重力作用于管道上的荷載壓力，N/m；q3為管道上方土層重力作用于管道上的荷載壓力，N/m。

其中

管道自重q1計算如公式（14）所示：

式中：ρp為管道材料的密度，kg/m3；g為重力加速度，m/s2。

輸送介質(zhì)重量q2計算如公式（15）所示：

管道上方土層重量q3計算如公式（16）所示：

式中：h為管道埋深，m；ρsoil為管道上方覆土密度，kg/m3。

黏聚力計算qc計算如公式（17）所示：

式中：qc為黏聚力作用力，N/m；c為土壤黏聚力，Pa。

摩擦力計算qf計算如公式（18）所示：

式中：qf為摩擦力作用力，N/m；k0為系數(shù)，計算如公式（19）所示：

式中：φ為土的內(nèi)摩擦角(°)。

1.3 管道應(yīng)力計算過程

土體沉降作用下埋地管道軸向應(yīng)力計算的整個過程如圖3所示。公式求解采用牛頓迭代法，首先確定計算所采用的荷載壓力q；然后假定N0，逐步調(diào)整N0，使管道材料彈性變形計算過程中計算出的撓度近似相等，求出N0、M0；最終求得管道的軸向應(yīng)力。計算過程復(fù)雜，可通過編制MATLAB程序?qū)崿F(xiàn)快速計算。

圖3 迭代流程Fig.3 Iterative process

2 數(shù)值模擬

FLAC 3D軟件基于有限差分拉格朗日元法，對于計算巖土工程領(lǐng)域常見的幾何非線性大變形問題具有很大的優(yōu)勢。使用FLAC3D對沉降過程進行模擬，探究沉降作用下管道應(yīng)力分布規(guī)律。

FLAC3D 殼型結(jié)構(gòu)單元提供了殼體與網(wǎng)格的非線性接觸，考慮了管土間的法向壓縮、拉伸作用和剪切摩擦相互作用，可模擬管道與土體的接觸，因此本模擬選用該結(jié)構(gòu)單元建立管道模型。管材軸向拉伸曲線簡化為三折線模型取值。摩爾-庫倫模型適合于模擬抗壓強度和抗拉強度不同的巖石以及砂土粒狀材料，并且用常規(guī)實驗可以很容易獲得摩爾-庫倫模型的各項參數(shù)，使用起來簡單便捷，被廣泛應(yīng)用于巖土工程的數(shù)值分析[17]，因此本模擬采用摩爾-庫倫模型建立土體模型。

根據(jù)王樂[9]的研究，當采用固定邊界模型時，管道相當于具有固定邊界的簡支梁，為了避免出現(xiàn)邊界效應(yīng)，整個沉降模型長度應(yīng)至少為管道直徑的60倍，即Ltotal=L+60D。式中Ltotal為模型的軸向長度，m；L為沉降區(qū)長度，m；D為管道的直徑，m。結(jié)合張杰[18]的研究，模型的橫向和縱向尺寸可以取為10 倍管徑，綜上確定沉降模型規(guī)模按()L+60D×13D×10D（長×寬×高）建立。

土體選用六面體網(wǎng)格劃分，分別對沉降區(qū)、沉降區(qū)與非沉降區(qū)交界處附近及管道周圍的網(wǎng)格進行加密處理，對土體的非沉降區(qū)下表面及各側(cè)面分別施加法向約束，對管道兩端施加全約束。D508 管道由4 960 個結(jié)構(gòu)單元組成，土體由37 448 個實體單元組成（圖4）；D813管道由5 280個結(jié)構(gòu)單元組成，土體由49 764個實體單元組成。管道及沉降土體的性質(zhì)參數(shù)見表1 和表2，其中管道參數(shù)根據(jù)某燃氣公司實際管道選取，土體參數(shù)根據(jù)易發(fā)生沉降地區(qū)的土體實測得到。

表1 實驗管道參數(shù)Tab.1 Experimental pipeline parameters

表2 沉降土體基本參數(shù)Tab.2 Basic parameters of settlement soil

圖4 FLAC 3D模型圖Fig.4 Figure of FLAC 3D model

數(shù)值模擬整體分為兩個步驟：①對模型施加重力以實現(xiàn)地應(yīng)力平衡；②沉降區(qū)域底面施加豎向沉降位移。為了驗證修正后的力學(xué)模型計算方法的可行性，分別討論了不同沉降范圍、土體黏聚力及管道壁厚下管道的受力特征，分析了力學(xué)模型計算結(jié)果與模擬結(jié)果間的誤差。

3 模擬結(jié)果分析

3.1 沉降范圍對管道的影響

不同的管道沉降災(zāi)害可能出現(xiàn)不同的沉降范圍，圖5 為不同沉降范圍下管道軸向應(yīng)力分布圖。由圖5可知，沉降范圍的變化對于管道受力具有顯著影響。沉降范圍為10 m 時，D813 管道并未發(fā)生塑性變形，隨著沉降范圍增加，管道中部受力達到屈服極限，管道進入塑性變形階段。相同沉降范圍下，D813 管道的中間塑性變形段管道長度＜D508管道，這表明管道的管徑越大，管道中部的應(yīng)力越??；若管道中部應(yīng)力達到管材屈服應(yīng)力，則管道中間部位達到屈服應(yīng)力后發(fā)生塑性變形的管段長度減小。對于D508 的管道，在沉降范圍為30 m 時，沉降中心管道應(yīng)力明顯降低，這表明沉降范圍的增加并不能無限增加管道中間塑性變形段的長度。由此可見，在管道沉降災(zāi)害中，沉降邊界附近的管道受到破壞的概率最大，需要重點關(guān)注。

圖5 不同沉降范圍下不同管徑管道的軸向應(yīng)力分布Fig.5 Axial stress distribution of different pipeline diameters under different settlement ranges

3.2 土體黏聚力對管道的影響

黏聚力是表征土體性質(zhì)的重要指標之一，不同類型的土體具有不同的黏聚力，純砂土黏聚力幾乎為0，而黏土黏聚力往往超過10 000 Pa，根據(jù)模擬中的黏聚力5 660 Pa，取1/2和2倍數(shù)關(guān)系，模擬得到不同土體黏聚力下的管道軸向應(yīng)力分布如圖6所示。對于D508 管道，隨著黏聚力的增加，管道中部塑性變形管道范圍略有增加，且管道中部軸向應(yīng)力出現(xiàn)略為減小的趨勢；對于D813 管道，管道中部塑性變形管道范圍同樣略有增加，但塑性變形外管道軸向應(yīng)力變化更加劇烈。研究結(jié)果表明，土體黏聚力變化對管道受力的影響較小，隨著黏聚力增加，沉降區(qū)中間管道發(fā)生塑性變形的范圍略有增加。

圖6 不同土體黏聚力下不同管徑管道的軸向應(yīng)力分布Fig.6 Axial stress distribution of different pipeline diameters under different soil cohesion

3.3 管道壁厚對管道的影響

不同沉降災(zāi)害中管道的壁厚存在差異，上述沉降過程中，D508 管道的壁厚為7.9 mm，D813 管道的壁厚為11.9 mm，根據(jù)GBT 21835《焊接鋼管尺寸及單位長度重量》，對D508 管道另選6.3 mm 和8.8 mm壁厚，對D813管道另選10 mm和14.2 mm壁厚進行研究，圖7為不同的壁厚下的管道軸向應(yīng)力分布圖。由圖7 可知，壁厚對管道受力的影響很小，隨著壁厚增加，管道的中間塑性區(qū)范圍略有降低。

圖7 不同管道壁厚下不同管徑管道的軸向應(yīng)力分布Fig.7 Axial stress distribution of different pipeline diameters under different pipeline wall thicknesses

4 模擬與計算結(jié)果的對比

4.1 不同沉降范圍下的對比

4.1.1 D508管道

圖8 為不同沉降范圍下D508 管道力學(xué)模型計算得到的管道應(yīng)力與模擬所得數(shù)據(jù)的對比結(jié)果。沉降范圍對結(jié)果影響較大，沉降區(qū)長度為10 m 時，D508 管道的力學(xué)模型計算所得應(yīng)力曲線與模擬所得應(yīng)力曲線貼合度較好；沉降區(qū)長度為20 m、30 m時，D508 管道中間部位發(fā)生塑性變形，管道中部應(yīng)力不再繼續(xù)增加，而力學(xué)計算按照管道彈性開展計算，故力學(xué)模型計算所得最大應(yīng)力與模擬所得應(yīng)力在此段誤差較大，力學(xué)模型計算所得應(yīng)力明顯偏大。因此，力學(xué)模型計算適用于彈性變形區(qū)應(yīng)力大小計算，對于超過屈服強度的部分建議在該區(qū)段采取基于應(yīng)變的管道破壞準則計算。

圖8 D508管道不同沉降范圍下應(yīng)力模擬結(jié)果與理論計算結(jié)果對比Fig.8 Comparison between stress simulation results and theoretical calculation results under different settlement ranges of D508 pipeline

4.1.2 D813管道

圖9 為不同沉降范圍下D813 管道力學(xué)模型計算得到的管道應(yīng)力與模擬所得數(shù)據(jù)的對比結(jié)果。沉降區(qū)長度較小時，D813 管道的力學(xué)模型計算所得應(yīng)力曲線與模擬所得應(yīng)力曲線貼合度較好；沉降區(qū)長度較大時，D813 管道中間出現(xiàn)塑性區(qū)，導(dǎo)致力學(xué)模型計算結(jié)果在此段與模擬結(jié)果出現(xiàn)差別，且隨著沉降區(qū)長度的增大，管道塑性區(qū)內(nèi)的計算結(jié)果誤差不斷增大。

圖9 D813管道不同沉降范圍下應(yīng)力模擬結(jié)果與理論計算結(jié)果對比Fig.9 Comparison between stress simulation results and theoretical calculation results under different settlement ranges of D813 pipeline

4.2 不同土體黏聚力下的對比

4.2.1 D508管道

圖10 為不同土體黏聚力下D508 管道力學(xué)模型計算得到的管道應(yīng)力與模擬所得數(shù)據(jù)的對比圖。雖然模擬結(jié)果受黏聚力的影響較小，但是力學(xué)模型計算結(jié)果受影響較大。隨著黏聚力增大，力學(xué)模型計算所得最大應(yīng)力不斷增加，模擬所得應(yīng)力基本不變，因此對于該直徑的管道，力學(xué)模型計算結(jié)果與模擬結(jié)果在塑性段的誤差不斷增大。

圖10 D508管道不同土體黏聚力下應(yīng)力模擬結(jié)果與理論計算結(jié)果對比Fig.10 Comparison between stress simulation results and theoretical calculation results under different soil cohesion of D508 pipeline

4.2.2 D813管道

圖11 為不同土體黏聚力下D813 管道力學(xué)模型計算得到的管道應(yīng)力與模擬所得數(shù)據(jù)的對比圖。隨著黏聚力增大，力學(xué)模型計算結(jié)果與模擬結(jié)果貼合度不斷增加，但是當黏聚力增大至11 320 Pa 時，力學(xué)模型計算結(jié)果顯著大于模擬結(jié)果，呈現(xiàn)出過于保守的情況，已無法準確體現(xiàn)模擬所得管道應(yīng)力分布。在各類土體中，黏土黏聚力很大，由此可知，管道處于黏土環(huán)境下，對于該直徑的管道，采用上述計算方法得到的管道應(yīng)力十分保守。

圖11 D813管道不同土體黏聚力下應(yīng)力模擬結(jié)果與理論計算結(jié)果對比Fig.11 Comparison between stress simulation results and theoretical calculation results under different soil cohesion of D813 pipeline

4.3 不同管道壁厚下的對比

4.3.1 D508管道

圖12 為不同壁厚下D508 管道力學(xué)模型計算得到的管道應(yīng)力與模擬所得數(shù)據(jù)的對比圖。管道應(yīng)力計算結(jié)果受壁厚的影響很小，不同壁厚下，力學(xué)模型計算所得結(jié)果與模擬結(jié)果的誤差基本不變，但管道塑性變形區(qū)受壁厚的影響而變化，壁厚越大，模擬所得管道塑性區(qū)范圍越小。

圖12 D508管道不同管道壁厚下應(yīng)力模擬結(jié)果與理論計算結(jié)果對比Fig.12 Comparison between stress simulation results and theoretical calculation results under different pipeline wall thicknesses of D508

4.3.2 D813管道

圖13為不同壁厚下D813 管道力學(xué)模型計算得到的管道應(yīng)力與模擬所得數(shù)據(jù)的對比圖。不同壁厚下，力學(xué)模型計算所得結(jié)果曲線與模擬結(jié)果曲線較為貼合，隨著壁厚越大，貼合度不斷增加。

圖13 D813管道不同管道壁厚下應(yīng)力模擬結(jié)果與理論計算結(jié)果對比Fig.13 Comparison between stress simulation results and theoretical calculation results under different pipeline wall thicknesses of D813

5 結(jié)論

基于燃氣公司實際管道情況，通過數(shù)值模擬分析了土體沉降作用下管道在不同工況時的受力特征，結(jié)合修正后的力學(xué)模型計算方法，討論了計算結(jié)果與模擬結(jié)果間的誤差，得到如下結(jié)論：

（1）沉降范圍對管道的受力結(jié)果影響較大。沉降范圍增加，管道中部達到屈服應(yīng)力后發(fā)生塑性變形的管段明顯增加。但沉降范圍的增加并不能無限增加管道中間塑性變形段的長度，沉降范圍過大時管道中間部位不會發(fā)生塑性變形，而沉降邊界附近的管道始終發(fā)生塑性變形。

（2）土體黏聚力變化對沉降災(zāi)害下管道受力影響較小，而管道壁厚的變化幾乎對其無影響。

（3）修正后的力學(xué)模型計算結(jié)果在管道受力未超過屈服強度時與數(shù)值模擬結(jié)果具有較好的吻合度，超過管道屈服應(yīng)力時其計算結(jié)果不再可靠，需采取基于應(yīng)變的管道破壞準則進行計算。