刘小川， 陈 琴， 周 亶， 储小雪， 颜家兴， 樊德强， 蒋东云， 农以宁

(1.桂林电子科技大学 电子工程与自动化学院, 广西 桂林 541010； 2.桂林电子科技大学 生命与环境科学学院, 广西 桂林 541010； 3.桂林市环境卫生管理处环境科学研究所, 广西 桂林 541010)

我国畜禽养殖场的规模化，随之带来养殖废弃物的不断增加，已逐渐成为主要农业污染源[1]。养殖过程产生的粪便、污水、病死尸体等一系列污染源若得不到有效处理，极有可能引发畜禽疾病，损害养殖者的经济利益[2]。畜禽养殖产生的污染还可能引发一系列社会问题[3]。厌氧消化作为一种污水处理手段成为目前国内外养殖废水处理技术的重要环节[4]。但是，相比其他处理方式，厌氧消化工艺过程更为复杂，目前缺乏完善的过程在线监测和自动控制技术[5]。

ADM1模型已经被广泛证明是厌氧消化工艺设计、工程运行结果模拟和预测的有效工具[7-10]。ADM1 是采用微分代数方程组(DAE)和微分方程组(DE)描述生化和物化过程的结构化模型：DAE模型包括26个动态浓度变量、19个生化动力学过程、3个气-液转换动力学过程以及8个隐式代数变量；DE模型包括32个动态浓度变量和6个附加酸-碱动力学过程。该模型对厌氧系统内的生化过程和物化过程进行了详细分析，明确划分了模型组分并建立了相应的反应动力学方程，从而实现了厌氧消化的可计量性[6]。目前，ADM1 在厌氧消化理论研究领域和厌氧处理实际应用方面引起了广泛关注，Bikash[7-8]等人利用ADM1模型模拟厌氧消化过程中的微量元素和产酸的过程，Andres[9-10]等人利用 ADM1分别模拟废水和固体废弃物的厌氧消化。Daniele和Fezzani[11-15]等人利用ADM1模型模拟不同混合基质的厌氧共消化过程。ADM1 模型能较好地模拟和预测不同厌氧工艺在不同运行工况下的运行状态和效果，但是在实际应用中，要经过简化、扩充或修正才能对不同实例进行模拟[16]。Xiaofei Zhao和Victor Rivera-Salvador[17-18]利用数学方法校准ADM1模型参数模拟厌氧消化，修改后的模型可以为厌氧工艺的设计、运行和优化控制提供理论指导和技术支持。

本项工作是基于ADM1建立农村有机废弃物(ROW)厌氧消化沼气产量模型，并通过中试设备的产气数据对进行校准。改进后的模型可以模拟不同投料方式的农村有机废弃物的厌氧消化产气规律。

1 材料与方法

1.1 中试规模厌氧消化反应器和发酵底物

厌氧消化反应器为LBP一体化中试系统[19]，主体为CSTR全混式搅拌反应器。发酵底物取自重庆璧山区喜观村示范点，主要构成为农业生产生活废弃物、杂草、牛粪和少量实验性碳源。LBP一体化中试反应系统的厌氧消化基质含固率为5%～10%。

1.2 序批式实验

中温条件下的序批式厌氧消化实验，使用2000 L带有定时搅拌的不锈钢反应器进行实验。使用自行研制的沼气监测模组测量反应器的沼气产量以及气体成分。厌氧消化在37°C±1°C下进行，按比例加入ROW浆液至反应器，有效容积2000 L,每3个小时搅拌30 min。依照沼气产气量和沼气成分的变化进行序批式投出料，每次出料约200 kg。

1.3 数据处理

Batstone[20]等人使用Xpr，Xch，Xli和XI作为主要ADM1的投入，但是对于混合基质的建模，使用Xc会使模型更容易，因为每一次调整基质都可以改变Kdis，fch，fpr，fli来实现，这个方法在Zaher和Espositoet[21，23]等人的研究中都有成功实施。本次实验使用Xc，Xpr，Xch，Xli作为ADM1主要的投料。实验使用的原料分为两个部分，第1部分是农村有机废弃物ROW(Rural Organic Waste)，对应的化学计量系数和动力学参数为Kdis_ROW，fch_ ROW，fpr_ ROW，fli_ ROW，第2部分是混合有实验性碳源等其他物质的混合物，对应的化学计量系数和动力学参数为Kdis_m，fch_ m，fpr_ m，fli_ m(具体定义见表1)。使用ADM1中微分方程组(DE)为本实验的平衡方程，其他参数的选择参考Batstone[6]等人的研究方法。ADM1微分方程使用MATLAB中的ODE23s求解器函数求解。

表1 模型中修改的变量和参数的名称

在平衡状态下，气体产量的主要是通过公式(1)方程进行计算，在平衡状态下，CH4，H2，CO2气液传输过程的动力学方程可以分别由亨利定律使用方程式(2)～(4)描述[6]：

(1)

ρT,CH4=kLaCH4(Sliq,CH4-16KH,CH4ρgas,CH4)

(2)

ρT,H2=kLaH2(Sliq,H2-16KH,H2ρgas,H2)

