全国数学建模竞赛获奖论文-长江水质的评价和预测

长江水质的评价和预测

摘要

水是生命之源，保护水就是保护我们自己，保护水的重中之重就是保护大江大河。

本文对近两年的水质分析，综合评价，得出了部分地区的水质污染情况，并根据十年的

数据，对未来十年水质污染发展趋势做了预测，本文可以得出结论：保护母亲河的行动

迫在眉睫！

对于问题一，为了便于综合评价，本文设出了综合水质标识指数P

i

和单因子水质标

识指数p

ik

（具体公式计算见模型建立与求解），我们通过对单个城市28个月的综合的

评价标识指数求平均值,数据如下（1.95222.1162.23012.41842.10192.2515

2.04483.54692.25092.75411.78032.8682.56282.3923.5888

2.44352.3802），综合的评价标识指数平均值越大，表示污染越严重。

对于问题二，为了判断主要污染源分布地区，本文采取判断本地排放主要污染物k

的量Q

ijk

，十三个月的Q

ijk

求和取平均值来断定主要污染源。计算数据用数列表示如下：

当为高锰酸盐指数时，（8.986，37.1748，50.907，70.4526，58.196，59.9114，58.259）

当为氨氮时，（0.4816，3.0496，4.1418，6.3864，5.0473，5.0276，2.4794）

取该数据较大的几个为污染源，为主要污染源分布地区，结果如下：高锰酸盐指数和氨

氮的污染源主要所在地分别为：湖南岳阳城陵矶，江西九江河西水厂，安徽安庆皖河

口,江苏南京林山四地；湖南岳阳城陵矶，江西九江河西水厂，安徽安庆皖河口三地。

对与问题三，对为来十年的排污量进行预测时，建立了灰色系统模型。对这十年的预测

值如下：（322．5221343．2881365.3912388.9175413.9585440.6118

468.9812303．123499.1772531.3174）

对于问题四，本文根据第三问对将来十年废水排放的预测值建立了废水排放与IV、

V类水的百分比之间的关系，Ⅳ，我们建立了百分比y与废水派放量x之间的关系y=f(x),

令y≤20,求出x的上限，则预测的废水排放量与x的上限的差值即为需要处理的污水，

从而将IV、V类水的百分比控制在20%，劣V类为0，求出了每年需要处理的污水量。

对于问题五，本文参考以上问题得出的数据，并参考一些文献资料，呼吁保护长江

人人有责，保护长江一定要采取行之有效的行动！

1

一、问题重述

水是人类赖以生存的资源，保护水资源就是保护我们自己，对于我国大江大河水资

源的保护和治理应是重中之重。专家们呼吁：“以人为本，建设文明和谐社会，改善人

与自然的环境，减少污染。”

长江是我国第一、世界第三大河流，长江水质的污染程度日趋严重，已引起了相关

政府部门和专家们的高度重视。2004年10月，由全国政协与中国发展研究院联合组成

“保护长江万里行”考察团，从长江上游宜宾到下游上海，对沿线21个重点城市做了

实地考察，揭示了一幅长江污染的真实画面，其污染程度让人触目惊心。为此，专家们

提出“若不及时拯救，长江生态10年内将濒临崩溃”（附件１），并发出了“拿什么拯

救癌变长江”的呼唤（附件2）。

附件3给出了长江沿线17个观测站（地区）近两年多主要水质指标的检测数据，

以及干流上７个观测站近一年多的基本数据（站点距离、水流量和水流速）。通常认为

一个观测站（地区）的水质污染主要来自于本地区的排污和上游的污水。一般说来，江

河自身对污染物都有一定的自然净化能力，即污染物在水环境中通过物理降解、化学降

解和生物降解等使水中污染物的浓度降低。反映江河自然净化能力的指标称为降解系

数。事实上，长江干流的自然净化能力可以认为是近似均匀的，根据检测可知，主要污

染物高锰酸盐指数和氨氮的降解系数通常介于0.1~0.5之间，比如可以考虑取0.2(单

位：1/天)。附件4是“1995~2004年长江流域水质报告”给出的主要统计数据。下面的

附表是国标(GB3838-2002)给出的《地表水环境质量标准》中4个主要项目标准限值，

其中Ⅰ、Ⅱ、Ⅲ类为可饮用水。

请你们研究下列问题：

（1）对长江近两年多的水质情况做出定量的综合评价，并分析各地区水质的污染

状况。

（2）研究、分析长江干流近一年多主要污染物高锰酸盐指数和氨氮的污染源主要

在哪些地区?

（3）假如不采取更有效的治理措施，依照过去10年的主要统计数据，对长江未来

水质污染的发展趋势做出预测分析，比如研究未来10年的情况。

（4）根据你的预测分析，如果未来10年内每年都要求长江干流的Ⅳ类和Ⅴ类水的

比例控制在20%以内，且没有劣Ⅴ类水,那么每年需要处理多少污水？

（5）你对解决长江水质污染问题有什么切实可行的建议和意见。

附表:《地表水环境质量标准》（GB3838—2002）中4个主要项目标准限值

单位：mg/L

分

序类

