1 引言(Introduction)

碳是组成生命体的核心元素, 也是物质循环和能量流动的枢纽元素, 在水体生态系统中发挥着重要作用.河口有机碳循环作为全球碳循环的组成部分, 是碳循环研究中的一个重要方向(Cole et al., 2007; Regnier et al., 2013); 另外, 有机碳复杂多样的生化过程影响着河口生态环境的稳定性, 例如, 有机碳降解时会消耗大量的溶解氧, 甚至会造成河口水体出现缺氧(林晶, 2007; 张恒, 2009), 而溶解氧是海洋生物和生物地球化学过程的基础, 其含量的变动能引起河口、海洋生产力、生物多样性和生物地球化学循环的重大变化(Breitburg et al., 2018).因此, 研究河口有机碳的分布、迁移及转化机理, 对深入研究河口有机碳循环和河口缺氧等生态环境问题有着十分重要的意义.

不同河口有机碳的分布特征及迁移转化规律不尽相同, 例如, 长江口和Amazon河口有机碳分布主要受陆源有机碳输送的影响(Wu et al., 2007; Bouchez et al., 2014), 美国东海岸的York和Parker河口有机碳分布主要受藻类有机碳(即海源有机碳)的影响(Brucepeterson, 1994; Raymond et al., 2001).珠江口位于中国南海北部, 所在区域经济发达、人口众多、环境问题突出.珠江口三江汇流, 八口入海, 年径流量高达3.3×10¹¹ m³, 年输沙量8.87×10⁷ t (Shi et al., 2017), 每年有9.2×10⁵ t左右的有机碳由珠江河网通过珠江口输入南海(Ni et al., 2008).关于珠江口有机碳的研究在20世纪80年代已有开展, 相关研究表明, 珠江口有机碳浓度分布特征主要表现为口门高于近海, 靠近口门的上游水域有机碳主要受陆源输送的影响(蔡艳雅等, 1990; 刘庆霞, 2012; 郭威等, 2016; 林建荣, 2007), 而近海水域则主要受藻类有机碳的影响(刘庆霞, 2012).近年来, 珠江口水域环境发生着巨大变化, 例如, 上游修建大型水库拦沙截流, 导致珠江口有机碳输入减少，特别是颗粒态有机碳(Dai et al., 2008; Zhang et al., 2008), 再加上工、农、养殖业的发展使得珠江口富营养化进一步加重, 促进了藻类的生长(Harrison et al., 2008; Lu et al., 2015), 海源有机碳贡献逐步提升(刘庆霞, 2012).这些环境因素的改变会影响珠江口有机碳的来源和分布, 进而影响珠江口有机碳的生物地球化学过程.以往研究对珠江口有机碳的分布特征及其与环境因子之间的关系进行了详细探究, 但缺乏有机碳物理、生物化学过程等方面的定量化研究, 同时也未对珠江口有机碳进行过收支研究, 有关珠江口物理和生物化学过程对有机碳分布的影响程度及作用机制仍不清晰.

针对现有研究存在的不足, 本文通过构建珠江口三维水动力-水质模型对该海域有机碳的分布特征及迁移转化过程进行模拟研究, 并结合有机碳收支分析, 定量化研究珠江口物理过程与各生物化学过程对有机碳分布的贡献, 以揭示其中的主导动力学过程.

2 模型介绍(Model introduction) 2.1 水动力模型

珠江是中国第二大河流, 由西江、北江、东江及珠江三角洲诸河汇聚而成的复合水系, 水动力过程受到河网径流和南海潮汐双重作用的影响(刘德洪等, 2016).本文将珠江河网、河口及近岸海域视为一个整体系统, 采用前期已经过验证和应用的一维河网与三维河口耦合的水动力模型(胡嘉镗等, 2008)对珠江口复杂的水动力过程进行动态模拟, 研究区域见图 1a, 计算范围见图 1b.

