喬 彤，周 建*，張?zhí)祢湥Y熠誠

(1.浙江大學(xué) 建筑工程學(xué)院，浙江 杭州 310012；2.濱海和城市巖土工程研究中心，浙江 杭州 310058)

對(duì)于水下隧道各類問題的研究，隧道滲流場(chǎng)是目前一個(gè)重要的方面[1-3]，絕大多數(shù)學(xué)者常將隧道截面假設(shè)為圓形以便求解。Zhang等[4]在無限承壓含水層井流模型的基礎(chǔ)上，利用鏡像法與疊加原理推導(dǎo)出了多水平地層條件下圓形隧道涌水量的等效計(jì)算公式，分析了隧道埋深、半徑和海水深度對(duì)涌水量的影響。Meye等[5]基于保角映射法，假定洞周土為各向同性多孔彈性介質(zhì)，推導(dǎo)了圓形隧道滲流量和周圍水頭分布公式，研究了不同防水排水類型下的滲流場(chǎng)。在等水壓和等水頭兩種水力邊界條件下，Zou等[6]采用坐標(biāo)變換法，將直角坐標(biāo)系中穩(wěn)定滲流控制方程求解問題轉(zhuǎn)化為雙極坐標(biāo)系中的求解問題，簡化了圓形隧道滲流量及襯砌外水壓力的求解過程。謝寒松[7]以復(fù)變函數(shù)的反演變換和地下水力學(xué)理論為基礎(chǔ)，推導(dǎo)了圓形隧道滲流場(chǎng)解析解，并對(duì)隧道與溶洞的間距和隧道半徑對(duì)隧道壓力水頭的影響進(jìn)行了分析。張丙強(qiáng)[8]采用鏡像法和滲流場(chǎng)疊加法則，推導(dǎo)了含注漿圈與襯砌作用下的圓形隧道滲流場(chǎng)解析解，且其對(duì)于非平行的雙孔水下隧道同樣適用。郭玉峰等[9]基于質(zhì)量守恒定律和達(dá)西定律，采用保角變換法和Schwartz迭代法推導(dǎo)了水下雙線平行圓形隧道穩(wěn)態(tài)滲流場(chǎng)的解析式，討論了隧道間距、埋深和相對(duì)大小對(duì)雙孔平行圓形隧道滲流場(chǎng)水頭分布、滲流量的影響?；贒upuit假設(shè)和積分法，Meng等[10]獲得了穩(wěn)定滲流條件下的各向同性垂直分層潛水層中圓形排水隧道滲流場(chǎng)解析解。Tang等[11]基于保角映射復(fù)變函數(shù)理論推導(dǎo)了考慮滲流各向異性的水下圓形隧洞滲流場(chǎng)解析解，但在映射過程中地層表面將發(fā)生變形，因此只能算是一種近似解?；诟飨蛲苑€(wěn)定滲流假定，Lin等[12]采用保角映射及最優(yōu)化技術(shù)推導(dǎo)了考慮滲流各向異性的水下圓形隧道滲流場(chǎng)的半解析式，并研究了等水壓和等水頭滲水模式下的隧道滲流量變化情況，但未考慮到襯砌及注漿圈的影響。朱成偉等[13]以保角映射法為基礎(chǔ)，結(jié)合地下水連續(xù)滲流理論及疊加原理，推導(dǎo)了水下雙線平行圓形隧道滲流場(chǎng)的解析解，獲得了水頭分布及隧道涌水量的解析表達(dá)式。趙建平等[14]結(jié)合雙極坐標(biāo)及鏡像法，推導(dǎo)了水下隧道滲流場(chǎng)解析解，并對(duì)支護(hù)參數(shù)取值問題進(jìn)行了合理分析。徐長節(jié)等[15]在滲流各向同性圓形隧洞滲流場(chǎng)解析解的基礎(chǔ)上，通過坐標(biāo)變換及保角映射法，推導(dǎo)了各向異性滲流下深埋圓形隧洞的滲流場(chǎng)解析解，但無限含水層豎井理論的假設(shè)與實(shí)際不符。