Ⅰ类Ⅱ类Ⅲ类Ⅳ类Ⅴ类劣Ⅴ类

号标准值

项目

7.5

溶解氧

1（或饱和65320

(DO)≥

率90%）

2

高锰酸盐指数

2461015∞

(CODMn)≤

氮

0.153

氨

0.51.01.52.0∞

（NH3-N）≤

4PH值（无量纲）6---9

2

二、符号说明

iiP

ij

第个观测站（地区）第j月的综合水质标识指数（=1、2、3、„、17；j=1、

2、3、„、28）

p

ijk

第i个观测站（地区）第j月的第k种检测项目的单因子水质标识指数（k1、

2、3；i=1、2、3、„、17；j=1、2、3、„、28）

（地区）第j月第k种主要检测项目指数（

k1,2,3

；i=1、2、3、„、x

ijk

第i个观测站

17；j=1、2、3、„、28）

t

kmmax

第m类第k种检测项目指数的上限（k1、2、3；m=Ⅰ、Ⅱ、Ⅲ、Ⅳ、Ⅴ、

劣Ⅴ）

t

kmmin

第m类第k种检测项目指数的下限（k1、2、3；m=Ⅰ、Ⅱ、Ⅲ、Ⅳ、Ⅴ、

劣Ⅴ）

A

i

第i个观测站（地区）28个月的综合的评价标识指数的平均值（i=1、2、3、„、

17）

降解系数（=0.1～0.5）

Q

ijk

第i个观测站（地区）第j月k种主要检测项目的单位时间内本地排放量（i=1、

2、3、„、17；k2、3；

j

=1、2、3、„、13）

（地区）k种主要检测项目的残留量（k2、3；i=1、U

ijk

第j月污染物在第i个观测站

2、3、„、17；

j

=1、2、3、„、13）

v

ij

第i个观测站（地区）第j月的水速（i=1、2、3、„、17；j=1、2、3、„、13）

S

i(i1)

第i个观测站（地区）与上游观测站（地区）的距离（i=1、2、3、„、17）

T

i(i1)j

在第j月水流从上游观测站（地区）到第i个观测站（地区）所需要时间（i=1、

2、3、„、17；j=1、2、3、„、13）

V

ij

第i个观测站（地区）第j月的水流量（i=1、2、3、„、17；j=1、2、3、„、

13）

（地区）k种主要检测项目的13个月的平均值（k2、W

ik

第j月污染物在第i个观测站

3

3；i=1、2、3、„、17；

j

=1、2、3、„、13）

X(0)数据统计数列

X(1)数据累加数列

（i=1，2，„，10）B

i

对未来第i年的Ⅳ，Ⅴ类水的所占河长比例的预测值。

（i=1，2，„，10）B

i

'劣Ⅴ类水的所占河长比例的预测值。

（i=1，2，„，10）X

i

对未来第i年总排污量值。

i个观测站依次是：四川攀枝花，重庆朱沱，湖北宜昌南津关，湖南岳阳城陵矶，江

西九江河西水厂，安徽安庆皖河口，江苏南京林山，四川乐山岷江大桥，四川宜宾

凉姜沟，四川泸州沱江二桥，湖北丹江口胡家岭，湖南长沙新港，湖南岳阳岳阳楼，

湖北武汉宗关，江西南昌滁槎，江西九江蛤蟆石，江苏扬州三江营。

三、模型假设

1、长江干流的自然净化能了近似均匀

2、降解系数在一定时间段固定不变

3、干流污染物的富集主要受上游影响

四、问题说明

4.1对问题一的分析理解

对附件3长江流域主要城市水质检测报告，我们对所给数据做了统计，有统计结果

我们可以看出水质分类的的依据是：单个项目超过所对应的标准限制，该项目所对应的

类型，即为水质类型。我们现在是以单个城市为研究对象，分析水质污染状况，如果我

们选取单个项目分析一个城市的污染状况，那么四个项目，就要对同一个城市做四次分

析，再经过比较分析，需要得到单个城市受主要项目指标的影响程度，对十七个城市，

每个城市做四次分析，我们可以看到，工作量是多么巨大。并且如此得到的结果是：单

个城市近两年多水环境质量标准类型主要是受那个污染物的影响。由于原问题要求给出

定量的综合评定，我们如果用上述方法，那么必然不符合定量分析的要求，所以我们设

出了一个综合评价的指数，于是设该指数为：综合水质标识指数。

在设出综合水质标识指数后，我们按月份计算出每个月的综合水质标识指数的数

值，这样可以以月份为横轴，综合水质标识指数为纵轴，画出单个的城市随月份不同，

综合水质标识指数变化的曲线。用这种方法，我们仅用十七个图就综合展现了各个城市

水质量污染状况的变化趋势，这样更加直观简洁，同时我们使用综合水质标识指数，也

确实做到了定量分析，综合来说这种方法是很具有说服力的。

4．2对问题二的分析理解

对长江干流近一年多主要污染物高锰酸盐指数和氨氮的主要污染源的判断，我们必

须考虑到上游对下游的的影响。上游污染物通过水流可实现对下游的富集，故必须通过

