寒区隧道渗流场和温度场耦合问题的三维非线性分析

第28卷 第9期 岩 土 工 程 学 报 Vol.28 No.9 2006年 9月 Chinese Journal of Geotechnical Engineering Sep., 2006

寒区隧道渗流场和温度场耦合问题的三维非线性分析

张学富

1，2

，喻文兵，刘志强

22

(1.重庆交通大学土木建筑学院，重庆 400074；2.中国科学院寒区旱区环境与工程研究所冻土工程国家重点实验室，甘肃 兰州 730000)

摘 要：根据传热学和渗流力学的基本理论，考虑冻土渗流场和温度场的耦合影响，建立冻土渗流场和温度场耦合问题的三维数学模型，应用Galerkin法进行有限元分析，并编制有限元分析程序。运用该数学模型和有限元分析程序，进行了青藏铁路昆仑山隧道DK977+578～DK977+682浅埋段的温度场以及渗流场和温度场耦合问题的三维非线性分析。分析结果表明：渗流场对隧道围岩的温度场影响较大，在寒区隧道工程设计时，应该考虑渗流场与温度场的耦合影响。

关键词：寒区隧道；渗流场；温度场；耦合问题

中图分类号：O241 文献标识码：A 文章编号：1000–4548(2006)09–1095–06

作者简介：张学富(1972– )，男，重庆合川人，博士，2004年毕业于中国科学院寒区旱区环境与工程研究所，现任副教授，主要从事桥梁结构、岩土与隧道工程以及寒区岩土工程数值计算等方面的科研教学工作。

Three-dimensional nonlinear analysis for coupled problem of seepage field and

temperature field of cold regions tunnels

ZHANG Xue-fu1,2, YU Wen-bing2, LIU Zhi-qiang1

(1. School of Civil Engineering and Architecture, Chongqing Jiaotong University, Chongqing 400074, China；2. State Key Laboratory of Frozen Soil Engineering, Cold & Arid Regions Environmental and Engineering Research Institute, Chinese Academy of Sciences, Lanzhou

730000, China)

Abstract: According to the essential theories of heat transfer and seepage, considering the coupled effects of seepage field and temperature field, three-dimensional computational model of the coupled problem was given, the finite element formulae of this problem were obtained by Galerkin's method, and the computer program of the finite element was written. Using this model and the program, three-dimensional nonlinear analysis for the coupled problem of seepage field and temperature field of the flat excavated part from DK977+578～DK977+682 in Kunlun Mountain tunnel on Qinghai-Tibet Railway were made. It was shown that the effect of seepage field on temperature field in the tunnels in cold regions was very large, and seepage field should be taken into account during the design.

Key words: tunnels in cold region; seepage field; temperature field; coupled problem

0 引 言

寒区隧道与其它地区的隧道不同，就是必须分析隧道围岩的温度场，这往往是修建寒区隧道成功与否的关键。寒区隧道由于冻害的影响，有的衬砌裂缝宽达5 cm，有的常年有8～9个月不能使用，严重影响正常交通。另外，在加拿大、俄罗斯和中国的东北及西北有许多寒区隧道，有的隧道由于渗流的影响，隧道围岩得不到及时冻结而出现隧道漏水、挂冰以及涎流冰等病害[1]。因此，为了解决寒区隧道这些普遍存在的冻害问题，进行寒区隧道渗流场和温度场耦合问题的分析研究是非常必要的。就此问题国内外学者进行过一些研究。在国外，Banacina et al为了解决带相变的热传导问题提出了一种数值解法[2]；Comini et al

对相变热传导温度场非线性问题进行有限元分析[3]；Bathe提出了没有单元网格迭代的渗流自由面的有限元计算法[4]；Zhao和Valliappan用瞬态无界元对无限介质中的渗流问题进行了研究[5]。在国内，赖远明等既进行了寒区隧道渗流场和温度场耦合问题的非线性分析，又进行了寒区隧道渗流场、温度场和应力场耦

───────

基金项目：重庆交通学院启动基金资助项目；交通部公路科学研究所预研基金资助项目；国家杰出青年科学基金项目(40225001)；国家自然科学基金资助项目(40401010)；中国科学院知识创新工程重要方向项目(KZCX3-SW-351)；国家基础学科人才培养资助项目(J0130084)；中国科学院知识创新工程重大资助项目(INF105-SDB-1-21)；国家科技攻关项目“西部开发科技行动”重大资助项目(2004BA901A02) 收稿日期: 2005–06–27

1096 岩 土 工 程 学 报 2006年

