Chapter 3 SHUD模型软件与数据

3.1 安装 SHUD and rSHUD

本章介绍SHUD模型和rSHUD的安装过程，以及必要的开发库的安装。

3.1.1 SUNDIALS/CVODE

SHUD模型的求解需要SUNDIAL中CVODE库的支持。 SUNDIALS（SUite of Nonlinear and Differential/ALgebraic equation Solvers）由六个求解器组成；其中的 CVODE 是一个用于求解刚性和非刚性常微分方程系统（初始值问题）的数值求解器，数学表达为： \(y' = f(t,y)\). CVODE支持可变阶精度、可变步长数值方法。CVODE字母组合含义为“C版本的可变系数常微分方程求解器(Variable-Coefficient Ordinary Differential Equation Solver)” 。使用SHUD模型只需要CVODE支持，因此无需安装SUNDIALS所有求解器，只需要安装CVODE即可。

因为SUNDIALS/CVODE持续更新中，代码和函数定义等等变化较大，因此SHUD模型通常仅支持CVODE较新版本，并非所有版本都可支持。

SUNDIALS/CVODE的安装方法有两种：（1）自行下载安装包； （2）使用SHUD源代码的安装脚本。

3.1.1.1 使用脚本安装CVODE

源代码文件夹中的configure.sh文件是自动安装CVODE的脚本——可以自动从github上下载CVODE最新源代码并通过cmake编译并安装在用户目录下。

3.1.1.2 自行下载并安装CVODE

SUNDIALS/CVODE 代码的下载地址为：https://computation.llnl.gov/projects/sundials/sundials-software。

CVODE 5.x安装步骤：

  1. 进入命令行界面，并解压下载好的cvode压缩包。
  2. 在CVODE的根目录下新建文件夹 builddir， 并进入该文件夹。
mkdir builddir
cd builddir/
  1. 测试 ccmake. 如果没有 cmake请安装。
ccmake
  1. 运行ccmake配置安装环境。
ccmake ../sundials/cvode-5.0.0
Screenshot of Step 1 in ccmake
Screenshot of Step 1 in ccmake

首先出现的是空白配置页面。按c开始配置。

Screen short of Step 2 in ccmake
Screen short of Step 2 in ccmake

上图为默认的配置信息图片上的四处配置信息需要注意，值需要改为下面的值。OPENMP值需要您的系统支持OpenMP，因此在确认您的系统支持OpenMP的情况下打开。

BUILD_CVODE = ON
CMAKE_INSTALL_PREFIX = ~/sundials
EXAMPLES_INSTALL_PATH = ~/sundials/examples

以上值修改完成后，按c确认配置。

Screenshot of Step 3 in ccmake
Screenshot of Step 3 in ccmake

ccmake开始自动配置安装参数。当配置完成后，按g生成配置文件并退出。

  1. 然后使用make编译并安装：
make
make install

3.1.2 编译SHUD

编译信息保存在文件Makefile，其中关键的编译参数:

  1. SUNDIALS_DIR的保存路径（非常关键）. 默认值为 ~/sundials，如果您已经安装到此路径，即可保持默认值。
  2. OpenMP的路径。当您的系统支持OpenMP时，可根据实际安装路径配置。
  3. SRC_DIR路径。默认值为 SRC_DIR = .，即当前源码所在路径。
  4. BUILT_DIR路径，默认值为BUILT_DIR = . ，即当前源码所在路径。

配置好Makefile文件即可编译，使用以下命令进行编译：

make clean
make shud

更多可选的编译参数：

  • make all - 清理旧的编译文件，并编译SHUD和SHUD_omp文件。
  • make help - 展示编译的帮助信息。
  • make shud - 编译shud模型的可执行文件。
  • make shud_omp - 编译OpenMP支持的并行版本，编译后的可执行文件名为shud_omp。

3.1.2.1 OpenMP（可选）

当需要并行计算时，需要安装OpenMP。以下为Mac OS上安装的参数：

brew install llvm clang
brew install libomp
compile flags for OpenMP: 
  -Xpreprocessor -fopenmp -lomp
Library/Include paths:
  -L/usr/local/opt/libomp/lib 
  -I/usr/local/opt/libomp/include

3.1.2.2 运行 SHUD 可执行文件

运行SHUD模型，请在命令行内执行一下命令：

./shud <projectname>

Default screen print of shud and command tips 命令的参数模式：

./shud [-0][-p projectfile] [-o output] [-n Num_Threads] <project_name>

  • -0 空计算。读取输入文件并写出结果文件，但是内部没有计算。仅仅用于验证输入文件格式的可靠性和可用性。

  • <project name> 项目名称。所有的项目文件必须以此开头。

  • [-p projectfile] 指定的项目配置文件，文件中包含了所有输入文件的地址。

  • [-o output_folder] 制定的结果文件输出地址。默认的输出地址为output/projname.out。

  • [-n Num_Threads] 并行计算的线程数量。此参数仅在并行计算打开时可用。

当shud 程序正常运行时，屏幕输出应该如下图所示： Screen print of SHUD running.

3.1.3 rSHUD

rSHUD是一个标准R包，但尚未向CRAN提交，因此需要通过github以源码方式安装。在R环境中执行以下命令即可自动安装rSHUD。

如果您的电脑上尚未在R环境中安装devtools，请先安装。

install.packages('devtools')

然后通过devtools从github上安装rSHUD。

devtools::install_github('SHUD-System/rSHUD')

rSHUD依赖的R工具包包括：

  • Rcpp
  • reshape2
  • ggplot2
  • gridExtra
  • grid
  • fields
  • xts
  • hydroGOF
  • zoo
  • raster (>= 2.1.0)
  • sp
  • rgeos
  • RTriangle
  • rgdal (>= 1.1.0)
  • proj4
  • abind
  • utils
  • lubridate
  • geometry
  • methods
  • ncdf4
  • GGally
  • doParallel

3.2 模型输入数据

本节从原始数据要求到输入文件内涵，向用户详细介绍模型输入数据。

3.2.1 原始数据

类型 数据名 必要性 要求(格式、步长、数量） 说明
空间数据 流域边界 高 矢量 连续、封闭、唯一多边形；
空间数据 河流网络 高 矢量 河流具有方向性；河流无向下分叉；
空间数据 高程 高 0 单位：米
空间数据 土地利用分类 中 矢量/栅格 \(n_{lc}\)个分类
空间数据 土壤分类 中 矢量/栅格 \(n_s\)个分类
空间数据 气象站点 中 矢量/栅格 \(n_{fc}\)个站点
参数 土地利用 中 \(n_{lc}\) 组 土地利用的水力学参数
参数 土壤 中 \(n_s\)组 土壤性质：Silt百分比、Clay百分比、有机质含量、Bulk Density
时间序列数据 气象驱动 高 小时～日 降雨、气温、湿度、辐射、风速、气压
时间序列数据 叶面积指数(LAI) 中 小时～月 每类土地利用一个LAI时间序列
时间序列数据 融化系数 低 小时～月 一个时间序列
时间序列数据 观测数据 高 小时～月 通常为径流数据、地下水、蒸散发等

3.2.1.1 数据实例

DWR
南苏丹Pongo流域原始数据实例：流域边界（绿色实线），河流（红色实线）、DEM(背景)和气象站点覆盖范围（黑色虚线）。气象数据站点为GLDAS数据，因此为$ 0.25 x 0.25 $度覆盖范围。
DWR DWR
(a) Percentage of silt (b) Percentage of Clay
DWR DWR
(c) Organic mater (Organic Carbon) (d) Bulk density
南苏丹Pongo流域土壤数据：Silt百分比，Clay百分比，有机物含量，Bulk Density。数据来源：ISRIC_SoilGrids。
DWR
南苏丹Pongo流域土地利用信息。土地利用数据来自USGS GLC数据。
DWR
FLDAS数据的覆盖范围，即气象站点位置信息。
DWR
FLDAS的气象驱动数据。

3.2.2 模型输入文件

输入文件列表

File Category 备注 Header # of column
.mesh sp 三角形单元定义 是
.att sp 单元属性表 是
.riv sp 河流 是
.rivseg sp 三角形单元与河道单元拓扑信息 是
.calib cfg 物理参数校准文件 是
.para cfg 模型运行设置参数文件 是
.ic cfg 模型初始条件 是
.geol para 地下水层的水力学参数 是
.soil para 浅层土壤的水力学参数 是
.lc para 土地利用的水力学和能量参数 是
.forc tsd 气象驱动数据文件的列表 是
.csv tsd 气象驱动时间序列数据 是
.lai tsd 叶面积指数 是
.obs tsd 用于校准观测时间序列数据 是
.mf tsd 融雪参数时间序列数据 是
.rl tsd 计算潜在蒸散发的粗糙度时间序列数据 是
gis/domain 三角形单元的矢量文件 x x
gis/river Shapefile 河段的矢量文件 x x
gis/seg Shapefile 被三角形分割的河段矢量文件 x x
The screenshot of input files for SHUD
The screenshot of input files for SHUD

在文件夹 gis 和 fig 中的文件与SHUD模拟无关，然后在数据分析和可视化中有很大作用，因此保留在输入文件当中。

3.2.3 空间数据

3.2.3.1 .sp.mesh 文件

Example of .sp.mesh file (block 1) Example of .sp.mesh file (block 2)

.mesh 文件中有两个表格，第一个表格为三角形的组成和拓扑关系；第二个表格为组成三角形的节点坐标信息。

  • 表格 1 (三角形单元的几何定义)

  • 表头

值1 值2
行数 ( \(N_{cell}\)) 列数 (\(8\))
  • 表
列 含义 取值范围 单位 备注
ID 单元序号 \(i\) 1 ~ \(N_{cell}\) -
Node1 三角形单元第1角点\(i\) 1 ~ \(N_{node}\) -
Node2 三角形单元第2角点 \(i\) 1 ~ \(N_{node}\) -
Node3 三角形单元第3角点 \(i\) 1 ~ \(N_{node}\) -
Nabr1 三角形单元第1邻居序号\(i\) 1 ~ \(N_{cell}\) -
Nabr2 三角形单元第2邻居序号\(i\) 1 ~ \(N_{cell}\) -
Nabr3 三角形单元第3邻居序号\(i\) 1 ~ \(N_{cell}\) -
Zmax 三角形单元z中心点地表高程 -9999 ~ +inf \(m\)

  • Block 2 (三角形的三个角点的信息)

  • 表头:

值1 值2
行数( \(N_{node}\)) 列数 (\(5\))
  • Table
列 含义 取值范围 单位 备注
ID 点序号 \(i\) 1 ~ \(N_{cell}\) -
X x坐标 \(i\) 1 ~ \(N_{node}\) -
Y y坐标 \(i\) 1 ~ \(N_{node}\) -
AqDepth 含水层厚度(地表至不透水层距离)\(i\) 0 ~ +inf \(m\)
Elevation 点的高程（z坐标） \(i\) -9999 ~ +inf \(m\)

3.2.3.2 .sp.att文件

Example of .sp.att file
Example of .sp.att file
  • 表头
值1 值2
行数( \(N_{cell}\)) 列数 (\(7\))
  • Table
列 含义 取值范围 单位 备注
ID 单元序号\(i\) 1 ~ \(N_{cell}\) -
SOIL 土壤类型序号 1 ~ \(N_{soil}\) -
GEOL 类型序号 1 ~ \(N_{geol}\) -
LC 类型序号 1 ~ \(N_{lc}\) - \(N_{lc}\) = \(N_{lai}\)
FORC 气象站点序号 1 ~ \(N_{forc}\) -
MF 融雪指数序号 1 ~ \(N_{mf}\) -
BC 边界条件序号 1 ~ \(N_{bc}\) -
SS 源汇序号 1 ~ \(N_{bc}\) -

3.2.3.3 .sp.riv文件

Example of .sp.riv file
Example of .sp.riv file
  • 表头
值1 值2
行数( \(N_{riv}\)) 列数 (\(5\))
  • 表
列 含义 取值范围 单位 备注
ID 河段序号\(i\) 1 ~ \(N_{river}\) -
DOWN 下游河段序号 1 ~ \(N_{river}\) - Negative vlaue indicates outlet
Type 河流参数序号 1 ~ \(N_{rivertype}\) -
Slope 河床底坡度 -10 ~ 10 \(m/m\) Height/Length
Length 河段长度 \(i\) 0 ~ inf \(m\)

3.2.3.4 .sp.rivseg文件

Example of .sp.rivseg file
Example of .sp.rivseg file
  • 表头
值1 值2
行数( \(N_{segment}\)) 列数 (\(4\))
  • Table
列 含义 取值范围 单位 备注
ID 片段序号 \(i\) 1 ~ \(N_{segment}\) -
iRiv 所属河段序号 1 ~ \(N_{river}\) -
iEle 相交单元序号 1 ~ \(N_{cell}\) -
Length 片段长度 \(i\) 0 ~ inf \(m\)

3.2.4 水力学参数

3.2.4.1 .para.soil

土壤表层的水力学参数

3.2.4.2 .para.geol

地下水层的水力学参数

3.2.4.3 .para.lc

与土地利用有关的参数

3.2.5 模型配置文件

3.2.5.1 .cfg.para文件

Example of .cfg.para file
Example of .cfg.para file
  • Table
列 含义 取值范围 单位 推荐值
VERBOSE 输出冗余信息 - - 0
INIT_MODE 初始条件模式 0~3 - 3 (0=Relief condition, 1=Dry condition, 2=Default guess, 3=Warm start)
ASCII_OUTPUT 是否输出文本格式结果 1/0 - 0
Binary_OUTPUT 是否输出二进制格式结果 1/0 - 1
CRYOSPHERE 冰冻圈模块是否开启 1/0 - 0
SPINUPDAY 模型预热天数 0 ~ inf \(day\) 0
SCR_INTV 屏幕输出间隔 0 ~ \(N_{threads}\) \(min\) 1440
ABSTOL CVODE绝对容差 1e-6 ~ 0.1 - 0.0001
RELTOL CVODE相对容差 1e-6 ~ 0.1 - 0.0001
INIT_SOLVER_STEP 初始迭代步长 - \(min\) 1
MAX_SOLVER_STEP 最大迭代步长 1~60 \(min\) 10
ET_STEP 正散发计算步长 1~360 \(min\) 60
START 模型开始时间（天数） 0 ~ inf \(day\) 0
END 模型结束时间（天数） - \(day\) -
dt_ye_snow 储量输出步长，积雪 0 ~ inf \(min\) 1440
dt_ye_surf 储量输出步长，地表水 0 ~ inf \(min\) 1440
dt_ye_unsat 储量输出步长，未饱和层 0 ~ inf \(min\) 1440
dt_ye_gw 储量输出步长，地下水层 0 ~ inf \(min\) 1440
dt_Qe_surf 单元流量输出步长，地表水 0 ~ inf \(min\) 1440
dt_Qe_sub 单元流量输出步长， 地下水 0 ~ inf \(min\) 1440
dt_qe_et0 单元流量输出步长，截流蒸发 0 ~ inf \(min\) 1440
dt_qe_et1 单元流量输出步长，蒸腾 0 ~ inf \(min\) 1440
dt_qe_et2 单元流量输出步长，蒸发 0 ~ inf \(min\) 1440
dt_qe_etp 单元流量输出步长，潜在蒸散发 0 ~ inf \(min\) 1440
dt_qe_prcp 单元流量输出步长，降雨 0 ~ inf \(min\) 1440
dt_qe_infil 单元流量输出步长，下渗 0 ~ inf \(min\) 1440
dt_qe_rech 单元流量输出步长，地下水补给 0 ~ inf \(min\) 1440
dt_yr_stage 河段储量输出步长，河道 0 ~ inf \(min\) 1440
dt_Qr_down 河段流量输出步长，向下游 0 ~ inf \(min\) 1440
dt_Qr_surf 河段流量输出步长，向坡面 0 ~ inf \(min\) 1440
dt_Qr_sub 河段流量输出步长，向地下水 0 ~ inf \(min\) 1440
dt_Qr_up 河段流量输出步长，向上游 0 ~ inf \(min\) 1440
dt_Lake 湖泊相关变量的输出步长 0 ~ inf \(min\) 1440

注：

  • INIT_MODE:
    • 0 = 含水层完全饱和作为初始条件,
    • 1 = 含水层完全干涸,
    • 2 = 未饱和层水量为含水层40%，饱和层水储量为含水层的40%,
    • 3 = 读取.cfg.ic作为初始条件

3.2.5.2 .cfg.calib文件

Example of .cfg.calib file
Example of .cfg.calib file
  • Table
列 含义 取值范围 单位 备注
GEOL_KSATH 水平水力传导度，地下水层 ? -
GEOL_KSATV 垂直水力传导度，地下水层 ? -
GEOL_KMACSATH 大孔隙水平水力传导度，地下水层 ? -
GEOL_DMAC 大孔隙深度，地下水层 -
GEOL_THETAS 地下水层孔隙度，地下水层 -
GEOL_THETAR 土壤残留含水量，地下水层 -
GEOL_MACVF 大孔隙面积比，地下水层 -
SOIL_KINF 饱和水力传导度，未饱和层 ? -
SOIL_KMACSATV 大孔隙垂向水力传导度，未饱和层 ? -
SOIL_DINF 下渗深度参数，未饱和层 ? -
SOIL_ALPHA \(\alpha\)值，van Genuchten公式 -
SOIL_BETA \(\beta\)值，van Genuchten公式 -
SOIL_MACHF 大孔隙面积比，未饱和层 -
LC_VEGFRAC 植被覆盖度 -
LC_ALBEDO 反照率 -
LC_ROUGH 地表曼宁粗糙度 -
LC_SOILDGD 土壤劣化系数 -
LC_IMPAF 不透水面积比 -
LC_ISMAX 最大截流系数 -
AQ_DEPTH+ 含水层深度 \(m\)
TS_PRCP 降雨 -
TS_SFCTMP+ 温度 \(C\)
ET_ETP 潜在蒸散发 -
ET_IC 冠层截流 -
ET_TR 植被蒸腾 -
ET_SOIL 直接蒸发 -
RIV_ROUGH 河道曼宁粗糙度 -
RIV_KH 河床水力传导度 -
RIV_DPTH+ 河道深度 \(m\)
RIV_WDTH+ 河道宽度 \(m\)
RIV_SINU 河道绵延度 -
RIV_CWR 谢才(Chezy)公式系数\(C_{wr}\) -
RIV_BSLOPE+ 河床坡度 \(m/m\)
IC_GW+ 地下水水位初始条件 \(m\)
IC_RIV+ 河道水量初始条件 \(m\)

注：

  • 表格中+代表该值的含义为该值与模型参数相加得到新参数。即\(\theta ^* = \theta + \zeta\), \(\theta\)为参数， \(\zeta\)为校准调节值。
  • 除+以外的所有值，都是与参数相乘关系。即\(\theta ^* = \theta \times \zeta\)。

3.2.5.3 .cfg.ic文件

Example of .cfg.ic file
Example of .cfg.ic file

  • 表1 (三角单元初始条件)

  • 表头

值1 值2
行数( \(N_{cell}\)) 列数 (\(6\))
  • Table
列 含义 取值范围 单位 备注
ID 单元序号\(i\) 1 ~ \(N_{cell}\) -
Canopy 冠层截流 \(i\) 0 ~ inf \(m\)
Snow 积雪 \(i\) 0 ~ inf \(m\)
Surface 地表水 \(i\) 0 ~ inf \(m\)
Unsat 未饱和层水 \(i\) 0 ~ inf \(m\)
GW 地下水 \(i\) 0 ~ inf \(m\)

  • 表2 (河段初始条件)

  • 表头:

值1 值2
行数( \(N_{riv}\)) 列数 (\(2\))
  • Table
列 含义 取值范围 单位 备注
ID 河段序号 \(i\) 1 ~ \(N_{riv}\) -
Stage 河段水位 \(i\) 0 ~ inf \(m\)

  • 表3 (湖泊初始条件)

  • 表头:

值1 值2
行数( \(N_{lake}\)) 列数 (\(2\))
  • Table
列 含义 取值范围 单位 备注
ID 湖泊序号 \(i\) 1 ~ \(N_{lake}\) -
Stage 湖泊水位 \(i\) 0 ~ inf \(m\)

3.2.6 时间序列数据

3.2.6.1 .tsd.forc文件

Example of .tsd.forc file
Example of .tsd.forc file
  • Line 1: 气象站点数量 | 开始日期 (YYYYMMDD)
  • Line 2: 驱动数据文件所在文件夹
  • Line 3~N: 时间序列数据文件名
Example of .csv forcing file
Example of .csv forcing file
  • 表头:
值1 值2
( \(0\)) 列数 (\(6\))
  • Table
列 含义 取值范围 单位 备注
Day 日数 0 ~ \(N_{day}\) \(day\)
PRCP 降雨 0 ~ 500 \(mm/day\) 降雨强度超过阈值模型会警告，但不妨碍运行
TEMP 气温 -100 ~ 70 \(C\)
RH 相对湿度 0 ~ 1 \(-\)
wind 风速 0 ~ inf \(m/s\) 非负
Rn 太阳净辐射 ? \(W/m^2\) 最大值为太阳常数
Press 气压 ? \(Pa\) 可缺省

3.2.6.2 .tsd.lai文件

Example of .tsd.lai file
Example of .tsd.lai file
  • 表头:
值1 值2 值3
天数 ( \(N_{time}\)) 列数 (\(N_{lc}\)) 开始日期 (YYYYMMDD)
  • Table
列 含义 取值范围 单位 备注
第一列 时间 0 ~ \(N_{time}\) \(day\)
第2列 叶面积指数1 0 ~ inf \(m^2/m^2\)
第i列 叶面积指数\(i-1\) 0 ~ inf \(m^2/m^2\)
… … … …

3.2.6.3 .tsd.mf文件

Example of .tsd.mf file
Example of .tsd.mf file

融雪度日模型(Degree-Day Model)的系数。

  • 表头:
值1 值2 值3
天数 ( \(N_{time}\)) 列数 (\(N_{mf}\)) 开始日期 (YYYYMMDD)
  • Table
列 含义 取值范围 单位 备注
第一列 时间 0 ~ \(N_{time}\) \(day\)
第2列 融雪因子1 0 ~ inf -
第i列 融雪因子\(i-1\) 0 ~ inf -
… … … …

3.2.6.4 .tsd.obs文件

Example of .tsd.obs file
Example of .tsd.obs file
  • 表头:
值1 值2 值3
天数 ( \(N_{time}\)) 列数 (\(N_{obs}\)) 开始日期 (YYYYMMDD)
  • Table
列 含义 取值范围 单位 备注
第一列 时间 0 ~ \(N_{time}\) \(day\)
第2列 观测1 ? ?
第i列 观测 \(i-1\) ? ?
… … … …

3.3 结果输出文件

3.3.1 文件名格式

为了方便结果的分类和用户理解，作者规则化定义了输出文件的文件名；用户通过文件名即可读懂文件所包含的变量类型、所属单元和读写格式。输出文件的文件名格式：

[项目名].[识别码].[格式]

或者：

[Project Name].[Identifier].[format]

  • [项目名]（[Project Name]） 是用户给模拟定的名称，通常为流域的名字简写。例如：黄河流域，可以写成hh或者huanghe.

  • [格式] （[format]）可以使用csv或者dat. csv为规则文本格式，可使用任意文本编辑器或者WPS表格打开。dat二进制格式，文件占用硬盘空间较小

  • [识别码] （[Identifier]） 定义了输出结果内涵。识别码基本结构为 [模型单元][变量类型][变量名]（[Model cell][Variable Type][Variable Name]）。 [Model cell] 是ele （element，模型三角单元）, riv （river，河段） or lak （lake，湖泊）. 变量类型包括 y, v and q，分别代表状态变量 （单位\(L\)）, 单位流量或速度（单位\(L^3 L^{-2} T^{-1}\)）和体积流量（in \(L^3 T^{-1}\)）——\(L\)为长度单位，\(T\)为时间单位。

下表为SHUD输出文件列表：

识别码 单元 类型 变量名 含义 单位
.eleyic. ele y ic 储量，冠层截流 \(m\)
.eleysnow. ele y snow 储量，雪水当量 \(m\)
.eleysurf. ele y surf 储量，地表水 \(m\)
.eleyunsat. ele y unsat 储量，未饱和层 \(m\)
.eleygw. ele y gw 储量，地下水水位 \(m\) .GW
.elevetp. ele v etp 通量，潜在蒸散发 \(\frac{m^3}{m^2 d}\)
.eleveta. ele v eta 通量，实际蒸散发 \(\frac{m^3}{m^2 d}\)
.elevetic. ele v etic 通量，冠层截流蒸发 \(\frac{m^3}{m^2 d}\)
.elevettr. ele v ettr 通量，植被蒸腾 \(\frac{m^3}{m^2 d}\)
.elevetev. ele v etev 通量，直接蒸发 \(\frac{m^3}{m^2 d}\)
.elevprcp. ele v prcp 通量，降雨 \(\frac{m^3}{m^2 d}\)
.elevnetprcp. ele v netprcp 通量，净降雨量 \(\frac{m^3}{m^2 d}\)
.elevinfil. ele v infil 通量，下渗 \(\frac{m^3}{m^2 d}\)
.elevexfil. ele v infil 通量，出渗（反下渗） \(\frac{m^3}{m^2 d}\)
.elevrech. ele v rech 通量，地下水补给 \(\frac{m^3}{m^2 d}\)
.eleqsurf. ele q surf 流量，地表径流（坡面流） \(m^3/d\)
.eleqsub. ele q sub 流量， 地下水水平流动 \(m^3/d\)
.rivystage. riv y stage 储量，河段水位 \(m\)
.rivqup. riv q up 流量，向上游 \(m^3/d\)
.rivqdown. riv q down 流量，向下游 \(m^3/d\)
.rivqsurf. riv q surf 流量，向坡面 \(m^3/d\)
.rivqsub. riv q sub 流量，向地下水 \(m^3/d\)

3.3.2 文本格式输出文件（.csv）

N - 输出文件列数（不含时间列） m - 步长总数 StartTime - 开始时间 （YYYYMMDD or YYYYMMDD.hhmmss）

N StartTime
\(T_1\) \(v_{1 \cdot 1}\) \(v_{1 \cdot 2}\) … \(v_{1 \cdot N}\)
\(T_2\) \(v_{2 \cdot 1}\) \(v_{2 \cdot 2}\) … \(v_{2 \cdot N}\)
\(T_3\) \(v_{3 \cdot 1}\) \(v_{3 \cdot 2}\) … \(v_{3 \cdot N}\)
… … … … …
\(T_{m}\) \(v_{m \cdot 1}\) \(v_{m \cdot 2}\) … \(v_{m \cdot N}\)

3.3.3 二进制输出文件（.dat）

二进制文件输出的结果核心矩阵与文本格式一致，但是使用了不同的存储结构。二进制文件格式紧凑，读写效率高，占用较小磁盘空间，因此作为SHUD模型的默认输出方式。 二进制文件格式不仅可被各类计算机语言（C/C++, Basic, Fortran，……）快速读取，并且也能在R、python、 Matlab、 Mathmetica等交互式语言中高效读取，如果读取有障碍，请联系模型作者获取帮助。

表格中，重要的变量是数据列数\(N_x\)。 \(N_x\)代表输出变量的列数，即该变量的输出计算单元数。 \(0 < N_x \le N_{max}\), 例如模型有\(N_{cell}\)个三角形单元，则.eleygw的地下水输出文件中\(N_{max} = N_{cell}\)；模型的河流径流结果中\(N_{max} = N_{riv}\)。

二进制输出文件现有v1.0, v2.0和v2.1三种格式，三种输出结果的格式见下表：

v 1.0 （\(N_x = N_{max}\)）

内容 字节数 类型 示例
数据列数 8 double 100
起始时间 8 double 20220101.120000
数据矩阵 8 double 1440, 1,2,…; 2880,2,4,…

v 2.0（\(N_x = N_{max}\)）

内容 字节数 类型 示例
文件描述 1024 char groundwater export, bla bla bla
起始时间 8 double 20220101.120000
数据列数 8 double 100
数据矩阵 8 double 1440, 1,2,…; 2880,2,4,…

v 2.1 以上版本（\(0 < N_x \le N_{max}\)）

内容 字节数 类型 示例
版本号 8 double 2.0
文件描述 2048 char groundwater export, bla bla bla
起始时间 8 double 20220101.120000
数据列数 8 double 100
索引ID 8 double 1,3,5,…;
数据矩阵 8 double 1440, 1,2,…; 2880,2,4,…
v2.1版本输出文件的逻辑和内存结构
v2.1版本输出文件的逻辑和内存结构

注：double为计算机的双精度数据格式，默认长度为8个比特位

数据为\(m \times N_x\)的二维矩阵。

  • 矩阵中的值全部为8字节double型。
  • \(m\)为结果输出的最终时间长度，\(m=100\)表明模型共输出100个时间步长的数据，则数据矩阵的行数为100。
  • \(N_x\)为输出结果的数量。 在v2.1之后的版本中，模型可灵活输出指定序号的结果，因此在v2.1结果中，结果矩阵的列数并不一定等于模拟单元的数量；例如模拟三角形单元数为1000(\(N_{max} = 1000\))，模型配置中指定值输出偶数序号的结果，则\(N_x = 500\)，索引ID则为\(i_{o} = 2,4,6,...,998,1000\)。
  • \(i\)为矩阵列序号，取值范围：\(0 \le i \le N_{x}\)，第一列(\(i=0\))为时间，\(i=1~N_x\)则为输出结果值。
  • \(j\)为时间序号，取值范围：$0 j m $。
  • \(N_{h}\)为文件头的长度，换算为8字节double型的数量。

数据矩阵的存储方式见下表：

文件存储位置 矩阵行号(\(j=\)) 数据列序号(\(i=\)) 值 值含义
\(N_{h} + 1\) 1 0 \(T_1\)
\(N_{h} + 2\) 1 1 \(v_{1 \cdot 1}\)
\(N_{h} + 3\) 1 2 \(v_{1 \cdot 2}\)
… 1 … …
\(N_{h} +N_x\) 1 \(N_x\) \(v_{1 \cdot N_x}\)
\(N_{h} +（N_x+1） * (2-1) + i\) 2 0 \(T_2\)
\(N_{h} +（N_x+1） * (2-1) + i\) 2 1 \(v_{2 \cdot 1}\)
\(N_{h} +（N_x+1） * (2-1) + i\) 2 2 \(v_{2 \cdot 2}\)
\(N_{h} +（N_x+1） * (2-1) + i\) 2 … …
\(N_{h} +（N_x+1） * (2-1) + i\) 2 \(N_x\) \(v_{2 \cdot N_x}\)
\(N_{h} +（N_x+1） * (3-1) + i\) 3 0 \(T_3\)
\(N_{h} +（N_x+1） * (3-1) + i\) 3 1 \(v_{3 \cdot 1}\)
\(N_{h} +（N_x+1） * (3-1) + i\) 3 2 \(v_{3 \cdot 2}\)
\(N_{h} +（N_x+1） * (3-1) + i\) 3 … …
\(N_{h} +（N_x+1） * (3-1) + i\) 3 \(N_x\) \(v_{3 \cdot N_x}\)
\(N_{h} +（N_x+1） * (j-1) + i\) \(j\) … …
\(N_{h} +（N_x+1） * (j-1) + i\) \(j\) … …
\(N_{h} +（N_x+1） * (j-1) + i\) \(j\) … …
\(N_{h} +（N_x+1） * (j-1) + i\) \(j\) … …
… … … …
\(N_{h} +（N_x+1） * （m-1） + i\) m 0 \(T_{m}\)
\(N_{h} +（N_x+1） * （m-1） + i\) m 1 \(v_{m \cdot 1}\)
\(N_{h} +（N_x+1） * （m-1） + i\) m 2 \(v_{m \cdot 2}\)
\(N_{h} +（N_x+1） * （m-1） + i\) m … …
\(N_{h} +（N_x+1） * （m-1） + i\) m \(N_x\) \(v_{m \cdot N_x}\)