基于数值方法的大型水库水温结构影响因素分析

邹 琳，刘 伟，蒋定国

（1.三峡大学水利与环境学院，湖北宜昌 443002；2.中国长江三峡集团有限公司，北京 100038）

0 引 言

水电作为我国的重要能源资源，对国家经济建设和社会发展具有重大作用。然而水电在实现发电效益、促进节能减排的同时，对河流生态环境也产生了环境效应［1，2］，如夏季下泄水温明显低于下游河道水温，从而对农作物产量［3］、鱼类的生长繁殖产生不利影响［4-6］。随着水库的建成蓄水，库区内水体流速减缓，与未建水库前的天然河道相比，水温过程发生明显改变，而作为水环境评价的重要因子之一，水温问题目前在水利工程生态效应和水库运行管理中是一个重要课题，也是备受关注的热点问题［7-9］。

目前，国内外针对水库水温研究方法大都采用现场监测和数值模拟［10-13］。学者们基于国内外多个水库的实测资料，提出了许多经验性水温估算方法如东勘院张大发法［14］、朱伯芳法［15］以及垂向水温分布公式［16］等。随着计算机技术的快速发展，数值模拟成为研究水库水温结构及演变规律的重要工具。数学模型的研究从一维［17］开始，但垂向一维数学模型适用范围有一定限制，一般用于回水长度短且深的水库，在入流影响较大的水库，垂向一维数学模型并不适用。为提高水库水温模拟精度，立面二维数学模型随之产生，在国内外得到了广泛应用［18，19］。由于水流运动具有三维特性，所以水温问题实际上是一个三维问题，因此，近年来国内外研究学者将研究方向转向三维数学模型［20-22］。任华堂［23］等建立三维水温数值模型对阿海水库进行模拟，研究发现取水口位置对下泄水温、均温层位置产生影响。张士杰［24］等利用MIKE3 建立二滩水库的三维水温模型，研究了上游来水水温对水库水温结构的影响，结果表明，随着上游来水温度升高，水库水温相应升高。黄廷林［25］等建立黑河金盆水库的三维水温模型，研究发现气温与金盆水库水体热分层的稳定性有显著的正相相关关系。

河道型水库的水动力条件、热力学性质复杂，自然条件下，水温结构的影响因素众多、共同作用、相互耦合，仅根据所有影响因素共同作用后的温度值进行分析，难以探明各个因素的作用程度和范围。故以西南地区某水库实测资料概化成的模型水库作为研究对象，采用数值模拟方法对重要影响因素的单独作用情况进行模拟，分析其对水库水温结构的影响程度和范围，为相关研究和工程设计提供资料参考。

1 计算方法

1.1 Delft 3D-Flow模块计算原理

水库主要在垂向上出现水温分层，同时受发电、灌溉以及泄洪的影响，水库水流具有明显的三维特征，因此本文采用Delft 3D-Flow 模块进行三维模拟。Delft 3D 是由荷兰代尔夫特研究所开发，能模拟二维和三维的水体流动、波浪运动、水质演变等，其中水动力（Flow）模块是模型的核心。为了减少由于边界不规则引入的离散化误差，水平方向采用正交曲线坐标系，垂向上采用σ坐标系。模型在浅水和Boussinesq 假定下求解不可压缩流体的纳维-斯托克斯方程，紊流模型使用k-ε模型，热输运采用三维的对流-扩散方程，方程中利用源项和汇项来模拟热量的进出，在水气界面考虑了风、辐射、热对流等因素作用，使得模型可以模拟复杂的热输运过程［26，27］。模型得到了董菁水库、香溪河流域、大宁河流域等实测资料验证［28，29］。

1.2 计算网格

根据西南地区某水库实测资料概化成模型水库，水库长40 km，库容2.000 8亿m3，模型从上游入口至800 m处为平底，入口河底高程98 m，从800 m 处开始变坡，底部高程逐渐降低，至坝前处为斜坡段，坝前河底高程0 m，河道宽100 m。水库划分网格时采用正交曲线网格，网格的余弦值小于0.02，满足计算要求，生成的网格数为75 000 个，沿河道方向网格数为1 000 个，横向网格数为5个，垂向分15层。

