当前水体富营养化现象已成为全球最受关注且亟待解决的水环境问题之一^[1-2].我国水库污染形势十分严峻，2000—2015年，水质劣于地表水Ⅲ类标准的水库比例由15 %上升至20 %，处于富营养状态的水库比例由30 %上升至37 % ^[3].经济发达地区的水库富营养化问题则更为突出^[4].处于富营养化状态的水库为藻类生长提供了适宜的水温、总氮(TN)、总磷(TP)等水质条件和缓慢的流场等水动力条件^[5-8]，常伴有水华暴发而严重影响周边人民的生活质量.因此研究水库水温、TN、TP的时空变化过程对科学预测水体富营养化及防控水华暴发有重要意义.

水库水质的时空分布特征受边界条件变化的直接影响^[9-11]，调水工程改变了其原有的进口水动力—水质条件，必然影响库区水质变化.我国调水工程众多，如“南水北调”、“引江济太”、“引江济巢”等，国内学者就输水对水库(湖泊)水质的影响进行了深入研究. Wang等^[12]、赵世新等^[13]研究了“南水北调”东线工程沿程水质变化，发现调水期上游骆马湖水质良好，下游输水河流遭受城市点源污染汇入水质变差进而影响南四湖入湖水质，经湖泊稀释净化后出湖水质则明显提升；Hu等^[14]、Li等^[15]及潘晓雪等^[16]就“引江济太”对太湖水质的影响研究表明，作为降低水华暴发风险的应急措施长江引水能增强太湖水体交换，有效降低部分湖区TN、叶绿素a浓度，对TP浓度却无明显影响，当引水水体中氮磷浓度未被降低至合理水平时，长期引水反而会增加太湖富营养化的风险；Chen等^[17]、谢兴勇等^[18-19]对“引江济巢”水动力—水质变化过程的数值模拟显示输水能有效削减巢湖TN、TP浓度峰值，降低湖区富营养化程度，且距江水入湖口越近水质改善效果越明显.

于桥水库是“引滦入津”工程中重要的调蓄水库，其水质优劣同样受进口边界-引滦输水影响严重.近年来，于桥水库水体富营养化程度加剧，且个别月份暴发轻微水华^[20-21]，严重威胁了天津市的供水安全.为研究于桥水库水质的时空变化特点，张庆强等^[22]对输水沿线水质进行调查，结果表明库区TN主要源于上游潘家口、大黑汀水库，TP受沿线污染物汇入影响较大；Chen等^[23]基于长时间序列的水质资料分析，发现于桥水库上游营养物浓度长期处于较高水平且水温具有下降趋势；葛宁等^[24]对于桥水库水质空间分布的数值模拟表明TN、TP由东南向西北扩散浓度逐渐降低；Liu等^[25]运用“水动力—水质”耦合数学模型发现果河入流及风生流是库区水质分布重要影响因素.然而果河在引滦输水期与非输水期水动力、水质条件差异明显，上述学者均未阐明不同时期果河入库水质差异及库区水质变化与入库流量的响应关系.

为研究引滦输水期与非输水期果河入库水质差异，本文基于2011—2015年实测资料(2016年未调水)，对入库水质(水温、TN、TP)与果河流量进行相关关系分析；为详尽地反映果河—库区水质变化随流量的动态响应过程，本文选取污染严重、流量变化剧烈的2012年作为典型年，运用平面二维“水动力—水质”数学模型模拟于桥水库水质的时空分布，对比分析了输水期与非输水期库区水质变化随果河流量的响应关系.研究成果可为改善于桥水库水质及降低水华暴发风险提供水质基础数据.

1 研究方法 1.1 流量监测与水质取样

于桥水库位于天津市蓟县城东(图 1a)，是天津市生活及工农业用水的重要水源地之一.其正常蓄水位21.16 m，蓄水面积113.8 km²，最大回水东西长30 km，南北宽8 km，平均水深4.6 m，总库容15.59亿m³.库区南岸地势陡峭，受原河道深泓控制水深大；北岸地势则相对平坦水深浅(图 1b).水库控制流域属温带大陆性季风性半湿润气候，年平均气温为10.4~11.5℃，多年平均降水量为748.5 mm^[23-24]，主要集中在6—9月；冬季气候寒冷干燥，水库进入结冰期.流域内主要入库河流为果河和淋河，分别位于水库东南岸和东北岸(图 1c).果河由沙河和黎河汇入而成，黎河为引滦输水通道，上游连接输水隧洞，受输水影响果河入库氮、磷浓度持续超标^[25]，库区水质远低于规划要求的Ⅲ类标准.

