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

1、2011高教社杯全国大学生数学建模竞赛承诺书我们仔细阅读了中国大学生数学建模竞赛的竞赛规则我们完全明白，在竞赛开始后参赛队员不能以任何方式（包括电话、电子邮 件、网上咨询等）与队外的任何人（包括指导教师）研究、讨论与赛题有关的问 题。我们知道，抄袭别人的成果是违反竞赛规则的，如果引用别人的成果或其他公 开的资料（包括网上查到的资料），必须按照规定的参考文献的表述方式在正文引 用处和参考文献中明确列出。我们郑重承诺，严格遵守竞赛规则，以保证竞赛的公正、公平性。如有违反竞赛规则的行为，我们将受到严肃处理。我们参赛选择的题号是（从 A/B/C/D中选择一项填写）：我们的参赛报名号为（如果赛区设置报名

2、号的话） ： 所属学校（请填写完整的全名）：参赛队员（打印并签名）：1.2. 3. 指导教师或指导教师组负责人（打印并签名）：日期：年月日赛区评阅编号（由赛区组委会评阅前进行编号）：编号专用页赛区评阅编号（由赛区组委会评阅前进行编号）：赛区评阅记录（可供赛区评阅时使用）：评阅人评分备注全国统一编号（由赛区组委会送交全国前编号）：全国评阅编号（由全国组委会评阅前进行编号）：城市表层土壤重金属污染分析摘要随着城市经济的快速发展和城市人口的不断增加，人类活动对城市环境质量 的影响日益突出。通过土壤地质的查证开展城市环境质量评价日益成为人们关注 的焦点。针对问题一，首先应用 Excel 软件筛选所给采

3、样点数据，将其整合到五个区 域中，利用 Matlab 软件绘出采样点的散点分布图， 然后应用 Matlab 软件中插值拟 合函数得到重金属元素在该城区的空间分布图。最后采用单因子指数法和内梅罗 综合指数法对土壤重金属污染进行评价，得到生活区中Cu、Zn污染相对严重，综合污染等级为警戒级。工业区中 Hg、Cu、Pb、Zn、Cd污染相对严重，综合 污染等级为严重污染。山区中综合污染等级为清洁级。主干道路区只有Hg为强度 污染，综合污染等级为轻度污染。公园绿地区大部分元素为轻度污染，综合污染 等级为清洁级。针对问题二，将重金属污染的主要原因进行分区域讨论，运用主成分分析法 确定各区域的主要污染原因，

4、根据其结果得到生活区受污染的主要原因是公路两 侧含铅汽油的燃烧和汽车轮胎的磨损，因此 Zn 、 Cr 的含量明显高于国家背景值； 工业区“三废”的排放和金属的冶炼导致 As、Cd、Cr、Ni、Pb、Zn元素的含量 偏高；山区主要由于金属矿山的开采导致 Cu元素含量偏高；主干道路区污染的主 要原因是由机动车尾气排放使得 Pb、Cd元素的含量偏高造成的，农药化肥的使 用和地下水的灌溉，致使Pb、Zn元素浓度过量。针对问题三，本文根据重金属污染物在土壤中的传播特征，依据污染物离子 迁移方程的理论建立重金属污染物浓度与迁移距离的微分方程，从而得到重金属 污染物的传播规律，并确定了污染源的位置。针对问题

5、四，本文重点讨论了土壤的酸碱度、组分、温度等因素对重金属污 染物的扩散和吸附所产生的重要影响，因此需要采集这些信息进行相关性分析， 并据此对模型进行改进以便更好地研究城市地质环境的演变模式。关键词： 单因子指数法 内梅罗综合指数法 主成分分析法 微分方程一 问题重述随着城市经济的快速发展和城市人口的不断增加，人类活动对城市环境质量 的影响日显突出。对城市土壤地质环境异常的查证，以及如何应用查证获得的海 量数据资料开展城市环境质量评价，研究人类活动影响下城市地质环境的演变模 式，日益成为人们关注的焦点。按照功能划分，城区一般可分为生活区、工业区、山区、主干道路区及公园 绿地区等，分别记为 1 类