合问题的非线性分析[6-7]，但他们的研究只是针对二维平面问题；张有天提出了具有自由面的渗流分析的初流量法[8]。综合国内外的研究状况，到目前为止，还没有人进行寒区隧道渗流场和温度场耦合问题的三维非线性分析研究。而不论渗流场还是温度场，隧道轴向剖面对隧道横截面是有影响的，因此进行寒区隧道渗流场和温度场耦合问题的三维分析研究是非常必要的。

本文就此问题，根据热传导和渗流力学的基本理论，建立了渗流场和温度场耦合问题的三维数学模型，应用Galerkin法进行了数值分析，并编制了有限元分析程序。然后，运用该模型和有限元分析程序，对青藏铁路昆仑山隧道的DK977+578～DK977+682浅埋段的温度场以及渗流场和温度场耦合问题进行了三维非线性分析研究。

Su

Hu u Hu u Hu

= Kx+ Ky +

yy t x x

u Hu

+Q ， (2b) Kz

z z

1 渗流场和温度场耦合问题的三维数

学模型及其有限元公式

本文所涉及的寒区隧道是指地处于年平均气温低于0 ℃的隧道，寒区隧道衬砌表面和地表年年发生冻融循环，这样一来，本文讨论的问题便成了带相变的热传导问题。考虑渗流对温度场的影响，根据热传导和渗流力学的基本理论，可导出温度场和渗流场耦合问题的三维数学模型。

在正冻区 f内有 Tf Tf Tf Tf C= f λf + λf + λf txxyyzz

f Tff Tff Tf Cwρw Vx x+Vy z+Vz z ，

S Hf= Kf Hf + Kf Hf +

y x

f t x x y y f Hf (1) Ky +Q，

zz

f Hf

Vxf= Kx，

x

f

f Hf

， Vy= Ky

y

Hf

。 Vzf= Kzf

z

Hu

， x Huu ，Vyu= Ky (2c) y Hu

。Vzu= Kzu z

式中 符号f和u分别表示冻结和融化状态；Tf、Cf和λf分别为正冻区 f内寒区隧道岩土的温度、体积

ff、Ky和Kzf分别为正热容和导热系数；Hf、Sf、Kx

冻区的水头、给水度、x方向、y方向和z方向的渗透系数；Q为渗流场的源或汇；Cw和ρw为水的比热和容重；带有符号“u”的符号分别为融区的相应物理量。

在移动界面s(t)上，必须满足连续条件和能量守恒条件，即

Tf(s(t),t)=Tu(s(t),t)=Tm ， (3)

Vxu= Kxu

ds(t) T Tf

λuu=L ， (4)

dt n n

式中 L为含水岩土的相变潜热；n为移动边界s(t)的方向矢量；Tm为发生相变时的温度。

固定边界上的边界条件为

λf

T=Tw， λ T=α(T Ta) ， (5)

n1

式中 α为对流换热系数；Ta为环境温度；n1为固定边界的方向矢量。

初始条件为

Tut=0=T0， ut=0=0 。 (6) 渗流场内水头H还必须满足以下两类边界条件： (1) 第一类边界条件，指水头边界

H=H0 。 (7)

(2) 第二类边界条件，指流量边界

H HKxcos(n2,x)+Kycos(n2,y)+

x y HKzcos(n2,z)= q0 ， (8)

z

式中 n2为渗流边界的方向矢量。

式（1）～（8）便是渗流场和温度场耦合问题的三维数学模型、初始条件及其边界条件。由于岩土的比热和导热系数是随着温度的变化而变化，加上相变s(t)界面的位置也不是固定的，因此，s(t)界面的能量守恒条件是非线性的。此外渗流场自由面也是待确定，故本问题在数学上是一个强非线性问题。

在未冻区 u内有 T T T T

Cuu= λuu + λuu + λuu

t x x y y z z

T T T

Cwρw Vxuu+Vyuu+Vzuu ， (2a)

x y z

第9期 张学富，等. 寒区隧道渗流场和温度场耦合问题的三维非线性分析

1097

用显热容法，假设相变是发生在Tm附近的一个温度范围(Tm± T)内。在构造等价热容时还必须考虑 T的影响，构造的热容表达式和导热系数的表达式分别为[9-10]

T<Tm T， Cf,

C+Cu L

,Tm T≤T≤Tm+ T，C*= +f (9)

22T

T>Tm+ T， Cu,和

T<Tm T， λf,

λu λf λ+ T (Tm T) ， (10) λ*= f2 T Tm T≤T≤Tm+ T，

T>Tm+ T。 λu,

由于此问题在数学上是一个强非线性问题，无法

