2011年全国大学生数学建模优秀论文.doc

2011高教社杯全国大学生数学建模竞赛承 诺 书我们仔细阅读了中国大学生数学建模竞赛的竞赛规则.我们完全明白，在竞赛开始后参赛队员不能以任何方式（包括电话、电子邮件、网上咨询等）与队外的任何人（包括指导教师）研究、讨论与赛题有关的问题。我们知道，抄袭别人的成果是违反竞赛规则的, 如果引用别人的成果或其他公开的资料（包括网上查到的资料），必须按照规定的参考文献的表述方式在正文引用处和参考文献中明确列出。我们郑重承诺，严格遵守竞赛规则，以保证竞赛的公正、公平性。如有违反竞赛规则的行为，我们将受到严肃处理。我们参赛选择的题号是（从A/B/C/D中选择一项填写）： A 我们的参赛报名号为（如果赛区设置报名号的话）： 所属学校（请填写完整的全名）： 参赛队员 (打印并签名) ：1. 2. 3. 指导教师或指导教师组负责人 (打印并签名)： 日期： 年 月 日赛区评阅编号（由赛区组委会评阅前进行编号）：2011高教社杯全国大学生数学建模竞赛编 号 专 用 页赛区评阅编号（由赛区组委会评阅前进行编号）：赛区评阅记录（可供赛区评阅时使用）：评阅人评分备注全国统一编号（由赛区组委会送交全国前编号）：全国评阅编号（由全国组委会评阅前进行编号）：城市表层土壤重金属污染分析摘要针对城市表层土壤重金属污染的问题，本文建立了相应的数学模型，解决了城市重金属污染的空间分布、污染程度、污染因素及污染源确定等问题，在建模的过程中运用了MATLAB和EXCEL软件。对于问题1，运用MATLAB软件对数据进行处理，给出8种重金属元素在该城区的空间分布图。再根据单项指数法和内梅罗综合指数法，求出该城区不同区域重金属的污染程度为山区土壤环境质量属于轻度污染，公园绿地区土壤环境质量属于中度污染的，生活区、工业区和主干道区土壤环境质量属于重度污染，而且工业区的污染程度要大于主干道区，更远远大于生活区。对于问题2，先求出土壤中元素含量的变异系数，分析变异系数得到土壤中重金属元素污染程度和分散情况，再求出单项污染指数的变异系数，根据该变异系数得到土壤各重金属元素单因子污染指数值的平均变异程度,综合上述两个步骤最终确定出该城区重金属污染的主要原因是工业生产中“三废”的排放和城市生活垃圾的堆积。对于问题3，首先分析重金属的传播特征，建立多元线性回归模型，再对每个因素线性拟合得到8种重金属元素所对应的多元线性回归函数，并对其进行误差分析，最后解出污染源位置。对于问题4，首先对模型进行评价，再确定影响重金属污染物的传播特征的其他因素，最后对所有因素进行综合考虑，优化多元线性回归模型，并加以推广。关键词：空间分布 单项污染指数法 内梅罗指数法 变异系数 传播特征 多元线性回归 23一 问题重述随着城市经济的快速发展和城市人口的不断增加，人类活动对城市环境质量的影响日显突出。对城市土壤地质环境异常的查证，以及如何应用查证获得的海量数据资料开展城市环境质量评价，研究人类活动影响下城市地质环境的演变模式，日益成为人们关注的焦点。按照功能划分，城区一般可分为生活区、工业区、山区、主干道路区及公园绿地区等，分别记为1类区、2类区、5类区，不同的区域环境受人类活动影响的程度不同。现对某城市城区土壤地质环境进行调查。为此，将所考察的城区划分为间距1公里左右的网格子区域，按照每平方公里1个采样点对表层土（010 厘米深度）进行取样、编号，并用GPS记录采样点的位置。应用专门仪器测试分析，获得了每个样本所含的多种化学元素的浓度数据。另一方面，按照2公里的间距在那些远离人群及工业活动的自然区取样，将其作为该城区表层土壤中元素的背景值。附件1列出了采样点的位置、海拔高度及其所属功能区等信息，附件2列出了8种主要重金属元素在采样点处的浓度，附件3列出了8种主要重金属元素的背景值。现要求通过数学建模完成以下任务：(1) 给出8种主要重金属元素在该城区的空间分布，并分析该城区内不同区域重金属的污染程度。(2) 通过数据分析，说明重金属污染的主要原因。(3) 分析重金属污染物的传播特征，由此建立模型，确定污染源的位置。(4) 分析你所建立模型的优缺点，为更好地研究城市地质环境的演变模式，还应收集什么信息？有了这些信息，如何建立模型解决问题？二 问题分析通过题目中所给的数据，建立浓度体系，即用坐标轴表示该城区的平面地形，做出每种重金属浓度的等高线图，从而直观的得出8种主要重金属元素在该城区的空间分布情况；运用单项指数法和内梅罗指数法对数据进行处理，可以得到8种主要元素在该城区不同区域重金属的污染程度。通过对问题1的结论进行分析，求出污染程度最为严重的元素，再根据元素找到相应的分区，最终找出重金属污染的主要原因。根据题中信息和实际情况，可以大致确定影响该地区污染程度的主要因素是浓度、海拔、轴和轴，因为不管重金属在土里还是在空中都要像分子扩散一样由高浓度向低浓度扩散，由于地球引力的影响，元素应从高海拔向低海拔扩散，不同功能区（采样点坐标）所排放的污染物不同。对四者综合考虑即是重金属污染物的传播特征。三 模型假设1. 采样时采集的是老成土，未采集填埋土和堆积时间不长的新成土；2. 在进行土壤地质调查过程中不考虑气候和人为因素，并且不发生地质活动和重大核事故；3. 不考虑农业活动对城区土壤中元素含量的影响；4. 空间分布图符合地图的方向原则（上北下南，左西右东）。四 符号说明:土壤中重金属元素的单项污染指数；:土壤中重金属元素的实测数据；:土壤中重金属元素的评价标准，即其背景值的平均值；:土壤综合污染指数；:第种元素在中的第个采样点的轴坐标；:第种元素在中第个采样点的轴坐标；:第种元素在中第个采样点的海拔；:第种元素对应函数的常数项；:第种元素在中第个采样点的浓度；:第种元素在中的采样点个数；五 模型的建立与求解5.1 问题1的求解5.1.1 8种主要重金属元素在该城区的空间分布用Matlab对附件1中各采样点的、坐标和附件2中各重金属元素在每一个采样点的浓度进行处理，得到8种元素浓度在该城区的空间分布（如图18），即在平面坐标系内用各采样点的、坐标表示所有采样点的空间分布，并用不同的图标加以区分各个分区的位置，在空间分布图中，最右侧的浓度梯度表示该元素城区的浓度水平，即图中的浓度等高线，由于浓度是非负的，所以不考虑负值浓度梯度；颜色越深表示元素浓度越高。图1 城区土壤中金属As的空间分布图2 城区土壤中金属Cd的空间分布图3 城区土壤中金属Cr的空间分布图4 城区土壤中金属Cu的空间分布图5 城区土壤中金属Hg的空间分布图6 城区土壤中金属Ni的空间分布图7 城区土壤中金属Pb的空间分布图8 城区土壤中金属Zn的空间分布该城区土壤中Cr和Ni 两种元素的含量分布规律比较相似（如图3、6），在城区西南位置出现一个较小峰值，并没有出现明显的富集，说明人类活动对Cr、Ni两元素的分布影响不大，另外它的含量和背景值非常接近，所以可以推断该城区土壤中的Ni元素主要是自然来源。 Cu和Pb 两种元素在该城区土壤中含量的分布比较相似（如图4、7），从整体上看在西南区域形成两个峰值，而两个峰值主要出现在生活区、工业区和主干道区。生活区Cu和Pb 两种元素含量高可能与生活垃圾的堆放有关，工业区Cu元素含量高可能与冶炼厂“三废”的排放有关，而主干道区Pb元素含量高则可能与含铅汽油的燃烧有关。 Hg和Zn 两种元素在该城区土壤中含量的分布比较相似（如图5、8），从整体上看在地理位置大致相同的地方形成3个峰值，而3 个峰值主要出现在生活区、工业区和主干道区。城市生活垃圾的堆放和工业生产废料的排放导致Hg和Zn 两种元素在生活区和工业区土壤中的含量增高，而含有添加剂Zn的机动车轮胎磨损产生的含Zn粉尘使Zn元素在主干道区土壤中含量上升。 As元素含量的分布在山区出现峰值（如图1），说明在山区的As元素的含量较高，而As元素含量较高可能与在山区内进行矿山冶炼有关。Cd元素在该城区的多处区域出现峰值（如图2），而峰值主要出现在生活区、工业区。燃煤、石油的燃烧和垃圾废弃物的焚烧导致Cd元素在生活区、工业区土壤中的含量增高。5.1.2 分析该城区内不同区域重金属的污染程度由于8种重金属元素可以看成独立的单项指标，从而可以通过单项污染指数法比较直观地反映出环境中各项污染指标的情况，即用公式 (1)对所有采样点的数据进行处理，其中当时表示土壤未受污染物的污染, 当时表示土壤受污染,越大受污染程度越严重，得到表格1如下：表格 1不同元素单项污染指标值区域As单项污染指数Cd单项污染指数Cr单项污染指数Cu单项污染指数Hg单项污染指数Ni单项污染指数Pb单项污染指数Zn单项污染指数11.74182.23052.22643.74272.65831.49122.22921.692922.01433.02391.72299.661818.3531.61073.00131.985231.12341.17171.25681.31191.17021.25641.17920.5235341.58562.76931.87274.713312.7661.43232.04951.734751.73992.1581.40762.28733.28551.24311.95831.1017通过表格1可以得到8种重金属元素分别在各个区域的单项污染指数，再通过内梅罗指数法对不同区域的单项污染指标进行综合计算，运用内梅罗指数法不仅考虑到了所有评价因子单项污染程度的平均水平，而且还考虑到了最大污染指数，因而最终能够科学、综合地反映评价区域内总体环境质量状况，即用公式 (2)得到各个区域的综合污染指数，通过Matlab编程，得到表格2；表格 2 不同区域的综合污染指数区域12345Pi3.088513.4831.22169.38222.6829根据GB15618-1995土壤环境质量标准2级标准值作为土壤污染评价标准并对土壤环境质量进行分级，得到图土壤污染评价分级标准，即表格3，如下：表格 3 土壤污染评价分级标准等级划分综合污染指数P污染等级污染水平P0.7安全清洁0.7P1警戒值尚清洁1P2轻度污染土壤轻度受污染23重度污染土壤重度受污染从表格3可以看出：3类区属于等级，5类区属于等级，1类区属于等级，1、2、4类区属于等级。综上可以得出该城区不同区域重金属的污染程度，得到以下结论：山区土壤环境质量是受到轻度污染的；公园绿地区土壤环境质量是受到中度污染的；生活区、工业区和主干道区土壤环境质量为重度污染，而它们的综合污染指数为分别为3.0885、13.483、9.3822，虽然从级别上看三个分区都是重度污染，但工业区的污染程度要大于主干道区，更远远大于生活区。 5.2 问题2的求解在对附件2所给数据进行处理的过程中，通过元素含量的变异系数来表示不同重金属元素在土壤中的含量变化程度，变异系数越大则该元素在土壤中的分布就越分散，其中元素含量变异值=(元素含量标准值/元素含量平均值)，通过对附件2进行数据计算后得到表格4，如下：表格 4城区土壤重金属元素含量统计表元素样本数（个）最大值（ug/g）最小值(ug/g)平均值(ug/g)标准差（ug/g）背景值平均值（ug/g）变异系数(%)As31930.131.615.683.020873.653.18Cd3191.61980.040.30240.224990.1374.4Cr319920.8415.3253.5170.001831130.82Cu3192528.482.2955.02162.91513.2296.11Hg319160.008570.299711.629540.035543.71Ni319142.54.2717.269.9414212.357.6Pb319472.4819.6861.7450.05783181.08Zn3193760.8232.86201.2339.23369168.6从表格4中可以看出，在该城区土壤中，8种重金属元素的平均含量都不同程度的超过了背景值，其中As、Ni、Cd、Pb的变异系数分别为53.18%、57.6%、74.4%、81.08%，在所有元素变异系数中相对较小；Cr、Zn、Cu的变异系数分别为130.82%、168.6%、296.11%，在所有元素变异系数中相对较大；而Hg的变异系数最大，为543.71%。反应了该城区土壤中As、Ni、Cd、Pb 4种元素的含量变化一般，Cr、Zn、Cu 3种元素的含量有较大差异，Hg元素的含量有很大的差异。通过以上分析可以看出，该城区土壤中As、Ni、Cd、Pb 4种元素污染一般、在土壤中不太分散；Cr、Zn、Cu 3种元素有较大程度的污染、在土壤中较分散；而Hg元素污染最严重、在土壤中最分散。即土壤中重金属污染程度和分散情况依次为：HgCuZnCrPbCdNiAs。在对元素含量进行变异系数进行处理后，进一步对单项污染指数进行变异系数处理，单项污染指数的变异系数表示该元素在城区的污染程度，变异系数越大则该元素在城区的污染程度越严重，通过对表格1进行整理计算后得到表格5，如下：表格 5 土壤重金属污染指数统计分析各元素单项污染指数P样本数（个）最大值最小值平均值标准差变异系数(%)As52.01431.12341.6410.32819.99Cd53.02391.17172.270680.71471.38Cr52.22641.25681.697280.38422.62Cu59.66181.31194.34343.24874.78Hg518.3531.17027.64667.52898.45Ni51.61071.24311.406740.15711.17Pb53.00131.17922.08530.65231.27Zn51.98520.523531.4076060.59141.99从表格5中可以看出，该城区土壤重金属单项污染指数值由大至小的顺序为Hg、Cu、Cd、Zn、Pb、Cr、As、Ni，同时根据土壤中重金属污染程度和分散情况可以得到土壤中不同重金属单项指数值变异系数规律与重金属污染程度规律基本相同，即导致城区土壤污染的重金属元素主要为Hg、Cu，而其中Hg元素在该城区土壤中污染最严重、分布最分散，为主要污染因子，有明显的外源污染的现象，Cu元素在该城区土壤中污染较严重、分布较分散，为次要污染因子，有较明显的外源污染的现象。结合相关文献可以知道土壤中Hg元素来自于工业生产废料的排放和城市生活垃圾的堆积，Cu元素的来自于含铜矿的开采和冶炼厂“三废”的排放，则该城区土壤重金属污染的原因是工业生产中“三废”的排放和城市生活垃圾的堆积。5.3对问题3的求解5.3.1筛选重污染区间首先对问题1所得的重金属元素在该城区的空间分布图进行分析，由于浓度分布杂乱不均匀，所以我们只对每一种重金属元素在该地区污染最严重的区域进行讨论，得到每一个重金属元素在该地区污染最严重区域的区间范围如表格6：表格6每一个重金属元素在该地区污染最严重区域的区间范围元素符号区间范围x轴y轴As (g/g)15000，230008000，12000Cd (ng/g)0，40002000，4000Cr (g/g)2000，50004000，7000Cu (g/g)0，40002500，4000Hg (ng/g)0，50002000，4000Ni (g/g)1000，50004000，7000Pb (g/g)2000，70004000，7000Zn (g/g)12000，150008000，12000然后给每一个重金属元素建立采样点集合，并根据每个重金属元素在该地区污染最严重区域的区间范围对原题采样点数据分别进行筛选，最后把每次筛选得到的采样点放入对应元素所建立的采样点集合中。5.3.2模型的建立8种重金属污染物的传播特征受浓度、海拔、轴和轴四个因素的影响，且每个因素对污染物的传播的影响程度不同，对此我们需要对离散点进行拟合，基于此我们综合考虑四个因素建立了多元线性回归模型： （3） 其中：，为与，相关的未知系数，称为回归系数；第个元素在中的各采样点的坐标矩阵如下：第个元素所对应的回归系数矩阵如下：第个元素在中的各采样点的浓度矩阵如下：正规方程组的矩阵形式为:当矩阵满秩时, 为可逆矩阵，可解出；把解出的各个重金属元素所对应的回归系数矩阵代入（3）式，得到各个重金属元素所对应的多元线性回归函数，然后把各个重金属元素所对应的中的各采样点的坐标代入，得到各个重金属元素在每个采样点的理论浓度值，用实际浓度值得到相对误差，用得到绝对误差。我们规定当有60%以上的采样点在内时，选用的系数符合多元线性回归模型的要求，能够用来解决具体问题，反之不能用来解决具体问题。5.3.3模型的解答根据上面对模型的分析，首先利用MATLAB（见程序3）得到集合和回归系数矩阵如下：接着根据系数矩阵得到各个重金属元素所对应的多元线性回归函数如下：把各个重金属元素所对应的中的各采样点的坐标代入，得到各个重金属元素在每个采样点的理论浓度值，再根据得到绝对误差，对进行分析大概有68%的采样点的绝对误差在内，说明该模型对于本问题可用。最后在浓度矩阵中选出每种重金属元素中浓度最高的浓度，将其代入对应的中，解出其对应的坐标如表格7：表格7污染源坐标元素XYZ浓度As10142166286.7152Cd1049212712934.5022Cr459246036552.3388Cu242739712447.4537Hg10492127125.32E+03Ni24865999254.5621Pb352643577183.1446Zn10398103600393.0065.4对问题4的求解对模型优缺点的分析如下：优点：操作简单，可行性好。对每个因素的影响曲线进行拟合得到综合影响曲线，这样能更全面更实际反应出传播特征曲线。缺点：具有一定的局限性，可能存在其他因素对传播特征的影响。 选取高浓度区间较少，得到的污染源较少。 模型得到的理论结果与实际结果的误差较大。为更好地研究城市地质环境的演变模式，经分析，影响金属污染物的传播特征的因素还存在自然降解力、气候和人为因素等。所以还应收集各个功能区的植被密度、植被类型和土壤中微生物类型及该地区的气候特点人类活动特点等信息。首先对属于同一类因素的不同因素进行权重处理并经综合加权得到该因素影响金属污染物的传播特征的能力，然后考虑所有因素建立多元线性回归模型，最后对得到的结果应用单层线性神经网络模型中的最小均方差LMS学习算法：使误差达到预选设定要求，由解出的多元线性回归方程解出污染源位置。六 模型的推广在本题中多元线性回归模型解决了重金属污染物的传播特征问题，在现实生活中该模型可以应用到分子扩散特征、病毒传播特征、人体温度分布特征和信号传播特征等问题。七 参考文献1（美）Gerald Recktenwald,数值方法和MATLAB实现与应用，北京：机械工业出版社，2004年9月；2 王济，王世杰，土壤中重金属环境污染元素的来源及作物效应，贵州师范大学学报，第23卷第2期：1-8，2005年；3 李玲，吴克宁，张雷，吕巧灵，郑州市郊区土壤重金属污染评价分析，土壤通报，第39卷第5期：1-5,2008年；4 杜廷松，沈艳军，覃太贵，数值分析及实验，北京：科学出版社，2006年12月；5 神经网络与Matlab仿真, /view/511c0886b9d528ea81c7795f.html, 2011年9月11号；6 多元线性回归模型,/view/1586224.htm, 2011年9月11号。八 附录程序1：%第一问clear allclcA=xlsread(cumcm2011.xls,1,B4:D322);B=xlsread(cumcm2011.xls,1,F4:M322); %导入excel文件 x=A(:,1)y=A(:,2)for i=1:8; z=B(:,i); figure X,Y,Z=griddata(x,y,z,linspace(min(x),max(x),linspace(min(y),max(y),v4); surf(X,Y,Z) contourf(X,Y,Z,15) hold onx1=A(1:44,1);y1=A(1:44,2);plot(x1,y1,k*);hold onx2=A(45:80,1);y2=A(45:80,2);plot(x2,y2 ,c.);hold onx3=A(81:146,1);y3=A(81:146,2);plot(x3,y3,b+);hold onx4=A(147:284,1);y4=A(147:284,2);plot(x4,y4,m);hold onx5=A(285:319,1);y5=A(285:319,2);plot(x5,y5,r);enda=rand(5,8);for i=1:8; a(1,i)=sum(B(1:44,i)/44;end %1区的八种元素的平均值for i=1:8; a(2,i)=sum(B(45:80,i)/36;end %2区的八种元素的平均值for i=1:8; a(3,i)=sum(B(81:146,i)/66;end %3区的八种元素的平均值for i=1:8; a(4,i)=sum(B(147:284,i)/138;end %4区的八种元素的平均值for i=1:8; a(5,i)=sum(B(285:319,i)/35;end %5区的八种元素的平均值a;a(:,1)=a(:,1)/3.6;a(:,2)=a(:,2)/130;a(:,3)=a(:,3)/31;a(:,4)=a(:,4)/13.2;a(:,5)=a(:,5)/35;a(:,6)=a(:,6)/12.3;a(:,7)=a(:,7)/31;a(:,8)=a(:,8)/140; a=a; %a 每个区的各元素平均值的标准后值b=max(a);ave=sum(a)/8;c=rand(1,5);for k=1:5; c(1,k)=sqrt(b(1,k).2+ave(1,k).2)/2); % 内梅罗综合指数end %对每个浓度标准化C=;C(:,1)=B(:,1)/3.6;C(:,2)=B(:,2)/130;C(:,3)=B(:,3)/31;C(:,4)=B(:,4)/13.2;C(:,5)=B(:,5)/35;C(:,6)=B(:,6)/12.3;C(:,7)=B(:,7)/31;C(:,8)=B(:,8)/140;程序2：%第二问的单值指数程序clear allclcA=xlsread(cumcm2011.xls,1,B4:D322);B=xlsread(cumcm2011.xls,1,F4:M322); %导入excel文件 D=B;D(:,1)=B(:,1)/3.6;D(:,2)=B(:,2)/130;D(:,3)=B(:,3)/31;D(:,4)=B(:,4)/13.2;D(:,5)=B(:,5)/35;D(:,6)=B(:,6)/12.3;D(:,7)=B(:,7)/31;D(:,8)=B(:,8)/140;d=;e1=D(1:44,:)1;d(1,:)=sum(e1);e2=D(45:80,:)1;d(2,:)=sum(e2);e3=D(81:146,:)1;d(3,:)=sum(e3);e4=D(147:284,:)1;d(4,:)=sum(e4);e5=D(285:319,:)1;d(5,:)=sum(e5);程序3：%第三问的金属污染源确定的程序clear allclcA=xlsread(cumcm2011.xls,1,B4:D322);B=xlsread(cumcm2011.xls,1,F4:M322);%1 AsD1=;e1=;j=1;for i=1:319; if A(i,1)10000&A(i,1)1000&A(i,2)6000 D1(j,:)=A(i,:); % 选定的点的坐标矩阵 e1(j,1)=B(i,1); % 浓度矩阵 j=j+1; %As选取的点的各值 endendk=j-1;D=D1,ones(k,1); % 坐标矩阵b1=inv(D*D)*D*e1; %As选取点的系数B1=D*b1; % B1 为确定参数回带x y z所确定的浓度v=abs(e1-B1);for i=1:k; v1(i,1)=v(i,1)/e1(i,1); %相对误差endl1=(sum(v1)-max(v1)/(k-1); % As的浓度高峰处总体误差a1=v10&A(i,1)2000&A(i,2)4000 D2(j,:)=A(i,:); e2(j,1)=B(i,2); j=j+1; %Cd选取的点的各值 endendk=j-1;D=D2,ones(k,1);b2=inv(D*D)*D*e2; %Cd选取点的系数B2=D*b2;v=abs(e2-B2);for i=1:k; v2(i,1)=v(i,1)/e2(i,1); %相对误差endl2=(sum(v2)-max(v2)/(k-1); % Cr的浓度高峰处总体误差a2=v22000&A(i,1)4000&A(i,2)7000 D3(j,:)=A(i,:); e3(j,1)=B(i,3); j=j+1; %Cr选取的点的各值 endendk=j-1;D=D3,ones(k,1);b3=inv(D*D)*D*e3; %Cr选取点的系数B3=D*b3;v=abs(e3-B3); for i=1:k; v3(i,1)=v(i,1)/e3(i,1); %相对误差endl3=(sum(v3)-max(v3)/(k-1); % Cr的浓度高峰处总体误差a3=v30&A(i,1)0&A(i,2)5000 D4(j,:)=A(i,:); e4(j,1)=B(i,4); j=j+1; %Cu选取的点的各值 endendk=j-1;D=D4,ones(k,1);b4=inv(D*D)*D*e4; %Cu选取点的系数B4=D*b4; % B1 为确定参数