水库调洪演算程序

版本 1.1
输入库容曲线 · 泄流曲线 · 多个洪水过程线 · 起调水位,可同时算设计洪水、校核洪水等多个工况: 逐工况给出最高库水位、最大泄流量、削峰流量,并生成多工况成果对比表,一键导出 Excel
梯形平均式隐式试算· 直线 / 四点三次插值· 多工况并列计算· 2026年10月23日· Copyright © longlong

1基本参数起调水位须落在库容曲线与泄流曲线的节点范围内

各工况共用同一调洪起始水位
默认梯形平均式,三种格式解同一水量平衡方程
库容、泄流、洪水过程线三条曲线同时生效

2基础资料输入直接从 Excel 复制两列粘贴即可,Tab / 空格 / 逗号均识别

库容曲线0 节点
水位 Z (m)库容 V (万 m³)
水位递增;可加第 3 列水面面积(万 m²)
泄流曲线0 节点
水位 Z (m)泄流量 q (m³/s)
即 q = q(Z),自由泄流关系
洪水过程线0 节点
工况名称
时间 t (h)入库流量 Q (m³/s)
时间递增,首末宜为 0;每个工况一条洪水过程线,末节点即该工况洪水历时
计算方法与公式说明

一、基本方程。水库水量平衡的微分形式为

ƒ(z)·dz/dt = Q(t) − q(z)

式中 ƒ = ƒ(z) 为库水面面积;z = z(t) 为库水位;Q = Q(t) 为入库流量;q = q(z) 为出库流量。 移项并定义调洪函数

F(t, z) = [Q(t) − q(z)] / ƒ(z) , dz/dt = F(t, z),z(t₀) = z₀

程序按所选方法对该初值问题逐步积分;ƒ(z) 由库容曲线对水位求导得到(也可在库容表第 3 列直接给水面面积)。

二、梯形平均式(隐式试算,默认)。工程上最常用的差分格式,对相邻两时段取平均:

ΔV = [ (Q₁+Q₂)/2 − (q₁+q₂)/2 ] · Δt

由于 q₂ 依赖 V₂(即依赖未知的 Z₂),方程为隐式。程序逐时段迭代:假定 V₂ → 由库容曲线反查 Z₂(二分法) → 插值求 q₂ → 回代上式更新 V₂,直到相邻两次迭代的 V₂ 之差小于 1×10⁻¹¹ 万 m³,收敛后再推进下一时段。

三、多工况计算。每个洪水过程线为一个独立工况(如 100 年一遇设计洪水、2000 年一遇校核洪水), 共用同一库容曲线、泄流曲线与起调水位,各自以自身洪水历时为计算终止时间并行计算, 成果汇总成对比表,可直接得到设计洪水位、校核洪水位及其相应下泄流量。

四、插值口径。下拉框选择后对库容曲线、泄流曲线、洪水过程线同时生效:

· 直线插值(默认):节点间线性连接,稳定、无过冲,与手工查图习惯一致,洪水总量守恒;
· 四点三次多项式:取目标点相邻 4 个节点做拉格朗日三点式(三次)插值,曲线更平滑, 首末区间固定取前 4 点或末 4 点。节点少于 4 个时自动退化为直线插值。
  注:节点稀疏时三次插值会在节点间产生轻微过冲,洪水总量与节点直线积分略有差异,比对已有成果时建议先用直线。

五、计算步长。Δt 不需人工输入,由程序按各工况洪水过程线的节点最小间隔自动确定 (取其 1/2 后就近取常用值 0.005 / 0.01 / 0.02 / 0.025 / 0.05 / 0.1 / 0.2 / 0.25 / 0.5 / 1 / 2 / 5 h); 末时段按实际剩余时长收尾,使成果表最后一行的绝对时间正好等于该工况的洪水历时。

六、单位与量纲。Q、q 为 m³/s,Δt 为 h,V 为万 m³,则

ΔV (万m³) = Q (m³/s) × Δt (h) × 3600 / 10000 = 0.36·Q·Δt

七、成果校核。表中"水量平衡差"列为该时段 ΔV 与方程右端之差,梯形隐式格式下应为 0(仅存舍入误差); 若用龙格—库塔或改进欧拉,该列给出该时段积分格式与梯形平均式的偏差,可用于判断时段步长是否足够小。 全过程校核:入库总量 − 出库总量 = 末库容 − 起调库容。

八、导出 Excel。工作表为:库容曲线、泄流曲线、洪水过程线(各工况并列)、 调洪成果(单工况时即该工况成果表;多工况时为汇总对比 + 每个工况一张明细表), 明细表列式为:时间、入流量、出流量、库水位、库容、库容变差(原《水库洪水调节计算成果表》中的 "下游控制断面流量"列在自由泄流条件下与出流量相同、恒为空白,已略去)。

说明:本程序按自由泄流(q = q(Z))计算,不含闸门控泄、机组过流与下游防洪补偿等约束; 若需按入库流量分级控泄,可另在 q(Z) 曲线上按控制条件分段给定。