6、区、2 类区、, 、 5 类区，不同的区域环境受人类活动 影响的程度不同。现对某城市城区土壤地质环境进行调查。为此，将所考察的城区划分为间距 1 公里左右的网格子区域， 按照每平方公里 1 个采样点对表层土 (010 厘米深度) 进行取样、编号，并用 GPS 记录采样点的位置。应用专门仪器测试分析，获得了 每个样本所含的多种化学元素的浓度数据。另一方面，按照 2 公里的间距在那些 远离人群及工业活动的自然区取样，将其作为该城区表层土壤中元素的背景值。附件 1 列出了采样点的位置、海拔高度及其所属功能区等信息，附件 2列出 了 8 种主要重金属元素在采样点处的浓度，附件 3 列出了 8 种主要重

7、金属元素的 背景值。现要求你们通过数学建模来完成以下任务：(1) 给出 8 种主要重金属元素在该城区的空间分布，并分析该城区内不同区 域重金属的污染程度。(2) 通过数据分析，说明重金属污染的主要原因。(3) 分析重金属污染物的传播特征，由此建立模型，确定污染源的位置。(4) 分析你所建立模型的优缺点，为更好地研究城市地质环境的演变模式， 还应收集什么信息？有了这些信息，如何建立模型解决问题？二 模型假设1. 问题一中国家评价标准采用二级标准。2. 假设全部采样点在自然情况下产生，排除其他因素的干扰。3. 假设该城区没有受到灾难性的污染。4. 假设污染源是均匀向外扩散的。三 符号说明符号名称

8、符号意义iIiCi第i种重金属元素污染元素 i 的分指数某区域 i 元素的平均浓度Si国家评价标准P综合污染指数c1重金属污染源以扩散形式传播的离子浓度Rf为污染源离子浓度延迟因数D为污染源离子浓度扩散系数c2重金属污染源以连续分布吸附传播的离子浓度c污染源处重金属污染物最大浓度aibi(t)土壤酸碱度、组分、温度、等因素对重金属污染物 传播的影响系数为土壤酸碱度、组分、温度、等因素随时间变化的 函数四 问题分析对于问题一已知城区由五个不同的区域组成，想要得到该城区重金属元素的 空间分布图，首先需要得到重金属在五个区域浓度和地势分布图，由于所给采样 点不是位于同一个区域，这就需要对所给采样点的

9、数据进行分区域筛选和整合， 然后将整合后的各个区域采样点的浓度和三维坐标值导入 Matlab 软件绘出八种 重金属在各个区域的浓度和地势分布图，即得到重金属元素在该城区的空间分布 图；由于还需要分析不同区域重金属的污染程度， 污染程度需要抽象为污染指数， 通过比较污染指数与国家背景值，确定污染程度较大的金属元素，针对这个问题 我们采用单因子指数法和内梅罗综合指数法，得到不同区域重金属的污染程度。对于问题二经过数据的分析，重金属污染的主要原因取决于重金属在各个区 域中的含量，由于每个区域的重金属含量不同，我们首先需要确定哪些重金属元 素为污染环境的主要因素，采用主成分分析法确定主要污染环境的元素

10、。对于问题三，重金属污染物在土壤中以离子形式存在，其迁移形式一般为对 流、弥散、扩散、吸附等形式。研究重金属污染物的传播特征即可以转化为研究 重金属污染物在土壤中离子形式迁移的运动规律，由运动规律分析污染物浓度在 整个考察的土壤区域的分布，进而可以确定其污染源的位置。对于问题四，土壤酸碱度、土壤组分、动植物的影响、人为污染和破坏、气 候变化等等因素都会使得城市地质环境发生重要的变化。研究这些因素对地质环 境的演变影响，需要了解相关方面的研究理论、采集相关的统计数据和实验数据 加以分析，结合上述的研究方法可以得到城市地质环境演变的规律。五模型的建立与求解5.1 问题一本文给出了一个城区的五个不同

