数学建模之长江水质监测问题

1、长江水质监测摘要本文解决的是长江水质的评价与监测问题，通过分析过去十年不同监测站收集到的长江水质数据，运用不同的理论建立不同的模型，对长江过去十年的水质情况作出评价，然后再预测未来十年长江水质的变化情况。针对问题一：考虑到问题一中需要对长江水质情况作出定量的评价，并分析各地区水质的污染状况，为此，建立模糊综合评价模型确定了其隶属度函数，建立评判因子的权重矩阵，求得最终结果为：水质最差的地方是江西南昌滁槎（15号），其次水质差的地方为四川乐山岷江大桥（8号）、湖南长沙新港（12号）以及四川泸州沱江二桥（10号），此四处水质污染严重；水质最好的地方是湖北丹江口胡家岭（11号）。针对问题二：根据长江

2、的降解系数，可得到污染物随时间的变化量。由于污染源的污染物排放量等于本地区污染物的流量与上游流下的污染物流量之差。因此，建立污染物流量随时间变化的微分方程模型。最后求得：高钮酸钾指数和氨氮的污染源主要集中在宜昌至岳阳之间。针对问题三：根据已知的过去10年的主要统计数据，建立了灰色预测模型。在相对误差较小的情况下对未来10年的水质情况作出了预测，分析得出结论：未来10年可饮用水所占的比例越来越低，排污量有明显的上升趋势。针对问题四：在问题四中建立多元线性回归方程，利用最小二乘法求解系数，在满足问题四要求的前提下，求出未来10年的允许最大相对排污量，继而求得未来10年每年的相应排污量，后者与前者的

3、差值与未来10年的长江水总流量的乘积，求得最终结果如下表：未来10年预处理的排污量年代2005200620072008200920102011201220132014预处理排污量（亿吨）71.2483.1194.98106.86118.73130.60142.48154.35166.22178.09针对问题五：分析总结前几个问题的结果，找出水质污染的根本原因。结合考察团的调查结果，给出合理的建议和意见。最后，对模型中运用的方法进行了优、缺点评价，在模型的推广中提出了可以建立类似模型解决生活中的一类问题。关键词：模糊评价微分方程灰色预测线性回归1问题重述1.1 问题背景长江是我国第一、世界第三大

4、河流，其水质的污染程度日趋严重，已引起了相关政府部门和专家们的高度重视。2004年10月，由全国政协与中国发展研究院联合组成“保护长江万里行”考察团，从长江上游宜宾到下游上海，对沿线21个重点城市做了实地考察，揭示了一幅长江污染的真实画面，其污染程度让人触目惊心。为此，专家们提出“若不及时拯救，长江生态10年内将濒临崩溃”。长江水是许多人赖以生存和发展的资源，保护水资源就是保护我们自己，就长江近年来的水质情况，采取合理的保护和治理措施刻不容缓。1.2 解决问题为了制定出合理的治理长江水质污染的方案，根据长江地区近两年多主要水质指标的检测数据，现讨论以下几个问题：(1)对长江近两年多的水质情况做

5、出定量的综合评价，并分析各地区水质的污染状况。(2)研究、分析长江干流近一年多主要污染物高钮酸盐指数和氨氮的污染源的主要分布地区。(3)假如不采取更有效的治理措施，依照过去10年的主要统计数据，对长江未来水质污染的发展趋势做出预测分析。(4)根据预测分析，如果未来10年内每年都要求长江干流的IV类和V类水的比例控制在20%以内，且没有劣V类水，那么每年需要处理多少污水？(5)发表对解决长江水质污染问题切实可行的建议和意见。2问题分析问题需要对长江水质作出评价和预测，根据附件中已知的数据，建立相应的评价预测模型，分析得出长江水质过去十年以及未来十年的污染状况。2.1 问题一的分析问题提供的数据主

