淺埋矩形頂管求解的復(fù)變函數(shù)實踐

李新源，劉國彬

(同濟(jì)大學(xué) 土木工程學(xué)院，上海 200092)

淺埋矩形頂管求解的復(fù)變函數(shù)實踐

李新源，劉國彬

(同濟(jì)大學(xué) 土木工程學(xué)院，上海 200092)

同時考慮洞周邊界與地表邊界，采用復(fù)變函數(shù)法求解淺埋矩形頂管施工引起的土體應(yīng)力場、位移場。根據(jù)黎曼存在定理、復(fù)變函數(shù)的三角插值理論，提出計算含矩形的半無限域到同心圓環(huán)域共形映射的計算方法。在此基礎(chǔ)上，將邊界條件等式兩邊都展成洛朗級數(shù)，采用冪級數(shù)解法求得了復(fù)應(yīng)力函數(shù)的系數(shù)。從該解法與有限元解的對比來看，在大部分點處結(jié)果相差不大，誤差在2%左右。結(jié)果表明：(1)提出的保角映射函數(shù)形式及系數(shù)的求解方法適用于淺埋矩形隧洞；(2)提出的復(fù)變函數(shù)解法具有步驟清晰、收斂快、操作簡單等特點。

半無限體；矩形頂管；復(fù)變函數(shù)；冪級數(shù)法

隨著社會的發(fā)展，城市建設(shè)活動中的頂管工程也越來越多?，F(xiàn)已有許多方法來求解頂管施工引起的土體附加應(yīng)力場、位移場[1]。Verruijt[2]首先將復(fù)變函數(shù)法用于求解淺埋圓形隧洞開挖問題；Strack[3]得到了隧洞受浮力作用時的解；童磊[4]給出了隧道任意變形下圍巖的應(yīng)力場、位移場的解答；江學(xué)良等[5]給出了地面集中荷載作用下淺埋圓形洞室的應(yīng)力解答；Lu等[6]得到了地表荷載及圍巖自重作用下淺埋圓形隧洞圍巖的應(yīng)力場。雖然已有許多淺埋隧洞開挖問題的復(fù)變函數(shù)解法，但目前都僅限于淺埋圓形隧洞，淺埋矩形隧洞開挖問題的復(fù)變函數(shù)解法鮮見文獻(xiàn)報道。本文提出了含矩形隧洞的半無限域到同心圓環(huán)域保角映射函數(shù)的形式，對淺埋矩形隧洞開挖問題進(jìn)行了求解，為掌握隧洞應(yīng)力分布，指導(dǎo)隧洞設(shè)計、施工提供了依據(jù)。

1 問題的描述

假設(shè)圍巖為線彈性體并在開挖前已沉降穩(wěn)定，不考慮重力作用時淺埋矩形隧洞開挖問題可簡化為如圖1所示平面應(yīng)變問題，圖中地平線邊界Γ1為零應(yīng)力邊界，隧洞邊界Γ2為已知位移邊界，h為隧洞中心埋深，θ為隧洞邊界上點的幅角。隧洞周圍巖體的應(yīng)力場、位移場是由Γ2的位移引起的。

設(shè)σij、ui(i，j=1，2)為圖1中圍巖在直角坐標(biāo)系下的應(yīng)力、位移分量，則他們可用兩個復(fù)應(yīng)力函數(shù)φ(z)、ψ(z)表示為[2，7]

圖1 淺埋矩形隧洞開挖問題Fig.1 Shallow cavern with rectangular excavation crosssection

其中，φ(z)、ψ(z)的表達(dá)式可根據(jù)圖1的邊界條件式(1)、(2)求得：

2 共形映射函數(shù)

映射函數(shù)的求取是采用復(fù)變函數(shù)求解的第一步。本文所求的映射函數(shù)應(yīng)能將圖1中的D映射為圖2中的同心圓環(huán)域D′。基于前人研究，作者提出了此映射函數(shù)的形式：

式中：K、A、Ci均為系數(shù)。系數(shù)可采用“兩步走”求解步驟，其具體過程如下：

(1)采用恰當(dāng)?shù)挠成浜瘮?shù)z=z(w)將z平面的半無限域D映射為w平面上的域D'。

(2)在第一步映射的基礎(chǔ)上，采用特定的映射函數(shù)w=w(ζ )將w平面上的外邊界Γ'1映射為ζ平面上的單位圓環(huán)Γ''1。

根據(jù)Schinzinger[8]的研究成果，可將A選為0.3，且根據(jù)ρ=K/2h=0.3來選取K的值，之后根據(jù)邊界對應(yīng)條件來求取Ci與r的值。

3 復(fù)應(yīng)力函數(shù)求解過程

由于φ(ζ )、ψ(ζ )可以展開成洛朗級數(shù)的形式：

圖2 保角映射后的圖形Fig.2 The region after conformal mapping

因此問題最終歸結(jié)為根據(jù)邊界條件求解ak、bk、ck、dk的值。設(shè)ζ =ρσ，其中σ =eiα，α為ζ平面上邊界點的幅角。對于地表邊界ζ =rσ，由式(5)得