在隧道建設(shè)中，矩形、直墻曲拱形、多圓心形等具有高效空間利用率的隧道越來越常見?，F(xiàn)有關(guān)于水下隧道滲流場(chǎng)的研究常將隧道截面形狀假設(shè)為圓形，對(duì)于非圓形隧道，常采用等效面積法和等效圓周法將其轉(zhuǎn)化為圓形截面的隧道[16-17]。少數(shù)關(guān)于非圓形隧道滲流場(chǎng)的研究也未考慮到滲流各向異性的影響，關(guān)于聯(lián)合考慮滲流各向異性及非圓形隧道的研究未見報(bào)道。

基于此，在地下水滲流力學(xué)的基礎(chǔ)上，采用坐標(biāo)變換、保角映射、最優(yōu)化技術(shù)，求解穩(wěn)定滲流連續(xù)方程，給出了在滲透各向異性半無限多孔介質(zhì)中水下非圓形隧道的滲流量及襯砌外水頭的解析推導(dǎo)式。本文提供的方法適用于考慮洞周土滲透各向異性的矩形、直墻曲拱形、多圓心形等任意非圓形隧道的滲流場(chǎng)解析研究。

1 計(jì)算模型

本文研究對(duì)象為考慮滲透各向異性的水下非圓形隧道，即：水下隧道襯砌橫截面并非圓形，而是類似于馬蹄形。地層表面以上水頭高度為Hw；襯砌內(nèi)邊界形狀控制方程為f1=L(x,y), 襯砌外邊界形狀控制方程為f2=S(x,y) ，滲透系數(shù)為k1；隧道中心到地表距離定義為埋深h；洞周土水平向滲透系數(shù)為kx，豎直向滲透系數(shù)為ky，與滲透主軸方向一致，且均為常數(shù)?；鶞?zhǔn)面建立在地層表面；將地表邊界與襯砌外邊界之間的區(qū)域定義為區(qū)域1，將襯砌單獨(dú)定義為區(qū)域2。模型示意圖具體如圖1所示。

圖1 非圓形水下隧道模型示意圖Fig.1 Schematic diagram of non-circular underwater tunnel model

圖2 坐標(biāo)轉(zhuǎn)換示意圖Fig.2 Schematic diagram of coordinate transformation

圖3 區(qū)域1保角映射Fig.3 Region 1 conformal mapping

為方便推導(dǎo)，模型做如下基本假設(shè)：

1）洞周土體為滲透各向異性介質(zhì)；

2）隧道處于穩(wěn)定滲流狀態(tài)；

3）水流服從達(dá)西定律；

4）襯砌與土體接觸面為密實(shí)接觸界面，在縱向上滿足平面應(yīng)變條件；

5）排水通過襯砌實(shí)現(xiàn)均勻滲水。

按照等水頭滲水模式模擬隧道襯砌的排水方式[18]，假設(shè)襯砌內(nèi)水壓力為0，則隧道襯砌內(nèi)緣水頭為：

式中，H為滲流域某點(diǎn)水頭，h為隧道埋深。

2 解析解推導(dǎo)過程

由于區(qū)域1為滲透各向異性介質(zhì)，區(qū)域2襯砌支護(hù)為滲透各向同性介質(zhì)，問題比較復(fù)雜。為方便求解，將區(qū)域1與區(qū)域2分開計(jì)算，再按照不同區(qū)域滲流量相等原則聯(lián)立求解[19]。

2.1 區(qū)域1滲流場(chǎng)求解

2.1.1 區(qū)域1坐標(biāo)轉(zhuǎn)換

根據(jù)地下水滲流的連續(xù)性原理，可求得均質(zhì)、滲透各向異性含水層的2維穩(wěn)定滲流基本微分方程為：

通過坐標(biāo)變換將式（2）轉(zhuǎn)化為Laplace方程形式，以方便利用現(xiàn)存Laplace方程求解方法求解區(qū)域1滲流場(chǎng)。本文假設(shè)滲透系數(shù)kx、ky為常數(shù)，不發(fā)生變化，式（2）兩邊同除以ky得：

令