11、的环境区域，生活区、工业区、山区、主干道路区及公园绿地区和大量采样点的浓度值和三维坐标，由于采样点分布的离散程 度大，所以需要利用Excel软件对大量的数据进行筛选和整合，得到五个区域采 样点的散点分布图（见图1）4003002001001.54X 100.54X 10图1筛选后采样点的散点分布图为了得到更加形象具体的重金属空间分布图，我们利用As、Cd、Cr、Cu、Hg、Ni、Pb、Zn八种重金属的浓度分布图和生活区、工业区、山区、主干道 路区及公园绿地区五个区域的地势分布三维图对其进行说明。首先通过引用Matlab软件中插值拟合函数得到该城区八种重金属的浓度分布图（见图 2,图3, 图4,

12、图5,图6,图7,图8,图9）3018000+1600025-+1400012000201000015800060001040002000000.5122.5x 10图2 As浓度分布图1800030160002514000120002010000158000600010400020000.51.52.5F S人 F .相落+$耳、斗.汁 _|_a t .g曲” i-知+ /+4 七V护 4-汽 2 47.1 护亠 七+工*-严7工业区 交通区 生活区 公园绿地区 山区 t*化咯i占I X ,1.50x 10图3 Cd浓度分布图3000.5工业区 交通区 生活区 公园绿地区 山区4252015

13、105x 10图4 Cr浓度分布图0.51.52.530252015104图5 Cu浓度分布图1800016000140001200010000800060004000200001800016000140001200010000800060004000200005x 1025201510图6 Hg浓度分布图工业区 交通区 生活区 公园绿地区 山区302520151018000160001400012000100008000600040002000018000160001400012000100008000600040002000000.511.522.54图7 Ni浓度分布图30工业区 交通区

14、 生活区 公园绿地区 山区0.51.54252015105x 10180001600014000120001000080006000400020000180001600014000120001000080006000400020000图8 Pb浓度分布图0.5工业区 交通区 生活区 公园绿地区 山区0。304 +2520+ 15101.52.5x 10n 44*4：-+ + 4-4- 4.4-+ + #+十4斗”，十+ K气4*图9 Zn浓度分布图接下来用Matlab软件中插值拟合函数得到该城五个区域的三维地势分布图（见图10）400300Ai200100-10021.54X 100.54 X

15、 10图10五个区域的地势分布图 单因子指数法接下来我们对该城区不同区域的污染程度进行评价，由于每个区域都存在八 种重金属，所以我们对每个区域的八种重金属含量分别进行评价，评价模式分为 单因子指数法和内梅罗综合指数法。首先建立单项污染指数式为：ii=C(1)式中：i取18分别对应 As、Cd、Cr、Cu、Hg、Ni、Pb、Zn八种元素，.为污染元素i的分指数，C为某区域中i元素的平均浓度，S为i元素的国家评价 标准。内梅罗综合指数法单项污染指数得到每个区域各种重金属污染程度，为了得到每个区域的综合 污染程度，我们利用内梅罗综合指数法求解。内梅罗综合指数式为：J卩.5“ - （平均值）+（I最大

16、值）(2)式中：P为综合污染指数，I平均值为各元素污染指数平均值，I最大值为各元素污染 指数中的最大值。表1各个城区8种主要重金属元素的平均值城市功能 区As平均值（卩 g/g）Cd平均 值（卩g/g）Cr平均值（卩 g/g）Cu平均值（卩 g/g）Hg平均 值（卩g/g）Ni平均值（卩 g/g）Pb平均值（卩 g/g）Zn平均值（卩 g/g）生活区6.270.2969.0249.400.9318.3469.11237.01工业区7.250.3953.41127.54P 0.6419.8193.04277.93山区4.040.1538.9617.320.0415.4536.5673.29交通区

17、5.710.3658.0562.21r 0.4517.62 :63.53 1242.85公园区6.260.2843.6430.190.1115.2960.71154.24表2 各个城区8种主要重金属元素的土壤单项污染指数及综合指数城市功 能区土壤单项污染指数土壤综合污 染指数As（卩 g/g）Cd（卩 g/g）Cr（卩 g/g）Cu（卩 g/g）Hg（卩 g/g）Ni（卩 g/g）Pb（卩 g/g）Zn（卩 g/g）生活区1.742.232.233.7426.571.492.233.4319.18工业区2.013.021.729.6618.291.613.004.0313.49山区1.121.