获得解析解，故只能采用数值计算方法获得数值解。本文采用Galerkin法，可得该问题的有限元计算公式为

T [M] +[K]{T}={F} ， (18)

t H [N] +[G]{H}={P} ， (19)

t 式中 Kij=∑∫

Ni Nj Ni Nj Ni Nj

+ + d + e

xxyyzz

Nj Nj Nj

CV+V+Vρ∑∫ eww x xy yz z Nid +

λ

∑∫

Γe2

αNiNjdΓ ， (20)

Mij=∑∫

e

C NiNjd ， (21)

对于渗透系数也可以写成一个统一的表达式

Wu β

Kxu，T≤Tm， Kx= W (11)

Ku， T>T，

mx

Wu β

u

，T≤Tm， Ky (12) Ky= W

Ku， T>T，

ym

Wu β

Kzu，T≤Tm， Kz= W (13)

Ku， T>T，

zm

式中 W和Wu分别为寒区隧道岩土的含水量和未冻

水含量；β为常数，由试验确定。

应用式（9）～（13），则式（1）～（4）可化简为

T T T T C = λ + λ + λ t x x y y z z T T T Cwρw Vx +Vy +Vz ， (14)

xyz

H H H H S = Kx + Kz + Ky +Q， t x x y y z z

(15)

式中

T≤Tm， S，

S = f (16)

S，T>T。um Vx 、Vy 和Vz 为x、y和z方向的渗流速度，其

表达式为

H H HVx = Kx，Vy = Ky，Vz = Kz 。(17)

x z y

Fi=∑∫eαTaNidΓ ， (22)

Γ2

Nj Ni Nj Ni

+Ky + Gij=∑∫e(K

x x y y

N Nj

)d ， (23) Kz i

z z

Nij=∑∫eS NiNjd ， (24)

x

Pi=∑∫

e

QNid 。 (25)

先确定渗流问题的自由面，然后应用Crank-Nicolson格式对式（18）、（19）进行逐时段求解，便得本问题的数值结果。

2 算例分析

2.1 分析模型、参数和条件

青藏铁路是国家重大工程建设项目，昆仑山隧道是青藏铁路的关键控制工程。昆仑山隧道的工程概况和地质纵断面图如文献[9]所述。在2003年5月以后，昆仑山隧道出现漏水现象，其漏水的主要是由于浅埋段冲沟的地表水渗入隧道所致，所以本文选择昆仑山隧道DK977+578～DK977+682浅埋段作为研究对象。其计算域的范围如图

1和图2所示。

图1 计算域的纵向面图

Fig. 1 The longitudinal section of the computation domain

1098 岩 土 工 程 学 报 2006年

2π 5.2 12sin th+arcsin 。

(26)

876012

图2 计算域的横剖面图

Fig. 2 The cross-section of the computation domain

昆仑山隧道DK977+578～DK977+682浅埋段位

于Ⅳ级和Ⅵ级围岩中。而隧道围岩的级别不同，其隧道衬砌尺寸也有所不同。但是，隧道衬砌尺寸差别不是很大，而且与温度场有主要关系的保温材料厚度是一样的，所以为了简化模型，在计算分析中隧道衬砌结构按Ⅵ级围岩中的衬砌结构来处理，如图3所示。

图4 单元网格划分

Fig. 4 The calculation model and element meshes

隧道衬砌表面按中铁西南研究院的实测资料计算。而隧道的山顶面AIJH是阴面，太阳辐射较少，并考虑植被的影响，根据中国科学院寒区旱区环境与工程研究所吴紫汪研究员的研究成果应为

2.6

Ta= 5.2 1.5+×th+

50×8760 2π 5.2 12sin th+arcsin 。 (27)

876012 由于两边山峰的影响，GHJK沟底表面的太阳辐

来计算。而KGFE射同样较少，所以同样按照式（26）

是阳面，太阳辐射较多，并考虑植被的影响，同样根据吴紫旺研究员的研究成果应为

2.6

×th+ Ta= 5.2+1.5+

50×8760 2π 5.2 12sin th+arcsin 。(28)

876012 空气与围岩的对流换热系数α＝15.0

(W/m2.k)[11]。根据勘测资料，该地区的平均地温为－2.5 ℃。然后冻土中按2.5%的地热梯度增加可得到围岩的初始温度。根据《铁路隧道设计规范》[12]和室内试验，可得隧道围岩、衬砌混凝土、抛石以及保温材料的热学参数和渗透系数，见表1所示。表中，Kxu、u

和Kzu相等。 Ky