式中，X、Y為轉(zhuǎn)換后G平面的橫縱坐標(biāo)。

可得在G平面上的滲流連續(xù)方程：

由于內(nèi)邊界形狀不規(guī)則，在G平面中直接求解式（5）難以實(shí)現(xiàn)，可借助保角變換公式將區(qū)域1轉(zhuǎn)換為同心圓環(huán)進(jìn)行求解。

2.1.2 區(qū)域1保角映射

區(qū)域1主要由兩個(gè)邊界條件控制，可采用如式（6）所示的復(fù)變函數(shù)式[21]，該式可將區(qū)域1映射為同心圓環(huán)域，將地表邊界映射為半徑為1的外圓，內(nèi)邊界=s′(X,Y) 映 射為一半徑為 α（α ＜1）的內(nèi)圓。

式中：z為原平面點(diǎn)； ?為像平面點(diǎn)， ?=ξ+iη=ρeiθ，ξ、η 為復(fù)變函數(shù)w(?) 的實(shí)部和虛部，i為虛數(shù)單位，ρ為極徑， θ 為極角； μ 為實(shí)常數(shù)，k=1,2,3,···,μ；a、βk為待定實(shí)常數(shù)。當(dāng)給出具體的地表及內(nèi)邊界控制方程=s′(X,Y)時(shí)，也可以確定上述值，具體求解過程如下。

1）在地表邊界：將地層表面映射為半徑為1的圓。

首先，將坐標(biāo)改寫：

將式（7）代入式（6），分離實(shí)部和虛部可得：

通過式（8）可知，坐標(biāo)原點(diǎn)被映射為 ξ=-1，η=0 ；無窮遠(yuǎn)點(diǎn)映射為 ξ=1 ， η=0；地層表面其余各點(diǎn)也在單位圓上一一對(duì)應(yīng)。

2）在區(qū)域1內(nèi)邊界：將內(nèi)邊界映射為半徑為 α的圓。

將式（9）代入式（6），分離實(shí)部和虛部可得在內(nèi)邊界：

2.1.3 區(qū)域1滲流連續(xù)方程求解

映射后的圓環(huán)域內(nèi)外半徑分別為α、1，邊界條件不變， ?平面滿足滲流連續(xù)方程：

其通解形式為：

式中，H為水頭，C1、C2為待定常數(shù)，C3、C4為待定常數(shù)族，t為常數(shù)。代入邊界條件：

式中，h1為襯砌外緣水頭。

最終得區(qū)域1內(nèi)的水頭表達(dá)式：

式中，H1為區(qū)域1內(nèi)任一點(diǎn)水頭。

在已知洞周土域水頭分布的情況下，可以對(duì)其進(jìn)行積分得到區(qū)域1的滲流量Q1：

2.2 區(qū)域2滲流場(chǎng)求解

2.2.1 區(qū)域2保角映射

區(qū)域2為馬蹄形非圓截面襯砌（圖4），包括內(nèi)外緣兩個(gè)邊界，控制方程為：內(nèi)邊界f1=L(x,y)，外邊界f2=S(x,y)。由于是非圓形截面，直接求解困難，本文采用呂愛鐘等[22]提出的復(fù)變函數(shù)式（17）將區(qū)域2映射為同心圓環(huán)，其中，將襯砌內(nèi)緣映射為半徑為1的內(nèi)圓，將外緣映射為半徑為rs的外圓。

圖4 區(qū)域2保角映射Fig.4 Region 2 conformal mapping

保角映射公式一般以勞倫級(jí)數(shù)展開式表示：

式（17）確定后，區(qū)域1在式（17）下可映射為像平面的同心圓環(huán)域，Laplace方程在像平面依舊成立，可直接進(jìn)行求解。

2.2.2 區(qū)域2滲流連續(xù)方程求解

在 λ平面內(nèi)，Laplace方程形式為：

代入邊界條件：

解得襯砌水頭H2和滲流量Q2：

式中，H2、Q2分別為區(qū)域2的水頭和流量。

2.3 區(qū)域1、2聯(lián)立求解

