BOD降解耗氧系数K1估算方法的程序实现

孟冲冲，杨欢，孙鹏，苏平定

BOD降解耗氧系数K1估算方法的程序实现

孟冲冲1，杨欢2，孙鹏1，苏平定2

（1.西安景天水利水电勘测设计咨询有限公司，陕西西安710000）；（2.陕西省渭南市水利水电勘测设计院，陕西渭南714000）

对BOD实验室测定值估算降解耗氧系数K1的两种方法最小二乘法和单变量求极值法进行了程序实现，并通过实例进行了验证。结果表明：两个程序均可以完成降解耗氧系数K1的估算，并讨论了最小二乘法不同初值k对结果的影响，k在一定的范围内，可以得出正确的结果，在这个范围之外则结果无效。采用最小二乘法进行降解系数推求时，初值k选定为0.5d-1结果较为合理，这两种方法对初值的要求不同，单变量极值法的0.618法，对于初值区间的选定需要进一步分析研究。

BOD降解耗氧系数K1；最小二乘法；单变量求极值法；0.618法

目前BOD降解耗氧系数的估算一般采用两种方法：即由BOD的实验室资料估算和由野外水质监测资料估算[1]。目前由BOD的实验室资料估算K1值的方法有很多：最小二乘法、差分法、斜率法、图解法、矩量法[2]、系列法[3]、单变量求极值法[4]等。差分法、斜率法、图解法计算简单，但精度不高，目前已很少使用；矩量法虽然计算简单，但只能对特定时间序列的数据进行估算，适用范围小；系列法估算过程比较繁琐，也很少使用。而最小二乘法和单变量求极值法精度较高，目前较为常用，但其计算复杂，必须借助于计算机，为便于应用，本文进行了这两种方法的程序实现，并进行了验证。

1 Theriaut的最小二乘法估算BOD降解耗氧系数K1

1.1 方法介绍[1]

降解耗氧系数K1的大小实际上体现了水中耗氧微生物对有机物降解转化的速度，因此影响微生物活性的许多因素都会直接、间接第影响K1值。如水中污染物的特性，包括污染物成分、组成比例、浓度大小等，还有微生物生存环境特性，包括水文、pH值、溶解氧浓度、水动力学因素、悬浮物等。

BOD降解过程符合反应动力学方程y=L0（1-e-kit），其中y为消耗掉的氧量，L0为初始时的BOD。

塞里奥特（Theriaut）提出把上式中的K1表达为：k1=k'1+h，K1是BOD的降解系数，是待估的K1的近似值，h是K1的允许偏差值，根据精度要求定出，因h甚小，exp(-ht)可近似为（1-ht），故反应动力学方程可以变形为：

其中，a=L0,d=L0h,f1=1-exp（-K'1t）,f2=t[ exp（-K'1t ）］

按最小二乘法求解a和b，得：

先假定一个K'1,求出a、d，并解得h=d/L0进行判别，一般精度为0.0001d-1，若h≥0.0001d-1，则令K'1=K'1.0+h，返回重新计算a、d、h，直至h＜0.0001d-1，此时的L0=a，K1=K'1+h，即所得结果。

最小二乘法的实质是对反应动力学方程y=L0（1-e-kit）求离差平方和的极小值,即:

图1 最小二乘法程序图

1.2 最小二乘法程序实现

2 单变量求极值法

2.1 方法介绍[4]

单变量求极值法是对最小二乘法的目标函数对L0求偏导数，并令其为零。可得：

再带回原目标函数得：然后用0.618法对其进行搜索。

选出K值可能范围(K1，K4)。

（1）在该区间中选出两点：K2=K1+(K4-K1)×0.382、K3=K1+ (K4-K1)×0.618；

（2）当目标函数值f(K2)＜f(K3)时,则区间缩小为(K1,K3)，回到第一步;若f(K2)＞f(K3)时，则区间为(K2,K4)，回到第一步。

如此重复迭代计算,直到缩小后区间(K1，K4)满足一定精度为止，即K4-K1＜ε,ε为允许误差,视计算精度要求而定,为了和最小二乘法比较，取ε=0.0001。

2.2 程序实现

单变量极值法的0.618法的程序图见最后页图2。

3 程序验证

本文选取文献[4]中的BOD实验数据(见表1)，用最小二乘法程序求解耗氧系数K1及初始BOD浓度。

表1 BOD实验数据

图20 .618法程序图

3.1 Theriaut的最小二乘法

赋值：k=0.5;t=[1,2,3,4,5,6,7,8,9,10];

y=[2.0567;3.6904;4.9881;6.0189;6.8377;

7.4881 ;8.0047;8.4151;8.7411;9.0]；

[L0,K1]=Zxec[k,t,y]

输出结果为：L0=10.00005,K1=0.23025

其中：k为BOD降解系数的初值，t为培养时间，y为与时间对应的BOD浓度mg/L，L0为初始时的BOD。（为便于比较，k值取与文献[4]相同）

