1、数学建模之长江水质监测问题数学建模之长江水质监测问题Document number： WTWYT-WYWY-BTGTT-YTTYU-2018GT长江水质监测摘要本文解决的是长江水质的评价与监测问题，通过分析过去十年不同监测站收集 到的长江水质数据，运用不同的理论建立不同的模型，对长江过去十年的水质 情况作出评价，然后再预测未来十年长江水质的变化情况。针对问题一：考虑到问题一中需要对长江水质情况作出定量的评价，并分析各 地区水质的污染状况，为此，建立模糊综合评价模型确定了其隶属度函数，建 立评判因子的权重矩阵，求得最终结果为：水质最差的地方是江西南吕滁槎（15号），其次水质差的地方为四川乐山岷江
2、大桥（8号）、湖南长沙新港（12号）以及四川泸州沱江二桥（10号），此四处水质污染严重；水质最好的 地方是湖北丹江口胡家岭（11号）。针对问题二：根据长江的降解系数，可得到污染物随时间的变化量。山于污染 源的污染物排放量等于本地区污染物的流量与上游流下的污染物流量之差。因 此，建立污染物流量随时间变化的微分方程模型。最后求得：高镭酸钾指数和 氨氮的污染源主要集中在宜昌至岳阳之间。针对问题三：根据已知的过去10年的主要统计数据，建立了灰色预测模 型。在相对误差较小的情况下对未来10年的水质情况作出了预测，分析得出结 论：未来10年可饮用水所占的比例越来越低，排污量有明显的上升趋势。针对问题四：在
3、问题四中建立多元线性回归方程，利用最小二乘法求解系 数，在满足问题四要求的前提下，求出未来10年的允许最大相对排污量，继而 求得未来10年每年的相应排污量，后者与前者的差值与未来10年的长江水总 流量的乘积，求得最终结果如下表：未来10年预处理的排污量年代2005200620072008200920102011201220132011预处理 排污童 （亿吨）针对问题五：分析总结前儿个问题的结果，找出水质污染的根本原因。结 合考察团的调查结果，给出合理的建议和意见。最后，对模型中运用的方法进行了优、缺点评价，在模型的推广中提出了 可以建立类似模型解决生活中的一类问题。关键词：模糊评价微分方程灰色
4、预测线性回归1问题重述问题背景长江是我国第一、世界第三大河流，其水质的污染程度日趋严重，已引起了相 关政府部门和专家们的高度重视。2004年10月，山全国政协与中国发展研究 院联合组成“保护长江万里行”考察团，从长江上游宜宾到下游上海，对沿线 21个重点城市做了实地考察，揭示了一幅长江污染的真实画面，其污染程度让 人触目惊心。为此，专家们提岀“若不及时拯救，长江生态10年内将濒临崩 溃” o长江水是许多人赖以生存和发展的资源，保护水资源就是保护我们自己，就长 江近年来的水质情况，采取合理的保护和治理措施刻不容缓。解决问题为了制定岀合理的治理长江水质污染的方案，根据长江地区近两年多主要水质 指标
5、的检测数据，现讨论以下儿个问题：(1)对长江近两年多的水质情况做岀定量的综合评价，并分析各地区水质的 污染状况。(2)研究、分析长江干流近一年多主要污染物高镭酸盐指数和氨氮的污染源 的主要分布地区。(3)假如不采取更有效的治理措施，依照过去10年的主要统讣数据，对长江 未来水质污染的发展趋势做岀预测分析。(4)根据预测分析，如果未来10年内每年都要求长江干流的IV类和V类水的 比例控制在20%以内，且没有劣V类水,那么每年需要处理多少污水(5)发表对解决长江水质污染问题切实可行的建议和意见。2问题分析问题需要对长江水质作出评价和预测，根据附件中已知的数据，建立相应的评 价预测模型，分析得出长江
6、水质过去十年以及未来十年的污染状况。问题一的分析问题提供的数据主要反映了三方面的内容即：水质划分等级标准、17个测站点 在28个月份中主要污染物的浓度和水质等级。而观测点水质等级的确定是山主 要污染物的浓度决定的，因此为了对长江及各观测点作岀定量的综合评价我们 需要将主要污染物的浓度归一到一个比较量中。这个比较量要能够快速而准确 的反映出各测站点的水质污染状况，因此我们评价整个长江流域的水质状况采 用模糊综合评价模型。虽然PH值是影响水质等级的因素之一，但是给出的所有PH值都在69这个范 圉内，即PH值对于水质等级儿乎没有影响，所以在以后的问题讨论中不再考 虑PH值的影响。问题二的分析因为污染
7、物在时间和空间上是动态变化的，为简便起见，不考虑支流的因素， 只从纵向的角度，结合本地与上游污水进行水质的污染分析。根据常识可知， 污染物的污染源就是新增污染物较多的地方。因此，考虑到降解系数，求得污 染物随时间的变化量，通过相邻两主要观测点之间每千米的相对排污量判定污染 源，即本地新增污染物的量等于本地区现有的污染物减去上游通过降解后到达 该处的污染物的总量，再比上两地之间的距离，得到相对排污量。建立污染物 的量随时间变化的微分方程模型，最后通过新增污染物的相对量大小关系得出 主要的污染源地区。问题三的分析测第三问根据过去10年长江的总体水污染状况的监测数据，可以看岀长江总体水 流量变化不大
8、，但年排污的总量在增加，这使得污染河段比例增加，污染的严 重程度呈现快速增长的趋势，即每年污染情况主要与当年的排污量和总水流量 等因素有关。为此，首先利用回归分析方法确定岀可饮用水的比例与总排污量 和总水流量的关系式，然后根据过去10年排污量，利用灰色预测方法对未来的 年排污量做岀预测，最后根据总排污量的增长趋势来推断岀可饮用水比例的变 化趋势，从而可以预测出未来10年长江水质的变化情况。问题四的分析在问题四中建立多元线性回归方程，利用最小二乘法求解系数，在满足未来10 年内每年都要求长江干流的IV类和V类水的比例控制在20%以内，且没有劣V 类水的前提下，求岀未来10年的允许最大相对排污量，
9、继而求得未来10年每 年的相应排污量，后者与前者的差值与未来10年的长江水总流量的乘积，即可 求得未来10年要处理的排污量。问题五的分析分析总结前儿个问题的结果，找出水质污染的原因。结合考察团的调查结果, 给岀合理的建议和意见。3模型假设与符号说明模型假设假设一：降解系数在一定时间段固定不变；假设二：长江干流的自然净化能了近似均匀；假设三：干流的相邻两观测站点之间的排污口主要集中在下游处；假设四：干流污染物的富集主要受上游影响，支流的影响忽略不计。符号说明符号符号意义氧气属于第i等级的隶属度高镒酸盐属于第i等级的隶属度氨氮属于第i等级的隶属度第i种指标的权重各观测点模糊关系矩阵第i种指标的实时
10、均值平均排污量相对排污量断面水流平均流速某组分子在兀处的浓度降解常数每年长江总水流量每年总排污量V类水所占的比例劣V类水所占比例4模型的建立与求解问题一的解答(1)确立评判指标：根据地表水环境质量标准(GB38382002)中4个 主要项LI标准限值设置表，结合28个月的各观测点的测量结果，选取溶氧量(DO)、高镒酸盐指数(CODMn)和氨氮(NH3-N)为评价因子，设置评价因 子集：“=&,吃宀)，求出17个观测点28个月的个评价因子的均值，得到17个点的评价因子集。(2)建立评价集序 号I类II类111类IV类V类劣V类1溶解氧(DO)(或饱和率90%)653202高猛酸盐(CODMn)2
11、461015OO3氨氮（NH3-N） OO4PH值(无量纲)694个指标将水质分为6个等级。(3)建立隶属函数，进行单因素评价：山于水质污染程度和水质分级标准都是模糊的，所以用隶属度来描述分级界限 较为合理，先根据各指标的6级标准，做岀6个级别的隶属函数，其中DO的 评价指标以数值大为优，其余两个指标以数值小为优，函数如下：DO (A)以升半梯形分布建立隶属函数COOMn(B)和NH3N (C)以升半梯形分布建立隶属函数：有了各指标的隶属函数，就可以进行单因素评价，将各观测点的3个指标带入 相应的隶属函数，计算出指标的隶属度，得到各观测点模糊关系矩阵心(j=h2, 3，17),具体数据见附录一
12、。(4)建立评价因素的权重集：山于DO.COOMn ffl NH3-N等污染指标对水质的影响不同，因此对各指标应赋 予不同的权重，根据污染物对水质的污染大、权重的的原则来决定权重的大 小。定义：对于DO越大越优型：q =玉，对于COOMn和NH3-N越小越优型: 你=仝，其中，为第i种指标的权重，“是第i种指标的实时均值，坷分 别为多级浓度标准的最大值和最小值。再将权重归一化处理，Q=$,评/-I价因素权重的集合为A=o(5)模糊综合评价：评价因素的权重集A乘以单因子矩阵得到模糊综合评价结果，B = AR o将6个等级看作一种相对位置使其连续化，设各等级分别用1, 2, 3, 4, 5, 6表
13、 示，并称为各等级的秩，然后用B中对应的量将各级的秩加权求积，得到被评 等级的相对位置，即：式中久为隶属于第k级的隶属度。具体计算结果见下表：表1:各观测站水质等级序号I23456789水质级序号1011121314151617水质级山表1可知，各观测站中水质等级最大的为江西南昌滁槎（15号），可知此处 水质情况最差，其次水质差的地方为四川乐山岷江大桥（8号）、湖南长沙新 港（12号）以及四川泸州沱江二桥（10号），此四处水质污染严重；水质最好 的地方是湖北丹江口胡家岭（11号）。问题二的解答因为污染物在时间和空间上是动态变化的，为简便起见，不考虑支流的因素， 只从纵向的角度，结合本地与上游污
14、水进行水质的污染分析。建立污染物流量 随时间变化的微分方程模型，最后通过新增污染物的流量大小关系得出主要的 污染源地区。（1） 某一污染物扩散所满足的微分方程是一个抛物线方程，结合实际问题的假 设，常假设其水流近似的处于稳定状态，断面均匀，普通对流扩散方程为 上式是解决污染物的一般方程，结合问题二的情况，假设观测站干流的污染物 浓度在一个月内保持不变，则=0, 乂因为D = 0,所以上式进一步简化为：ot其中，“为断面水流平均流速，C为某组分子在x处的浓度，k为降解常数取。 山此解出长江干流污染物浓度C,与据观测站x米处的函数关系：（2） 对于长江上任一段干流AB, A为笫i个观测点，B为第i
15、 + 1个观测点，两 观测点的距离为閔，假设该段干流之间有”个排污口（包括支流入口和直排 口），为计算的方便，假设排污口聚集在一起，其总的排污口的流量为平 均流量为“，污染物的浓度为c, 即为AB段的排污量，则对于任意一个江段AB,污染物的浓度CrQ,水流量流速Ug,距离d刖和降解系数k都已知，则可以算出AB段污染物的排放量旳o（3）根据所给数据可以算出，对于每个月每个江段的排污区间W,丄工吧，表示平均排污量,13 i.i（j=l,2,.,13个月）,取13个月的平均值祈=则该江段每千米的排污量v = 为相对排污量，每个江段都有一个相对排污 “AB量，以此为指标可以判断长江干流的主要污染源。按
16、步骤运用MATLAB编求解（程序见附录3）,结果如下表: 表2： COOMn的排污量江段46002相对排污 量排序651423表3： NH3-N的排污量江段相对排污 量排序541326问题二用长江干流7个观测站点将长江分为6个江段，通过污染物流量随时间 的变化关系建立了微分方程模型，先计算出每月每段的高镭酸盐和氨氮的量， 再对六个污染源近一年所排放的污染物的数量求期望，分析得出污染物高镭酸 盐指数和氨氮的污染源主要集中在宜吕至岳阳之间。问题三的解答当系统内部信息和特性是部分已知的，另一部分是未知时，人们往往难以建立 客观的物理原形，内部因素难以辨识，以及相互之间的关系较为隐蔽，难已准 确了解这
17、类系统的行为特征，因此，对于这类问题进行定量描述，即建立模型 难度较大，便选择建立灰色系统，对问题求解。1 灰色动态模型附件4中给出了 1995-2004年长江各等级在丰水期，枯水期，水文年的百分 比，为了预测此趋势之下未来十年长江水质的变化，通过对水文年全流域水质 的变化研究来对长江未来水质污染的发展趋势作出预测分析。建立灰色预测模 型，首先舍去异常数据:1995年IV类、劣V类数据，然后从附件4中提取岀 水文年全流域随时间变化的原始数据。对于原始数据进行累加处理，即加生成数列为x=.v(1|(l),x(.x(“)” =10采用一阶单变量微分方程进行拟合，得到白化方程的GM (1,1)模型+
18、 6M(,)(r) = dt 7其中式中的a上为待定系数。综上所述，灰色动态模型为:2 模型求解取水文年全流域的饮用水预测为例，可通过对其变化研究对长江未来水质污 染的发展趋势作出预测分析。求解步骤：笫一步：从附表中提取出数据第二步：利用公式对兀做累加生成得 第三步：构造矩阵B和数据向量，工)与川满足关系K(/I) = ha其中笫四步：求解常系数0.0231 88.8436第五步：得岀表达式: 令牟)詁爲-冷)所以冷)表示还原后的值。(2)结果求解按步骤运用MATLAB编求解(程序见附录4),运用此程序根据过去十年的 数据可以预测出未来10年的排污量和可饮用水如下表：表4：未來10年排污虽年代
19、2005200620072008200920102011201220132014排污址年份2005200620072008200920102011201220132014可饮水从表5可以看岀，可饮用水所占的比例越来越低，所以，可以看岀水质越来越 差，从2012年起，可饮用水的白分比其至低于60%,如此严重的水污染，应该 引起足够重视，采用合理的措施来减少对长江水的污染。问题四的解答(1)长江干流的IV类和V类水以及劣V类水的比例与每年长江总水流量“(/) 和总排污量兀2(/)有关，因此分别以IV类和V类水所占的比例儿和劣V类水所占 比例$2为响应变量，以历年来的Xi(f) = 9205,951
20、3,9171 .26,9513,9924,8892 .8,10210,9980,9405 x2(f)= 174,179,183,207,234,220.5,256,270,285 为解释变量，做出多元线性回归，得到IV类和V类水的总和及劣V类水与长江 总水流量及总排污量的关系，做出一般的多元线性回归模型：运用mat lab用最小二乘法求解得到回归系数就求得了 IV类和V类水的总和及劣V类水与长江总水流量及总排污量的关系。问题四要解决问题是：未来10年内每年都要求长江干流的IV类和V类水的 比例控制在20%以内，且没有劣V类水。则约束条件为：综上所述，所建立的模型为：对上述回归方程进行回归检验时
21、，发现回归效果并不是太好，考虑到长江干流 的W类和V类水的比例与江总水流量及总排污量不是直接的线性关系，因此对 模型进行改进。(2)优化模型长江干流的IV类和V类水以及劣V类水的比例与每年2(/)总排污量和长江 总水流量的比值Z(0=y相对排污量有关，因此分别以IV类和V类水 所占的比例比和劣V类水所占比例儿为解释变量，z(/)为响应变量，做多元线 性回归。对表四的数据进行观测，得到Z)与时间/成线性关系，因此列出如下 线性回归方程：运用mat lab用最小二乘法进行求解，得到回归系数分别为计算时剔除了 2003年这个突变点，且检验到回归性较好。得到的回归方程如 下：将yi(r)20,y2(r
22、) = 0带入Z/f)的表达式中，得到最大允许相对排污量 Z2(r) = 0.0225 ,将( = 11,12,13,14，20带入乙(f)的表达式中得到未来10年每年预测的相对排污量，则污水处理量Hr）=nnxzi（r）-Z2（r）.0xL,因为每年的总水流量（除特殊年份如1998年）变化不大，故取前10年的总水流量均值为每1 10 ,、年的总水流量厶=0）,得到每年的污水处理量如下表:i z表6：未来10年预处理的排污量年代2005200620072008200920102011201220132014预处理 排污量（亿吨）从表6可以看出，未刁尺10年所需处士里的污水越来越多，所以，可以看
23、出力K质越来越差，如此严重的水污染，应该引起足够重视，采用合理的措施来减少对长 江水的污染。通过问题一的讣算可以看出，长江水质的枯水期、丰水期与平水期都山水流量 的不同导致水质情况也不同。从问题三、四的求解结果可以看出，在未来的十 年里，长江水质不断恶化，废水排放量逐年递增，需要处理的污水也越来越 多。根据模型假设，在不发生特大洪涝灾害时长江的水流量是一定的，总结得出解 决长江水质问题的关键，在于减少废水的排放量的同时增强废水的处理能力。鉴于此，我们认为可以从以下儿个方面着手解决长江水质问题：（1） 加大处罚力度。（2） 对污染严重地区进行重点治理。（3） 加强宣传，增强民众保护长江的意识。（
24、4） 严格进行检测，处理违规企业，严惩腐败。（5） 保护长江立法，引起政府和人民的高度重视。5模型的评价、改进与推广模型的评价（1） 问题一中建立了模糊综合评价模型，通过精确的数字手段处理模糊的评价 对象，对蕴藏信息呈现模糊性的资料作岀了比较科学、合理、贴近实际的量化 评价；（2） 问题二的模型，简化了许多因素，忽略了支流的影响，但这种理想化还是 合理的，较符合实际情况；（3）问题三采用灰色预测模型，具有少数据性、良好的时效性、较强的系统性 和关联性等特征，可以合理的对数据做岀预测。在求解模型时，为了使问题得到方便的解决，往往采用简化的手段进行求解， 因此求出的结果与真实值会存在一定的偏差。模
25、型的改进在问题三中运用了灰色动态模型对未来10年长江的水质作出了预测，但并未对 预测的结果进行检验，不知预测结果的准确度有多高，在此，最好能将模型三 的基础上对灰色模型的预测结果进行检验，以保证预测结果有较高的可信度。模型的推广模型一中，除了用模糊算法评价外，还可以用BP神经网络来求长江水质的在 综合评价而；而在求解问题三和问题四时采用的灰色预测，还可以用回归分析 预测、时间序列预测和神经网络预测等代替。6参考文献1宋晓秋编着，模糊数学原理与方法，北京，中国矿业大学岀版社,19992张俊福等编着，应用模糊数学，北京，地质出版社，20013岳超源，决策理论与方法，北京：科学出版社，20034邓聚
26、龙，灰理论基础,武汉：华中科技大学出版社，20025韩中庚，数学建模方法及其应用，北京：高等教育出版社，2005附录各观测点模糊关系矩阵:1000001000001号0.21610.216100002号0.04820.048200000.0940.09400000.51940.51940000_ 100000_ 100000_3号0.43750.437500004号0.89290.892900000.32660.32660000_0.90.90000-10000o0.97030.9703000o5号0.21430.214300006号0.28750.287500000.02970.029700
27、000.2280.2280000_10 0 0 0 0_0.0733 0.07330 0 0 0_11号10 0 0 0 012号0.2429 0.24290 0 0 010 0 0 0 00 0.83280.8328 0 0 0 1 000 0 0_0.9476 0.9476 0 0 0 0_13号0 0.0965 0.09650 0 014号0.6625 0.6625 0 0 0 00.6734 0.673400 0 00.1357 0.1357 0 0 0 0_ 0 0.6982 0.69820 0 0_1 0 0 0 0 0_15号0.1620 0.162000 0 016号0.871
28、5 0.8715 0 0 0 00 000 1 00.3897 0.3897 0 0 0 0附录二问题一的程序(以17号观测站为例)： clc;cleara二100000;0000;0000;b 二;c=b*a;d二 0;e二 0;fori=l:6d=d+c (1, i)*i :endfori=l:6e二e+c(l, i);endB 二 d/e附录三问题1的程序： clc;clear numl 1J =xlsread ;num29J =xlsread; fori=l:6Ca=numll(i,6):Cb=numll(i+1, 6);Qa=num29(2, i+2);Qb=num29(2, i+3); dl=num29(l, i+1); d2=num29(l, i+2);d二(d2-dl)*1000;Ua=num29(3, i+2);Ub=num2