6、要反映了三方面的内容即：水质划分等级标准、17个测站点在28个月份中主要污染物的浓度和水质等级。而观测点水质等级的确定是由主要污染物的浓度决定的，因此为了对长江及各观测点作出定量的综合评价我们需要将主要污染物的浓度归一到一个比较量中。这个比较量要能够快速而准确的反映出各测站点的水质污染状况，因此我们评价整个长江流域的水质状况采用模糊综合评价模型。虽然PH值是影响水质等级的因素之一，但是给出的所有PH值都在69这个范围内，即PH值对于水质等级几乎没有影响，所以在以后的问题讨论中不再考虑PH值的影响。2.2 问题二的分析因为污染物在时间和空间上是动态变化的，为简便起见，不考虑支流的因素，只从纵向的

7、角度，结合本地与上游污水进行水质的污染分析。根据常识可知，污染物的污染源就是新增污染物较多的地方。因此，考虑到降解系数，求得污染物随时间的变化量，通过相邻两主要观测点之间每千米的相对排污量判定污染源，即本地新增污染物的量等于本地区现有的污染物减去上游通过降解后到达该处的污染物的总量，再比上两地之间的距离，得到相对排污量。建立污染物的量随时间变化的微分方程模型，最后通过新增污染物的相对量大小关系得出主要的污染源地区。2.3 问题三的分析测第三问根据过去10年长江的总体水污染状况的监测数据，可以看出长江总体水流量变化不大，但年排污的总量在增加，这使得污染河段比例增加，污染的严重程度呈现快速增长的趋

8、势，即每年污染情况主要与当年的排污量和总水流量等因素有关。为此，首先利用回归分析方法确定出可饮用水的比例与总排污量和总水流量的关系式，然后根据过去10年排污量，利用灰色预测方法对未来的年排污量做出预测，最后根据总排污量的增长趋势来推断出可饮用水比例的变化趋势，从而可以预测出未来10年长江水质的变化情况。2.4 问题四的分析在问题四中建立多元线性回归方程，利用最小二乘法求解系数，在满足未来10年内每年都要求长江干流的IV类和V类水的比例控制在20%以内，且没有劣V类水的前提下，求出未来10年的允许最大相对排污量，继而求得未来10年每年的相应排污量，后者与前者的差值与未来10年的长江水总流量的乘积

9、，即可求得未来10年要处理的排污量。2.5 问题五的分析分析总结前几个问题的结果，找出水质污染的原因。结合考察团的调查结果,给出合理的建议和意见。3模型假设与符号说明3.1 模型假设假设一：降解系数在一定时间段固定不变；假设二：长江干流的自然净化能了近似均匀；假设三：干流的相邻两观测站点之间的排污口主要集中在下游处；假设四：干流污染物的富集主要受上游影响，支流的影响忽略不计3.2 符号说明符号符号意义Aui氧气属于第i等级的隶属度Bui高钮酸盐属于第i等级的隶属度Cui氨氮属于第i等级的隶属度a第i种指标的权重Rj各观测点模糊关系矩阵Xi第i种指标的实时均值W平均排污量w相对排污量u断面水流平