3.2 单变量求极值法程序计算

赋值：m=0.2;n=0.3;t=[1,2,3,4,5,6,7,8,9,10];

y=[2.0567;3.6904;4.9881;6.0189;6.8377;7.4881;8.0047;8.4151;8.7411;9.0];

输出结果为：L0=10.00004,K1=0.23026。

表2 方法结果比较s

表4 BOD实验数据

为便于区分降解系数K1计算结果与初值，程序中将K值可能范围(K1，K4)设定为（m，n）。t，y含义同上。

3.3 程序计算结果与参考文献结果比较

通过上表可以看出，程序实现的结果与参考文献一致，说明两种方法的计算程序完全满足计算要求，可以进行降解系数的推求。

4 问题讨论

4.1 初值K1的选取对结果的影响

试讨论将初值分别设为0.1～0.9，间隔为0.1进行计算，结果见表3。

通过表3可以看出，当初值为0.1～0.6这六个数时，得出的结果都是L0=10.000，K1=0.230，与前面的计算结果一致。但当初值为0.7、0.8、0.9这三个数时，计算结果显示为负值，显然违背客观现实。因此当采用Theriaut的最小二乘法的程序进行计算时，初值最大不得超过0.6。为进一步证实此结论，选取参考文献[1]中的一组实验数据再次进行验证，实验数据及结果见表4、表5。

计算结果显示，当初值为0.2～0.8时，推求的结果一致，并与参考文献相同，但初值k选定为0.1和0.9时，计算结果为负值，结果不合理。分析表1和表5的BOD实验数据，其值差别比较大，出现计算结果为负值的初值设定也不一样，因此建议，采用最小二乘法进行降解系数推求时，k初值选定为0.5d-1结果较为合理。

4.2 降解系数K1值的讨论

不同河流的特征不同、影响耗氧速率的因素又及其复杂，致使耗氧系数K1值千差万别，小者仅为0.01d-1，大者可达5.20d-1[5]。同时K1值的大小还与污水初始浓度L0有关系。日本南部等学者实验结果表明[6]：在0＜L0＜40mg/L的范围内，K1随L0增大而有所增大；大过40mg/L后，K1值基本为一固定数。因此计算结果的准确性可依此做参考。

5 结论

（1）分别用最小二乘法和单变量极值法的0.618法MATLAB程序验证了实验数据，结果与参考文献一致，最小二乘法和单变量极值法的0.618法MATLAB程序可以用来进行降解系数K1的估算；

（2）采用最小二乘法进行降解系数推求时，通过选取不同初值k的计算结果对比，得出初值k选定为0.5d-1结果较为合理；

（3）用最小二乘法MATLAB程序推求出降解系数K1值后，可以用污水初始BOD浓度L0作参考来判断其准确性。

[1]雒文生,李怀恩.水环境保护[M].北京:中国水利水电出版社,2009.

[2]Moore E W,Thomas H A,Snow W B.Simplified method for analysis of BOD data[J].Sewage Ind.Wastes,1950,22(10):1343-1355.

[3]Marske D M,Polkowski L B.Evaluation of methods for estimating biochemical oxygen demand parameters[J].Water Pollution Control Federation,1972,44,(10):1987-2000.

[4]裘庆芗,施诊.降解系数K1的计算方法[J].水资源保护,1987(1):28-31.

[5]马巍,雒文生,张金存.水体耗氧系数与CBOD/NBOD相关关系的研究[J].水电能源科学,2000,18(4):49-50.

[6]雒文生,宋星原.水环境分析及预测[M].武汉:武汉大学出版社,2000.

Program Implementation on Estimating Method of the Degradation of Oxygen Consumption Coefficient K1

Meng Chongchong,Yang Huan,Sun Peng,Su Pingding
（Xi'an Jingtian Water Conservancy and Hydropower Survey and Design Consulting Co.,Ltd Xi’an 710000,Shaanxi；Weinan Water Conservancy and Hydropower Survey and Design Institute,Weinan 714000,Shaanxi）

Two methods of estimating the oxygen consumption coefficient K1in the BOD laboratory,the least square method and the single variable extreme value method are tobe implemented and tested by the actual example.The results show both procedures can complete the estimation of the degradation oxygen consumption coefficient K1and discuss the effect of the least squares initial value k on the results.K In a certain range,the correct result can be obtained,and the result is invalid outside this range The When the least squares method is used to derive the degradation coefficient,the initial value k is chosen to be 0.5d-1.The two methods have different requirements on the initial value.The 0.618 method of the univariate extreme value method is chosen for the initial value interval.Need tofurther analysis ofthe study.

BOD the degradation of oxygen consumption coefficient K1；he least square method and the single variable extreme value method and the method of0.618

表3 不同初值计算结果(1)

表5 不同初值计算结果表(2)

X8；TP311

A

1673-9000（2017）02-0100-03

2016-01-16

孟冲冲（1990-），男，河北邢台人，助理工程师，主要从事水利水电设计工作。