图1 模型水库网格图Fig.1 Model reservoir grid diagram

1.3 模型参数、率定及验证

采用西南某水库2012年实测水温数据进行对比验证，并对模型参数进行率定，率定后的模拟结果见图2，可以看出，计算得到的坝前水库水温值与实测值变化趋势基本相同，绝对误差最大值0.32 ℃，模拟结果较好，模拟精度可满足要求。此时曼宁系数为0.037，水平涡流黏滞系数为1.0 m2/s，斯坦顿数为0.001 3，重力加速度为9.81 m/s2，水的密度为1 000 kg/m3。

图2 坝前水温验证结果Fig.2 Water temperature verification results in front of dam

2 计算方案与边界条件

水库水温的影响因素有水文水力因素、气象因素以及水库几何特征等，现有研究表明入流条件、气象条件以及取水口位置对水温分布影响显著［24］，因此，重点分析入流边界和水气界面交换对水库水温结构的影响，同时考虑取水口位置的影响。模型上游入流边界采用流量边界，下游设置取水口，为简化计算，入流边界主要考虑入流流量和入流水温，水气界面考虑气温的周期性变化，其余因素如云量、相对湿度、辐射量、风速等参考西南地区实际情况进行输入。

入流流量是水库水温结构的重要影响因素，为分析不同入流流量对水库水温结构的影响程度和范围，考虑取水口位置的不同，设置不同的入流流量进行模拟分析，具体方案见表1。水库初始水温设为15 ℃，初始水位设为100 m，采用流量入流边界；出库流量设定值与入库流量相同，入流水温设为20 ℃，取水口高7 m，宽20 m。

表1 入流流量方案Tab.1 Incoming flow scheme

入流水温的周期性变化也会对水库水温结构产生影响，为分析不同周期的入流水温对水库水温结构的影响程度和范围，设置不同周期的入流水温进行模拟分析，具体方案见表2。水库初始水温设为15 ℃，初始水位设为100 m，采用流量入流边界，出库流量设定值与入库流量相同，为100 m3/s，来流水温简化为正弦变化，均值15 ℃，振幅5 ℃，即入流水温在10～20 ℃之间进行周期性变化，周期分别为12、24、168 h，取水口高7 m，宽20 m，取水口底板高程75 m。

表2 周期性入流水温方案Tab.2 Periodic inflow water temperature scheme

不同时间尺度的水气界面边界条件对水库水温结构的影响效果不同，主要考虑气温，其他因素像云量、相对湿度、辐射量、风速等参考西南地区实际情况进行输入，为分析不同周期的气温对水库水温结构的影响程度和范围，设置不同周期的气温进行模拟分析，具体方案见表3。水库初始水温设为15 ℃，初始水位设为100 m，采用流量入流边界，出库流量设定值与入库流量相同，为100 m3/s，入流水温为15 ℃，气温设为正弦变化，均值15 ℃，振幅5 ℃，即气温在10 ℃至20 ℃之间进行周期性变化，周期分别为7、30、90 d，取水口高7 m，宽20 m，取水口底板高程75 m。

表3 周期性气温方案Tab.3 Periodic temperature scheme

3 模拟结果与分析

3.1 入流流量对水库水温结构影响分析

水库初始水温设为15 ℃，来流水温设为20 ℃，可以通过水库水温的变化过程来分析入流流量对水库水温结构的影响。图3 表示在入流流量为100 m3/s 的情况下，取水口底板高程为85 m、75 m 时的水库水温历时图。由图可知，在同一来流流量下，随着水深的增加，水库水温达到与入流水温同温所需要的时间增长，同时对比图3（a）和（b）还可以看出，以取水口底板高程为分界线，在分界线以上，水库水温升温速度较快，而在分界线以下，水体的升温速度明显减缓。这是受取水口位置影响，取水口底板高程以上水体，流速大，掺混强，因而可在一定时间内使这部分水体水温趋于均匀，而取水口底板高程以下水体，流速小，掺混弱，升温极其缓慢。