流量监测站及水质取样点如图 1c所示.沙河逐日流量监测站位于水平口，黎河引滦输水逐日流量监测站位于前毛庄，果河流量为黎河沙河流量之和. 2011—2015年果河多年平均径流量为6.77亿m³，淋河为0.21亿m³且年内常处于断流状态，研究时忽略淋河对库区水质的影响.非输水期果河流量主要来自于沙河天然径流，在2~15 m³/s范围内；输水期果河流量则主要由黎河引滦输水贡献，在55~75 m³/s范围内，汛期为确保防洪安全，引滦输水流量控制在30~40 m³/s，若遭遇洪水果河流量可达100~321 m³/s.沙河水质取样点位于沙河桥；黎河水质取样点位于黎河桥；入库水质取样点位于果河桥；库区水质取样点有3个：峰山南、库中心、放水洞.

除峰山南、放水洞每年仅监测5—10月水质外，其他各测点水质取样基本维持每月月初一次；2012—2014年引滦输水频繁，水质取样月中再增加一次.取样时天气状况良好，冬季封冻期取样于冰盖之下.水温采用JM222L型手持式数字温度计现场测定，河流水温监测点位于水下20 cm，库区水温监测点位于水下40 cm. TN、TP水质样品使用有机玻璃深水取样器按照《水质采样技术指导(HJ494—2009)》、《水质采样方案设计技术规定(HJ495—2009)》进行采集，并在水质样品中加入硫酸使得pH≤2固定保存；取样后用SKALAR SAN⁺⁺ SYSTEM连续流动分光光度仪进行分析，TN浓度采用连续流动-盐酸萘乙二胺分光光度法^[26]测定，TP浓度采用连续流动-钼酸铵分光光度法^[27]测定.

1.2 “水动力—水质”数学模型的建立 1.2.1 基本方程

基于MIKE21 FM，采用Boussinesq假定和流体静压假定下的平面二维不可压缩雷诺时均N-S方程，采用温度影响下的密度流模式，建立于桥水库“水动力—水质”数学模型.模型的控制方程由水流连续性方程、运动方程和水质对流扩散方程组成^[28-29]，在笛卡尔坐标系下，如下所述.

水流连续方程：

$ \frac{{\partial h}}{{\partial t}} + \frac{{\partial \left( {h\bar u} \right)}}{{\partial x}} + \frac{{\partial \left( {h\bar v} \right)}}{{\partial y}} = hS $

(1)

水流运动方程：

x方向：

$ \begin{array}{c} \frac{{\partial h\bar u}}{{\partial t}} + \frac{{\partial (h{{\bar u}^2})}}{{\partial x}} + \frac{{\partial \left( {h\overline {uv} } \right)}}{{\partial y}} = f\;\bar vh - gh\frac{{\partial \eta }}{{\partial x}} - \frac{h}{{{\rho _0}}}\frac{{\partial {p_a}}}{{\partial x}} - \frac{{g{h^2}}}{{2{\rho _0}}}\frac{{\partial \rho }}{{\partial x}} + \\ \frac{{{\tau _{{\rm{s}}x}}}}{{{\rho _0}}} - \frac{{{\tau _{{\rm{b}}x}}}}{{{\rho _0}}} + \frac{{\partial (h{T_{xx}})}}{{\partial x}} + \frac{{\partial (h{T_{xy}})}}{{\partial y}} + h{u_{\rm{s}}}S \end{array} $

(2)

y方向：