18、171.26M.311.141.261.181.061.25交通区1.592.771.874.7112.771.432.053.529.43公园区1.742.161.412.293.141.241.962.242.64表3污染等级划分标准指数范围污染等级指数范围污染等级3.0极严重污染由单因子污染指数分析可知，该城市各功能区表层土壤中As、Cd、Cr、Cu、Hg、Ni、Pb、Zn的单因子污染指数均大于1,说明各功能区表层土壤均存在一定 程度的 As、Cd、Cr、Cu、Hg、Ni、Pb、Zn 污染。As、Cd、Cr、Cu、Hg、Ni、 Pb、Zn在生活区、工业区、山区、交通区和公园绿地区的单因子

19、污染指数分别为1.12 2.01、1.15 3.00、1.29 2.33、1.31 9.66、1.14 26.57、1.24 1.61、1.18 3.02、1.064.03, Hg的污染指数最高，其次为Cu, As的污染指数最低。在该城 市的五个区域中Hg污染最严重，在山区属于轻度污染，其余四个区域都属于重 度污染；其次为Cu，在山区属于轻度污染，在公园绿地区属于中度污染，其余三 区都属于重度污染；As的污染较轻，除了工业区为中度污染外，其余各城区均属于轻度污染。由综合污染指数分析可知，该城市的生活区、工业区、交通区表层土壤As、Cd、Cr、Cu、Hg、Ni、Pb、Zn的综合污染指数均大于3,

20、综合污染指数分别为 19.18、13.49、9.43,均属于重度污染；公园绿地区的表层土壤重金属综合污染指 数为2.64，属于中度污染；山区的表层土壤重金属综合污染指数为1.25,属于轻度污染。由上述分析可见，该城市的生活区、工业区表层土壤重金属污染最为严 重，其次为交通区，山区污染较轻。该城市各功能区表层土壤重金属的污染状况 基本上反映了该市土地利用格局，与实际情况相符合。5.2问题二主成分分析法主成分分析也称主分量分析，旨在利用降维的思想，把多指标转化为少数几 个综合指标。在实证问题研究中，为了全面、系统地分析问题，我们必须考虑众 多影响因素。这些涉及的因素一般称为指标，在多元统计分析中也

21、称为变量。因 为每个变量都在不同程度上反映了所研究问题的某些信息，并且指标之间彼此有 一定的相关性，因而所得的统计数据反映的信息在一定程度上有重叠。在用统计 方法研究多变量问题时，变量太多会增加计算量和增加分析问题的复杂性，人们 希望在进行定量分析的过程中，涉及的变量较少，得到的信息量较多。首先我们通过将生活区的重金属元素原始数据导入 SPSS软件进行主成分分 析处理，得到生活区各种重金属元素的主成分荷载矩阵（见表 4）及方差分析主 成分分析（见表5）。表4生活区重金属元素的成分矩阵成份矩阵成份123As.669-.646-.010Cd.784.171-.417Cr.643.234.493Cu

22、.729-.246.024Hg.492.130-.437Ni.686-.253.523Pb.803.112-.348Zn.501.691.267提取方法：主成份。a.已提取了 3个成份表5生活区重金属元素的特征值解释的总方差成份初始特征值提取平方和载入合计方差的%累积%合计方差的%累积%13.61645.19945.1993.61645.19945.19921.13314.16559.3651.13314.16559.36531.07513.43272.7971.07513.43272.7974.80710.08382.8805.5246.54489.4246.4455.55694.9807.

23、2362.95497.9348.1652.066100.000提取方法：主成份分析。表6工业区重金属元素的成分矩阵成份矩阵a成份12As.518.758Cd.786.074Cr.916-.206Cu.868-.463Hg.845-.459Ni.767.421Pb.858.049Zn.859.188提取方法：主成份。a.已提取了 2个成份。表7山区重金属元素的成分矩阵成份矩阵成份123As-.009.668.478Cd.601-.678.141Cr.761.480-.365Cu.517.469.615Hg.324-.074.747Ni.737.499-.409Pb.605-.630.173Zn.

