重金属污染国赛题

我们仔细阅读了中国大学生数学建模竞赛的竞赛规则.

我们完全明白，在竞赛开始后参赛队员不能以任何方式（包括电话、电子邮件、网上咨询等）与队外的任何人（包括指导教师）研究、讨论与赛题有关的问题。

我们知道，抄袭别人的成果是违反竞赛规则的, 如果引用别人的成果或其他公开的资料（包括网上查到的资料），必须按照规定的参考文献的表述方式在正文引用处和参考文献中明确列出。

我们郑重承诺，严格遵守竞赛规则，以保证竞赛的公正、公平性。如有违反竞赛规则的行为，我们将受到严肃处理。

我们参赛选择的题号是（从A/B/C/D中选择一项填写）：

我们的参赛报名号为（如果赛区设置报名号的话）：

所属学校（请填写完整的全名）：

参赛队员(打印并签名) ：1.

2.

3.

指导教师或指导教师组负责人(打印并签名)：

日期：年月日赛区评阅编号（由赛区组委会评阅前进行编号）：

2010高教社杯全国大学生数学建模竞赛

编号专用页

赛区评阅编号（由赛区组委会评阅前进行编号）：

全国统一编号（由赛区组委会送交全国前编号）：全国评阅编号（由全国组委会评阅前进行编号）：

城市表层土壤重金属污染分析

摘要

本文运用Matlab、SPSS及Excel软件，综合分析该城区各区重金属污染程度及其重金属污染的主要原因，并根据重金属的传播特征建立数学模型，寻找出污染源，同时更进一步修改模型以更好研究城市地质环境的演变模式。

针对问题一，利用Matlab作出各重金属浓度空间分布图，可看出重金属浓度大多呈高浓度向外逐步扩散现象。同时，利用单因素污染指数公式及污染等级规则，通过模糊综合评价模型，得出各类区域重金属污染程度如下：

针对问题二，利用SPSS软件对各重金属变量进行R型聚类分析，将重金属污染分为四类：Cr、Ni、Cu为一类，Cd、Pb、Zn为一类，而As、Hg各为一类。并计算出各区域各单金属污染指数进行对照分析，说明该城区重金属污染主要来源于工业活动及其生产品和交通主干道污染。

针对问题三，我们根据重金属物质在土壤中扩散的特点，利用有衰减的扩散模型，构造出浓度与坐标的关系，并利用matlab编程，多次试验，合理筛选数据，通过回归分析的方法，求解参数，从而得到不同重金属的不同污染源。

针对问题四，由于问题三求解的局限性，我们在模型三的基础上，增加了时间、水文、人类活动等因素，得出一个更为具体并符合实际的模型。有利于计算土壤中重金属的浓度。

关键词单因素污染指数模糊综合评价 R型聚类分析有衰减的扩散模型扩散规律

一、问题重述

随着城市经济快速发展和城市人口不断增加，人类活动对城市环境质量的影响日显突出。对城市土壤地质环境异常的查证，以及展开城市环境质量评价，研究人类活动影响下城市地质环境的演变模式，日益成为人们关注的焦点。

按照功能划分，城区一般可分为生活区、工业区、山区、主干道路区及公园绿地区等，分别记为1类区、2类区、……、5类区，不同的区域环境受人类活动影响的程度不同。

现对某城市城区土壤地质环境进行调查，将城区划分为间距为1公里左右的网格子区域，按照每平方公里1个采样点对表层土（0~10 厘米深度）进行取样、编号，并用GPS记录采样点的位置。应用专门仪器测试分析，获得了每个样本所含的多种化学元素的浓度数据。另一方面，按照2公里的间距在那些远离人群及工业活动的自然区取样，将其作为该城区表层土壤中元素的背景值。

现根据采样点的位置、海拔高度及其所属功能区，以及其8种主要重金属元素在采样点处的浓度和背景值等信息，建立数学模型，解决以下问题：

(1) 给出8种主要重金属元素在该城区的空间分布，并分析该城区内不同区域重金属的污染程度。

(2) 通过数据分析，说明重金属污染的主要原因。

(3) 分析重金属污染物的传播特征，由此建立模型，确定污染源的位置。

(4) 分析所建立模型的优缺点，并收集更多信息，建立数学模型，以更好地研究城市地质环境的演变模式。

二、问题分析