水流量，水流速，站点距离对污染物富集的影响。

4

如果我们单纯的只确定污染源，根据主要检测项目值就可简单判断，由于污染物从

上游到下游，要经过降解过程，如果下游污染物的检测项目值大于上游的值，那么就可

简单判断此地是污染源。但根据题意我们必须判断的是主要污染源，由此不能把污染源

当作主要污染源，我们必须比较的是当地的排污能力，才能断定那一地点是主要污染源。

综合上述针对问题三，我们主要的目的就是建立模型计算出单个地点的排污能力，

再综合比较，得出主要污染源，另外由于降解系数（0.1—0.5）不同，我们可以得到

不同的污染物富集结果，由此也得到不同的单个地点的排污能力。

通过这种方法，我们就可方便解决判断主要污染源的地点问题。

4．3对问题三的分析理解

首先我们对评价河长，这一概念进行了分析，通过概率论的有关知识，我们可以对

评价河长，理解为随机取样，有部分研究全体的问题。由此有评价河长，得出的各类水

质所占河长的百分比，可认为是全流域的该类水质所占河长的比例。这样就达到了对评

价河长的研究，实现了对全流域的研究。

由于该问题是通过过去十年的数据实现对未来十年的预测，所以有必要建立一个各

类水质类型河长所占全流域比例，随年度变化的函数模型。但通过对各类水质类型河长

所占全流域比例，随年度变化的图形，我们看到它们的在坐标图（以年度为横轴，比例

为纵轴）上的点是很离散的，我们很难得到残差比较小，并且拟合程度很好的函数。所

以我们又考虑到对后三类水质类型河长所占全流域比例之和，与年度变化建立函数模型

（该模型如果建立，也可对污染程度预测，因为后三类为非饮用水，是重污染水域），

但结果同上，依然得到很离散的图形（该图见图（一）），我们同样通过线形回归很难得

到拟合很好的函数。这时我们

图（一）图（二）

又考虑用设一个水质年度级别来判断，这样每年对应一个年度级别，综合预测未来十年

的水质级别，但这个级别很难用一个有说服力的标准来确定，另外问题四，也要用到未

来十年的Ⅳ类，Ⅴ类，劣Ⅴ水的比例数据，故设出水质年度级别，并不方便解题，并无

端增加数据统计上的麻烦。另外我们试图建立各个类水的比例与排污量的关系图（以单

个类的水十年所占百分比为横轴上的点（单位为：%），年排污量为纵轴上的点（单位为：

亿吨）），我们依然看到它们的图形是很离散的（图二为Ⅰ类水的图形）。我们不优先考

虑此类模型。最后我们不得不考虑最初所设想的模型的可行性，建立每年排污量与年度

的对应关系函数时，我们选择了建立了灰色系统模型，这样比较精确的对问题求解。同

样年度排污量与年度的变化关系也有灰色系统模型求出。

5

4．4对问题四的分析理解

只要得出Ⅳ，Ⅴ类水的所占河长比例的预测值、劣Ⅴ类水的所占河长比例的预测值

和未来第i年总排污量值就可以简单的算出排污量。我们便可以得到所需模型。

五、模型的建立与求解

5．1根据所给两年多主要检测项目的数据，对长江两年多的水质情况做出定量的综合

评定，并分析各个地区水质的污染状况。

依照《地表水环境质量标准》（GB3838-2002）中的四个主要项目标准限值，我们根

据两年多的水质情况表，我们判断水质类别是有单个项目超标决定的，但为了综合评定，

我们有必要把各个项目都考虑在内，由此得到综合的评价标识指数。我们首先建立单因

子水质标识指数函数模型。我们可将模型分为两类，这是因为容解氧(DO)的单因子水质

标识指数判断是有该类型的下限值比较所决定的，而高锰酸盐指数（KMn）与氨氮（NH3-N）

的单因子水质标识指数是有该类型的上限值比较所确定的。由于PH值在判断水质类别

时，没有明显作用，我们在考虑综合的评价标识指数时不在考虑它对该指数的影响。最

后我们可以对单个城市28个月的综合的评价标识指数求平均值，作为最后水质污染状

况的综合评价依据，该数据越大表示污染程度越深，并设该平均值为A

i

：

在这里我们求出第i个观测站（地区）第j月容解氧(DO)的单因子水质标识指数p

ijk

，

我们为了与水质类别建立联系，故得到如下的计算公式：

p

ijk



t

kmmax

x

ijk

t

kmmax

t

kmmin

+f(m)，（k=1）

上式中当m=Ⅰ时f(m)=1，当m=Ⅱ时f(m)=2，当m=Ⅲ时f(m)=3，当m=Ⅳ时

f(m)=4，当m=Ⅴ时f(m)=5，当m=劣Ⅴ时，f(m)=6。

有上式我们先判断k种（k=1表示容解氧(DO)，k=2表示高锰酸盐指数，k=3表示

氨氮）所属类型m对应值的t

kmmax

和t

kmmin

，再进行计算。其中当k=1，m=Ⅰ时，根据附

录3给出的过去两年多的数据，经统计得出t

kmmax

=14.4，t

kmmin

