长江水质的量化评价与预测分析建模.doc

高教社杯全国大学生数学建模竞赛承 诺 书我们仔细阅读了中国大学生数学建模竞赛的竞赛规则.我们完全明白，在竞赛开始后参赛队员不能以任何方式（包括电话、电子邮件、网上咨询等）与队外的任何人（包括指导教师）研究、讨论与赛题有关的问题。我们知道，抄袭别人的成果是违反竞赛规则的, 如果引用别人的成果或其他公开的资料（包括网上查到的资料），必须按照规定的参考文献的表述方式在正文引用处和参考文献中明确列出。我们郑重承诺，严格遵守竞赛规则，以保证竞赛的公正、公平性。如有违反竞赛规则的行为，我们将受到严肃处理。我们参赛选择的题号是（从a/b/c/d中选择一项填写）： b 我们的参赛报名号为（如果赛区设置报名号的话）： 所属学校（请填写完整的全名）： 参赛队员 (打印并签名) ：1. 2. 3. 指导教师或指导教师组负责人 (打印并签名)： 日期： 年 月 日赛区评阅编号（由赛区组委会评阅前进行编号）：高教社杯全国大学生数学建模竞赛编 号 专 用 页赛区评阅编号（由赛区组委会评阅前进行编号）：赛区评阅记录（可供赛区评阅时使用）：评阅人评分备注全国统一编号（由赛区组委会送交全国前编号）：全国评阅编号（由全国组委会评阅前进行编号）：长江水质的量化评价与预测分析摘要本文分析评价了长江近年的水质情况，确定出长江主要污染源，并分析预测了未来10年长江的污染情况。建立了模糊综合评价模型来评价长江水质。由所给17个观测站28个月的水质数据，分别求出每个观测站水质处于各类污染的隶属度，建立单因子模糊评判矩阵，结合评价指标的权系数向量，求出反映17个观测站水质状况的模糊综合评判矩阵，并进行归一化处理。求解得长江全流域1类水质断面占17.65%，2类水质断面47.06%，3类水质断面23.53%，4类水质断面5.88%，5类水质断面5.88%，并得到各断面的水质情况。改进稳态一维对流扩散水质模型，分别求出长江干流上六个江段高锰酸盐（codmn）和氨氮(nh3-n)的污染量。再结合支流的地理位置及支流观测站的污染浓度数据，分析相关图象，得出长江干流近一年多主要污染物高锰酸盐和氨氮的污染源均主要分布在：湖北宜昌至湖南岳阳江段沿岸、重庆朱沱至湖北宜昌江段沿岸，以及岷江流域的四川乐山地区。根据过去10年的长江流域水质数据，利用灰色预测理论中的gm(1,1)模型，对长江未来水质污染的发展趋势做出预测并进行了残差检验，预测结果精度优良。用matlab编程求解得出10年后可饮用水的河长预测比例。结果如下（单位：%）：年份2005200620072008200920102011201220132014比例67.5566.2264.7163.2461.8160.4059.0357.6956.3855.10若要求未来10年内每年长江干流的类和类水的比例都控制在20%以内，且没有劣类水，则每年需要处理的污水量为（单位：亿吨）： 20052006200720082009污水年处理量24.089728.606433.065137.616742.3310 20102011201220132014污水年处理量47.205152.256257.559763.084668.8667最后针对所作的水质评价与预测，提出解决长江水质污染问题的切实可行的建议和意见。一、问题重述长江是我国第一、世界第三大河流，长江水质的污染程度日趋严重，已引起了相关政府部门和专家们的高度重视。现已知长江沿线17个观测站（地区）近两年多主要水质指标的检测数据、干流上个观测站近一年多的基本数据（站点距离、水流量和水流速）以及“19952004年长江流域水质报告”给出的主要统计数据。一般说来，江河自身对污染物都有一定的自然净化能力，反映这种能力的指标称为降解系数。根据检测，主要污染物高锰酸盐指数和氨氮的降解系数可以考虑取0.2(单位：1/天)。给出国标(gb3838-2002)地表水环境质量标准中的4个主要项目标准限值（见附件8），其中、类为可饮用水。现需要做以下工作：（1）对长江近两年多的水质情况做定量综合评价，并分析各地区水质的污染状况。（2）研究、分析长江干流近一年多主要污染物高锰酸盐指数和氨氮的污染源主要在哪些地区?（3）假如不采取更有效的治理措施，依照过去10年的主要统计数据，对长江未来水质污染的发展趋势做出预测分析，比如研究未来10年的情况。（4）根据预测分析，如果未来10年内每年都要求长江干流的类和类水的比例控制在20%以内，且没有劣类水,那么每年需要处理多少污水？ （5）对解决长江水质污染问题提出切实可行的建议和意见。二、问题分析水质的评价与预测的难点在于某观测点的水质状况一直处于动态变化之中，且在一定程度上受到上游水质的影响。对于问题一，水质情况是一个模糊概念，在进行长江水质综合评价时很难给出确切的表达，可以考虑利用模糊理论进行水质评价。模糊综合评价就是由给出的评价标准和实测值，经过模糊变换，对待评价对象给出总的评价的一种方法。建模时需要找出影响水质的各主要因素，确定评价因子集、评价集、隶属函数，然后通过计算各因素的权重和隶属度，得到综合隶属度，确定水质级别。对于问题二，要研究分析长江干流近一年多高锰酸盐指数和氨氮的主要污染源地区，可考虑分两步的策略第一步：确定七个观测站之间受污染最严重的江段。污染物在河流中符合一级反应动力学，可以建立稳定一维均匀河流水质模型，分别求出长江干流上六个江段高锰酸盐和氨氮的污染量（包含了该江段支流所汇入的污染、该江段沿岸产生的污染、上游污染经自然净化后剩余的污染）；第二步：对于上一步中求得的污染严重的江段，进一步缩小范围确定主要污染源地区。再结合支流、湖泊的地理位置及支流观测站的污染浓度数据，作出相关图象，进行分析。对于问题三，根据过去10年的主要统计数据，要对未来10年的水质污染作出科学的预测，这也是治理长江的迫切需要。准确的水质预测能让人们及早地发现潜伏着的自然灾害。然而，由于只有过去十年的统计数据，甚至要预测10年的污染状况。信息不完整，预测范围过大，故用曲线拟合法、多元回归模型等可能得不到令人满意的结果，所以应寻求其他方法。对于问题四，长江水体的污染日益严重，急待解决的是如何使未来10年内每年的干流、类水的比例控制在20%以内，且没有劣类水。由于污水的定义自身不能明确表明各种污染指标的浓度，因此考虑自然降解的情况下，预测未来10年内每年的污水处理量将非常困难，但是可以考虑从数据上来挖掘信息。三、模型假设1 长江干流的自然净化能力近似均匀，在任何情况下都保持恒定0.2.2 所有的数据都来自科学采集3 入河排污口水量与水质变化稳定4 长江各观测河段的水流水质状态稳定5 四川攀枝花上流的污染忽略。四、变量与符号说明符号说明单位q水流量m3/sc污染物浓度mg/lw污染指标权重u水流速m/sx河道长度kmt时间sk降解系数 k=0.21/天五、模型建立与求解（一）问题一水质污染情况是一个模糊概念，在进行长江水质综合评价时很难给出确切的表达，可以考虑利用模糊理论进行水质评价1。模糊综合评价是由给出的评价标准和实测值,经过模糊变换,对待评价对象给出总的评价的一种方法。水质污染程度是一个模糊概念，在进行评价时很难给出确切的表达，应用模糊理论处理问题，评价结果比较合理、更加接近客观实际。下面建立模糊综合评价模型来解决两年多的水质综合评价问题。1.建立模糊评判矩阵记模糊评判矩阵为，其中表示在第i个观测站测得的处于第j级污染程度的隶属度，隶属度是通过对隶属函数的计算来确定的，隶属函数一般采用“降半梯形”的函数。由于劣类污染物指标较极端，我们根据标准将水质分成5级。以溶解氧（do）为例，即溶解氧应有对应于5个级别的隶属函数。以do的监测值为自变量x,对第j级别的隶属度为（水质标准见附录8） 将各监测断面的监测数据代入前面确定的隶属函数中，就可以计算其隶属度，进而建立每个断面的单因子模糊评价矩阵（程序见附录2）。2.确定评价指标的权系数向量当对多个（m个）目标进行综合模糊评价时，还要对各个目标分别加权。权重是衡量因子集中某一因子对水质污染程度影响相对大小的量,权重系数越大,则该因子对水质的影响程度越大。设第i个目标权系数为wi，则可得权系数向量a=(w1,w2,wm)，满足。这里仅需要考虑do、codmn，nh3-n这三种污染物对水质的影响，对它们赋予不同的权重wi。我们以污染物的超标情况决定权重，各因素的监测值相对于水质标准的超标越大,对污染的贡献越大,从而权重越大，可用下面的公式求权重系数：对于codmn和nh3-n指标：对于do指标，因为于do与其它因素性质相反,实测do浓度大,说明水质污染不严重,水质好。所以do的权重赋值取c i/ si的倒数，即其中，ii （i =1,2,3）表示第i种污染物的权重；ci 表示第i种污染物的浓度实测均值；si为第i种污染物的六个级别浓度标准限值的均值。为了进行模糊复合运算，还必须对各因子权重进行归一化处理,即:根据各观测站中三种污染因素的实测浓度ci 和平均浓度si ，得到了因素权重分配集x的归一处理化结果（程序见附录3）,如下表： 表1：三种因素权重的的归一化结果观测站溶解氧do高锰酸盐指数codmn氨氮nh3-n10.50350.32230.174120.46500.25030.284630.46140.32440.214240.39420.37250.233350.55610.30110.142860.52510.28980.185170.60660.27340.120080.34490.28900.366190.39940.28200.3187100.35580.23450.4096110.58870.30720.1041120.35040.17800.4716130.37520.37610.2486140.49710.35260.1505150.14630.05570.7980160.43120.36700.2018170.45670.32290.2204 3.求模糊综合评判矩阵b利用矩阵的模糊乘法得到模糊综合评判矩阵b ：经过归一化处理后得：同理可以得出剩余16个观测站水质数据经归一化处理后所得的模糊综合评判矩阵b2b17（程序见附件4），见下表：表2：17个观测站的归一化处理结果 1类2类3类45所属类别四川攀枝花（干流）0.38030.24380.24380.131701重庆朱沱 （干流）0.44960.27520.2752001湖北宜昌（干流）0.36990.36990.2601002湖南岳阳（干流）0.33950.33950.3209002江西九江（干流）0.50000.50000002安徽安庆（干流）0.40610.40610.1878002江苏南京（干流）0.44470.49470.0606002四川乐山（支流）0.11200.21490.22340.22540.22444四川宜宾（支流）0.24290.24290.19380.19380.12652四川泸州（支流）0.20700.20900.19530.18070.20802湖北丹江口（支流）0.65710.34290001湖南长沙（支流）0.21310.21310.28690.286903湖南岳阳（支流）0.27140.27200.27200.184603湖北武汉（支流）0.30490.30490.30490.085303江西南昌（支流）0.09170.09170.18720.12930.50015江西九江（支流）0.27010.27010.22990.229902江苏扬州（支流）0.30710.30710.30710.078703注：类别中的15级分别表示水质标准中的、类 4.评价结果图1：长江全流域水质比例 表3：长江全流域水质 1类2类3类4类5类水质比例17.65%47.06%23.53%5.88%5.88% （一）问题二由问题分析，污染源地区的确定分两步走：1确定七个观测站之间受污染最严重的江段 稳态一维河流水质模型由于江河自身对污染物都有一定的自然净化能力，所以污染物在江河中迁移时，随着迁移距离的增大，浓度会不断下降。其次水流速度对污染物的浓度也有明显的影响。例如：水流速度很大时，许多污染物到下游时还不及降解，浓度变化缓慢。基于上述分析，可利用污染物质在河流中运动变化的基本模型“稳态一维对流扩散模型” 2：其中，c为整个断面的平均浓度u为流速，e为离散系数 x为河道长度t为时间s1、s2分别反映内部降解变化和外部影响 由于长江水流状态一般为稳定、均匀流动状态，与江长相比，水面宽度很小，可认为污染物浓度在横断面上分布比较均匀，即对同一地点进行水质采样在无外界环境的影响下，可认为污染物浓度不随时间变化。其次，因为长江水流状态稳定、均匀，可认为污染物浓度随迁移距离线性变化，且不受外部影响。s1与江水的自然降解和污染物浓度有关。所以一维对流扩散的水质模型可简化得到这里的稳定一维河流水质模型：其中x为河段长度 u为河段平均流速 c为污染物浓度 k为降解解数对（2）式分离变量积分得：其中c0为初始位置的污染物浓度。 利用（3）式，建立反映六江段污染物总量的模型考查干流上相邻两点间的污染情况，即在考虑自然降解的情况下这段流域上支流和干流的污染状况。干流上各江段污染物总量可以表达为：其中，表示干流上第i+1个观察站（地区）的污染量，表示水流量，为污染物浓度。表示干流第i个观察站（地区）的污染经过自然降解到下一个观测站（地区）时的剩余污染量。 模型求解由于所给数据为13个月的各观测点（干流上）的离散数据，考虑用13个月的污染物总合来表征各江段的污染状况。利用（4）式分别对高锰酸盐指数和氨氮进行计算（程序见附件5），计算结果：表4：长江六个江段的污染总量（单位：万吨）江段1-22-33-44-55-66-7codmn109.5481118.8481155.0080100.601656.5220105.8197江段1-22-33-44-55-66-7nh3-n9.510910.753715.00068.62364.01400.5338由上表的数据可知：对于高锰酸盐（codmn）和氨氮（nh3-n）两种污染，湖北宜昌南津关至湖南岳阳城陵矶江段的污染最严重，其次为重庆朱沱至湖北宜昌南津关，和四川攀枝花重庆朱沱。2缩小范围确定污染源上文已经确定了长江上高锰酸盐（codmn）和氨氮（nh3-n）两种污染最严重的江段。为了确定具体的污染源地区，我们采用作图法： 两湖区进出口处污染比较codmn浓度比较：图2：两湖进出口codmn污染浓度比较结论：a. 洞庭湖出口处codmn浓度增大，从类水变成了类水。但比照地表水环境质量标准（gb3838-2002）可知，其仍属于可饮用水范围。可见洞庭湖湖畔工厂对洞庭湖水质中的codmn指标有清度污染。由上文结果可分析得知污染最严重的湖北宜昌南津关至湖南岳阳城陵矶江段的codmn主要污染源为该江段沿岸地区，而不是洞庭湖。b. 鄱阳湖出入口的水质按照国家质量标准，均属于类，为可饮用水。即鄱阳湖湖畔基本无codmn污染。 nh3-n浓度比较：图3：两湖进出口nh3-n污染浓度比较 结论：a. 洞庭湖出口处nh3-n浓度减小，从类水变成了类水。但比照地表水环境质量标准（gb3838-2002）可知，其仍属于可饮用水范围。可见洞庭湖湖畔工厂对洞庭湖水质基本没有nh3-n污染，湖水湖畔入口流域的nh3-n污染起到了净化作用。由上文结果可分析得知污染最严重的湖北宜昌南津关至湖南岳阳城陵矶江段的nh3-n主要污染源仍为该江段沿岸地区，而不是洞庭湖。b. 鄱阳湖对上游nh3-n污染的净化作用更加明显。分析可知鄱阳湖入口江西南昌滁槎河段为nh3-n的主要污染源。 岷江流域两观测点污染比较 codmn浓度比较：图4：四川乐山与四川宜宾codmn污染浓度比较结论：近一年中，四川乐山的codmn污染浓度普遍高于四川宜宾，其中第1、11月属于污染状态。nh3-n浓度比较：图5：四川乐山与四川宜宾nh3-n污染浓度比较注：“交汇前”指四川乐山观测站，“交汇后”指四川宜宾观测站结论：近一年中，四川乐山的nh3-n污染浓度普遍高于四川宜宾，且第1、4、10、11月这四个月份均属于污染状态。综上所述：长江干流近一年多高锰酸盐（codmn）和氨氮（nh3-n）的主要污染源是：湖北宜昌至湖南岳阳江段沿岸、岷江流域的四川乐山地区、重庆朱沱至湖北宜昌江段沿岸。（二）问题三对未来进行水质预测，是治理长江的迫切需要。准确的水质预测能让人们及早地发现潜伏着的自然灾害。然而，由于只有过去十年的统计数据，信息不完整，故用曲线拟合法、多元回归模型等皆不能得到令人满意的结果，所以用对信息质量要求不高的灰色系统分析法3来作预测，建立gm（1，1）模型。由于1998年长江流域遭遇罕见的洪水灾害，数据变化很大，因此我们先将长江1998年的数据作为灾变数据剔除。记其中表示第i年长江流域水文年可饮用水（、类）的河长比例。模型建立：1令为gm（1，1）建模序列，表示灰导数其中，2令为的ago序列，3令为的均值（mean）序列，表示白化背景值则得到gm（1，1）的灰微分方程模型为其中，经变换后得到模型求解：1.利用matlab编程求出当k=11,1220时的（程序见附件6）2.因为，变换后得到最后预测出长江流域未来十年水文年可饮用水河长比例。图6：可饮用水河长比例预测值进行残差检验，结果如下：表5：长江全流域可饮用水的河长比例实测值与预测值比较（单位：%）年份1995199619971998199920002001200220032004预测93.183.3181.4279.5777.7776.0074.2772.5970.9469.33实测93.185.380.788.480.27473.776.77768残差0-1.970.72-8.82-2.432.000.57-4.11-6.561.33相对误差0-2.330.899.993.032.700.785.36-8.461.95从上表可以看出，除1998年因为长江发生罕见自然灾害，预测值误差较大外。其余相对误差都小于9%，可以认为使用gm（1，1）模型预测效果是令人满意的。表6：对长江流域未来十年可饮用水河长比例的预测（单位：%）年份2005200620072008200920102011201220132014比例67.5566.2264.7163.2461.8160.4059.0357.6956.3855.10 由上表分析得出：未来十年内，长江流域可饮用水的河长比例逐年下降。到2014年，全流域可饮用水的河长比例下降至55.10%，即近一半为非饮用水河段。污染如此严重，长江生态濒临崩溃。（三）问题四根据问题三的预测，长江水体的污染将日益严重，急待解决的是如何使未来10年内每年的干流、类水的比例控制在20%以内，且没有劣类水。利用解决问题三的方法，基于过去10年的统计数据，分别预测未来10年（2005=i7.5, do(j,1)=do(j,1)+1; elseif x1=6&x15&x13&x12&x1=3, do(j,4)=do(j,4)+x1-2; do(j,5)=do(j,5)-x1+3; elseif x12, do(j,5)=do(j,5)+1; end if x22&x24&x26&x210&x215, mn(j,5)=mn(j,5)+1; end if x30.15&x30.5&x31&x31.5&x32, nh3(j,5)=nh3(j,5)+1; end endend3.确定评价指标的权系数向量function myfun2load datafor i=1:17, for j=1:5,do(i,j)=do(i,j)./sum(do(i,:);mn(i,j)=mn(i,j)./sum(mn(i,:);nh3(i,j)=nh3(i,j)./sum(nh3(i,:);endends=zeros(1,3);s(1)=mean(7.5 6 5 3 2);s(2)=mean(2 4 6 10 15);s(3)=mean(0.15 0.5 1 1.5 2);c=zeros(17,3);for i=1:17, for k=1:28, c(i,:)=d(k-1)*17+i,2:4)+c(i,:); endendc=c./28;a=zeros(17,3);a(:,1)=s(1)./c(:,1);a(:,2)=c(:,2)./s(2);a(:,3)=c(:,3)./s(3);w=zeros(17,3);%w(i,j)表示第i个地区第j项指标的权重for i=1:3, for j=1:17, w(j,i)=a(j,i)/sum(a(j,:); endend 4.求综合模糊评判矩阵function myfun3load datafor i=1:17,%17个地区r=do(i,:);mn(i,:);nh3(i,:);a=w(i,:);%a(i)表示第i项指标的权重 for j=1:5, for k=1:3, tj(k)=min(a(k),r(k,j); end b(i,j)=max(tj); end endfor i=1:17, for j=1:5, c(i,j)=b(i,j)/sum(b(i,:); endenddisp(综合模糊评判矩阵:)c计算结果:综合模糊评判矩阵:c = 0.9339 0.0661 0 0 0.8676 0.1324 0 0 0.6836 0.3164 0 0 0.8471 0.1529 0 0 0.4348 0.5652 0 0 0.4486 0.5514 0 0 0.3333 0.6297 0 0.0370 0.2820 0.2745 0.1690 0.2745 0.6168 0.1642 0.1642 0.0547 0.4280 0.3050 0.0381 0.2288 0.8918 0.1082 0 0 0.4474 0.3453 0 0.2072 0.6811 0.2551 0 0.0638 0.4988 0.2506 0 0.2506 0.2500 0.2143 0.0357 0.5000 0.4780 0.4822 0 0.0398 0.4611 0.3952 0 0.14375.分别求出长江干流各江段污染物总量function wenti2%分别求出长江干流各江段污染物总量load datav=3.7 2.1 0.9 0.9 1.0 1.1 1.2 3.7 1.9 0.8 0.9 1.1 1.1 1.2 3.9 2.1 1.2 1.3 1.5 1.5 1.6 4.1 2.3 1.4 1.5 1.5 1.6 1.7 3.8 2.1 1.4 1.4 1.5 1.7 1.7 5.1 4.8 1.7 1.9 2.1 3.4 3.4 3.1 2.3 1.5 1.6 1.6 1.7 1.9 2.7 1.9 0.7 0.8 0.9 0.9 1.0 3.1 1.5 0.7 0.8 0.8 0.8 0.9 2.1 1.5 0.5 0.6 0.7 0.7 0.8 2.0 1.0 0.4 0.6 0.7 0.7 0.8 1.9 0.9 0.4 0.6 0.8 0.8 0.9 2.1 1.2 0.4 0.5 0.7 0.8 0.8;%v(i,j)为i月第j个观测点的流速for i=1:6, for j=1:13, jun_v(i,j)=(v(j,i+1)+v(j,i)/2; endendjun_v=jun_v;s=950 (1728-950) (2123-1728) (2623-2123) (2787-2623) (3251-2787);s=s*1000%两相临观测点间距离,单位换算为米；for i=1:13, t(i,:)=s./(jun_v(i,:);endt=t./3600./24;%时间换为天k=3600*24*30;%一个月时间换算成秒，这里每个月按30天计算l_codmn=zeros(1,7);%l_codmn(i)表示第i个站点排放codmn总量l_nh3=zeros(1,7);%l_nh3(i)表示第i个站点排放nh3_n总量a=0.2%降解系数，单位为（1/天）q=3690 13800 21000 25600 28100 29500 29800 3720 13100 19800 20500 29800 34000 34500 4010 14200 20300 22600 29500 32100 33100 4660 16400 22700 24100 27000 31900 32100 3740 10600 24000 25900 32100 33400 35100 6280 47600 53500 53800 72800 74200 81000 3260 16200 19100 22300 24800 31000 38400 1500 8170 10600 12000 14600 17000 19600 951 6550 7400 10700 13200 14100 14900 712 4020 4570 8190 10900 12300 14400 612 3603 4510 7980 10300 13700 15100 623 4740 5180 7040 14300 21400 21500 642 3650 5400 7240 15100 20200 22100;%q(i,j)表示第i月j站点测定流量(单位：立方米/秒)l_codmn(1)=sum(codmn(:,1).*q(:,1)*k;l_nh3(1)=sum(nh3(:,1).*q(:,1)*k;for i=2:7, l_codmn(i)=k*sum(codmn(:,i).*q(:,i)-codmn(:,i-1).*q(:,i-1).*exp(-a.*t(:,i-1); l_nh3(i)=k*sum(nh3(:,i).*q(:,i)-nh3(:,i-1).*q(:,i-1).*exp(-a.*t(:,i-1);endl_codmn=l_codmn/1000/1000/10000%单位换成万吨l_nh3=l_nh3/1000/1000/10000%单位换成万吨disp(各江段分别排放codmn总量（单位：万吨）)l_codmndisp(各江段分别排放nh3-n总量（单位：万吨）)l_nh3计算结果：各江段分别排放codmn总量（单位：万吨）l_codmn = 30.2793 109.5481 118.8481 155.0080 100.6016 56.5220 105.8197 各江段分别排放nh3-n总量（单位：万吨）l_nh3 =1.6227 9.5109 10.7537 15.0006 8.6236 4.0140 0.53386.利用灰度预测对长江全流域未来10年可饮用水河长比例预测function wenti3%利用灰度预测对长江全流域未来10年可饮用水河长比例预测n=9;%除98年外共9年g0=25.8+42.6+24.7,15.3+20.2+49.8,12.2+24.9+43.6,5.2+39.8+35.2,5.6+32.8+35.6,5.9+33.1+34.7,4.4+44+28.3,4.7+41.5+31.3,1.2+26.9+39.9;%x0(i)表示第i年可饮用水的比例x0=sqrt(g0);x1=zeros(1,n);x1(1)=x0(1);for i=2:n, x1(i)=x1(i-1)+x0(i);endfor i=1:n-1, z1(i)=(x