基于藻类活动影响的湖库pH值建模

黄熠丽,崔婧嫄,陈思绪,张海平

(同济大学环境科学与工程学院,上海 200092)

pH是水溶液的一个重要特性参数,是影响水生生物繁殖,反映水质状况的一个综合指标。pH过低可能引发水体酸化等环境问题,过高则可能伴随水体富营养化等风险[1-2]。在饮用水处理和输送过程中,pH被认为是水处理化学中最重要的指标之一[3],其值高低和波动程度对管网输送、原水预处理、混凝沉淀、饮用水消毒等过程都有很大的影响。因此,建立湖库pH机理性模型,有助于深刻认识和描述pH与水质因子间的关系,探讨湖库pH的变化特征与影响因素;也有助于预测湖库pH变化趋势,指导后续水处理工艺中混凝剂的投加量,以节约经济成本;在湖库实测点位有限的情况下,还有助于进一步了解pH的空间分布,对于湖库特别是饮用水水源型水库的管理有着重要意义。

然而,由于pH的非保守性质[4],通过模拟水域生物、化学及物理过程以模拟pH的变化并不容易,目前应用于自然水体pH值建模与预测的成熟模型还相对较少。已有的模型可归纳为经验模型[5-7]、pH控制模型[8-10]、生物地球化学模型[11-14]3类。经验模型多基于pH与其他水质指标间的定性与定量关系建立经验公式,未能揭示pH的影响机理;pH控制模型则多用于养殖体系和实验室培养等人为控制下的水体,以监测、调节并预测水体pH,此类模型较少考虑水动力的对流扩散影响,无法模拟较大流域中pH的时空异质性;生物地球化学模型通常具有更复杂的模拟结构和较高的模拟精度[15],但此类模型通常需要大量的参数,在实际应用中,部分参数(如各类菌落数、磷酸盐浓度、铝离子浓度等)难以获得,使得其在一些场合不能适用,且目前的研究主要集中于海洋、酸水排放区等水域,其生态系统与湖库有较大区别。

在湖库条件下,pH的变化与藻类活动呈现出了显著的相关性,被认为是叶绿素a(Chl-a)浓度的被动因子[16]。但已有的研究[17]中,较少有进一步研究湖库净生产力并量化其对水体pH的作用,在模拟过程中或是简化了浮游植物的生化影响,或是通过复杂的模板模拟藻类生长、死亡、排泄等过程,而藻类实际的生物活动更为复杂,参数精度受限,从而影响pH模拟的结果。鉴于以上背景,本研究基于藻类活动对湖库pH的显著影响,利用MIKE 21建立水动力模型,以长期实测Chl-a浓度数据为基础,旨在提出一套实用型的湖库pH模拟方法,并具备以下特征:(1)模型参数易在湖库的日常管理中获得;(2)模型方法能够反映湖库存在的生物活动、水化学平衡及流体动力学过程;(3)模拟结果能够体现藻类活动和湖库流态对pH时空波动的影响。最后,将模拟方法置于华东地区两座重要饮用水水源型水库中进行验证,可为湖库的pH预警预报提供技术支撑,促进湖库的科学管理和供水安全保障。

1 模型方法与构建

1.1 湖库碳循环

研究[5]发现,湖库水体pH的动态平衡主要受游离CO2及碳酸盐的平衡系统影响,即碳循环是影响淡水湖库pH的主要框架,水体溶解无机碳(DIC)浓度与总碱度(TA)共同决定了水体pH的大小,这一过程如图1所示,涉及到浮游植物、浮游动物、沉积物等对象,它们通过不同的生物化学作用影响水中碳酸盐的含量。浮游植物在水体中的生命活动主要由藻类生长、新陈代谢,以及被浮游动物摄食、沉降和死亡这几部分组成;浮游动物的生命过程则主要包括摄食藻类、呼吸、排泄及死亡等过程;湖泊和底泥中的碎屑有机物、粪便排泄物在生物及微生物作用下会发生矿化降解[18];而在水气表面,气体受水体与大气的分压差溶解或逸出。这些过程都伴随着CO2的产生与消耗,引起水体DIC浓度的变化。