一维河网水动力模型基于一维圣维南方程组, 采用三级联合解法进行求解, 共概化了299条河道及1726个断面.三维河口水动力模型采用ECOMSED水动力模型(Blumberg, 2002), 模拟变量主要包括水位、温度、盐度、流场等, 模型基本方程主要包括连续性方程、Navier-Stokes运动方程、温度和盐度守恒方程及海水状态方程等, 采用模态分裂技术进行求解.三维模型计算范围如图 1b所示, 西起镇海湾, 东至红海湾, 南至外海75 m等深线附近水域, 水平方向划分正交曲线网格186×183(南北向×东西向), 垂向上分为17个sigma层, 网格平面空间分辨率从河口内400~500 m逐渐过渡到外海3~4 km.一维河网模型和三维河口模型以8个口门作为连接断面, 通过互赋边界进行耦合计算(胡嘉镗等, 2008).

2.2 水质模型

采用RCA三维水质模型对珠江口的水质变化过程进行模拟研究.RCA水质模型由美国Hydroqual公司开发, 是由AESOP箱式水质模型及美国环保署WASP模型发展而来的先进水质模型.RCA模型中设有氮、磷、硅等主要营养盐及藻类、叶绿素a、颗粒态有机碳(POC)、溶解态有机碳(DOC)等多个状态变量, 详细介绍见文献(Fitzpatrick, 2004).本文重点关注珠江口有机碳过程, 涉及的主要状态变量见表 1.模型中包含C、N、P、Si、溶解氧(DO)等多变量循环交互系统, 并耦合了SFM底泥通量模型(Di Toro, 2001)(图 2).SFM中底泥分为耗氧层和厌氧层, 水体中颗粒态有机物沉降到底泥后, 一部分被埋藏, 一部分通过再悬浮返回水体, 一部分经过复杂的生化反应后以无机盐或气体的形式返还水体, SFM将底泥和水体两大系统动态关联, 更好地还原出底层水环境的真实性.RCA三维水质模型中水质因子控制方程如下所示：

(1)

式中, C为水质因子浓度(mg·L^-1), t为时间(s), x、y、z分别代表纵向、横向、垂向3个坐标方向, E为扩散系数(m²·s^-1), U为流速(m·s^-1), S为水质因子源汇项(mg·L^-1·s^-1), W为外源输入项(mg·L^-1·s^-1).

RCA中碳循环模块共包括7种形态的有机碳, 通过藻类和整个系统耦合, 相关转化关系如下：藻类通过光合作用合成并分泌ExDOC(变量相关描述见表 1), 并在被捕食过程产生RPOC、LPOC、RDOC、LDOC; POC在重力作用下发生沉降, 同时经过水解作用转化为DOC, DOC有两种消耗途径, 一种是反硝化作用消耗部分LDOC, 另一种是氧化作用消耗部分RDOC、LDOC、ReDOC、ExDOC.各变量具体反应方程可见文献(Fitzpatrick, 2004).

3 模型设置与验证(Model settings and validation) 3.1 模型设置

本文的模拟时间段为2006年1—8月, 其中, 1—6月为模型预热阶段, 以使模型进入正式模拟分析时处于相对稳定的状态, 7—8月为模型验证时段; 模型验证后, 基于7月的模拟结果对珠江口夏季有机碳的分布特征及收支进行深入分析.一维河网模型上游边界条件采用高要、石角、博罗水文站的逐日实测流量及老鸦岗和石咀水文站的逐时水位观测数据, 三维河口模型外海水位边界采用全球潮汐预测模型(OTIS)在南海北部的优化子集数据, 外海温度、盐度边界参考大万山站的历史监测数据进行设定, 风速、风向等海表面驱动条件采用欧洲中期天气预报中心(ECMWF)的再分析气象数据.河流水质边界和外海部分水质边界主要参考南海分局在八大口门的逐月监测数据和相关的研究文献(Wang et al., 2017), 外海POC、DOC、藻类边界设置主要参考相关的观测研究(Hung et al., 2007; Chen et al., 2008; 吴凯, 2017), 底泥模块的设置主要参考刘德洪等(2016)的研究成果.另外, 水质模型中各个主要反应参数取值参考目前珠江口已有的应用成果, 具体可参见相关文献(刘德洪等, 2016; Wang et al., 2017), 此处不再赘述.一维模型的时间步长为5 s, 三维水动力模型的外模时间步长为5 s, 内模时间步长为30 s, 水质模型的计算步长设为60 s.