$ \begin{array}{c} \frac{{\partial h\bar v}}{{\partial t}} + \frac{{\partial \left( {h\overline {uv} } \right)}}{{\partial x}} + \frac{{\partial (h{{\bar v}^2})}}{{\partial y}} = f\;\bar uh - gh\frac{{\partial \eta }}{{\partial y}} - \frac{h}{{{\rho _0}}}\frac{{\partial {p_a}}}{{\partial y}} - \frac{{g{h^2}}}{{2{\rho _0}}}\frac{{\partial \rho }}{{\partial y}} + \\ \frac{{{\tau _{{\rm{s}}y}}}}{{{\rho _0}}} - \frac{{{\tau _{{\rm{b}}y}}}}{{{\rho _0}}} + \frac{{\partial (h{T_{xy}})}}{{\partial x}} + \frac{{\partial (h{T_{yy}})}}{{\partial y}} + h{v_{\rm{s}}}S \end{array} $

(3)

物质浓度对流扩散方程：

$ \frac{{\partial h\bar C}}{{\partial t}} + \frac{{\partial \left( {h\bar u\bar C} \right)}}{{\partial x}} + \frac{{\partial \left( {h\bar v\bar C} \right)}}{{\partial y}} = h{F_{\rm{C}}} - h{k_{\rm{p}}}\bar C + h{C_{\rm{S}}}S $

(4)

温度对流扩散方程：

$ \frac{{\partial h\bar T}}{{\partial t}} + \frac{{\partial \left( {h\bar u\bar T} \right)}}{{\partial x}} + \frac{{\partial \left( {h\bar v\bar T} \right)}}{{\partial y}} = h{F_{\rm{T}}} - h\dot H + h{T_{\rm{S}}}S $

(5)

式中，h为水深，m；u、v分别为x、y方向垂线平均流速，m/s；S为源汇项，s^-1；f为地转偏向参数，s^-1；η为水面高程，m；ρ为水体密度，kg/m³；ρ₀为水体参照密度，kg/m³；g为重力加速度，m/s²；p_a为当地大气压强，Pa；τ_s为表面切应力，N/m²；τ_b为底部切应力，N/m²；T_ij为横向应力项，m²/s²；u_s、v_s分别为源汇项x、y方向流速，m/s；C为垂线平均物质浓度，mg/L；F_C为浓度水平扩散项，mg/(L ·s)；k_p为衰减系数，s^-1；C_S为浓度源汇项，mg/L；T为垂线平均温度，℃；F_T为温度水平扩散项，℃/s；Ḣ为大气热交换项，℃/s；T_S为源汇项温度，℃.

MIKE21 FM采用非结构化网格，求解方法为单元中心的有限体积法，控制方程离散时，结果变量u、v位于单元中心，跨边界通量垂直于单元边.计算出每个控制体边界沿法向输入(出)的流量和动量通量后，对每个控制体分别进行水量和动量平衡计算，得到计算时段末各控制体的平均水深和流速，再由多个控制体的方程联合求解节点数据，具体参见MIKE21 FM科学手册^[28].

1.2.2 模型构建

参照于桥水库遥感影像，基于库区1:25000地形图，将果河桥作为进口边界，坝前放水洞作为出口边界，以正常蓄水位21.16 m淹没区域作为模型计算范围.整个区域采用三角形网格进行划分，如图 1b所示，并对果河段及进出口地形高程变化剧烈区域进行网格加密.经测试，果河段网格边长控制在20~50 m，库区网格边长控制在80~150 m时模型能准确刻画计算区域地形变化并且计算耗时可以接受，此时共有44085个网格单元，22739个网格节点.模型进口边界采用逐日实测流量数据及半月一次实测水温、TN浓度、TP浓度变化过程；出口边界采用逐日实测水位过程及零梯度水质出口条件.考虑水体与大气的热交换过程，大气温度、相对湿度采用遵化站实测数据，辐射强度采用经验模式^[29]根据于桥水库经纬度及气象条件进行计算，显热及潜热过程采用默认参数设置；并考虑风对库区水质变化的影响，采用蓟县站2012年实测风速数据，以西向风为主，平均风速1.6 m/s.初始水位为3月31日坝前实测水位19.3 m；初始地形采用北京54坐标、大沽高程下的于桥水库实测地形图；初始温度采用4月初实测水温；采用试算15天后的TN、TP浓度分布结果作为水质初始条件.

1.2.3 模型的率定与验证

