关于城市表层土壤重金属污染的研究

承 诺 书我们仔细阅读了中国大学生数学建模竞赛的竞赛规则.我们完全明白，在竞赛开始后参赛队员不能以任何方式（包括电话、电子邮件、网上咨询等）与队外的任何人（包括指导教师）研究、讨论与赛题有关的问题。我们知道，抄袭别人的成果是违反竞赛规则的, 如果引用别人的成果或其他公开的资料（包括网上查到的资料） ，必须按照规定的参考文献的表述方式在正文引用处和参考文献中明确列出。我们郑重承诺，严格遵守竞赛规则，以保证竞赛的公正、公平性。如有违反竞赛规则的行为，我们将受到严肃处理。我们参赛选择的题号是（从 A/B/C/D 中选择一项填写）： 我们的参赛报名号为（如果赛区设置报名号的话）： 所属学校（请填写完整的全名）： 参赛队员 (打印并签名) ：1. 2. 3 指导教师或指导教师组负责人 (打印并签名)： 日期：2013 年 8 月 28 日赛区评阅编号（由赛区组委会评阅前进行编号）：2013 全国大学生数学建模竞赛编 号 专 用 页赛区评阅编号（由赛区组委会评阅前进行编号）：赛区评阅记录（可供赛区评阅时使用）：评阅人评分备注全国统一编号（由赛区组委会送交全国前编号）：全国评阅编号（由全国组委会评阅前进行编号）：关于城市表层土壤重金属污染的研究摘 要本文通过主要重金属元素在该城区的空间分布，分析该城区内不同区域重金属的污染程度和原因，进一步由重金属污染物的传播特征确定污染源，从而为减少污染提出对应的解决措施。针对问题一，首先运用 软件绘制出该区域地势分布图和 5 个功能区在城区matlb的确切分布，并通过每种重金属污染等值线图可以大致得出每种重金属污染在各个功能区污染的程度，由图可知生活区、工业区、主干道路区所受重金属污染较重，并且(砷)、 (镉)、 (镍) 在城市中污染范围分布较广。为了进一步对重金属污染AsCdNi进行量化分析，采用污染负荷指数分析法，将重金属的污染程度分为无污染、中等污染、强污染、极强污染四个等级，最后根据污染负荷指数，计算出五个功能区的污染负荷指数分别为：2.348，1.921，1.834，1.578，1.060，因此受污染程度较重的前三个功能区依次为工业区，主干道路区，生活区。针对问题二，首先通过表中重金属污染负荷指数与污染程度的关系将5个区的受污染程度进行分类，结合问题一中重金属元素的空间分布图及五个功能区的地理位置，从而得出各地区受重金属污染的主要原因。分析结果表明：工业区的重金属污染是最严重的，其次是生活区和山区，受污染程度相对较轻的是主干道区和公园绿地区。针对问题三，为了确定污染源的位置，首先针对 8 种重金属元素的污染特征进行分析，发现重金属污染主要是通过扩散进行传播。然后运用 软件计算所给数据matlb中海拔与污染浓度的相关系数，结果表明污染浓度与海拔相关性较弱，因此可以建立采样点的坐标和重金属浓度的二维传播模型，根据距离污染源越近重金属浓度越大的原理，通过求解偏微分方程，从而求得各种重金属元素的污染源的个数以及坐标。以(砷)为例，计算出它的污染源有三个，它们的坐标分别是（18134，10046，41） 、As（12696，3024，27） 、 （4742，7293，9） 。 针对问题四，首先针对二维传播模型进行了评价，肯定了该模型能够快速准确的确定污染源，但是并未联系土壤的渗透以及水源流动对污染造成的影响。为更好地研究城市地质环境的演变模式，在原有模型基础上，结合渗流速度和梯度方向，得到改善后的二维对流方程，进而能够确定更为精确的污染源位置。以 (砷)为例得到污染As源坐标为(4507.5，7495.7)、(12506，2994)、(18500，10491.4)。最后结合该区域的5 个功能区对城市环境的影响，并分析污染源所处的功能区，因此可以将污染源缩小为两个点(4507.5,7495.7)、(18500,10491.4)。关键词：等值线分布图；污染负荷指数法；二维传播模型；二维对流方程一、问题重述随着城市经济的快速发展和城市人口的不断增加，人类活动对城市环境质量的影响日显突出。对城市土壤地质环境异常的查证，以及如何应用查证获得的海量数据资料开展城市环境质量评价，研究人类活动影响下城市地质环境的演变模式，日益成为人们关注的焦点。按照功能划分，城区一般可分为生活区、工业区、山区、主干道路区及公园绿地区等，分别记为 1 类区、2 类区、5 类区，不同的区域环境受人类活动影响的程度不同。本题中给出某城市五大区重金属含量的相关资料，要求建立相应的建立数学模型分别讨论下列问题：问题一：给出所需研究重金属的空间分布，并分析该城区不同区域重金属的污染程度。问题二：分析数据并说明重金属污染的主要原因。问题三：分析重金属污染物的传播特征，并建立模型确定污染源的位置。问题四：分析所建模型的优缺点，讨论为更好的研究城市地质环境演变模式所需搜集的信息，已经应该如何利用这些信息建立模型解决问题。二、问题分析对于问题一，我们首先分析所给区域的地貌情况，然后运用 软件描绘城区matlb的五种功能分布图和八种重金属元素浓度的等值线图，将两种分析图结合在一起，能够清晰地得到重金属元素的空间分布状况，进一步利用污染负荷指数法，对整个区域每个功能区的各种重金属进行了定量评价，根据污染负荷指数与污染程度关系分级标准，最后可以得到知该区域内不同功能区的重金属污染程度。对于问题二，首先我们利用表 2 中重金属污染负荷指数与污染程度的关系将 5 个区的受污染程度分别进行分类，这样就能清楚地了解到各地区受重金属污染的情况。再结合该地区的实际地理位置分布，查阅关于重金属危害的相关资料，进而得出这五个区受重金属污染的主要原因。 对于问题三，首先我们对重金属污染物传播特征的分析是否能够得出污染物浓度大小与距污染源的距离和所处的地理特征等因素有关。我们建立模型确定污染源位置之前，先分析海拔高度与污染物浓度的相关性，确定海拔和污染浓度无关时，然后建立二维传播模型，根据距离污染源越近污染浓度越大的规律，求拟合函数的极值，即可得到污染源的坐标对于问题四，通过问题三的求解，我们可以简便准确的找到各种重金属的污染源。我们认为重金属的污染传播与其位置关系有关，重金属的金属含量浓度越高，这个点为即为污染源。地表水的流动性作用是重金属在土壤表层中传播的主要原因，故假设重金属传播以对流方式为主，由此设二维对流方程，并根据等高线估算其梯度方向的余弦，得出的其他七种重金属元素的污染源坐标，结合等高线和这种元素的等值线分布图确定重金属污染源所在的功能区结合重金属的污染的特征规律，得到更为精确地污染源的位置。三、模型假设(1) 由于考虑土壤和水流对重金属污染物传播的影响，所以不考虑风力和重力的影响；(2)该城区在长时间内没有发生重大污染事故； (3)假设污染物是在一个时刻以点源的形式传播；(4)不考虑取样时间差的对结果所造成的影响； (5)取样点的数据较好的反映了该地区的污染物浓度。四、符号说明：元素 的实测含量iCi：元素 的评价标准，即背景值i：元素 的高级污染系数iFi：某一点的污染负荷指数PLI：区域污染负荷指数N：不同的金属元素i：相关系数i：浓度的标准差：拟合点与污染源的距离l：污染源坐标yx,：重金属的坐标及浓度tC：渗流速度：渗流速度沿 坐标轴的分量ax：渗流速度沿 坐标轴的分量by：采样的时刻0t五、模型的建立与求解5.1 问题一的模型建立与求解随着城市经济的快速发展和城市人口的不断增加，人类活动对城市环境质量的影响日显突出。对城市土壤地质环境异常的查证，研究人类活动影响下城市地质环境的演变模式，日益成为人们关注的焦点。为了更清楚的反映重金属污染物同功能区分布的关系，采用地势图来体现该城区的地势分布情况，以经纬为纵横坐标，海拔为竖坐标，可利用 软件 （详见附matlb4录程序 1）真实的还原该城区的地貌环境，以便能够进行对污染原因进行全面分析。等值线图能够反映出某种元素在城市的分布情况和污染程度。 5.1.1 重金属元素的空间分布由附件所给的数据，在已知不同功能区检测点坐标值的情况下，以经纬为纵横坐标，海拔为竖坐标，我们可以得到该地区的地势分布图： 图 1 地势分布图由地势分布图 1 可以看到，该地区总体东北方向地势较高，多为山丘和丘陵地带，延西南方向，地势逐渐放缓，西部多广阔平原，中部为小盆地。运用 软件（详见附录程序 1）做出该城区的生活区、工业区、公园绿地区、matlb交通区、公园绿地五个功能区的分布图，如图 2。从功能区分布图中我们可以清楚地看到每个区域大致所在的范围，以便于我们进行分析，从而更清楚的了解重金属污染发生的区域。我们结合图 1 和图 2 进行初步地理分析，由于东部多山丘和丘陵，不利于生活和工业制造，因此多为绿地，而西部地势平坦广阔为平原地带，进行生产活动条件便利，因此主要分布着生活区、工业区以及主干道路区，随着人类居住以及工业生产该地区产生了大量的工业废气、废渣以及生活垃圾等等，环境污染情况比较严重。中部的盆地地区，由于地势较低，空气流动性较差，在此周围进行的各种人类活动产生的“三废” ，污染物很难扩散，极易引起重金属污染物的沉淀。0 0.51 1.52 2.53x 10400.511.52x 104-2000200400x(m)y(m)业业业m业-10001002003000 0.5 1 1.5 2 2.5 3x 10400.511.52x 104x(m)y(m)业业业业业业业业业业业业业业业业图 2 功能区分布图5 个功能区 8 种重金属元素使用 软件（详见附录程序 1）描绘出的 8 种重金atlb属元素在该区的平面分布图（图 3-图 10） ，以浓度作为等值线，颜色越深表示浓度越大，也就表明该城市地区的受重金属污染的程度越重：图 As 的污染程度 图 3 的等值线分布 图 4 的等值线分布 As Cd从图 3 我们可以清楚的发现， （砷）在该区域分布较广，公园绿地区的As（砷）的等值线比较稀疏，表明 的浓度较小，所以在公园绿地区 的含量较少。s As而在生活区、主干道路区和工业区的颜色较深，表明 在生活区、主干道路区和工业As区的含量很多。从图 4 的重金属含量浓度的等值线可以看出，该地区 （镉）污染十分严重。山d区的 （镉）的等值线颜色较浅，分布比较分散表明 （镉）的浓度较小，所以在Cd C山区 （镉）的含量较少。而生活区、主干道路区、工业区的颜色较深分布密集，表x(m)y(m)0 0.5 1 1.5 2 2.5 3x 1045001001500x(m)y(m)0 0.5 1 1.5 2 2.5 3x 1045001001500明 （镉）在主干道路区和工业区的含量很多。Cd图 5 的等值线分布 图 6 的等值线分布Cr Cu从图 5 中我们可以看到，山区以及公园绿地的 （铬）的等值线分布稀疏， 的r Cr浓度较小，所以在山区以及公园绿地的 （铬）的含量较少。而在生活区和主干道路Cr区的颜色较深，等值线密集，表明 (铬）在生活区和主干道路区的含量很多，浓度很大。从图 6 的 （铜）含量浓度的等值线可以看出，山区和公园绿地的 （铜）的Cu u等值线几乎没有，这两个地区的 （铜）含量几乎没有。而生活区和工业区的颜色较u深，尤其是生活区， （铜）分布格外密集，表明 （铜）在生活区和工业区的含u量很多。图 7 的等值线分布 图 8 Ni 的等值线分布Hg从图 7 中我们可以看到，山区、主干道路区以及公园绿地的 （汞）几乎没有分Hg布，所以在山区、主干道路区以及公园绿地的 （汞）的含量较少或者几乎没有。而Hg在生活区和工业区的颜色较深，等值线密集，表明 （汞）在生活区和工业区的含量很多，浓度很大。x(m)y(m)0 0.5 1 1.5 2 2.5 3x 1045001001500x(m)y(m)0 0.5 1 1.5 2 2.5 3x 1045001001500x(m)y(m)0 0.5 1 1.5 2 2.5 3x 1045001001500x(m)y(m)0 0.5 1 1.5 2 2.5 3x 1042004006008001001200140016001800x(m)y(m)0 0.5 1 1.5 2 2.5 3x 1045001001500从图 8 的 （镍）含量浓度的等值线可以看出，该地区基本都被 （镍）所污染，Ni Ni但公园绿地区的 （镍）的等值线稀疏，公园绿地区的 （镍）浓度较小。而生活区、i山区、工业区和主干道路区的颜色较深，尤其是生活区 （镍）分布格外密集，表明（镍）在生活区和工业区的含量很多。i图 9 的等值线分布 图 10 的等值线分布 Pb Zn从图 9 的 （铅）含量浓度的等值线可以看出，该地区基本都被 （铅）所污Pb染，但公园绿地区和山区的 （铅）的等值线稀疏，这两个地区的 （铅）浓度较Pb小。而生活区、工业区和主干道路区的颜色较深，尤其是生活区 （铅）分布格外密集，颜色较深，表明 （铅）在生活区含量最多，污染最严重。从图 10 的 （锌）含量浓度的等值线可以看出，山区的 （锌）的等值线几乎Zn Zn没有，山区的 （锌）基本未被 （锌）所污染。而生活区、工业区和主干道路区Zn的颜色较深，尤其是生活区和工业区， （锌）分布格外密集，表明 （锌）在生活区和工业区的浓度很大，污染十分严重。由图 3-图 10 我们看到 基本覆盖了整个地区，该铅 ）（ 镍 ）（ 镉 ）砷 ） (,(PbNiCdAs受这四种重金属污染范围较广，污染程度较深； （汞）的分布较少，主要集中于生Hg活区、工业区，其他地区基本未受到污染。但是工业区、生活区、山区、主干道路区都被这八种重金属所污染，并且污染程度较深，考虑到金属污染的难降解性，这个地区重金属污染的治理迫在眉睫，否则污染程度会越加深刻，环境会急剧恶化，给居民生活带来更为严重的伤害，因此我们需要结合污染分布图和功能分区图，分析污染产生的原因，从源头上治理，才能取得较好的治理效果。通过上述分析，该地区的重金属污染在 5 种功能区的空间分布情况大致可由表 1表示：x(m)y(m)0 0.5 1 1.5 2 2.5 3x 10450001000015000表 1 重金属空间分布情况重金属种类 空间分布（砷）As 主干道路区、生活区、山区（镉）Cd主干道路区、工业区、生活区（铬）r 工业区、主干道路区（铜）u 工业区、生活区（汞）Hg 工业区、生活区（镍）Ni 主干道路区、工业区、生活区（铅）Pb主干道路区、工业区、生活区（锌）Zn主干道路区、工业区、生活区通过表 1 我们可以发现，该地区的主干道路区、生活区、工业区的重金属污染较重，应当加强治理，以防止污染的进一步恶化和扩散。5.1.2 污染负荷指数法土壤重金属污染程度的评价分别采用研究土壤中重金属含量的自然背景值和国家一级标准值为保护区域自然生态，维持自然背景的土壤环境质量的限制值为参考标准，当评价标准为背景值时，相应的评价结果称之为元素富集，当以国家一级标准为评价标准时，称之为元素污染 。污染负荷指数法是在从事重金属污染水平的分级研究中6提出来的一种评价方法，它被广泛应用于土壤重金属污染的评价，某一点的污染负荷指数的公式如下：/iiCF(1)nnPLI321其中， 为元素 的高级污染系数， 元素 的实测含量，单位为 ； 为元素iFi i gu/iC的评价标准，即背景值，一般选取该城区自然表层土壤中元素的平均值作为背景值。i为某一点的污染负荷指数； 为评价元素的个数，而某一区域的污染负荷指数为：PLI nPLPLPL nN III321其中 为区域污染负荷指数； 为评价点的个数（即采集样点的个数） ，污染负荷指NPLI数法的优点是该指数有评价区域所包含的多种金属成分共同构成。通过这种方法对整个区域各个点位各种重金属进行了定量评价，并对各点的污染程度进行分级，以重金属污染负荷指数与污染程度关系作为评价标准，见表 2：表 2 重金属污染负荷指数与污染程度的关系污染程度 无污染 中等污染 强污染 极强污染值PL1133污染等级 0 根据重金属污染负荷指数与污染程度关系，根据某一点的污染负荷指数不同分为无污染，中等污染，强污染，极强污染四个等级，分别用符号 0，来表示，根据污染负荷指数法计算得到该区域内不同区域的污染负荷指数值，参照污染负荷指数与污染程度关系分级标准，可知该区域内不同区域重金属污染程度见表：表 3 不同地区的重金属污染程度元素 指数值PL污染等级 污染程度生活区 1.833561244 中等污染工业区 2.348448408 强污染山区 1.060221 中等污染主干道路区 1.920920706 中等污染公园绿地区 1.57799 中等污染该城市的污染十分严重，不仅空间分布较广涉及到生活区、工业区以及主干道路区，而且污染的程度很深，金属污染不同于其他的污染，金属污染较难解决，污染造成的后果严重，应当引起当地部门的重视，加强治理和预防，以避免更深程度的污染和污染区域的进一步扩大。5.2问题二的模型建立与求解为说明重金属污染的主要原因，我们需要利用问题一中统计出各地区的重金属污染情况，及表 2 中重金属污染负荷指数与污染程度的关系来对生活区、工业区、山区、主干道道路区、公园绿地区这五个地区进行逐个分析。详细分析如下:表 4 生活区各种金属污染情况表重金属 AsCdrCuHgNiPbZn样本量 44 44 44 44 44 44 44 44无污染（个） 5 7 4 2 13 7 6 5无污染（%） 11.36 15.91 9.10 4.55 29.55 15.91 13.64 11.36中等污染（个） 25 17 31 16 12 31 24 20中等污染（%） 56.82 38.64 70.45 36.36 27.27 70.45 54.55 45.45强污染（个） 12 13 5 8 8 6 8 6强污染（%） 27.27 29.55 11.36 18.18 18.18 13.64 18.18 13.64极强污染（个） 2 7 4 18 11 0 6 13极强污染（%） 4.55 15.91 9.10 40.91 0.25 0 13.64 29.55由表 4 统计的结果可知，生活区除了 Ni 元素之外，其它元素都有让部分地区受到极强污染，其中 Cn 污染最严重，受到极强污染的地区占 40.91%。As，Hg 较少些。每种重金属都有 10%左右地区处于强污染，处于无污染与污染临界区所占范围较大，因此在生活区这些重金属的含量都需要引起重视。表 5 工业区各种金属污染情况表重金属 AsCdrCuHgNiPbZn样本量 36 36 36 36 36 36 36 36无污染（个） 5 2 7 1 5 6 0 1无污染（%） 13.89 5.56 19.44 2.28 13.89 16.67 0 2.78中等污染（个） 16 9 21 7 9 21 16 13中等污染（%） 44.44 25 58.33 19.44 25 58.33 44.44 36.11强污染（个） 12 11 6 9 4 8 10 14强污染（%） 33.33 30.56 16.67 25 11.11 22.22 27.78 38.89极强污染（个） 3 14 2 19 18 1 10 14极强污染（%） 8.33 38.89 5.56 52.78 50 2.78 27.78 38.89工业区有大概一半地方处于极强污染区，这些地方主要集中 Cu，Cd，Pb，Zn 这些重金属，其中 Cn 和 Hg 的极强污染最为严重，都达到 50%以上，Ni 的极强污染程度最轻。Pb 污染了整个工业区，该地区 Cn，Hg，Pb 的含量较高需引起高度重视。另外，该地区整体上来看，工业区的污染比较严重，需要引起人们的极大关注，应尽快采取措施，否则污染程度会进一步加大。表 6 山区各种金属污染情况表重金属 AsCdrCuHgNiPbZn样本量 138 138 138 138 138 138 138 138无污染（个） 23 16 16 2 45 24 17 16无污染（%） 16.67 11.59 11.59 1.45 33.33 17.39 12.32 11.59中等污染（个） 93 41 99 45 35 104 59 42中等污染（%） 67.39 29.71 71.74 33.33 25.36 75.36 42.75 30.43强污染（个） 19 37 14 20 14 8 36 36强污染（%） 13.77 26.81 10.14 14.49 10.14 5.80 26.09 26.09极强污染（个） 3 44 9 71 44 2 26 44极强污染（%） 2.17 31.88 6.52 51.45 31.88 1.45 18.84 31.88山区处于无污染和污染临界的面积占大多数，应对这些地区采取一定的防范措施，避免污染情况加重。其中受到 元素极强污染的区域达到 50%以上，应采取适当措施。CuNi 的污染力度最小。表 7 主干道道路区各种金属污染情况表重金属 AsdruHgNiPbZn样本量 66 66 66 66 66 66 66 66无污染（个） 36 33 33 31 30 33 34 36无污染（%） 54.55 50 50 46.97 45.45 50 51.51 54.54中等污染（个） 24 26 24 28 31 25 29 28中等污染（%） 36.36 3.03 36.36 42.42 46.97 37.88 43.94 42.42强污染（个） 5 5 8 5 4 6 1 1强污染（%） 7.58 7.58 12.12 7.58 6.06 9.09 1.52 1.52极强污染（个） 1 2 1 2 1 2 2 1极强污染（%） 1.52 3.03 1.52 3.03 1.52 3.03 3.03 1.52主干道路区处于无污染状况的区域相对较大，都占到了 30%左右。处于污染临界值的区域也占了大多数，但是处于强污染，极强污染的区域较小，应对没有受到污染的地区继续保持，其他地区采取一定的措施。表 8 公园绿地区各种金属污染情况表重金属 AsCdruCHgNiPbZn样本量 35 35 35 35 35 35 35 35无污染（个） 4 6 2 2 11 9 2 10无污染（%） 11.43 17.14 5.71 5.71 31.43 25.71 5.71 28.57中等污染（个） 22 18 30 16 10 23 24 17中等污染（%） 62.86 51.43 85.71 45.71 28.57 65.71 68.57 48.57强污染（个） 8 4 2 15 4 3 3 4强污染（%） 22.85 11.43 5.71 42.86 11.43 8.51 8.51 11.43极强污染（个） 1 7 1 2 10 0 6 4极强污染（%） 2.86 0.2 2.86 5.71 28.57 0 17.14 11.43公园绿地区除了 Ni 元素以外，均有地区达到极强污染，特别是 Hg，Pb 对该地区造成严重污染。处于无污染和污染临界的面积相对较大，应继续保持。表 9 各区污染程度综合评价表生活区 工业区 山区 主干道路区 公园绿地区中度污染 极强污染 强污染 中度污染 中度污染表 10 各区污染情况区域 生活区 工业区 山区 主干道路区 公园绿地区重金属污染情况描述As、Cr、Ni、Pb 中度污染，Cu、Zn重度污染，其中 Cu 的污染最严重。Ni 中度污染，其他全部重度污染，Cu，Hg，Zn的极强污染程度也非常严重全部中度污染，Cu、Zn极强污染严重未受污染的地区很多，其他基本在未污染和污染的临界值上大部分区域处在临界区，Hg,Pb,Zn 有极强污染分析生活区的情况：生活区是人们主要居住地，而周围有大量工厂，这些工厂排放产生大量有害重金属，这些重金属向四周环境逐渐扩散，对环境造成恶劣污染，同时，废弃物在雨水的冲洗下，也会使重金属的迁移能力增强，因此，城市垃圾堆的聚集地成了重金属污染的潜在来源。此外，随着城市化不断发展，这些污染企业应该如何可持续发展也成了亟待解决的问题，当然，现在生活中也有很多生活用品含有大量重金属，例如染发剂，化妆品，电池等等。分析工业区的情况：工厂作业产生的工业垃圾均含有大量重金属，该地区可能有重金属冶炼厂，化工厂，印刷厂等，所以该地区重金属均重度污染。分析山区的情况：山区虽然远离市区，可是厂商可能会将工业垃圾存放在这里，因此，垃圾中的重金属会向四周环境扩散，垃圾焚烧也会产生大量粉尘，粉尘中将会含有大量的重金属颗粒，所以山区也不是想象中那样重金属含量很低，相反，山区也慢慢地被污染了。分析主干道区的情况：这里的重金属主要是 Pb，Zn，Cd，Cu，Cr 等，他们主要来自含铅汽油的燃烧、汽车轮胎磨损产生的含锌粉尘等。所以会有一些重金属污染，不是过于严重是因为主干道的绿化相对较好一些。分析公园绿地区的情况：由分析可知公园绿地也受到了一定成的重金属影响。因为工厂产生的工业废水废气未经彻底处理就直接排放进附近的河流和大气中，这些重金属随着大气和水流流动，最终进入土壤。此外，农药化肥中也含有大量重金属，这也加剧了公园绿地区重金属污染的情况。综合上述分析，工业区的重金属污染是最严重的，其次是山区和生活区，主干道区和公园绿地区好一点。工业区重金属含量超标可以理解因为工厂作业产生的各种垃圾均含有大量重金属，生活区污染比较严重是因为一方面离工业区较近，也没有采取一定的防范措施；另一方面现在越来越多的生活用品中含有一定量的重金属。而山区则是因为很多污染严重的工厂将垃圾堆放在这里，日积月累，重金属污染日益严重。主干道和公园绿地区稍好一些，因为有很多绿色植物，特别是远离工厂区的公园绿地区，相对来说是受重金属污染最轻的。5.3问题三的模型建立与求解目前工厂排放重金属污染物基本是通过废渣废水废气三种途径进行排放，因此在考虑金属污染物的污染特征时，我们认为重金属元素的传播途径可以分为三类，一类是在土壤中传播，一类是通过河流传播，还有一类是在大气中传播，分别对应工厂排放废物的废渣、废水和废气。通过以上对重金属污染物传播特征的分析可知：污染物浓度大小与距污染源的距离和所处的地理特征等因素有关。我们在建立模型确定污染源位置之前，先分析海拔高度与污染物浓度的相关性。5.3.1 重金属污染物的污染特征(砷)的污染特征 :砷主要来自工业生产及含砷农药的使用、煤的燃烧，都可As3产生含砷废水、废气和废渣，对环境造成污染。(镉)的污染特征：镉的污染主要来源于铅、锌、铜的矿山和冶炼厂的废水。Cd由于土壤的吸附作用很强，镉几乎不发生向下的再迁移，因而累积于土壤表层，但是在降水的影响下，土壤表层的镉的可溶态部分随水流动就可能发生水平迁移，进入界面土壤和附近的河流或湖泊而造成次生污染。（铬）的污染特征:土壤中铬的污染主要来源于铁、铬、电镀、金属酸洗、皮r革鞣制、耐火材料、铬酸盐和三氧化铬工业的“三废”排放进入土壤。铬主要残留积累于土壤表层，由于铬在土壤中多以难溶性且不能被植物所吸收利用的形式存在，因而铬的生物迁移作用较小。(铜)的污染特征：主要污染来源是铜锌矿的开采和冶炼、金属加工、机械制Cu造、钢铁生产等。在靠近铜冶炼厂附近的土壤，含有高浓度的铜。岩石风化和含铜废水灌溉均可使铜在土壤中积累并长期保留。(汞)的污染特征：土壤中汞的氧化还原，可从土壤中挥发进入大气环境，被植Hg物的根系和叶片吸收，并随着土壤温度的升高，其挥发的速度加快，而被土壤胶体吸附的能力相对较弱，因此二甲基汞较易进行气迁移和水迁移。（镍）的污染特征：冶炼镍矿石及冶炼钢铁时，部分矿粉会随气流进入大气。Ni在焙烧过程中也有镍及其化合物排出，镍可在土壤中富集。（铅）的污染特征：土壤中铅的污染主要来自大气污染中的铅沉降和铅应用工Pb业的“三废”排放，主要是通过空气、水等介质形成的二次污染（锌）的污染特征：主要污染源有锌矿开采、冶炼加工、机械制造以及镀锌、Zn仪器仪表、有机会合成和造纸等工业的排放。在酸性土壤中容易迁移，具有很强的吸附性。5.3.2 海拔高度与污染物浓度的相关性设高度的各个数据量为 ,对应的浓度值为 ，相关性系数nh321, nQ321,函数为： )var()r(,coyxQ把每种重金属元素的浓度以及各取样点的海拔高度等数据带入上式，通过软件编程（详见附录程序2）计算可以得到：matlb表11 相关性分析i12345678相关性 -0.289 -0.284 -0.152 -0.138 -0.084 -0.163 -0.235 -0.178根据此表可知，海拔高度与污染物浓度的相关性很小，即在分析污染物的传播时可忽略海拔高度的影响。5.3.3二维传播模型根据问题一做出的各种重金属污染物浓度在空间分布的状态，找出等值线分布中比较密集颜色较深的区域（且此点不应为变异点） ，将此点假想为污染源，那么在此污染源附近的区域内，污染物浓度值的大小将随距污染源距离的变化而变化，因此设： ),(),(yxFQ因为浓度还受到其他各因素的影响，修正后的不考虑海拔高度 拟合函数形式为：h2)(leQ其中 为浓度的标准差。通过 编程（详见附录程序3）可求得标准差 （ 表示不同的重金属元素，matlb i，分别表示 ）8,21i ZnPbNiHgCurdAs,表12 标准差i12345678标准差 3.0243 224.9876 70.0018 162.9151 1629.5 9.9414 50.0578 339.2325从等值线分布图我们可以发现，颜色越深，污染浓度越高，其周围的点的污染浓度也很高，由此我们假设这个点为污染源，假设污染源的平面坐标为 。因为该),(0yx区域的地形式已知的，所以只需知道平面坐标 ，即可以得到 点与污染源的),(0yx距离： 2020l由此拟合函数形式可以变为： 200),(yxeyQ求拟合函数的浓度的极值，可以利用 软件（详见附录程序4）求解偏微分方程。matlb所以令： 0yx以 为例,已知 ，假设污染源的坐标为 所以浓度可以表示砷As0243.11046,83为：(2)2043.1618),(yxeyQ求解得到 ，所以我们可以近似认为浓度越大的样本点就是重金属元素的yx0,污染源。这样我们就得到重金属元素的污染源位置如下表所示：表13 污染源位置表重金属元素 污染源个数 污染源坐标 所属区域（18134，10046，41） 主干道路区（12696，3024，27） 主干道路区As3（4742，7293，9） 工业区（21439，11383，45） 主干道路区Cd2（3299，6018，4） 主干道路区r1 （3299，6018，4） 主干道路区u1 （2385，3692，7） 工业区（2708，2295，22） 主干道路区（13694，2357，33） 主干道路区Hg3（15248，9106，16） 主干道路区（3299，6018，4） 主干道路区Ni2 （22193，12185，79） 山区（4777，4897，8） 生活区Pb2（1647，2728，6） 工业区（13797，9621，18） 主干道路区Zn2（9328，4311，24） 生活区根据上表我们可以看出重金属元素的污染源主要为主干道路区和工业区，特别是主干道路区。分析原因，一方面，主干道路区汽车较多，汽车尾气的排放是重金属气体传播的重要来源之一；另一方面，主干道路区遍布全市，污染点的坐标可能有一些误差，结果使其落入主干道路区。5.4问题四的建模与求解通过问题三的求解，我们可以简便准确的找到各种重金属的污染源。我们认为重金属的污染传播与其位置关系有关，重金属的金属含量浓度越高，这个点为即为污染源，符合我们生活中的认知情况。以浓度拟合函数为定解条件，解偏微分方程，即可以得到污染源的坐标。但是在建立模型时：(1)假设污染物是在一个时刻以点源的形式排放；(2)未考虑重金属在土壤中的吸附转化等化学过程；(3)未区分城市功能区、河流等诸多因素的影响；从计算结果来看，计算出的污染源位置和采样得到的浓度最大值点相距较近，其原因是没有考虑到河流、生活污水等承载金属污染物的流体的流动所产生的迁移，所以模型有待完善。5.4.1结合地质演变模式要准确的分析该城区的地质演变模式，首先数据需要在一定时间间隔做一次采样，也即得到动态的数据，观察同一采样点在不同时刻重金属浓度的变化情况，其次需要考虑降雨量的数据、水文地质、地下水埋藏和分布规律、土壤中的含水量等数据。地表水的流动性作用是重金属在土壤表层中传播的主要原因，故假设重金属传播以对流方式为主，所以设二维对流方程为：0yCxtC其中， 为重金属的坐标及浓度， 为渗流速度沿两个坐标轴的分量。),(tyxCyx,设地下水的渗流速率是常数，变化的只是方向。经查阅资料一般可设渗流速度为设 ，把采样的时刻设为 ，将问题三的函数 做为sm7105.bayx, 0t ),(yxQ定解条件： ，那么传播模型为一阶线性齐次偏微分方程的初值问题),(,0Qty为： ),(),(0yxQtyxCbao所以方程的通解为 ，令 ：, 000tyx得t(3)0,byaxy是初始时刻的重金属浓度的分布，此时函数的极大值点为污染源的位置。由于本问题所给的数据中没有涉及采样的时间，我们设采样时刻为污染物排放后的半年，即 。st15200由问题三我们大致得出了重金属污染的污染源位置坐标，通过 编程（详见matlb附录程序5）我们可以得到该地区的等高线分布图：x(m)y(m)0 0.5 1 1.5 2 2.5x 104050001000015000-150-100-50050100150200250300图11 等高线分布以 为例， 的浓度极大值点为 ，由图估算其梯度方向的余弦为Ass750,4则可设21,3 77105.2,108.4)23( ba带入方程(2)，可以得出这个污染源为 ，同理可以求出 的另外两个污.95,. As染源的位置为 。同理我们可以得到其他七种污染源的位置。410,85,9412506我们结合问题一中的功能区分布图2、等值线分布图3和等高线分布图11，可以确定从图中我们可以确定每种污染源所在的城区进而能够更加准确的确定最源头的污染位置，例如刚求的 的污染源的三个坐标为As7.495,.0，由此我们可以根据图11判断出他们所在的功能区分别是4.109,85,2941506生活区、主干道路区和山区，由于 （砷）主要产生于煤炭的燃烧、对土壤具有较强As的吸附作用，因此我们可以将污染源缩小为两个点 和 。7.495,.04.109,8同样的我们根据公式(3)得出的其他七种重金属元素的污染源坐标，结合等高线图11 和这种元素的等值线分布图确定重金属污染源所在的功能区结合重金属的污染的特征规律，得到更为精确地污染源的位置。六、模型的评价与推广(1)问题一将数据与图形相结合，通过两者的比对进行分析更具有说服力。本文在考虑空间分布时，忽略了海拔高度对重金属污染浓度的影响，虽然这样做在一定程度上简化了求解的繁琐性，但同时也降低了结果的说服力。(2)问题二中所建立的模型与实际紧密联系，通过利用简单的模型就能达到很好的效果，具有很好的通用性和推广性；缺点在于处理数据和求解过程中不可避免的会出现各种误差，会在一定程度上影响到模型求解的精确度。 (3)问题三的模型计算方便，与软件相结合，预测数据更为准确，且具有可信度，基本符合事实规律。分析后我们发现，这一模型有一定的偏差，但是总体上是符合我们的要求，模型通过第四问的改进，增加了流动性因素，更加贴切实际，但是模型还应该增加其他对环境污染的因素进行完善，改进后的模型结果会更加准确。本模型的建立对于研究减少城市污染问题和保护环境具有重要意义，尤其在当今高能耗高污染的生产模式为主的工业时代，该模型的建立对于研究城市规划和工厂位置选择以及交通建设时具有重要的参考价值，同时利用该模型也可以研究物质扩散规律。参考文献1顾以秋，沈洁洁，姚文旭.重金属污染程度分析J.西昌学院报.27(1)：25-27,2013.2尹君.基于绿色食品基地土壤环境质量评价方法研究J.农业环境保护.20(6):8-9,2001.3蔡杨阳，童伟娟.某城市表层土壤重金属污染的特征分布及评价J.河南科学，31(3): 363-364,2013.4刘卫国.MATLAB 程序设计与应用M.高等教育出版社.2006.5董臻圃.数学建模方法与实践M.国防工业出版社.2006.6中国环境监测总站.中国土壤元素背景值M.中国科学出版社.1990.附录程序 1：x=data(:,2); y=data(:,3); As=data(:,6); Cd=data(:,7); Cr=data(:,8); Cu=data(:,9); Hg=data(:,10); Ni=data(:,11); Pb=data(:,12); Zn=data(:,13); f1=44;f2=80;f3=146;f4=284;f5=319; figure(1) plot(x(1:f1),y(1:f1),*) hold on plot(x(f1+1:f2),y(f1+1:f2),ko) plot(x(f2+1:f3),y(f2+1:f3),k.) plot(x(f3+1:f4),y(f3+1:f4),k+) plot(x(f4+1:f5),y(f4+1:f5),ks) legend(生活区,工业区,山区 ,交通区,公园绿地区) hold off X,Y=meshgrid(0:100:30000,1:100:20000); figure(2) ZAs=griddata(x,y,As,X,Y); surf(X,Y,ZAs,50) contour(X,Y,ZAs,30 20 15 10 8 6 5 4 3 2 1) title(contour of As) figure(3) ZCd=griddata(x,y,Cd,X,Y); surf(X,Y,ZCd) contour(X,Y,ZCd,50) title(contour of Cd) figure(4) ZCr=griddata(x,y,Cr,X,Y); surf(X,Y,ZCr) contour(X,Y,ZCr,50) title(contour of Cr) figure(5) ZCu=griddata(x,y,Cu,X,Y); surf(X,Y,ZCu) contour(X,Y,ZCu,100) title(contour of Cu) figure(6) ZHg=griddata(x,y,Hg,X,Y); surf(X,Y,ZHg) contour(X,Y,ZHg,50) title(contour of Hg) figure(7) ZNi=griddata(x,y,Ni,X,Y); surf(X,Y,ZNi) contour(X,Y,ZNi,50) title(contour of Ni) figure(8) ZPb=griddata(x,y,Pb,X,Y); surf(X,Y,ZPb) contour(X,Y,ZPb,50) title(contour of Pb) figure(9) ZZn=griddata(x,y,Zn,X,Y); surf(X,Y,ZZn) contour(X,Y,ZZn,50) title(contour of Zn) 程序 2：相关性分析h=data(:,4);y1=data(:,6);y2=data(:,7);y3=data(:,8);y4=data(:,9);y5=data(:,10);y6=data(:,11);y7=data(:,12);y7=data(:,12);y8=data(:,13);p1=corrcoef(h,y1)p2=corrcoef(h,y2)p3=corrcoef(h,y3)p4=corrcoef(h,y4)p5=corrcoef(h,y5)p6=corrcoef(h,y6)p7=corrcoef(h,y7)p8=corrcoef(h,y8)程序 3：标准差a1=std(data(:,6)a2=std(data(:,7)a3=std(data(:,8)a4=std(data(:,9)a5=std(data(:,10)a6=std(data(:,11)a7=std(data(:,12)a8=std(data(:,13)程序 4：解方程syms x y x0 y0 bz=exp(x-x0)2+(y-y0)2)/b)diff(z,x0)diff(z,y0)x=74;1373;1321;0;1049;1647;2883;2383;2708;2933;4233;4043;2427;3526;5062;4777;5868;6534;5481;4592;2486;3299;3573;4741;5375;5635;5394;5291;4742;4948;5567;7004;7304;7048;8180;9328;9090;8049;8077;8017;6869;7056;7747;8457;9460;9062;9319;10631;10685;10643;11702;11730;11482;10700;10630;11678;11902;13244;12746;12855;13797;14325;15467;12442;13093;13920;14844;16569;16387;16061;15658;14298;14177;15092;12778;17044;17087;17075;17962;18413;19007;18738;17814;18134;17198;17144;18393;19767;21006;21091;22846;23664;22304;21418;21439;20554;20101;21072;20215;18993;19968;21766;22674;22535;25221;26453;2