24、905-.173-.183提取方法：主成份。a.已提取了 3个成份。表8主干道路区重金属元素的成分矩阵成份矩阵成份12As.235-.131Cd.621.458Cr.874-.319Cu.906-.250Hg.170.743Ni.888-.320Pb.703.477Zn.643.115提取方法：主成份。a.已提取了 2个成份。表9公园绿地区重金属元素的成分矩阵成份矩阵成份123As.636-.569.349Cd.811.103-.232Cr.809-.426-.003Cu.679.457-.151Hg.203.354.889Ni.663-.622.012Pb.782.544.097Zn.798

25、.286-.242提取方法：主成份。a.已提取了 3个成份主成分个数提取原则为主成分对应的特征值大于 1的前m个主成分，通过表 21,可知提取3个主成分。用表16中数据除以主成分对应特征值的开平方便得到三 个成分中每个指标对应的系数，即特征向量 A1，A2，A3，将特征向量与标准化 后数据项城变得主成分表达式：F1=A1*ZX1;F2=A2*ZX2;F3=A3*ZX3。以每个主成分所对应的特征值占所提取主成分总的特征值之和的比例作为权 重计算主成分综合模型，即用第一主成分 F1中每个指标所对应的系数乘上第一主 成分F1所对应的贡献率再除以所提取的三个主成分贡献率之和，然后加上第二主 成分F2对

26、应系数乘以其贡献率百分比，再加上第三主成分F3对应系数乘以其贡献率百分比，就得到了综合主成分F中每个指标所对应的权重系数A。根据权重系数A值大小降序排列，系数大的对应的重金属元素就是主要污染的主要原因。这里 只将生活区列举出来，其它区域分析方法与上述方法相同，就不做过多叙述，其 他四个区域的各种重金属元素的主成分荷载矩阵（见表 6表9）。表10各个区域主要污染元素分布生活区权重系 数工业区权重系 数山区权重系 数Cr0.34Ni0.34Cu0.35Zn0.34Zn0.33As0.23Ni0.27As0.31:Cr0.23 :Pb0.22Pb0.31Ni0.22Cd0.21Cd0.29Hg0.2

27、1Cu0.2Cr0.29Zn0.17Hg0.11Cu0.23Pb0.06As0.1Hg0.22Cd0.04 J主干道区权重系 数公园绿 地权重系 数Pb 0.38Pb0.36Cd0.34Cu0.27Cu0.29Hg0.27Ni0.27Zn0.26Zn0.27Cd0.23Cr0.26Cr0.16Hg0.23As0.14As0.06Ni0.08结果评价对该城区不同区域采集的土壤 As、Cd、Cr、Cu、Hg、Ni、Pb、Zn元 素进行主成分分析，主成分分析结果（见表10）。通过主成分分析得到：（1）生活区元素Zn、Cr在该城市土壤中的含量较高。这是由于生活区存在 大量的公路，公路两侧含铅汽油的燃烧

28、和汽车轮胎的磨损的粉尘会增加土壤中Zn、Cr等元素的含量，而且由于生活区人类活动比较频繁，日常生活用品丢弃 后成为垃圾也会导致这三种元素的增加。（2）工业区主要污染元素为 As、Cd、Cr、Ni、Pb、Zn。因为工业区“三 废”排放，采矿和冶炼会增加这些元素在土壤中的含量。（3）山区只有Cu元素为主要污染元素。这是由于山区远离城市和工业区的污 染源，只有一部分金属矿山的开采会导致一些金属元素的外露。（4）主干道路区Pb、Cd元素较高。机动车尾气排放既是城市大气的主要污染源，也显著引起公路两侧土壤的重金属污染，汽车汽油、发动机、轮胎、润滑 油和镀金部分都能燃烧或磨损而释放出Hg、Zn元素。（5）

29、公园绿地区中Pb为主要重金属元素，Hg、Cu、Zn为次要重金属元素。公园绿地区农药化肥的使用会不同程度的影响重金属的污染，而且地下水灌溉和 塑料薄膜会增加、Cu、Zn元素的土壤含量。5.3问题三土壤中重金属的传播一般为对流、弥散、扩散、吸附四种形式。由于重金属 的性质决定了对流和弥散在重金属传播过程中只占有很少的比例，在此情况下， 假设土壤性质为饱和粘土且不考虑重金属污染物的蜕变，故重金属污染物的传播 主要是两种扩散和吸附形式，并且这两种形式对重金属污染物的传播影响可以取 叠加和的关系。首先，对重金属污染物扩散形式建立模型得出其重金属污染物浓度，再对重 金属吸附形式建立模型求得重金属污染浓度，