(3)

ρT,CO2=kLaCO2(Sliq,CO2-16KH,CO2ρgas,CO2)

(4)

式中：KH,CH4，KH,H2，KH,CO2为CH4，H2，CO2亨利定律平衡常数，kmol·m-310-5Pa-1；ρgas,CH4，ρgas,H2，ρgas,CO2为CH4，H2，CO2气体压力，105Pa；Sliq,CH4，Sliq,H2，Sliq,CO2为液相中CH4，H2和CO2的液体浓度，kgCOD·m-3；kLaCH4，kLaH2，kLaCO2为CH4，H2和CO2总质量传递系数与比传递面积的乘积，d-1；kLaCH4，kLaH2，kLaCO2取决于反应器尺寸，液体流量，气体流量和扩散率值等[24-25]。3个系数kLaCH4，kLaH2，kLaCO2分别用于CH4，H2和CO2。因为介质中气体的传质系数为与扩散率的平方根成正比，kLaCO2和kLaH2的计算公式如下来自kLaCH4：

(5)

(6)

式中：DCO2，DH2，DCH4从Stockle[24]等人的文章中得出分别是4.65，1.98，1.57。通过调整kLaCH4可以得到相应的kLaCO2和kLaH2。

1.4 模型评估

确定系数R2和均方根误差RMSE用于评估：

(1)R2也称为拟合优度指数，使用MATLAB2017软件确定的，R2值越接近1代表拟合效果越好[8]：

(2)RMSE是预测值和测量值之间的均方根误差RMSE值低表示拟合度更好[8]：

2 结果与讨论

2.1 实验结果与模型分析

1.2中提到的实验设备运行80天的CH4，O2，CO2浓度变化以及沼气产气量累积的曲线如图1～4所示，图3中的O2浓度曲线表明了整个阶段的氧气浓度变化，氧气浓度出现剧烈变化的点既是投料点。1～10天是厌氧消化的启动阶段，这时罐内的氧气浓度缓慢下降。在第8，15，22天进行第1，2，3次投料，从图1和图3中可以看出投料之后甲烷浓度迅速增加进入产甲烷阶段。这几次投料时虽然氧浓度变化剧烈，但是甲烷浓度均迅速恢复至正常产气的70%。厌氧消化处于ADM1这3个阶段(水解、酸化、产甲烷)中的产甲烷阶段，产甲烷菌的数量上升甲烷产量增加，在第32天进行第4次投料，从图1～图3中可以看出投料之后甲烷浓度也在增加，但是浓度不高，一直在40%,说明此时产甲烷菌的活性受到了抑制，可能是因为酸化导致产甲烷菌活性降低，从而导致产甲烷浓度下降。在43～50天时，进行了第5次投料，此后CH4的浓度升到70%后逐渐下降到50%，这时产甲烷菌属于活性最高的时候，处于大量产气阶段。

图1 CH4浓度占比曲线

图2 CO2浓度占比曲线

图3 O2浓度-投料操作关系曲线

图4 沼气累积产气量

Cruz[26，28]等人研究了不同的克服甲烷生成抑制策略，例如添加添加剂(例如痕量金属、生物炭、颗粒状活性炭)和厌氧共消化(AcoD)。这些研究说明共消化可以提高甲烷的产生效率。本实验采用混合投料的方式，农村有机废弃物含有农业生产生活废弃物、杂草、牛粪等，混合基质以农村有机废弃物为基础加入实验性碳源等物质提高甲烷产气效率。进行动力学分析时，Souza[29]等人认为动力学参数可以在默认值的25%和400%之内变化，本实验研究的动力学参数是与Xc，Xpr，Xch，Xli相关的农村有机废弃物的Kdis_ROW和fch_ ROW，fpr_ ROW，fli_ ROW以及与混合基质相关的Kdis_m和fch_ m，fpr_ m，fli_ m。

2.2 模型校准

以CH4，CO2产气量为模型目标产物，根据Batstone[6]等人的文章确定除Kdis，fch_ Xc，fpr_Xc，fli_ Xc外的其他化学计量数和动力学参数，在1～4次投料时使用的是农村有机废弃物，农村有机废弃物的Kdis_ROW和fch_ ROW，fpr_ ROW，fli_ ROW取自Batstone[6]等人的研究分别是0.4，0.2，0.2，0.25。第5次投料投入的是混合基质，根据Daniele[11]等人的研究，进行共发酵时化学计量参数和动力学的参数均会下降，参考Daniele[11]等人研究，设置混合基质的参数Kdis_m和fch_ m，fpr_ m，fli_ m为0.3，0.1，0.15，0.1。由En Shi和Poggio[8,10]等人的研究可知当kLaCH4的初始值为1.00 d-1，根据公式(5)(6)得出kLaH2和kLaCO2的值为1.72和1.12 d-1，由于kLa值并不固定并且受搅拌、温度和液体性质的影响变化非常大[24]，通过改变kLa值可以更好地拟合产气模型。根据CH4，CO2的单日产量调整kLaCH4，kLaH2，kLaCO2的值得出合适的模拟曲线。

2.3 模型验证