对于问题一，可通过Matlab画图函数分别直观描述出该城区各类重金属的空间分布图，包括各采样点的平面位置、海拔以及其浓度分布，并大致观测各类区重金属浓度的变化。同时，对该城区内不同区域重金属的污染程度进行分析，首先需引入单因素污染指数及其计算公式，并给出重金属元素的污染等级评判标准。基于此，本文建立因素集和评价集，以简单模糊统计方法构造隶属函数，并在各重金属元素本身对土壤地质的作用是一样的假设下，即不考虑各重金属本身的毒性，通过模糊综合评价模型对该城区各区域的各类重金属总污染程度进行综合评定。

对于问题二，需对数据进行分析，来说明重金属污染的主要原因。考虑到重金属种类有8种，比较多，可通过R型聚类分析对相关性较高的重金属元素聚分成几类，基于此法，可大致判断出哪几类重金属元素很可能源于同一性质或同一类污染物或者相似人为性活动，同时计算出各区域各重金属的污染指数和污染等级，由此进行对照分析，综合推测出该城区重金属污染的主要原因。

对于问题三，考虑到土壤中重金属的传播特征与分子运动扩散性质很相似，即一般以污染源为起点向周围辐射状传播，并且浓度逐渐地衰减，传播方向是从高浓度到低浓度。针对这些传播特点，建立衰减扩散偏微分方程。首先根据问题一所做的图进行初步判断，观察高浓度的大约范围，由此选择数据。然后利用有衰减的扩散数学模型，和回归分析，由此进行参数估计，得到污染源的坐标。

对于问题四，我们考虑到问题三所建立的有衰减的扩散数学模型具有一定的局限性，

因此在问题三的基础上，研究考虑更多能够影响土壤重金属浓度的因素，并将其加入问题三的问题中，由此得到一个更全面的新模型。

三、模型假设

（1）问题一中，各重金属元素本身对土壤地质的作用是一样的，即综合模糊评价模型中的各重金属元素所占的权重一样。

（2）问题三中，土壤中重金属的传播服从扩散定律，并不考虑其它大气或水流等因素的影响。

（3）问题三中，假设重金属在水平方向上的扩散速度是一样的，土壤表面顺滑，不存在重金属物在非污染源处堆积现象。

四、符号说明

五、模型的建立与求解

5.1 问题一的求解

5.1.1 各金属元素空间分布

利用matlab 软件，以x 表示横坐标，y 表示纵坐标，z 表示平面上点(x,y)对应的海拔高度，利用色卡，其蓝色至红色渐变表示重金属浓度由低到高变化，得到如下图（其余见附录）；

图5.1.1-1：As 空间分布图

由以各分布图可直观看出，各类金属在空间分布上大多呈现出局部集中的现象，浓度较高，再从四周逐步扩散，浓度逐渐降低。

5.1.2 模糊综合评价模型的建立

针对问题一，可引入单因素污染指数，作为评价各类重金属污染程度的标准。因此，若评价一个区域的污染程度，可进行单因素污染指数的评价，再次对该区域的各类重金属进行综合评价。单因素污染指数的计算公式如下：

i p =

i

i

s c ，i =1,2,3…,8 其中，i p 是第i 类重金属污染指数，i c 是第i 类重金属的实际浓度，i s 是第i 类金属的背景浓度值。