由于模拟耗时长，模型率定时段取2012年4—7月，该时段包含两次输水期与两次非输水期，有利于提高模型在不同水动力条件下参数的率定效率.模型验证时段选为同年8—10月，该时段包含输水期与非输水期各一次.当糙率在0.015~0.025之间，水平涡黏度采用Smagorinsky公式计算，温度水平扩散系采用水平涡黏度比拟格式，比拟系数为2.5，TN、TP水平扩散系数为0.015 m²/s，TN衰减系数为0.006~0.015 d^-1，TP衰减系数为0.008~0.010 d^-1时，模型能很好地还原于桥水库水质的时空变化过程.下文给出了模型参数率定后，2012年4—10月模拟结果与实测数据的对比.

根据模拟结果果河段流速在0.06 m/s左右；进口段过水断面扩大流速放缓为0.04 m/s；库区流速均小于0.02 m/s.无论输水期还是非输水期，于桥水库整个库区流速均较为缓慢，且不存在明显的水面比降，因此可仅采用库中心实测水位对模型的水动力模块进行校核. 图 2给出了库中心水位模拟值与实测值之间的对比，平均绝对误差为0.005 m，平均相对误差为0.02 %，均方根误差(Root Mean Squared Error, RMSE)为0.02 m，纳什效率系数(Nash-Sutcliffe Efficiency coefficient, NSE)为0.99，模拟结果良好.

采用峰山南、库中心、放水洞实测水温、TN、TP浓度过程对模型的水质模块进行校核.水温模拟值与实测值之间具有一致的变化趋势(图 3)，三监测点平均绝对误差为1.2~1.5℃，平均相对误差为6 % ~7 %，RMSE为1.4~1.7℃，NSE为0.89~0.93(表 1)；由于气温存在日内差异，库区水温监测并非日内的同一时刻，模拟值与实测值出现上述误差是可以接受的. TN、TP浓度模拟值与实测值之间变化趋势一致(图 4)，TN浓度平均绝对误差0.32~0.45 mg/L，平均相对误差14 % ~18 %，RMSE为0.37~0.58 mg/L，NSE为0.73-0.78；TP浓度平均绝对误差0.005~0.008 mg/L，平均相对误差10 % ~20 %，RMSE为0.006~0.009 mg/L，NSE为0.63~0.65(表 1).水质进口边界数据频率为半月一次，不能完全真实地反映入库水质状态，数值模拟结果出现上述误差是可以接受的. 8月模拟误差最大，这与模型未考虑淋河入流对库区水质的影响有关，期间淋河遭遇了流量达66.2 m³/s的洪水过程必然影响库区水质分布.总体来看模拟误差在可接受的范围之内，数学模型能弥补实测资料空间跨度大采样频率低的不足，更详尽地反映于桥水库水质的时空变化；后文研究输水对库区水质分布的影响也只针对2012年，此时段内数值模拟结果具有足够的可靠性.

2 结果 2.1 输水对入库水质的影响 2.1.1 入库水温变化

从果河桥及库中心水温随果河流量的时间变化过程(图 5)可以看出：1)11月—次年4月，果河桥及库中心水温均受当地气温控制，年际差异不大，随季节变化呈先降低后升高的趋势；冬季封冻期水体取样于冰盖之下，受冰盖隔绝作用^[30]，水温明显高于当地气温；2)5—10月，非输水期，果河流量在2~15 m³/s范围内，果河桥水温略低于库中心水温；输水期，果河流量骤增为55~75 m³/s，引滦输水取代气温成为果河桥水温的重要影响因素，使得果河桥水温明显降低，与库中心水温差异增大，遭遇洪峰时水温差异则更大；3)输水期库中心水温略有下降，洪峰流量达300 m³/s以上时，库中心水温明显降低；且2013—2015年输水期变长使得库中心年均水温较2011—2012年下降1.6℃. Chen等^[23]基于1989—2009年的实测资料，同样发现库区上游水温具有下降趋势.

为进一步研究5—10月输水期果河桥与库中心的水温差异，定义入库温差ΔT=库中心水温—果河桥水温，根据2011—2015年的实测数据，对入库温差ΔT与果河流量进行回归分析(图 6a)，发现入库温差ΔT随果河流量增加呈显著的线性增大趋势.为解释这一现象，对黎河上游引滦输水隧洞5—10月的水温进行调查(图 6b，部分月份数据缺失)，得知：输水期，隧洞出流水温在平均值13.8℃上下波动，其远低于库中心水温(平均值25.7℃).因此5—10月受黎河上游隧洞引水影响，输水使得果河桥水温降低，果河流量越大则降温效果愈明显，入库温差ΔT愈大.

