《论文2011模型.pdf》由会员分享，可在线阅读，更多相关《论文2011模型.pdf（17页珍藏版）》请在taowenge.com淘文阁网|工程机械CAD图纸|机械工程制图|CAD装配图下载|SolidWorks_CaTia_CAD_UG_PROE_设计图分享下载上搜索。
1、-1-A A 题题 城城 市市 表表 层层 土土 壤壤 重重 金金 属属 污污 染染 分分 析析 摘摘 要要 本文利用采样获得的数据对该城区 8 种重金属污染的空间分布进行分析，并对多个功能区的污染程度做了定量分析，通过建立克里金(Kriging)模型，确定污染源位置。鉴于对样本进行的数据处理，使用 Mathematica 软件绘制了城区地形图与各元素在城区的分布图，初步说明了污染物的空间分布状况。我们进一步借助单因子指数模型和内梅罗指数评价模型定量分析不同功能区重金属的污染程度，结果如下表所示.各个功能区综合污染指数 区域 生活区 工业区 山区 交通区 公园绿地区 综合综合 0.72 1.0
2、1 0.25 0.77 0.49 结果表明，总体说来，该城市工业区受到轻度污染，交通区接近污染，其它区无污染。接下来，我们再通过现状值与背景值比对进行检验。结果表明汞(Hg)、铜(Cu)、铅(Pb)和锌(Zn)元素在生活区、工业区、交通区积累明显，说明了这些重金属污染的主要来源有三个方面：居民生活污染源、交通污染源、工矿企业污染源。元素砷(As)、镉(Cr)、镍(Ni)的浓度在各功能区相差不大，说明了这些元素主要是地质起源，受人类活动影响较小。进一步，通过对污染物的传播特征的分析，我们建立了克里金模型计算某一区域重金属的最高浓度，进而预测出污染源的位置。采用普通克里金法，在约束条件下根据采样点
3、的含量对城区土壤重金属含量进行最优、无偏估值。由于矿物质含量在一定的空间区域内成结构性变化的这一性质由变差函数进行刻画，而各污染物元素的浓度同样在一定的空间区域内成结构性的变化，因而我们可以类比矿物质的含量变化。通过对样本数据进行函数拟合，得到所需的变差函数，运用克里金方法，预测出指定区域当中特定 9个点的浓度值，作为这一区域的预测，进而通过比较确定出污染源，结果如下表所示：8 种重金属污染源预测位置 坐标 As Cd Cr Cu Hg Ni Pb Zn x 18134 16439 3799 1883 1708 2380 4670 2083 y 10546 10383 6518 4192 27
4、95 3590 4396 3292 z 41 45 4 7 22 7 8 7 通过克里金模型不仅能够预测结果，还能获得预测误差，这样有利于评估预测结果的不确定性，并且它是在有限区域对区域内变化量的取值进行无偏，最优估计的一种方法。最后，我们提出，如果能收集到重金属随时间变化的数据，再收集到河流流向信息、风向信息、交通流量等重要数据信息，就可以建立重金属迁移的对流扩散微分方程模型，得到地质环境的演变模式。关键词：关键词：城区土壤，重金属污染，指数评价，变差函数，克里金法.-2-一、一、问题重述问题重述 随着城市经济的发展和城市人口的不断增加，人类对城市环境质量的影响日显突出。现在对某城市城区土壤
5、地质环境进行调查，通过对不同功能区的调查所获得的数据建立模型。我们将要解决如下问题：1.给出 8 种主要重金属元素在该城区的空间分布，并分析该城区内不同区域重金属的污染程度；2.通过数据分析，说明重金属污染的主要原因；3.分析重金属污染物的传播特征，由此建立模型，确定污染源的位置；4.分析所建立模型的优缺点，为更好地研究城市地质环境的演变模式，提出一些改进方案。二、二、模型假设模型假设 1.忽略金属元素之间的污染相关性。2.忽略动植物对重金属元素吸收含量。3.污染物的评价标准采用国家土壤质量二级标准。4.假设重金属元素浓度的局部平均值与总体平均值无关。5.假设区域化变量存在空间相关性。6.假设
6、区域化变量是平稳的，不存在偏差。三、三、符号说明符号说明 :土壤中污染物 i 的环境质量指数；:土壤中污染物 i 的实测值；:土壤中污染物 i 的评价标准；综合综合 :各区域污染综合指数；:污染物 i 的标准偏差；:污染物 i 的现状值；：表示样本点 i 到样本点 j 的距离;()：未知点 p 的估值;()：p 点的真值.四、四、问题分析问题分析 环境质量是各个环境要素优劣的总和概念。衡量环境优质优劣的因素很多，通常用环境中污染物质的含量来表达。环境质量指数就是这样一个有代表性的无量纲数，是质量好坏的表征，表示污染物在环境中实际浓度超过评价标准的程度及超标倍数。分析重金属的污染程度，就是需要建
7、立一个重金属污染物对环境影响程度的评价标-3-准。在参考大量文献后，我们选择了单因子指数模型以及内梅尔综合指数模型作为评价标准。并且，由于样品数据的随机性，可能造成结果出现偏差。因此，我们运用了现状值检验来消除一些极端数据的影响。在预测污染源的问题中有大量的样本数据，我们考虑运用经典统计学的方法，找出样本点数据的统计规律，进而预测出各种重金属的污染源。但经典统计方法在地质预测问题中又存在着一定的局限性。于是我们找到了一种既能保持概率统计的有效性，又考虑到地学变量特点的方法，即克里金(Kriging)方法。最后，通过所建立的模型，我们研究了城市地质环境的演变模式，并为地质环境污染的预测，提供了进
8、一步的改进方案。五、五、模型分析与求解模型分析与求解 首先我们根据数据，使用 Mathematica 软件做出了样本点及其所在区域的地貌三维图，如图 1 所示：图 1：城市地形图 由图 1，我们对抽象的样品数据有一个直观上的认识，进而可以分析出更详尽的地形特点：东西长 30000 米，南北宽 20000 米，海拔高度为 0 到 150 米。整体地势成东高西低，南北两侧隆起，中间低凹。-4-问题一问题一分析不同区域元素的污染程度分析不同区域元素的污染程度 根据假设 1 和假设 2，在无其他因素的影响下，重金属的迁移与积累处于动态平衡的状态。研究方法：1.1.数据来源：数据来源：本研究主要采用城市
9、土壤单点样采集。按功能将城区分为五个地区，将城区划分为间距 1 公里左右的网络子区域，按照每平方公里 1 个采样点取样，共取 319个表层土（010 厘米深度）样品，并编号，记录取样点并用 GPS 记录取样点的位置。2.2.数据分析：数据分析：将数据整理汇总，用 Mathmatica 软件将研究地区数字化作出采样点分布图，利用软件做出各元素的分布图，并根据土壤污染国家评定标准进行分析，按评价标准有：单因子指数公式：=;综合污染指数法：综合综合=()+().表 1：土壤环境质量标准值二级标准 As(g/g)Cd(ng/g)Cr(g/g)Cu(g/g)Hg(ng/g)Ni(g/g)Pb(g/g)Z
10、n(g/g)二级标准=30=600=200=200=500=50=300 250 其公式中含有评价参数中最大的单项污染指数，突出了污染指数最大的污染物对环境的影响和作用，克服了平均值法各个污染物分担的缺陷，但是其没有考虑土壤中各污染物对作物毒害的差别，而且最大值对所的结果的影响很大。根据内梅罗综合指数可将土壤分为 5 个等级，如表 2 所列.表 2：内梅罗法土壤污染分级标准 等级划分 P综合 污染等级 污染水平 1 P综合 0.85 非污染 清洁 2 0.85 P综合 1.7 轻度污染 土壤轻污染，作物已受损 3 1.7 P综合 2.56 中度污染 土壤、作物均受重污染 4 2.56 P综合
11、重度污染 土壤、作物均受污染相当严重 根据综合污染指数法，可计算出各个区域的污染指数，结果如表 3 所示.-5-表 3：城市不同功能区的综合污染指数 区域 生活区 工业区 山区 交通区 公园绿地区 P综合 0.72 1.01 0.25 0.77 0.49 通过表 3 的数据与表 2 作比对可知：工业区为中度污染，而生活区、山区、公园绿地、交通区为非污染区。3.3.计算与讨论计算与讨论 各元素选取的样本数据不成正态分布,必须对其进行逐步剔除，使数据接近正态分布：先计算所有样本的平均值和标准差。然后进行标准差检验，即：将超出算术平均值加减2倍的标准差的实测值舍弃，不参加现状值的计算，直到全部样本无
12、法被剔除为止。根据这些数据计算出来的平均值均 加减偏差值 就作为研究区的土壤现状值。现状值公式:iiiXXS均 标准偏差公式:1 22111niijijSXXN均 计算结果列于表 4表 8.表 4：生活区各元素统计结果 元素 平均值 标准差 现状值（上限）现状值（下限）As(g/g)6.04 1.90 7.93 4.14 Cd(ng/g)260.89 123.80 384.69 137.09 Cr(g/g)53.31 28.33 81.64 24.98 Cu(g/g)41.95 31.46 73.41 10.48 Hg(ng/g)70.80 55.38 126.18 15.42 Ni(g/g)
13、17.68 4.88 22.57 12.80 Pb(g/g)55.70 26.70 82.40 29.01 Zn(g/g)175.23 171.97 347.20 3.26 表 5：工业区各元素统计结果 元素 平均值 标准差 现状值（上限）现状值（下限）As(g/g)6.49 2.88 9.37 3.62 Cd(ng/g)334.10 135.51 469.61 198.58 Cr(g/g)46.78 19.04 65.81 27.74 Cu(g/g)58.94 53.40 112.33 5.54 Hg(ng/g)274.99 427.48 702.47 152.49 -6-Ni(g/g)19
14、.19 7.59 26.78 11.60 Pb(g/g)74.50 36.22 110.72 38.28 Zn(g/g)204.75 174.86 379.61 29.89 表 6：山区各元素统计结果 元素 平均值 标准差 现状值（上限）现状值（下限）As(g/g)3.79 1.37 5.16 2.42 Cd(ng/g)141.20 60.39 201.59 80.81 Cr(g/g)36.02 16.81 52.83 19.20 Cu(g/g)15.85 6.81 22.66 9.05 Hg(ng/g)38.40 18.75 57.16 19.65 Ni(g/g)14.07 6.44 20.
15、51 7.63 Pb(g/g)33.38 9.71 43.09 23.67 Zn(g/g)69.30 20.65 89.95 48.65 表 7：交通区各元素统计结果 元素 平均值 标准差 现状值（上限）现状值（下限）As(g/g)5.31 1.68 6.99 3.64 Cd(ng/g)322.58 175.52 498.10 147.06 Cr(g/g)48.34 20.02 68.36 28.33 Cu(g/g)52.71 44.63 97.33 8.08 Hg(ng/g)125.63 247.09 372.72 121.45 Ni(g/g)16.71 4.94 21.65 11.76 P
16、b(g/g)59.69 26.02 85.72 33.67 Zn(g/g)194.81 153.88 348.69 40.92 表 8：公园绿地区各元素统计结果 元素 平均值 标准差 现状值（上限）现状值（下限）As(g/g)5.96 1.65 7.61 4.32 Cd(ng/g)219.27 123.97 343.24 95.30 Cr(g/g)41.01 10.21 51.22 30.80 Cu(g/g)26.86 11.44 38.31 15.42 Hg(ng/g)78.98 71.19 150.18 7.79 Ni(g/g)14.22 3.61 17.82 10.61 Pb(g/g)5
17、2.04 29.02 81.06 23.02 Zn(g/g)117.91 85.73 203.65 32.18 问题二问题二说明重金属污染的主要原因说明重金属污染的主要原因 经过上述问题的大量数据计算分析，将各个区域中各个元素的现状值与在那些远离人群及工业活动的自然区取样得到的背景值进行比对，进而研究污染原因。-7-图 2：城市功能区图 表 9：8 种主要重金属元素的背景值(供分析以下问题所用)元素 平均值 标准偏差 范围 As(g/g)3.6 0.9 1.85.4 Cd(ng/g)130 30 70190 Cr(g/g)31 9 1349 Cu(g/g)13.2 3.6 6.020.4 Hg
18、(ng/g)35 8 1951 Ni(g/g)12.3 3.8 4.719.9 Pb(g/g)31 6 1943 Zn(g/g)69 14 4197 由表 4表 8 以及表 9 的分析：（1）研究区砷(As)现状值为土壤背景值的 0.961.74 倍。由图 2 与图 3(左)，研究区内的 As 没有地质积累，更不存在 As 的污染。As 的含量与汽车尾气排放和汽车轮胎磨损所产生的气体有关，得出在工业区与交通区含量较高。-8-图 3：左：砷(As)的浓度，右：镉(Cd)的浓度 （2）Cd 现状值为土壤背景值的 1.062.62 倍。研究区 Cd 元素的分布见图 3（右），Cd也常出现在工业与交通
19、区。由图可以看出，研究区的西部 Cd 含量较高，主要来源是道路交通区。（3）Cr 现状值为土壤背景值的 1.051.67 倍，接近背景值，所以研究区内 Cr 元素含量受人为影响不大。图 4：左：铬(Cr)的浓度，右：铜(Cu)的浓度 （4）Cu 现状值为土壤背景值的 1.15.5 倍，土壤中的 Cu 主要来源于含铜矿的开采和冶炼厂“三废”的排放，有图 4(右)可看出西部工业区附近的道路土壤 Cu 含量较为集中，因此此地污染来源可能来自矿产的运输。-9-图 5：左：汞(Hg)的浓度，右：镍(Ni)的浓度 （5）Hg 现状值为土壤背景值的 1.1213.77 倍，大气汞的干湿沉降也可以引起土壤中汞
20、的含量增高。大气汞通过干湿沉降进入土壤后，被土壤中的粘土矿物和有机物的吸附或固定，富集于土壤表层，造成土壤汞的浓度的升高，造成 Hg 在五个区域内的含量都较高。（6）Ni 现状值为土壤背景值的 0.91.35 倍，基本都在背景值范围内，人类活动对 Ni在土壤中的含量影响不大，在少数居民区 Ni 含量较高，与生活污水及工业废水有关。图 6：左：铅(Pb)的浓度，右：锌(Zn)的浓度 （7）土壤中的 Pb 主要来源于大气沉降，电池、燃料和焊接等工业中的应用。Pb 现状值为土壤背景值的 1.892.57 倍，西部城区边缘工业区与道路区影响 Pb 含量较高。东部山区及南部部分绿化带铅(Pb)含量在背景
21、值范围内。（8）土壤中 Zn 主要来源于含锌工业如合金、橡胶、玻璃等排放。Zn 现状值为土壤背景-10-值的 0.933.91 倍。由图 6(右)，西部交通、工业集中地区及南部居民区锌(Zn)的含量相对较高。可见含锌(Zn)高的地区是工业、交通密集处及人类活动频繁处。问题三问题三建模求污染源位置建模求污染源位置 首先，搜索出各元素物质浓度最大的样本点的坐标。然后分别以这些点为中心做出一个小区域，运用克里金(Kriging)算法预测出每个小区域中均匀分布的的 9 个点的浓度值。进而通过比较，将每一个小区域的最高浓度值所对应的点，作为每种重金属污染源。为了运用 Kriging 算法预测出每个小区域
22、的未知点的浓度值，必须通过样本点的数据拟合出变差函数。首先，计算出所有样本点中任意两个点的距离，形成集合 d（集合 d 中不允许出现相同的元素）。然后，做出以元素对,i,j作为元素的集合 pt。其中，i，j 表示样本 点 的 编 号，表 示 样 本 点 i 到 样 本 点 j 的 距 离。做 出 滞 后 距 离 集 合h=*1,2,3,+,用以对集合 pt 中元素进行分类，做出变差函数的离散点，即：=12()()()2/23/2 其中()表示 i 点的浓度值，()表示 pt 中满足/2 3/2 的元素个数。再次，通过变差函数的离散点，可以拟合出我们所需要的变差函数f(x)。最后，通过f(x)就
23、可以解一个，预测点 p 的浓度值的 Kriging 方法的方程：A w=b 其中，A=(11)(1)(1)1()1 1 10)，b=(1)()1)，w为点 p 的权重向量。得到w为后，将w的最后一个分量剔除。最后：c(p)=1 求得预测点 p 的浓度值c(p)。按照上述 Kriging 算法，我们详细计算了 As,Cd,Cr 等 8 种重金属污染物的区域预测点。然后经过比较，筛选出了最高浓度点的坐标值，并将其作为重金属污染源。其 8 种重金属污染源的地理位置如下：表 10：预测重金属污染源位置 As Cd Cr Cu Hg Ni Pb Zn X 18134 16439 3799 1883 17
24、08 2380 4670 2083 -11-y 10546 10383 6518 4192 2795 3590 4396 3292 z 41 45 4 7 22 7 8 7 为了测试 Kriging 算法的有效性并分析这个方法的误差，我们利用误差公式：=z(p)z(p)以重金属污染物 As 为例，分别将 As 的样本点当做预测点进行计算，然后与实际样本点的浓度值进行对比，得到下表：As 数据 1 数据 2 数据 3 预测值 7.32 9.15 7.22 实测值 6.77 9.39 7.44 相对误差 8.1%2.5%3.0%从表中可以明显看出，通过 Kriging 算法建立的模型，相对误差较小
25、，符合我们的预测要求，可以对污染源进行预测。六六、模型模型评价评价 通过 Kriging 模型不仅能够预测结果，还能获得预测误差，这样有利于评估预测结果的不确定性，并且它是在有限区域对区域内变化量的取值进行无偏，最优估计的一种方法。对于单因子模型和内梅罗综合指数模型，影响土壤环境的因素很多，而每一种因素的影响程度各不相同，而又会因所处的自然环境条件、时间等的差异而表现出不同的特征。这些因素互相联系，关系极为复杂，不易确定，成相关的模糊性，应在土壤质量评价中，需尽量多选择有代表性的影响因子，在单因子指数法和综合污染指数法的评价基础上，选用一些其他的评价方法来评价土壤质量。对于 Kriging 模
26、型，在拟合变差函数时，由于需要数据具有一定的正态性，所以对样本点的统计性要求较高。并且在拟合时，对函数的选取也有一定的偏差。在解 Kriging方程时，往往会出现坏条件数的系数矩阵，对求解方程的准确性有一定的影响。问题四问题四地质环境的演变模式地质环境的演变模式 在以上我们所建立的模型基础上，如果要分析出地质环境的演变模式，还需获得河流流向信息、风向信息、交通流量等重要数据信息。通过新获得的数据，建立出重金属迁移的对流扩散微分方程模型，得到地质环境的演变模式，为长期预测提供方便。-12-参考文献参考文献 1 姚 德 等 青 岛 城 区 土 壤 重 金 属 环 境 地 球 化 学 研 究 J,中
27、 国 地 质，135(3):539-549,2008.6 2王军，徐春晓，陈芳 铜陵林冲尾矿库复垦土壤的重金属污染评价J,合肥工业大学学报（自然科学版）,28(2):142-144,2005.2 3周爱国、周建伟、梁合诚等，地质环境评价M，中国地质大学出版社，2008.5 4阳正熙，吴堑虹，地学数据分析教程M，科学出版社，2008.9 5张长波，骆永明，吴龙华，土壤污染物源解析方法及其应用研究进展J,土壤,39(2)：190-195,2007 6张长波，李志博等，污染场地土壤重金属含量的空间变异特征及其污染源识别指示意义J,土壤,38(5)：525-533,2006 7张蕾，周启星，城市地表径
28、流污染来源的分类与特征J,生态学杂志,29(11)：2272-2279,2010.9 8李保杰，于法展，纪亚洲，徐州市九里矿区土壤重金属插值分析及污染评价J,测绘科学，35(6)：166-169,2010 9孙梅，杨富贵，蒋恒毅等，青岛城区土壤 Hg 元素环境地球化研究J,山东理工大学学报，22(2)：37-40,2008.3 10国家环境保护局，土壤环境质量标准 GB15618-1995S,北京：中国标准出版社，1995:84-86 -13-附附 录录 ListPlot3DdataLocal,ColorFunction(Hue#&),ClippingStyleAutomatic orData
29、=ImportorData.txt,Table;metal=Importmetal.txt,Table;asData=Tablemetali,1,i,Lengthmetal;cdData=Tablemetali,2,i,Lengthmetal;crData=Tablemetali,3,i,Lengthmetal;cuData=Tablemetali,4,i,Lengthmetal;hgData=Tablemetali,5,i,Lengthmetal;niData=Tablemetali,6,i,Lengthmetal;pbData=Tablemetali,7,i,Lengthmetal;znD
30、ata=Tablemetali,8,i,Lengthmetal;-14-locData1=TableTakeorDatai,3,i,LengthorData;p=FlattenPositionasData,MaxasData 84 x0,y0=TakeorDatap,2 18134,10046 dq=2000;xMin,xMax=x0-dq,x0+dq 16134,20134 yMin,yMax=y0-dq,y0+dq 8046,12046 asS=;Do IfxMin=locData1i,1xMax&yMin=locData1i,2yMax,AppendToasS,locData1i,asD
31、atai ,i,LengthlocData1 asS asN=TableasSi,2,i,LengthasS locData=TableasSi,1,i,LengthasS 17044,10691,93,17087,11933,43,18413,11721,88,19007,11488,84,18738,10921,53,17814,10707,64,18134,10046,41,17198,9810,37,17144,9081,20,18393,9183,26,19767,8810,46,20101,10774,40,16301,8299,24,17904,8287,25,19072,851
32、9,36,16428,9069,20,16289,10072,43,16267,11058,60 pt=;Do DoAppendTopt,NNormlocDatai-locDataj,j,i ,j,i-1 ,i,LengthlocData ppt=Sortpt d=Tableppti,1,i,Lengthppt Lengthd 153 -15-bcH=200;h=TablebcH*i,i,RoundMind/bcH+1,RoundMaxd/bcH-1 828,1028,1228,1428,1628,1828,2028,2228,2428,2628,2828,3028,3228,3428,362
33、8,3828,4028,4228 vAs=Table k=0;sum=0;sum1=0;Do Ifhj-bcH/2ppti,1hj+bcH/2,sum=sum+(asNppti,2,1-asNppti,2,2)2;sum1=sum1+ppti,1;k+,i,Lengthppt;sum1/k,sum/(2k),j,Lengthh;pvAs=ListPlotvAs fvCd=FitvAs,1,x,x2,x3,x-16.7282+0.149207 x-0.0000777673 x2+1.0389?10-8 x3 fx_:=-16.72818567912944+0.149206893916502 x-
34、0.00007776727809825729 x2+1.0388981351333464*-8 x3 fcd=PlotfvCd,x,1,5000 15002000250030003500400020406080100-16-ShowpvAs,fcd asA=Table Whichij,0,iLengthlocData1+1|j LengthlocData1+1,1,True,fNNormlocData1i-locData1j ,i,LengthlocData1+1,j,LengthlocData1+1 bp_:=TablefNNormlocData1i-p,i,LengthlocData1 label=FlattenPositionasData,MaxasData 1000200030004000500020406080100020003000400050002020406080100-17-84 label=96 96 pLabel=locData1label 20554,11228,43 bv=bpLabel;AppendTobv,1;w=LinearSolveasA,bv w=Dropw,-1;w.asData 7.31951 MaxasData 30.13 asDatalabel 6.77