图3 入流流量为100 m3/s时水库水温历时图Fig.3 Diagram of water temperature duration when inflow flow is 100 m3/s

图4 比较了当入流流量分别为50、100、150、200 和300 m3/s时，在水库坝前断面处取水口底板高程位置，水库水体的升温过程。由图4 可知，在水库同一位置处，随着入流流量的增大，水温升温所需要的时间减少。这是因为随着入口断面处入流流量的增大，水体流速增大，置换作用增强，因而可以加速水体混合时间。

图4 不同入流流量下取水口底板高程处水温升温过程对比图Fig.4 Comparison diagram of water temperature warming process at intake floor elevation under different inflow flows

为分析入流流量、水库水温变化时间和水库库容间的关系，拟采用公式（1）对其进行分析。

式中：q为入流流量，m3/s；t为入流流量作用时间，s；V为入流流量影响的库容，m3。

由于入流流量对取水口底板高程以下水体影响不大，仅考虑取水口底板高程以上水体，V=V0-V1，其中V0为水库全部库容，V1为取水口底板高程以下库容。

表4 为根据公式（1）计算所得的值，其中t0为水库受入流影响水温开始发生变化的时间，以水温变化率大于10-4（℃/s）的点为准，此时计算所得的为F0，t1为水库受入流影响水温变化至稳定状态所需时间，以水温变化率小于10-4（℃/s）的点为准，此时计算所得的为F1，可见入流流量、水库水温变化时间和水库库容存在对应关系。在遇到F小于F0的洪量，对水库水温结构一般没有影响，在遇到F大于F0且小于F1的洪量，水库水温结构发生变化，在遇到F等于F1的洪量，水库水温结构变化至另一稳定状态，此后的来水不会再对水温结构产生影响。根据F值，在已知水库取水口位置和入流流量的情况下，可以大致预估水库水温结构。

表4 入流流量、水库水温变化时间和库容关系表Tab.4 Table of relation between inflow discharge，water temperature change time and storage capacity

传统的判别水库水温结构的β法定义为一次洪水总量比总库容，分层型水库在遇到β大于1 的洪量，会变成临时混合型，在遇到β小于0.5 洪量，对水库水温结构一般没有影响，在遇到β在0.5 至1 之间的洪量，对水库水温结构的影响一般介于两者之间。表5 以入流流量100 m3/s 为例，在水库取水口底板高程为85、75 m 的情况下，β值和F值的对比表。由表5 可知，在取水口底板高程为75 m 时，当t=300 h，洪量为1.08 亿m3时，β值为0.54，F值为1.50，β法和F值法均表明此次洪量已经对坝前水库水温结构产生影响；当t=496 h，洪量达到1.785 6 亿m3时，β值为0.89，表明坝前水温结构已经接近混合，而F值为2.48，表明来水对水库水温结构影响已经结束，此后的来水不会再对水温结构产生影响；当t=700 h，洪量达到2.52 亿m3时，β值为1.26，表明坝前水温结构已经混合。由图5 可知，F值判断更为准确。

图5 不同洪量下水库水温分布图Fig.5 Distribution diagram of reservoir water temperature at different flood volumes

表5 β值与F值对比表Tab.5 Table of β and F values

3.2 周期性入流水温对水库水温结构影响分析

