搞水文预报的人基本都有过这种体验:手里拿到气象模式输出的降水预报,图看着挺漂亮,但心里还是没底——不知道这些雨落到具体的流域里,到底会形成多大的洪峰、什么时候到。反过来,做水文模型的人又常年被降水输入坑,站点太稀疏、雷达估测偏差大,模型参数调得再好,雨没给对,流量照样对不上。
WRF-Hydro 就是冲着这个痛点来的一个开源社区模型,核心价值在于「气象水文耦合」:它把气象模式输出的降水、辐射、风、温湿等大气强迫,和坡面产汇流、土壤水过程、河道汇流、地下水过程放进同一套框架里模拟,让“天上下雨”和“河道涨水”这两个原本分属不同学科的问题,在一条代码流水线里闭环解决。这篇文章我从建模实操的角度,把这套模型的原理、编译、数据准备、参数配置、案例评估和踩坑经验完整串一遍,适合刚接触 WRF-Hydro、打算在自己流域里跑通第一个算例的人参考。
1. 先把两个系统之间的“翻译层”讲清楚
1.1 单跑水文模型时你在猜什么
传统水文模型(SWAT、HEC-HMS、新安江这些都算)的基本输入是降水,输出是流量过程。但这里有个长期没解决好的问题:降水到底从哪里来?如果是历史事件模拟,你可以用站点插值或雷达栅格;如果是预报场景,你只能拿气象模式的输出。可是气象模式本身也是模拟结果,空间分辨率通常在几公里到几十公里,落到一个几平方公里的小流域上,一场暴雨可能就劈成两半,一半在流域内、一半在流域外。
更麻烦的是热力过程。一场雨落下来,有多少入渗、多少蒸发、多少变成壤中流,取决于土壤湿度、植被状态、辐射和风速,这些恰恰是传统水文模型最弱的部分——很多模型干脆把蒸散发简化成一个经验系数。结果是:你往往不是在模拟,而是在猜。这也是为什么同样一场雨,放在不同的前期土壤条件下,实际洪水差异能有一倍甚至更大。
1.2 WRF-Hydro把哪些过程接进了闭环
WRF-Hydro 的基本思路是把陆面过程模式(Noah-MP)和水文汇流模块串起来。Noah-MP 负责算每个网格上的能量平衡和水量平衡:降水落到地表后,先过植被截留,再算入渗、土壤水再分配、蒸发蒸腾,最后算出每个网格的产流量。这些产流量并不直接进河道,而是先进入坡面汇流层,在网格之间的地形梯度驱动下,通过坡面汇流和饱和地表径流,慢慢聚集到河道格点,再由河道汇流模块往下游推进。
模型里明确区分了几个空间过程和尺度:
- 陆面网格(Land Surface Grid):通常 1km 左右,跑 Noah-MP 的地表过程。
- 汇流网格(Routing Grid):比陆面网格更细(100m、250m 常见),在这个尺度上算坡面汇流、地下水和河道之间的转换。
- 河网(Reach/Channel):通过 Route_Link 文件描述每一段河道的几何和水力特征。
另外还有一个地下水蓄水池(bucket)简化模块,用来模拟基流过程。模式可以离线跑(用再分析或观测强迫数据驱动),也可以在线耦合进 WRF 气象模式,实现大气-陆面-水文全耦合预报。实际业务里离线模式用得多,因为计算量可控、调试方便。
2. 编译安装最容易被依赖库耗死
2.1 版本和编译器怎么搭配
WRF-Hydro 当前主流是 5.x 系列,官网社区版代码从 GitHub 上拿,更新很快。我的建议是不要一味追新,5.3.0、5.5.x 这类发布过一段时间、社区里踩坑文档多的版本更省心。
编译这块核心矛盾集中在编译器套件和 NetCDF 库的匹配。实践下来,Intel ifort 编译器和英特尔的 NetCDF/HDF5 库兼容性最好,踩坑最少;如果你只有 gfortran,也能编,但要注意 NetCDF-Fortran 必须和 gcc/gfortran 是同一套配置编译出来的,否则 Fortran 接口经常对不上。
依赖库的编译顺序基本是固定的:zlib → HDF5 → NetCDF-C → NetCDF-Fortran。MPI 库(MPICH 或 OpenMPI)可以放在最前面装好,后面所有库都用它来编,保证运行的时候 MPI 环境一致。
2.2 依赖库编译顺序与常见报错
我整理了一个相对稳妥的编译路径:
# 1. 设置环境变量 export CC=mpicc export CXX=mpicxx export FC=mpifort export F77=mpifort export NETCDF=/path/to/netcdf export HDF5=/path/to/hdf5 # 2. 编译 HDF5(并行版) ./configure --prefix=$HDF5 --enable-parallel --enable-fortran make -j8 && make install # 3. 编译 NetCDF-C ./configure --prefix=$NETCDF --disable-dap make -j8 && make install # 4. 编译 NetCDF-Fortran export LD_LIBRARY_PATH=$HDF5/lib:$NETCDF/lib:$LD_LIBRARY_PATH export CPPFLAGS="-I$HDF5/include -I$NETCDF/include" export LDFLAGS="-L$HDF5/lib -L$NETCDF/lib" ./configure --prefix=$NETCDF make -j8 && make install # 5. 编译 WRF-Hydro cd WRF_HYDRO # 5.x 用 CMake cmake -DNETCDF_DIR=$NETCDF -DHDF5_DIR=$HDF5 -DCMAKE_INSTALL_PREFIX=$HOME/wrfhydro_run .. make -j8新手最容易在 NetCDF-Fortran 这里卡住。它 configure 的时候要找到 NetCDF-C,经常因为环境变量没配对,报 “C compiler cannot create executables” 这类误导性错误。遇到这种问题别慌,先检查 LD_LIBRARY_PATH,再确认 NetCDF 版本位数是否一致(32 位 / 64 位混用必挂)。
编译完成后,确认一下wrf_hydro.exe已经生成,然后跑一次性感测试:用官方自带的 small test case 跑通一遍,看输出目录里有没有生成水文输出文件。这个步骤能帮你把“模型本身的问题”和“你自己流域数据的问题”分开。
2.3 功能裁剪:serial 跑小流域其实够用
很多人一上来就想着上 MPI 并行,其实没必要。对几百平方公里的小流域,汇流网格撑死几万个,单核跑完全能接受。WRF-Hydro 编译时可以选择不带 MPI 的 serial 模式,好处是少一个变量,排查问题简单很多。
我个人的流程是:第一次跑通,用 serial;确认物理过程和参数没问题以后,再编译 MPI 版本做正式实验。这样即使后面报了内存或并行通信相关的错,也知道不是数据或配置的问题。
3. 让模型认识流域:空间离散与GIS数据预处理
3.1 陆面网格和汇流网格的关系
WRF-Hydro 一个让新人困惑的点是:为什么同一套模拟里有两套网格,而且网格大小还不一样。原因很简单,陆面过程需要相对大的网格来稳定地刻画能量平衡,而产汇流对地形细节敏感,河道位置偏一个格点,模拟的洪峰可能差出半小时。
这两套网格之间通过“聚合权重”建立映射:每个汇流网格根据高程、坡度被分配到它对应的陆面网格,产流量从陆面网格传递到精细汇流网格时,按权重分配。这个关系写在一个叫 Fulldom 的文件里,全称通常是Fulldom_hires.nc,这个文件是模型能不能正确跑起来的关键。
Fulldom 里记录了几样核心信息:流向(flow direction)、坡度、每个格点是不是河道格点,以及和陆面网格的对应关系。流向计算错了,水就往山上走,模拟出来的流量过程基本不能用。
3.2 Fulldom和Route_Link是怎么生成的
生成这两个文件的标准路径是用 NCAR 提供的 GIS 预处理工具(现在常见的是基于 QGIS 的插件集,也有 Python 脚本版本)。大致流程:
- 准备 DEM 数据,分辨率和你要生成的汇流网格一致(比如 100m)。
- 对 DEM 做填洼处理,避免洼地格点困住水流。
- 计算流向(D8 算法)和累积流量。
- 按累积流量阈值提取河网,并对河网做分级(Strahler 分级)。
- 划分流域边界,生成每段河道(reach)的起止点。
- 运行预处理工具,输出
Fulldom_hires.nc和Route_Link.nc。
这个流程里最容易被忽略的是投影坐标系统。DEM 的投影和模型运行时用的坐标必须一致,否则网格对齐直接错位,运行报错还好,怕的是不报错但结果很怪。建议全程统一用投影坐标(比如 UTM),只在最后生成文件时保留 WGS84 的经纬度信息供模型读取。
Route_Link 文件描述河道的拓扑和水力属性:每一段河道的长度、坡度、糙率(Manning's n)、上下游连接关系。这个文件直接驱动河道汇流模拟,它的河段顺序和检查点的关系如果搞错,模型会把水送到错误的河段。
3.3 气象强迫数据准备中最容易忽略的事
离线跑 WRF-Hydro,你需要准备一套大气强迫数据,核心变量包括:
- 降水(precipitation)
- 近地面气温
- 比湿
- 气压
- 风速
- 向下短波辐射
- 向下长波辐射
时间分辨率至少小时级。ERA5 和 NLDAS 是常用的公开数据源。这里我要特别提三个坑:
第一,单位。降水的单位经常是 kg/m²/s,换算成 mm/h 要乘 3600,很多人在这一步差了三到四个数量级还没发现,模拟出来流量大得离谱。
第二,坐标。ERA5 的格点和 WRF-Hydro 的陆面网格一般不重合,需要重采样。重采样时要注意插值方法,降水变量不要用平滑的样条插值,它会把你好不容易保留的降水峰值给磨平掉;用最近邻或者距离反比加权更合适。
第三,时间戳。模型对强迫数据的时间连续性很敏感,文件命名里的时间戳和 namelist 里的开始结束时间必须对得上,否则模型在某个时刻找不到数据,直接崩溃。
4. 两个namelist文件决定模型行为
4.1 namelist.hrldas的离线开关与输出
WRF-Hydro 的配置核心是两个文本文件。第一个是namelist.hrldas,控制陆面部分(Noah-MP)和强迫数据读取。离线模式下关键配置大致长这样:
&NOAHLSM_OFFLINE START_TIME = "2016-06-01_00:00:00", END_TIME = "2016-07-01_00:00:00", NOAHLSM_OFFLINE = .true., FORCING_TYPE = 1, DX = 1000, DT = 300, RESTART = .false. / &OUTPUT OUTPUT_TIMESTEP = "3600", RESTART_FREQUENCY = "1440", OUTPUT_TYPE = "netcdf" /这里FORCING_TYPE = 1表示从 NetCDF 格式格点强迫文件读取;DX是陆面网格大小,单位通常为米;DT是陆面步长。我刚接触时有个误区,觉得步长越好越低精度越高,实际跑下来,Noah-MP 在 1km 网格上用 300 秒到 600 秒的步长已经足够稳定,太小只会徒增计算量。
4.2 hydro.namelist里的每个关键参数
第二个文件是hydro.namelist,控制产汇流过程。它决定模型的行为走向,也是率定阶段主要动刀的地方。核心参数包括:
&HYDRO HYDRO_REAL = 1, DTRT_OPT = 3, DTRT_SUBOPT = 1, DTRT_CH_OPT = 3, MANNING_OV = 0.15, MANNING_CH = 0.035, GWBASEMODEL = 1, ROUTE_TIMESTEP = "60", OUTPUT_TIMESTEP = "3600", CHANRTSWCRT = 1, TERRAIN_ROUTING = 1, /逐个解释:
HYDRO_REAL:是否启用水文路由模块,离线跑必须设 1。DTRT_OPT和DTRT_SUBOPT:坡面产流和地下水采用哪种数值方案,我建议先用文件附带的默认方案,不要乱改。MANNING_OV:坡面漫流曼宁糙率,影响洪峰形态和汇流速度。值越大,流速越慢,洪峰越平缓。MANNING_CH:河道曼宁糙率,直接决定河道里的流速。默认 0.035 可作为起点,实际要根据河床质和断面情况在 0.02~0.08 之间试。GWBASEMODEL:地下水简化模型开关,开启后能模拟基流过程,对枯季流量和退水段影响非常大。ROUTE_TIMESTEP:汇流步长,这个参数直接关系到数值稳定性,细网格下太大很容易算发散。
另一个经常被忽略的是GWBUCKPARM.nc,它是地下水蓄水池参数的栅格文件,决定基流衰减系数。很多新手调了半天退水段还是对不上,问题往往不在曼宁糙率,而在 GWBUCKPARM 里的slope和exponent参数。
4.3 时间步长搭配和运行稳定性
时间步长组合是运行稳定性的关键。我通常会采用“陆面步长 + 汇流步长”的搭配来做测试。汇流步长受数值稳定性条件约束,太大会出现流量锯齿状振荡甚至溢出。判断标准很简单:看输出流量过程线有没有不自然的上下抖动。
一个相对稳的起点:1km 陆面网格配 600 秒陆面步长,汇流网格 250m 配 60 秒汇流步长。如果你把汇流网格细化到 100m,又是山区陡坡地形,60 秒可能就不稳了,试试压到 30 秒或者 20 秒。跑一次大流域这种步长差距带来的计算量差异不是线性增长,而是指数级的,所以不要为了追求精细而无脑加密。
5. 一个丘陵小流域洪水模拟的完整案例
5.1 流域情况和数据来源
我用一个典型的南方丘陵小流域作为例子说明:集水面积约 520 平方公里,主河道长 35 公里,主河槽平均比降约 0.003,植被以次生林和耕地为主,土壤偏黏性。DEM 原始分辨率 30m,重采样到 100m 汇流网格,陆面网格设为 1km,河网提取阈值设在 20 平方公里以上,最终生成了 40 段河道。
强迫数据用 ERA5-Land 小时数据,重采样到模型格点后用当地 4 个雨量站的逐小时降雨做偏差校正。模拟时段选 2016 年 6 月中旬一次锋面降雨过程,过程总降雨 128mm,最大小时雨强 22mm/h。
5.2 从冷启动到spin-up的运行路径
第一次跑建议不要直接从洪水过程当天开始,因为土壤初始含水量是未知的。标准做法是先冷启动跑一段预热期(spin-up),让模型自行把土壤水调整到相对合理状态,再从重启文件继续跑目标事件。
具体操作:
- 从 2016 年 4 月 1 日开始冷启动,跑两个月。
- 到 6 月 1 日输出
RESTART文件。 - 把 namelist 里
RESTART改成.true.,START_TIME改为 6 月 1 日,引用上一步的 restart 文件。 - 正式跑 6 月的洪水过程。
运行命令很简单:
mpirun -n 16 ./wrf_hydro.exe这个 520 平方公里的流域,16 核跑 30 天模拟,大约 1.5 到 2 小时。计算负担主要不在格点数量,而在汇流步长:汇流步长越小,总步数越大,计算时间成倍增加。
5.3 模拟结果怎么评估才算数
模拟完成后,拿流域出口断面的模拟流量和实测流量对比。评估指标不要只看纳什系数(NSE),还要看 KGE(Kling-Gupta Efficiency)和水量偏差 PBIAS。我通常还会画一张包括实测流量、模拟流量和降雨的过程线图,因为指标只能告诉你“好坏”,过程线才能告诉你“哪里不对”。
调参之前,我的第一版模拟 NSE 只有 0.66,KGE 0.51,PBIAS 高达 28%,说明模拟径流总量明显偏大。看过程线发现,降雨开始后模拟洪峰比实测早了将近 3 个小时,退水段又跌得太快,基流偏小。这一眼就能判断出问题出在哪里:峰现时间早,大概率是坡面糙率太小、汇流太快;基流不足,是地下水参数没调对。
调整思路如下:
- 把
MANNING_OV从 0.10 提高到 0.18,推迟坡面汇流速度。 - 修改地下水蓄水池参数,增加基流衰减时间常数。
- 对河道曼宁糙率先不动,因为河道部分的不确定性相对小,放在最后微调。
调整后 NSE 提升到 0.83,KGE 0.76,PBIAS 降到 7%,峰现时差控制在 30 分钟以内。这个结果对于丘陵区中小流域已经算不错了。
6. 率定顺序、敏感参数和踩坑实录
6.1 先水量再过程:率定参数的顺序逻辑
率定顺序有一套通用逻辑,核心原则是“先水量平衡,再流量过程,最后峰形细节”。我见过很多新手上来就调曼宁糙率,因为觉得洪峰对不上就是糙率的事,结果调了半天峰形好了,总量还是差得很远。
正确的顺序应该是:
- 检查降雨总量和径流系数,如果 PBIAS 大于 20%,先修正产流端:入渗参数、土壤饱和导水率。
- 再看基流和退水段,这时候动地下水蓄水池参数。
- 然后看峰现时间和洪峰形态,这时候调坡面糙率和河道糙率。
- 最后微调峰形,调整河道糙率的下游分布,或者初始土壤含水量的空间分布。
每一步只动一个参数,记录结果和影响方向。我自己的做法是每次调参后保存一个表:
| 目标问题 | 优先调整参数 | 调整方向 | 影响 |
|---|---|---|---|
| 径流总量偏大 | 入渗相关参数 | 增大入渗量 | PBIAS 下降 |
| 基流偏低、退水太快 | 地下水蓄水池系数 | 增大衰减时间常数 | 退水段变缓 |
| 洪峰偏早 | 坡面曼宁糙率 | 增大糙率 | 峰现时间推迟 |
| 洪峰偏高 | 河道曼宁糙率 | 增大糙率 | 洪峰削减、展平 |
6.2 高频报错对照表
把我在实际运行中遇到频率最高的几类问题整理成一张表:
| 现象 | 可能原因 | 处理方法 |
|---|---|---|
| 输出全为 NaN 或 -9999 | 初始土壤含水量为负;GWBUCKPARM 缺失或参数异常 | 检查土壤初始场,重建地下水参数文件 |
| 河道流量锯齿状振荡 | 汇流步长太大,数值不稳定 | 把 ROUTE_TIMESTEP 减半,逐次测试 |
| 运行到一半直接崩掉 | MPI 通信配置问题或内存不足 | 先用 serial 模式单核跑,排除并行问题 |
| NetCDF 打开变量失败 | geo_em 或 Fulldom 文件与模型版本不匹配 | 确认 GIS 预处理工具和模型版本配套 |
| 修正后结果完全不变 | 参数文件被缓存重建方式不对 | 确认 namelist 修改后重新完整初始化,不要沿用旧 restart 文件 |
第 5 条是隐蔽性最强的一个坑。我印象很深的一次,改了半天参数结果完全不变,最后发现是启动时自动读取了旧 restart 文件,模型根本没有重新跑陆面初始化,所有参数改动都被旧状态覆盖了。遇到这种情况,先检查 restart 文件的时间戳,再决定是否需要删掉重跑。
6.3 我的一点个人体会
WRF-Hydro 这套系统,不能说难,但确实复杂。它横跨了气象、水文、计算机并行、GIS 数据处理四门功夫,任何一个环节不熟都会卡住很久。我个人实际接触下来,最大的体会是:数据质量永远比参数率定重要。降水强迫差 20%,你靠调参数把 NSE 调上去,也只是把错误“补偿”掉了,换一场雨就露馅。参数率定只能校正模型结构误差,不能弥补输入数据错误。
另外一个很实用的建议是,新流域启动时先别急着追求高分辨率。1km 大气强迫配 250m 汇流网格,对大部分中小流域洪水模拟已经够用。你真把一个 500 平方公里的流域加密到 100m 汇流网格,计算量翻三四倍,收益却相当有限,因为上游的气象强迫数据本身也就那点分辨率,水文模型单方面加密,瓶颈还是在输入。
最后分享一个平时容易忽略却很有用的操作:把输出结果跟遥感蒸散发产品、土壤湿度产品做交叉验证。流量只有出口一个断面能对,但土壤湿度分布和蒸散发能反映陆面过程是否合理。我后来在不少项目里靠这个办法识别出了数据预处理阶段的系统误差——那些在流量过程线上看不太出来、但真实存在的问题,往往藏在这些辅助变量里。