2011数学建模国家二等奖

2011高教社杯全国大学生数学建模竞赛 承 诺 书 我们仔细阅读了中国大学生数学建模竞赛的竞赛规则. 我们完全明白，在竞赛开始后参赛队员不能以任何方式（包括电话、电子邮件、网上咨询等）与队外的任何人（包括指导教师）研究、讨论与赛题有关的问题。 我们知道，抄袭别人的成果是违反竞赛规则的, 如果引用别人的成果或其他公开的资料（包括网上查到的资料），必须按照规定的参考文献的 表 关于同志近三年现实表现材料材料类招标技术评分表图表与交易pdf视力表打印pdf用图表说话 pdf 述方式在正文引用处和参考文献中明确列出。 我们郑重承诺，严格遵守竞赛规则，以保证竞赛的公正、公平性。如有违反竞赛规则的行为，我们将受到严肃处理。 我们参赛选择的题号是（从A/B/C/D中选择一项填写）： A 我们的参赛报名号为（如果赛区设置报名号的话）： A 甲01720 所属学校（请填写完整的全名）： 中国石油大学（华东） 参赛队员 (打印并签名) ：1. 欧阳剑锋 2. 马洋洋 3. 李顺丽 指导教师或指导教师组负责人 (打印并签名)： 乔田田 日期： 2011 年 9 月 12 日 赛区评阅编号（由赛区组委会评阅前进行编号）： 2011高教社杯全国大学生数学建模竞赛 编 号 专 用 页 赛区评阅编号（由赛区组委会评阅前进行编号）： 赛区评阅记录（可供赛区评阅时使用）： 评 阅 人 评 分 备 注 全国统一编号（由赛区组委会送交全国前编号）： 全国评阅编号（由全国组委会评阅前进行编号）： 城市表层土壤重金属污染分析 摘要 本文研究了城区不同区域污染物的空间分布及污染物的传播特征问题，同时，分析了不同区域的污染程度、土壤金属污染的主要原因和城市地质环境的演变模式。 文中采用Matlab绘制重金属元素在城区的空间分布图，通过对图形的整体观察，可知城区海拔右高左低，污染物分布山区少平原多，各功能区内重金属元素呈局部积聚小范围传播的特点。结合平均浓度值和相应背景值分析计算地积累指数又称Muller指数和土壤的富集系数，并分别与土壤环境质量 标准 excel标准偏差excel标准偏差函数exl标准差函数国标检验抽样标准表免费下载红头文件格式标准下载 对照来表征不同区域重金属的污染程度，这样可以从自然因素中分离出人为因素，更直观地表征区域污染程度。 利用各区8种元素相关性、主成分分析法、各区的污染程度分析重金属污染的主要原因并结合物理、化学、生物、社会生活等因素分析出该地区主要的污染元素和污染原因，最终得到1区即生活区Ni属于无污染，As、Cd、Hg、Pb和Cr属于轻度-中等污染，Cu和Zn属于中等污染，主要污染元素是Zn、Cu、Hg和Pb，污染原因主要是城市生活垃圾和交通污染；2区即工业区As，Cr，Ni属于轻度-中等污染，Cd、Pb和Zn属于中等污染，Cu属于中等-强污染，Hg属于强污染，主要污染元素是Cu、Zn、Pb和Cd，污染原因在于工业“三废”、农业污染和土壤本身来源；3区即山区污染很小，因为山区本身人类活动少；4区即主干道区Ni属于无污染，As，Cd，Cr，Pb属于轻度-中等污染，Cu和Zn属于中等污染，Hg属于强污染，主要污染元素是Hg、Cu、Zn和Cd，污染原因是工业、农业污染和土壤本身含量；5区即公园绿地区Cr、Ni属于无污染，As、Cd、Cu、Pb、Zn属于轻度-中等污染，Hg属于中等污染，主要污染元素是Hg，污染原因在于工业污染。 根据文献资料并结合第一问中重金属元素在城区的空间分布图，了解到重金属污染物的传播特征可认为是以污染源为中心的径向传播。同时，污染物浓度随径向距离成指数递减，且时间为定常量即浓度分布与时间无关。选取分布相对规则、具有代表性的Hg元素分析，根据重金属离子在土壤中扩散的规律，并将已知数据与经插值得到的数据共同拟合出污染物扩散的稳态方程： 该问题是一个反问题，因此，对于反复尝试拟合出的模型要进行回判。本文对于该数学模型的准确性，利用回归分析计算判定系数 等于0.877，说明该方程是收敛、稳态的，可以作为土壤中重金属污染物的传播特征方程。利用径向距离r与x,y的关系，使用Matlab工具箱将三维方程转化为二维进行拟合，用拟合出的多项式求解极值，得出污染源位置： 。考虑到实际取点距离过大、重金属元素实际分布较集中、浓度衰减迅速等因素，将已建方程指数一次项修正为二次项，再次拟合，得到与实际点拟合程度非常好的拟合曲线。 从模型结构可以看出局部小范围内污染物的传播特征，而且方程形式简单易于分析。同时，也存在一定的局限性。首先，模型不是完全地反映整个城市的污染物浓度分布；其次，模型忽略了海拔、水流、风等自然因素的影响，并不具备一般性。由于模型本身的不足，需再收集有关城市的地形地貌、土壤成分、气候特征、主导风向等信息，并分析这些因素与城市地质环境的关系，结合已建的传播模型，预测城市地质环境的演变模式。 关键词：重金属污染分布 Muller指数 相关性 主成分分析 插值拟合 回归分析 1、 问题重述 由于经济增长、人口增加，城区环境质量成为人们关注的焦点。因此，要学会城市土壤地质环境异常的查证，同时利用已得数据进行环境质量评价，进而研究城市地质环境的演变模式。 城区按功能一般可分为生活区、工业区、山区、主干道路区及公园绿地区等，分别记为1类区、2类区、……、5类区，不同的区域环境受人类活动影响的程度不同。 为对某城市城区土壤地质环境进行调查，将所考察的城区划分为间距1公里左右的网格子区域，按照每平方公里1个采样点对表层土（0~10 厘米深度）进行取样、编号，并用GPS记录采样点的位置。应用专门仪器测试分析，获得了每个样本所含的多种化学元素的浓度数据。另一方面，按照2公里的间距在那些远离人群及工业活动的自然区取样，将其作为该城区表层土壤中元素的背景值。 采样点的位置、海拔及所属功能区8种主要元素在采样点处的浓度已知，同时8种主要重金属元素背景值已测得。 要求通过数学建模来完成以下任务： (1) 给出8种主要重金属元素在该城区的空间分布，并分析该城区内不同区域重金属的污染程度。 (2) 通过数据分析，说明重金属污染的主要原因。 (3) 分析重金属污染物的传播特征，由此建立模型，确定污染源的位置。 (4) 对所建模型进行评价，分析优缺点，为更好地研究城市地质环境的演变模式，还应收集的信息，及利用信息建立模型。 二、问题分析 本题是一个城市表层土壤重金属污染分析问题，就是需要建立重金属污染物扩散模型得出城市污染物的传播特征方程。而模型的建立源于对城区重金属污染物的分布和不同元素对不同功能区污染程度的分析，污染物由污染源向各个方向扩散，由此可以判断污染物浓度与污染源的距离有关。 第一问对城市土壤地质环境调查得到的数据，用Excel按功能区重新统计，用Matlab编程画出X、Y坐标与一种金属浓度及X、Y坐标与海拔的关系图，综合两种分布图可分别得到8种重金属元素在该城区不同功能区的分布。用算术平均值求得不同区域各种浓度的平均值，结合背景值，计算地积累指数并与土壤环境质量标准对照分析污染程度。 第二问中将不同区域8种污染物浓度进行相关性分析、主成分分析，同时结合第一问中的不同区重金属污染程度结果，综合考虑地理环境、化学、生物等相关知识，可以说明重金属污染的主要原因。 第三问中，分析第一问中城区重金属空间分布图和假设条件，可知重金属浓度在离污染源越近的地方越大，越远越小，且衰减很快。根据一维扩散方程，将已知数据和Matlab插值所得数据统计分析，则得到重金属浓度与据污染源距离成指数关系。再用Matlab将数据拟合可得重金属污染物传播模型，又根据假设条件各种金属传播规律相近，由此可以根据某一功能区一种重金属浓度得到重金属污染物的传播方程。再将实际浓度值与方程得到浓度值进行方差计算，则可得所建模型的收敛性、稳定性，若模型收敛，则建立模型正确；反之，模型建立不恰当。对于污染源位置确定问题，利用径向距离r与x,y的关系，使用Matlab工具箱将三维方程转化为二维进行拟合，对拟合出的多项式求解极值，得出污染源位置： 。 第四问中评价模型优缺点，由于在建模型时，仅考虑了主要因素（重金属污染物在土壤中依靠水土流失自然扩散），且为稳定状态。但实际上污染物传播还与重金属的种类、位置、地理条件等有关，那么模型在一定程度上可以说明重金属污染物的传播，但存在局限性，需要修正。因此，针对实际取点距离过大、重金属元素实际分布较集中、浓度衰减迅速等因素，将已建方程指数一次项修正为二次项，再次拟合。 针对该城市重金属污染问题说明该城市地质环境变化，但不足以体现地质环境的演化，因此要收集该城市的地形地貌、河流分布、主导风向以及环境气候等信息，利用这些信息，分析这些信息与土壤中重金属污染物浓度的关系得到新的数学模型来表征城市地质环境的演化过程。 三、模型假设 （1）土壤中重金属污染物浓度不随时间变化，处于稳态； （2）土壤中重金属元素是沿径向传播、圆形扩散的； （3）土壤中重金属污染物的传播途径在径向满足扩散方程； （4）土壤中污染物由载体带动运动，且载体在一定时间范围内可视为以恒定的速度运动； （5）同种元素在不同地域的传播特性相同，即某种元素局部传播特性可以反映整体规律； （6）海拔高度对重金属污染物的传播影响可以忽略； （7）重金属污染物浓度最大的位置为污染源。 四、符号说明 c：土壤中重金属污染物浓度值； r：被检测点到污染源的距离； t：时间； x：被测点的x轴坐标值； y：被测点的y轴坐标值； D：扩散系数； ν：为水的通量密度； s：吸附在土壤上重金属离子的浓度； ρ：土壤介质的密度； θ：空隙率。 五、模型的建立与求解 5.1重金属元素空间分布及污染程度 5.1.1重金属元素空间分布 As (μg/g) Cd (ng/g) Cr (μg/g) Cu (μg/g) Hg (ng/g) Ni (μg/g) Pb (μg/g) Zn (μg/g) (a) (c) (d) (e) (f) (g) (h) (i) 图一 重金属元素空间分布图 5.1.2各区域重金属元素污染程度 对于重金属评价方法，目前主要有综合评价法、污染指数评价法、聚类法、因子分析法。这些方法能对测试区污染程度进行较为全面的评价，但无法从自然异常中分离人为异常，获得人为因素所造成的影响，这里采取地积累指数法和富集系数法对土壤的污 染进行评价，这两种方法均引入背景值作为标准对土壤重金属含量进行归一化处理，表示过程中重金属元素的人为污染情况。其结果不仅有助于了解重金属分布的自然变化特征，而且可以判别人类活动对环境的影响，使土壤中重金属的污染得到更好的评价。 地积累指数法是德国海德堡大学沉积物研究所科学家Muller在1969年提出的，它是研究水环境沉积物中重金属污染程度的定量指标，在欧洲被广泛采用，近年来被国内外的学者专家广泛用于人为活动产生的重金属对土壤污染的评价。其计算公式如下： （5.1.2-1）； （5.1.2-1）式中， 表示某一样品中元素n的浓度； 表示一粘质沉积岩(普通页岩)中该元素的地球化学背景值。常量1．5是为消除各地岩石差异可能引起背景值的变动转换系数。 (1) 沉积物重金属地积累指数分级与污染程度之间的相互关系列于表一： 表一：Muller地积累指数分级标准 地积累指数Igeo 分级 污染程度 (5,10] 6 极严重污染 (4,5] 5 强-极严重污染 (3,4] 4 强污染 (2,3] 3 中等-强污染 (1,2] 2 中等污染 (0,1] 1 轻度-中等污染 (-∞,0] 0 无污染 (2)土壤环境质量标准是土壤中污染物的最高容许含量。根据土壤应用功能和保护目标，划分为三类。城市土壤重金属评价通常采用中华人民共和国《土壤环境质量标准》，见表二： 表二：土壤环境质量标准（二级） 项目 PH<6.5 PH 6.5-7.5 PH>7.5 Cd 0.30 0.30 0.60 Cr 150 200 250 Cu 50 100 100 Ni 40 50 Pb 250 300 350 Zn 200 250 300 As 30 25 20 Hg 0.30 0.50 1.0 (3)根据对城区不同功能区表层土壤样品重金属浓度计算平均值，得到各个功能区土壤重金属的平均浓度，如表三： 表三：城区不同功能区表层土壤重金属平均浓度 功能区 As（ug/g) Cd（ng/g) Cr(ug/g) Cu(ug/g) Hg(ng/g) Ni（ug/g) Pb(ug/g) Zn(ug/g) 1 6.27 289.96 69.02 49.4 93.04 18.34 69.11 237.01 2 7.25 393.11 53.41 127.54 642.36 19.81 93.04 277.93 3 4.04 152.32 38.96 17.32 40.96 15.45 36.56 73.29 4 5.71 360.01 58.05 62.21 446.82 17.62 63.53 242.85 5 6.26 280.54 43.64 30.19 114.99 15.29 60.71 154.24 (4)城区土壤重金属含量背景值见表四： 表四：城区土壤元素背景值 元素 平均值 As (μg/g) 3.6 Cd (ng/g) 130 Cr (μg/g) 31 Cu (μg/g) 13.2 Hg (ng/g) 35 Ni (μg/g) 12.3 Pb (μg/g) 31 Zn (μg/g) 69 (5)评价结果：将不同功能区表层土壤重金属平均浓度值及城区土壤元素背景值代入Mul1er地积累指数计算公式，得到如下结果，详见表五： 表五：城区不同功能区表层土壤重金属Muller指数评价 功能区 As（ug/g) Cd（ng/g) Cr(ug/g) Cu(ug/g) 分级 分级 分级 分级 1 0.215506 1 0.57238 1 0.569784 1 1.319011 2 2 0.425022 1 1.011459 2 0.199879 1 2.687377 3 3 -0.4186 0 -0.35637 0 -0.25524 0 -0.19306 0 4 0.080531 1 0.884563 1 0.320065 1 1.651646 2 5 0.213203 1 0.524732 1 -0.09158 0 0.60857 1 功能区 Hg(ng/g) Ni（ug/g) Pb(ug/g) Zn(ug/g) 分级 分级 分级 分级 1 0.825534 1 -0.00863 0 0.577166 1 1.195317 2 2 3.612993 4 0.102608 1 1.00062 2 1.425091 2 3 -0.3581 0 -0.25601 0 -0.34696 0 -0.49794 0 4 3.089304 4 -0.06641 0 0.450207 1 1.230435 2 5 1.133112 2 -0.27103 0 0.384703 1 0.575546 1 根据表五得到不同功能区各重金属污染程度的评价结果： （1）一类区（生活区）内的Ni属于无污染；As,Cd，Hg,Pb和Cr属于轻度-中等污染；Cu和Zn属于中等污染。 （2）二类区（工业区）内的As，Cr，Ni属于轻度-中等污染；Cd，Pb和Zn属于中等污染；Cu属于中等-强污染；Hg属于强污染。 （3）三类区（山区）内八种主要重金属元素都属于无污染。 （4）四类区（主干道路区）内的Ni属于无污染；As，Cd，Cr，Pb属于轻度-中等污染；Cu和Zn属于中等污染；Hg属于强污染。 （5）五类区（公园绿地区）内的Cr，Ni属于无污染；As，Cd，Cu，Pb，Zn属于轻度-中等污染；Hg属于中等污染。 5.2重金属污染主要原因 通过对不同功能区各元素相关性分析及主成分提取，并结合不同功能区主要人类活动会导致的污染情况，分析重金属污染的主要原因。 图二 城区内不同功能区分布 （1） 一类区（生活区）： 表六 1类区（生活区）相关性 As Cd Cr Cu Hg Ni Pb Zn As 1 .381* .238 .531** .293 .605** .450** -.017 Cd 1 .349* .499** .397** .283 .802** .346* Cr 1 .376* .150 .527** .416** .412** Cu 1 .198 .434** .502** .238 Hg 1 .211 .340* .242 Ni 1 .300* .334* Pb 1 .328* Zn 1 相关性分析： *. 在 0.05 水平（双侧）上显著相关。 **. 在 .01 水平（双侧）上显著相关。 主成分提取结果： 表七 成份矩阵a 成份 1 2 3 As .669 -.646 -.010 Cd .784 .171 -.417 Cr .643 .234 .493 Cu .729 -.246 .024 Hg .492 .130 -.437 Ni .686 -.253 .523 Pb .803 .112 -.348 Zn .501 .691 .267 提取方法 :主成份。 城市土壤重金属来源于成土母质和人类活动，同一来源的重金属之间存在着相关性，根据相关性可以判断土壤重金属污染来源是否相同 。如果重金属之间存在显著正相关，则其来源可能相同，否则来源可能不止一个，根据相关性分析结果只有Pb与Cd及Ni与As之间呈显著的相关性，相关系数分别为：0．802和0．605，表明上述几种重金属可能具有相同的污染源，而其他没有显著相关性，说明生活区的重金属来源相对复杂，很难就某种污染源予以确认。为了进一步明确城区不同功能区土壤中重金属的污染来源，利用SPSS软件对上述不同功能区土壤重金属的8项指标进行主成分分析，通过主成分分析计算，生活区的8个变量的全部信息可由3个主成分表示。即对前3个主成分进行分析已经能够反映全部数据的大部分信息。第一主成分的特点表现为因子变量在元素Pb与Cd上有较高载荷，而Pb主要在汽车尾气颗粒物含量较多，因此可能是由于生活区进出车辆尾气所致，对于Cd该种重金属在农业肥料，农药中较常见，在塑料薄膜中也有含量，因此其污染可能与工业区附近不合理使用农肥农药的农业生产有关。第二主成分主要在Zn和As上有较高载荷，因大部分居民生活区设有临时停车场，而轮胎与地面的磨擦是Zn产生的一个重要途径，当然，由于居民区生活垃圾的堆放也可能是Zn污染的一个重要原因之一。第三主成分贡献不多，主要在Cr上有载荷，反映了工业污染的原因。综合分析 可以认为该区重金属污染的主要原因为工业污染，农业污染，交通污染及土壤来源。 （2）二类区（工业区）： 表八 2类区（工业区）相关性 As Cd Cr Cu Hg Ni Pb Zn As 1 .329 .380* .153 .181 .690** .395* .518** Cd 1 .541** .566** .533** .489** .829** .754** Cr 1 .920** .902** .698** .675** .695** Cu 1 .983** .503** .670** .622** Hg 1 .479** .612** .590** Ni 1 .578** .634** Pb 1 .739** Zn 1 表九 成份矩阵a 成份 1 2 As .518 .758 Cd .786 .074 Cr .916 -.206 Cu .868 -.463 Hg .845 -.459 Ni .767 .421 Pb .858 .049 Zn .859 .188 提取方法 :主成份。 根据相关性分析结果，Hg和Cu与Cr以及Hg和Cu以及Pb和Cd之间均呈显著的相关性，相关系数分别为：0．920，0．902，0.983,0.829；表明上述几种重金属可能来自于同种污染源。通过主成分分析计算，生活区的8个变量的全部信息可由2个主成分表示。第一主成分的特点表现为因子变量在元素Cr与Cu，Pb,Zn上有较高载荷， Pb在汽车尾气颗粒物种含量较多，因此可能是由于生活区进出车辆尾气所致，而且Pb污染受冶炼厂、化工厂工业废气影响较为突出，尤其是冶炼厂排出的废气中Pb的含量较高，其排放量达2.68t/a，因此以上反映了工业污染可能是其主要的污染原因之一，对于Cu，Zn可能源于居民生活垃圾堆放和工业生产。第二主成分的特点表现为因子变量在元素As上有较高载荷；综合分析，可以认为该区重金属污染的主要原因为工业污染，农业污染，交通污染及土壤来源。 （3）三类区（山区）： 表十 3类区（山区）相关性 As Cd Cr Cu Hg Ni Pb Zn As 1 -.291* .113 .527** .075 .078 -.205 -.176 Cd 1 .066 .090 .246* .049 .766** .606** Cr 1 .364** -.006 .945** .107 .627** Cu 1 .505** .358** .122 .252* Hg 1 -.045 .226 .170 Ni 1 .028 .629** Pb 1 .590** Zn 1 *. 在 0.05 水平（双侧）上显著相关。 **. 在 .01 水平（双侧）上显著相关。 表十一 成份矩阵a 成份 1 2 3 As -.009 .668 .478 Cd .601 -.678 .141 Cr .761 .480 -.365 Cu .517 .469 .615 Hg .324 -.074 .747 Ni .737 .499 -.409 Pb .605 -.630 .173 Zn .905 -.173 -.183 提取方法 :主成份。 结合相关性分析及主成分提取，可以明显观察到只有Ni和Cr以及Pb和Cd有显著相关性，相关系数分别为0.945,0.766；表明上述几种重金属可能来自于同种污染源。通过主成分分析计算，生活区的8个变量的全部信息可由3个主成分，第一主成分的特点表现为因子变量在元素Cr与Zn上有较高载荷，Zn可能源于居民生活垃圾堆放和工业生产以及交通污染。第二主成分的特点表现为因子变量在元素Cr与As上有较高载荷，第三主成分的特点表现为因子变量在元素Cu有较高载荷综合分析，可以认为该区重金属污染的主要原因为工业污染，生活污染，交通污染及土壤来源。 （4）四类区（主干道区）： 表十二 4类区（主干道区）相关性 As Cd Cr Cu Hg Ni Pb Zn As 1 .121 .139 .092 -.004 .228** .060 .188* Cd 1 .373** .424** .211* .351** .615** .294** Cr 1 .894** .012 .869** .428** .395** Cu 1 .032 .886** .506** .432** Hg 1 .040 .266** .118 Ni 1 .396** .503** Pb 1 -.030 Zn 1 **. 在 .01 水平（双侧）上显著相关。 *. 在 0.05 水平（双侧）上显著相关。 表十三 成份矩阵a 成份 1 2 As .235 -.131 Cd .621 .458 Cr .874 -.319 Cu .906 -.250 Hg .170 .743 Ni .888 -.320 Pb .703 .477 Zn .643 .115 提取方法 :主成份。 结合相关性分析及主成分提取，可以观察到Cr与Cu及Ni，Ni与Cu有显著相关性，相关系数分别为0.894,0.896,0.886；表明上述几种重金属可能来自于同种污染源。通过主成分分析计算，生活区的8个变量的全部信息可由2个主成分，第一主成分的特点表现为因子变量在元素Cr，Cu，Ni上有较高载荷，Cr和Cu的污染可能来自于工业生产，Zn,Ni存在于汽车轮胎老化及工业污染。第二主成分主要是Hg上有载荷，可能源于水银制碱工业污染。综合分析，可以认为该区重金属污染的主要原因为工业污染，交通污染及土壤来源。 （5）五类区（公园绿地区） 表十四 5类区（公园绿地区）相关性 As Cd Cr Cu Hg Ni Pb Zn As 1 .358* .689** .107 .176 .691** .265 .285 Cd 1 .564** .500** .054 .433** .598** .712** Cr 1 .357* .023 .739** .397* .509** Cu 1 .136 .267 .756** .521** Hg 1 -.048 .389* .063 Ni 1 .168 .298 Pb 1 .748** Zn 1 表十五 成份矩阵a 成份 1 2 3 As .636 -.569 .349 Cd .811 .103 -.232 Cr .809 -.426 -.003 Cu .679 .457 -.151 Hg .203 .354 .889 Ni .663 -.622 .012 Pb .782 .544 .097 Zn .798 .286 -.242 提取方法 :主成份。 根据相关性分析结果，Ni与Cr以及Pb和Cu以及Pb和Zn之间均呈较显著的相关性，相关系数分别为：0.739，0.756，0.748；表明上述几种重金属可能来自于同种污染源。通过主成分分析计算，生活区的8个变量的全部信息可由3个主分表示。第一主成分的特点表现为因子变量在元素Cr与Cd，Pb,Zn上有较高载荷， Pb在汽车尾气颗粒物种含量较多，因此可能是由于生活区进出车辆尾气所致，而且Pb污染受冶炼厂、化工厂工业废气影响较为突出，尤其是冶炼厂排出的废气中Pb的含量较高，其排放量达2.68t/a，因此以上反映了工业污染可能是其主要的污染原因之一，对于Zn可能源于工业生产。第二主成分的特点表现为因子变量在元素Ni上有较高载荷，可能是交通污染导致。第三主成分的特点表现为因子变量在元素Hg上有较高载荷,应该是工业污染所致。综合分析，可以认为该区重金属污染的主要原因为工业污染，农业污染，交通污染及土壤来源。 5.3建立重金属污染物传播模型 5.3.1 重金属污染物的传播特征 经过查阅相关资料，并结合第一问中重金属污染物在城区的空间分布图，得出金属离子在多孔介质中迁移取决于水相中的分子扩散，固相上的吸附、压缩及固相中的弯曲度，所以重金属污染物以污染源为中心径向传播，且浓度随径向距离增大而急速衰减。因此，采用一维传播方程计算重金属元素在土壤中的浓度分布如下式: (5.3.1-1) 其中c为重金属离子在土壤中的浓度，D为扩散系数，ν为水的通量密度，s为吸附在土壤上重金属离子的浓度,t为时间,r为被考察点到重金属离子污染源的距离，ρ为土壤介质的密度，θ为空隙率。 在本文中，时间对重金属浓度的影响微乎其微，可以忽略，因此， 、 等于0。 所以可得： （5.3.1-2） 解这个二阶线性偏微分方程得： （5.3.1-3），其中 、 为常系数对此式经过反复修正，加一个一次项和常数项，最终得到的结果是： （5.3.1-4），其中 EMBED Equation.DSMT4 、 、 为常系数，此式即重金属离子污染物传播特征方程。 5.3.2 重金属污染物传播模型的建立及验证 建立：由5.3.1-3式建立重金属传播模型： ， （5.3.2-1） 其中x,y为被测点的坐标，则重金属污染物浓度与考察点到污染源的距离成指数关系。接下来要求得常系数 、 、 和k的值。 选取元素分布比较标准的Hg元素分析，利用附件中给出的数据，选择重金属分布集中的点，记录相应的坐标值。同时，用随机点插值的方法，在Matlab中按一定距离插分出一系列点。将插值所得点和真实测量值在Excel中统计分析，取出合理的点组成集合，利用程序在Matlab中拟合可以得到三维立体曲面图如下所示： 图三 拟合所得三维立体曲面图 由此，可以得到常系数 EMBED Equation.DSMT4 、 、 、k的值分别为： 所以重金属污染物传播模型为： （5.3.2-2） 验证：为表明所得重金属污染物传播模型的准确性即能否符合附件中所给不同区域重金属浓度分布，采用将Hg元素被测点坐标值代入模型（5.3.2-2）中，比较真实浓度值与模型计算所得浓度值，得到方差分析如下所示： 根据方差分析可得：所得重金属污染物传播模型是稳定、收敛的，可以作为重金属污染物在土壤中传播的数学模型。 5.3.3 污染源位置的确定 考虑到r和x,y的关系 ，将三维模型在形式上转化为关于r的二维方程， 在matlab工具箱中利用二维cftool进行拟合，通过对拟合函数的选取，可以找到用五次多项式拟合该系列数据吻合很好，该多项式方程为： （5.3.3-1） 拟合图： 图四 五次方多项式拟合图 对（5.3.3-1）式进行极值求解，可得 , 这里所研究的是Hg元素在某一局部区域的传播数据，故找到了在该局部区域内Hg元素的扩散传播和污染源的位置。在该局域的污染源位置图： 图五 污染源位置 5.4模型评价及城市地质环境的演变模式 5.4.1模型的评价与修正 （1）土壤中重金属污染物传播模型： 优点：在误差允许的范围内，较好地解决了小范围内重金属污染物的扩散问题；而且模型方程形式简单，便于理解并分析城市不同区域污染的程度及主要原因；根据方程可以求得污染源的位置坐标，确定污染源的位置，有利于城市对污染问题进行处理。 局限性：首先，由于采用网格取样法，而污染物浓度与距离成指数关系，使得取点不能很好地反映污染物分布；其次，采用局部范围内某一种元素分析，使建立的模型不能充分表征整个城市同区域不同污染物和不同区同种污染物的分布；再有，忽略了海拔高度、水体、风等自然因素对污染物浓度的影响；还有，污染物不可能严格按照一维、径向扩散。 （2）对方程的进一步修正： 考虑到实际测得的点在相对较短的距离上过于分散，及重金属元素在实际扩散中可能遇到的吸附、流失、富集等因素，建立在较短距离上衰减更为强烈的方程，即将指数的一次项改为两次项，并加入调和因子m,运用Matlab进行拟合，得到了更加接近实际变换规律的拟合曲线。方程具体形式及拟合结果如下： General model: Coefficients (with 95% confidence bounds): Goodness of fit: SSE: 8.04e+004 R-square: 0.9996 Adjusted R-square: 0.9995 RMSE: 94.52 得到进一步修正方程： 图六 拟合曲线图 由图六可知，金属元素分布在0m至40m内衰减剧烈，因此在实际取点时不易取到衰减过程中的点。而在50m以后曲线几乎接近0点，说明由于重金属元素受富集、吸附、流失、地势等其他环境因素的影响损失严重。该修正结果对真实数据的吻合度很高。 5.4.2 城市地质环境的演变模式 地质环境评价是指地质环境对人类生存与发展适宜性的综合性评价，为区域经济开发和环境建设与保护提供决策依据。地质环境评价主要是依据环境地质问题与地质灾害对人类生存与发展的不利影响，按照“无问题（灾害）即优良”的基本原则，作出安全意义上的好坏评判。 本题中，城市表层土壤已经受到多种重金属元素的污染，如果得不到及时解决将很有可能危及到城市居民今后的正常生活及生产。因此，研究地质环境演变模式，对于城市地质灾害预测、城市生产居住规划产生很好的引导作用。 为研究城市地质环境的演变模式，还应收集有关该城市的地形地貌、河流分布、主导风向、气候特征等信息；此外，经取样所获数据对于建立演变模型明显不够，故需有更细致的取样方式来获得足够的土壤中重金属浓度数据。 在已经建立的传播方程的基础上，再考虑到风、水的搬运迁移作用、城市地势、气候等多种因素的影响，需要对已建方程进行进一步修正，以趋近于实际地质环境。分析重金属元素空间分布图可知，本城东北部地势偏高、西南部地势较低，由东北向西南地势逐级降低。此种地势条件产生的海拔、风向、水文特征影响了某些元素的分布及扩散。因此可将海拔看做影响重金属分布的环境参量之一代入方程。另外，本地区常年的风速值会通过元素转移升降改变分布状态，因此也应为一个重要的参量。考虑到元素迁移过程中金属离子在多孔介质中迁移取决于水相中的分子扩散，固相上的吸附、压缩及固相中的弯曲度，因此也应将这些参数考虑到方程中。 综上所述，修正方程不仅应该考虑元素迁移扩散的基本规律，还应考虑实际环境参数。由此可见，重金属元素的迁移是一个多种因素并存的问题，仅靠单一因素不能解决实际问题。 由于现代地质环境科学研究的问题是多个有关学科的前沿问题，特别是多学科交叉的边缘科学领域问题。如今随着科技的发展，人类社会也产生了一些系统科学的理论和方法去解决城市地质问题，如网络信息技术、卫星通讯技术和高速信息传输技术、遥感技术、地理信息系统（GIS)、全球定位系统(Glob～l positioning systems，GPS)、技术生物工程技术 、高精度的分析测试技术和方法、高分辨率的年代学测年技术，相信通过这些高新技术的大量推广及应用，城市地质环境问题、以及和地址有关的污染问题将会得到更好的解决。 参考文献 [1] 段雪梅，蔡焕兴，巢文军，南京市表层土壤重金属污染特征及污染来源，环境科学与管理，第35卷第10期，2010。 [2] 钱翌，赵世刚，青岛市不同生态功能区表层土壤重金属污染初步评价，中国农学通报，2010。 [3] 阴雷鹏，赵景波，西安市主要功能区表层土壤重金属污染现状评价，陕西师范大学学报（自然科学版），第34卷第3期，2006。 [4] 刘志钧，刘小华，程建光，应用统计学方法建立矿山表层土壤重金属污染基线模型，能源环境保护，第22卷第5期，2008。 [5] 胡恭任，于瑞莲，应用地积累指数和富集因子法评价324国道塘头段两侧土壤的重金属污染，中国矿业，第17卷第4期，2008。 [6] 袁淑君，孟庆茂，数据统计分析—SPSS/PC+原理及其应用，北京师范大学出版社，1995。 [7] 王正东，数学软件与数学实验（第二版），科学出版社，2010。 [8] 袁震东，数学建模方法，华东师范大学出版社，2003。 [9] 同登科，周生田，张高民，计算方法，中国石油大学出版社，2009。 [10] 章毛连，王祥科，陈磊，重金属离子在土壤中迁移的模型研究，吉林大学报（理学版），第44卷第4期，2006。 附录 各区域内元素相关性和主成分分析 1. 公因子方差 初始 提取 VAR00001 1.000 .865 VAR00002 1.000 .818 VAR00003 1.000 .712 VAR00004 1.000 .592 VAR00005 1.000 .450 VAR00006 1.000 .809 VAR00007 1.000 .778 VAR00008 1.000 .800 提取方法：主成份分析。 解释的总方差 成份 初始特征值 提取平方和载入 合计 方差的 % 累积 % 合计 方差的 % 累积 % 1 3.616 45.199 45.199 3.616 45.199 45.199 2 1.133 14.165 59.365 1.133 14.165 59.365 3 1.075 13.432 72.797 1.075 13.432 72.797 4 .807 10.083 82.880 5 .524 6.544 89.424 6 .445 5.556 94.980 7 .236 2.954 97.934 8 .165 2.066 100.000 提取方法：主成份分析。 2. 公因子方差 初始 提取 VAR00009 1.000 .844 VAR00010 1.000 .623 VAR00011 1.000 .880 VAR00012 1.000 .967 VAR00013 1.000 .924 VAR00014 1.000 .765 VAR00015 1.000 .739 VAR00016 1.000 .774 提取方法：主成份分析。 解释的总方差 成份 初始特征值 提取平方和载入 合计 方差的 % 累积 % 合计 方差的 % 累积 % 1 5.254 65.670 65.670 5.254 65.670 65.670 2 1.263 15.786 81.455 1.263 15.786 81.455 3 .781 9.763 91.218 4 .267 3.336 94.554 5 .227 2.834 97.388 6 .149 1.868 99.256 7 .051 .634 99.890 8 .009 .110 100.000 提取方法：主成份分析。 3. 公因子方差 初始 提取 VAR00001 1.000 .674 VAR00002 1.000 .840 VAR00003 1.000 .943 VAR00004 1.000 .865 VAR00005 1.000 .668 VAR00006 1.000 .959 VAR00007 1.000 .793 VAR00008 1.000 .883 提取方法：主成份分析。 解释的总方差 成份 初始特征值 提取平方和载入 合计 方差的 % 累积 % 合计 方差的 % 累积 % 1 3.042 38.022 38.022 3.042 38.022 38.022 2 2.036 25.445 63.467 2.036 25.445 63.467 3 1.549 19.362 82.829 1.549 19.362 82.829 4 .689 8.616 91.445 5 .249 3.116 94.562 6 .231 2.883 97.445 7 .161 2.017 99.463 8 .043 .537 100.000 提取方法：主成份分析。 4. 公因子方差 初始 提取 VAR00009 1.000 .072 VAR00010 1.000 .595 VAR00011 1.000 .865 VAR00012 1.000 .884 VAR00013 1.000 .581 VAR00014 1.000 .891 VAR00015 1.000 .723 VAR00016 1.000 .427 提取方法：主成份分析。 解释的总方差 成份 初始特征值 提取平方和载入 合计 方差的 % 累积 % 合计 方差的 % 累积 % 1 3.751 46.883 46.883 3.751 46.883 46.883 2 1.287 16.084 62.967 1.287 16.084 62.967 3 .994 12.424 75.391 4 .749 9.359 84.751 5 .685 8.560 93.311 6 .336 4.195 97.505 7 .116 1.450 98.955 8 .084 1.045 100.000 提取方法：主成份分析。 5. 公因子方差 初始 提取 VAR00017 1.000 .851 VAR00018 1.000 .722 VAR00019 1.000 .835 VAR00020 1.000 .693 VAR00021 1.000 .956 VAR00022 1.000 .827 VAR00023 1.000 .917 VAR00024 1.000 .777 提取方法：主成份分析。 解释的总方差 成份 初始特征值 提取平方和载入 合计 方差的 % 累积 % 合计 方差的 % 累积 % 1 3.907 48.840 48.840 3.907 48.840 48.840 2 1.615 20.182 69.023 1.615 20.182 69.023 3 1.056 13.202 82.225 1.056 13.202 82.225 4 .546 6.829 89.054 5 .336 4.203 93.257 6 .242 3.030 96.288 7 .213 2.659 98.947 8 .084 1.053 100.000 提取方法：主成份分析。 6. 插值找点程序matlab： xi=12400:500:15810; yi=1063:500:3512; ci=griddata(x,y,c,xi,yi','cubic') subplot(1,2,1), plot(x,y,'*') subplot(1,2,2), mesh(xi,yi,ci) ci = 1.0e+004 * NaN 0.0002 -0.0237 -0.0414 -0.0286 -0.0116 -0.0039 NaN 0.0429 0.0992 0.2525 0.1333 -0.0103 -0.0277 NaN 0.5032 1.0691 1.2053 0.7050 0.0619 0.0006 NaN 0.5215 1.2201 1.2331 0.5101 -0.0245 -0.0449 NaN 0.0743 0.2902 0.3193 0.0299 -0.0300 NaN >> xi xi = Columns 1 through 6 12400 12900 13400 13900 14400 14900 Column 7 15400 >> yi yi = 1063 1563