2.1.2 入库TN、TP浓度变化

从果河桥及库中心TN浓度随果河流量的时间变化过程(图 7)可知：1)无论输水期还是非输水期，果河桥TN浓度均远高于库中心TN浓度，果河入流必然使得库区水质变差；2)非输水期，果河桥TN浓度维持较高水平，平均值为8.16 mg/L；输水期，果河流量增大使得TN浓度急剧降低，平均值降为5.21 mg/L，并且流量越大，果河桥TN浓度降低越明显；3)非输水期，库中心TN浓度受降解作用而呈衰减势态；输水期，库中心TN浓度受果河入流影响而升高，由于库区流速放缓使得库中心TN浓度变化明显滞后于果河流量变化.

为进一步研究入库TN浓度变化与果河流量的关系，考虑到水质随流量变化的滞后效应，并排除汛期降雨带来面源污染的情况，对果河桥当日TN浓度与果河前七日平均流量进行回归分析(图 8a)，发现果河桥TN浓度随流量增加呈显著的幂函数降低趋势.为解释这一现象，对果河上游来流沙河、黎河的TN浓度进行分析(图 8b)，得知：无论输水期还是非输水期，沙河TN浓度皆在平均值($\overline {{C_{沙}}} $ =11.33 mg/L)上下波动；而输水期黎河TN浓度在平均值($\overline {{C_{黎}}} $ =4.78 mg/L)上下波动，远低于沙河TN浓度.假定沙河来流处于稳态，流量为2011—2015年平均流量$\overline {{Q_{沙}}} $ =3.83 m³/s，TN浓度为$\overline {{C_{沙}}} $ =11.33 mg/L，黎河来流TN浓度为$\overline {{C_{黎}}} $ =4.78 mg/L；在不考虑营养物降解且果河水体由沙河、黎河充分混合而成的条件下，果河TN浓度${C_{{\rm{TN}}}} = \frac{{\overline {{C_{沙}}} \cdot \overline {{Q_{沙}}} + \overline {{C_{黎}}} \cdot {Q_{黎}}}}{{{Q_{果}}}} = \frac{{25.09}}{{{Q_{果}}}} + 4.78$，其计算值如图 8a所示，与实测值变化趋势一致.由于实际沙河、黎河TN浓度处于波动状态，所以果河桥实测TN浓度值也在理论计算值两侧波动.因此非输水期果河TN浓度主要受沙河来流影响处于较高水平；引滦输水期TN浓度低的黎河来流稀释了TN浓度高的沙河来流，导致果河入库TN浓度随流量增加呈幂函数降低趋势.

从2011—2015年果河桥及库中心实测TP浓度随果河流量的时间变化过程(图 9)可知：1)总体上，果河桥TP浓度明显高于库中心TP浓度，果河TP入流同TN一样也使得库区水质变差；2)与TN浓度变化趋势相反，非输水期，果河桥TP浓度维持较低水平，平均值为0.08 mg/L；输水期，果河流量增大使得TP浓度急剧升高，平均值增加到0.28 mg/L，并且引滦输水流量越大，果河TP浓度越高；3)非输水期，库中心TP浓度略有降低；输水期，库中心TP浓度受果河入流影响浓度明显升高，库区流速放缓同样使得库中心TP浓度变化明显滞后于果河流量变化.