结合参考文献]1[，得到各类重金属污染的等级判定规则，如表5.1.2-1：

该城区各区域的污染程度分析主要依据其各类重金属元素浓度高低，因此可用模糊数学模型中的模糊综合评价法对重金属浓度进行综合评价。其中，重金属的污染指数构成因素集U 。并记为：

{=U 0u ，1u ，2u ，3u ，4u ，5u ，6u ，7u }

评价集中含有四个元素，分别为非污染，轻度污染，中度污染和重度污染，并记为：

{=V 0v ，1v ，2v ，3v }

即各因素集和评价集元素表示如表5.1.2-2：

5.1.3模糊综合评价模型的求解

由于题设要求对五个区域进行污染程度评价，因此以1类区（生活区）为例，进行模糊评价。

（1） 计算各重金属污染指数，作第一次定量。数据处理如表5.1.2-3：

（2） 用模糊统计方法构造隶属函数，同时构造模糊评判矩阵A 。隶属函数即为：

ij a =

n

u ij ，j i , =1,2,3,4,5,6,7,8

得：

?????????????

??

???????????????=4413223115445223221911622302232219447411121134413229112

11422

111

144544311114474413441744744311344

24111A （3） 确定权向量e ，作第二次量化并求解D

设每种金属对土地污染程度的作用是相同的，并保证归一化，取权向量为：

E={1/8,1/8,1/8,1/8,1/8,1/8,1/8,1/8}

令D=EA,得：

D ={0.14，0.52，0.27,0.18}

（4） 取值，得出评价结果。

取α=0.8，1β=0.8，2β=0.9，3β=0.75. 则1d =0.14<α,

1d +2d ×1β=0.556<α, 1d +2d ×1β+3d ×2β=0.82>α,

因此，1类区（生活区）属于中度污染区。

以同样的方法，可求出其余四个区域的重金属污染程度，结果如表5.1.2-4：

5.1.4 结果分析

由表5.1.2-4可知，2类区（工业区）与4类区（主干道路区）的重金属污染程度属

重度污染，该城区应加强对其重金属污染进行重点控制；其次，1类区（生活区）与5类区（公园绿地区）属中度污染，可见人们生活环境及城市绿化水平正趋于严重化；而3类区（山区）属轻度污染，需做好防范措施防止其进一步恶化。

5.2 问题二的求解

5.2.1 R 型聚类模型的建立与求解 R 型聚类分析模型主要通过研究指标变量间的相似关系将其聚合成若干类，进而找出影响系统的主要因素。基于此模型，针对该城区，分析各重金属间的相关系数，从而推断出重金属污染的具体主要因素。

（1）对数据进行标准化处理

记重金属样本矩阵为n m ij x X ?=)(，-

j x ,j S 分别为第j 类重金属的均值和标准差，有：

j

j

ij ij S x x x -

-=

~

（2）确定相关性度量

记各重金属浓度j x T nj j j x x x ),,,(21 ，j =1,2,…,m 。则以两重金属浓度j x 与k x 的样本相关系数作为其相似性度量为：

2

/11

1

2

1])()([)

)((∑∑∑==-

-

=-

-----=

n i n

i k ik

j iij n

i k ik j ij

jk x x

x x x x x x

r ，

因此，可得到各重金属两两相关系数矩阵A 。 （3）对各重金属浓度进行相关性显著检验 在进行聚类分析之前，运用SPSS 软件，对各重金属浓度进行Pearson 相关性显著检验，得到表5.2-1：

由上表可知，在显著水平 =0.01上，各重金属间浓度显著相关，即可利用Pearson 相关性度量标准来进行聚类分析。初步判断，Cr和Ni、Pb和Cd的相关系数分别为0.717和0.665，相关性比较大，很有可能来源于同一污染物或同一污染性活动。

（4）进行聚类分析

针对上述，运用SPSS软件的分析工具，进行变量分群，使用组间连接聚类方法和Pearson相关性度量标准，选择Z-Scores标准化，可得到图5.2-2：

图5.2-2

由树状图可看出，这8种重金属可大体分为四类， Cr、Ni、Cu为一类，Cd、Pb、Zn为一类，而As、Hg各为一类。这四类，各类内重金属的相关性相对较高，其污染可能会来源于同一性质的或关联性较高的活动和污染因素，从而推测出重金属污染的主要原因。

5.2.2 对比分析重金属污染的主要原因

基于单因素污染指数计算公式和等级划分规则，做出该城区各类区关于各重金属的污染指数和污染指数等级表，如表5.2-3：

由上表分析看出，该城区的各金属浓度均高于自然区背景值，因此其重金属污染主要是由于一些人为活动造成的。

Cr 在生活区属于中度污染，主要源于劣质化妆品、皮革及镀铬金属部件等废物堆放，而Cu 污染主要集中于生活区、工业区和主干道路区，特别处于工业区，铜污染指数高达9.66，主要来源是铜锌矿的开采和冶炼、金属加工、机械制造、钢铁生产等。Ni 污染虽高于母质，但几乎属于轻度污染，一般来源于镀镍、铸币、制造催化剂和玻璃陶瓷和不锈钢和抗腐蚀合金的制造与堆放。并且，Cr 与Ni 的相关性很高，很可能大多源于性质较为相同的工业制造品。因此，该城区Cr 、Cu 和Ni 污染主要由于工业活动以及工业类制造品的废物堆放。

Cd 、Pb 和Zn 主要在工业区和主干道路区污染严重，多分布于工矿的周围和公路、铁路的两侧，源于工业生产、汽车尾气排放及汽车轮胎磨损产生的有害气体和粉尘。

As 在工业区属于中度污染，主要源于化学制药、玻璃工业脱色剂、硬质合金、杀虫剂、砷酸盐药物、化肥以及农药等，从而对公园绿地区等起到轻度污染影响。

Hg 在工业区、主干道路区和公园绿地区的污染属于重污染，特别在工业区和主干道路区，污染指数分别达到了18.35和12.93。

由此看来，该城区的重金属污染主要来源于工业活动及其生产品和交通主干道污染。

5.3问题三的求解

5.3.1 有衰减的扩散数学模型的建立 （1）建立基本模型

根据题目所给的条件，寻找土壤中的污染源，我们需要忽略大气、水文等影响，仅考虑在土壤中的有衰减扩散。以污染源为中心，以辐射的形式向四周扩散。因此，我们根据有衰减的扩散数学模型和扩散定律，可以建立以下的微分方程：

??

??

?---=-??+??+??=??)

()()()0,,,(0002

2

2

2222222z z y y x x M z y x u c k z c c y c b x c a t c δδδ ① 其中，M 为扩散源的质量，根据质量守恒定律，我们取待求金属的最大质量。2a ，2b ，

2c 分别为x 轴，y 轴和z 轴上的扩散系数。2k 为衰减系数。

（2）化简模型的形式，以便求解

得到②式后，我们将它进一步化简，以便由matlab 程序，得到污染源的坐标。同时令t=1，则将②式化简为：

??

???---???--=

22

202202202

-4)-4)(4)(exp 8),,(k c z z b y y a x x a M z y x c （ππ ② 两边区对数，得：

??

????++-+---=22

202202202

34)-(4)(4)()ln()