10、均流速C某组分子在X处的浓度k降解常数Xi（t）每年长江总水流量X2（t）每年总排污量y)V类水所占的比例y2（t）劣V类水所占比例4模型的建立与求解4.1 问题一的解答4.1.1 模型一的建立与求解(1)确立评判指标：根据地表水环境质量标准(GB3838-2002)中4个主要项目标准限值设置表，结合28个月的各观测点的测量结果，选取溶氧量(DO)、高钮酸盐指数(CODMn和氨氮(NH3-N)为评价因子，设置评价因子集：u=(Xi,X2,X3),求出17个观测点28个月的个评价因子的均值，得到17个点的评价因子集。(2)建立评价集顷目标准值分类I类II类田类IV类V类劣V类1溶解氧(DO)&g

11、t;7.5(或饱和率90%)653202:高钮酸盐(CODMn)<2461015OO3氨氮（NH3-N）<0.150.51.01.52.0OO4PH值(无量纲)694个指标将水质分为6个等级。(3)建立隶属函数，进行单因素评价：由于水质污染程度和水质分级标准都是模糊的，所以用隶属度来描述分级界限较为合理，先根据各指标的6级标准，做出6个级别的隶属函数，其中DO的评价指标以数值大为优，其余两个指标以数值小为优，函数如下：尸1ax-6%x=二0Au3（x）=T0x-5x-3.5-3x-7.56二x二7.5x<6x-6orx-35x:63：x<50（x>7.5orx&l

12、t;5）x-6Au2x=7.5-6x-5，0.x-3Amx=5-3Jx-126二x二7.55二x<6x-5orx-23二x二52*30Ausfx）=?x-2x2x_3orxm02二x二30:x<20人6x=；21x-20：x二2x-0DO(A)以开半梯形分布建立隶属函数COOMn(B)口NH3-N(C)以升半梯形分布建立隶属函数:1Bu1X=7-Q4-200_/x4Bu3（X）=，6-4x-610-6，0,x10Bu5x二1x<22二x:二4x-4x_10orx<44：x<66二x二10x£010二x<15x150x-2Bu2（x）=42x-4.6

13、4Bu4x=0x-310-6x-1015-10x-0.150.5-0.15x<0.150.15二x二0.5x-0.5x-0.5品3仪）=1-0.5x-11.5-10cx.x-1.5Cu5x=1x-1.5orx-0.50.5二x<11:x：1.5x三1.51.5:x<2x2x_2orx_62：x<44:二x<6x-15orx<66<x<1010二x二15Bu6x=1x，二x<0.15orx_1Cu2x=Cu4x=x-0.150.5-0.15x-0.5J-0.51.5-1x-1.52-1.50.15：x<0.50.5二x三1x一2orx三

14、11：x<1.51.5<x二2Cu6x=1x>二有了各指标的隶属函数，就可以进行单因素评价，将各观测点的3个指标带入相应的隶属函数，计算出指标的隶属度，得到各观测点模糊关系矩阵Rj(j=1,2,3.,17),具体数据见附录一。（4）建立评价因素的权重集：由于DO,COOM吊口NH3-N等污染指标对水质的影响不同，因此对各指标应赋予不同的权重，根据污染物对水质的污染大、权重的的原则来决定权重的大小。x定义：对于DO越大越优型：a=，对于COOMn和NH3-N越小越优型：minai=,其中，ai为弟i种指标的权重，Xi是用i种指标的头时均值，mi,n分X别为多级浓度标准的最大值和

15、最小值。再将权重归一化处理，Qi=a-,评价、ai因素权重的集合为A=0.145,0.226,0.629。（5）模糊综合评价：评价因素的权重集A乘以单因子矩P$Rj得到模糊综合评价结果,Bj=ARj。将6个等级看作一种相对位置使其连续化，设各等级分别用1,2,3,B中对应的量将各级的秩加权求积,4,5,6表示，并称为各等级的秩，然后用得到被评等级的相对位置，即：6、bkkB*二-'bkk1式中bk为隶属于第k级的隶属度。具体计算结果见下表:表1：各观测站水质等级厅P123456789水质1.2991.4111.4031.4561.2401.5001.1642.5001.445级1689

16、40_800序号1011121314151617水质级2.12541.00002.38881.91561.50004.15721.43241.41654.1.2 问题一的结果分析由表1可知，各观测站中水质等级最大的为江西南昌滁槎（15号），可知此处水质情况最差，其次水质差的地方为四川乐山岷江大桥（8号）、湖南长沙新港（12号）以及四川泸州沱江二桥（10号），此四处水质污染严重；水质最好的地方是湖北丹江口胡家岭(11号)4.2 问题二的解答因为污染物在时间和空间上是动态变化的，为简便起见，不考虑支流的因素,只从纵向的角度，结合本地与上游污水进行水质的污染分析。建立污染物流量随时间变化的微分方程模

17、型，最后通过新增污染物的流量大小关系得出主要的污染源地区。4.2.1 模型二的建立与求解(1)某一污染物扩散所满足的微分方程是一个抛物线方程，结合实际问题的假设，常假设其水流近似的处于稳定状态，断面均匀，普通对流扩散方程为2c_ccu=D-2kC.:t：:xx上式是解决污染物的一般方程，结合问题二的情况，假设观测站干流的污染物浓度在一个月内保持不变，则生=0,又因为D=0,所以上式进一步简化为:ft：C0u=-kCfx其中，u为断面水流平均流速，C为某组分子在x处的浓度，k为降解常数取0.2。由此解出长江干流污染物浓度Cx与据观测站x米处的函数关系：kxCxGeui(2)对于长江上任一段干流A

18、B,A为第i个观测点,B为第i+1个观测点,两观测点的距离为心书，假设该段干流之间有n个排污口(包括支流入口和直排口)，为计算的方便，假设排污口聚集在一起，其总的排污口的流量为q,平均流量为u,污染物的浓度为c,qc即为AB段的排污量,则,xi1-kcqcQi,Cieci;qcQACAeUCbQiQbB=qc=CbQb-CaQaSuUaUb2对于任意一个江段AB,污染物的浓度Ca,Cb,水流量Qa,Qb,流速Ua,Ub,距离dAB和降解系数k者B已知，则可以算出AB段污染物的排放量WjWj（3）根据所给数据可以算出，对于每个月每个江段的排污区间113一（j=1,2,.,13个月），取13个月的

19、平均值W=£Wj,表示平均排污量，则13j4J该江段每千米的排污量w=也为相对排污量，每个江段都有一个相对排污量,dAB以此为指标可以判断长江干流的主要污染源。4.2.2 模型二的求解按步骤运用MATLAB编求解（程序见附录3）,结果如下表：表2:COOMn的排污量江段S1S20S3S3S4S4S5S5S6S6S7W32510.5435266.774600229855.6211104.1531404.31相对排污量34.2245.33116.4659.71105.7560.05排序651423表3:NH3-N的排污量江段S1S2S2S3S3S4S4S5S5S6S6S7W2822.62

20、3191.464451.772559.311191.38158.46相对排污量2.974.1011.275.127.260.34排序5413264.2.3 问题二的结果分析问题二用长江干流7个观测站点将长江分为6个江段，通过污染物流量随时间的变化关系建立了微分方程模型，先计算出每月每段的高钮酸盐和氨氮的量，再对六个污染源近一年所排放的污染物的数量求期望，分析得出污染物高钮酸盐指数和氨氮的污染源主要集中在宜昌至岳阳之间。4.3问题三的解答4.3.1 模型三的建立与求解当系统内部信息和特性是部分已知的，另一部分是未知时，人们往往难以建立客观的物理原形，内部因素难以辨识，以及相互之间的关系较为隐蔽，

21、难已准确了解这类系统的行为特征，因此，对于这类问题进行定量描述，即建立模型难度较大，便选择建立灰色系统，对问题求解。1 .灰色动态模型附件4中给出了1995-2004年长江各等级在丰水期，枯水期，水文年的百分比，为了预测此趋势之下未来十年长江水质的变化，通过对水文年全流域水质的变化研究来对长江未来水质污染的发展趋势作出预测分析。建立灰色预测模型，首先舍去异常数据：1995年IV类、劣V类数据，然后从附件4中提取出水文年全流域随时间变化的原始数据。xf)=葭匕)xf*2).xC%n'n=10对于原始数据进行累加处理，即kx1k!>x0ik=1,2.10i4加生成数列为x1-k11,