采用与前文TN浓度同样的分析方式，对入库TP浓度与果河流量的关系进行回归分析(图 10a)，发现果河桥TP浓度随流量增加呈显著的线性升高趋势.对上游沙河、黎河来流的TP浓度进行调查(图 10b)，得知：无论输水期还是非输水期，沙河桥TP浓度均在平均值($\overline {{C_{沙}}} $ =0.10 mg/L)上下波动，远低于黎河桥TP浓度.而黎河桥TP浓度与沙河桥不同，其受引滦输水影响，随流量呈显著的线性升高关系(C_TP=0.004Q_黎+0.14).这与输水流量增大、河床底质扰动增强和TP沿程汇入增加有关^[22].在与TN相同的假定条件下，同样将沙河视为$\overline {{Q_{沙}}} $ =3.83 m³/s，$\overline {{C_{沙}}} $ =0.10 mg/L的稳态污染源，则果河桥TP浓度${C_{{\rm{TP}}}} = \frac{{\overline {{C_{沙}}} \cdot \overline {{Q_{沙}}} + {C_{黎}} \cdot {Q_{黎}}}}{{{Q_{果}}}} = 0.004{Q_{果}} + 0.109 - \frac{{0.09}}{{{Q_{果}}}}$.由于输水期果河流量均在20 m³/s以上，忽略掉尾项小量，则${C_{{\rm{TP}}}} = \frac{{\overline {{C_{沙}}} \cdot \overline {{Q_{沙}}} + {C_{黎}} \cdot {Q_{黎}}}}{{{Q_{果}}}} = 0.004{Q_{果}} + 0.109$(图 10a).果河桥实测TP浓度值同样也在理论计算值两侧波动.因此非输水期果河TP浓度主要由沙河控制，维持较低水平(平均0.08 mg/L)；引滦输水期TP浓度高的黎河来流成为影响果河TP浓度的主要因素，使得果河入库TP浓度随流量呈线性升高趋势.

2.2 输水对库区水质分布的影响 2.2.1 库区水温变化

对南岸深水区及北岸浅水区的水温进行分析，以研究输水对库区水温沿程分布的影响.从南北两岸(图 1b，起于果河桥讫于坝前)2012年6月输水期、7月非输水期、8月输水期水温模拟结果的平均值(图 11)可知：1)无论输水期还是非输水期，于桥水库南北两岸水温沿程(自果河桥至坝前)均呈先升高后降低的趋势；果河段受引滦输水隧洞低温水体影响水温最低；库区下游水深大，受水库放水影响水体混合更充分，而库区上游水深小，水体也趋近静止，接受相同的热量，库区上游则更易升温从而水温最高；2)非输水期7月，果河段与库区水温差异小，入库水温梯度平缓；输水期6月，果河平均流量为60.2 m³/s，引滦输水使得果河段水温迅速下降，入库温差增大，温度梯度陡增；输水期8月遭遇洪水，果河平均流量达92.4 m³/s，导致果河段水温降低更明显，甚至低于7月水温，并且库区上游也出现了降温现象，流量骤增使得入库温差更大、温度梯度更陡；3)相对于库区北岸，库区南岸地势低，更易遭受果河来流影响，因此输水期库区南岸降温效果更明显，6月和8月南岸库区上游平均水温较北岸分别低0.6和1.8℃.

以南岸深水区为例，研究库区水温变化与入库流量的响应关系，从沿南岸距果河桥15 km(位于进口段)、20 km(位于库区中间)、25 km(位于坝前)处水温模拟结果随果河流量的变化过程(图 12)可知：1)按照设计流量输水时(Q_果=55~75 m³/s)，距果河桥15 km处水温明显降低，水温谷值滞后于流量变化3天，而距果河桥20和25 km处水温几乎不受果河流量影响；2)汛期引滦输水流量控制在30~40 m³/s，果河遭遇洪水时(Q_果=100~321 m³/s)，距果河桥15和20 km处水温均明显降低且流量越大降温效果越明显，水温谷值仍滞后于果河流量变化且距离越远滞后时间越长，而25 km处水温受果河影响不明显.分析北岸水温与果河流量响应关系，发现具有同样的规律，这里不再赘述.

2.2.2 库区TN、TP浓度变化