图1 湖库无机碳循环基本框架及其对pH值的影响Fig.1 Basic Framework of Inorganic Carbon Cycle in Lakes and/or Reservoirs and the Effect on pH Values

1.2 藻类活动对pH的影响

如前所述,CO2含量直接与水中浮游生物特别是水生植物的数量、生命活动的旺盛程度有关,而Chl-a浓度可用于表征水中藻类生物量[22]。因此,湖库pH的动态变化理论上和水体Chl-a浓度具备相关性。对华东地区某水库2020年的pH与Chl-a浓度建立回归关系,如图2所示。总体而言,pH随Chl-a浓度的增加而升高。

图2 华东地区某水库输水口pH值与Chl-a浓度相关性Fig.2 Correlation between pH Values and Chl-a Concentration at a Reservoir Outlet in East China

同时,在任一Chl-a浓度下,pH仍有较大幅度波动,表明pH同时还受其他因子影响,例如受水温和光照强度的规律性影响及水生植物初级生产的综合性影响,pH表现出明显的昼夜变化、季节循环以及与许多其他因素相关的强烈不规则波动[23-25],总体呈现为白天升高、夜间降低,春夏季高、秋冬季低的趋势[26-28]。水温主要从CO2溶解度对pH产生影响,随水温升高,CO2溶解度大幅降低,逸出量减少,pH变低,且温度升高同样提高了藻类的生命活性;光照强度则主要通过藻类光合作用对pH产生影响,白天藻类光合作用占据优势,水体pH升高,到傍晚达到最高值,夜间水生生物呼吸产生CO2,pH逐渐降低,到次日早晨达到最低值。因此,本文通过模拟藻类在不同光照、温度等情景下由于初级净生产力产生的CO2,同时考虑水-气交换、有机物、浮游动物等作用,模拟湖库pH变化过程。

1.3 CO2的模拟与计算

CO2的模拟是构建湖库pH模型的基础,本文应用MIKE 21软件中的ECOLab模板平台,编制CO2计算程序,并与水动力模型耦合。

CO2的模拟分为对流扩散、水-气交换和生化反应3个方面。对流扩散作用以水动力模型为基础,前期已通过对氯化物浓度的率定获得扩散系数取值;CO2在水-气界面的转移量采用大气碳通量(Cflux)进行计算;生化反应则主要考虑了藻类呼吸作用(CChla_respt)、藻类光合作用(CChla_phtsyn)、浮游动物呼吸作用(Czoop_respt)、有机物降解(CCODMn-CO2)4个过程,上述过程可表达为式(1)。

(1)

其中:CCO2——CO2质量浓度,mg/L;

τ——单位时间,d;

Cflux——CO2传输通量,mg/(L·d);

CChla_respt——藻类呼吸作用产生的CO2,mg/(L·d);

CChla_phtsyn——藻类光合作用消耗的CO2,mg/(L·d);

Czoop_respt——浮游动物呼吸作用产生的CO2,mg/(L·d);

CCODMn-CO2——有机物降解产生的CO2,mg/(L·d)。

1.3.1 水-气交换

基于1992年Wanninkhof[29]理论和2014年更新的Wanninkhof[30]理论,在水-大气接触面的CO2传输通量采用式(2)计算,其中0.010 56为单位转换因子。

Cflux=0.010 56kCO2K0(pCO2Air-pCO2Water)/d

(2)

其中:kCO2——CO2转移速率,cm/h;

K0——CO2的溶解度,mol/(L·atm),1 atm=101 325 Pa;

pCO2Air——大气的CO2达到平衡状态时的分压,设为恒定值380 μatm[31];

pCO2Water——水表面的CO2达到平衡状态时的分压,μatm,在给定盐度、温度、压力等环境条件下,可通过输入DIC和pH,在CO2SYS-Program程序中求解得到(返回单位为Pa);

d——网格水深,m。