在實際計算時，設(shè)Ek的項數(shù)為s，當(dāng)k＞s時，Ek、E-k均為0。Ek的值可根據(jù)式(8)左右兩側(cè)σ同冪次系數(shù)相等來求得。

式(3)中x=h邊界條件可表示為

比較上式等號左右兩邊σ的同冪次系數(shù)可得一組求解ak、bk、ck、dk的方程

考慮洞口邊界條件，將得到另一組求解方程。(4)式中隧洞邊界條件可表示為

對于地表邊界ζ =σ有

在實際計算時取Hk的項數(shù)為s。Hk的值可根據(jù)式(13)左右兩側(cè)σ同冪次系數(shù)相等來求得。

將式(12)、(13)代入式(11)并化簡后可得

比較上式等號左右兩邊關(guān)于的同冪次系數(shù)得

聯(lián)立式(10)、(15)即可求得ak、bk、ck、dk的值。需要注意的是，根據(jù)Verruijt的研究，僅靠式(10)、(15)是無法求得唯一的ak、bk、ck、dk，在實際計算中也驗證了這個結(jié)論。需根據(jù)ak、bk、ck、dk的收斂性來補(bǔ)充一個方程[2]，首先令a0=0計算當(dāng)k足夠大時ak的值(設(shè)此時ak=a)，再令a0=1計算當(dāng)k足夠大時ak的值(設(shè)此時ak=b)，則合適的a0取值為

之后再根據(jù)式(10)、(15)計算ak、bk、ck、dk的值，將這些求得的值代回式(6)、(7)即可得到所要求的復(fù)應(yīng)力函數(shù)。

4 應(yīng)力及位移的求解

根據(jù)求得的φ(z)、ψ(z)及式(1)、(2)即可求出某一點處的σij、ui，需要注意的是，在點∞處由式(2)計算得到的ux不一定為0，設(shè)此點處ux=，由于本文將參考點取在此處，為滿足改點ux=0的條件，將每一點處所求的ux都減去，即：

由于增加一個剛體位移并不會改變應(yīng)力，顯然這樣處理后依然滿足應(yīng)力邊界條件。

5 算例及與其驗證

設(shè)某一矩形頂管其隧洞幾何尺寸及坐標(biāo)如圖1所示，圖中a=2 m，b=2 m，h=10 m。根據(jù)地質(zhì)報告，取泊松比v=0.3，彈性模量E=10 MPa。矩形頂管存在均勻與非均勻收斂模式，在實際工程中其位移模式一般為非均勻變形，故可設(shè)=設(shè)計要求地層損失率控制在5‰以內(nèi)，可得

在將所有的邊界條件都表示成的級數(shù)形式后，利用Matlab將4節(jié)求解過程編寫為程序，其求解由計算機(jī)自動完成。在實際求解過程中，φ(ζ )、ψ(ζ )的項數(shù)、Ek的項數(shù)、Gk的項數(shù)、Hk的項數(shù)不可能取無限，一般可取一個較大的數(shù)，本文求解過程中分別取為100、50、50、50。同時為了驗證本文求解的正確性，將所求結(jié)果與有限元結(jié)果相比較。有限元采用線彈性模型，采用與解析解相同的參數(shù)，并采用平面應(yīng)變模型進(jìn)行模擬?？紤]到頂管的埋深約為10 m，且寬度較窄而長度較長，空間效應(yīng)明顯，故本模型寬度方向取距頂管邊約5倍開挖深度即總寬度為10 m×5×2=100 m，在深度方向上取10 m×5=50 m。

5.1 位移分析

圖3為頂管施工引起的地表沉降對比圖。從圖中可以看出：復(fù)變函數(shù)計算所得的沉降與有限元結(jié)果吻合較好，經(jīng)計算兩者的誤差在1%左右；矩形頂管開挖引起的地層位移也近似呈正態(tài)分布，反彎點大致位于1.5倍的開挖深度處，地表沉降主要發(fā)生在4倍的開挖深度范圍內(nèi)；埋深越大沉降越大，地表沉降最小，其值約為洞室邊界沉降的70%左右。

5.2 應(yīng)力分析

圖3 地表位移分布Fig.3 Deformation of the surface

圖4為不同位置處的σx、σv應(yīng)力分布圖，圖中的負(fù)值表示壓應(yīng)力。從圖中可以看出：除部分角點處，復(fù)變函數(shù)計算所得的應(yīng)力結(jié)果與有限元結(jié)果吻合較好，經(jīng)計算兩者的誤差在1%左右；在角點處有限元所得的應(yīng)力集中系數(shù)比本文結(jié)果大得多，這是因為有限元在建模時采用了直角，引起了較大的應(yīng)力集中；對σx而言，最大的壓應(yīng)力發(fā)生在頂板與側(cè)幫的角點處，最大壓應(yīng)力值約為40 kPa，最大拉應(yīng)力發(fā)生在側(cè)幫與底板的角點處，最大拉應(yīng)力約為35 kPa；對σv而言，最大的壓應(yīng)力發(fā)生在頂板與側(cè)幫的角點處，最大壓應(yīng)力值約為30 kPa，最大拉應(yīng)力發(fā)生在側(cè)幫與頂板的角點處，最大拉應(yīng)力約為25 kPa。