=7.5。对于当k=1，m=

Ⅱ、Ⅲ、Ⅳ、Ⅴ、劣Ⅴ时，参照《地表水环境质量标准》上的数据便可得到t

kmmax

和t

kmmin

的值。

用同样方法我们可求出第i个观测站（地区）第j月高锰酸盐指数（KMn）与氨氮

（NH3-N）的单因子水质标识指数p

ijk

，故得到如下的计算公式：

p

ijk



x

ijk

t

kmmin

t

kmmax

t

kmmin

+f(m)，（k=2、3）

6

上式中同样当m=Ⅰ时

f(m)

=1，当m=Ⅱ时

f(m)

=2，当m=Ⅲ时

f(m)

=3，当m=

Ⅳ时

f(m)

=4，当m=Ⅴ时

f(m)

=5，当m=劣Ⅴ时，

f(m)

=6。

有上式我们先判断k种（k=1表示容解氧(DO)，k=2表示高锰酸盐指数，k=3表示

氨氮）所属类型m对应值的t

kmmax

和t

kmmin

，当k=2，m=Ⅳ时，t

kmmax

=10，t

kmmin

=6；当k=2，

m=Ⅴ、劣Ⅴ时，由于根据附录3给出的过去两年多的数据，经统计得出高锰酸盐指数

的上限值为9.9，所以在这里t

kmmax

和t

kmmin

不需要再给出。当k=2，m=Ⅰ、Ⅱ、Ⅲ时参

照《地表水环境质量标准》上的数据便可得到t

kmmax

和t

kmmin

的值。

当k=3，m=劣Ⅴ时，根据附录3给出的过去两年多的数据，经统计得出氨氮的

t

kmmax

=24.2，t

kmmin

=15。当k=2，m=Ⅰ、Ⅱ、Ⅲ、Ⅳ、Ⅴ时，参照《地表水环境质量标

准》上的数据便可得到t

kmmax

和t

kmmin

的值。

有以上单因子水质标识指数p

ijk

确定后，我们可得到综合的评价标识指数P

ij

，并有

计算公式可求出，公式如下：

13

P

ij

p

ijk3

k1

有以上公式我们用matlab编程（见附录一）分别得出如下各个观测站（地区）第

28个月的综合的评价标识指数，我们列出部分城市数据及根据月份，单因子水质标识指

数和综合的评价标识指数变化的曲线图（其余数据见附录一），如下：

四川攀枝花：

（1.74442.45991.72101.69832.14921.63941.81432.16202.4505

1.93451.97562.12871.77661.77583.01082.44891.53721.8444

1.64251.55851.72581.85331.66501.86301.93162.12371.8976

2.1295）

重庆朱沱：

7

（2.27031.97831.96361.95022.06702.07341.96702.06002.19032.2765

2.35162.13462.17782.23462.00702.36511.92881.99882.08322.0694

2.17872.20322.36732.23671.84342.11342.13092.0273）

根据matlab编程，绘图所得结果，我们通过对单个城市28个月的综合的评价标识指

数求平均值A

i

，公式如下：

128

A

i

P

ij28

j1

我们分析各个地区的水质污染状况如下：i个地区的28个月的综合的评价标识指数

求平均值A

i

依次如下：（1.95222.1162.23012.41842.10192.25152.0448

3.54692.25092.75411.78032.8682.56282.3923.58882.4435

2.3802）。A

i

越大表示该地区污染越严重，有该数据我们判断四川乐山岷江大桥，江西

南昌滁槎两地污染最严重，重庆朱沱、湖北宜昌南津关、湖南岳阳城陵矶、江西九江

河西水厂、安徽安庆皖河口、江苏南京林山、四川宜宾凉姜沟、四川泸州沱江二桥、

西九江蛤蟆石、江苏扬州三江营、湖南长沙新港、湖南岳阳岳阳楼、湖北武汉宗关属

于中度污染，其余两地污染不是很严重。

5.2根据近一年多的干流统计数据，研究、分析主要污染物高锰酸盐指数（KMn）与氨

氮（NH3-N）的主要污染源所在地。

由于水从上游到下游流动的过程中，必然会产生污染物的富集。又由于水系自身有

自净化能力，我们可以通过上下游主要监测项目的指标值的对比，简单判断某一个地区

是否存在污染源。但我们的目的是为了得到主要污染源，所以我们必须建立模型，通过

对比各地污染物的排放量，确定污染源所在的主要地区。

我们首先确定污染物第j月水从上游观测站（地区）流到第i个观测站（地区）经

过时间T

i(i1)j

后，第j月污染物在第i个观测站（地区）k种主要检测项目的残留量U

ijk

（k=2表示高锰酸盐指数，k=3表示氨氮）的计算公式，计算过程如下：

由于

dQ

(i1)jk

dT

i(i1)j

Q

(i1)jk

（Q

(i1)jk

为上游观测站（地区）第j月k种主要检测项目的单

位时间内本地排放量，为降解系数）

于是可以得到在在T

i(i1)j

时间内第i个观测站（地区）的上游污染物的降解量，再有

上游观测站（地区）第j月k种主要检测项目的本地排放量减去降解量，便可得到残留