在淡水条件下,CO2转移速率用式(3)～式(4)计算,对于风速为3～15 m/s有效,适用温度为-2～40 ℃,总体上不确定性约为20%[30]。

(3)

Sc=1 923.6-125.06t+4.377 3t2-
0.085 681t3+0.000 702 8t4

(4)

其中:u10——水面以上10 m处的风速,m/s;

Sc——施密特系数;

600——淡水在20 ℃下的施密特系数;

t——水温,℃。

在淡水条件下,CO2的溶解度K0通过式(5)计算。

(5)

其中:T——水温,K。

1.3.2 藻类呼吸作用

藻类呼吸作用产生的CO2通过式(6)进行计算。

(6)

其中:rresp——藻类呼吸速率,d-1;

tresp——呼吸作用的温度系数;

CChla——Chl-a质量浓度,mg/L,作为实测值输入到模型中;

CDO——DO质量浓度,mg/L,作为实测值输入到模型中;

Dresp——DO半饱和常数,mg/L;

pC/pChla——C/Chl-a,取33;

pCO2/pC——CO2/C,为3.67。

研究[32]表明,水温对控制藻类呼吸作用具有重要意义,被认为是控制季节性生长的主要因素。水温对于呼吸作用的影响通过引入温度校正系数实现,根据参考温度为20 ℃的Arrhenius表达式返回温度相关性。DO作为藻类新陈代谢的限制因子之一,可通过米氏方程量化其限制作用。

1.3.3 藻类光合作用

藻类光合作用消耗的CO2通过式(7)～式(8)进行计算。

(7)

(8)

其中:mave——日均藻类生长系数,d-1,通过式(8)进行计算;

tphtsyn——光合作用的温度系数;

Cphtsyn——CO2半饱和常数,mg/L;

mday——白天的藻类生长系数,d-1;

Tsunrise——当日日出时刻,h;

Tsunset——当日日落时刻,h。

光照对藻类光合作用有决定性的影响。本模型对光照强度做了简化处理,认为白天的光照强度为一恒定值,夜间为0。水温对于光合作用的影响同样通过温度校正系数实现,并通过基于CO2半饱和常数的米氏方程计算水体CO2浓度对于藻类光合作用的限制作用。

1.3.4 浮游动物呼吸作用

浮游动物呼吸作用产生的CO2通过式(9)进行计算。

Czoop_respt=rzoop×tzoopt-20×z

(9)

其中:rzoop——浮游动物呼吸作用产生CO2的速率,d-1;

tzoop——浮游动物呼吸作用的温度系数;

z——浮游动物与细菌生物量,mg/L,设为常数。

1.3.5 有机物降解

有机物降解产生的CO2通过式(10)进行计算。

(10)

其中:KCODMn——CODMn降解系数,d-1;

tCODMn——CODMn降解的温度系数;

CCODMn——CODMn质量浓度,mg/L,作为实测值输入到模型中,(CCODMn+1)表示将实测CODMn数据转化为TOC/(mg·L-1)的经验值;

DCODMn——CODMn降解的半速率反应常数,mg/L。

1.4 DIC与pH的模拟与计算

在平衡状态中,CO2进入水中的反应可通过式(11)～式(14)进行描述。

(11)

(12)

(13)

(14)

其中:K1——碳酸的一级解离常数,在淡水条件下受温度影响;

K2——碳酸的二级解离常数,在淡水条件下受温度影响;

[H+]——H+摩尔浓度,mol/L;

NDIC——DIC摩尔浓度,mol/L。

(15)

在淡水条件下,不考虑磷碱度及硅碱度,淡水碱度可通过式(16)～式(17)。

NTA=AC+[OH-]-[H+]

(16)

(17)

其中:NTA——TA摩尔浓度,mol/L;

AC——碳碱度的摩尔浓度,mol/L;

[OH-]——OH-摩尔浓度,mol/L。

结合上述公式,可以得到TA、DIC、[H+]与之间存在的关系式,如式(18)。

(18)

其中:Kw——水的离子积常数,在淡水条件下受温度影响。

2 模型应用

2.1 模拟水体