图6 表示入流水温周期分为12、24、168 h（7 d）时，库尾水温历时对比图。取水口底板高程75 m，由图6可知，当入流水温周期为12 h，选取的库尾断面处总水深约为25 m，表中层受来流水温周期变化基本一致，底层振幅略有衰减且升温过程发生比较明显变形；当入流水温周期为24 h，选取的库尾断面处总水深约为26 m，表中层受来流水温周期变化基本一致，底层升温过程发生比较明显变形；当入流水温周期为7 d，选取的库尾断面处总水深约为30 m，表中底层受入流水温周期变化基本一致。随着入流水温周期的增大，影响的库尾水深范围变大，使得小于等于此水深的水库库尾范围内整体变化趋势与入流水温变化趋势基本一致，但仍可以看出入流水温对库尾的最大影响水深仍在取水口底板高程分界线附近。

图6 不同入流水温周期下库尾水温历时对比图Fig.6 Comparison diagram of reservoir tail water temperature duration under different periodic inflow water temperature

表6列出了不同入流水温周期在坝前的水温波动范围。由表可知，不论入流周期如何变化，随着水深的增加，水温逐渐降低，这是因为当入流水温以周期性进行变化时，高温水从表层流入水库，低温水逐渐下沉。水温波动范围的差值先增大后减小，大约在取水口底板高程分界线处差值达到最大。当处于同一水深，随着入流周期的增大，水温波动范围也随之增大。

表6 不同入流水温周期下水温波动范围表Tab.6 Water temperature fluctuation range under different periodic inflow water temperature

3.3 周期性气温边界条件对水库水温结构影响分析

由于水面是和空气相接触的，水气界面需要进行热交换，气温在一定程度上会对水库水温产生影响。周期性气温边界条件对水库水温结构影响分析考虑了气温的周期性变化，其余因素如云量、相对湿度、辐射量、风速等参考西南地区实际情况进行输入。图7～9表示气温周期为7、30、90 d时，选取坝前断面的水库水温随时间变化图。由图可知，以水温变化大于0.2 ℃的区域设为气温影响区域，当气温变幅一致时，随着气温周期的增加，对水库水温的影响范围逐渐增大，当气温周期为7 d时，经过一个周期后，影响水深约为20 m，表现为影响区域以内水温升高，当气温周期为30 d 时，经过一个周期后，影响水深约为25 m，表现为影响区域以内水温升高，当气温周期为90 d时，经过一个周期后，影响整个影响水库水深，表现为随着水深的增加，水温逐渐降低。同时，随着气温周期的增加，表层水温的波动范围增加，当气温周期为7 d 时，表层水温波动范围为15.0～15.4 ℃，当气温周期为30 d 时，表层水温波动范围为15.0～16.0 ℃，当气温周期为90 d 时，表层水温波动范围为15.0～17.1 ℃。

图7 气温周期为7 d时坝前水库水温随时间变化图Fig.7 Diagram of water temperature variation with time in front of dam when temperature periodicity is 7 d

图8 气温周期为30 d时坝前水库水温随时间变化图Fig.8 Diagram of water temperature variation with time in front of dam when temperature periodicity is 30 d

图9 气温周期为90 d时坝前水库水温随时间变化图Fig.9 Diagram of water temperature variation with time in front of dam when temperature periodicity is 90 d

4 结 论

针对概化的模型水库，采用数值方法分析入流流量、入流水温和气温对大型水库水温结构的影响程度和范围。研究表明，受取水口位置影响，取水口底板高程以上水体，流速大，掺混强，而取水口底板高程以下水体，流速小，掺混弱，因此对于大型深水库，入流流量、入流水温对水库水温结构的影响区域与取水口高程有一定相关性；随着入流流量的增大，水库水温结构变化至另一稳定状态所需要的时间减少；随着入流水温周期的增大，影响的水深范围增大，但最大影响水深仍在以取水口底板高程为分界线附近；随着气温周期的增大，对水库水温的影响范围逐渐增大，经过一个周期后，当气温周期为7 d 时，影响水深约为20 m，当气温周期为30 d时，影响水深约为25 m，当气温周期为90 d 时，可影响整个影响水库水深；同时改进后的F值法，仅考虑取水口底板高程以上水体，相比较于β法，可更为准确地判断坝前水温结构的变化。