30、最后将两种形式相结合可得重金属 污染物浓度与污染源的模型。5.3.1重金属污染物扩散形式模型的建立由生物和化学知识可知，土壤中重金属污染物以离子形式存在和传播。现假 设重金属污染物离子的传播只进行扩散形式而不考虑吸附形式，则污染源离子浓 度应该是随时间t和传播距离变化的函数r，可以表示为d =d(t,r)。：t通过了解相关土壤研究的知识和资料查询可以获得t和r的具体关系式:(3)其中Rf为延迟因数，D为扩散系数，假定污染物离子浓度为一个与时间变化f(t)无关的函数，故可用求解微分方程的方法求解3式，解为：1 Rf 2ci(t,r)f(t) -r C( 4)2 D上式中，C为常数，f(t)为一个

31、常函数，该式即为重金属扩散形式模型。5.3.2重金属污染物吸附形式模型的建立重金属进入水体或土壤之后大多都结合在天然有机物的表面，天然有机物包 括腐植酸、富里酸、沉积物、悬浮物等多种组成和和形式，其吸附、解析、沉积 和迁移决定了重金属的去向和归宿。目前，重金属吸附普遍采用的模型为 Langmuir和Freundlich等模型，这些 模型为经验模型一般不能确切的表达吸附机理，并且这些模型中的结合常数还要依赖外部条件，当土壤酸碱性、离子浓度发生变化后，模型不能很好地反映重金 属浓度变化规律。于是，在此背景下英法等国学者提出了NICA-D onnan模型。我们利用此模型的建立原理为依据，建立重金属污

32、染物吸附模型。I局部吸附模型将土壤当作吸附剂，由于土壤对不同的污染物吸附能力不同，即对不同的污染物有不同的吸附常数。假设土壤对第i种确定的污染物其吸附常数是确定的，记为Ki，于是单位土壤面积上吸附物质的浓度计算式为G = KiCin，式中n表示温度对土壤吸附能力影响参数，假设温度在局部保持恒定，故在此处设定n为1。n连续分布吸附模型在建立了局部的吸附模型后，需要将局部模型扩大到整个区域的吸附模型。 由于所得局部吸附模型是在离散情形下建立的，故利用微积分学中的微元法和求 定积分的方法得到如下模型：c2 = Gd；丁二 KjCidn(5)在式5中Ki、Ci随着土壤性质的不同而不同，其中 G变化趋势

33、应该是第i种重金属污染物元素浓度随着离污染源距离增大而减小。故应有当G二g(r)时g(r) _ 0，将此式代入上面式5中得出：Q * M g(r)d二(6)川重金属污染物浓度与污染源的模型求出重金属污染物在土壤中的扩散模型和吸附模型之后，联立式4和式6得到重金属污染物浓度与污染源的模型，即1 Rf 2c=C!(t,r) C2f (t) - rKj!g(r)d二 C(7)D其中，c = Cmax，即污染源处重金属污染物浓度最大。由最小二乘法原理利用Matlab编程对数据进行二次曲性拟合，拟合曲线如图 11所示50040030020010080060040020012015010040010008

34、00600400200300025002000150010005000-50-100-150-20060050040030020012008006004002001000100不同区域之间重金属污染物的污染源是不同的，即每个区有各自的污染源， 同一个区域内某几个元素变化趋势相似， 而污染物的浓度却有明显差异。基于此 特点，结合题目所给数据及所建立模型，确定该城市污染源的位置如下表所示：表11污染源的位置标号污染源坐标污染源包含的重金属污染源靠近的功能区1(2589,10730)Pb ,Cd工业区，交通区2(2226,3422)Ni ,Cu,Cr,Hg,Pb,Cd工业区，交通区3(4812,7213)As工业区4(18240,10080)As交通区5(21620,11420)As交通区6(13720,2260)Hg生活区7(15300,9294)Hg交通区8(13720,9618)Zn工业交通5.4问题四在5.3中我们综合考虑了影