量的计算公式：U

ijk

Q

(i1)jk

Q

(i1)jk

eT

i(i1)j。

有以上公式，我们便可求出第i个观测站（地区）第j月k种主要检测项目的本地排

8

单位时间内放量Q

ijk

（k=2表示高锰酸盐指数，k=3表示氨氮）的计算公式：

Q

ijk

x

ijk

V

ij

U

ijk

（V

ij

表示第i个观测站（地区）第j月的水流量）

113

为了方便分析，我们对十三个月的Q

ijk

求和取平均值W

ik

Q

ijk

，来确定污染源

13

j1

的主要所在地。

经过matlab编程（见附录二），我们可以计算出当=0.2（=0.1～0.5）时W

ik

的

的数值为：

当k=2时（k=2表示高锰酸盐指数）求W

ik

（单位为：吨）用数列表示为：（8.986，

37.1748，，50.907，70.4526，58.196，59.9114，58.259）

得到长江七个干流城市排污量对比图如下：

对比0.1NH3-N

8.986

58.25937.1748

50.907

59.9114

58.196

70.4526

1

2

3

4

5

6

7

8

由此我们得出高锰酸盐的主要污染源所在地为i=4，5，6,7三地。即湖南岳阳城陵矶，

江西九江河西水厂，安徽安庆皖河口,江苏南京林山三地。

当k=3时（k=3表示氨氮）求W

ik

（单位为：吨）用数列表示为：（0.4816，3.0496，

4.1418，6.3864，5.0473，5.0276，2.4794）

得到长江七个干流城市排污量对比图如下：

9

长江七个干流城市排污量对比表

2.4794

5.0276

0.4816

3.0496

4.1418

5.0473

6.3864

1

2

3

4

5

6

7

由此我们得出氨氮污染源的主要所在地为i=4，5，6三地。即湖南岳阳城陵矶，江西

九江河西水厂，安徽安庆皖河口三地。

另外根据取值不同，所得每个城市各个月的Q

ijk

的大小也会不同，就可以断定出

污染源的主要所在地也不同，我们列表表示如下（具体计算数据见附录）：

当=0.1k=2时

NH3-N对比图

777

1

2

3

4

5

6

7

有以上图形我们可以看出降解系数对污染源断定影响并不是很大，经过分析我们可

以断定污染源的主要所在地排量比积累量大的多，更加说明了模型的正确性。

5.3在不采取有效治理的措施的情况下，依照过去是十年的主要统计数据，对长江未来

水质污染发展趋势作出预测分析，比如研究未来十年的情况。

通过对《附件4：1995年～2004年长江流域水质报告》的分析，我们决定对长江

未来水质污染做两个方面的预测，这两个方面是：一、对未来十年每年排污量进行预测。

二、对未来十年Ⅳ、Ⅴ、劣Ⅴ类水所占全流域比例作出预测。

在对为来十年的排污量进行预测时，我们建立了灰色系统模型。当系统内部信息和

特性是部分已知的，另一部分是未知时，人们往往难以建立客观的物理原形，内部因素

难以辨识，以及相互之间的关系较为隐蔽，难已准确了解这类系统的行为特征，因此，

对于这类问题进行定量描述，即建立模型难度较大，我们便选择建立灰色系统(Grey

System)，对问题求解。灰色系统分析建模方法主要是根据具体灰色系统的行为特征数

10

据，充分利用数量不多的数据和信息寻求相关各因素自身和因素之间的数学关系，我们

对第一方面排污量进行预测，建立相应的数学模型模型建立过程如下：

由已知数据，对于1995年～2004年每年废水排污量统计为以下数据：

X(0)(X(0)(1),X(0)(2),,X(0)(10))

(174,179,183,189,207,234,220.5,256,270,285)

对

X(0)作一次累加，则：

X(1)X(1)

，

X(l)X(0)(i)(l1,2,3,,10)；(0)(1)(1)

i1

l

于是得到类加数列为：

X(1)(X(1)(1),X(1)(2),,X(1)(10))

取

X(1)的加权均值，则

Z(1)(k)X(1)(l)(1)X(1)(l1)(k2,3,,10)

，为确定参数，设

0.5，记

Z(1)(Z(1)(2),Z(1)(3),,Z(1)(10))（5.3--1）

于是GM(1,1)的白化微分方程模型为

dX(1)

aX(1)b

（5.3--2）

dt

其中

a

是发展灰度，b是内生控制灰度.

由于

X(1)(l)X(1)(l1)x(0)(l)

，取

X(0)(l)

为灰导数，

Z(1)(l)

为背景值，则将方程对

（5.3--2）相应的灰微分方程为

X(0)(l)aZ(1)(l)b(l2,3,,10)

或

X(0)(l)aZ(1)(l)b(l2,3,,10)

即矩阵形式为

Y(0)B(a,b)T



Z(1)(2)Z(1)(3)Z(1)(10)



(0)(0)(0)(0)其中

Y(X(2),X(3),,X(10))

,

B



.

11



1



T

用最小二乘法求得参数的估计值为