3.2 模型验证

水质验证数据为908专项在2006年7—8月珠江口夏季航次的水质监测数据, 包括温度、盐度、有机碳浓度(均为瞬时值).其中, 研究区域内的有机碳采样点分布见图 1a, 模型验证结果见图 3.

从验证结果来看, 温度、盐度的模拟验证效果较好(图 3a、3b), 模拟值与观测值相关系数均高于0.83, 归一化的均方根误差均小于0.56, 归一化的标准差均接近于1(图 3c); 珠江口有机碳模拟(除GC15、16、17站点外, 见图 3d)的均方根误差为1.18 mg·L^-1, 绝对误差为-0.27 mg·L^-1, 相对误差为-5.5%, 相对平均偏差为31.55%;GC15、16、17站点位于香港南侧的岛屿附近, 模拟结果低于观测值, 这可能是由于模型没有全面考虑香港南侧岛屿有机碳的点源排放所造成的.另外, 根据后续有机碳收支的研究需要, 本文结合前人调查结果, 补充了相关指标的对比验证, 包括珠江口水体DOC消耗速率(本文模拟值为0.026~0.039 mg·L^-1·d^-1, 文献(He et al., 2010)报道结果为0.012~0.036 mg·L^-1·d^-1)、初级生产力(PP)和叶绿素a、DOC、POC等水质指标(表 2).由表 2可知, 模拟结果与以往观测数据相近.总的来看, 水质模型可以反映出珠江口水质的动态变化过程, 模型验证效果良好, 模拟结果可用于后续的研究分析.需要说明的是, 本文采用的水动力模型在前期已经过验证, 本文不再赘述, 具体验证结果可参见文献(Wang et al., 2017).

4 模拟结果与分析(Results and discussion) 4.1 珠江口夏季有机碳浓度的分布特征 4.1.1 有机碳浓度的水平分布特征

珠江口夏季表、底层有机碳浓度分布如图 4a1、4a2所示.从图中可以看出, 珠江口夏季表、底层有机碳浓度在水平分布上较为相似, 均呈明显的空间变化, 与盐度表现出显著的负相关, 即从口门到外海随着盐度升高, 有机碳浓度逐渐降低, 这一规律与刘庆霞(2012)和郭威等(2016)的研究一致.整体来看, 表层有机碳平均浓度为2.42 mg·L^-1, 高于底层有机碳平均浓度(1.91 mg·L^-1).有机碳浓度高值区均集中在各大口门附近水域, 而这些水域由于悬浮泥沙浓度高, 不利于藻类生长(Lu et al., 2015), 相应的藻类生化产碳(即藻类分泌及被捕食过程产生的海源有机碳)也较少, 这说明珠江口夏季有机碳浓度高值区的水平分布主要受陆源有机碳输送的影响.

在表层, 冲淡水携带大量的陆源有机碳从东四口门流入, 而藻类由于高浓度悬浮泥沙的影响生长缓慢(Lu et al., 2015), 叶绿素a浓度较低(图 4c1), 生化产碳较少, 陆源有机碳的输送是造成内伶仃洋水域有机碳浓度明显较高的首要原因.随后冲淡水携带丰富的有机碳继续向近海岸扩展, 在科氏力的作用下, 冲淡水向西扩展明显, 并与西四口门输送的有机碳汇合, 在西四口门近海域也形成了有机碳浓度高值区; 与此同时, 冲淡水也提供了充足的营养盐, 该水域光照条件改善, 藻类开始迅速生长(Lu et al., 2015), 在磨刀门口海域和黄茅海及其西南沿岸形成了叶绿素a高值区(图 4c1), 藻类生化产碳增多, 因此, 与内伶仃洋相比海源有机碳贡献有所增加.另外, 受西南季风的影响, 冲淡水也会向东扩展, 部分有机碳及营养盐经冲淡水携带向香港岛、大亚湾以南海域输运, 其中, 营养盐输送促进了藻类的生长, 使该区域的叶绿素a浓度高于外海, 为海源有机碳的生产提供了物质基础.