同样对南北两岸的TN、TP浓度进行分析，以研究输水对库区营养物浓度沿程分布的影响.从2012年输水期(4、6、8月)和非输水期(5、7、9、10月)TN、TP浓度模拟结果的平均值(图 13)可知：1)无论输水期还是非输水期，于桥水库TN、TP浓度均沿程(自果河桥至坝前)递减，并且在水库进口段存在较大的浓度梯度，这是由于水域面积突然扩大，流速放缓，TN、TP浓度被迅速稀释所致；2)非输水期，果河段水动力强度不足TN汇集浓度高；输水期，引滦输水TN浓度低使得果河段TN浓度迅速下降，进口段TN浓度梯度变小，TN被输移至库区上游导致库区TN浓度梯度增大；非输水期，库区上游TN继续向下游扩散，扩散速度大于衰减速度，使得下游TN浓度升高，整个库区TN浓度梯度减小；3)非输水期果河段TP浓度低；输水期，引滦输水TP浓度高使得果河段TP浓度急剧升高，进口段TP浓度梯度增大，TP被输移至库区上游导致库区TP浓度梯度增大；非输水期，同TN一样，库区上游TP继续向下游扩散使得下游TP浓度升高，整个库区TP浓度梯度减小；4)与水温分布类似，库区南岸比北岸更易遭受果河来流影响，输水期南岸库区上游TN和TP平均浓度分别比北岸高0.72和0.01 mg/L，从而南岸营养物浓度梯度较北岸更平缓.

从沿南岸距果河桥15、20和25 km处TN、TP浓度模拟结果随果河流量的变化过程，以研究库区水质变化与入库流量的响应关系(图 14)可知：1)非输水期，库区TN、TP均衰减而浓度降低，输水期距果河桥15和20 km处TN、TP浓度升高而出现峰值，且峰值出现时间明显滞后于果河流量变化过程；受库区水体稀释及降解作用影响，距果河桥越远TN和TP浓度峰值越小，滞后时间也越长；2)按照设计流量输水时(Q_果=55~75 m³/s)，距果河桥15和20 km处的TN浓度峰值分别为4.5和3.2 mg/L，TP浓度峰值分别为0.08和0.06 mg/L，营养物浓度峰值滞后时间分别为8和15天；25 km处未出现明显的浓度峰值；3)汛期引滦输水流量控制在30~40 m³/s，果河遭遇洪水时(Q_果=100~321 m³/s)，距果河桥15、20和25 km处TN、TP浓度峰值均明显增大，TN浓度峰值分别升高为5.9、4.2和3.4 mg/L，TP浓度峰值分别升高为0.13、0.09和0.07 mg/L；且峰值滞后时间也变短，分别为7、13和20天.分析北岸TN、TP浓度与果河流量响应关系，发现具有同样的规律，这里不再赘述.

3 结论

于桥水库水质(水温、TN和TP)变化受进口边界-果河来流影响显著，在引滦输水期与非输水期，果河入库及于桥库区水质随果河流量变化呈现出不同的时空分布特征及响应关系，具体如下：

1) 11月—次年4月，果河桥、库中心水温主要受当地气温控制先降低后升高；5—10月，非输水期，果河桥入库水温略低于库中心水温；输水期，果河流量增大使得果河桥水温明显低于库中心水温，且入库温差随流量增大呈线性升高趋势.

2) 相对于沙河，黎河引滦输水TN浓度低、TP浓度高；非输水期，果河入库水质受沙河来流控制TN浓度高、TP浓度低；输水期，引滦输水导致果河入库TN浓度随流量增大呈幂函数降低趋势，而TP浓度呈线性升高趋势.

3) 对于桥水库水温分布的数值模拟表明，夏季(6—8月)非输水期库区上游水温高，果河段、进口段及库区下游水温略低；输水期，流量增大能明显降低果河段、进口段及库区上游水温，且水温谷值滞后于果河流量变化；相对于北岸，果河来流对库区南岸的降温效果更明显.

4) 对于桥水库TN、TP分布的数值模拟表明，非输水期，TN、TP在果河段汇集浓度高，进口段浓度梯度大，库区浓度梯度小且TN、TP浓度呈衰减势态；输水期，库区TN、TP浓度升高出现峰值且滞后于果河流量变化，营养物被输移至库区上游导致库区浓度梯度增大；相对于北岸，库区南岸更易受果河来流影响从而TN、TP浓度更高.

5) 无论输水期还是非输水期，果河段TN、TP浓度均远高于库区，并且汛期输水果河流量增大使得TN、TP浓度峰值更高、滞后时间更短、库区污染更严重.为改善于桥水库水质、降低水华暴发风险，在果河段采取一定的净化措施是十分必要的，并且应减少汛期输水的情况，以削减洪峰及库区TN、TP浓度峰值.