昆仑山隧道的排水是按照“防、排、截和堵，多道防线，综合治理”，以防和堵为主，排和截为辅的原

图3 昆仑山隧道浅埋段衬砌结构图

Fig. 3 The lining structures of Kunlun Mountain tunnel

计算时，采用8节点六面体块体单元，整个模型

共划分43 426个单元，其中围岩划分33 694个，内衬划分2 744个，保温材料划分1 372个，外衬划分4 116个，其单元网格划分如图4所示。

由于模型的边界较远，所以ABCI、DEF、ABHGF和CDEKJI面为绝热边界。计算模型底面BCD的地热热流密度q=0.06W/m2。考虑气候变暖，该地气温变化曲线为

2.6

Ta= 5.2+×th+

50×8760

表1 隧道围岩、衬砌混凝土、抛石以及保温材料的热学参数和渗透系数

Table 1 Thermal and seepage parameters of the insulation layer, concrete, ripped rock and the rock surrounding the tunnel

物理量 围岩 混凝土 碎石 保温 材料

w/% 13.5 1.0 8.1 1.2

Kxu/(m·s-1)

5×10-6 2×10-7 1×104

-

λu/(W·m-1·℃-1)

3.00 1.00 0.36

Cu/(J·m-3·℃-1)

2.52×106 1.49×106 1.01×106

λf/(W·m-1·℃-1) Cf/(J·m-3·℃-1)

3.50 1.96 1.00 1.49 0.69 1.01

1×10-10 0.03 0 0.03 0

第9期 张学富，等. 寒区隧道渗流场和温度场耦合问题的三维非线性分析 1099

则进行设计的。在设计中，在隧道衬砌之间的保温材料内外都铺设了防水层，隧道墙脚设双侧保温沟，墙

通过“三通”及φ50 脚纵向设φ100 mm的PVC盲沟，

mm的PVC泄水管与隧道内保温侧沟连接，侧沟水通过检查井与格尔木端洞外保温暗沟连接引排。因此，

CDEKJI面以及渗流计算时的边界条件为ABFGH面、

隧道衬砌外表面为不渗水边界，AI线上的水头为H0＝56 m，EF线上的水头为H1＝68 m，GH线上给定一个水源41.67 m3/d，JK线上给定一个汇25 m3/d，在隧道围岩墙脚线上也给定一个汇0.96 m3/d。 2.2 计算结果分析

将上述的各种热学参数、边界条件和初始条件，代入有限元公式中，进行有限元计算。本文首先进行了昆仑山隧道DK977+578～DK 977+682浅埋段的渗流场和温度场耦合问题的三维非线性分析；然后，为了研究渗流对温度场的耦合影响，还进行了该段温度场的三维非线性分析。

由于此计算模型是三维空间的，如果将计算结果取成三维空间图，则不易看出隧道衬砌附近围岩的冻融状态，所以仅绘出了两个具有代表意义的剖面图：一个是沿隧道中心线的纵剖面，一个是沿隧道浅埋段冲沟中心的横剖面，即 DK977+618处的横剖面。两

种条件下的计算结果如图5～6所示，具体分析如下。

岩出现了极少季节性融化；在冲沟下的土体、山体阴

坡下端土体以及山体阳坡的土体都存在季节性融化；在冲沟处围岩已经没有冻土层，而围岩全部融化，这时如果冲沟中有地表水通过，则地表水就会通过融化层到达隧道壁后。从计算结果来分析，如果昆仑山隧道在浅埋段还会出现渗漏水现象。

图6 当不考虑渗流场和温度场的耦合影响时，第50年的9月

30日的温度分布图

Fig. 6 The temperature distributions of the tunnel without considering the coupled effects of seepage field and

temperature field on Sept. 30 in the 50th year

图5 当考虑渗流场和温度场的耦合影响时第50年的9月30

日的温度分布图

Fig. 5 The temperature distributions of the tunnel after considering the coupled effects of seepage field and temperature field on

Sept. 30 in the 50th year

由图6可以看出，在考虑气候变暖即50 a升高2.6℃的条件下，昆仑山隧道在不考虑渗流场和温度场的耦合影响条件下，在第50年的9月30日，浅埋段衬砌壁后围岩出现了极少季节性融化，格尔木端也出现了极少季节性融化；在冲沟下的土体、山体阴坡下端土体以及山体阳坡的土体都存在季节性融化；在冲沟处围岩还有1.2 m左右厚的冻土层。因为冻土是很好的隔水层，所以这1.2 m左右厚的冻土层可以阻碍冲沟中的地表水渗入隧道衬砌壁厚。从计算结果来看：如果昆仑山隧道将不会出现渗漏水现象。