在底层, 冲淡水携带大量有机碳向外输送, 在外海高盐水的顶托作用下(图 4b2), 输送范围有所收缩, 于内伶仃洋及西四口门附近水域形成有机碳浓度高值区(图 4a2).在内伶仃洋及西四口门以外水域, 由于有机碳和营养盐输运受阻, 该水域藻类少、叶绿素a浓度低(图 4c2), 不能产生更多的海源有机碳, 同时, 陆源有机碳又无法及时补充, 有机碳浓度普遍较低.

底表层有机碳浓度特征差异显著(图 4a3), 两者单因子方差分析的p值为0, 这一差异与底表层盐度差密切相关(图 4b3).Pearson相关性分析显示, 底、表层有机碳浓度差与盐度差相关性系数达-0.90, 呈显著负相关(p < 0.01).冲淡水携带陆源有机碳和藻类产生的海源有机碳由口门向近海扩展, 沿途会有部分有机碳在其自身重力及盐度锋面等复杂动力的作用下沉降(张恒, 2009; Wang et al., 2017)至底层水体(图 4d, 图中所展示的有机碳沉降速率是指单位时间单位面积上有机碳沉降至水体底部的质量), 引起底表层盐度差高值区的上游水域(即沉降速率较大的水域), 底层有机碳浓度稍高于表层; 在底表层盐度差高值区及下游水域, 由于表层冲淡水扩展范围大且表层叶绿素a浓度高, 造成底层有机碳浓度低于表层.

4.1.2 有机碳浓度的垂向分布特征

为了进一步了解珠江口夏季有机碳的垂向分布特征, 图 5展示了珠江口夏季两个代表剖面(位置见图 1a, 其中, 剖面1从口门沿西槽向外海延伸, 剖面2则横跨外伶仃洋东、西两侧海域, 并横穿磨刀门海域盐度锋区域)的有机碳、盐度和叶绿素a的垂向分布.从图 5a1可以看出, 珠江口向外有机碳浓度逐渐降低, 在22.5°N以北的水域, 主要由有机碳含量较高的冲淡水占据, 盐度垂向分布均匀(图 5a2), 水体混合条件较好, 有机碳垂向分布也较为均匀(图 5a1); 在22.5°N以南水域, 由于受高盐海水的顶托, 冲淡水扩展范围从上到下逐渐收缩, 相应的有机碳浓度从表层到底层逐渐降低; 从图 5a3中可以看出, 叶绿素a浓度表层高于底层, 加剧了表层有机碳浓度高于底层的现象; 另外, 22.5°N北侧叶绿素a浓度低于南侧, 这也说明在内伶仃洋水域, 泥沙遮光限制了藻类生长, 海源有机碳贡献较弱.

从图 5b1中可以看出, 有机碳高值区出现在113.0°~114.0° E之间的6 m以上水层, 结合剖面2所在的位置和图 5b3可以发现, 一方面有机碳浓度高值区正对西四口门, 口门径流携带有机碳含量较高的低盐水漂浮在有机碳含量较低的海水之上; 另一方面, 叶绿素a浓度表层高于底层, 相应的表层海源有机碳的产量也较高, 这都导致了表层有机碳浓度高于底层.此外, 从剖面2的盐度分布图(图 5b2)来看, 113.0°~114.0°E之间的水体层化不利于上下层水体有机碳及时、充分混合(Dong et al., 2004), 更加剧了表层有机碳浓度高于底层的现象.有机碳浓度分层和水体层化情况较为契合, 这与刘庆霞(2012)的研究结论相一致.