2(ln

),(ln k c z z b y y a x x abc M y x c π ③

5.3.1 有衰减的扩散数学模型的求解

（1）将③试中除x,y,z 外的数据看做参数，即化简为：

z z y y x x y x c 6254232210),(ln βββββββ++++++= ④ 其中：

0β为常数项，2β=-20x 1β，4β=-20y 3β，6β=-20z 5β。

也就是说，利用matlab 软件对方程进行回归分析，得出参数的值，即可得到污染源的坐标（0x ，0y ，0z ）。

由于不同元素的污染源头不止一个，因此我们根据图中所示的污染源，取值计算它的污染源具体坐标。

（2）求解As 的污染源。

由图可知，As 的污染源共有三个。因此，我们自左而右将其分别记为污染源1、污染源2、污染源3. 如图所示：

（3）计算污染源1所取的数据如下：

由此计算得：0x =6190.089；0y =11152.52.且参数可信度为95%。 （4）利用相同的方法，计算出As 污染源2、3的坐标如下：

由此得：0x =12780.85；0y =3388.884.且参数可信度为：93%。

由此得：0x =18217.61；0y =9094.62.且参数可信度为99.9%。

（5） 利用同样的方法，计算出其余7种重金属的污染源。结果如下表：

5.4 问题四的求解

5.4.1 对于问题三所建立的模型的评价 （1）模型的优点

针对问题三所建立的模型，我们通过有衰减的扩散数学模型，同时利用回归分析，得出了不同金属的几个污染源的准确位置；模型三中，由于取样所得的点是离散的，因此在取样的时候容易出现偏差较大的点。我们选择缩小范围，取具有代表性的点，拟合出符合实际要求的浓度与位置的曲线，最后由浓度最大的点视为污染源。 （2）模型的缺点

按照实际来说，不同的金属元素的传播方式不同，如Cd 、Pb 和Zn 主要通过大气传播，As 主要通过污水传播。我们在建立模型和求解的过程中，均忽略了大气传播、金属元素沉淀等因素，单纯考虑金属在土壤中的四周辐射式扩散，因此模型具有局限性。 5.4.2 关于研究城市地质环境的演变模式 （1）信息采集

从实际出发，若要探究城市地质环境的演变模式，必须加入时间为条件。在时间演变下，探究地质环境的演变模式。

除此之外，必须从深层的传播条件出发。土壤的重金属含量与当地水文条件、气候、所属功能区以及人类活动有关，因此在采集信息的时候，需要考虑土壤中的湿度，当地的水文状况，这两个因素都能促进重金属元素扩散的速度，由此计算出重金属的扩散系数。

不同的功能分区对土壤中的重金属含量也有影响。因此需要给功能区赋值，主要考察该区排放的重金属元素对土壤含量的影响。

5.4.3模型的准备

（1）在模型三的基础上，加入时间段，记为T 。在T 时间内，测量不同数据的土壤重金属浓度和经纬度，得到大量数据。

（2）考察收集样本地点的水文条件，计算河流的流速为s 。

（3）根据不同元素，和不同功能区的活动活跃情况，给不同功能区赋值，以表示此种元素在土壤中的浓度与人类活动的关系，记为P 。如：在工业区中，Cu 的污染指数最高，在研究Cu 的污染时，令有以下计算方法：

废水

的含该区年均工业生产排除废水

的含收集样本时间段内排除Cu Cu P =

5.4.4 模型的建立

已知模型三中的有衰减扩散数学模型，可知：

??

???---???--=22

202204)(4)(exp 8),(k b y y a x x ab M

y x c ππ 修改后得：

??

???+---???--=T s k T b y y T a x x T Tab PM

T y x c )(4)(4)(exp 8),,(2

2

20220ππ 由此，便可计算和预测某城市的某种重金属在土壤中的含量，即某种重金属的污染

程度。

六、模型优缺点及推广

模型的优缺点：

（1）模糊综合评价模型能较合理地综合分析各类区重金属污染程度的，但对于权重和隶属函数的选择比较模糊，因此本文基于各金属元素权重相等的前提下，使用简单统计方法构造隶属函数，虽在本文涉及的问题和大量数据中是合理的，但现实中各类重金属对土壤的污染作用各不相同，因此该模型具有一定的局限性。

（2）通过SPSS 软件，R 型聚类分析较容易实现多变量的聚合分类，推断出哪几类重金属比较相对相关，但其本身没有函数模型，结果不如主成分分析详细精确，而基于问题重金属污染的主要原因说明， R 型聚类分析就足以解决问题。

（3）有衰减的扩散数学模型非常符合重金属在土壤中的传播特征，但不考虑其它水流、风力的影响，并需结合Matlab 软件，对多元非线性回归方程进行多次试验后，选择置信度较高的参数初值，从而推出污染源的位置，因此该模型的求解有一定的难度和误差。

模型的推广：

模糊综合评价模型可以推广到涉及多个因素或多个指标的事物或系统的综合评价，评价其优劣好坏等程度；而R 型聚类分析可以针对于相互之间有一定相关的多变量指标分析问题，可将其分成聚分成几个大类，从而找出主要指标；有衰减的扩散数学模型适合一切满足分子扩散规则的问题，通过数据，可找出点源，也可实现从点源出发的扩散

过程问题解决。

七、参考文献

［1］魏志远.海南省文昌水稻土重金属污染评价［J］.安徽农业科学，2009-08-10. ［2］李少荣.关于确定表层土壤重金属污染源的研究［J］.河南科学，2012,3.

［3］韩中庚.《数学建模方法及其应用》,高等教育出版社，2005

［4］姜启源.《数学模型》（第四版）,高等教育出版社，2010

［5］魏秀国.广州市菜园土和蔬菜中镉含量水平及污染评价［J］.土壤与环境，2002,11(2):129-132.

［6］金忠青.用脉冲谱-优化法求解对流-扩散方程源项控制反问题［J］.河南大学学报，1992,3(2).

【附录】

（1）A=[74 781 5 7.84

………

5985 2567 44 7.56

7653 1952 48 9.35 ]; %X、Y、海拔和元素浓度

x=A(:,1);y=A(:,2);z=A(:,3);c=A(:,4)

[X,Y,Z]=griddata(x,y,z,linspace(min(x),max(x))',linspace(min(y),max(y)),'v4');%插值

[X,Y,C]=griddata(x,y,c,linspace(min(x),max(x))',linspace(min(y),max(y)),'v4');%插值

figure,meshc(X,Y,Z,C)

caxis([min(c) max(c)]);

xlabel('X(m)');

ylabel('Y(m)');

zlabel('海拔(m)');

title('As立体分布图');

clear all;close all;clc;

x=[1 7776 10613 9 2.329227

1 7458 8920 36 1.900614

1 7405 10981 6 1.521699

1 6423 8831 40 1.214913

1 6395 10443 4 1.962908

1 5567 678

2 7 2.354228

1 5438 3994 10 2.235376

1 5006 8846 6 2.159869

1 4948 7293 6 2.911263

1 474

2 729

3 9 3.085116];%元素数据

xx=x(:,1:4)

yy=x(:,5)

beta=[ 9.38648061866336 -0.00402303999797 0.00000014548075 0.00068845595413 -0.00000005250534

0.02252737015359 -0.00108482795115];%初值

[beta,R,J]=nlinfit(xx,yy,@user_function,beta()) %求解系数

betaci=nlparci(beta,R,J) %置信区间

（2）其他金属元素空间分布图

图5.1.1-2：Cd空间分布图

图5.1.1-3：Cu空间分布图图5.1.1-4：Ni空间分布图

图5.1.1-5：Pb 空间分布图 5.1.1-6：Zn 空间分布图

图5.1.1-7：Cr 空间分布图 图5.1.1-8：Hg 空间分布图