有限元結(jié)果可以認(rèn)為是精確解，從本文解法與有限元解的對比來看，在大部分點處本文結(jié)果與有限元結(jié)果相差不大，因此總體來說本文的解法是精確的，可以滿足工程精度要求。

圖4 不同位置處應(yīng)力分布圖Fig.4 Stresses curves at different position

6 結(jié)論

1)半無限體矩形隧洞開挖問題復(fù)變函數(shù)法的難點在于保角映射函數(shù)的求取。采用本文提出的先將含矩形隧洞的半無限域映射為中間域，再將中間域映射為同心圓環(huán)域的“兩步走”的求解方法可較方便地求解出該函數(shù)，且便于應(yīng)用現(xiàn)有研究成果。

2)從本文解法與有限元解的對比來看，在大部分點處本文結(jié)果與有限元結(jié)果相差不大，誤差在2%左右。本文給出的求解方法求解過程步驟清晰、收斂快、精度高，易于編程實現(xiàn)，有很強(qiáng)的可操作性。

3)對于類似本文的工程，矩形頂管開挖引起的地層位移近似呈正態(tài)分布，反彎點大致位于1.5倍的開挖深度處，地表沉降主要發(fā)生在4倍的開挖深度范圍內(nèi)；頂管施工引起的地表變形約為洞室邊界位移的70%左右。

4)對σx而言，最大的壓應(yīng)力發(fā)生在頂板與側(cè)幫的角點處，最大壓應(yīng)力值約為40 kPa，最大拉應(yīng)力發(fā)生在側(cè)幫與底板的角點處，最大拉應(yīng)力約為35 kPa；對σv而言，最大的壓應(yīng)力發(fā)生在頂板與側(cè)幫的角點處，最大壓應(yīng)力值約為30 kPa，最大拉應(yīng)力發(fā)生在側(cè)幫與頂板的角點處，最大拉應(yīng)力約為25 kPa。

[1] 陳子蔭.圍巖力學(xué)分析中的解析方法[M]. 北京：煤炭工業(yè)出版社，1994：46-73.

[2] VERRUIJT A. A complex variable solution for a deforming circular tunnel in an elastic half-plane[J].International Journal for Numerical & Analytical Methods in Geomechanics，1997，21(21)：77-89.

[3] STRACK O E，VERRUIJT A. A complex variable solution for a deforming buoyant tunnel in a heavy elastic halfplane[J]. International Journal for Numerical & Analytical Methods in Geomechanics，2002，26(12)：1235-1252.

[4] 童 磊，謝康和，盧萌盟，等.盾構(gòu)任意襯砌變形邊界條件下復(fù)變函數(shù)彈性解[J].浙江大學(xué)學(xué)報：工學(xué)版，2010 (9)：1825-1830.

[5] 江學(xué)良，楊 慧，曹 平.邊坡下伏地下圓形洞室的彈性應(yīng)力解析[J].計算力學(xué)學(xué)報，2012，29(1)：62-68.

[6] LU Aizhong，ZENG Xiantai，XU Zhen. Solution for a circular cavity in an elastic half plane under gravity and arbitrary lateral stress[J].International Journal of Rock Mechanics & Mining Sciences，2016，89：34-42.

[7] 陳行威，宋振森.加勁十字形軸壓桿考慮初始扭轉(zhuǎn)缺陷的扭轉(zhuǎn)位移函數(shù)[J].河北工程大學(xué)學(xué)報：自然科學(xué)版，2016，33(3)：8-12.

[8] SCHINZINGER R，LAURA PAA. Conformal mapping：methods and applications[M]. New York City：Dover Publications，2003.

A complex variable solution for rectangle pipejacking in elastic half-plane

LI Xinyuan，LIU Guobin
(School of Civil Engineering，Tongji University，Shanghai 200032，China)

Considering the boundary and the surface conditions，the stress fi eld and displacement fi eld caused by the construction of the shallow rectangular pipe jacking are solved by the complex function method. According to the Riemann’s existence theorem and basic complex variable theory，a conformal mapping function which can transform the half-plane with a rectangle cavity into the concentric ring is established. The both sides of boundary conditions equation are expanded into Laurent series，and then the coef fi cients of complex stress function are solved by power series method. The derived solution is applied to an example and a comparison is made using FEM method to show the accuracy of the methods，the results of this paper are almost the same as those of the FEM method，and the error is about 2%. The result shows：(1) The method presented in this paper is applicable to a shallow buried rectangular tunnel；(2) The complex function method proposed in this paper is characterized by clear steps，fast convergence and simple operation.

half-plane；rectangle pipejacking；complex variable solution；power series method.

TU4

A

1673-9469（2017）04-0001-04

10.3969/j.issn.1673-9469.2017.04.001

2017-08-11

國家自然科學(xué)基金資助項目(51378389)；國家重點基礎(chǔ)研究發(fā)展計劃(“973”計劃)項目(2015CB057806)

李新源(1986-)，男，山東濰坊人，博士，主要從事地下結(jié)構(gòu)受力研究。