选取A、B两座华东地区饮用水源水库对本文所建立的pH模型进行验证。两座水库的水文水质特征有较大差别。A水库总面积达66.26 km2,设计有效库容为4.35亿m3,供水规模逾800万m3/d,库区水力停留时间较长,约为20 d,水流缓慢,局部存在藻类过度繁殖和聚集等问题,库内Chl-a质量浓度在2～41 μg/L,输水区pH值则明显大于上游来水,为8.0～9.0,pH偏高且波动剧烈,对于自然水体而言变化明显。B水库总面积为1.35 km2,设计有效库容为953万m3,供水能力逾200万m3/d,水体交换时间短,水力停留时间仅为3～5 d,整体流通性强,大型藻类、沉水植物难以生长,大部分区域浮游藻类质量浓度低,常年为1～5 μg/L,输水区pH值为7.8～8.4,处于相对较低的水平。

利用MIKE 21软件分别建立两座水库的二维水动力水质模型。A水库模型地形网格数量共102 826个,平均每个网格分辨率为644 m2,B水库模型地形网格数量共26 042个,平均每个网格分辨率为52 m2。利用模型计算每个网格的垂向平均水动力参数(水位、流量和流速等)和水质参数(CO2、DIC、pH等)。模型计算时间为2020年5月1日—2021年5月1日。

利用库区实测水位资料对模型进行率定,实测数据来自5 min/次的在线数据。A水库平均相对误差(MRE)为4.6%,均方根误差(RMSE)为0.16 m;B水库的MRE为0.73%,RMSE为0.06 m,表明所构建的水动力模型较好地模拟了水文水动力要素的主要变化特征,为后续pH建模提供了可靠基础。

图3 模拟与实测pH值比较(2020年5月—2021年5月)Fig.3 Comparison of Simulated and Measured pH Values (from May 2020 to May 2021)

2.2 模型条件设置与参数取值

根据水库上游库外的实测pH、TA、水温等数据,利用CO2SYS-Program程序计算上游库外的CO2浓度作为模型边界条件。库区风速、水深、水温、Chl-a浓度、CODMn、DO浓度、大气CO2分压、盐度和TA是模型输入条件。其中,水深数据从水动力模型读取;风速、水温、Chl-a浓度、CODMn、DO浓度为实测值;大气CO2分压设为恒定380 μatm。通过将实测日碱度输入和平均值输入进行测试比较发现,模型对这两座水库的碱度输入不敏感,结果相差较小,为简化起见,采用全年平均值输入到模型中参与计算,但这点是否适用于其他水库,还需进一步观察。

参数取值的合理性直接影响到模型对水体系统表征的可靠性与合理性[33-34],进而影响到模型精度。在本研究中,通过文献查阅及过往运行经验得出参数取值范围后,用实测数据对模型进行率定,得到的参数值如表1所示,与其他淡水湖库的研究结果总体一致[35-36]。A水库藻类活动强烈,实测得到的CODMn包含了藻类生长而产生的有机物,因此,在表观上几乎未发生降解,即综合降解系数为0。

表1 模型参数设置Tab.1 Setup of Model Parameter

2.3 模型结果

将模拟得到的输水口pH与实测数据进行比较,如图3所示。实测数据的监测频率为每天一次。A水库的MRE为1.02%,平均绝对误差(MAE)为0.086,RMSE为0.111;B水库的MRE为0.77%,MAE为0.062,RMSE为0.082。

本模型考虑了对流扩散、水-气交换、藻类活动等过程,模拟的总体结果显示,模型对于两座水库pH的高低及全年变化趋势都能有较好地重现与反演,且模拟趋势与实测趋势具有良好一致性,表明模型的设置总体合理。

2.4 模型验证

基于表1中的率定参数,选取2022年4月1日—6月30日的数据对模型进行验证,结果如图4所示。A水库的MRE为1.24%,MAE为0.106,RMSE为0.133;B水库的MRE为0.97%,MAE为0.081,RMSE为0.106,验证精度与率定模型相当。