第1次投料，甲烷的日产量峰值为8.8 L·d-1。CO2的日产量峰值为3.3 L ·d-1。根据甲烷的测量值调整ADM1模型的kLaCH4值，调整时发现kLaCH4在0.2附近时接近测量值，所以调整kLaCH4在0.1～0.3之间间隔0.01取值，计算每1次取值后的模型产气，与实际产气量一起计算拟合优度指数R2，得出kLaCH4与R2关系如图5。

图5 不同kLaCH4的产气量与测量值的拟合度与拟合曲线

再通过MATLAB软件得出3阶拟合曲线，计算拟合曲线极大值得出R2最大的点的kLaCH4值为0.1885，根据Poggio[10]等人得出的kLaCH4，kLaH2，kLaCO2之间的关系，计算出kLaCH4，kLaH2，kLaCO2分别为0.1885，0.3242，0.2111。其他投料的kLaCH4，kLaH2，kLaCO2根据第1次计算规则计算，得出的结果见表2。前3次投料，甲烷的浓度均能上升到70%，说明厌氧消化进程顺利。第4次投料，甲烷的浓度一直保持在40%，产生了酸化问题，所以在调整pH值之后进行了第5次投料。第五次投料使用的是混合基质，改变ADM1模型中的Kdis_ROW,fch_ ROW,fpr_ ROW,fli_ ROW为Kdis_m,fch_m,fpr_ m,fli_ m。根据计算第1次投料时的计算规则，计算其他投料时刻的kLa值，见表2。从图6、图7、图8、图9可知,调整后ADM1模型的每日产气量与实际值变化趋势基本符合,R2在0.8341到0.8674之间。

图6 第5～40天CH4日产量和ADM1模拟的结果

图7 第60～80天CH4日产量和ADM1模拟的结果

图8 第5～40天CO2日产量和ADM1模拟的结果

图9 第60～80天CO2日产量和ADM1模拟的结果

表2 第2，3，4，5次投料的kLa值

由于发酵底物种类复杂,其中化学性质差异较大。不同物质的水解时间不同，存在一些物质溶胀、水解滞后的现象，导致发酵和产气的滞后。ADM1模型尽管涵盖了众多的物化和生化过程，但是仍然不能囊括所有实际物质在厌氧环境下的表现。对分散性较大的溶胀、水解和厌氧发酵滞后现象没有很好的调整手段，使得曲线拟合度仍然不够理想。图10～图13是CH4和CO2实际测量浓度和基于ADM1模型拟合结果的计算浓度比值。

图10 第5～40天CH4浓度比和ADM1模拟的结果

利用CH4和CO2在总产气量中各自所占比值计算，比值浓度的拟合效果一般使用偏离实际值的误差评估[9]。本实验使用均方根误差RMSE表示模拟值与实际值的偏差，图10和图12的RMSE在10左右，图11和图13的RMSE在14左右。后两个曲线的误差较大，其原因主要是在75天左右未新增发酵基质而出现了个每日产气量的峰值。由于投加的发酵底物中含有干杂草，其中纤维素厌氧发酵前所需的水解时间较长，约需20～30天。此时，纤维素的水解产物参与到产甲烷的生化反应中，这就产生了75天的产气峰值。由于此时没有投料和出料等操作，ADM1模型对此种情况没有应对的变量，模型中也不能体现关联，导致计算值曲线在此处没有出现相应的峰值，图11和图13的RMSE值较大。

图11 第60～80天CH4浓度比和ADM1模拟的结果

图12 第5～40天CO2浓度比和ADM1模拟的结果

图13 第60～80天CO2浓度比和ADM1模拟的结果

表3 图6～9的每日产气量拟合度

表4 图10～13的每日产气量浓度比的均方根误差

2.4 累积产气量的ADM1拟合

通过校准沼气中的kLaCH4和CO2调整ADM1，图14～16显示了实际测量的累积CH4产气量、CO2产气量以及沼气产气量与ADM1模拟累积CH4产气量、CO2产气量以及沼气产气量。测量值与ADM1模拟曲线R2分别为0.9942，0.9854，0.9912，3条曲线的拟合效果令人满意。CO2的R2较低是因为投料时空气中的CO2进入导致拟合出现偏差，致使拟合度有下降，产甲烷效率增加，产气量大幅度提升。实际工程中，选择适当的投料策略，选取适当的ADM1参数，可以得到可靠的沼气产量模拟数据。

图14 CH4累积产量曲线和ADM1模拟曲线

图15 CO2累积产量曲线和ADM1模拟曲线

图16 沼气累积产量曲线和ADM1模拟曲线

3 结论

通过在中试规模LBP一体化设备上进行序批式农村生产生活废弃物混合厌氧消化中试实验，结合投料方式，通过校准kLa值模拟CH4,CO2产气量。得到的累积产气量模型拟合度较好，均高于0.98。可见校准后的ADM1模型能够反映序批式混合有机质厌氧消化的产气规律。但是，如果参与厌氧发酵的有机质底物成分复杂，在初始的溶胀水解阶段，不同生物质的降解速率各不相同，ADM1对此造成的发酵动力学差异并无校准方法，有待进一步研究。