高教社杯全国大学生数学建模竞赛A题城市表层土壤重金属污染分析

2011高教社杯全国大学生数学建模竞赛承诺书我们仔细阅读了中国大学生数学建模竞赛的竞赛规则.我们完全明白，在竞赛开始后参赛队员不能以任何方式（包括电话、电子邮件、网上咨询等）与队外的任何人（包括指导教师）研究、讨论与赛题有关的问题。我们知道，抄袭别人的成果是违反竞赛规则的,如果引用别人的成果或其他公开的资料（包括网上查到的资料），必须按照规定的参考文献的表述方式在正文引用处和参考文献中明确列出。我们郑重承诺，严格遵守竞赛规则，以保证竞赛的公正、公平性。如有违反竞赛规则的行为，我们将受到严肃处理。我们参赛选择的题号是（从A/B/C/D中选择一项填写）：A 我们的参赛报名号为（如果赛区设置报名号的话）：所属学校（请填写完整的全名）：浙江理工大学参赛队员(打印并签名)：1.冯世扣2.李漪3.陆天栋指导教师或指导教师组负责人(打印并签名)：数模组日期：2011年9赛区评阅编号（由赛区组委会评阅前进行编号）：2011高教社杯全国大学生数学建模竞赛 编号专用页 赛区评阅编号（由赛区组委会评阅前进行编号）：赛区评阅记录（可供赛区评阅时使用）：评阅人评分备注全国统一编号（由赛区组委会送交全国前编号）：全国评阅编号（由全国组委会评阅前进行编号）：城市表层土壤重金属污染分析摘要随着城市经济的快速发展和城市人口的不断增加，人类活动对城市环境质量的影响日显突出，因此研究城市地质环境的变化十分重要。对于问题一，我们首先对数据进行预处理，然后利用地统计学软件surfer运用反距离加权插值法（InverseDistancetoaPower）做出8种重金属的空间等值线分布图，紧接着我们用C语言编程算出五个功能区中各种重金属含量的平均值。在分析该城区内不同区域重金属的污染程度时，我们利用已经统计的结果采用了单因子指数法，内梅罗综合污染指数法以及均方根综合污染指数法进行评价。用单因子评价方法得出了，该城区和Hg污染较为严重。用内梅罗综合污染指数法以及均方根综合污染指数法得出五个功能区污染程度依次为工业区、主干道路区、生活区、公园绿地区和山区，其中工业区和主干道路区重度污染，生活区属于中度污染，公园绿地区中度污染，山区为轻微污染。对问题二，根据计算出来的单因子数值以及8张重金属空间等值线分布图，发现图中有较多的地方含有高浓度的，且这些区域较为分散。这些区域可能是主干道路区，主要是因为汽车尾气中含有大量的。其他7张重金属空间等值线分布图分别只有几个的区域含有高浓度的，Hg，，，，，，且这些区域集中是在工业区。因为工业区排出来的“三废”经过大气和水体最终污染土壤，使土壤中重金属的浓度增加。问题三我们把主成份分析法和扩散问题的偏微分方程结合来来建立模型。我们先用主成份分析法得到各重金属间的相关性，并画出主因子的空间等值线分布图，从而先确定污染源的大概位置。然后再建立扩散问题的抛物型偏微分方程模型，为了确定源位置，通过微分方程的解析解和测量数据的关系，构造了一个目标函数优化问题，通过使用源附近浓度最大的点作为基准点，巧妙地把非线性问题转化为线性问题，用线性回归确定了解析解的参数，并用SPSS软件求解污染源的准确位置。得出和来自同一个污染源（3863，7725，36），和来自共同污染源（2029，2717，34）（21524，11778，47），的污染源在（2016，5221，18）。的三个污染源分别在（6559，6947，34）（1356，5201，54）（17770，7515，54），的污染源（11209，9257，17）Hg污染源在（2894，2345，21）。对于问题四由于扩散问题的偏微分方程模型是一个静态的模型，因为地下水的作用，重金属在土壤里的传播是一个流动的模型。这时就要引入对流项，此时就会变成一个非常系数的偏微分方程，这种方程求解十分的复杂，可以采用反演的算法来求解。模型的创新性在于引入一个基准点，把非线性的偏微分方程转化为线性的偏微分方程。但是此模型是通过SPSS内部的拟合的函数求得污染源位置，因此会存在一定的误差。关键词：重金属土壤污染内梅罗法主成分分析法扩散问题的偏微分方程一问题重述1.1基本情况为了研究城市表层土壤重金属污染程度，现对某城区地质环境进行勘察调查。按照功能划分，城区一般可分为生活区、工业区、山区、主干道路区及公园绿地区等，分别记为1类区、2类区、……、5类区。我们将城区划分为间距1公里左右的网格子区域，按照每平方公里1个采样点对表层土（0～10厘米深度）进行取样、编号，并用GPS记录采样点的位置。应用专门仪器测试分析，获得了每个样本所含的多种化学元素的浓度数据。另一方面，按照2公里的间距在那些远离人群及工业活动的自然区取样，将其作为该城区表层土壤中元素的背景值。1.2需要解决的问题(1)给出8种主要重金属元素在该城区的空间分布，并分析该城区内不同区域重金属的污染程度。(2)通过数据分析，说明重金属污染的主要原因(3)分析重金属污染物的传播特征，由此建立模型，确定污染源的位置。(4)分析你所建立模型的优缺点，为更好地研究城市地质环境的演变模式，还应收集什么信息？有了这些信息，如何建立模型解决问题？二符号说明及基本假设2.1符号说明（1）是指土壤中污染物的环境质量指数；（2）是指污染物的实测值，单位为；（3）是指污染物的平均值，单位为；（4）则为采样点处的浓度，为点源G的质量；（5），，，为重金属扩散速度的算数平方根；（6）k为植物等对重金属衰减系数的算术平方根；（7）表示采样点的位置；（8）表示拟合出来的线性方程2.2基本假设（1）平稳性假设，即假设在每个采样点区域范围内重金属污染物的浓度是均匀的，不考虑因某种特殊原因下某一点的某种重金属含量显著高于周围的情况，如土壤中埋入了一节废旧电池等。（2）假设同一种重金属在不同的功能区的土壤传播特征是一样的（3）不考虑其他重金属对该城区的污染（4）不考虑8种重金属之间的化学物理反应（5）假设附件中给出的数据都是准确的三问题分析3.1问题1的分析为了给出8种主要重金属元素在该城区的空间分布并且分析该城区内不同区域重金属的污染程度。题目中附录一和附录二提供了大量的数据，我们的任务是从这些庞大的调查数据中找出重金属元素在不同区域的变化规律，这可以归结为一个数理统计的问题。只要我们把各个重金属的位置和测量值统一标在同一幅图上，就可以得出重金属在该城区的空间分布。对于分析不同城区重金属的污染程度，我们引进了新的评判指标——均方根综合污染指数法和内梅罗综合污染指数法，并且用编程求出城区内不同区域的均方根综合污染指数法和内梅罗综合污染指数，然后根据各自方法的评判指标，得出该城区内不同功能区重金属的污染程度。3.2问题2的分析对于问题2，要分析重金属的污染的主要原因，我们要结合重金属在该城区的空间分布以及内梅罗综合污染指数法和均方根综合污染指数法的结果，再查阅这几年重金属污染的相关文献，找到了污染的主要原因3.3问题3的分析要根据附录中数据，确定污染源的位置。重金属在土壤中的传播属于一个典型的扩散问题，但是由于该城区重金属的采样点较为的多，给我们用扩散问题的偏微分方程模型来寻找污染源造成了很大的困难。因此我们需要缩小污染源的范围。另外由于不同重金属之间可能存在相关性，即不同重金属可能来自同一个污染源，因此我们想到了用主成份分析法。因为主成份分析法不但可以减少变量的个数，大大的减少了计算量，还可以画出主因子的空间分布图，这样我们就可以找出重金属较为集中的地方，也就是能够寻到污染源的大概位置。在确定污染源大概位置后，我们就可以用扩散问题的偏微分方程模型的方法来求出具体污染源的具体位置。四模型的建立与求解4.1基于内梅罗综合污染指数法的模型建立要建立8种重金属在区域的空间分布，我们采用用反距离加权插值法（InverseDistancetoaPower）画出每种重金属的空间等值分布图。在分析该城区内不同区域重金属的污染程度过程中，我们先建立了内梅罗综合污染指数法评判模型。其中，单因子分指数的计算公式为：式中是指土壤中污染物的环境质量指数；是指污染物的实测值，单位为；是指污染物的平均值，单位为。内梅罗综合污染指数的计算公式为：式中，为土壤污染中污染指数的最大值；是指土壤污染中污染指数的平均值。我们再建立均方根综合污染指数法，其公式为：式中，P表示污染综合指数；表示参加评价的污染物种类；表示污染物的实测值，；污染物的评价标准，；利用上述的公式算出各自的内梅罗综合指数和均方根指数后，再根据污染划分指标，分析该城区内不同区域重金属的污染程度。4.2内梅罗综合污染指数法评判模型求解4.2.1对8种主要重金属元素在该城区的空间分布情况的分析图4-1的空间等值线分布图图4-2的空间等值线分布图从图中可以看出，只有五个区域的浓度较为集中，其他地方的浓度都是较低的。因此可能是这个区域是加工工业区。对于表，图中的有较多的地方含有高浓度的，且这些区域较为分散，这些可能是主干道路区。图4-3的空间等值线分布图图4-4的空间等值线分布图从图中可以看出，只有两个区域和的浓度较为集中，其他地方的浓度都是较低的。因此可能是这个区域是加工和工业区图4-5Hg的空间等值线分布图图4-6的空间等值线分布图从图中可以看出，只有三个区域Hg的浓度较为集中，其他地方Hg的浓度都是较低的。因此可能是这个区域是Hg加工工业区。对于表，图中的只有两个区域的浓度较为集中，其他地方的浓度都是较低的，因此可能是这个区域是加工工业区图4-7的空间等值线分布图图4-8的空间等值线分布图对于表，图中和的有较多的地方含有高浓度的，且这些区域较为分散，这些区域可能是主干道路区。4.2.2内梅罗法和均方根法分析污染程度在模型的求解过程中，我们先将五个功能区的各种重金属的数值带入单因子分指数的计算公式：和内梅罗公式：算出各个重金属在不同功能区的单因子指数值和五个功能区的内梅罗指数，数据见下表：表4-1各种金属污染物在不同功能区的单因子指数值AsCdCrCuHgNiPbZn内梅罗指数生活区1.742.661.492.233.433.17工业区2.003.021.729.6618.351.613.004.0313.533山区61.311.171.261.181.061.253主干道路区1.592.771.874.7112.771.432.053.529.427公园绿地区1.742.161.412.293.291.241.962.242.734表4-2土壤单因子指数分级标准污染等级未污染轻度污染中度污染重度污染极度污染研究单因子指数值可以发现，山区的污染最小，所有重金属的单因子指数均在1和2之间，都属于轻度污染。研究地土壤污染以Hg污染最为严重，其中在工业区和主干道路区Hg的单因子指数分别达到了18.35和12.77，都是属于重度污染。生活区和公园绿地区分别则属于中度污染，只有山地是轻度污染。五个功能区的和的单因子数值均在1到2之间，都是属于轻度污染。五个功能区的除了生活区是中度污染外，其他地区均小于2，都属于轻度污染。工业区的的单因子指数为9.66，已经属于极度污染，主干道路区属于重度污染，生活区和公园绿地区为中度污染。工业区的的单因子指数为4.03,属于重度污染，生活区，主干道路区和公寓绿地区均为中度污染。表4-3土壤综合污染指数分级标准污染等级污染指数污染等级污染水平1安全清洁2警戒线尚清洁3轻污染土壤轻污染农作物开始受到污染4中污染土壤作物均受中度污染5重污染土壤作物均受污染且相当严重通过内梅罗指数我们可以很清楚的看到，工业区的污染是最严重的，内梅罗指数为13.533，属于重污染。主干道路区内梅罗指数为9.427，也是属于重污染。生活区的内梅罗指数为3.170，也属于重污染。公园绿地区的内梅罗指数为2.734，属于中污染。山地的污染最小，内梅罗指数为1.253，属于轻污染。然后我们将五个功能区的各种重金属的数值带入均方根公式：得出五个功能区的均方根指数，如下表:表4-4各功能区的均方根指数表功能区生活区工业区山区主干道路区公园绿地区均方根指数2.5747.6981.1945.2142.123表4-5土壤综合污染指数分级标准污染等级污染指数污染等级污染水平1未污染清洁2轻微污染尚清洁3偏中度污染土壤轻污染农作物开始受到污染4中度污染土壤作物均受中度污染5偏重度污染土壤作物均受污染且相当严重6重污染土壤作物均受污染且相当严重7极度污染土壤作物均受污染且相当严重通过均方根指数我们可以很清楚的看到，工业区的污染是最严重的，均方根指数为9.698，属于极度污染。主干道路区均方根指数为5.214，也是属于重污染。生活区的均方根指数为2.574，也属于偏重度污染。公园绿地区的均方根指数为2.734，属于偏重度污染。山地的污染最小，均方根指数为1.194，属于轻微污染。4.3分析重金属污染的主要原因土壤中的重金属主要来源于成土母质和人类活动[1]。根据题目给出的土壤各重金属的背景值和各个监测点的实测值，再结合我们计算出来的单因子数值，我们可以判定，人类的活动是造成这些功能区重金属超标的主要原因。这五个功能区污染程度依次为:工业区﹥主干道路区﹥生活区﹥公园绿地区﹥山区。工业区中，Hg，，，已经严重超标，另外三种也超标，主要由于工业区有较多的金属冶炼厂，蓄电池加工厂，印刷厂和化工厂等，这些工厂都需要大量的金属原材料，煤炭和石油等，这些都会产生大量的重金属废渣废气。这些废渣废气排到土壤，河流和大气里，会造成工业区土壤重金属含量严重超标，另外这些废气废水会随着大气水体影响附近的地区的土壤。比如蓄电池加工厂周围的土壤中的含量会明显超标。不同的金属冶炼厂对各重金属污染在周边土壤中的积累贡献各异，如冶炼厂附近土壤中不仅的含量会增加，同时的含量也会有所增加[4]。主干道路区的,,,污染也较为严重，主要是由于机动车会产生大量的尾气。汽车尾气中含有大量的，，，等重金属元素，这些重金属会严重影响主干道路区的土壤。通过对汽车尾气颗粒物中重金属含量的测定发现，的含量为37%，，，的含量分别为34.5%，22.6%，2.6%[2]。土壤中的元素是不能随着时间被土壤降解的，这就导致元素会随着时间慢慢积累。另外元素是汽车硬度的添加剂，在汽车轮胎摩擦产生的粉尘会引起土壤，含量的增多。生活区的污染主要是由于生活垃圾和生活污水的排放以及生活用煤。垃圾堆积在土地上，垃圾中的重金属会渗透到土壤里，这就会导致土壤重金属含量的增加。生活污水的排放，重金属就会流入河流，最终影响生活区的土壤和其他地区的土壤。居民用煤也是产生重金属的重要原因。煤炭燃烧过程中产生的易于富集在干灰中，且和更易于富集于细灰粒中，然后随着细灰粒排放到大气中，这些重金属就会进入大气，最终沉降到土壤中。公园绿地区的污染较轻，其中，Hg，都是中度超标。这有可能是受工业区的废气和汽车尾气影响。因为大气中的重金属经过自然沉淀和雨淋会进入土壤中，这就有可能影响到公园绿地区的土壤。另外旅客产生的垃圾也会影响公园绿地去的土壤。山区的污染在这几个功能区里面是最轻的。山区的重金属来源可能来自大气，附近的河流和农药的使用。但是植物能够降解重金属，这会使得山区的土壤重金属含量大大的降低。4.4基于扩散问题的偏微分方程建立首先我们先建立主成分分析法模型，来找出污染源的大概位置。第一步，我们先对数据进行预处理，把数据根据功能区重新排列，依次为1区，2区，3区，4区，5区。第二步我们再利用软件计算出标准化数据的相关系数矩阵，并且算出相关系数矩阵的特征值，特征向量，方差贡献率，累计方差贡献率。正交的方差最大旋转法来求解旋转后的特征向量，方差贡献率，累计方差贡献率。通过旋转后的累计方差贡献率，我们就可以确定哪几个因子为主因子，哪几个因子为次因子。第三步我们再进行正交互换，使用方差最大法，求出旋转因子的载荷矩阵。第四步我们通过分析旋转因子的载荷矩阵，就可以知道哪些重金属为主因子，哪些重金属属于次因子。为了更加精确的确定污染源的位置，我们做出了主因子在5个功能区的空间分布等直线图，这样能够更加直观的看出污染源大概在空间上的分布情况。然后我们就可以建立了扩散问题的偏微分方程模型，即（1）其中点G（，，）为假设的污染源扩散点，则为G点处的浓度，，，，为重金属扩散速度的算数平方根，k为植物等对重金属衰减系数的算术平方根，为点源G的质量。（，，）为取样点的位置。4.5扩散问题的偏微分模型的求解4.5在用主成份分析法的过程中，我们先用SPSS软件得出了各种重金属之间的相关性。相关性分析就是研究两个或两个以上变量之间的相互关系的统计分析方法。在环境研究中常用来定性的分析研究各环境变量之间的相互关系，而变量之间的相互关系用相关系数和相关系数的显著水平来表示[3]。如果相关系数为正，表示两个变量正相关，取值越小相关性越强。相反，如果相关系数为负值，表示两个变量负相关，值越小则表示相关性越强。土壤中重金属含量的相关性可以推测出重金属是不是来自同一个污染源[4]。如果重金属含量有着很强的相关性，表明重金属的来源可能相同，如果重金属的相关性不是很强，则重金属的来源可能不止一个。下表列出了各种重金属之间的相关性：表4-6重金属的相关性表AsCdCrCuHgNiPbZnAs10.2550.1890.160.0640.3170.290.247Cd0.25510.3520.3970.2650.3290.660.431Cr0.1890.35210.5320.1030.7160.3830.424Cu0.160.3970.53210.4170.4950.520.387Hg0.0640.2650.1030.41710.1030.2980.196Ni0.3170.3290.7160.4950.10310.3070.436Pb0.290.660.3830.520.2980.30710.494Zn0.2470.4310.4240.3870.1960.4360.4941从图中我们可以看出各个重金属之间的相关性不是很强。和的相关度为0.532，为中度极显著正相关，和的相关度为0.716，为高度极显著正相关，和均为中度极显著正相关，和的相关度为0.660，为中度极显著正相关，和的相关度为0.494，也为中度极显著正相关。4.5由于该城区有八种重金属的污染，而且有些金属之间还存在一定的相关性，这就会导致放映的信息重叠。因此我们要降低因子数。而因子分析则是通过降维方法进行简化得到综合指标，综合指标之间既互不相关，又能反映原来的观察指标的信息。我们用SPSS软件算出了特征值，贡献率和累计的贡献率。但是主成份分析法只强调离差而不是强调变量在样品中的百分含量，因此我们用正交的方差最大旋转法使每一个主因子只与最少个数的变量有相关的关系[6]。算出来的旋转变化后的特征值，贡献率和累计贡献率，如表所示：表4-7特征值和累计贡献率表成份旋转变化前旋转变化后特征值贡献率%累计贡献率%特征值贡献率%累积贡献率%13.5644.544.51.82722.83422.83421.1514.37758.8771.52519.06341.89730.96512.06370.9411.04513.06854.96540.7689.59680.5371.03212.89867.86350.5787.2287.7561.0212.74480.60760.4325.39993.1561.00412.54993.156从表格里可以看出，累计贡献率为93.156%的前提下，我们得到了6个主因子，也就是意味着6个主因子提供了93.156%的信息，满足主成份分析法的条件。这样大大地简化问题，且能偶提供93.156%的信息，有较高的准确性。另外我们发现在旋转前后累计贡献率都是为93.156%，表明总量信息没有丢失，我们的模型具有较高的准确性。