按照層間滲流量相等原則Q1=Q2， 聯(lián)立式（16）和（21）解得襯砌外水頭h1和 隧道滲流量Q：

2.4 存在注漿圈時(shí)穩(wěn)定滲流方程的求解

當(dāng)注漿圈存在時(shí)，其求解思路與第2.2節(jié)區(qū)域2滲流場(chǎng)的求解思路一致，也即利用復(fù)變函數(shù)式（17）將注漿圈域轉(zhuǎn)換為同心圓環(huán)域，結(jié)合邊界條件求解滲流連續(xù)方程，再由滲流連續(xù)性原則與洞周土域、襯砌域聯(lián)立，不再重復(fù)贅述。

3 算例驗(yàn)證

3.1 模型設(shè)置

假設(shè)水下隧道斷面為橢圓形，襯砌內(nèi)緣長軸l1=4.8 m ，短軸l2=3.0 m ， 厚 0.2 m，埋深20.0 m，水位高于地表20.0 m。以滲透各向異性比（ky為常量）表示洞周土滲透各向異性，取值1、2、3；n越大，表示土體滲透各向異性越強(qiáng)。此外，設(shè)襯砌內(nèi)緣水壓力為0，主要參數(shù)見表1。

表1 水下隧道計(jì)算參數(shù)Tab.1 Calculation parameters of an underwater tunnel

表2 區(qū)域1參數(shù)計(jì)算值Tab.2 Area 1 parameter calculation values

3.2 數(shù)值模擬

采用COMSOL數(shù)值軟件進(jìn)行對(duì)比驗(yàn)證，參數(shù)取值與表1一致。為消除邊界條件的影響，模型寬度取300 m，高度取120 m；采用穩(wěn)態(tài)滲流模式，流體密度取ρ=1 000 kg/m3，孔隙率取0.3；水力邊界條件設(shè)模型底部及左右側(cè)為不透水邊界，上表面為自由透水邊界；此外，將有限元網(wǎng)格劃分設(shè)為加密模式，保證計(jì)算結(jié)果的收斂性，劃分2 971個(gè)有限元單元，模型如圖5所示。

3.3 解析公式參數(shù)計(jì)算

3.3.1 區(qū)域1參數(shù)計(jì)算

由本文水下隧道模型參數(shù)可知，區(qū)域1的內(nèi)邊界形狀為長軸為5.0 m，短軸為3.2 m的橢圓，埋深20.0 m，根據(jù)第2.2節(jié)可計(jì)算a、 α 、 βk。經(jīng)試算， βk取3項(xiàng)即可獲得較高的精確度。具體參數(shù)計(jì)算值見表2。

3.3.2 區(qū)域2參數(shù)計(jì)算

區(qū)域2是外邊界長短軸分別為5.0、3.2 m，內(nèi)邊界長短軸分別為4.8、3.0 m的橢圓形襯砌。通過式（17）將其轉(zhuǎn)換為同心圓環(huán)。其中，對(duì)于一般的隧道截面形狀只需要取3～4項(xiàng)即可滿足精度要求[24]。取v=4，代入計(jì)算程序得到：rs=1.051 3，R=3.999 8，Ck值見表3。

表3 v=4時(shí)Ck計(jì)算結(jié)果Tab.3 Ck calculation results when v=4

最終得到映射函數(shù)為：

經(jīng)過上述映射函數(shù)式，可將內(nèi)、外邊界映射為半徑為1.000 0、1.051 3 m的圓環(huán)，且可保證內(nèi)、外邊界獲得較高的精確度，最大誤差不超過0.263 8 m。

3.4 驗(yàn)證與分析

