气沉降。在此，我们将模型理想化，假设只有这三种方式，分别对其建立模型，之后再将三种模型进行时间和空间上的几何叠加。
1.水体传播模型的建立及求解：
在此模型中，我们将问题理想化，认为水体传播仅仅是靠河流传播的。假设河流没有多大的改变。
在河流中，重金属污染物部分以可溶状态存在，随水流运动，部分以沙粒为载体，随沙运动，两者不断地发生化学和生物化学反应，一般情况下，泥沙吸附—解吸是重金属污染物从不饱和液相转移到固相的主要途径。因此，在考虑重金属在水环境中的物理迁移过程下，我们用水力学和泥沙运动学的原理处理。
假设，泥沙浓度和水相金属离子浓度反应分布均匀。在天然河流中, 水环境化学条件对某一河段、某一时段基本不变,变化较大并难处理的是水流泥沙因素。
泥沙颗粒在水溶液中对重金属污染物的吸附遵循兰格茂尔吸附运动学方程
dN?K1c(b?N)?K2N，将其对时间的偏微分 ，即对均匀沙有： dt?N?k1c(b?N)?k2N （11） ?t水流泥沙运动数学模拟研究已有较长的历史，在水流泥沙运动确定后, 根据建立的重金属迁移转化方程, 吸附动力学方程及边界、初始条件就可确定河流重金属污染物迁移转化过程, 构成河流一维重金属污染物迁移转化数学模型的基本方程如下：
??y水流连续方程 (Bhu)?B?0 （12）
?x?t?u?u?yu2水流运动方程 ?u?g?g2?0 （13）
?t?x?tc1k??y(Bhu)?B(hs)???B?(s?s0) （14） ?x?t?y??y河床变形方程 (Bhu)?B(hs)???B0?0 （15）
?x?t?t泥沙连续方程
有A?Bh,A3?By0及L=2h?B，则
?A1?N?c?c1??c1?s?N1???N?u?(AE1')???N13??N3L???s1?us1?E1'?（16） ?t?xA?x?xA?tA?x?x?x???t 假定悬浮泥沙的吸附量N1?N1?N1?N1'其中，N=N1'=0，对（11）式全断面积分得：
?N1?k1sc(bs?N1)?k2sN1 （17） ?t 解决重金属迁移转换模型可将（12）至（17）式简化成：
2?y0?N?Nv3m?c?c'?ch=常数，u?常数，?0,s??k(),?u?E12??s1?us1,?tgkw?t?x?x?t?x ?N1?k1c(b?N1)?k2N1?x 11
根据静态吸附的结果，取b=0.543g/kg,k1=0.0076 1/(ppm·s)，
k2=0.00084 1/s ,重金属污染物的初始条件为：
cx?0?c0?(x)式中：c0?1ppm?1?0x?0,Ntx?0?(x)=?t?0?0
边界条件为：
cx?0,?1?c0?(t),?(t)=??0t?0?2c,t?0?x2?0,N1x?Lx?0?0
在相同条件下，污染物的运动遵循以下方程：
2?c?c'?c?u?E12?0,ct?0?c0?(x),cx?0?c0?(t),?2c?x2?t?x?xx?L?0。用向后显式差分求解，结
果如下图2:
由图2可知：在不同时刻, 水相重金属浓度由于泥沙吸附, 峰值明显减少, 而且由于泥沙吸附需要时间, 因此河段上泥沙吸附量的峰值总是滞后于水相浓度的峰值。
从3中可以看出：污染物不受泥沙运动影响, 在输移过程中, 由于纵向综合紊动扩散作用, 峰型变化相同条件下重金属污染物输移, 尽管图形相似, 但重金属污染物的峰值因泥沙吸附而下降很多, 特别是后半部减小更为明显, 且随时间增加明显加大。
所以，重金属在水体中的传播特征是随着离污染源的距离增加而减少。速度也随之减小。
2.土壤传播模型的建立及求解：
重金属在土壤中的传播的影响因素也很多，不能单纯的看作是线性变化，在此，我们基于重金属在土壤中的质量分数，根据在水中传播的方式类似的定义了一个非线性沉淀模型。并推导出非线性沉淀情况下重金属速度的控制方程。
对于饱和土壤空隙内充填的溶液质量分数为ｗ，而土粒表面会吸附溶质，因此，
12
当溶质在土／水界面达到平衡之后，会存在以下关系：
Ft(j)?F?(j)?e?(j)/?0 （18）
式中：Ft(j) — 土壤中溶质的总量在单位质量土粒上的质量分数，包括土粒
上吸附溶质和孔隙溶液中的溶质；
F?(j) — 土颗粒表面所吸附的溶质（j）的质量分数；
?(j) — 孔隙中的溶质的质量浓度（ｇ／Ｌ）； ｅ — 土壤孔隙比；
?0 — 土样干密度（g/cm3）
当溶质在土粒上的吸附呈线性状况时，沉淀满足Henry定律：
Ft(j)?vd?(j)． （19）
式中：vd — 溶质在土水界面的分配系数．此时，满足下式：
Ft(j)?(A?B)?(j) （20）
式中：Ａ、Ｂ 均为模型参数，可表示为A?(e?b??s)/s?b,B???b
? — 土壤空隙率；
S — 土壤容重． ．
对于溶质在土壤介质上的吸附，其分配系数通常情况下并非常数，而是随溶质浓度变化而变化的参数．为了便于计算，本文对溶质在土柱中土水间的分配关系做出如下简化假定：
Ft(j)??(?(j))n， （21）
传统的运移控制方程如下：
??(j)???(j)??(j) （22） Rd?(D(?(j)))?v?t?x?x?x
式中：ｔ为时间，Rd为阻滞因子，D(?(j)）表示随土壤孔隙中溶质的质量浓度 ?(j)而变化的非线性的扩散系数，v为孔隙流体的流速．
式（22）是以土壤孔隙中溶质的平衡浓度为基础建立起来的．假定扩散系数Ｄ 为常数，考虑吸附的非线性特性，将式（21）代入式（22），可得改进的非线性重金属运移控制方程为
1?n?F(j)1?n?F(j)?F(j)?(?bn(?)1n/e)t?(DFt(j)nt)?vFt(j)nt （23）
?t?x?x?x
13
?Ft0(j)??(?0(j))n其中：???Ft(j)?0??xFt0(j)?(上边界)(下边界)
基于土粒上溶质质量分数的运移参数提取方法，根据有限差分法，可将运移方程（22）
写成如下差分格式：
RddFi(j)F(j)?Fi(j)F?Fi?1vFi?1(j)?Fi(j)2D?(Fi?1/2(j)i?1?Fi?1/2(j)i)?(dtxi?1?xi?1xi?1?xixi?xi?12xi?1?xiF(j)?Fi?1(j)?i)xi?xi?1 （23）
式中：xi — Ｘ 方向被差分之后的差分网格节点的空间位置；
Fi(j)≈F(j)，ｉ＝１，…，Ｎ－１；F0(j)、FN(j)为上下边界处溶质的质量分
数；Fi0(j)为0时刻土柱内溶质的质量分数．
上述差分格式可以利用Matlab软件进行求解，结合最小二乘法，利用内嵌的最优化工具箱求解运移参数．对比拟合结果与试验得到的浓度剖面数据，当方差最小时即认为此时的参数即为有代表性的典型参数．
由模型可得重金属在土壤中的传播模型如下图所示：
在图中可以清晰地看出重金属在土壤中的传播随时间的增减而减缓，而且速度较水体传播慢。
3. 大气扩散模型的建立及求解：
大气中重金属主要来源于工业生产、汽车尾气排放、汽车轮胎磨损产生的大量含有重金属的有害气体和粉尘。它们主要分布在工矿的周围及公路、铁路两侧。
重金属在大气中的传播一部分通过自然沉降和雨淋沉降进入土壤圈和水体中，一部
14
分，直接污染空间。
在这里因为数据不足更为了简化模型，我们忽略雨淋沉降，将通过大气沉降地表的重金属看做由自然沉降产生的。而自然沉降的主要作用是重力作用。
在垂直方向，考虑粒子受到重力、浮力、流体阻力和随机力的作用，建立粒子的运动方程如下：
dwG?F浮?F阻???(t) （24） dtm其中w为粒子垂直方向速度，F浮为粒子受到的浮力，F阻为粒子受到的流体阻力，
?(t)为随机力及大气流动引起的随机加速度，满足如下条件:
?(t)?0,?(t)?(t')?K?(t?t?) (25)
粒子受到的阻力为F阻?6????，其中?为动力粘性系数，r为粒子半径。 由式（24）得：
?gdw9??g(1?)?w??(t) 其中?g和?p分别为气体和粒2dt?p2?pr2?pr29?，由上式得到：
?t子密度，通常?g《?p。令??dw1 ?g?w??(t) 解此式可得：
dt?ts?tw(t)?g?(1?e?)?e?t?w(0)??e0??(s)ds （26）
将 w(t)分解成为平均量和脉动量之和： w(t)?w(t)?w?(t)
由这两式可得：
w(t)?g?(1?e?t?)?e??t?w(0) （27）
ts?tw?(t)?e?tw(0)??e0??(s)ds （28）
由上式求出w(t)的方差w?2(t)，由于w(t)和?(t)不相关，则，
w?2(t)?w?2(0)e?2t?)?2K?(1?e?t?2t?) （29）
由（28）式可得： w?(0)w?(t)?w?2(0)e?
2假设w(t)是平稳随机过程，则：w?2(0)?w?2(t)??w 2?w由上面三式可得：2K??? 即：K?
2?2w由（26）式可得：w(t??t)?g?(1?e?(t??t)?)?e?(t??t)?w(0)??(t??t)0es?t??t??(s)ds
15