地下工程渗流场和应力场耦合研究

∫ ∫ ∫ δ εT dσdV − δ uT dfdV − δ uT dtdS = 0
V
V
S
式中， t 为面力， f 为体力；δ ε 、δ u 分别为虚位移和虚应变。
以增量形式表示的本构关系为：
（6-19）
dσ′ = Dep (dε − dεl )
（6-20）
式中， Dep 为弹塑性矩阵， dεl 为孔隙流体压力引起的颗粒压缩，仅在正方向上产生变形，无剪方
z 是某指定参考面之上的高度， ρw 为流体密度， gradH 为水力梯度。
可以发现，式（6-1）表示的渗流速度是一种假想流速，即认为过水断面中任意一点均存在着水 流（包括固体骨架上也存在着水流），实际上，过水断面有两部分，一是孔隙空间，另一部分是固体
骨架。但是，水流仅能通过孔隙空间流动，可供水流通过的面积为 nA （ n 为孔隙率）；有时孔隙空
岩土体的流-固耦合是水体流动和介质变形相互作用、相互影响的结果，完全流-固耦合模型的 相互作用机理如图 6-1 所示，由变形和孔隙水相互作用而产生可称为“直接耦合”，如图中的过程Ⅰ 和Ⅱ；而“间接耦合”是孔隙度的改变引起渗透系数的变化，孔隙度和渗透系数的变化由有效应力
的变化而引起，孔隙的减小会引起介质截面积和所含水量的减小，进一步使得材料刚度增大，如图 中的Ⅲ和Ⅳ，间接耦合使得耦合系统呈非线性，导致多孔介质的渗透系数呈各向异性。
为是 Forchheimer 定律的特例。 Darcy 定律用于表述为层流条件下通过多孔介质的渗流速度与水力梯度满足线性关系，在一维
条件下有：
v′ = Q = −kgradH = kJ A
（6-1）
式中，v′ 为平均渗流速度；Q 为流量，A 为过水面积；k 为渗透系数；H 为测压水头，H = z + pw ， gρw
pw 与 pa 分别为孔隙水压力和孔隙气压力， χ 与饱和度和表面张力有关，由于试验数据较难获
得，一般可假定 χ = s ，显然 sw + sa = 1 。于是，有效应力的表达式又可表示为：
( ) σ i′j = σ ij + δij ⎡⎣sw pw + 1− sw pa ⎤⎦ = σ ij + δij p
独立的渗透系数为 3 个。
在多孔介质渗流中，雷诺数 Re 可表达为：
Re
=
vd ν
（6-7）
式中， d 为多孔介质中固相颗粒的平均粒径； v 为渗流速度；ν 为流体的运动粘滞系数。
Darcy 定律适合在低雷诺数（ Re < 1 ~ 10 ）条件下流体通过多孔介质的层流运动。当雷诺数大
于 1~10 之间的某个值时，而应当采用非线性的 Forchheimer 定律，一维条件下，水力梯度和渗流速 度满足如下非线性关系：
6.1 岩土介质渗流-应力耦合理论
流固耦合理论是渗流力学与固体力学交叉而生成的一个力学分支，它是研究地质环境中流体与 岩体相互作用的一门科学，其研究与应用已涉及到了水力水电工程的渗流与控制、水库诱发地震、 有害核废料处理、煤矿瓦斯泄漏及石油开发等领域。近年来针对岩石低渗透特性的研究已成为国内 外岩土工程界关注的热点之一，我国实施与规划中的石油/天然气地下能源储存、低渗透油气田开发、 高瓦斯矿井瓦斯抽放、放射性废料地质深埋处置等工程，都涉及到在复杂的地质结构中建造地下工 程。
度 v 和水力梯度 gradH 是矢量，若对于各向同性介质而言，k 是一个标量，对于各向异性介质来说，
k 是一个二阶张量。 在非饱和条件下，需考虑饱和度 s 的影响，Darcy 定律可表达式为：
各向异性渗透介质的达西定律可表示为：
snev = −kgradH = kJ
（6-4）
sne vi = −kij J j
+
m
dp 3KS
⎞ ⎟ dV ⎠
−
δ εT mdpdV −
V
δ uT dfdV −
V
δ uT dtdS = 0
S
（6-23）
由于渗流连续性方程含有时间项，为了将应力和渗流进行耦合，需对虚功方程（6-23）进行时
间的求导，具体表示式为：
∫ ∫ ∫ ∫ V
δ
εT
Dep
⎛ ⎜ ⎝
dε dt
+
m
1 3KS
dp dt
σi′j = σ ij + αδij ⎡⎣χ pw + (1− χ ) pa ⎤⎦
（6-15）
式中，σ i′j 为有效应力，σ ij 为总应力， δij 为 Kronecker 符号，α 为 Biot 系数，定义为：
α = 1− KV KS
（6-16）
式中， KV 、 KS 分别为岩石的体积压缩模量和固体颗粒的压缩模量。通常α 接近 1，取α = 1 。
按照渗透性是否随渗流方向变化，可将岩土介质分为各向同性介质和各向异性介质，而将渗透
性是否随空间位置的改变而变化，又可将其区分为均质和非均质两类。渗透张量各渗透系数 kij 随坐
标轴的转动而变化，其转换关系为：
krs = αriαsj kij
（6-10）
式中，αri 为 r 轴与 i 轴夹角的余弦。该式即为由 ij 坐标系转换为 rs 坐标系的渗透系数转化关系式。
⎞ ⎟ ⎠
dV
−
δ εT m dp dV −
V
dt
δ uT df dV − V dt
δ uT dt dS = 0 S dt
（24）
由于假定气压在整个非饱和区域是恒定的，即 dpa = 0 ，故平均孔隙流体压力对时间的导数可 dt
简化为：
( ) dp = d
dt
sw pw + (1− sw ) pa
间中一部分为死端孔隙，水流不能通过，可用有效孔隙率 ne 来描述，可供水流通过的面积为 ne A 。
因此，实际渗流速度 v 可表述为：
v = Q = v′ ne A ne
（6-2）
饱和条件下，考虑有效孔隙率 ne 情况下的 Darcy 定律表达式为：
nev = −kgradH = kJ
（6-3）
在一维条件下，渗透速度 v 、水力梯度 gradH 和渗透系数 k 是一个标量；但在三维条件下，渗透速
[ ] 根据张量坐标系统的变换法则，总有三个彼此正交的方向建立坐标系统，使得渗透系数 K 变
成主对角阵：
⎡k1
[K
]
=
⎢ ⎢
0
⎢⎣ 0
0 k2 0
0⎤
0 k3
⎥ ⎥ ⎥⎦
，
[
K
]
=
⎡k1
⎢ ⎣
0
0⎤
k2
⎥ ⎦
主对角阵上的元素称为主渗透系数，随对应的方向为渗透主方向。
任意渗透方向 vr 上渗透系数可表示为：
6.1.1 渗流-应力相互耦合的力学机理
在诸多岩土工程中，如石油天然气开采、地下水抽放、地热开采、核废料处理、煤层瓦斯突出、 水库诱发地震等渗流问题，均涉及到流固耦合渗流问题，此项研究已成为目前科学研究的热门课题。 传统的流-固耦合模拟中，孔隙度和渗透率是保持不变的，而在实际的渗流过程中，由于孔隙流体压 力的变化，一方面要引起多孔介质骨架有效应力变化，由此导致渗透率、孔隙度等的变化；另一方 面，这些变化又反过来影响孔隙流体的流动和压力的分布。因此，需考虑孔隙流体在多孔介质中的 流动规律及其对多孔介质本身的变形或者强度造成的影响，即考虑多孔介质内应力场与渗流场之间 的相互耦合作用。
（6-11）
1 = cos2 β1 + cos2 β2 + cos2 β3
kv
k1
k2
k3
（6-12）
式中，β1 、β2 、β3 分别为渗透方向 vr 与渗透主坐标轴 1、2、3 的夹角；cos β1 、cos β2 、cos β3 分
别为渗透方向 vr 的方向余弦。渗透系数 kv 的几何意义如图 6-2 所示，渗透系数张量对应的几何为一
应力场 应力、应变
(Ⅳ)
(Ⅱ)
力学参数
孔隙体积的改变
水力参数
(Ⅰ)
(Ⅲ)
渗流场 孔隙水压力、流量
图 6-1 完全流-固耦合系统相互作用的力学机理
6.1.2 多孔介质中流体渗流规律
孔隙流体的渗流行为遵循 Darcy 定律或 Forchheimer 定律，Darcy 定律一般适用于低渗流流速，
是线性关系；而 Forchheimer 定律是非线性定律，主要模拟更高流动速度的情况，Darcy 定律可以认
向上的变形，具体表达式为：
其中， m = [1,1,1, 0, 0, 0]T 。
dεl
=
−m
dp 3KS
（6-21）
式（6-17）可表示为：
σ′ = σ +αmp
（6-22）
取α = 1 ，结合式（6-22）和式（6-20），式（6-19）可表示为：
∫ ∫ ∫ ∫ V
δ εT Dep
⎛ ⎜ dLeabharlann ⎝⎧⎪kx x ⎪=
kx
+ kz 2
−
kx
− kz 2
cos 2θ
⎪⎨kz z ⎪
=
kx
+ kz 2
+
kx
− kz 2
cos 2θ
⎪⎪⎩kx z
=
kx
− kz 2
sin 2θ
（6-14）
z
kx z
θ
kz
kv
kz z
x
kx x
kx
图 6-2
渗透椭圆及渗透圆示意图
6.1.3 孔隙介质的有效应力原理
岩石是由固体骨架和相互连通的孔隙以及储存于骨架孔隙中的流体（油、气和水）三者组成的 多孔介质。岩石中的流体能够承担或传递压力，将其定义为孔隙压力，而通过岩石颗粒间的接触面 传递的应力称为有效应力。目前一般采用 Biot 有效应力，本文定义应力以拉力为正，孔压在饱和区 以压应力为正，在非饱和区域孔压为负。Biot 有效应力的表达式为：