圖6為隧道滲流量隨襯砌滲透系數(shù)的變化曲線。由圖6可知，取有限項(xiàng)v的解析解會(huì)與數(shù)值模擬解產(chǎn)生一定的誤差，但整體上理論解和數(shù)值解變化趨勢(shì)一致，誤差很小，擬合度較高，證明了本文解的正確性。觀察圖6可知，隧道滲流量Q隨襯砌滲透系數(shù)的增大而增大，但增大趨勢(shì)漸漸放緩。此外，考慮洞周土滲透各向異性對(duì)Q也會(huì)產(chǎn)生較大的影響，且襯砌滲透性越好，這種影響越明顯。例如：在k1=1×10-7m/s條件下，當(dāng)n=1時(shí)，解析值Q為8.14 m3·d-1·m-1；當(dāng)n=2、3時(shí)，Q為12.67、16.00 m3·d-1·m-1，相比n=1分別增加 55.65%、96.56%。在k1=5×10-8m/s條件下，當(dāng)n=1時(shí)，解析值Q為6.86 m3·d-1·m-1；n=2、3時(shí)，Q為9.81、11.70 m3·d-1·m-1，相比n=1分別增加了43.00%、70.55%。從計(jì)算結(jié)果看，滲透各向異性對(duì)隧道滲流量的影響顯而易見。與此同時(shí)，當(dāng)襯砌滲透性很低時(shí)，n值對(duì)Q的影響大大降低，當(dāng)k1≤1×10-9m/s時(shí)，n值的影響已經(jīng)很小了。由此可見，是否考慮洞周土滲透各向異性的影響要根據(jù)水下隧道洞周土及支護(hù)結(jié)構(gòu)的滲透性具體分析，不能一概而論。

圖6 隧道滲流量隨襯砌滲透系數(shù)的變化曲線Fig.6 Variation curves of tunnel seepage flow with lining permeability coefficient

圖7為水下隧道滲流場(chǎng)水頭整體分布規(guī)律和隧道襯砌外水頭隨襯砌滲透性的變化規(guī)律。圖7（a）為當(dāng)n=1，k1= 1×10-7m/s時(shí)，水下隧道滲流場(chǎng)水頭分布。觀察圖7（a）發(fā)現(xiàn)，本文計(jì)算的水頭分布規(guī)律與數(shù)值模擬結(jié)果基本符合，驗(yàn)證了本文解析法的正確性。圖7（b）為襯砌外水頭隨襯砌滲透系數(shù)變化曲線，縱坐標(biāo)為襯砌外水頭與全水頭的比（h1/H全，H全為全水頭）。由圖7（b）可知，襯砌外水頭的解析解與數(shù)值解基本吻合，且襯砌外水頭隨著k1的變化規(guī)律與隧道滲流量隨k1的變化規(guī)律一致。當(dāng)襯砌滲透性很低時(shí)，襯砌外緣承擔(dān)的水頭很高，接近于全水頭值，這對(duì)襯砌的抗壓強(qiáng)度提出了很高的要求，必要時(shí)需進(jìn)行注漿加固以減小襯砌外水頭帶來的水壓力。

圖7 本文解析法與數(shù)值法水頭計(jì)算結(jié)果對(duì)比Fig.7 Comparison of analytical and numerical methods for hydraulic head calculation

滲透各向異性比n也會(huì)對(duì)襯砌外水頭產(chǎn)生一定的影響。在圖7（b）中，當(dāng)k1=5×10-8m/s時(shí)，n=1的襯砌外水頭為全水頭的32%，n=2、3的襯砌外水頭為全水頭的45%、54%?？梢姡舨豢紤]洞周土滲透各向異性的影響有可能低估襯砌外水頭值，將給襯砌的強(qiáng)度設(shè)計(jì)帶來隱患。

目前，在隧道工程中對(duì)于非圓形隧道常采用等效面積或者效周長法將非圓形隧道轉(zhuǎn)化為圓形隧道，等效圓半徑計(jì)算公式為[25]：

等效周長法：

等效面積法：

式中，r3、r4為等效周長法和等效面積法得到的隧道等效半徑，l、S為隧道周長和面積。

本文隧道的內(nèi)緣長短軸分別為4.8 m、3.0 m，計(jì)算周長為26.05 m，面積為45.24 m2，按照等效周長法和等效面積法計(jì)算得到r3=4.15 m、r4=3.80 m。同樣地，對(duì)于襯砌外緣，長短軸分別為5.0、3.2 m。按照等效周長法和等效面積法計(jì)算得到的等效圓半徑分別為r5=4.35 m、r6=4.00 m。滲流量計(jì)算公式采用文獻(xiàn)[25]中提供的公式，具體不再展開。

