目 录
第一章 渗流理论基础 .............................................................................................................................. 1
第一章 渗流理论基础
1.1 渗流的基本概念
1.1.1 多孔介质及其特性
1.1.1.1多孔介质的概念
多孔介质(Porous medium):地下水动力学中具有空隙的岩石。广义上包括孔隙介质、裂隙介质和岩溶不十分发育的由石灰岩和白云岩组成的介质,统称为多孔介质。 孔隙介质:含有孔隙的岩层,砂层、疏松砂岩等;
裂隙介质:含有裂隙的岩层,裂隙发育的花岗岩、石灰岩等。
1.1.1.2 多孔介质的性质
(1) 孔隙性:有效孔隙和死端孔隙。
孔隙度(Porosity )是多孔介质中孔隙体积与多孔介质总体积之比(符号为n ),可表示为小数或百分数,n=Vv/V。
有效孔隙(Effective pores )是多孔介质中相互连通的、不为结合水所占据的那一部分孔隙。
有效孔隙度(Effective Porosity)是多孔介质中有效孔隙体积与多孔介质总体积之比(符号为n e ),可表示为小数或百分数,n e =Ve /V。
死端孔隙(Dead-end pores )是多孔介质中一端与其它孔隙连通、另一端是封闭的孔隙。
(2) 连通性:封闭和畅通,有效和无效。
(3) 压缩性:固体颗粒和孔隙的压缩系数推导。
(4) 多相性:固、液、气三相可共存。其中固相的成为骨架,气相主要分布在非饱和带中,液相的地下水可以吸着水、薄膜水、毛管水和重力水等形式存在。
固相—骨架matrix
气相—空气,非饱和带中
液相—水:吸着水 Hygroscopic water
薄膜水 pellicular water
毛管水 capillary water
重力水 gravitational water
1.1.1.3多孔介质中的地下水运动
比较复杂,包括两大类,运动特点各不相同,分别满足于孔隙水和裂隙岩溶水的特点。
(1)第一类为地下水在多孔介质的孔隙或遍布于介质中的裂隙运动,具有统一的流场,运动方向基本一致;
(2)另一类为地下水沿大裂隙和管道的运动,方向没有规律,分属不同的地下水流动系统。
1.1.2渗透与渗流
1.1.2.1渗透
渗透是地下水在岩石空隙或多孔介质中的运动,这种运动是在弯曲的通道中,运动轨迹在各点处不等。为了研究地下水的整体运动特征,引入渗流的概念。
图1-2 岩石中的渗流
(a )实际渗透 (b)假想渗流
1.1.2.2渗流
渗流(seepage flow ):具有实际水流的运动特点(流量、水头、压力、渗透阻力),并连续充满整个含水层空间的一种虚拟水流;是用以代替真实地下水流的一种假想水流。其特点是:
(1)假想水流的性质与真实地下水流相同;
(2)充满含水层空隙空间和岩石颗粒所占据的空间;
(3)运动时所受的阻力与实际水流所受阻力相等;
(4)通过任一断面的流量及任一点的压力或水头与实际水流相同。
渗流场(flow domain):假想水流所占据的空间区域,包括空隙和岩石颗粒所占的全部空间。
1.1.2.3典型单元体
典型单元体(REV ,Representative Elementary V olume )又称代表性单元体,是渗流场中其物理量的平均值能够近似代替整个渗流场的特征值的代表性单元体积。
REV 具备两个性质:
(1) 其体积和面积,大于个别空隙而小于渗流场,其中的渗流可以从一点连续运动到另一点;
(2) 通过单元体的运动要素(流量Q 、水头h 、压力p 、实际水头受到的阻力R )与真实水流相等,运动要素是连续变化的。
REV 的作用:
(1) 把物理性质看作是坐标的函数,孔隙度n 、导水系数T 、给水度μ和渗透系数均连续。
(2) 渗流的要素可以微分、积分,可以用微分方程来描述渗流要素。
1.1.2.4渗流速度
(1)过水断面(Cross-sectional area)是渗流场中垂直于渗流方向的任意一个岩石截面,包括空隙面积(A v )和固体颗粒所占据的面积(As ) ,A= Av + As 。渗流平行流动时为平面,弯曲流动时为曲面。
图1-3 渗流过水断面
(a )实际渗透 (b)假想渗流
(2)渗流量(Seepage discharge)是单位时间内通过过水断面的水体积,用Q 表示,单位m 3/d。
(3)渗流速度(Specific discharge/seepage velocity)又称渗透速度、比流量,是渗流在过水断面上的平均流速。它不代表任何真实水流的速度,只是一种假想速度。它描述的是渗流具有的平均速度,是渗流场空间坐标的连续函数,是一个虚拟的矢量。单位m/d,表示为
υ=Q
A (1-1)
(4)实际平均流速(Mean actual velocity)是多孔介质中地下水通过空隙面积的平均速度;地下水流通过含水层过水断面的平均流速,其值等于流量除以过水断面上的空隙面积,量纲为L/T。记为。它描述地下水锋面在单位时间内运移的距离,是渗流场空间坐标的离散函数。表示为:
=Q
A ⋅n (1-1a )
υ与之间存在以下关系:
υ= n (1-2)
若确定渗流场中任一点的渗流速度,可以按以下方法进行讨论:
设以P 点为中心的REV 的平均渗流速度矢量为v ,令REV 的体积为∆V 0,其中空隙体积为n ∆V 0,在空隙中的不同地点,流速u 不同,将u 在全部空隙体积n ∆V 0中求积分,再除以REV 体积∆V 0,即为渗流速度,表示为
v =1⎰(∆V 0) udVv ∆V 0(1-3)
1⎰(∆V v ) 0udVv (∆V v ) 0(1-4) =
可得V= n
1.1.3地下水的水头与水力坡度
(1)地下水水头(hydraulic head):渗流场中任意一点的总水头近似等于测压水头(piezometric head),即
H =Z +P
γ(1-5)
通常称为渗流水头。
在水力学中定义总水头(total head):
u 2
H =Z ++γ2g (1-6) P
式中右端三项分别称为位头(potential head )、压头(pressure head) 和速头(velocity head) 。
总水头(Total head )为测压管水头和流速水头之和。
测压管水头(Piezometric head)为位置水头与压力水头之和,H n =Z +P
γ。
压力水头(pressure head):含水层中某点的压力水头(h )指以水柱高度表示的该点水的压强,量纲为L ,即:h =P/ γ,式中 P 为该点水的压强;γ为水的容重,γ=ρg 。
速度水头(velocity head):在含水层中的某点水所具有的动能转变为势能时所达到的高
u 2
h v =2g ,式中u 为地下水在该点流动的速度;g 为重力加速度。 度,量纲为L ,即
u 2
由于在地下水中水流的运动速度很小,故速头2g 可以忽略,所以h 近似等于H ,即
H ≈H n =Z +P
γ=Z +P
ρg (1-7)
意义:渗流场中任意一点的水头实际上反映该点单位质量液体具有的总机械能,地下水在运动过程中不断克服阻力,消耗总机械能,因此沿地下水流程,水头线是一条降落曲线。
(2) 水力坡度[水力梯度](hydraulic gradient):在渗流场中大小等于梯度值,方向沿等水头面的法线并指向水头下降方向的矢量,用J 表示。
J =-dH n dn (1-8)
n 式中——法线方向单位矢量。在空间直角坐标系中,其三个分量分别为:
J x =-∂H ∂H ∂H , J y =-, J z =-∂x ∂y ∂z (1-9)
(3)等水头面与等水头线
等水头面:渗流场中水头值相同的各点相互连接所形成的一个面。可以是平面也可为曲面。
等水头线(groundwater contour):等水头面与某一平面的交线。
等水头面上任意一条线上的水头都相等。等水头面(线)在渗流场中是连续的,不同大小的等水头面(线)不能相交。
1.1.4 地下水运动特征分类
(1)渗流运动要素(Seepage elements)是表征渗流运动特征的物理量,主要有渗流量Q 、渗流速度V 、压强P 、水头H 等。
地下水运动方向(Groundwater flow direction)为渗透流速矢量的方向。
(2)层流与紊流
层流(laminar flow):水流流束彼此不相混杂、运动迹线呈近似平行的流动。
紊流(turbulent flow):水流流束相互混杂、运动迹线呈不规则的流动。
图1-4 空隙岩石中地下水的层流和紊流
根据Reynolds number判别地下水流态,通常
R e =νd 0νd =γ(0. 75n +0. 23) γ2(1-10)
式中:ν—地下水的渗流速度;
d —含水层颗粒的平均粒径;
d 0—含水层颗粒的有效粒径;
ν—地下水的运动粘度(粘滞系数)。
通常,确定d 的方法有:(1)d=d10;(2)Collins(1961):d =K n ;(3)Ward(1964):
d K , 其中n 为孔隙度。
若ReRe临界,则地下水处于紊流状态,此时液体质点无秩序地相互混杂地流动。 Re 临界≈ 150~300。天然地下水多处于层流状态。
(2)稳定流与非稳定流
根据渗流运动要素是否与时间有关而进行的划分。
稳定流(steady flow ):渗流运动要素不随时间变化;在一定的观测时间内水头、渗流速度等渗透要素不随时间变化的地下水运动。
非稳定流(unsteady flow ):渗流运动要素随时间变化;水头、渗透速度等任一渗透要素随时间变化的地下水运动。
(3)一、二、三维流
根据渗流方向与所选坐标轴方向之间的关系来划分。
一维流运动:当地下水沿一个方向运动,将该方向取为坐标轴,此时地下水的渗透速度只有沿该坐标轴的方向有分速度,其余坐标轴方向的分速度为0。
一维流(one-dimensional flow ),也称单向运动,指渗流场中水头、流速等渗流要素仅随一个坐标变化的水流,其速度向量仅有一个分量、流线呈平行的水流。
图1-5 承压水的一维流动
二维流运动:若地下水的渗透速度沿两个坐标轴方向都有分速度,仅一个坐标轴方向的分速度为0。
二维流(two-dimensional flow ),也称平面运动,地下水的渗透流速沿空间二个坐标轴方向都有分速度、仅仅一个坐标轴方向的分速度为零的渗流;水头、流速等渗流要素随两个坐标变化的水流,其速度向量可分为两个分量,流线与某一固定平面呈平行的水流。
平面二维流(Two-dimensional flow in plane),由两个水平速度分量所组成的二维流。 剖面二维流(two-dimensional flow in section),由一个垂直速度分量和一个水平速度分量组成的二维流。
单宽流量(Discharge per unit width):渗流场中过水断面单位宽度的渗流量,等于总流量Q 与宽度B 之比。即
q=Q/B。 (1-11)
总渗流量Q 为单宽流量q 与宽度B 的乘积,Q=qB。
图1-6 渠道向河流渗漏的地下水二维流动
(a )平面图 (b)剖面图
三维流运动:地下水的渗透流速沿空间三个坐标轴的分量均不为0。
三维流(three-dimensional flow),也称空间运动,地下水的渗透流速沿空间三个坐标轴的分量均不等于零的渗流;水头、流速等渗流要素随空间三个坐标而变化的水流。
图1-7 河弯处潜水的三维流动
(a )平面图 (b)剖面图
图1-8 均质各向同性含水层中潜水井抽水时的地下水运动
(a ) 平面图 (b)剖面图
(b )
1.2 渗流基本定律
1.2.1 达西定律(线性渗透定律)
图1-9 Darcy 实验装置
(1)达西定律表达式
实验条件:定水头、定流量、均质砂。
此时地下水做一维均匀运动,渗流速度与水力坡度的大小和方向沿流程不变。
达西定律(1856年) 表达式:
Q =KAJ =KA H 1-H 2
L (1-12)
V =Q =KJ A (1-13)
其中Q ——渗透流量(出口处流量),亦即通过过水断面(砂柱各断面)A 的流量(m3/d);volumetric flow rate.
K ——多孔介质的渗透系数(m/d);
2A ——过水断面面积(m) ;cross-sectional area of flow.
H 1、H 2——上、下游过水断面的水头(m);
L ——渗透途径 (m);
J ——水力梯度(J = (H1-H 2)/L),等于两个计算断面之间的水头差除以渗透途径,亦即渗透路径中单位长度上的水头损失。
达西定律的微分形式:
v =KJ =-K dH
dn (1-14)
v x =-K
J =-dH dH dH , v y =-K , v z =-K dx dy dz (1-15) dH =-gradH dn
达西定律的矢量形式:
v =v x i +v y j +v z k (1-16)
(2)达西公式讨论
达西定律反映了能量转化与守恒。
V 与I 的一次方成正比;当K 一定时,当V 增大时,水头差增大,表明单位渗透途径上被转化成热能的机械能损失越多,即V 与机械能的损失成正比关系;当V 一定时,K 越小,水头差越大,即K 与机械能的损失成反比关系。
(3) 达西公式适用范围
Re
Re>10-100,层流,不适用,地下水流速增大,为过渡带,由粘滞力占优势的层流转变为以惯性力占优势的层流运动;
Re>100,紊流,不适用。
达西定律的下限:地下水在粘性土中运动时存在一个起始水力坡度J 0。当时及水力坡度J
图1-10 渗透系数与水力坡度的实验关系(J. Bear)
1.2.2渗透系数、渗透率与导水系数
(1)渗透系数(K )(hydraulic conductivity)
V=KI ,当I=1时,V=K,即K 在数值上等于渗流速度,具有速度的单位,它又可以称为水力传导系数,反映含水介质对渗流阻力大小的系数。常用单位:m/d,cm/s。
渗透系数是反映岩石透水性的指标,可以根据渗透系数的大小进行岩石透水性分级(表1-1)。
表1-1 岩石透水性划分
K 的影响因素:
① 岩石的性质:粒度、成分、颗粒排列、充填状况、裂隙性质及其发育程度等,空隙大小起主导作用;
② 流体的物理性质:容重、粘滞性等。 (2)渗透率(intrinsic permeability)
k ρg dH 达西定律可表示为: *v =-
μ
ds (1-17)
式中ρ——液体的密度(density ) ;
g ——重力加速度;
μ——动力粘度;P the dynamic viscosity of the fluid H *=z +
γ——测压水头;
k ——描述多孔介质本身的渗透性能的常数,表示介质能使流体通过其本身的性能,它
不随渗透液体的物理、力学性质而变化。表征岩层渗透性能的参数;渗透率只取决于岩石的性质,而与液体的性质无关,记为k 。单位为cm 2或D ,1D=9.8697×10-9cm 2。
2ρgd 2ρgb u =J u =J
32μ12μ管流公式:; 缝流公式:
nd 2ρg nb 2ρg v =J v =J
32μ; 裂隙介质:12μ 多孔介质:
nb 2ρg ρg nd 2ρg ρg
K ==k K ==k
12μμ 32μμ 或于是有:ρg g
K =k =k
μμ(1-18)
达西(D )的定义:当液体的动力粘滞度为0.001Pa·s ,压强差为101325Pa 的情况下,通过面积为1cm 2、长度为1cm 岩样的流量为1cm 3/s时岩样的渗透率,记为D 。
尺度效应是指渗透系数与试验范围有关,随着试验范围的增大而增大的现象,K=K(x)。亦即抽水时间t 长、降深s 大的群孔抽水试验所得K 较抽水时间t 短、降深s 小的抽水试验所得K 大。
(3)导水系数(transmisivity )
Q=KMBJ
Q=Q/B=KMJ=TJ (1-19)
式中T=KM,称为导水系数,反映含水层出水能力的水文地质参数,其物理意义是水力坡度为1时,通过整个含水层厚度上的单宽流量。它是定义在一维或二维流中的水文地质参数。
单位:m 2/d。
图1-11 导水系数的概念
1.2.3非线性运动方程
(1) Re>1-10,P. Forchheimer(1901)公式:
J =av +bv 2(1-20)
m
或J =av +bv (1-20a )
式中的a ,b 由实验确定的常数, 1.6≤m ≤2。
1(2)当a=0时,有Chezy 公式:
v =k c J 2(1-21)
gk g k (1-22) d 2
k =
360,其中d 2是颗粒直径。 式中
1.3 岩层透水特征及水流折射定律
1.3.1岩层透水特征分类
(1)均质与非均质
根据岩层透水性随空间坐标的变化情况划分,若渗流场中,任意点都具有相同的渗透系数,或渗透系数不随空间坐标的变化而变化,则该岩层是均质的,反之则为非均质。岩石的非均质分两类,一类是渐变的,另一类是突变的。
均质岩层(Homogeneous strata):渗流场中所有点都具有相同参数的岩层。
非均质含水层(inhomogeneous strata ):渗流场中所有点不都具有相同参数的岩层,渗透系数K=K(x,y,z),为坐标的函数。
非均质分为两类,即渐变的和突变的。 (2)各向同性与各向异性
根据岩层透水性与渗流方向的关系划分,若渗流场中,某一点的K 与渗流方向无关,则该岩层是各向同性的,反之则为各向异性。
各向同性岩层(Isotropic strata):渗流场中某一点的渗透系数不取决于方向,即不管渗流方向如何都具有相同渗透系数的岩层。
各向异性岩层(anisotropic strata ):渗流场中某一点的渗透系数取决于方向,渗透系数随渗流方向不同而不同的岩层。
v 0. 55 2
(3)Ward J (=υ+υ
1.3.2渗透系数张量
岩石的透水性是用渗透系数来衡量的。渗透系数实际上是个张量。
(1)对于各向同性介质,其中任一点的渗透系数值与渗流方向无关,是一个标量,水力坡度与渗流方向是一致的。此时,υ可以表示为如下表达式:
υx =-K
∂H ∂H ∂H
, υy =-K , υz =-K ∂x ∂y ∂z
(2)对于各向异性介质,K 与渗流方向有关,K 不再是标量,水力坡度与渗流方向一般是
不一致的。此时,υ可以表示为如下表达式:
∂H ∂H ∂H
-K xy -K xz ∂x ∂y ∂z ∂H ∂H ∂H
υy =-K yx -K yy -K yz
∂x ∂y ∂z ∂H ∂H ∂H
υz =-K zx -K zy -K zz
∂x ∂y ∂z (1-23)
υ=-K xx
即可写成
⎡K xx
⎢K =⎢K yx
⎢K zx ⎣
K xy K yy K zy
K xz ⎤⎥K yz ⎥K zz ⎥⎦ (1-24a )
⎡K xx
K =⎢
⎣K yx 在二维空间中,
即有 υ=K·J (1-25)
K xy ⎤
K yy ⎥⎦ (1-24b )
渗透系数是对称张量,即K xy =Kyx ,K xz =Kzx ,K yz =Kzy 。
在各向异性介质中,水力坡度与渗流方向不一致,但在三个方向上两者是平行的,而且这三个方向称为主方向。
主渗透系数(主值)是指沿主方向测得的渗透系数,用 K 1、K 2、K 3表示,有K xx =K1,K yy =K2,K zz =K3,此时:
⎡K 1
K =⎢⎢0
⎢⎣0
0K 2
0⎤0⎥⎥K 3⎥⎦
(1-26)
1.3.3 层状岩层的等效渗透系数
在自然界中很常见的非均质岩层多是由许多透水性各不相同的薄层相互交替组成的层状岩层。当每一分层的渗透系数K i 和厚度面M i 已知时,可求出平行于层面的渗透系数K p 和垂直于层面的渗透系数K v 。
(1)水流平行层面
图1-11 层状岩层中平行于层面的渗流
特点:水流为稳定流,岩层水平分布,各段流量之和等于各部分流量之和,且各段具有统一的水头,各段具有相同的水力坡度。
根据达西定律有:故
q =∑q i =∑K i M i
i =1
i =1
n n
∆H ∆H
q =K p M
L ,若把其视为整体时,有L ,
n
∆H ∆H K p M =∑K i M i
L L 。水平岩层的等效渗透系数为: i =1
K p =
∑M
i =1
n
i
K i
M
(1-27)
等效导水系数为 T p =∑T i =∑M i K i
i =1
i =1
n
n
(1-28)
垂直方向岩性渐变时,有
1-12 层状岩层中垂直于层面的渗流
1
K p =
M
M 0
⎰
M
K (z ) dz
T p =⎰K (z ) dz
(1-29)
(2)水流垂直层面
特点:水流垂直层面运动,每段水流具有相同的单宽流量,且每段水力坡度不同。
q i =K i b
由
M ∆H ∆H q M 1
, q 2=K 2b ∆H 1=, HH 2=b 2M 1M 2,由此推导出,b K 1K 2
q n M i ∆H =∆H 1+∆H 2+... +∆H n =∑
b i =1K i 依次类推,有
K v =
∑M
i =1
n
n
i
M i ∑i =1K i (1-30)
可见,K v 取决于K i 最小的分层(阻力最大),K i =0,则K v =0。 另外,总是有K p >K v 。 1.3.4突变界面的水流折射定律
根据水流连续性条件,当水流斜向由一种介质进入另一种介质时,会发生折射。 如图所示:水流由K 1介质进入 K 2介质中,二者交界面上某一点的渗流速度和水头在两介质中的值依次为V 1、V 2和H 1、H 2。对于界面上的任一点应满足以下条件:
⎧H =H 2⎨
⎩V 1n =V 2n (1-31)
由图中几何条件有
tg θ1=
v 1τv , tg θ2=2τv 2n v 2n ,则有
tg θ1v 1τ
=tg θ2v 2τ
∂H 1=
∂H 2
-K 2
∂x -K 1
∂H 1∂H 2
=∂x ∂x ,则得到水流折射定律(渗流折射时必须满足的方程)因为:
tg θ1K 1
=
tg θ2K 2(1-32)
图1-12 渗透水流的折射
讨论上式可以的出以下结论:
(1)若K 1=K2,则θ1=θ2,表明在均质介质中水流不发生折射。
(2)若K 1≠K 2,且K 1,K 2均不为0,若θ1=0,则θ2=0,表明水流垂直通过界面时水流不发生折射。
(3)若K 1≠K 2,且K 1,K 2均不为0,若θ1=90,则θ2=90,表明水流平行于界面
时水流不发生折射。
(4)当水流斜向通过界面时,介质的渗透系数越大,θ值也越大,流线也越靠近界面。介质相差越大,两角的差值也越大。
根据水流折射原理和达西定律,可以帮助分析流场的水动力条件的变化。
1.4 流网及其应用
1.4.1 流网的概念
(1)流网(flow net ):渗流场中由一组流线与由一组等势线(当容重不变时为一组等水头线)相交组成的网格。对各向同性介质组成正交网。
流线(Streamline )渗流场内处处与渗流速度矢量相切的曲线。
地下水动力学中流线的概念和水力学中的概念是完全一致的。流线应是一根处处和渗流速度矢量相切的曲线。因此,流线簇就代表渗流区内每一个点的水流方向。 1.4.2流函数方程
(1) 流线的方程
根据上述定义,没有水流穿越流线。现在来研究描述流线的方程式。如图,在任一流线上取任意两点M(x, y)和M' (x+dx, y+dy)。M 点的渗流速度矢量为v ,它与它的两个分量Vx ,Vy 构成一个三角形MAB 。自M' 点作垂线Mb ,并延长至a ,见图。
图1-13 流线
当M 与M' 无限逼近时,弧线MM ' 可用切线Ma 来代替,故有Mb= dx,ab=dy。因为∆MAB ≈∆Mab ,有以下等式成立: 流线方程
dx dy
=
v x v y (1-33a)
v x dy -v y dx =0
(1-33b)
M 和M’是任意流线上任选的两点。因此,上式对流线上的任一点都是正确的,可以把它看成是流线的方程,用它来描述流线。上面的流线方程无论对各向同性和各向异性介质都是适用的。在各向异性介质中,如果选取的坐标轴(直角坐标系) 的方向分别与渗透系数的主方向一致,则上式变为:
dy dx
=∂H ∂H K xx K yy
∂x ∂y
对于各向同性介质,则式中的K xx =K yy =K 。由于(1-33b )式只涉及一个点的水流情况,故也适用于非均质介质。
(2) 流函数方程
d ψ=
设有二元ψ函数Ψ(z,y),其全微分为:
∂ψ∂ψ
dx +dy ∂x ∂y 。若取这样一种函数,使
∂ψ∂ψ
=-v y , =-v x ∂x ∂y (1-34)
d ψ=
则
∂ψ∂ψ
dx +dy =v x dy -v y dx =0∂x ∂y (1-35)
对其积分得:ψ=常数。表明沿同一流线,函数ψ为常数,不同的流线则有不同的函数值。称函数ψ为流函数,又称Lagrange 流函数,量刚为[L 2T -1]。
(3)流函数的物理意义
在无限接近的两条流线ψ和ψ+d ψ上沿某等水头线取两个点a(x,y)和b(x+dx,y+dy)。自a 、b 分别做垂线和水平线,相交于c 。见图1-14。
通过两流线之间的单宽流量dq 可看作是通过ac 和bc 的流量的代数和。将渗流速度分解则有:dq =v x ac+vy bc ,但ac=dy,bc=-dx,所以有dq=vx dy-v y dx
把式(1-34)和(1-35)代入上式,则得到:
dq =
将(35)式在ψ1和ψ1区间积分得:
∂ψ∂ψ
dy +dx =d ψ∂y ∂x (1-36)
图1-14 流函数与流量之间的关系
q =⎰d ψ=ψ2-ψ1
ψ1
ψ2
(1-37)
由(1-37)可以得出:在平面运动中,两流线之间的单宽流量等于和这两条流线相应的流函数之差。在同一条流线上,d ψ=0,q=0,ψ=常数。
由达西定律和(1-48)式,有:
v x =-K
∂H ∂ψ∂H ∂ψ
=, v y =-K =-∂x ∂y ∂y ∂x (1-38)
将(1-38)中第一式对y 求导,第二式对 x求导,得到:
∂2H ∂2ψ∂2H ∂2ψ
-K =, -K =-2
∂x ∂y ∂y 2∂y ∂x ∂x ∂2ψ∂2ψ
=-2, 2
∂y ∂x 整理得:
满足该方程。
(4) 流函数的特性
① 对于一给定的流线,流函数是常数。不同的流线有不同的常数值。流函数决定于流线。ψ=c
∂2ψ∂2ψ
+=0∂y 2∂x 2
(1-39)
表明在均质各向同性介质中,流函数满足Laplace 方程,而在其他情况下,流函数均不
② 在平面运动中,两流线之间的单宽流量等于和这两条流线相应的流函数之差。q=ψ
1-ψ2
③ 在均质各向同性介质中,流函数满足Laplace 方程;而在其他情况下,流函数均不满足该方程。
④ 在非稳定流中,流线不断地变化,只能给出某一瞬时的流线图。只有对不可压缩的液体的稳定流动,流线才有实际意义。 1.4.2流网的性质
(1)在各向同性介质中,流线与等水头线处处垂直,流网为正交网格。 由(1-38)式,得:
-K
消去K ,得:
∂H ∂ψ∂H ∂ψ
=K
∂y ∂y ∂x ∂x (1-40)
∂H ∂ψ∂H ∂ψ
+=0
∂x ∂x ∂y ∂y (1-41)
gradH =∇H =
等水头线
∂H ∂H
i +j ∂x ∂y
grad ψ=∇ψ=
流线
式中i ,j ——单位矢量。
∂ψ∂ψi +j ∂x ∂y
∇H ⋅∇ψ=
∂H ∂ψ∂H ∂ψ
+=0
∂x ∂x ∂y ∂y (1-42)
∇ϕ⋅∇ψ=0(1-43)
在非均质各向同性介质中,上式亦成立。
(2)在均质各向同性介质中,流网中每一网格的边长比为常数。
dH =
∂H ∂H 1dx +dy =-υds ∂x ∂y K
∂H ∂H 1d ψ=dx +dy =υds
∂x ∂y K
ds dH =-K dl d ψ(1-44)
式中dl ——相邻流线的间距;
ds ——等势线的间距。
通常取ds/dl=1,流网为曲边正方形。
(3)若流网中各相邻流线的流函数差值相同,且每个网格的水头差值相等时,通过每个网格的流量相同。
∆q =KJ ∆l =K
∆H ∆l
∆l =K ∆H ∆s ∆s (1-45)
∆l
当∆s =1时,∆q =K ∆H (1-46)
式中∆s ——网格相邻两等势线间的平均长度;
∆l ——网格相邻两流线间的平均宽度。
若上下游总水头差H r =H1-H 2,则m 个水头带中每一网格的水头差为
∆H =
H r
m (1-47)
(4)若两个透水性不同的介质相邻时,在一个介质中为曲边正方形的流网,越过界面进入另一个介时则变成曲边矩形。
∆q =K 1
∆H ∆H
∆l 1=K 2∆l 2∆s 1∆s 2
∆l 1K 2∆l 2
=
∆s 1K 1∆s 2
∆l 1∆l 2∆l 1K 1
=1=≠1∆s ∆s ∆s K 12当K 1≠K 2,且1时,2。
1.4.3流网的绘制与应用 1.4.3.1 流网的绘制
可采用解析法、各种模型试验法、徒手绘渐进法绘制流网。 (1)确定边界条件
① 河渠的湿周为一条等水头线。 ② 平行于隔水边界可绘制流线。
③ 无入渗补给及蒸发排泄,有侧向补给,做稳定运动时,地下水面是一条流线。 ④ 有入渗补给时,地下水面既不是流线也不是等水头线。
(2)根据已经确定的边界条件,根据流网的性质可以判定另一条件,作出流网。 (3)根据流网的性质绘制,各向同性含水层中,流线与等水头线处处正交,网格边长比为常数。
(4)在同一渗流区内,除奇点外,流线与等水头线各自不能相交;如遇透水性大的透
镜体时,则流线向该点汇集,反之则绕行。流线穿越突变界面时,应用水流折射定律绘制。 1.4.3.2流网的应用
(1)定量计算渗流区中的渗流运动要素
①水头H 、渗透压强P
P
γ
=H ±Z , P =γ(H ±Z )
(1-48)
②水力梯度J 、渗流速度v
J =
∆H
, ∆s
v =KJ
(1-49)
③流量q
q =K ∆Hn =
Kn (H 1-H 2)
m (1-50)
式中m 、n ——水头带数目、流带的数目。
(2)定性分析渗流区的水文地质条件及其变化。 (3)主要用于解决稳定渗流问题。
1.5 渗流连续方程
1.5.1含水层的状态方程
含水层的状态方程主要包括地下水的状态方程和多孔介质的状态方程。 1.5.1.1 地下水的状态方程
Hooke 定律:
E =-V
dp dV
式中E ——体积弹性系数(体积弹性模量),20℃时,E=2.1×105N/cm2。其倒数为压缩系数。
等温条件下,水的压缩系数(coef. of compressibility)为
β=-
1dV
V dp
(1-51)
积分(p →p 0,V →V 0)改写得:
-β(p -p 0)
体积:V =V 0e
(1-52)
β(p -p )
ρ=ρe 0密度:
(1-53)
(1-54)
按Taylor 级数展开,得到近似方程: V =V 0[1-β(p -p 0)] 和
ρ=ρ0[1+β(p -p 0)]
(1-55)
因 d (ρV ) =ρdV +Vd ρ=0(质量守恒),故有
d ρ=-ρ
dV
=ρβdp V
(1-56)
1.5.1.2 多孔介质的状态方程
多孔介质压缩系数(Coefficient of compressibility)表示多孔介质在压强变化时的压缩性的指标,用α 表示。
多孔介质压缩系数α的表达式为
α=-
1dV b
V b d δ(1-57)
式中,V b =V s +Vv ——多孔介质中所取单元体的总体积;
V s ——单元体中固体骨架(solid matrix)体积; V v ——为其中的孔隙(voids )体积。 δ——介质表面压强; V v =nV b ;V s =(1-n)Vb
α=-
1dV s 1dV v 1-n dV s n dV v
-=---
V b d δV b d δV s d δV v d δ
(1-58)
(1-59)
α=(1-n ) αs +n αp
αs =-
式中
1dV s
V s d δ——多孔介质固体颗粒压缩系数,表示多孔介质中固体颗粒本身的压缩性
的指标, αs
αp =-
1dV v
V v d δ——多孔介质中孔隙压缩系数(Compressibility of the pores of a porous
medium ),表示多孔介质中孔隙的压缩性的指标。
n ——多孔介质的孔隙度。
因(1-n ) αs
考虑承压含水层受力情况,取一水平横截面AB ,按Terzaghi (1883~1963)观点:
σ=λσs +(1-λ) p
(1-60)
式中σ——上覆荷重引起的总应力(total stress);
σs ——作用在固体颗粒上的粒间应力 (intergranular stress);
λ——横截面面积中颗粒与颗粒接触面积所占的水平面积比;
p ——水的压强。
图1—1 一个可压缩的承压含水层(J. Bear)
Terzaghi 令λσs =σ' ,称为有效应力(effective stress)。λ很小,(1-λ)p ≈ p,因此有
σ=σ' +p
(1-61)
在水位下降为∆H 时,有σ=(p -γ∆H ) +(σ' +γ∆H ) 。即作用于固体骨架上的力增加了γ∆H 。作用于骨架上力的增加会引起含水层的压缩,而水压力的减少将导致水的膨胀。含水层本来就充满了水,骨架的压缩和水的膨胀都会引起水从含水层中释出,前者就象用手挤压充满了水的海绵会挤出水—样。
dV b dn
=
1-n 。 因Vs=constant,故V b
dV b d (∆z ) =
dz ,故 只在垂直方向上有压缩,V b
d (∆z ) =∆z αdp
dn =(1-n ) αdp
(1-62) (1-63)
上两式表示垂直厚度变化、孔隙度变化与水的压强变化的关系。
为了讨论水头降低时含水层释出水的特征,取面积为1 m2、厚度为l m (即体积为l m3) 的含水层,考察当水头下降1m 时释放的水量。此时,有效应力增加了γ∆H =ρg ×1 = ρg 。介质压缩体积减少所释放出的水量(dV b )与水体积膨胀所释放出的水量(dV )之和为
-dV b =αV b dp =α⋅1⋅ρg =αρg
dV =-βdp =-βn (-ρg ) =n βρg
μs =-dV b +dV =αρg +n βρg
或μs =ρg (α+n β)
(1-64)
式中μs ——贮水率[释水率](specific storativity),量纲[L-1],为弹性释水[贮水]。
μ* = μs M
式中M ——含水层厚度(m);
μ*贮水系数(storativity )。
贮水系数μ* 和贮水率 μs 都是表示含水层弹性释水能力的参数,在地下水动力学计算中具有重要的意义。
贮水率
μs =(α+n β) ρg ,表示含水层水头变化一个单位时,从单位体积含水层中,
因水体积膨胀(压缩)以及骨架的压缩(或伸长)而释放(或储存)的弹性水量。单位1/L。
贮水率是描述地下水三维非稳定流或剖面二维流中的水文地质参数,既适用于承压水也适用于潜水。对于平面二维非稳定流地下水运动,当研究整个含水层厚度上的释水情况时,用贮水系数来体现。
贮水系数又称释水系数或储水系数,为含水层水头变化一个单位时,从底面积为一个单位,高度等于含水层厚度的柱体中所释放(或贮存)的水量;指面积为一个单位、厚度为含水层全厚度M 的含水层柱体中,当水头改变一个单位时弹性释放或贮存的水量,无量纲。既适用于承压含水层,也适用于潜水含水层。
⎧μs M =μ*.......... (承压含水层)
E =⎨
⎩μs (H -Z ) =μ...(潜水含水层)
μ*范围值:n ×10-3~ n×10-5;μ范围值:0.05~0.3。实际测出的值往往小于理论值。
弹性释水与重力给水:对于含水层而言,由于受埋藏条件的限制,抽水时,水的给出存在着不同。
潜水含水层在抽水过程中,大部分水在重力作用下排出,疏干作用于水位变动带(饱水带)和包气带两部分,由于包气带的存在,使得饱水带中水的释放存在延滞和滞后现象。当水头下降时,可引起二部分水的排出。在上部潜水面下降部位引起重力排水,用给水度μ 表示重力排水的能力;在下部饱水部分则引起弹性释水,用贮水率μ* 表示这一部分的释水能力。必须区分两者之间的不同,潜水含水层还存在滞后疏干现象。
承压含水层抽水时,水的释放是由于压力减少造成的,这一过程是瞬时完成的。只要水头下降不低到隔水顶板以下,水头降低只引起含水层的弹性释水,可用贮水系数μ* 表示这种释水的能力。 1.5.1.4导压系数
描述含水层水头变化的传导速度的参数,其数值等于含水层的导水系数与贮水系数之比或渗透系数与贮水率之比。
a =
T K =
E μs ,量刚为L 2T -1。
1.5.2 渗流连续方程
由于渗流场中各点的渗流速度大小、方向都不同,为了反映液体运动的质量守恒关系,需要在三维空间中建立微分方程形式表达的连续性方程。
在渗流场中任意取一点P(x, y, z),以P 为中心沿直角坐标轴取一微小的六面体,体积为
∆x ∆y ∆z ,称为特征单元体,设单元体无限小,但保证单元体穿过介质骨架和空隙。
图1-15 渗流区中的单元体
设v x ,v y ,v z 分别为该点在X 、Y 、Z 方向上的渗流速度。Abcd 面中点沿X 轴方向流入:
P 1=(x -
∆x
, y , z ) 2。
M xi =ρυx 1=ρυx (x -
流出:
M xo =ρυx 2=ρυx (x +
∆x
, y , z ) ∆y ∆z ∆t 2;
∆x
, y , z ) ∆y ∆z ∆t 2。
利用Taylor 级数展开,略去二阶导数以上的高次项,有:
∆M x =M xi -M xo =[ρυx -
1∂(ρυx ) 1∂(ρυx )
∆x ]∆y ∆z ∆t -[ρυx +∆x ]∆y ∆z ∆t
2∂x 2∂x
∂(ρυx )
∆x ∆y ∆z ∆t ∂x = -
∂(υz ρ) ∂(υx ρ)
∆M =-∆x ∆y ∆z ∆t ∆M y =-∆x ∆y ∆z ∆t z
∂z ∂x 同理 ;
单元体本身水质量在Δt 时间内的变化量由质量守恒定律,得到渗流的连续性方程:
∆M =
∂
(ρn ∆x ∆y ∆z ) ∆t ∂t ,ρ为液体密度。
∆M x +∆M y +∆M z =∆M
⎡∂(ρυx ) ∂(ρυy ) ∂(ρυz ) ⎤∂
-⎢++∆x ∆y ∆z =(ρn ∆x ∆y ∆z ) ∆t ⎥
∂x ∂y ∂z ∂t ⎣⎦(1-65)
或
-div (ρυ) ∆x ∆y ∆z =
∂
(ρn ∆x ∆y ∆z ) ∂t
上式即为非稳定流的渗流连续方程,表明渗流场中任意体积含水层流入、流出该体积含水层中水质量之差永远恒等于该体积中水质量的变化量。它表达了渗流区内任何一个“局部”所必须满足的质量守恒定律。
若把含水层看作刚体,ρ=constant,n 不变,即水和介质没有弹性变形或渗流为稳定流,
则渗流连续性方程为
div (υ) =
∂v x ∂v y ∂v z
++=0∂x ∂y ∂z (1-66)
此式表明,在同一时间内流入单元体的水体积等于流出的水体积,即体积守恒。
连续性方程是研究地下水运动的基本方程,各种研究地下水运动的微分方程都是根据连续性方程和反映质量守恒定律的方程建立起来的。
1.6 渗流基本微分方程
1.6.1承压水运动的基本微分方程
基本假设:
(1)单元体体积无限小,为承压含水层;
(2)含水层侧向受到限制,∆x 、∆y 为常量,∆z 为变量,存在垂向压缩,水的密度ρ、孔隙度n 和∆z 随压力p 而变化;
(3)由ρ引起的变化
(υx
∂ρ∂ρ∂ρ+υy +υz ) ∆x ∆y ∆z ∂x ∂y ∂z ,远小于单元体内液体质量的
∂ρ
∂c ,可忽略不计;
变化量(含ρ)
υgrad υ
(4)水流服从Darcy’s Law;
(5)K 不因ρ=ρ(p ) 的变化而变化;
(6)μs 和K 也不受n 变化(由于骨架变形)的影响。 流体的质量:M =ρV =ρn ∆x ∆y ∆z
由于含水层的侧向受到限制,可假设∆x 、∆y 为常量,只考虑垂向压缩。于是,只有水的密度ρ、孔隙度n 和单元体高度∆z 三个量随压力而变化,于是有:
⎡∂(∆z ) ∂p ∂n ∂p ∂ρ∂p ⎤∂(ρn ∆z ) ∂M ∂
n ρ+ρ∆z +n ∆z =(ρn ∆x ∆y ∆z ) =∆x ∆y ⎢⎥∆x ∆y ∂p ∂t ∂p ∂t ∂p ∂t ⎦∂t ∂t ∂t =⎣
⎡∂(∆z ) ∂n ∂ρ⎤∂p
n ρ+ρ∆z +n ∆z ⎢⎥∂t ∆x ∆y ∂p ∂p ∂p ⎦=⎣ (1-67)
d ρ
=βρ, dp 由含水层状态方程,
d (∆z ) dn p
=α(1-n ), =α∆z H =Z +dp dp γ,所以有,因为
p =γ(H -Z ) ,Z 为定值,则dp =γdh ,则可得到:
∂M ∂H ∂H
=[n ρα∆z +ρ∆z α(1-n ) +n ∆z βρ]γ∆x ∆y =(α+n β) ρ2g ∆x ∆y ∆z ∂t ∂t ∂t (1-68)
于是连续性方程(1-65)变为:
⎡∂(ρυx ) ∂(ρυy ) ∂(ρυz ) ⎤2∂H -⎢++∆x ∆y ∆z =(α+n β) ρg ∆x ∆y ∆z ⎥
∂y ∂z ⎦∂t ⎣∂x (1-69)
∂p ∂H ∂ρ∂H p ∂ρ
=ρg +(H -Z ) g =ρg +
∂t ∂t ∂t ρ∂t 又 ∂t
∂p ρg ∂H ∂H
=≈ρg
∂t (1-70) 则 ∂t 1-βp ∂t
∂M ∂H
=ρμs ∆x ∆y ∆z
∂t 令 μs =(α+n β) ρg ,则 ∂t
根据连续性原理有:
⎡∂(ρυx ) ∂(ρυy ) ∂(ρυz ) ⎤∂M
=-⎢++⎥∆x ∆y ∆z ∂t ∂x ∂y ∂z ⎣⎦
⎡∂(ρυx ) ∂(ρυy ) ∂(ρυz ) ⎤∂M ∂H
=ρμs ∆x ∆y ∆z =-⎢++⎥∆x ∆y ∆z ∂t ∂t ∂x ∂y ∂z ⎣⎦则有:
⎡∂(ρυx ) ∂(ρυy ) ∂(ρυz ) ⎤∂H -⎢++∆x ∆y ∆z =ρμ∆x ∆y ∆z ⎥s
∂x ∂y ∂z ∂t ⎦即:⎣
将
υx =-K
∂h ∂h ∂h
, υy =-K , υz =-K ∂x ∂y ∂z 代入,整理得:
ρ⎢
⎡∂∂H ∂∂H ∂∂H ⎤∂H (K ) +(K ) +(K ) ⎥∆x ∆y ∆z =ρμs ∆x ∆y ∆z
∂x ∂y ∂y ∂z ∂z ⎦∂t ⎣∂x ,所以有 ∂∂H ∂∂H ∂∂H ∂H (K ) +(K ) +(K ) =μs ∂x ∂x ∂y ∂y ∂z ∂z ∂t (1-71)
上式为三维流微分方程,也可写成:
div K grad (h ) =μs
]
∂H
∂t
物理意义:渗流空间内任一单位体积含水层在单位时间内流入与流出该体积含水层中的弹性水量的变化量,即单位体积含水层的水量均衡方程。
基本微分方程是研究承压含水层中地下水运动的基础。它反映了承压含水层中地下水运动的质量守恒关系,表明单位时间内流入、流出单位体积含水层的水量差(左端)等于同一时间内单位体积含水层弹件释放(或弹性贮存) 的水量(右端)。它还通过应用Darcy 定律反映了地下水运动中的能量守恒与转化关系。可见,基本微分方程表达了渗流区中任何一个“局部”都必须满足质量守恒和能量守恒定律。
数学意义:表示渗流空间内任一点任一时刻的渗流规律。
在柱坐标系中,有
1∂∂H 1∂2H ∂2H μs ∂H
(r ) +2+2=r ∂r ∂r K ∂t (1-72a ) r ∂θ2∂z
∂2H 1∂H 1∂2H ∂2H μs ∂H ++2+2=22r ∂r r ∂θK ∂t (1-72b ) ∂z 或 ∂r
由地下水流基本微分方程(1-71),在均质各向同性介质中,方程简化为:
∂2H ∂2H ∂2H μs ∂H
+2+2=
K ∂t (1-73) ∂x 2∂y ∂z
对于各向异性介质,若把坐标轴方向和各向异性介质的主方向定为一致,则有
∂∂H ∂∂H ∂∂H ∂H
(K xx ) +(K yy ) +(K zz ) =μs ∂x ∂x ∂y ∂y ∂z ∂z ∂t (1-74)
在二维流情况下,基本微分方程可表示为:
∂∂H ∂∂H ∂H (T ) +(T ) =μ*∂x ∂x ∂y ∂y ∂t (1-75)
上式即为承压水平面二维流微分方程,该方程是研究承压水含水层中地下水运动的基础,反映了承压水含水层中地下水运动的质量守恒关系,表明单位时间流入、流出单位体积含水层的水量差等于同一时间内单位体积含水层弹性释放(或贮存)的水量。
在实际渗流问题中若存在抽、注水及越流影响,只要在微分方程中的左端中通过加、减W 项,通常把该项称为源汇项。所谓的源项表示在垂直方向上有水流入含水层,此时W 为正;汇指在垂直方向上有水流出含水层,此时W 为负。
此时(1-71)式变成:
∂∂H ∂∂H ∂∂H ∂H (K ) +(K ) +(K ) +w =μs ∂x ∂x ∂y ∂y ∂z ∂z ∂t (1-76)
二维流情况下:
∂∂H ∂∂H ∂H (T ) +(T ) +w =μ*∂x ∂x ∂y ∂y ∂t (1-77)
*
a =T /μ在二维流情况下,令压力传导系数(导压系数),则均质各向同性含水层
基本微分方程为:
∂2H ∂2H 1∂H
++w =
a ∂t (1-78) ∂x 2∂y 2
非均质各向同性含水层中的稳定流运动:
∂∂H ∂∂H ∂∂H (K ) +(K ) +(K ) =0∂x ∂x ∂y ∂y ∂z ∂z (1-79)
均质各向同性含水层中的稳定流运动:
∂2H ∂2H ∂2H
+2+2=02
∂x ∂y ∂z (1-80)
上式也称Laplace 方程。稳定运动方程的右端都等于零,意味着同一时间内流入单元体的水量等于流出的水量。这个结论不仅适用于承压含水层,也适用于潜水含水层和越流含水层。
1.6.2 越流含水层中地下水运动的基本微分方程
在自然界中,存在以下情况,承压含水层的上、下岩层并不是绝对隔水的,其中一个或两个可能是弱透水层。虽然含水层会通过弱透水层和相邻含水层发生水力联系,但它还是处于承压状态,将其称为半承压含水层。当该含水层和相邻含水层间存在水头差时,地下水就会从高水头含水层通过弱透水层流向低水头含水层。这种现象称为越流。半承压含水层称为越流含水层。
假设:主含水层渗透系数K 远远大于若透水层的渗透系数K 1;主含水层弹性释放的水量、弱透水层的越流量远远大于弱透水层弹性释放的水量。主含水层中的水流近似地看作二维流
1M
H =H(x,y, t) =H(x,y, z, t)dz ⎰0M 问题, 。
对于均衡单元体,根据水均衡原理可以写出下列形式的连续性方程:
[(Q x -
∂Q y ∆y ∂Q y ∆y ∂Q x ∆x ∂Q ∆x
) -(Q x +x )]∆t +[(Q y -) -(Q y +)]∆t ∂x 2∂x 2∂y 2∂y 2
+(v 2-v 1) ∆x ∆y ∆t =μ*
∂H
∆x ∆y ∆t ∂t (1-81)
式中,v 1,v 2分别为通过上部和下部弱透水层的垂直越流速率或越流强度,即
v 1=-K 1
∂H 1H -H 1∂H 2H -H 1
=K 1, v 2=-K 2=K 22∂z m 1∂z m 2(1-82)
其中,H 1(x, y, t)和H 2(x, y, t)分别为上含水层和下含水层中的水头,如T 表示主含水层的导水系数,则得到不考虑弱透水层弹性释水条件下非均质各向同性越流含水层中非稳定运动的基本微分方程:
H -H 1H -H ∂∂H ∂∂H ∂H
(T ) +(T ) +K 1+K 22=μ*∂x ∂x ∂y ∂y m 1m 2∂t (1-83)
对于均质各向同性介质来说,有:
∂2H ∂2H H -H 1H 2-H μ*∂H
+++=2222
T ∂t (1-84) ∂x ∂y B 1B 2
B 1=
式中
Tm 1
, K 1
B 2=
Tm 2
K 2(1-85)
分别称为上、下两个弱透水层的越流因素。
越流因素B 的量纲为[L]。弱透水层的渗透性愈小,厚度越大,则B 越大,越流量越小。在自然界中,越流因素值的变化很大,可以从几米到若干公里。对于一个完全隔水的覆盖层来说,B 为无穷大。
另一个反映越流能力的参数是越流系数σ' 。其定义为:当主含水层和供给越流的含水层间的水头差为一个长度单位时,通过主含水层和弱透水层间单位面积界面上的水流量。因此,
σ' =
K 1
m 1(1-86)
K 1、m 1分别为弱透水层的渗透系数和厚度。σ' 越大,相同水头差下的越流量越多。 1.6.3 潜水运动的基本微分方程 1.6.3.1Dupuit 假设
图1-16 Dupuit假设
在潜水面上任意取一点P ,有:
dH dz J =-=-=-sin θ
ds ds (1-87)
该点的流速v 方向与潜水面相切,则由达西定律有:v s =-KJ=-Ksin θ。
当θ很小时,tg θ=sinθ。此时,(1)潜水面比较平缓,等水头面呈铅直,水流基本水平,可忽略渗流速度的垂直分量v Z ;(2)隔水底板水平,铅垂剖面上各点的水头都相等,各点的水力坡度和渗流速度都相等,H (x , y , z , t ) 可以近似地用H (x , y , t ) 代替,此即著名的Dupuit 假设。
渗流速度:
v x =-K
dH
dx ,H=H(x) (1-88)
通过宽度B 的铅直平面的流量为
Q x =-Kh
dH
dx , H=H(x) (1-89)
式中Q x ——x 方向的流量;
h ——潜水含水层厚度;h=H(隔水层水平时)。 对于更一般情况,H=H(x,y) 有:
V x =-K
∂h ∂H
, V y =-K , ∂x ∂y (1-90) ∂H ∂H , Q y =-KhB ∂x ∂y (1-91)
则得:
Q x =-KhB
由于Dupuit 假设的引入,将垂直方向的水流速度忽略,减少了z 变量,简化了计算,但会
2
产生一定的误差,经验证明当i
h 2⎡h 2⎤-⎢h -⎥2⎣2⎦i 2dh 0≤
dx h 1+i
2(1-92)
表1- 用tg θ 代替sin θ的误差
θ(ο) 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22
tg θ 0.01746 0.03492 0.05241 0.06993 0.08749 0.10510 0.12278 0.14054 0.15838 0.17633 0.19438 0.21256 0.23087 0.24933 0.26795 0.28675 0.30573 0.32492 0.34433 0.36397 0.38386 0.40403
sin θ 0.01745 0.03490 0.05234 0.06976 0.08716 0.10453 0.12187 0.13917 0.15643 0.17365 0.19081 0.20791 0.22495 0.24192 0.25882 0.27564 0.29237 0.30902 0.32557 0.34202 0.35837 0.37461
D -0.0000027 -0.0000213 -0.0000718 -0.000170 -0.000333 -0.000576 -0.000915 -0.001368 -0.001950 -0.002679 -0.003571 -0.004645 -0.005917 -0.007406 -0.009130 -0.011108 -0.013359 -0.015903 -0.018759 -0.021950 -0.025496 -0.029420
% -0.01523 -0.06095 -0.1372 -0.2442 -0.3820 -0.5508 -0.7510 -0.9828 -1.247 -1.543 -1.872 -2.234 -2.630 -3.061 -3.528 -4.030 -4.569 -5.146 -5.762 -6.418 -7.114 -7.853
23 24 25 26 27 28 29 30
0.42447 0.44523 0.46631 0.48773 0.50953 0.53171 0.55431 0.57735
0.39073 0.40674 0.42262 0.43837 0.45399 0.46947 0.48481 0.50000
-0.033744 -0.038492 -0.043689 -0.049361 -0.055535 -0.062238 -0.069499 -0.077350
-8.636 -9.464 -10.338 -11.260 -12.233 -13.257 -14.335 -15.470
对于各向异性介质,上式中的i 代以(K xx /K zz ) i 。
2
2
Dupuit 假设无效的地区: (1)存在入渗的潜水分水岭地段;
(2)渗出面附近。渗出面(seepage surface) 是在下游边界面上,潜水面以下、下游水面以上的地段。渗出面上潜水面往往和边界面相切,有较大的垂向分速度。 (3)垂直的隔水边界附近。
图1-17 Dupuit假设无效的地区
1.63.2 Boussinesq方程
图1-18 潜水的非稳定运动
根据Dupuit 假设,可建立有关潜水含水层中的地下水流方程。
(1)潜水一维流方程(沿x 方向运动)
在Δt 时间内,上、下游流入、流出单元体的水量差为:
(q -
∂(υx H ) ∂q ∆x ∂q ∆x ∂q
⋅) ∆t -(q +⋅) ∆t =-⋅∆x ∆t =-∆x ∆t ∂x 2∂x 2∂x ∂x
在该段时间内,垂直方向的补给量为W ∆x ∆t ,故Δt 时段水量总的变化量为
⎡∂(υx H ) ⎤
-+W ⎢⎥∆x ∆t ∂x ⎣⎦。
∂H
由于水量的变化引起潜水面的升降,设其变化的速率为∂t ,则在Δt 时段,由于潜水
面的变化而引起的小土体内的水体积的增量为
μ
∂H
∆x ∆t ∂t ,则有
⎡∂(υx h ) ⎤dh -+W ∆x ∆t =μ∂H ∆x ∆t υx =-K ⎢⎥∂x ⎣⎦dx (均质各向同性)代入上式,可以得∂t 。将
到有入渗补给的潜水含水层中地下水非稳定运动的一维流方程,又称为Boussinesq 方程:
∂∂H W μ∂H
(h ) +=∂x ∂x K K ∂t
(1-93)
式中K 、μ—— 潜水含水层的渗透系数、给水度;
W ——含水层单位时间、单位面积上的垂向补排量,补给为正,排泄为负。 (2)潜水二维流方程
均质各向同性含水层,Boussinesq 方程为:
∂∂H ∂∂H W μ∂H
(h ) +(h ) +=∂x ∂x ∂y ∂y K K ∂t (1-94)
式中h=H-Z,Z=0时,h=H。
非均质含水层,Boussinesq 方程为:
∂∂H ∂∂H ∂H (Kh ) +(Kh ) +W =μ∂x ∂x ∂y ∂y ∂t (1-95)
在推导潜水基本微分方程时应用了Dupuit 假设,忽略了弹性储存,所选的单元体是一个包括了整个含水层厚度在内的土柱,这与承压水非稳定运动时选取的无限小的单元体不同。所以,应用潜水运动基本方程得到的H(x,y ,t) 只能代表该点整个含水层厚度上平均水头的近似值,不能用来计算同一垂直剖面上不同点的水头变化。
(3)潜水三维流方程
若不用Dupuit 假设,Boussinesq 方程的一般形式:
∂∂H ∂∂H ∂∂H ∂H (K ) +(K ) +(K ) =μs ∂x ∂x ∂y ∂y ∂z ∂z ∂t (1-96)
在上面的潜水基本运动微分方程中右端项为贮水率而不是给水度,其原因在于,当不考虑Dupuit 假设时,单元体位于渗流区内部,其贮存量的变化只能是弹性释水而不是疏干排水,因此推导出的潜水非稳定运动方程和承压水非稳定运动方程形式一样。在这种情况下,地下水非稳定运动的特点由边界条件来反映。
对于各向异性介质,坐标轴方向同主方向,有]
∂∂H ∂∂H ∂∂H ∂H (K xx ) +(K yy ) +(K zz ) =μs ∂x ∂x ∂y ∂y ∂z ∂z ∂t (1-97)
假设固体骨架是不可压缩的,μs =0,同时假设忽略水的压缩性, 即ρ=常数,有:
∂∂H ∂∂H ∂∂H (K ) +(K ) +(K ) =0∂x ∂x ∂y ∂y ∂z ∂z (1-98)
∂∂H ∂∂H ∂∂H (K xx ) +(K yy ) +(K zz ) =0∂x ∂x ∂y ∂y ∂z ∂z 或 (1-99)
(4)潜水稳定运动的微分方程
没有入渗和蒸发时,潜水稳定运动的方程式为:
∂∂H ∂∂H
(Kh ) +(Kh ) =0
∂x ∂y ∂y 非均质 ∂x (1-100)
∂∂H ∂∂H
(h ) +(h ) =0∂x ∂x ∂y ∂y 或 均质 (1-101)
∂∂H ∂∂H ∂H (5)地下水运动基本微分方程的统一形式:(F ) +(F ) +W =E ∂∂x ∂y 在承压含水层区∂y ∂t (1-102) T =KM ⎧x F =⎨
-Z ) 潜水含水层区 μKh *=K (H 在承压含水层区⎧⎩式中
E =⎨
潜水含水层区 ⎩μ
Z ——含水层底板标高。
1.7 数学模型的建立及求解
1.7.1数学模型的有关概念
同一形式的偏微分方程代表了整个一大类的地下水流的运动规律,而对于不同边界性质、不同边界形状的含水层,水头的分布是不同的。而且对于偏微分方程而言,方程本身并不包含反映特定渗流区条件的全部信息,方程可能存在无数个解,如需要从大量的可能解中求得与特定区域条件相对应的唯一特解,就必须提供反映特定区域特征的信息。这些信息包括:
*
T , μ, W ,当这些参数确定后,微分方程才能被确定下来。(1)微分方程中的有关参数
(2)渗流区范围和形状,当微分方程所对应的区域被确定之后才能对方程求解。 (3)边界条件(boundary conditions):表示渗流区边界所处的条件,用以表示水头H (或渗流量q )在渗流区边界上所应满足的条件,也就是渗流区内水流与其周围环境相互制约的关系。
(4)初始条件(initial conditions) :表示渗流区的初始状态,某一选定的初始时刻(t=0)渗流区内水头H 的分布情况。
将边界条件和初始条件并称为定解条件(definite solution condition) ,微分方程和定解条件一起构成渗流场的数学模型。
数学模型:描述某一研究区地下水流运动的数学方程与其定解条件共同构成的表示某一实际问题的数学结构。亦即从物理模型出发,用简洁的数学语言,即一组数学关系式来刻画它的数量关系和空间形式,从而反映所研究地质体的地质、水文地质条件和地下水运动的基本特征,达到复制或再现一个实际水流系统基本状态的目的的一种数学结构。
其中微分方程表示地下水的流动规律,定解条件表明研究对象所处的特定环境条件,即所研究的地下水流的真实状态。
定解问题是给定了方程(或方程组)和相应定解条件的数学物理问题。 建立模型是指建立数学模型的过程。 1.7.2 定解条件 1.7.2.1定解条件
定解条件指水头、流量等渗流运动要素在流场边界上的已知变化规律,这种变化规律是由流场外部条件引起的,但它不断地影响流场内部的渗流过程,并在整个期间一直起作用。
定解条件包括边界条件和初始条件。 1.7.2.2边界条件
边界条件是渗流区边界所处的条件,用以表示水头H (或渗流量q )在渗流区边界上所应满足的条件,也就是渗流区内水流与其周围环境相互制约的关系。
(1) 第一类边界条件(Dirichlet条件) :如果在某一部分边界(设为S l 或Γ1) 上,各点在每一时刻的水头都是已知的,则这部分边界就称为第一类边界或给定水头的边界,表示为:
H (x , y , z , t ) |s 1=ϕ(x , y , z , t ),
或
(x , y , z ) ∈s 1
(1-103)
H (x , y , z , t ) |Γ1=ϕ(x , y , t ),
(x , y ) ∈Γ1
(1-103)
给定水头边界不一定就是定水头边界。
可以作为第一类边界条件来处理的情况:
① 河流或湖泊切割含水层,两者有直接水力联系时,这部分边界就可以作为第一类边界处理。此时,水头φ是一个由河湖水位的统计资料得到的关于t 的函数。但要注意,某些
河、湖底部及两侧沉积有一些粉砂、亚粘土和粘土,使地下水和地表水的直接水力联系受阻,就不能作为第一类边界条件来处理。
在自然界,这种情况很少见。就是附近有河流、湖泊,也不一定能处理为定水头边界,还要视河流、湖泊与地下水水力联系的情况,以及这些地表水体本身的径流特征而定。在没有充分依据的情况下,不要随意把某段边界确定为定水头边界,以免造成很大误差。
② 区域内部的抽水井、注水井或疏干巷道也可以作为给定水头的内边界来处理。此时,水头通常是按某种要求事先给定,例如给定抽水井的允许降深等。上面介绍的都只是给定水头的边界。注意,给定水头边界不一定是定水头边界。
③ 排泄地下水的溢出带、冲沟或排水渠的边界也可近似看作给定水头边界。 (2)第二类边界条件(Neumam条件) :当知道某一部分边界(设为S 2或Γ2) 单位面积(二维空间为单位宽度) 上流入(流出时用负值) 的流量q 时,称为第二类边界或给定流量的边界。相应的边界条件表示为:
K ∂H
∂n
∂H ∂n
S 2
=q 1(x , y , z , t ), (x , y , z , t ) ∈S 2
(1-105)
或
T
Γ2
=q 2(x , y , t ), (x , y ) ∈Γ2
(1-106)
式中,n 为边界S 2或Γ2的外法线方向。q 1和q 2则为已知函数,分别表示S 2上单位面积和Γ2上单位宽度的侧向补给量。
常见的这类边界条件: ① 隔水边界(流线、分水岭):
∂H
∂n
Γ
=0
(1-107)
② 抽水井或注水井:
T
∂H ∂n
Γw
=-
Q
2πr w (1-108)
③ 补给或排泄地下水的河渠边界上,如已知补给量。
∂H
(3)第三类边界条件:某边界上H 和∂n 的线性组合是已知的,即有:
∂H
+αH =β∂n (1-109)
又称混合边界条件,α, β为已知函数。
σ' =
边界为弱透水层(渗透系数为K 1,厚度或宽度为m 1),
m 1
K 1,有
K
∂H ∂n ∂H ∂n
=
K 1
(H n -H ) =q (x , y , z , t ) m 1
H n -H
=0σ' (1-110) H n -H
=0σ' (1-111)
在s 3上,
K -
S 3
∂H
在Γ3上,∂n
T
-M
Γ3
浸润曲线的边界条件:
K
∂H ∂n
=q
c 2
(1-112)
当浸润曲线下降时,从浸润曲线边界流入渗流区的单位面积流量q 为:
∂H *q =μcos θ
∂t (1-113) 式中,μ为给水度,θ为浸润曲线外法线与铅垂线间的夹角。
1.7.2.3初始条件
初始条件:某一选定的初始时刻(t=0)渗流区内水头H 的分布情况。
H (x , y , z , t ) |t =0=H 0(x , y , z ),
或
(x , y , z ) ∈D (1-114)
H (x , y , t ) |t =0=H 0(x , y ),
(x , y , z ) ∈D (1-115)
其中,H 0为D 上的已知函数。 1.7.3 渗流数学模型的分类
(1)线性、非线性模型
模型由线性方程所组成,称为线性模型,如均质各向同性承压二维流方程。模型由非线性方程所组成,称为非线性模型,如潜水模型方程。
(2)静态、动态模型
根据模型中未知变量与时间的关系进行划分,若未知变量与时间无关,如稳定流模型,称为静态模型,反之,则为动态模型。
(3)集中、分布参数模型
模型中不含有空间坐标变量的模型,称为集中参数模型,如抽水井流量与降深之间的经验公式。模型中含有空间坐标变量的模型,称为分布参数模型。
(4)确定性与随机性模型
确定性模型:数学模型中各变量之间有严格的数学关系的模型。 随机性模型:数学关系式中含有一个或多个随机变量的模型。 用确定性模型来描述实际地下水流时,如前述,必须具备下列条件:
① 有一个(或一组) 能描述这类地下水运动规律的偏微分方程;同时,确定了相应渗流区的范围、形状和方程中出现的各种参数值。
② 给出相应的定解条件。对所建立的模型进行检验,即把模型预测的结果与通过抽水试验或其它试验对含水层施加某种影响后所得到的实际观测结果或一个地区地下水动态长期观测资料进行比较,看两者是否一致。若不一致,就要对模型进行校正,即修正条件(1)和(2)直至满意拟合为止。这一步骤称为识别模型或校正模型。
经过校正后的模型,能代表所研究的地质体,或者说是实际水流系统的复制品了,因而可以根据需要,用这个模型进行计算或预测,例如预测矿床疏干时的涌水量及地下水污染情况预测等。
③ 解(即满足条件①和②的解) 是存在的(存在性) ;
④ 解是唯一的(唯一性) ;要求所提问题的解存在和唯一是不言而喻的。
⑤ 解对原始数据是连续依赖的(稳定性) 。即稳定性的要求,意味着当参数或定解条件发生微小变化时,所引起的解的变化也是很微小的。只有有了这条保证,当参数和定解条件的数据有某些误差时,所求得的解才能仍然接近于真解;否则,解是不可信的,并应该认为此时的数学模型是有毛病的。在实际工作中,原始数据有某种误差,在所难免,所以这个条件很重要。
适定问题(Well –posed problem )是指数学模型满足(1)解是存在的(存在性),(2)解是唯一的(唯一性),(3)解对原始数据是连续依赖的(稳定性)这三个条件的问题。只要有一条不满足就是不适定问题。
正问题是根据数学模型、给定的含水层水文地质参数和定解条件求解水头的问题,又称水头预报问题。
逆问题(inverse problem)是根据数学模型、动态观测资料或抽水试验资料反过来确定含水层水文地质参数的问题。 1.7.4 建立数学模型的基本要点
(1)确定研究区的范围及渗流区的边界;
(2)确定渗流区的水力特征(包括埋藏条件、渗流状态、介质特征); (3)确定渗流区的边界条件; (4)确定渗流区的源汇项; (5)选择微分方程; (6)确定渗流区的初始条件。 1.7.5渗流数学模型的解法
(1)解析法(analytic method)
用参数分析及积分变换等方法直接求解数学模型解的方法。其解为精确解,使用简单,但该方法存在一定的局限性,只适用于含水层几何形状规则、方程式简单、边界条件单一的情况。
解析解(Analytic solution)又称精确解,是用解析方法求解数学问题所得到的解析表达式。
(2)数值法(numerical method)
用数值方法(离散化方法)求解数学模型的方法,其解为近似解。该方法是求解大型地下水流问题的主要方法。它把整个渗流区分割成若干个形状规则的小单元,每个小单元近似处理成均质的,然后建立每个单元地下水流动的关系式。把形状不规则的、非均质的问题转化为形状规则的均质问题。根据研究需要,确定单元划分数量,对于非稳定流还要对时段进行划分。最后,把局部整合起来,加上定解条件。
数值解(Numerical solution)是 用数值方法求得的数值解,是一种近似解。 (3)模拟法:利用物理现象与水流的相似性,在实验室内采用模拟的方法求解。
目 录
第一章 渗流理论基础 .............................................................................................................................. 1
第一章 渗流理论基础
1.1 渗流的基本概念
1.1.1 多孔介质及其特性
1.1.1.1多孔介质的概念
多孔介质(Porous medium):地下水动力学中具有空隙的岩石。广义上包括孔隙介质、裂隙介质和岩溶不十分发育的由石灰岩和白云岩组成的介质,统称为多孔介质。 孔隙介质:含有孔隙的岩层,砂层、疏松砂岩等;
裂隙介质:含有裂隙的岩层,裂隙发育的花岗岩、石灰岩等。
1.1.1.2 多孔介质的性质
(1) 孔隙性:有效孔隙和死端孔隙。
孔隙度(Porosity )是多孔介质中孔隙体积与多孔介质总体积之比(符号为n ),可表示为小数或百分数,n=Vv/V。
有效孔隙(Effective pores )是多孔介质中相互连通的、不为结合水所占据的那一部分孔隙。
有效孔隙度(Effective Porosity)是多孔介质中有效孔隙体积与多孔介质总体积之比(符号为n e ),可表示为小数或百分数,n e =Ve /V。
死端孔隙(Dead-end pores )是多孔介质中一端与其它孔隙连通、另一端是封闭的孔隙。
(2) 连通性:封闭和畅通,有效和无效。
(3) 压缩性:固体颗粒和孔隙的压缩系数推导。
(4) 多相性:固、液、气三相可共存。其中固相的成为骨架,气相主要分布在非饱和带中,液相的地下水可以吸着水、薄膜水、毛管水和重力水等形式存在。
固相—骨架matrix
气相—空气,非饱和带中
液相—水:吸着水 Hygroscopic water
薄膜水 pellicular water
毛管水 capillary water
重力水 gravitational water
1.1.1.3多孔介质中的地下水运动
比较复杂,包括两大类,运动特点各不相同,分别满足于孔隙水和裂隙岩溶水的特点。
(1)第一类为地下水在多孔介质的孔隙或遍布于介质中的裂隙运动,具有统一的流场,运动方向基本一致;
(2)另一类为地下水沿大裂隙和管道的运动,方向没有规律,分属不同的地下水流动系统。
1.1.2渗透与渗流
1.1.2.1渗透
渗透是地下水在岩石空隙或多孔介质中的运动,这种运动是在弯曲的通道中,运动轨迹在各点处不等。为了研究地下水的整体运动特征,引入渗流的概念。
图1-2 岩石中的渗流
(a )实际渗透 (b)假想渗流
1.1.2.2渗流
渗流(seepage flow ):具有实际水流的运动特点(流量、水头、压力、渗透阻力),并连续充满整个含水层空间的一种虚拟水流;是用以代替真实地下水流的一种假想水流。其特点是:
(1)假想水流的性质与真实地下水流相同;
(2)充满含水层空隙空间和岩石颗粒所占据的空间;
(3)运动时所受的阻力与实际水流所受阻力相等;
(4)通过任一断面的流量及任一点的压力或水头与实际水流相同。
渗流场(flow domain):假想水流所占据的空间区域,包括空隙和岩石颗粒所占的全部空间。
1.1.2.3典型单元体
典型单元体(REV ,Representative Elementary V olume )又称代表性单元体,是渗流场中其物理量的平均值能够近似代替整个渗流场的特征值的代表性单元体积。
REV 具备两个性质:
(1) 其体积和面积,大于个别空隙而小于渗流场,其中的渗流可以从一点连续运动到另一点;
(2) 通过单元体的运动要素(流量Q 、水头h 、压力p 、实际水头受到的阻力R )与真实水流相等,运动要素是连续变化的。
REV 的作用:
(1) 把物理性质看作是坐标的函数,孔隙度n 、导水系数T 、给水度μ和渗透系数均连续。
(2) 渗流的要素可以微分、积分,可以用微分方程来描述渗流要素。
1.1.2.4渗流速度
(1)过水断面(Cross-sectional area)是渗流场中垂直于渗流方向的任意一个岩石截面,包括空隙面积(A v )和固体颗粒所占据的面积(As ) ,A= Av + As 。渗流平行流动时为平面,弯曲流动时为曲面。
图1-3 渗流过水断面
(a )实际渗透 (b)假想渗流
(2)渗流量(Seepage discharge)是单位时间内通过过水断面的水体积,用Q 表示,单位m 3/d。
(3)渗流速度(Specific discharge/seepage velocity)又称渗透速度、比流量,是渗流在过水断面上的平均流速。它不代表任何真实水流的速度,只是一种假想速度。它描述的是渗流具有的平均速度,是渗流场空间坐标的连续函数,是一个虚拟的矢量。单位m/d,表示为
υ=Q
A (1-1)
(4)实际平均流速(Mean actual velocity)是多孔介质中地下水通过空隙面积的平均速度;地下水流通过含水层过水断面的平均流速,其值等于流量除以过水断面上的空隙面积,量纲为L/T。记为。它描述地下水锋面在单位时间内运移的距离,是渗流场空间坐标的离散函数。表示为:
=Q
A ⋅n (1-1a )
υ与之间存在以下关系:
υ= n (1-2)
若确定渗流场中任一点的渗流速度,可以按以下方法进行讨论:
设以P 点为中心的REV 的平均渗流速度矢量为v ,令REV 的体积为∆V 0,其中空隙体积为n ∆V 0,在空隙中的不同地点,流速u 不同,将u 在全部空隙体积n ∆V 0中求积分,再除以REV 体积∆V 0,即为渗流速度,表示为
v =1⎰(∆V 0) udVv ∆V 0(1-3)
1⎰(∆V v ) 0udVv (∆V v ) 0(1-4) =
可得V= n
1.1.3地下水的水头与水力坡度
(1)地下水水头(hydraulic head):渗流场中任意一点的总水头近似等于测压水头(piezometric head),即
H =Z +P
γ(1-5)
通常称为渗流水头。
在水力学中定义总水头(total head):
u 2
H =Z ++γ2g (1-6) P
式中右端三项分别称为位头(potential head )、压头(pressure head) 和速头(velocity head) 。
总水头(Total head )为测压管水头和流速水头之和。
测压管水头(Piezometric head)为位置水头与压力水头之和,H n =Z +P
γ。
压力水头(pressure head):含水层中某点的压力水头(h )指以水柱高度表示的该点水的压强,量纲为L ,即:h =P/ γ,式中 P 为该点水的压强;γ为水的容重,γ=ρg 。
速度水头(velocity head):在含水层中的某点水所具有的动能转变为势能时所达到的高
u 2
h v =2g ,式中u 为地下水在该点流动的速度;g 为重力加速度。 度,量纲为L ,即
u 2
由于在地下水中水流的运动速度很小,故速头2g 可以忽略,所以h 近似等于H ,即
H ≈H n =Z +P
γ=Z +P
ρg (1-7)
意义:渗流场中任意一点的水头实际上反映该点单位质量液体具有的总机械能,地下水在运动过程中不断克服阻力,消耗总机械能,因此沿地下水流程,水头线是一条降落曲线。
(2) 水力坡度[水力梯度](hydraulic gradient):在渗流场中大小等于梯度值,方向沿等水头面的法线并指向水头下降方向的矢量,用J 表示。
J =-dH n dn (1-8)
n 式中——法线方向单位矢量。在空间直角坐标系中,其三个分量分别为:
J x =-∂H ∂H ∂H , J y =-, J z =-∂x ∂y ∂z (1-9)
(3)等水头面与等水头线
等水头面:渗流场中水头值相同的各点相互连接所形成的一个面。可以是平面也可为曲面。
等水头线(groundwater contour):等水头面与某一平面的交线。
等水头面上任意一条线上的水头都相等。等水头面(线)在渗流场中是连续的,不同大小的等水头面(线)不能相交。
1.1.4 地下水运动特征分类
(1)渗流运动要素(Seepage elements)是表征渗流运动特征的物理量,主要有渗流量Q 、渗流速度V 、压强P 、水头H 等。
地下水运动方向(Groundwater flow direction)为渗透流速矢量的方向。
(2)层流与紊流
层流(laminar flow):水流流束彼此不相混杂、运动迹线呈近似平行的流动。
紊流(turbulent flow):水流流束相互混杂、运动迹线呈不规则的流动。
图1-4 空隙岩石中地下水的层流和紊流
根据Reynolds number判别地下水流态,通常
R e =νd 0νd =γ(0. 75n +0. 23) γ2(1-10)
式中:ν—地下水的渗流速度;
d —含水层颗粒的平均粒径;
d 0—含水层颗粒的有效粒径;
ν—地下水的运动粘度(粘滞系数)。
通常,确定d 的方法有:(1)d=d10;(2)Collins(1961):d =K n ;(3)Ward(1964):
d K , 其中n 为孔隙度。
若ReRe临界,则地下水处于紊流状态,此时液体质点无秩序地相互混杂地流动。 Re 临界≈ 150~300。天然地下水多处于层流状态。
(2)稳定流与非稳定流
根据渗流运动要素是否与时间有关而进行的划分。
稳定流(steady flow ):渗流运动要素不随时间变化;在一定的观测时间内水头、渗流速度等渗透要素不随时间变化的地下水运动。
非稳定流(unsteady flow ):渗流运动要素随时间变化;水头、渗透速度等任一渗透要素随时间变化的地下水运动。
(3)一、二、三维流
根据渗流方向与所选坐标轴方向之间的关系来划分。
一维流运动:当地下水沿一个方向运动,将该方向取为坐标轴,此时地下水的渗透速度只有沿该坐标轴的方向有分速度,其余坐标轴方向的分速度为0。
一维流(one-dimensional flow ),也称单向运动,指渗流场中水头、流速等渗流要素仅随一个坐标变化的水流,其速度向量仅有一个分量、流线呈平行的水流。
图1-5 承压水的一维流动
二维流运动:若地下水的渗透速度沿两个坐标轴方向都有分速度,仅一个坐标轴方向的分速度为0。
二维流(two-dimensional flow ),也称平面运动,地下水的渗透流速沿空间二个坐标轴方向都有分速度、仅仅一个坐标轴方向的分速度为零的渗流;水头、流速等渗流要素随两个坐标变化的水流,其速度向量可分为两个分量,流线与某一固定平面呈平行的水流。
平面二维流(Two-dimensional flow in plane),由两个水平速度分量所组成的二维流。 剖面二维流(two-dimensional flow in section),由一个垂直速度分量和一个水平速度分量组成的二维流。
单宽流量(Discharge per unit width):渗流场中过水断面单位宽度的渗流量,等于总流量Q 与宽度B 之比。即
q=Q/B。 (1-11)
总渗流量Q 为单宽流量q 与宽度B 的乘积,Q=qB。
图1-6 渠道向河流渗漏的地下水二维流动
(a )平面图 (b)剖面图
三维流运动:地下水的渗透流速沿空间三个坐标轴的分量均不为0。
三维流(three-dimensional flow),也称空间运动,地下水的渗透流速沿空间三个坐标轴的分量均不等于零的渗流;水头、流速等渗流要素随空间三个坐标而变化的水流。
图1-7 河弯处潜水的三维流动
(a )平面图 (b)剖面图
图1-8 均质各向同性含水层中潜水井抽水时的地下水运动
(a ) 平面图 (b)剖面图
(b )
1.2 渗流基本定律
1.2.1 达西定律(线性渗透定律)
图1-9 Darcy 实验装置
(1)达西定律表达式
实验条件:定水头、定流量、均质砂。
此时地下水做一维均匀运动,渗流速度与水力坡度的大小和方向沿流程不变。
达西定律(1856年) 表达式:
Q =KAJ =KA H 1-H 2
L (1-12)
V =Q =KJ A (1-13)
其中Q ——渗透流量(出口处流量),亦即通过过水断面(砂柱各断面)A 的流量(m3/d);volumetric flow rate.
K ——多孔介质的渗透系数(m/d);
2A ——过水断面面积(m) ;cross-sectional area of flow.
H 1、H 2——上、下游过水断面的水头(m);
L ——渗透途径 (m);
J ——水力梯度(J = (H1-H 2)/L),等于两个计算断面之间的水头差除以渗透途径,亦即渗透路径中单位长度上的水头损失。
达西定律的微分形式:
v =KJ =-K dH
dn (1-14)
v x =-K
J =-dH dH dH , v y =-K , v z =-K dx dy dz (1-15) dH =-gradH dn
达西定律的矢量形式:
v =v x i +v y j +v z k (1-16)
(2)达西公式讨论
达西定律反映了能量转化与守恒。
V 与I 的一次方成正比;当K 一定时,当V 增大时,水头差增大,表明单位渗透途径上被转化成热能的机械能损失越多,即V 与机械能的损失成正比关系;当V 一定时,K 越小,水头差越大,即K 与机械能的损失成反比关系。
(3) 达西公式适用范围
Re
Re>10-100,层流,不适用,地下水流速增大,为过渡带,由粘滞力占优势的层流转变为以惯性力占优势的层流运动;
Re>100,紊流,不适用。
达西定律的下限:地下水在粘性土中运动时存在一个起始水力坡度J 0。当时及水力坡度J
图1-10 渗透系数与水力坡度的实验关系(J. Bear)
1.2.2渗透系数、渗透率与导水系数
(1)渗透系数(K )(hydraulic conductivity)
V=KI ,当I=1时,V=K,即K 在数值上等于渗流速度,具有速度的单位,它又可以称为水力传导系数,反映含水介质对渗流阻力大小的系数。常用单位:m/d,cm/s。
渗透系数是反映岩石透水性的指标,可以根据渗透系数的大小进行岩石透水性分级(表1-1)。
表1-1 岩石透水性划分
K 的影响因素:
① 岩石的性质:粒度、成分、颗粒排列、充填状况、裂隙性质及其发育程度等,空隙大小起主导作用;
② 流体的物理性质:容重、粘滞性等。 (2)渗透率(intrinsic permeability)
k ρg dH 达西定律可表示为: *v =-
μ
ds (1-17)
式中ρ——液体的密度(density ) ;
g ——重力加速度;
μ——动力粘度;P the dynamic viscosity of the fluid H *=z +
γ——测压水头;
k ——描述多孔介质本身的渗透性能的常数,表示介质能使流体通过其本身的性能,它
不随渗透液体的物理、力学性质而变化。表征岩层渗透性能的参数;渗透率只取决于岩石的性质,而与液体的性质无关,记为k 。单位为cm 2或D ,1D=9.8697×10-9cm 2。
2ρgd 2ρgb u =J u =J
32μ12μ管流公式:; 缝流公式:
nd 2ρg nb 2ρg v =J v =J
32μ; 裂隙介质:12μ 多孔介质:
nb 2ρg ρg nd 2ρg ρg
K ==k K ==k
12μμ 32μμ 或于是有:ρg g
K =k =k
μμ(1-18)
达西(D )的定义:当液体的动力粘滞度为0.001Pa·s ,压强差为101325Pa 的情况下,通过面积为1cm 2、长度为1cm 岩样的流量为1cm 3/s时岩样的渗透率,记为D 。
尺度效应是指渗透系数与试验范围有关,随着试验范围的增大而增大的现象,K=K(x)。亦即抽水时间t 长、降深s 大的群孔抽水试验所得K 较抽水时间t 短、降深s 小的抽水试验所得K 大。
(3)导水系数(transmisivity )
Q=KMBJ
Q=Q/B=KMJ=TJ (1-19)
式中T=KM,称为导水系数,反映含水层出水能力的水文地质参数,其物理意义是水力坡度为1时,通过整个含水层厚度上的单宽流量。它是定义在一维或二维流中的水文地质参数。
单位:m 2/d。
图1-11 导水系数的概念
1.2.3非线性运动方程
(1) Re>1-10,P. Forchheimer(1901)公式:
J =av +bv 2(1-20)
m
或J =av +bv (1-20a )
式中的a ,b 由实验确定的常数, 1.6≤m ≤2。
1(2)当a=0时,有Chezy 公式:
v =k c J 2(1-21)
gk g k (1-22) d 2
k =
360,其中d 2是颗粒直径。 式中
1.3 岩层透水特征及水流折射定律
1.3.1岩层透水特征分类
(1)均质与非均质
根据岩层透水性随空间坐标的变化情况划分,若渗流场中,任意点都具有相同的渗透系数,或渗透系数不随空间坐标的变化而变化,则该岩层是均质的,反之则为非均质。岩石的非均质分两类,一类是渐变的,另一类是突变的。
均质岩层(Homogeneous strata):渗流场中所有点都具有相同参数的岩层。
非均质含水层(inhomogeneous strata ):渗流场中所有点不都具有相同参数的岩层,渗透系数K=K(x,y,z),为坐标的函数。
非均质分为两类,即渐变的和突变的。 (2)各向同性与各向异性
根据岩层透水性与渗流方向的关系划分,若渗流场中,某一点的K 与渗流方向无关,则该岩层是各向同性的,反之则为各向异性。
各向同性岩层(Isotropic strata):渗流场中某一点的渗透系数不取决于方向,即不管渗流方向如何都具有相同渗透系数的岩层。
各向异性岩层(anisotropic strata ):渗流场中某一点的渗透系数取决于方向,渗透系数随渗流方向不同而不同的岩层。
v 0. 55 2
(3)Ward J (=υ+υ
1.3.2渗透系数张量
岩石的透水性是用渗透系数来衡量的。渗透系数实际上是个张量。
(1)对于各向同性介质,其中任一点的渗透系数值与渗流方向无关,是一个标量,水力坡度与渗流方向是一致的。此时,υ可以表示为如下表达式:
υx =-K
∂H ∂H ∂H
, υy =-K , υz =-K ∂x ∂y ∂z
(2)对于各向异性介质,K 与渗流方向有关,K 不再是标量,水力坡度与渗流方向一般是
不一致的。此时,υ可以表示为如下表达式:
∂H ∂H ∂H
-K xy -K xz ∂x ∂y ∂z ∂H ∂H ∂H
υy =-K yx -K yy -K yz
∂x ∂y ∂z ∂H ∂H ∂H
υz =-K zx -K zy -K zz
∂x ∂y ∂z (1-23)
υ=-K xx
即可写成
⎡K xx
⎢K =⎢K yx
⎢K zx ⎣
K xy K yy K zy
K xz ⎤⎥K yz ⎥K zz ⎥⎦ (1-24a )
⎡K xx
K =⎢
⎣K yx 在二维空间中,
即有 υ=K·J (1-25)
K xy ⎤
K yy ⎥⎦ (1-24b )
渗透系数是对称张量,即K xy =Kyx ,K xz =Kzx ,K yz =Kzy 。
在各向异性介质中,水力坡度与渗流方向不一致,但在三个方向上两者是平行的,而且这三个方向称为主方向。
主渗透系数(主值)是指沿主方向测得的渗透系数,用 K 1、K 2、K 3表示,有K xx =K1,K yy =K2,K zz =K3,此时:
⎡K 1
K =⎢⎢0
⎢⎣0
0K 2
0⎤0⎥⎥K 3⎥⎦
(1-26)
1.3.3 层状岩层的等效渗透系数
在自然界中很常见的非均质岩层多是由许多透水性各不相同的薄层相互交替组成的层状岩层。当每一分层的渗透系数K i 和厚度面M i 已知时,可求出平行于层面的渗透系数K p 和垂直于层面的渗透系数K v 。
(1)水流平行层面
图1-11 层状岩层中平行于层面的渗流
特点:水流为稳定流,岩层水平分布,各段流量之和等于各部分流量之和,且各段具有统一的水头,各段具有相同的水力坡度。
根据达西定律有:故
q =∑q i =∑K i M i
i =1
i =1
n n
∆H ∆H
q =K p M
L ,若把其视为整体时,有L ,
n
∆H ∆H K p M =∑K i M i
L L 。水平岩层的等效渗透系数为: i =1
K p =
∑M
i =1
n
i
K i
M
(1-27)
等效导水系数为 T p =∑T i =∑M i K i
i =1
i =1
n
n
(1-28)
垂直方向岩性渐变时,有
1-12 层状岩层中垂直于层面的渗流
1
K p =
M
M 0
⎰
M
K (z ) dz
T p =⎰K (z ) dz
(1-29)
(2)水流垂直层面
特点:水流垂直层面运动,每段水流具有相同的单宽流量,且每段水力坡度不同。
q i =K i b
由
M ∆H ∆H q M 1
, q 2=K 2b ∆H 1=, HH 2=b 2M 1M 2,由此推导出,b K 1K 2
q n M i ∆H =∆H 1+∆H 2+... +∆H n =∑
b i =1K i 依次类推,有
K v =
∑M
i =1
n
n
i
M i ∑i =1K i (1-30)
可见,K v 取决于K i 最小的分层(阻力最大),K i =0,则K v =0。 另外,总是有K p >K v 。 1.3.4突变界面的水流折射定律
根据水流连续性条件,当水流斜向由一种介质进入另一种介质时,会发生折射。 如图所示:水流由K 1介质进入 K 2介质中,二者交界面上某一点的渗流速度和水头在两介质中的值依次为V 1、V 2和H 1、H 2。对于界面上的任一点应满足以下条件:
⎧H =H 2⎨
⎩V 1n =V 2n (1-31)
由图中几何条件有
tg θ1=
v 1τv , tg θ2=2τv 2n v 2n ,则有
tg θ1v 1τ
=tg θ2v 2τ
∂H 1=
∂H 2
-K 2
∂x -K 1
∂H 1∂H 2
=∂x ∂x ,则得到水流折射定律(渗流折射时必须满足的方程)因为:
tg θ1K 1
=
tg θ2K 2(1-32)
图1-12 渗透水流的折射
讨论上式可以的出以下结论:
(1)若K 1=K2,则θ1=θ2,表明在均质介质中水流不发生折射。
(2)若K 1≠K 2,且K 1,K 2均不为0,若θ1=0,则θ2=0,表明水流垂直通过界面时水流不发生折射。
(3)若K 1≠K 2,且K 1,K 2均不为0,若θ1=90,则θ2=90,表明水流平行于界面
时水流不发生折射。
(4)当水流斜向通过界面时,介质的渗透系数越大,θ值也越大,流线也越靠近界面。介质相差越大,两角的差值也越大。
根据水流折射原理和达西定律,可以帮助分析流场的水动力条件的变化。
1.4 流网及其应用
1.4.1 流网的概念
(1)流网(flow net ):渗流场中由一组流线与由一组等势线(当容重不变时为一组等水头线)相交组成的网格。对各向同性介质组成正交网。
流线(Streamline )渗流场内处处与渗流速度矢量相切的曲线。
地下水动力学中流线的概念和水力学中的概念是完全一致的。流线应是一根处处和渗流速度矢量相切的曲线。因此,流线簇就代表渗流区内每一个点的水流方向。 1.4.2流函数方程
(1) 流线的方程
根据上述定义,没有水流穿越流线。现在来研究描述流线的方程式。如图,在任一流线上取任意两点M(x, y)和M' (x+dx, y+dy)。M 点的渗流速度矢量为v ,它与它的两个分量Vx ,Vy 构成一个三角形MAB 。自M' 点作垂线Mb ,并延长至a ,见图。
图1-13 流线
当M 与M' 无限逼近时,弧线MM ' 可用切线Ma 来代替,故有Mb= dx,ab=dy。因为∆MAB ≈∆Mab ,有以下等式成立: 流线方程
dx dy
=
v x v y (1-33a)
v x dy -v y dx =0
(1-33b)
M 和M’是任意流线上任选的两点。因此,上式对流线上的任一点都是正确的,可以把它看成是流线的方程,用它来描述流线。上面的流线方程无论对各向同性和各向异性介质都是适用的。在各向异性介质中,如果选取的坐标轴(直角坐标系) 的方向分别与渗透系数的主方向一致,则上式变为:
dy dx
=∂H ∂H K xx K yy
∂x ∂y
对于各向同性介质,则式中的K xx =K yy =K 。由于(1-33b )式只涉及一个点的水流情况,故也适用于非均质介质。
(2) 流函数方程
d ψ=
设有二元ψ函数Ψ(z,y),其全微分为:
∂ψ∂ψ
dx +dy ∂x ∂y 。若取这样一种函数,使
∂ψ∂ψ
=-v y , =-v x ∂x ∂y (1-34)
d ψ=
则
∂ψ∂ψ
dx +dy =v x dy -v y dx =0∂x ∂y (1-35)
对其积分得:ψ=常数。表明沿同一流线,函数ψ为常数,不同的流线则有不同的函数值。称函数ψ为流函数,又称Lagrange 流函数,量刚为[L 2T -1]。
(3)流函数的物理意义
在无限接近的两条流线ψ和ψ+d ψ上沿某等水头线取两个点a(x,y)和b(x+dx,y+dy)。自a 、b 分别做垂线和水平线,相交于c 。见图1-14。
通过两流线之间的单宽流量dq 可看作是通过ac 和bc 的流量的代数和。将渗流速度分解则有:dq =v x ac+vy bc ,但ac=dy,bc=-dx,所以有dq=vx dy-v y dx
把式(1-34)和(1-35)代入上式,则得到:
dq =
将(35)式在ψ1和ψ1区间积分得:
∂ψ∂ψ
dy +dx =d ψ∂y ∂x (1-36)
图1-14 流函数与流量之间的关系
q =⎰d ψ=ψ2-ψ1
ψ1
ψ2
(1-37)
由(1-37)可以得出:在平面运动中,两流线之间的单宽流量等于和这两条流线相应的流函数之差。在同一条流线上,d ψ=0,q=0,ψ=常数。
由达西定律和(1-48)式,有:
v x =-K
∂H ∂ψ∂H ∂ψ
=, v y =-K =-∂x ∂y ∂y ∂x (1-38)
将(1-38)中第一式对y 求导,第二式对 x求导,得到:
∂2H ∂2ψ∂2H ∂2ψ
-K =, -K =-2
∂x ∂y ∂y 2∂y ∂x ∂x ∂2ψ∂2ψ
=-2, 2
∂y ∂x 整理得:
满足该方程。
(4) 流函数的特性
① 对于一给定的流线,流函数是常数。不同的流线有不同的常数值。流函数决定于流线。ψ=c
∂2ψ∂2ψ
+=0∂y 2∂x 2
(1-39)
表明在均质各向同性介质中,流函数满足Laplace 方程,而在其他情况下,流函数均不
② 在平面运动中,两流线之间的单宽流量等于和这两条流线相应的流函数之差。q=ψ
1-ψ2
③ 在均质各向同性介质中,流函数满足Laplace 方程;而在其他情况下,流函数均不满足该方程。
④ 在非稳定流中,流线不断地变化,只能给出某一瞬时的流线图。只有对不可压缩的液体的稳定流动,流线才有实际意义。 1.4.2流网的性质
(1)在各向同性介质中,流线与等水头线处处垂直,流网为正交网格。 由(1-38)式,得:
-K
消去K ,得:
∂H ∂ψ∂H ∂ψ
=K
∂y ∂y ∂x ∂x (1-40)
∂H ∂ψ∂H ∂ψ
+=0
∂x ∂x ∂y ∂y (1-41)
gradH =∇H =
等水头线
∂H ∂H
i +j ∂x ∂y
grad ψ=∇ψ=
流线
式中i ,j ——单位矢量。
∂ψ∂ψi +j ∂x ∂y
∇H ⋅∇ψ=
∂H ∂ψ∂H ∂ψ
+=0
∂x ∂x ∂y ∂y (1-42)
∇ϕ⋅∇ψ=0(1-43)
在非均质各向同性介质中,上式亦成立。
(2)在均质各向同性介质中,流网中每一网格的边长比为常数。
dH =
∂H ∂H 1dx +dy =-υds ∂x ∂y K
∂H ∂H 1d ψ=dx +dy =υds
∂x ∂y K
ds dH =-K dl d ψ(1-44)
式中dl ——相邻流线的间距;
ds ——等势线的间距。
通常取ds/dl=1,流网为曲边正方形。
(3)若流网中各相邻流线的流函数差值相同,且每个网格的水头差值相等时,通过每个网格的流量相同。
∆q =KJ ∆l =K
∆H ∆l
∆l =K ∆H ∆s ∆s (1-45)
∆l
当∆s =1时,∆q =K ∆H (1-46)
式中∆s ——网格相邻两等势线间的平均长度;
∆l ——网格相邻两流线间的平均宽度。
若上下游总水头差H r =H1-H 2,则m 个水头带中每一网格的水头差为
∆H =
H r
m (1-47)
(4)若两个透水性不同的介质相邻时,在一个介质中为曲边正方形的流网,越过界面进入另一个介时则变成曲边矩形。
∆q =K 1
∆H ∆H
∆l 1=K 2∆l 2∆s 1∆s 2
∆l 1K 2∆l 2
=
∆s 1K 1∆s 2
∆l 1∆l 2∆l 1K 1
=1=≠1∆s ∆s ∆s K 12当K 1≠K 2,且1时,2。
1.4.3流网的绘制与应用 1.4.3.1 流网的绘制
可采用解析法、各种模型试验法、徒手绘渐进法绘制流网。 (1)确定边界条件
① 河渠的湿周为一条等水头线。 ② 平行于隔水边界可绘制流线。
③ 无入渗补给及蒸发排泄,有侧向补给,做稳定运动时,地下水面是一条流线。 ④ 有入渗补给时,地下水面既不是流线也不是等水头线。
(2)根据已经确定的边界条件,根据流网的性质可以判定另一条件,作出流网。 (3)根据流网的性质绘制,各向同性含水层中,流线与等水头线处处正交,网格边长比为常数。
(4)在同一渗流区内,除奇点外,流线与等水头线各自不能相交;如遇透水性大的透
镜体时,则流线向该点汇集,反之则绕行。流线穿越突变界面时,应用水流折射定律绘制。 1.4.3.2流网的应用
(1)定量计算渗流区中的渗流运动要素
①水头H 、渗透压强P
P
γ
=H ±Z , P =γ(H ±Z )
(1-48)
②水力梯度J 、渗流速度v
J =
∆H
, ∆s
v =KJ
(1-49)
③流量q
q =K ∆Hn =
Kn (H 1-H 2)
m (1-50)
式中m 、n ——水头带数目、流带的数目。
(2)定性分析渗流区的水文地质条件及其变化。 (3)主要用于解决稳定渗流问题。
1.5 渗流连续方程
1.5.1含水层的状态方程
含水层的状态方程主要包括地下水的状态方程和多孔介质的状态方程。 1.5.1.1 地下水的状态方程
Hooke 定律:
E =-V
dp dV
式中E ——体积弹性系数(体积弹性模量),20℃时,E=2.1×105N/cm2。其倒数为压缩系数。
等温条件下,水的压缩系数(coef. of compressibility)为
β=-
1dV
V dp
(1-51)
积分(p →p 0,V →V 0)改写得:
-β(p -p 0)
体积:V =V 0e
(1-52)
β(p -p )
ρ=ρe 0密度:
(1-53)
(1-54)
按Taylor 级数展开,得到近似方程: V =V 0[1-β(p -p 0)] 和
ρ=ρ0[1+β(p -p 0)]
(1-55)
因 d (ρV ) =ρdV +Vd ρ=0(质量守恒),故有
d ρ=-ρ
dV
=ρβdp V
(1-56)
1.5.1.2 多孔介质的状态方程
多孔介质压缩系数(Coefficient of compressibility)表示多孔介质在压强变化时的压缩性的指标,用α 表示。
多孔介质压缩系数α的表达式为
α=-
1dV b
V b d δ(1-57)
式中,V b =V s +Vv ——多孔介质中所取单元体的总体积;
V s ——单元体中固体骨架(solid matrix)体积; V v ——为其中的孔隙(voids )体积。 δ——介质表面压强; V v =nV b ;V s =(1-n)Vb
α=-
1dV s 1dV v 1-n dV s n dV v
-=---
V b d δV b d δV s d δV v d δ
(1-58)
(1-59)
α=(1-n ) αs +n αp
αs =-
式中
1dV s
V s d δ——多孔介质固体颗粒压缩系数,表示多孔介质中固体颗粒本身的压缩性
的指标, αs
αp =-
1dV v
V v d δ——多孔介质中孔隙压缩系数(Compressibility of the pores of a porous
medium ),表示多孔介质中孔隙的压缩性的指标。
n ——多孔介质的孔隙度。
因(1-n ) αs
考虑承压含水层受力情况,取一水平横截面AB ,按Terzaghi (1883~1963)观点:
σ=λσs +(1-λ) p
(1-60)
式中σ——上覆荷重引起的总应力(total stress);
σs ——作用在固体颗粒上的粒间应力 (intergranular stress);
λ——横截面面积中颗粒与颗粒接触面积所占的水平面积比;
p ——水的压强。
图1—1 一个可压缩的承压含水层(J. Bear)
Terzaghi 令λσs =σ' ,称为有效应力(effective stress)。λ很小,(1-λ)p ≈ p,因此有
σ=σ' +p
(1-61)
在水位下降为∆H 时,有σ=(p -γ∆H ) +(σ' +γ∆H ) 。即作用于固体骨架上的力增加了γ∆H 。作用于骨架上力的增加会引起含水层的压缩,而水压力的减少将导致水的膨胀。含水层本来就充满了水,骨架的压缩和水的膨胀都会引起水从含水层中释出,前者就象用手挤压充满了水的海绵会挤出水—样。
dV b dn
=
1-n 。 因Vs=constant,故V b
dV b d (∆z ) =
dz ,故 只在垂直方向上有压缩,V b
d (∆z ) =∆z αdp
dn =(1-n ) αdp
(1-62) (1-63)
上两式表示垂直厚度变化、孔隙度变化与水的压强变化的关系。
为了讨论水头降低时含水层释出水的特征,取面积为1 m2、厚度为l m (即体积为l m3) 的含水层,考察当水头下降1m 时释放的水量。此时,有效应力增加了γ∆H =ρg ×1 = ρg 。介质压缩体积减少所释放出的水量(dV b )与水体积膨胀所释放出的水量(dV )之和为
-dV b =αV b dp =α⋅1⋅ρg =αρg
dV =-βdp =-βn (-ρg ) =n βρg
μs =-dV b +dV =αρg +n βρg
或μs =ρg (α+n β)
(1-64)
式中μs ——贮水率[释水率](specific storativity),量纲[L-1],为弹性释水[贮水]。
μ* = μs M
式中M ——含水层厚度(m);
μ*贮水系数(storativity )。
贮水系数μ* 和贮水率 μs 都是表示含水层弹性释水能力的参数,在地下水动力学计算中具有重要的意义。
贮水率
μs =(α+n β) ρg ,表示含水层水头变化一个单位时,从单位体积含水层中,
因水体积膨胀(压缩)以及骨架的压缩(或伸长)而释放(或储存)的弹性水量。单位1/L。
贮水率是描述地下水三维非稳定流或剖面二维流中的水文地质参数,既适用于承压水也适用于潜水。对于平面二维非稳定流地下水运动,当研究整个含水层厚度上的释水情况时,用贮水系数来体现。
贮水系数又称释水系数或储水系数,为含水层水头变化一个单位时,从底面积为一个单位,高度等于含水层厚度的柱体中所释放(或贮存)的水量;指面积为一个单位、厚度为含水层全厚度M 的含水层柱体中,当水头改变一个单位时弹性释放或贮存的水量,无量纲。既适用于承压含水层,也适用于潜水含水层。
⎧μs M =μ*.......... (承压含水层)
E =⎨
⎩μs (H -Z ) =μ...(潜水含水层)
μ*范围值:n ×10-3~ n×10-5;μ范围值:0.05~0.3。实际测出的值往往小于理论值。
弹性释水与重力给水:对于含水层而言,由于受埋藏条件的限制,抽水时,水的给出存在着不同。
潜水含水层在抽水过程中,大部分水在重力作用下排出,疏干作用于水位变动带(饱水带)和包气带两部分,由于包气带的存在,使得饱水带中水的释放存在延滞和滞后现象。当水头下降时,可引起二部分水的排出。在上部潜水面下降部位引起重力排水,用给水度μ 表示重力排水的能力;在下部饱水部分则引起弹性释水,用贮水率μ* 表示这一部分的释水能力。必须区分两者之间的不同,潜水含水层还存在滞后疏干现象。
承压含水层抽水时,水的释放是由于压力减少造成的,这一过程是瞬时完成的。只要水头下降不低到隔水顶板以下,水头降低只引起含水层的弹性释水,可用贮水系数μ* 表示这种释水的能力。 1.5.1.4导压系数
描述含水层水头变化的传导速度的参数,其数值等于含水层的导水系数与贮水系数之比或渗透系数与贮水率之比。
a =
T K =
E μs ,量刚为L 2T -1。
1.5.2 渗流连续方程
由于渗流场中各点的渗流速度大小、方向都不同,为了反映液体运动的质量守恒关系,需要在三维空间中建立微分方程形式表达的连续性方程。
在渗流场中任意取一点P(x, y, z),以P 为中心沿直角坐标轴取一微小的六面体,体积为
∆x ∆y ∆z ,称为特征单元体,设单元体无限小,但保证单元体穿过介质骨架和空隙。
图1-15 渗流区中的单元体
设v x ,v y ,v z 分别为该点在X 、Y 、Z 方向上的渗流速度。Abcd 面中点沿X 轴方向流入:
P 1=(x -
∆x
, y , z ) 2。
M xi =ρυx 1=ρυx (x -
流出:
M xo =ρυx 2=ρυx (x +
∆x
, y , z ) ∆y ∆z ∆t 2;
∆x
, y , z ) ∆y ∆z ∆t 2。
利用Taylor 级数展开,略去二阶导数以上的高次项,有:
∆M x =M xi -M xo =[ρυx -
1∂(ρυx ) 1∂(ρυx )
∆x ]∆y ∆z ∆t -[ρυx +∆x ]∆y ∆z ∆t
2∂x 2∂x
∂(ρυx )
∆x ∆y ∆z ∆t ∂x = -
∂(υz ρ) ∂(υx ρ)
∆M =-∆x ∆y ∆z ∆t ∆M y =-∆x ∆y ∆z ∆t z
∂z ∂x 同理 ;
单元体本身水质量在Δt 时间内的变化量由质量守恒定律,得到渗流的连续性方程:
∆M =
∂
(ρn ∆x ∆y ∆z ) ∆t ∂t ,ρ为液体密度。
∆M x +∆M y +∆M z =∆M
⎡∂(ρυx ) ∂(ρυy ) ∂(ρυz ) ⎤∂
-⎢++∆x ∆y ∆z =(ρn ∆x ∆y ∆z ) ∆t ⎥
∂x ∂y ∂z ∂t ⎣⎦(1-65)
或
-div (ρυ) ∆x ∆y ∆z =
∂
(ρn ∆x ∆y ∆z ) ∂t
上式即为非稳定流的渗流连续方程,表明渗流场中任意体积含水层流入、流出该体积含水层中水质量之差永远恒等于该体积中水质量的变化量。它表达了渗流区内任何一个“局部”所必须满足的质量守恒定律。
若把含水层看作刚体,ρ=constant,n 不变,即水和介质没有弹性变形或渗流为稳定流,
则渗流连续性方程为
div (υ) =
∂v x ∂v y ∂v z
++=0∂x ∂y ∂z (1-66)
此式表明,在同一时间内流入单元体的水体积等于流出的水体积,即体积守恒。
连续性方程是研究地下水运动的基本方程,各种研究地下水运动的微分方程都是根据连续性方程和反映质量守恒定律的方程建立起来的。
1.6 渗流基本微分方程
1.6.1承压水运动的基本微分方程
基本假设:
(1)单元体体积无限小,为承压含水层;
(2)含水层侧向受到限制,∆x 、∆y 为常量,∆z 为变量,存在垂向压缩,水的密度ρ、孔隙度n 和∆z 随压力p 而变化;
(3)由ρ引起的变化
(υx
∂ρ∂ρ∂ρ+υy +υz ) ∆x ∆y ∆z ∂x ∂y ∂z ,远小于单元体内液体质量的
∂ρ
∂c ,可忽略不计;
变化量(含ρ)
υgrad υ
(4)水流服从Darcy’s Law;
(5)K 不因ρ=ρ(p ) 的变化而变化;
(6)μs 和K 也不受n 变化(由于骨架变形)的影响。 流体的质量:M =ρV =ρn ∆x ∆y ∆z
由于含水层的侧向受到限制,可假设∆x 、∆y 为常量,只考虑垂向压缩。于是,只有水的密度ρ、孔隙度n 和单元体高度∆z 三个量随压力而变化,于是有:
⎡∂(∆z ) ∂p ∂n ∂p ∂ρ∂p ⎤∂(ρn ∆z ) ∂M ∂
n ρ+ρ∆z +n ∆z =(ρn ∆x ∆y ∆z ) =∆x ∆y ⎢⎥∆x ∆y ∂p ∂t ∂p ∂t ∂p ∂t ⎦∂t ∂t ∂t =⎣
⎡∂(∆z ) ∂n ∂ρ⎤∂p
n ρ+ρ∆z +n ∆z ⎢⎥∂t ∆x ∆y ∂p ∂p ∂p ⎦=⎣ (1-67)
d ρ
=βρ, dp 由含水层状态方程,
d (∆z ) dn p
=α(1-n ), =α∆z H =Z +dp dp γ,所以有,因为
p =γ(H -Z ) ,Z 为定值,则dp =γdh ,则可得到:
∂M ∂H ∂H
=[n ρα∆z +ρ∆z α(1-n ) +n ∆z βρ]γ∆x ∆y =(α+n β) ρ2g ∆x ∆y ∆z ∂t ∂t ∂t (1-68)
于是连续性方程(1-65)变为:
⎡∂(ρυx ) ∂(ρυy ) ∂(ρυz ) ⎤2∂H -⎢++∆x ∆y ∆z =(α+n β) ρg ∆x ∆y ∆z ⎥
∂y ∂z ⎦∂t ⎣∂x (1-69)
∂p ∂H ∂ρ∂H p ∂ρ
=ρg +(H -Z ) g =ρg +
∂t ∂t ∂t ρ∂t 又 ∂t
∂p ρg ∂H ∂H
=≈ρg
∂t (1-70) 则 ∂t 1-βp ∂t
∂M ∂H
=ρμs ∆x ∆y ∆z
∂t 令 μs =(α+n β) ρg ,则 ∂t
根据连续性原理有:
⎡∂(ρυx ) ∂(ρυy ) ∂(ρυz ) ⎤∂M
=-⎢++⎥∆x ∆y ∆z ∂t ∂x ∂y ∂z ⎣⎦
⎡∂(ρυx ) ∂(ρυy ) ∂(ρυz ) ⎤∂M ∂H
=ρμs ∆x ∆y ∆z =-⎢++⎥∆x ∆y ∆z ∂t ∂t ∂x ∂y ∂z ⎣⎦则有:
⎡∂(ρυx ) ∂(ρυy ) ∂(ρυz ) ⎤∂H -⎢++∆x ∆y ∆z =ρμ∆x ∆y ∆z ⎥s
∂x ∂y ∂z ∂t ⎦即:⎣
将
υx =-K
∂h ∂h ∂h
, υy =-K , υz =-K ∂x ∂y ∂z 代入,整理得:
ρ⎢
⎡∂∂H ∂∂H ∂∂H ⎤∂H (K ) +(K ) +(K ) ⎥∆x ∆y ∆z =ρμs ∆x ∆y ∆z
∂x ∂y ∂y ∂z ∂z ⎦∂t ⎣∂x ,所以有 ∂∂H ∂∂H ∂∂H ∂H (K ) +(K ) +(K ) =μs ∂x ∂x ∂y ∂y ∂z ∂z ∂t (1-71)
上式为三维流微分方程,也可写成:
div K grad (h ) =μs
]
∂H
∂t
物理意义:渗流空间内任一单位体积含水层在单位时间内流入与流出该体积含水层中的弹性水量的变化量,即单位体积含水层的水量均衡方程。
基本微分方程是研究承压含水层中地下水运动的基础。它反映了承压含水层中地下水运动的质量守恒关系,表明单位时间内流入、流出单位体积含水层的水量差(左端)等于同一时间内单位体积含水层弹件释放(或弹性贮存) 的水量(右端)。它还通过应用Darcy 定律反映了地下水运动中的能量守恒与转化关系。可见,基本微分方程表达了渗流区中任何一个“局部”都必须满足质量守恒和能量守恒定律。
数学意义:表示渗流空间内任一点任一时刻的渗流规律。
在柱坐标系中,有
1∂∂H 1∂2H ∂2H μs ∂H
(r ) +2+2=r ∂r ∂r K ∂t (1-72a ) r ∂θ2∂z
∂2H 1∂H 1∂2H ∂2H μs ∂H ++2+2=22r ∂r r ∂θK ∂t (1-72b ) ∂z 或 ∂r
由地下水流基本微分方程(1-71),在均质各向同性介质中,方程简化为:
∂2H ∂2H ∂2H μs ∂H
+2+2=
K ∂t (1-73) ∂x 2∂y ∂z
对于各向异性介质,若把坐标轴方向和各向异性介质的主方向定为一致,则有
∂∂H ∂∂H ∂∂H ∂H
(K xx ) +(K yy ) +(K zz ) =μs ∂x ∂x ∂y ∂y ∂z ∂z ∂t (1-74)
在二维流情况下,基本微分方程可表示为:
∂∂H ∂∂H ∂H (T ) +(T ) =μ*∂x ∂x ∂y ∂y ∂t (1-75)
上式即为承压水平面二维流微分方程,该方程是研究承压水含水层中地下水运动的基础,反映了承压水含水层中地下水运动的质量守恒关系,表明单位时间流入、流出单位体积含水层的水量差等于同一时间内单位体积含水层弹性释放(或贮存)的水量。
在实际渗流问题中若存在抽、注水及越流影响,只要在微分方程中的左端中通过加、减W 项,通常把该项称为源汇项。所谓的源项表示在垂直方向上有水流入含水层,此时W 为正;汇指在垂直方向上有水流出含水层,此时W 为负。
此时(1-71)式变成:
∂∂H ∂∂H ∂∂H ∂H (K ) +(K ) +(K ) +w =μs ∂x ∂x ∂y ∂y ∂z ∂z ∂t (1-76)
二维流情况下:
∂∂H ∂∂H ∂H (T ) +(T ) +w =μ*∂x ∂x ∂y ∂y ∂t (1-77)
*
a =T /μ在二维流情况下,令压力传导系数(导压系数),则均质各向同性含水层
基本微分方程为:
∂2H ∂2H 1∂H
++w =
a ∂t (1-78) ∂x 2∂y 2
非均质各向同性含水层中的稳定流运动:
∂∂H ∂∂H ∂∂H (K ) +(K ) +(K ) =0∂x ∂x ∂y ∂y ∂z ∂z (1-79)
均质各向同性含水层中的稳定流运动:
∂2H ∂2H ∂2H
+2+2=02
∂x ∂y ∂z (1-80)
上式也称Laplace 方程。稳定运动方程的右端都等于零,意味着同一时间内流入单元体的水量等于流出的水量。这个结论不仅适用于承压含水层,也适用于潜水含水层和越流含水层。
1.6.2 越流含水层中地下水运动的基本微分方程
在自然界中,存在以下情况,承压含水层的上、下岩层并不是绝对隔水的,其中一个或两个可能是弱透水层。虽然含水层会通过弱透水层和相邻含水层发生水力联系,但它还是处于承压状态,将其称为半承压含水层。当该含水层和相邻含水层间存在水头差时,地下水就会从高水头含水层通过弱透水层流向低水头含水层。这种现象称为越流。半承压含水层称为越流含水层。
假设:主含水层渗透系数K 远远大于若透水层的渗透系数K 1;主含水层弹性释放的水量、弱透水层的越流量远远大于弱透水层弹性释放的水量。主含水层中的水流近似地看作二维流
1M
H =H(x,y, t) =H(x,y, z, t)dz ⎰0M 问题, 。
对于均衡单元体,根据水均衡原理可以写出下列形式的连续性方程:
[(Q x -
∂Q y ∆y ∂Q y ∆y ∂Q x ∆x ∂Q ∆x
) -(Q x +x )]∆t +[(Q y -) -(Q y +)]∆t ∂x 2∂x 2∂y 2∂y 2
+(v 2-v 1) ∆x ∆y ∆t =μ*
∂H
∆x ∆y ∆t ∂t (1-81)
式中,v 1,v 2分别为通过上部和下部弱透水层的垂直越流速率或越流强度,即
v 1=-K 1
∂H 1H -H 1∂H 2H -H 1
=K 1, v 2=-K 2=K 22∂z m 1∂z m 2(1-82)
其中,H 1(x, y, t)和H 2(x, y, t)分别为上含水层和下含水层中的水头,如T 表示主含水层的导水系数,则得到不考虑弱透水层弹性释水条件下非均质各向同性越流含水层中非稳定运动的基本微分方程:
H -H 1H -H ∂∂H ∂∂H ∂H
(T ) +(T ) +K 1+K 22=μ*∂x ∂x ∂y ∂y m 1m 2∂t (1-83)
对于均质各向同性介质来说,有:
∂2H ∂2H H -H 1H 2-H μ*∂H
+++=2222
T ∂t (1-84) ∂x ∂y B 1B 2
B 1=
式中
Tm 1
, K 1
B 2=
Tm 2
K 2(1-85)
分别称为上、下两个弱透水层的越流因素。
越流因素B 的量纲为[L]。弱透水层的渗透性愈小,厚度越大,则B 越大,越流量越小。在自然界中,越流因素值的变化很大,可以从几米到若干公里。对于一个完全隔水的覆盖层来说,B 为无穷大。
另一个反映越流能力的参数是越流系数σ' 。其定义为:当主含水层和供给越流的含水层间的水头差为一个长度单位时,通过主含水层和弱透水层间单位面积界面上的水流量。因此,
σ' =
K 1
m 1(1-86)
K 1、m 1分别为弱透水层的渗透系数和厚度。σ' 越大,相同水头差下的越流量越多。 1.6.3 潜水运动的基本微分方程 1.6.3.1Dupuit 假设
图1-16 Dupuit假设
在潜水面上任意取一点P ,有:
dH dz J =-=-=-sin θ
ds ds (1-87)
该点的流速v 方向与潜水面相切,则由达西定律有:v s =-KJ=-Ksin θ。
当θ很小时,tg θ=sinθ。此时,(1)潜水面比较平缓,等水头面呈铅直,水流基本水平,可忽略渗流速度的垂直分量v Z ;(2)隔水底板水平,铅垂剖面上各点的水头都相等,各点的水力坡度和渗流速度都相等,H (x , y , z , t ) 可以近似地用H (x , y , t ) 代替,此即著名的Dupuit 假设。
渗流速度:
v x =-K
dH
dx ,H=H(x) (1-88)
通过宽度B 的铅直平面的流量为
Q x =-Kh
dH
dx , H=H(x) (1-89)
式中Q x ——x 方向的流量;
h ——潜水含水层厚度;h=H(隔水层水平时)。 对于更一般情况,H=H(x,y) 有:
V x =-K
∂h ∂H
, V y =-K , ∂x ∂y (1-90) ∂H ∂H , Q y =-KhB ∂x ∂y (1-91)
则得:
Q x =-KhB
由于Dupuit 假设的引入,将垂直方向的水流速度忽略,减少了z 变量,简化了计算,但会
2
产生一定的误差,经验证明当i
h 2⎡h 2⎤-⎢h -⎥2⎣2⎦i 2dh 0≤
dx h 1+i
2(1-92)
表1- 用tg θ 代替sin θ的误差
θ(ο) 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22
tg θ 0.01746 0.03492 0.05241 0.06993 0.08749 0.10510 0.12278 0.14054 0.15838 0.17633 0.19438 0.21256 0.23087 0.24933 0.26795 0.28675 0.30573 0.32492 0.34433 0.36397 0.38386 0.40403
sin θ 0.01745 0.03490 0.05234 0.06976 0.08716 0.10453 0.12187 0.13917 0.15643 0.17365 0.19081 0.20791 0.22495 0.24192 0.25882 0.27564 0.29237 0.30902 0.32557 0.34202 0.35837 0.37461
D -0.0000027 -0.0000213 -0.0000718 -0.000170 -0.000333 -0.000576 -0.000915 -0.001368 -0.001950 -0.002679 -0.003571 -0.004645 -0.005917 -0.007406 -0.009130 -0.011108 -0.013359 -0.015903 -0.018759 -0.021950 -0.025496 -0.029420
% -0.01523 -0.06095 -0.1372 -0.2442 -0.3820 -0.5508 -0.7510 -0.9828 -1.247 -1.543 -1.872 -2.234 -2.630 -3.061 -3.528 -4.030 -4.569 -5.146 -5.762 -6.418 -7.114 -7.853
23 24 25 26 27 28 29 30
0.42447 0.44523 0.46631 0.48773 0.50953 0.53171 0.55431 0.57735
0.39073 0.40674 0.42262 0.43837 0.45399 0.46947 0.48481 0.50000
-0.033744 -0.038492 -0.043689 -0.049361 -0.055535 -0.062238 -0.069499 -0.077350
-8.636 -9.464 -10.338 -11.260 -12.233 -13.257 -14.335 -15.470
对于各向异性介质,上式中的i 代以(K xx /K zz ) i 。
2
2
Dupuit 假设无效的地区: (1)存在入渗的潜水分水岭地段;
(2)渗出面附近。渗出面(seepage surface) 是在下游边界面上,潜水面以下、下游水面以上的地段。渗出面上潜水面往往和边界面相切,有较大的垂向分速度。 (3)垂直的隔水边界附近。
图1-17 Dupuit假设无效的地区
1.63.2 Boussinesq方程
图1-18 潜水的非稳定运动
根据Dupuit 假设,可建立有关潜水含水层中的地下水流方程。
(1)潜水一维流方程(沿x 方向运动)
在Δt 时间内,上、下游流入、流出单元体的水量差为:
(q -
∂(υx H ) ∂q ∆x ∂q ∆x ∂q
⋅) ∆t -(q +⋅) ∆t =-⋅∆x ∆t =-∆x ∆t ∂x 2∂x 2∂x ∂x
在该段时间内,垂直方向的补给量为W ∆x ∆t ,故Δt 时段水量总的变化量为
⎡∂(υx H ) ⎤
-+W ⎢⎥∆x ∆t ∂x ⎣⎦。
∂H
由于水量的变化引起潜水面的升降,设其变化的速率为∂t ,则在Δt 时段,由于潜水
面的变化而引起的小土体内的水体积的增量为
μ
∂H
∆x ∆t ∂t ,则有
⎡∂(υx h ) ⎤dh -+W ∆x ∆t =μ∂H ∆x ∆t υx =-K ⎢⎥∂x ⎣⎦dx (均质各向同性)代入上式,可以得∂t 。将
到有入渗补给的潜水含水层中地下水非稳定运动的一维流方程,又称为Boussinesq 方程:
∂∂H W μ∂H
(h ) +=∂x ∂x K K ∂t
(1-93)
式中K 、μ—— 潜水含水层的渗透系数、给水度;
W ——含水层单位时间、单位面积上的垂向补排量,补给为正,排泄为负。 (2)潜水二维流方程
均质各向同性含水层,Boussinesq 方程为:
∂∂H ∂∂H W μ∂H
(h ) +(h ) +=∂x ∂x ∂y ∂y K K ∂t (1-94)
式中h=H-Z,Z=0时,h=H。
非均质含水层,Boussinesq 方程为:
∂∂H ∂∂H ∂H (Kh ) +(Kh ) +W =μ∂x ∂x ∂y ∂y ∂t (1-95)
在推导潜水基本微分方程时应用了Dupuit 假设,忽略了弹性储存,所选的单元体是一个包括了整个含水层厚度在内的土柱,这与承压水非稳定运动时选取的无限小的单元体不同。所以,应用潜水运动基本方程得到的H(x,y ,t) 只能代表该点整个含水层厚度上平均水头的近似值,不能用来计算同一垂直剖面上不同点的水头变化。
(3)潜水三维流方程
若不用Dupuit 假设,Boussinesq 方程的一般形式:
∂∂H ∂∂H ∂∂H ∂H (K ) +(K ) +(K ) =μs ∂x ∂x ∂y ∂y ∂z ∂z ∂t (1-96)
在上面的潜水基本运动微分方程中右端项为贮水率而不是给水度,其原因在于,当不考虑Dupuit 假设时,单元体位于渗流区内部,其贮存量的变化只能是弹性释水而不是疏干排水,因此推导出的潜水非稳定运动方程和承压水非稳定运动方程形式一样。在这种情况下,地下水非稳定运动的特点由边界条件来反映。
对于各向异性介质,坐标轴方向同主方向,有]
∂∂H ∂∂H ∂∂H ∂H (K xx ) +(K yy ) +(K zz ) =μs ∂x ∂x ∂y ∂y ∂z ∂z ∂t (1-97)
假设固体骨架是不可压缩的,μs =0,同时假设忽略水的压缩性, 即ρ=常数,有:
∂∂H ∂∂H ∂∂H (K ) +(K ) +(K ) =0∂x ∂x ∂y ∂y ∂z ∂z (1-98)
∂∂H ∂∂H ∂∂H (K xx ) +(K yy ) +(K zz ) =0∂x ∂x ∂y ∂y ∂z ∂z 或 (1-99)
(4)潜水稳定运动的微分方程
没有入渗和蒸发时,潜水稳定运动的方程式为:
∂∂H ∂∂H
(Kh ) +(Kh ) =0
∂x ∂y ∂y 非均质 ∂x (1-100)
∂∂H ∂∂H
(h ) +(h ) =0∂x ∂x ∂y ∂y 或 均质 (1-101)
∂∂H ∂∂H ∂H (5)地下水运动基本微分方程的统一形式:(F ) +(F ) +W =E ∂∂x ∂y 在承压含水层区∂y ∂t (1-102) T =KM ⎧x F =⎨
-Z ) 潜水含水层区 μKh *=K (H 在承压含水层区⎧⎩式中
E =⎨
潜水含水层区 ⎩μ
Z ——含水层底板标高。
1.7 数学模型的建立及求解
1.7.1数学模型的有关概念
同一形式的偏微分方程代表了整个一大类的地下水流的运动规律,而对于不同边界性质、不同边界形状的含水层,水头的分布是不同的。而且对于偏微分方程而言,方程本身并不包含反映特定渗流区条件的全部信息,方程可能存在无数个解,如需要从大量的可能解中求得与特定区域条件相对应的唯一特解,就必须提供反映特定区域特征的信息。这些信息包括:
*
T , μ, W ,当这些参数确定后,微分方程才能被确定下来。(1)微分方程中的有关参数
(2)渗流区范围和形状,当微分方程所对应的区域被确定之后才能对方程求解。 (3)边界条件(boundary conditions):表示渗流区边界所处的条件,用以表示水头H (或渗流量q )在渗流区边界上所应满足的条件,也就是渗流区内水流与其周围环境相互制约的关系。
(4)初始条件(initial conditions) :表示渗流区的初始状态,某一选定的初始时刻(t=0)渗流区内水头H 的分布情况。
将边界条件和初始条件并称为定解条件(definite solution condition) ,微分方程和定解条件一起构成渗流场的数学模型。
数学模型:描述某一研究区地下水流运动的数学方程与其定解条件共同构成的表示某一实际问题的数学结构。亦即从物理模型出发,用简洁的数学语言,即一组数学关系式来刻画它的数量关系和空间形式,从而反映所研究地质体的地质、水文地质条件和地下水运动的基本特征,达到复制或再现一个实际水流系统基本状态的目的的一种数学结构。
其中微分方程表示地下水的流动规律,定解条件表明研究对象所处的特定环境条件,即所研究的地下水流的真实状态。
定解问题是给定了方程(或方程组)和相应定解条件的数学物理问题。 建立模型是指建立数学模型的过程。 1.7.2 定解条件 1.7.2.1定解条件
定解条件指水头、流量等渗流运动要素在流场边界上的已知变化规律,这种变化规律是由流场外部条件引起的,但它不断地影响流场内部的渗流过程,并在整个期间一直起作用。
定解条件包括边界条件和初始条件。 1.7.2.2边界条件
边界条件是渗流区边界所处的条件,用以表示水头H (或渗流量q )在渗流区边界上所应满足的条件,也就是渗流区内水流与其周围环境相互制约的关系。
(1) 第一类边界条件(Dirichlet条件) :如果在某一部分边界(设为S l 或Γ1) 上,各点在每一时刻的水头都是已知的,则这部分边界就称为第一类边界或给定水头的边界,表示为:
H (x , y , z , t ) |s 1=ϕ(x , y , z , t ),
或
(x , y , z ) ∈s 1
(1-103)
H (x , y , z , t ) |Γ1=ϕ(x , y , t ),
(x , y ) ∈Γ1
(1-103)
给定水头边界不一定就是定水头边界。
可以作为第一类边界条件来处理的情况:
① 河流或湖泊切割含水层,两者有直接水力联系时,这部分边界就可以作为第一类边界处理。此时,水头φ是一个由河湖水位的统计资料得到的关于t 的函数。但要注意,某些
河、湖底部及两侧沉积有一些粉砂、亚粘土和粘土,使地下水和地表水的直接水力联系受阻,就不能作为第一类边界条件来处理。
在自然界,这种情况很少见。就是附近有河流、湖泊,也不一定能处理为定水头边界,还要视河流、湖泊与地下水水力联系的情况,以及这些地表水体本身的径流特征而定。在没有充分依据的情况下,不要随意把某段边界确定为定水头边界,以免造成很大误差。
② 区域内部的抽水井、注水井或疏干巷道也可以作为给定水头的内边界来处理。此时,水头通常是按某种要求事先给定,例如给定抽水井的允许降深等。上面介绍的都只是给定水头的边界。注意,给定水头边界不一定是定水头边界。
③ 排泄地下水的溢出带、冲沟或排水渠的边界也可近似看作给定水头边界。 (2)第二类边界条件(Neumam条件) :当知道某一部分边界(设为S 2或Γ2) 单位面积(二维空间为单位宽度) 上流入(流出时用负值) 的流量q 时,称为第二类边界或给定流量的边界。相应的边界条件表示为:
K ∂H
∂n
∂H ∂n
S 2
=q 1(x , y , z , t ), (x , y , z , t ) ∈S 2
(1-105)
或
T
Γ2
=q 2(x , y , t ), (x , y ) ∈Γ2
(1-106)
式中,n 为边界S 2或Γ2的外法线方向。q 1和q 2则为已知函数,分别表示S 2上单位面积和Γ2上单位宽度的侧向补给量。
常见的这类边界条件: ① 隔水边界(流线、分水岭):
∂H
∂n
Γ
=0
(1-107)
② 抽水井或注水井:
T
∂H ∂n
Γw
=-
Q
2πr w (1-108)
③ 补给或排泄地下水的河渠边界上,如已知补给量。
∂H
(3)第三类边界条件:某边界上H 和∂n 的线性组合是已知的,即有:
∂H
+αH =β∂n (1-109)
又称混合边界条件,α, β为已知函数。
σ' =
边界为弱透水层(渗透系数为K 1,厚度或宽度为m 1),
m 1
K 1,有
K
∂H ∂n ∂H ∂n
=
K 1
(H n -H ) =q (x , y , z , t ) m 1
H n -H
=0σ' (1-110) H n -H
=0σ' (1-111)
在s 3上,
K -
S 3
∂H
在Γ3上,∂n
T
-M
Γ3
浸润曲线的边界条件:
K
∂H ∂n
=q
c 2
(1-112)
当浸润曲线下降时,从浸润曲线边界流入渗流区的单位面积流量q 为:
∂H *q =μcos θ
∂t (1-113) 式中,μ为给水度,θ为浸润曲线外法线与铅垂线间的夹角。
1.7.2.3初始条件
初始条件:某一选定的初始时刻(t=0)渗流区内水头H 的分布情况。
H (x , y , z , t ) |t =0=H 0(x , y , z ),
或
(x , y , z ) ∈D (1-114)
H (x , y , t ) |t =0=H 0(x , y ),
(x , y , z ) ∈D (1-115)
其中,H 0为D 上的已知函数。 1.7.3 渗流数学模型的分类
(1)线性、非线性模型
模型由线性方程所组成,称为线性模型,如均质各向同性承压二维流方程。模型由非线性方程所组成,称为非线性模型,如潜水模型方程。
(2)静态、动态模型
根据模型中未知变量与时间的关系进行划分,若未知变量与时间无关,如稳定流模型,称为静态模型,反之,则为动态模型。
(3)集中、分布参数模型
模型中不含有空间坐标变量的模型,称为集中参数模型,如抽水井流量与降深之间的经验公式。模型中含有空间坐标变量的模型,称为分布参数模型。
(4)确定性与随机性模型
确定性模型:数学模型中各变量之间有严格的数学关系的模型。 随机性模型:数学关系式中含有一个或多个随机变量的模型。 用确定性模型来描述实际地下水流时,如前述,必须具备下列条件:
① 有一个(或一组) 能描述这类地下水运动规律的偏微分方程;同时,确定了相应渗流区的范围、形状和方程中出现的各种参数值。
② 给出相应的定解条件。对所建立的模型进行检验,即把模型预测的结果与通过抽水试验或其它试验对含水层施加某种影响后所得到的实际观测结果或一个地区地下水动态长期观测资料进行比较,看两者是否一致。若不一致,就要对模型进行校正,即修正条件(1)和(2)直至满意拟合为止。这一步骤称为识别模型或校正模型。
经过校正后的模型,能代表所研究的地质体,或者说是实际水流系统的复制品了,因而可以根据需要,用这个模型进行计算或预测,例如预测矿床疏干时的涌水量及地下水污染情况预测等。
③ 解(即满足条件①和②的解) 是存在的(存在性) ;
④ 解是唯一的(唯一性) ;要求所提问题的解存在和唯一是不言而喻的。
⑤ 解对原始数据是连续依赖的(稳定性) 。即稳定性的要求,意味着当参数或定解条件发生微小变化时,所引起的解的变化也是很微小的。只有有了这条保证,当参数和定解条件的数据有某些误差时,所求得的解才能仍然接近于真解;否则,解是不可信的,并应该认为此时的数学模型是有毛病的。在实际工作中,原始数据有某种误差,在所难免,所以这个条件很重要。
适定问题(Well –posed problem )是指数学模型满足(1)解是存在的(存在性),(2)解是唯一的(唯一性),(3)解对原始数据是连续依赖的(稳定性)这三个条件的问题。只要有一条不满足就是不适定问题。
正问题是根据数学模型、给定的含水层水文地质参数和定解条件求解水头的问题,又称水头预报问题。
逆问题(inverse problem)是根据数学模型、动态观测资料或抽水试验资料反过来确定含水层水文地质参数的问题。 1.7.4 建立数学模型的基本要点
(1)确定研究区的范围及渗流区的边界;
(2)确定渗流区的水力特征(包括埋藏条件、渗流状态、介质特征); (3)确定渗流区的边界条件; (4)确定渗流区的源汇项; (5)选择微分方程; (6)确定渗流区的初始条件。 1.7.5渗流数学模型的解法
(1)解析法(analytic method)
用参数分析及积分变换等方法直接求解数学模型解的方法。其解为精确解,使用简单,但该方法存在一定的局限性,只适用于含水层几何形状规则、方程式简单、边界条件单一的情况。
解析解(Analytic solution)又称精确解,是用解析方法求解数学问题所得到的解析表达式。
(2)数值法(numerical method)
用数值方法(离散化方法)求解数学模型的方法,其解为近似解。该方法是求解大型地下水流问题的主要方法。它把整个渗流区分割成若干个形状规则的小单元,每个小单元近似处理成均质的,然后建立每个单元地下水流动的关系式。把形状不规则的、非均质的问题转化为形状规则的均质问题。根据研究需要,确定单元划分数量,对于非稳定流还要对时段进行划分。最后,把局部整合起来,加上定解条件。
数值解(Numerical solution)是 用数值方法求得的数值解,是一种近似解。 (3)模拟法:利用物理现象与水流的相似性,在实验室内采用模拟的方法求解。