4.2 珠江口夏季有机碳的收支研究

在口门有机碳输入和藻类生化产碳的共同作用下, 有机碳浓度从口门到近海逐渐降低, 空间分布差异显著.为了进一步量化物理、生物化学过程对珠江口夏季有机碳分布的影响程度, 下面对珠江口夏季有机碳的收支进行探究.

珠江口水环境复杂, 不同区域的有机碳受物理、生化过程的影响有所不同, 故将研究区域划分为4个区域进行研究, 具体划分见图 6a.其中, 区域1(内伶仃洋)靠近口门, 主要受东四口门物质输送影响; 区域2(中伶仃洋)是河口到海洋过渡的混合区域; 区域3(外伶仃洋)主要受西四口门输入及珠江口近海环境影响, 由于区域3水域较广, 东、西两侧物理与生化因素差异较大, 故将区域3划分为东、西两个子区域进行有机碳收支的分析.珠江口有机碳收支可按物理过程、生物化学过程进行划分, 包括有机碳在水动力作用下的输送过程(动力输送)、有机碳的点源排放(外源输入)和有机碳的沉降过程(沉降)3个物理过程, 以及藻类被捕食产生有机碳的过程(捕食)、藻类分泌有机碳的过程(分泌)、有机碳的氧化降解过程(氧化)和有机碳的反硝化过程(反硝化)4个生物化学过程.

通过数值模拟, 计算出珠江口夏季有机碳收支情况如图 6所示.图 6a为珠江口夏季水量收支平衡图.由图可知, 区域1的水量输送由东四口门汇入的径流主导; 到了区域2, 在科氏力的作用下大部分自内伶仃洋下泄的水量向西南方向流入区域3西侧; 在区域3, 除了西四口门的径流作用之外, 潮流、季风等作用加强(胡嘉镗等, 2009), 区域3西侧约有62.83%的水量在径流和科氏力的影响下沿西南方向流入南海, 剩余37.17%的水量受西南季风的影响流入区域3东侧, 与自东南向西北侵入的外海水流一同补偿了区域3东侧巨大的入海流.

口门输入的有机碳在水流的携带下主要沿西南方向输入外海, 少部分有机碳在径流、潮流、季风的影响下向东北输运, 并一同和外海输入的有机碳通过断面6进入南海(图 6b).在断面4上, 与水通量输送不同的是, 西侧输出的水量远小于东侧, 但西侧的有机碳输送通量却大于东侧, 这主要是因为西侧的有机碳浓度高于东侧(图 5b1), 西侧水流携碳量更高造成的.

区域1临近东四口门, 强劲的径流输送大量有机碳进入珠江口, 而该区域营养盐丰富但光限制严重(Yin et al., 2008), 不利于藻类生长(图 4c1、4c2), 藻类生化产碳和口门输入的陆源有机碳相比微不足道.因此, 有机碳主要来源于口门有机碳的动力输送(图 6c), 占区域1有机碳输入总量的83.80%.有机碳在该区域也发生着复杂的物理、生化耗碳过程, 其中以物理沉降为主, 口门向外, 流速渐缓, 在重力及盐度锋的作用下有机碳大量沉降, 沉降通量高达977.55 t·d^-1, 占区域1有机碳输入总量的83.18%.另外, 区域1生化耗碳量为128.48 t·d^-1(氧化+反硝化), 生化产碳量为130.17 t·d^-1(捕食+分泌), 分别占该区域有机碳输入总量的10.93%和11.08%.由上可知, 区域1内各生物化学过程较为平衡, 有机碳的分布主要受物理因素调控.