比较以上两种计算结果可以看出，在不考虑渗流场和温度场的耦合影响条件下，昆仑山隧道不会出现渗漏水现象，计算结果偏于安全。在考虑渗流场和温度场的耦合影响条件下，昆仑山隧道还会出现渗漏水现象，这与现场实际情况符合。

3 结 论

(1) 昆仑山隧道如果不采取防水措施，将还会出现漏水现象。

(2) 渗流场对隧道围岩的温度场影响较大，在寒区隧道工程设计时，应该考虑渗流场与温度场的耦合

由图5可以看出，在考虑气候变暖即50 a升高2.6℃的条件下，昆仑山隧道在考虑渗流场和温度场的耦合影响条件下，在第50年的9月30日，浅埋段衬砌壁后围岩出现了极少季节性融化，格尔木端仰拱底围

1100 岩 土 工 程 学 报 2006年

regions tunnels [J]. Cold Regions Science and Technology, 1999,29(1):89–96.

[7] LAI Y M, WU Z, ZHU Y, ZHU L. Nonlinear analysis for the

coupled problem of temperature, seepage and stress fields in cold-region tunnels [J]. Tunnelling and Underground Space Technology, 1998:435–440.

[8] 张有天,陈 平,王 镭.有自由面渗流问题分析的初流量法

[J].水利学报,1988,12(1):18–26.(ZHANG You-tian, CHEN Ping, WANG Lei. The original flux method for the seepage analysis with free surface [J]. Journal of Hydraulic Engineering, 1988,12(1):18–26.)

[9] ZHANG Xue-fu, LAI Yuan-ming, YU Wen-bing, et al.

Forecast analysis for the re-frozen of Kunlunshan permafrost tunnel on Qing-Tibet Railway in China[J]. Cold Regions Science and Technology, 2004,38(1):12–22.

[10] ZHANG Xue-fu, LAI Yuan-ming, YU Wen-bing, ZHANG

Shu-juan. Nonlinear analysis for the three-dimensional temperature fields in cold regions tunnels [J]. Cold Regions Science and Technology, 2002,35(3):123–221.

[11] ZHANG Z R. Heat transfer[M]. Beijing: Higher Educational

Press, 1989:133–134.

[12] 铁道部第二勘测设计院.铁路隧道设计规范[M].北京:中国

铁道出版社,1999:11–16.(The First Survey & Design Institute of China Railways. Code for design of railway tunnel [M]. Beijing: China Railway Publishing House, 1999:11–16.)

影响。

(3) 本文提供了寒区隧道渗流场和温度场耦合问题的三维计算模型和有限元方法，能够为寒区隧道工程的设计提供理论计算依据。 参考文献：

[1] LAI Yuan-ming, LIU Song-yu, WU Zi-wang, YU Wen-bing.

Approximate analytical solution for temperature fields in cold regions circular tunnels [J]. Cold Regions Science and Technology, 2002,34(1):43–49.

[2] BONACINA C, COMINI G, FASANO A, et al. Numerical

solution of phase-change problems[J]. Int J Heat Mass Transfer, 1973,16(6):1852–1832.

[3] COMINI G, GUIDICE S del, LEWIS R W, et al. Finite element

solution of nonlinear heat conduction problems with special reference to phase change[J]. Inter J for Numerical Methods in Engineering, 1974,8(6):613–624.

[4] BATHE K J, KHOSHGOFTAAR M R. Finite element free

surface seepage analysis without mesh iteration [J]. Inter J Numerical Analytical Methods in Geomechanics, 1979, 3(1):13–22.

[5] ZHAO C B, VALLIAPPAN S. Transient infinite elements for

seepage problems in infinite media [J]. Inter J Numerical Analytical Methods in Geomechanics, 1993,17:323–341. [6] LAI Y M, WU Z, ZHU Y, ZHU L. Nonlinear analysis for the

coupled problem of temperature and seepage fields in cold

岩土工程博士学位介绍启事

本刊今年第12期继续开辟专门栏目，报导我国岩土工程界博士学位的概况，内容包括博士学位获得者的简况，师从导师与学位授予学校和博士论文概要等。已在海内外获得岩土工程博士学位的我国学者均属报导范围。欢迎岩土工程博士学位

获得者向本栏目投来文稿，来稿请控制在1000字以内，在寄来书面文稿的同时发来文稿的电子文件，编辑部信箱ge@。文稿刊登不收版面费，不付稿费。

（本刊编辑部）

本文来源：https://www.bwwdw.com/article/cd7e.html