2.数值法 数值法
数值法是把刻划地下水运动的数学模型离散化，把定界问 题化成代数方程，解出渗流区域内有限个结点上的数值解 题化成代数方程，解出渗流区域内有限个结点上的数值解 （近似解）。 （近似解）。 数值法适用性广（复杂的含水层、定解条件等），通用性强 （计算机模拟）、并可程序化，修改模型方便。目前，在求解 大型地下水流问题时被广泛应用。 数值法模拟计算过程没有物理模拟法逼真、直观，而且计 算工作量大，需要借助于计算机进行模拟计算。 数值法中，最常用的是有限差分法（ 数值法中，最常用的是有限差分法（FDM）和有限单元法 ） （FEM)。有限差分法是建立在用差商代替导数的基础上；而有 。有限差分法是建立在用差商代替导数的基础上；而有 限单元法是建立在直接求函数的近似解的基础上。 限单元法是建立在直接求函数的近似解的基础上。
k i
T
∆t
+ （∆t ) 0
（8-5） ）
略去0(△ 和 △ 略去 △t )和 0 (△x)2 ，可得（8-5）式的对应的差分 可得（ ） 方程： 方程：
hk i −1 − 2hk i + hk i +1 (∆x)
2
=
µ ∗ hi k +1 − hi k
T ∆t
（8-6） ）
上式也可变为：
( hik−1 − 2 hik + hik+ 1 ) = hik + 1 − hik µ * (∆x)2 T∆t
若定义
λ=
T∆t µ (∆x)
∗
2
，则（8 ，则（8-6）式可变为：
k k k i+1 i−1 i
h
k +1
i+1
= λh + 1− 2λ）h + λh （
（8-7） ）
（8-7）式表明：只要知道了 时段初始时刻tk各 ）式表明：只要知道了k时段初始时刻 结点的h 值，便可计算出k时段末了时刻 结点的 ik值，便可计算出 时段末了时刻tk+1的hik+1 值（l≤i ≤l-1, 1≤k ≤m-1)，各方程可独立求解，因此， ， 这种方程称为显式有限差分方程。 这种方程称为显式有限差分方程。
0
图8-1 河间地块承压含水层示意图
在上述条件下，地下水有由高水位一边向低水位一 边流动，构成承压水一维非稳定流问题。取如上图所示 边流动，构成承压水一维非稳定流问题。取如上图所示 的坐标系，则该问题的数学模型为：
∂ 2H
2
Ⅰ
T ∂t ∂x H ( X ,0) = H 0 ( x) H (0，t） φ0 (t ) = H (l , t ) = φl (t )
k +1 l
k +1 l
= φ (t
l
k +1
)
2. 显式差分方程的求解 差分方程（ ）与离散后的定界条件（ ）、 差分方程（8-7）与离散后的定界条件（8-8）、 ）、(8-9） ） 构成了数学模型Ⅰ 构成了数学模型Ⅰ的显式差分方程问题。其求解步骤如下：
由（8-8）式所示的初始条件给出 0时刻各结点的水头值 ）式所示的初始条件给出t h00， h10，…hl0；再根据（8-7）式，在k=0时，分别取 再根据（8-7） k=0时 i=1，i=2,…i=l-1,便可求得 1时刻各内结点的水头值 便可求得t ， h11，…hl-10。 由以上计算的h 由以上计算的 11， h21 …hl-11值及由边界条件（8-9）式计 ） 算的h01和hl1,再次利用（8-7）式（取k=1， i=1， 再次利用（ ， ， i=2,…i=l-1，便可计算得 2时刻各结点的水头值。如此重复， ，便可计算得t 便可计算出t 便可计算出 3 ，t3 ，…各时刻的水头分布值。 各时刻的水头分布值。
hk +1i −1 − 2hk +1i + hk +1i +1 (∆x)
λ 若定义： =
∗
2
=
µ∗ hi k +1 − hi k
T ∆t
（8-9）
T ∆t µ ( ∆x )
2
，则（8-9）式可变为：
k +1 i
− λh
k +1 i −1
+ 1 + 2λ） h （
− λh
k +1 i +1
=h
k i
i=1,2,…,l-1 （8-10） K=1,2,… 有限差分法
有限差分法的基本思想
把渗流区域按一定的方式剖分成许多小区域 均衡域），用该区域中心点（结点） ），用该区域中心点 （均衡域），用该区域中心点（结点）的集 集合代替连续的渗流区域， 集合代替连续的渗流区域，在这些点上用差 商近似地代替导数， 商近似地代替导数，将描述地下水流问题的 数学模型化为一组以有限个未知函数值为未 知量的差分方程（代数方程） 知量的差分方程（代数方程）组，通过求 解差分方程组， 解差分方程组，得到所求解在结点上的近似 值。
(3) 定解条件离散化 在t=0的初始时刻，各结点水头值由初始条件给出： t=0的初始时刻，各结点水头值由初始条件给出：
h = H
0 i
0 i
= H (i∆ x )
0
（8-8） ）
边界结点上各时刻的水头值由边界条件给出：
h h
k +1 0
= H = H
k +1 0
= φ (t
0
k +1
)
（8-9） ）
(2) 微分方程的差分化
在t=tk时刻，（8-1）式左端用二阶导数的近似表达式， 时刻， ） 右端用时间导数的向前差分近似代替导数，有： 右端用时间导数的向前差分近似代替导数，有：
H − 2H + H (∆x)
k k i −1 i 2 k i +1
+ （[∆x] ) = 0
2
µ H −H
∗ k +1 i
差分的概念
用差商近似代替导数
一阶导数的近似表达式：
dH H ( x + ∆x ) − H ( x ) = lim dx ∆x →0 ∆x H ( x + ∆x) − H ( x) ∆x H ( x) − H ( x − ∆x) dH ≈ dx ∆x H ( x + ∆x) − H ( x − ∆x) 2 ∆x
=
µ ∗ ∂H
0≤x≤L, 0≤t≤Tsum (8-1) 0≤x≤L 0≤t≤Tsum 0≤t≤Tsum (8-2) (8-3) (8-4)
由于导数有不同的差分格式，所以差分方程也有不 同的形式。
二、一维显式有限差分格式
1. 显式差分方程的建立
（1） 离散化 ）
m
空间离散：首先将研究区 空间离散：首先将研究区 域[0，L]用直线等分为l [0，L]用直线等分为l 份，空间步长△ 份，空间步长△x=L/l;
四、一维中心有限差分格式
1. 中心差分方程的建立
这种差分格式是在t 这种差分格式是在 k与tk+1时刻之间，取中间过渡时 刻tk+△t/2时刻，以结点（i,k+1/2)处的微分方程为基础 △ 时刻，以结点（ ) 建立差分方程。 在t=tk+1/2时刻，（8-1）式右端用时间导数的中心差 时刻，（8 分近似代替导数，则有：
向前差分 向后差分 中心差分
二阶导数的近似表达式：
d H H ( x + ∆ x ) − 2 H ( x ) + H ( x − ∆x ) ≈ dx ( ∆x )
2 2 2
一、承压水一维非稳定流数学模型的概化
设在两条平行的河流之间 有一均质、各向同性、等 厚、无越流补给的承压含 水层； 两河流间距离为L，两河 ， 水位的变化规律为Φ0(t) 和Φl(t) ； 初始时刻含水层内的水头 分布为：H0(x)(0≤x≤L)。 。
∂H （ ) ∂t
k +1 / 2
h −h = ∆t
k +1 i
k i
（8-11） ）
3. 显式差分方程的收敛性和稳定性
差分方程的解h 是否逼近原微分方程的解H 差分方程的解 ik+1是否逼近原微分方程的解 ik+1?必须 必须 从差分方程的收敛性和稳定性两个方面回答此问题。 从差分方程的收敛性和稳定性两个方面回答此问题。 收敛性：如果在△ 收敛性：如果在△x, △t取得充分小时，差分方程的解和微 取得充分小时，差分方程的解和微 分方程的解析解很接近，便说这种差分格式是收敛的。 稳定性： 稳定性：差分计算时，每一步都有舍入误差，随着计算 时间或计算次数的增加，累积误差逐渐减小，以至于不 影响计算结果，那么这种差分格式便是稳定的。
显式差分格式收 敛和稳定的条件
1 0p ≤ λ 2
λ =
T∆t µ (∆ x)
∗
2
只有收敛和稳定的差分格式，才具有实用价值。因此， 只有收敛和稳定的差分格式，才具有实用价值。因此， 合理选取△ 合理选取△x, △t是很重要的。 是很重要的。
三、一维隐式有限差分格式
1. 隐式差分方程的建立 在t=tk+1时刻，（8-1）式左端用二阶导数的近似表达 时刻， ） 式，右端用时间导数的向后差分近似代替导数，略去0(△ 式，右端用时间导数的向后差分近似代替导数，略去 △t ) 和 0(△x)2 ，则有： △
有限差分法和有限单元法优缺点比较
数值法 优 点 缺 点
有限差 分 法
1.简单问题的数学表达式和计算的执 1.简单问题的数学表达式和计算的执 行过程比较直观、易懂； 行过程比较直观、易懂； 对自然边界处 2.算法效率比较高 算法效率比较高； 2.算法效率比较高； 理的灵活性较 3.计算精度高 计算精度高； 3.计算精度高； 差。 4.有可使用的商用软件 有可使用的商用软件。 4.有可使用的商用软件。 1.计算程序的通用性强； 1.计算程序的通用性强； 计算程序的通用性强 2.对不规则边界处理方便 对不规则边界处理方便； 2.对不规则边界处理方便； 占用计算内 3.计算单元划分灵活 计算单元划分灵活； 3.计算单元划分灵活； 存大， 存大，计工 4.水流问题 水流问题、 4.水流问题、物质输运问题解的精度 作量大。 作量大。 一般较FDM的精度高。 FDM的精度高 一般较FDM的精度高。