区域2距口门稍远, 有机碳输送通量受径流的影响相对区域1有所减弱, 动力输送通量下降至656.05 t·d^-1, 占区域2有机碳输入总量的40.54%.另外, 泥沙在区域1大量沉降后使得水体透明度逐渐增大, 光限制减弱, 再加上珠江口营养盐丰富, 藻类开始迅速生长(Harrison et al., 2008), 图 4c1、4c2也能表明区域2中藻类生化产碳作用显著提高.由图 6d可知, 捕食产碳量为406.13 t·d^-1, 分泌产碳量为346.75 t·d^-1, 生化产碳总量为752.88 t·d^-1, 占区域2有机碳输入总量的46.52%, 比动力输送作用高出5.98%.生化耗碳过程远小于生化产碳过程, 占区域2有机碳输入总量的24.53%, 有机碳的物理沉降过程依然是最大的耗碳过程, 沉降通量为967.77 t·d^-1, 占区域2有机碳输入总量的59.80%.总的来说, 区域2有机碳的分布受物理和生化因素共同调控.

区域3西侧受西四口门及上游动力输送的双重影响, 有机碳动力输送通量有所提高(图 6e), 占区域3西侧有机碳输入总量的53.23%.区域3西侧光限制进一步减弱, 藻类产碳加强, 生化产碳量高达1036.20 t·d^-1, 占区域3西侧有机碳输入总量的46.77%, 但受西四口门有机碳输送的影响, 动力输送比生化产碳作用高出6.46%.同时, 生化耗碳过程仍小于生化产碳过程, 占区域3西侧有机碳输入总量的24.38%, 有机碳的物理沉降作用占区域3西侧有机碳输入总量的60.92%.区域3西侧有机碳的分布受物理和生化因素共同调控, 其中物理因素占优势.

在区域3东侧, 冲淡水携带的营养盐受外海海水稀释作用的影响, 浓度较低, 藻类生长减缓, 使其生化产碳与区域3西侧相比明显降低(图 6f).然而有趣的是, 由于区域3东侧有机碳浓度较低, 水流携碳量不高, 动力输送作用迅速减弱, 动力输送通量仅为154.84 t·d^-1, 占区域3东侧有机碳输入总量的14.89%, 反而生化产碳以76.71%的绝对优势控制着区域3东侧有机碳的输入.另外, 物理沉降和生化耗碳作用相当, 分别占区域3东侧有机碳输入总量的31.65%和38.86%.综合而言, 区域3东侧有机碳的分布主要受生化因素调控.

5 结论(Conclusions)

1) 珠江口夏季有机碳的空间分布特征较为明显, 水平方向上, 从口门到近海有机碳浓度逐渐降低且与盐度呈负相关, 有机碳表、底层平均浓度分别为2.42和1.91 mg·L^-1.在有机碳沉降速率大的近口门水域, 表层有机碳浓度低于底层, 其它水域表层有机碳平均浓度高于底层.垂直方向上, 有机碳的浓度分布和水体层化密切相关, 在水体层化严重的水域, 有机碳浓度从上到下逐步降低, 其它水域上、下层水体中有机碳的浓度相差不大.

2) 珠江口水环境复杂, 不同研究区域有机碳的物理、生化过程差异显著.内伶仃洋水域, 有机碳主要来源于东四口门的径流输送, 主要耗碳项为物理沉降, 有机碳分布特征受物理因素调控.中伶仃洋水域, 有机碳来源以生化产碳为主, 以动力输送为辅, 主要耗碳项为物理沉降, 有机碳分布特征受物理和生化因素共同调控.在外伶仃洋西侧水域, 有机碳主要来源于径流输送, 主要耗碳项仍为物理沉降, 有机碳分布特征受物理和生化因素共同调控, 其中, 物理因素占优势.在外伶仃洋东侧水域, 有机碳主要来源于生化产碳, 主要耗碳项为生化耗碳, 有机碳分布特征受生化因素调控.

3) 复杂的生化过程是珠江口有机碳重要的生物地球化学过程, 从口门到近海生化过程均为产碳过程, 并且产碳量逐渐提高.其中, 捕食产碳过程是珠江口各生化过程中最重要的产碳过程, 氧化耗碳过程是珠江口各生化过程中最重要的耗碳过程.