ˆ

)T(BTB)1BTY(0)(a

ˆ

,b

于是方程（5.3--2）有响应（特解）

b



b



X(t1)



X(0)(1)



eat

a



a



则

b



X(l1)X(l1)X(l)



X(0)(1)



ealea(l1)（5.3—3）

a





(0)



(1)



(1)



(1)

11

对未来第i年的Ⅳ，Ⅴ类水的所占河长比例的，我们建立了百分比y与废水派放量

x之间的关系y=f(x)

我们可有公式（5.3—3）求出1995～2004十年的计算数值，经过残差检验分析（数

据见附录），我们断定所建模型是比较合适用于解决问题的。利用matlab编程（程序见

附录三）我们得出对未来十年预测的求解数据如下（322．5221，343．2881，365.3912，

388.9175，413.9585，440.6118，468.9812，303．123，499.1772，531.3174）

我们对残差分析如下图：

对未来第i年的Ⅳ，Ⅴ类水的所占河长比例的，我们建立了百分比y与废水派放量

x之间的关系y=f(x)

由于劣V类的趋势不明显，因而直接利用1减去第一大类与第二大类的比例即得劣

V类的比例。

关于对未来第i年的Ⅳ，Ⅴ类水的所占河长比例的预测值与对未来第i年总排污量值

的关系函数模型建立也使用该灰色系统模型。

枯水期干流

支流

丰水期

水文年

干流

支流

干流

支流

IV、V类水比例之和

Y=(x-170)-21.02.

y=-0.003257log(x-170)+1767.

Y=6.076log(x-170)-3.07

y=4.1178log(x-170)+1.31

y=8.4197log(x-170)-13.03.

Y=0.01177log(x-170)+17.22

5.4根据解题过程5.3中，对未来水质污染发展趋势做出的预测，做出未来十年内每年

长江干流的Ⅳ类和Ⅴ类水的比例控制在20%，并且没有劣Ⅴ类水的污水处理模型。

12

有5.3中所求出数据，即Ⅳ，Ⅴ类水的所占河长比例的预测值对未来第i年的Ⅳ，

Ⅴ类水的所占河长比例的，我们建立了百分比y与废水派放量x之间的关系y=f(x),令

y≤20,求出x的上限，则预测的废水排放量与x的上限的差值即为需要处理的污水量，

便得到了年处理量。

未来十年每年废水排放量和需处理量表注：单位：亿吨

年2014

份

排308.9331.8356.3382.5410.4440.0471.3504.2538.7575.0

放

量

处98.9121.8146.3172.5200.4230.0261.3294.2328.7365.0

理

量

5.5根据求解过程5.1和5.3所得结果和对长江污染问题的了解，提出切实可行的建议

和意见。

近些年来水质污染状况越来约严重，我们赖以生存的水的现状是———淮河浊臭、

黄河断流、长江告急。水是万物之灵，生命之源，如若对水质污染置之不理，任其发展，

那么子孙后代将会喝不到放心的饮用水，我们这代人将有愧与我们的民族！

通过我们自己对长江水系污染状况的研究，我们深深感受到拯救长江的时代责任感和紧

迫感！我们深信专家所言：长江已陷入深度危机，若不及时拯救，10年之内，长江水系

生态将濒临崩溃。所以拯救长江必须行动了，我们在下文对解决长江污染问题，提出了

自己的看法和意见。

第一，长江沿岸地区要改变GDP增长模式，提升对招商引资的技术要求。一些地区

只重视GDP的增长，对于招商引资项目缺乏管理。如江苏某市招商引资来的都是钢铁厂、

化工厂、造纸厂、造船厂、拆船厂这样的重污染企业，沿几十公里长的江岸一字排开。

在我国，消耗了世界1/4钢铁，1/3煤炭，1/2水泥才造就出世界1/30GDP。如此资源成

本造成长江上游严重水土流失，中游生态脆弱，下游污染惊人。而管理者更是缺乏对环

境保护的关注，甚至有人认为：沿江建厂，长江运输成本低、排污方便，直排长江，直

冲就东海。还有写地方重复建设严重，对经济只重视量的增长，而不重视质的提高！我

国提倡大力发展经济，重视招商引资，但决不是不分层次，只重视量，不重视质。所以

要改变经济增长方式，改资源密集型，劳动密集型的企业增长方式为技术密集型的企业

增长方式，扭转高消耗、高污染、低产出的状况，全面转变经济增长方式。同时对招商

引资进行管理，利用外资开发本地资源的同时，也要根据环保要求限制部分企业在沿江

流域的建设。GDP增长不应以生态赤字为交换条件，否则，不是福音而是灾难。

第二，完善法律，完善管理方式和制度，加强对长江水系的管理。例如：对于POP′S，

我国因无系统调研，也就无POP′S指标，所以排污就无人管理。所以我国应加快完善

检查指标，另外一些地区的干部缺乏管理意识，忽视环保，我们应该加强对对干部的环

抱教育。增强干部环保观念，并把环保作为考核干部指标之一。另外就如国家防汛指挥

部一样，我们应该成立一个同意的监管部门，统一管理长江水系。法不在立，而在必定，

所以必须完善环保法律，绝不让不法分子专法律的漏洞。

第三，增强舆论向导，提高全民环保意识。很多地区除了领导干部环保意识差外，

13

广大群众也是极其缺乏环保意识。严重污染成黑色的小河边，有的农民在那里钓鱼。脏

水里面都是剧毒物质，鱼闻起来有很重的腐臭味道，怎么能吃呢！另外根据有关资料，

沿江的生活污水很多就直接排入长江，长江变成了随意堆放、倾倒垃圾的“垃圾河”，

进入三峡后，垃圾带白茫茫一片覆盖整个江面，生活与工业垃圾、油类液体、死牛、死

羊都有。如果能够唤起群众意识，那么必然有一个异常强大的环保联盟出现。

第四，保护长江源。嘉陵江的水还有点蓝，长江的水就是蜡黄蜡黄的，像黄河一样，

含有大量泥沙。在重庆市，两江交汇的地方形成泾渭分明的一道奇景：一蓝一黄，分界

线特别明显，这说明长江上游的水土流失越来越严重。我们应当采取有效措施保护长江

源。

第五，加大对治污企业的扶持力度。有关数据显示：沿江351座污水处理厂，绝大

多数因缺乏资金或网管滞后而闲置。相当数量污水处理厂心理趋同“我花很多钱治污，

对岸的、隔壁的都不治污，我不污染，白不污染。”某市一先进的污水处理厂，该厂实

话实说，仅治理5%，就花了4.4亿元，95%便作为一般性处理了。先进处理，1千吨需

300元，一般性处理60元。所以政府部门要出台相应政策，引导该类企业发展，同时给

予优惠发展或优先发展的措施。另外银行部门向治理污水的企业提供优惠的贷款政策，

提高企业治理污水的积极性。

第六，严调狠打，加重对排污企业的处罚力度。《新民周刊》：对污染排放超标企业

最多处罚10万的力度是否过低？有这样一个案例：有个工厂污染排放超标，环保局为

了处罚它必须取样、化验„„一系列取证工作就花了50万。投入如此巨大，罚款却寥

寥无几，环保局哪里来积极性？排污企业窃喜：有本事你就罚我！违法成本低、守法执

法成本高。致使环保局有法不依、执法不严、违法不究。《水环境保护法》严重滞后，

最多罚10万元额度，不能有效遏制违法企业偷排、漏排、超标排放。利益驱动，其瞒

和骗到了无顾老百姓死活程度。所以应当增加罚款额，对违法企业重罚！

第七，因地制宜，重点城市重点治理。我们分析长江近些年的水质情况，可以知道，

不同地区，水质的污染程度不同，主要受的污染物也不同。长江岸边形成了许多污染带，

在干流21个城市中，重庆、岳阳、武汉、南京、镇江、上海六大城市累计污染带长度

占长江干流污染带总长的73%。我们要根据不同城市的状况出台有针对性的治理方案。

俗话：一方水土养一方人。眼下长江陷入一方人害一方水土，一方水土害一方人

的恶性循环。“水能载舟，亦能覆舟”，千古训诫。科学发展观扬载抑覆。否则，还要

吞噬多少人健康与生命？我们用以上提议，七剑合壁，实现教育治污、法律治污、科技

治污，投入治污。这样就能有效的保护长江，不让环境污染危害人类的悲剧上演！

六、误差分析

综合水质标识指数较为完整的标识了总体水质类别，水质情况，它急能定量评价，

也能定性评价；即不会因为个别水质指标较差就否定整体水质，又能对河流整体水质作

出合理评价；既可以在同一类别中比较水质优劣，也可以比较不同劣Ⅴ类水的污染严重

程度。但是单个项目标准限值决定水质分类，可以有效的解决水分类的饮用问题。综合

评价指数，对于Ⅱ、Ⅲ、Ⅳ类的计算数据上很难断定是否为饮用水。

通过对主要检测项目的本地单位时间内排放量十三个月的数值求和取平均做对比，

我们分析污染源的主要所在地，此方法比较综合的评价了一年多的某地区的排污情况，

经过对比很有说服力，但是这种方法缺少了对单个月份的变动情况考虑，对污染源的主

要所在地有可能在单个月改变缺少预测。

灰色系统模型，对于不太离散的问题情况分析较为精确，但对于离散的问题，使用

14

该模型时，就会出现失真的情况。

七、模型改进与推广

对于综合的评价标识指数和单因子水质标识指数，我们可以看到它的综合评价作用

是很优越的，我们可以把这两类指数推广到其他污染评价上，甚至可以对用一定方法对

作出经济发展指数，用于评价经济状况。该指数计算方法经过简单的变换，我们便很容

易的得到社会很多方面的评价指标。

灰色系统分析建模方法主要是根据具体灰色系统的行为特征数据，充分利用数量不

多的数据和信息寻求相关各因素自身和因素之间的数学关系，即建立相应的数学模型.

灰色系统可以推广到工程技术、经济管理、气象预报，以及政治、社会、工业、农业等

领域统计预测的方面，对于不是很离散的情况，该模型求解很有说服力。

对于根据长江水系的污染状况，我们提出的治污意见，对其他水系同样适用。

八、参考文献

[1]何寿强周萍张保平，嘉陵江广元段水质状况分析与评价，《四川环境》，

第14卷第3期：57～61页，1995年。

[2]幸治国钟成华，长江、嘉陵江重庆干流江段－维水质数学模型，《重庆环境

科学》，第16卷第3期：19～23页，1994年。

[3]黄树声高冬兰，南流江水资源保护问题不容忽视，《广西水利水电》，第四

期：24～28页，1995年。

15

附录一

functionzohn2

X=[174,179,183,189,207,234,220.5,256,270,285];

x=[174,179,183,189,207,234,220.5,256,270,285];

Y=[174,179,183,189,207,234,220.5,256,270,285];

y=[];

Z=[];

z=[];

r=[];

fori=2:10

X(i)=X(i)+X(i-1);

end

forj=2:10

Z(j)=0.5*X(j)+0.5*X(j-1);

end

Z(1)=[];

Y(1)=[];

B=[-Z;ones(1,9)]';

C=inv(B'*B)*B'*Y';

a=C(1);

b=C(2);

fork=1:9

y(k+1)=(174-b/a)*exp(-a*k)+b/a;

end

y(1)=174;

fort=1:9

z(t+1)=y(t+1)-y(t);

end

z(1)=174;

fors=1:10

r(s)=(x(s)-z(s))/x(s);

end

m=[732004];

errorbar(m,z,r);

title('模型误差棒棒图');

z

以上是对单个城市一个月的求解程序，如果城市改变，仅需要改变输入的X，x，Y值即

可。

附录二

16

functionzohn2

X=[174,179,183,189,207,234,220.5,256,270,285];

x=[174,179,183,189,207,234,220.5,256,270,285];

Y=[174,179,183,189,207,234,220.5,256,270,285];

y=[];

Z=[];

z=[];

r=[];

fori=2:10

X(i)=X(i)+X(i-1);

end

forj=2:10

Z(j)=0.5*X(j)+0.5*X(j-1);

end

Z(1)=[];

Y(1)=[];

B=[-Z;ones(1,9)]';

C=inv(B'*B)*B'*Y';

a=C(1);

b=C(2);

fork=1:19

y(k+1)=(174-b/a)*exp(-a*k)+b/a;

end

y(1)=174;

fort=1:19

z(t+1)=y(t+1)-y(t);

end

z(1)=174;

s=[79

20132014];

plot(s,z,'r');

z

以上是对单个城市一个月的求解程序，如果城市改变，仅需要改变输入的X，x，Y值即

可。

附录三

functionzohn33

A=[5.184.855.354.046.29.237.777.170.885.835.175.555.065.544.88

4.979.187.369.37.26.656.186.184.794.43.383.873.39];

B=[1.111.70.771.832.63.191.51.90.60.80.42.73.427.14.34.51.8.7

0.51.71.62.22.21.1];

C=[0.921.012.183.066.512.91414.54.22.090.731.482.531.19

1.620.637.055.3711.36.023.010.990.721.181.180.97

0.92];

17

a=[];

b=[];

c=[];

s=[];

x=[];

fori=1:28

ifA(i)<=14.4&7.5

a(i)=(14.4-A(i))/6.9+1;

elifA(i)<=7.5&6

a(i)=(7.5-A(i))/1.5+2;

elifA(i)<=6&5

a(i)=6-A(i)+3;

elifA(i)<=5&3

a(i)=(5-A(i))/2+4;

elifA(i)<=3&2

a(i)=3-A(i)+5;

elifA(i)<=2&0

a(i)=(2-A(i))/2+6;

end

end

a;

forj=1:28

ifB(j)<2&B(j)>=0,

b(j)=B(j)/2+1;

elifB(j)<4&B(j)>=2,

b(j)=(B(j)-2)/2+2;

elifB(j)<6&B(j)>=4,

b(j)=(B(j)-4)/2+3;

elifB(j)<10&B(j)>=6,

b(j)=(B(j)-6)/4+4;

end

end

b;

fork=1:28

ifC(k)<=0.15&C(k)>0,

c(k)=C(k)/0.15+1;

elifC(k)<=0.5&C(k)>0.15,

c(k)=(C(k)-0.15)/0.35+2;

elifC(k)<=1&C(k)>0.5,

c(k)=(C(k)-0.5)/0.5+3;

elifC(k)<=1.5&C(k)>1,

c(k)=(C(k)-1)/0.5+4;

elifC(k)<=2&C(k)>1.5,

c(k)=(C(k)-1.5)/0.5+5;

18

elifC(k)<=24.2&C(k)>2,

c(k)=(C(k)-2)/22.2+6;

end

end

c;

fort=1:28

s(t)=1/3*(a(t)+b(t)+c(t));

end

mean(s);

mean(a);

mean(b);

mean(c);

w=1:28;

subplot(4,1,1);

plot(w,a,'r');

subplot(4,1,2);

plot(w,b,'b');

subplot(4,1,3);

plot(w,c,'m');

subplot(4,1,4);

plot(w,s,'g');

不同情况对应不同输入数据。

19