22、x12.x1n1n=10采用一阶单变量微分方程进彳T拟合，得到白化方程的GM(1,1)模型dx。)八axt=udt其中式中的a,u为待定系数。综上所述，灰色动态模型为：式1T)=0.5x。k)+0.5xC*k-1)=X0/)+az，*k)=u2 .模型求解取水文年全流域的饮用水预测为例，可通过对其变化研究对长江未来水质污染的发展趋势作出预测分析。(1)求解步骤：第一步：从附表中提取出数据x0=93.1,85.3,80.7,88.4,80.2,74,73.7,76.7,77.5,681第二步：利用公式对xS撇累加生成得x1=93.1178.4259.1347.5427.7501.7575.465

23、2.1729.6797.61第三步：构造矩阵B和数据向量，x?匕x(1精足关系飞。厂b?其中Yn=1"2).x©&Jn=1085.380.788.480.27473.776.777.5681r|-135.751-218.751R(xH1)+xQ1)1-(x(u)+x。11-2(x(M1)+x。n)1=BTB，BTYn第五步：得出表达式:1-303.31- 387.61- 464.711-538.551- 613.7511-690.851- 690.851-763.61_0.023188.8436U-3753e,.0231k3846a所以Y(0/示还原后的值。(2)结

24、果求解按步骤运用MATLAB编求解(程序见附录4),运用此程序根据过去十年的数据可以预测出未来10年的排污量和可饮用水如下表：表4:未来10年排污量年代2005200620072008200920102011201220132014排污量303.01322.52343.29365.39388.92413.96440.61468.98499.18531.32表5:未来10年可饮用水所占比例年份2005200620072008200920102011201220132014可饮水69.6168.0266.4764.9563.4762.0260.658.2256.8754.554.3.2问题三的结果

25、分析从表5可以看出，可饮用水所占的比例越来越低，所以，可以看出水质越来越差，从2012年起，可饮用水的百分比甚至低于60%如此严重的水污染，应该引起足够重视，采用合理的措施来减少对长江水的污染。4.4问题四的解答4.4.1模型的建立与求解(1)长江干流的IV类和V类水以及劣V类水的比例与每年长江总水流量Xi(t)和总排污量x2(t次关，因此分别以IV类和V类水所占的比例%和劣V类水所占比例y2为响应变量,以历年来的Xt=9205,9513,9171.26,9513,9924,8892.8,102109980,94051x2t=174,179,183,207,234,220.5,256,270,

26、285】为解释变量，做出多元线性回归，得到IV类和V类水的总和及劣V类水与长江总水流量及总排污量的关系，做出一般的多元线性回归模型：必二aX1thx2tGy2=a2Xtb2X2tC2运用matlab用最小二乘法求解得到回归系数a1,b1,C1)=0.0032,0.079,15.9276a2,b2,C2=-0.0077,0.148456.2805就求得了IV类和V类水的总和及劣V类水与长江总水流量及总排污量的关系。问题四要解决问题是：未来10年内每年都要求长江干流的IV类和V类水的比例控制在20%Z内，且没有劣V类水。则约束条件为：必=-0.0077xi(t)+0.1485x2(t)+56.28

27、05<20V2=-0.0032xi(tV0.079x2(t)+15.9276=0综上所述，所建立的模型为：，i=-0.0032xi(t)+0.079x2(t)+15.9276y=-0.0077xi(t)+0.1845x2(t)+c，i=-0.0077xi(t计0.1485x2(t计56.2805<20y=-0.0032xi(t)+0.079x2(t)+15.9276=0对上述回归方程进行回归检验时，发现回归效果并不是太好，考虑到长江干流的IV类和V类水的比例与江总水流量及总排污量不是直接的线性关系，因此对模型进行改进。(2)优化模型长江干流的IV类和V类水以及劣V类水的比例与每年x

28、2(t)总排污量和长江总水流量的比值Z(t上黑相对排污量有关，因此分别以IV类和V类水所占的比例%和劣V类水所占比例y2为解释变量，Z(t)为响应变量，做多元线性回归对表四的数据进行观测，得到Z(t)与时间t成线性关系，因此列出如下线性回归方程：Z(t)=at+PZ2(t)=ayi(t)+by2(t)+c运用matlab用最小二乘法进行求解，得到回归系数分别为二；=0.0012,0.0165a,b,c=0.0002,0.00040.0185计算时剔除了2003年这个突变点，且检验到回归性较好。得到的回归方程如下：'Zi(t)=0.0012t+0.0165Z2(t)=0.0002y1(t

29、广0.0005y2(t)+0.0165将y1(t20,y2(t)=0带入Z2(t)的表达式中，得到最大允许相对排污量Z2ft)=0.0225,将t=11,12,13,14,.,20带入乙(t)的表达式中得到未来10年每年预测的相对排污量，则污水处理量w(t)=maxN(t)-Z2(t)0及L,因为每年的总水流量(除特殊年份如1998年)变化不大，故取前10年的总水流量均值为每年的110总水流量L=1£x(t),得到每年的污水处理量如下表：10tw表6:未来10年预处理的排污量年代2005200620072008200920102011201220132014预处理排污量(亿吨)71.

30、2483.1194.98106.86118.73130.60142.48154.35166.22178.094.4.2问题四的结果分析从表6可以看出，未来10年所需处理的污水越来越多，所以，可以看出水质越来越差，如此严重的水污染，应该引起足够重视，采用合理的措施来减少对长江水的污染。4.5解决长江水质问题的建议通过问题一的计算可以看出，长江水质的枯水期、丰水期与平水期都由水流量的不同导致水质情况也不同。从问题三、四的求解结果可以看出，在未来的十年里，长江水质不断恶化，废水排放量逐年递增，需要处理的污水也越来越多。根据模型假设，在不发生特大洪涝灾害时长江的水流量是一定的，总结得出解决长江水质问题

31、的关键，在于减少废水的排放量的同时增强废水的处理能力。鉴于此，我们认为可以从以下几个方面着手解决长江水质问题：(1)加大处罚力度。(2)对污染严重地区进行重点治理。(3)加强宣传，增强民众保护长江的意识。(4)严格进行检测，处理违规企业，严惩腐败。(5)保护长江立法，引起政府和人民的高度重视。5模型的评价、改进与推广5.1 模型的评价5.1.1 模型的优点(1)问题一中建立了模糊综合评价模型，通过精确的数字手段处理模糊的评价对象，对蕴藏信息呈现模糊性的资料作出了比较科学、合理、贴近实际的量化评价；(2)问题二的模型，简化了许多因素，忽略了支流的影响，但这种理想化还是合理的，较符合实际情况；(3

32、)问题三采用灰色预测模型，具有少数据性、良好的时效性、较强的系统性和关联性等特征，可以合理的对数据做出预测。5.1.2 模型的缺点在求解模型时，为了使问题得到方便的解决，往往采用简化的手段进行求解，因此求出的结果与真实值会存在一定的偏差。5.2 模型的改进在问题三中运用了灰色动态模型对未来10年长江的水质作出了预测，但并未对预测的结果进行检验，不知预测结果的准确度有多高，在此，最好能将模型三的基础上对灰色模型的预测结果进行检验，以保证预测结果有较高的可信度。5.3 模型的推广模型一中，除了用模糊算法评价外，还可以用BP神经网络来求长江水质的在综合评价而；而在求解问题三和问题四时采用的灰色预测，

33、还可以用回归分析预测、时间序列预测和神经网络预测等代替。6参考文献1宋晓秋编著，模糊数学原理与方法，北京，中国矿业大学出版社，19992张俊福等编著，应用模糊数学，北京，地质出版社，20013岳超源，决策理论与方法，北京：科学出版社，20034邓聚龙，灰理论基础，武汉：华中科技大学出版社，20025韩中庚，数学建模方法及其应用，北京：高等教育出版社，2005各观测点模糊关系矩阵:1号0,21610.2161.0.0940.094附录0.04820.0482-0.51940.5194一13号0,43750.4375。一04号I0.32660.32660.89290.892900.90.90。10

34、一15号0,21430.02970.99417号0.0465一19号0,36799.801110011号1000.21430.02970.99410.04650.367900.801101000113号0.0965-0.67340.6734一015号0,16200.69820.1620一117号0,51070.5107-0.39170.3917010.97030.97030.28750.28750.2280.22800.558608号00,6215010.55860.6215_00.84860.5767010号0,66970.84860.57670.66970.62360.62360.07330.0733012号0,24290.24290.096500.6982010.83280.83280101010.94760.947601014号0,66250.662501。13570.135701016号0,87150.8715附录二问题一的程序（以17号观测站为例）：p.38970.389701clc;cleara=100000;0.51070.51070000;0.39170.39170000；b=0.1450.2