图4 模拟与实测pH值比较(2022年4月—2022年6月)Fig.4 Comparison of Simulated and Measured pH Values (from April 2022 to June 2022)

3 结果分析与讨论

3.1 结果分析

从模拟结果来看,本文所构建的模型能够较好地重现pH在自然水体中的变化趋势,模型模拟结果与实测数据的MAE和RMSE均较低,表现出良好的模拟性能。所选择的模拟对象水体A、B水库的水文水质情况差异较大,而该模型均能较好地对其pH进行模拟,呈现出一定的普适性。

A水库的实测pH波动相对较大,模型对于pH起峰回落的过程及其所对应的时间响应良好,特别是pH谷值的模拟精度较高,在大多数情况下均较接近于实测值;pH峰值受藻类活性、气候条件、光照强度等环境因子的叠加影响,实际驱动作用可能比模型更为强烈,模拟值往往相较于实测值偏低,这与Cerco等[37]在富营养浅水海湾的模拟结果相符。

模型对相关数据获取困难的一些生化过程进行了简化描述,如将水体浮游动物量设为定值、不考虑优势藻种的季节性变化、底部沉积物溶解等,Liang等[17]的研究同样表明上述因素引起的pH变化并不明显,甚至可能增加模拟难度并降低精度。尽管在两座水库的模拟中未发现由此产生明显误差,将来在获得更多观测资料后,可进一步检验模型方法,优化模型结构。

3.2 参数敏感性分析

本文所构建的模型考虑了Chl-a浓度、DO浓度、水温、CODMn、碱度、风速、水深、光照、CO2浓度、pCO2等实测值或计算值,利用局部敏感性分析中的扰动法,量化所选参数对模拟结果的影响,其计算过程如式(19)。

(19)

其中:SS——参数敏感性数值;

ΔS——模拟结果的变化量,文中指ΔpH;

S——模拟结果,文中指pH;

ΔE——参数变化量,变化幅度设为合理范围内的±20%;

E——参数值。

计算得两座水库的敏感性数值及其绝对值排序如表2所示。

表2 参数敏感性及排序Tab.2 Parameter Sensitivity and Ranking

筛选出敏感性数值大于1%的参数,在A水库中表现为mday=光照时长>tresp>tphtsyn>tzoop>rzoop=zoop>Chl-a浓度>rresp,在B水库中表现为CO2浓度>tphtsyn>风速>水深>水温。总体而言,模型对Chl-a浓度、水温较为敏感,且部分因子呈现出异参同效的特征。Chl-a浓度及其相关参数(mday、rresp等)是影响湖库特别是富营养化湖库pH变化的重要驱动因子,其敏感性数值受到藻类指标浓度的高度影响[38];水温相关参数(tphtsyn、tresp等)敏感性也较高,且随温度季节性变化呈现出时间变异性特征[39];CO2作为边界条件,对B水库等面积较小、水动力较强的湖库影响较大。

3.3 结论

(1)本研究基于藻类活动对于湖库pH的显著影响,利用Chl-a浓度、水温、光照、风速、DO浓度、CODMn、CO2浓度、pCO2、TA等实测值或计算值,通过耦合水动力模型,考虑对流扩散作用、水-气交换、藻类呼吸作用、藻类光合作用、浮游动物呼吸作用、有机物降解等生态动力学过程,建立了适用于湖库的pH模型。

(2)所构建的模型应用于华东地区两座饮用水源水库。率定结果表明,pH模拟的MRE分别为1.02%和0.77%,MAE分别为0.086和0.062,RMSE分别为0.111和0.082,且模拟趋势与实测趋势具有良好一致性。模型验证精度基本相当,模拟精度较高,可较好反映湖库的pH动态变化过程。

(3)所模拟的两座水库存在较为不同的水文水质特征,而本研究所建立的模型均能较好地对其pH进行模拟,初步表明该模型具有华东地区湖库pH值建模的普适性。模型参数对于模型结果的影响较大,且参数敏感性随湖库水文水质特征表现出差异。同时,模型所需输入数据基本可从日常监测资料获取,模型的实用性较强。