n=1時(shí)，使用本文方法與等效周長法和等效面積法所得隧道滲流量計(jì)算值隨襯砌滲透系數(shù)的變化規(guī)律如圖8所示。從圖8可以看出，等效周長法滲流量解整體上大于本文解，等效面積法滲流量解整體上小于本文解，本文解與數(shù)值解更吻合，說明本文方法計(jì)算的準(zhǔn)確性相比等效法更高。此外，當(dāng)襯砌滲透性增強(qiáng)時(shí)，等效法與本文方法計(jì)算的滲流量差異將擴(kuò)大，說明在襯砌滲透性較好或者無襯砌支護(hù)的毛洞隧道工況下采用本文方法求解隧道滲流量更有意義。

圖8 本文方法與等效法下滲流量對(duì)比Fig.8 Comparison of seepage flow between the proposed method and the equivalent methods

襯砌外水壓力的準(zhǔn)確計(jì)算對(duì)于支護(hù)結(jié)構(gòu)設(shè)計(jì)具有重要意義。按照式（22）計(jì)算出襯砌外水頭值后，可以按式（27）計(jì)算出襯砌外水壓力值：

式中，p為水壓力，H總為總水頭，y為位置水頭， γw為水的重度。

表4為本文方法與等效法在n=1（滲透各向同性）情況下，拱腰和拱頂點(diǎn)處襯砌外水壓力計(jì)算值。

表4 本文方法與等效法襯砌外水壓力比較Tab.4 Comparison of external water pressure between the proposed method and the equivalent method kРa

從表4可以看出，本文方法和等效法在拱腰點(diǎn)處襯砌水壓力差異較小，而在拱頂點(diǎn)差異較大。主要原因是等效法將非圓形隧道等效為圓形隧道，改變了襯砌上各點(diǎn)的位置坐標(biāo)，進(jìn)而改變了各點(diǎn)的位置水頭，導(dǎo)致襯砌外水壓力出現(xiàn)較大差異。

綜合來看，本文的計(jì)算方法在隧道滲流量、襯砌上各點(diǎn)水壓力的計(jì)算方面相對(duì)于等效法更準(zhǔn)確，更符合實(shí)際，在隧道排水系統(tǒng)設(shè)置、襯砌抗壓強(qiáng)度設(shè)計(jì)方面具有重要工程應(yīng)用價(jià)值。

4 結(jié) 論

1）通過坐標(biāo)變換及兩個(gè)保角映射復(fù)變函數(shù)式，推導(dǎo)了考慮洞周土滲透各向異性的水下非圓形隧道滲流場(chǎng)解析解，并通過連續(xù)性原則聯(lián)立求解隧道滲流量及襯砌外水頭。本文提供的方法適用于考慮洞周土滲透各向異性的任意形狀水下隧道滲流場(chǎng)的解析研究。

2）通過對(duì)某水下隧道進(jìn)行算例驗(yàn)證，利用COMSOL數(shù)值模擬軟件建立模型，與解析解值進(jìn)行對(duì)比，從隧道滲流量及襯砌外水頭兩方面驗(yàn)證了本文解的正確性。此外，通過與等效周長法和等效面積法計(jì)算的隧道滲流量和襯砌外水壓力對(duì)比發(fā)現(xiàn)，本文方法的計(jì)算值更準(zhǔn)確，更符合實(shí)際情況。

3）洞周土滲透各向異性對(duì)水下隧道滲流場(chǎng)有重要影響。通過對(duì)隧道滲流量及襯砌外水頭的分析可知，不考慮洞周土滲透各向異性將會(huì)低估上述兩者的值，給隧道工程的設(shè)計(jì)帶來安全隱患，因此有必要將洞周土滲透各向異性的影響納入隧道排水系統(tǒng)及襯砌強(qiáng)度設(shè)計(jì)之中。

下一步將分析注漿圈對(duì)滲透各向異性下非圓形隧道滲流場(chǎng)的影響。