表方差极大正交旋转前的因子载荷矩阵。对旋转后的数据进行更加相信的分析，我们发现，1和2的贡献率大约在20%左右，这可以说明这两个因子可能是这个地区的最大的污染源，对这个地区的环境影响最大。而3，4，5，6四个因子的贡献率大概在13%左右，表明这四个因子对这个地区的污染没有1，2两个因子那么重要，但是也是也对土壤有着较大的影响。表4-8旋转前因子载荷矩阵123456As0.426-0.20.6810.551-0.026-0.065Cd0.7110.2810.282-0.322-0.2540.325Cr0.735-0.444-0.303-0.046-0.110.098Cu0.7560.125-0.3650.137-0.155-0.408Hg0.4080.673-0.2970.4490.1540.236Ni0.723-0.515-0.190.137-0.0140.2Pb0.7640.3140.237-0.248-0.158-0.217Zn0.699-0.0370.123-0.2410.654-0.06表4-9旋转后的因子载荷矩阵123456As0.1350.130.0160.9740.0430.084Cd0.2230.9180.1480.0860.0310.126Cr0.8590.187-0.0160.0040.2450.141Cu0.3950.1770.2730.0310.810.11Hg0.0170.1370.9670.0160.1730.069Ni0.890.1010.0510.1960.1030.17Pb0.0730.7180.0610.1620.4960.275Zn0.2610.2350.0820.0970.1220.916结合上面两幅表格一起分析可以得出：旋转前后因子载荷矩阵数值较为接近。然后我们继续分析旋转后的因子载荷矩阵，发现重金属与因子的载荷值越大，则表示重金属和因子的关系越密切。例如和与因子1的载荷分别为0.859和0.890。我们就认为因子1为和的组合。以此类推，因子2为和的组合，因子3为Hg，因子4为，因子5为，因子6为。因此和，和极有可能来自同一污染源。在前面的重金属相关性分析中，和，和都具有较好的。另外的重金属则来自不同的污染源。4.5.3首先我们运用傅立叶变换出求出偏微分方程的通解，求出了的通解为（2）如果认为经过相当长的时间后，扩散已经终止，物质分布处于平稳状态，则方程1中，于是方程1就可以有线性椭圆型方程的边值问题其中是区域D的边界。由于点源G的质量与位置，（，，）是已知的，前面我们用主成份分析法已经找出来的大概污染源的位置。为了确定方程（2）中的7个参数，我们给出下面一个目标函数：取样本点，使得最小,考虑到解的表达式中含有指数取对数，得：其中包括了项，是非线性的。所以我们引入一个浓度最大的点作为基准点，这样可以消掉一个参数M，并把非线性问题转换为线性问题，即求下面函数的最小值：现在我们先在大概污染源的中心取一点H，然后在大概污染源附近取一系列测量点：···。令t=1,把一系列的测量点和H点带入偏微分方程的通解中，再把测量点的物质密度与通解公式都取对数，得出两条等式，两条等式相减，则其中系数，，。我们把点源，点H，一些列的测量点···进行多元一次线性方程拟合，可以拟合出一条依次线性方程。然后根据系数就可以求出，，的值，即确定了污染源的位置。表4-10污染源的位置坐标位置AsCd和PbCr和NiCuHgZnx6559.792029.583863.542016289411209.1y6947.712717.967725.265221.1823459257.45z34.7384.1736.19182117.02x13560.5321524.8y5201.2111778.28z5447.64x17770.5y7515.97z54.284.6问题4求解4.6.1（1）我们先用主成份分析法确定污染源的大概位置后，在用扩散问题的偏微分方程模型来求解，这样大大的减少了建模过程中的计算量，并且能过很快的算出污染源的准确位置。（2）此模型还考虑了重金属在土壤中的衰减，比较符合实际，具有较好的有科学性。4.6.2（1）在建立模型三的过程中，我们考虑的是粉尘模型，没有考虑水流的作用，这样建立的数学模型与实际可能有一定的误差。（2）由于本模型在求解过程中，通过SPSS软件来拟合大量数据，从而求出污染源的位置。用拟合算法会存在一定的误差。4.6.3在知道了模型三的缺点后，我们需要采集相关的信息加以完善。比如采集多个不同时刻重金属在相同样本点的浓度。在模型中加上对流系数后，模型就不再是重金属在的静止的土壤中传播了，而是在相对流动的土壤中传播，这样更加的符合实际，也能更加准确的求出污染源的位置。模型三除了需要考虑土壤的对流系数，还可以考虑重金属在空间的变异性和自流系数。我们考虑到地下水流的影响，此时原函数就变为：我们甚至还应该考虑到重金属在各个方向上扩散系数的变化，即a、b、c是关于时间的函数。此时该方程是个非常难解的非常系数的偏微分方程，很难把解析解求出来，已知样本点的一些测量数据要求出源位置，这是偏微分方程反问题研究的范畴。已知源项和边值条件，由初边值条件求解即是人们熟知的正问题。而源项反问题的提法是，已知在某些点的测量数据或某个时刻T的测量数据反求源项，且该反问题在一般意义下是不适定的。由于该偏微分方程反问题含有很多参数，直接得到解析解是很困难的，必须使用数值微分的方法做反演计算（见[8,9]）。五模型的结果分析及评价5.1内梅罗综合污染指数法评判模型的结果分析对于模型一，我们用内梅罗综合污染指数法和均方根综合指数法来建立模型，虽然用内梅罗法求出来的数值稍大于均方根法求出来的数值，但是整体上的污染等级还是较为相近的。用两种方法建立模型求解的结果中，污染程度都是工业区>主干道路区>生活区>公园绿地区>山区。其中工业区和主干道区都是属于重度污染！因此两个模型都具有相当好的可行性。5.2扩散问题的偏微分方程模型的结果分析我们用模型3求解出来的污染源的位置和第一问中的重金属的空间等值线分布图进行了比较，发现污染源附近重金属浓度相当的高。这也说明模型三具有一定的准确性。5.3模型的优点（1）内梅罗综合污染指数法操作简单，运算也较为方便，只需要算出内梅罗综合污染指数，然后再跟评判标准进行比较，就可以得出地区的污染程度。另外内梅罗综合污染指数在加权过程中避免了权系数中主观因素的影响，兼顾了单因子污染指数的平均值和最高值。具有一定的可靠性。（2）模型三的优点在问题四中已有详尽的说明。5.4模型的缺点（1）内梅罗综合污染法在计算过程中，过分的看重了极大值对水体的污染，如果有一项指标的值很大，而其他值较小，这样也会使内梅罗综合指数偏大[7]。另外内梅罗综合污染法还忽略了重金属的毒性的差别。因此得到的结果往往比较机械。（2）模型三的缺点在问题四中已有详尽的说明。在改进内梅罗综合污染法的过程中，我们还考虑了不同重金属的毒性对的影响，得到改进后的内梅罗综合污染指数法，其公式为：其中为权重值最大的污染因子的与的比值。表示第种污染因子的相关性比值。令，其中为第种污染因子的权重值，。改进之后的内梅罗综合污染法一定程度上消除了极大值的影响。并提高了重金属毒性指标的影响，能够更加客观的评价地区的污染状况。参考文献[1]张嘉颖，大宝山矿区横石河流域重金属污染及其生态评价3/WFknowledgeServer_Mirror/D/Thesis_Y1476990.aspx，2011.9.9[2]张志红，杨文敏，汽油车排除颗粒物化学组分分析。中国公共卫生。2001，17（7）：623-624[3]林艳，基于地统计学与GIS的土壤重金属污染评价与预测，/p-68758645.html，2011.9.10[4]周燕，崇明岛土壤重金属污染的空间分析及污染评价，/kcms/detail/detail.aspx?dbCode=CMFD&dbName=CMFD2009&FileName=2009187292.nh&filetitle=，2011.9.10[5]张乃明，李保国，胡克林，太原污灌区土壤重金属与盐分含量的空间分布特征【J】．环境科学学报，21(3)：349～353，2001[6]王雄军，赖健青，鲁艳红，李德胜，周继红，王建武，基于因子分析法研究太原市土壤重金属污染的主要来源，17（2）：671-676，2008[7]李亚松，张兆吉，费宇红，王昭，内梅罗指数评价法的修正及其应用，25（6）：48-50，2009[8]YAMAMOTOM.Conditionalstabilityindeterminationofforcetermsofheatequationsinarectangle[J].Math.Comput.Modelling,18:79-88,1993[9]CHOULLIM,YAMAMOTOM.Conditionalstabilityindeterminingaheatsource[J].JournalofInverseandIll-PosedProblems,12:233-336,2004附录求解内梅罗指数的程序#include<stdio.h>#include<math.h>intmain(){doublep[5][8];doublec[5][8]={6.27,289.96,69.02,49.4,93.04,18.34,69.11,237.01, 7.25,393.11,53.41,127.54,642.36,19.81,93.04,277.93, 4.04,152.32,38.96,17.32,40.96,15.45,36.56,73.29,5.71,360.01,58.05,62.21,446.82,17.62,63.53,242.85,6.26,280.54,43.64,30.19,114.99,15.29,60.71,154.24};doubles[8]={3.6,130,31,13.2,35,12.3,31,69};inti,j;for(i=0;i<5;i++) for(j=0;j<8;j++) p[i][j]=c[i][j]/s[j];doublesum,ave,max,temp;doubleans[5];for(i=0;i<5;i++){sum=p[i][0];max=p[i][0];for(j=1;j<8;j++) { if(p[i][j]>max) max=p[i][j]; sum=sum+p[i][j]; } ave=sum/8.0; ans[i]=sqrt((max*max+ave*ave)/2.0);}for(i=0;i<5;i++)printf("%lf\n",ans[i]);return0;}求解单因子指数的程序#include<stdio.h>#include<string.h>intmain(){doubles[319][9];doublesum[6][9];doubless[8]={3.6,130,31,13.2,35,12.3,31,69};intcount[6]={0,0,0,0,0,0};FILE*fp;fp=fopen("1.txt","r");inti,j;memset(sum,'0',sizeof(sum));for(i=0;i<319;i++){ fscanf(fp,"%lf%lf%lf%lf%lf%lf%lf%lf%lf",&s[i][0],&s[i][1],&