SWAT 调参方法思路整理【3】

封面图ID

书接上文1

书接上文2

时隔一个月，回来继续整理编程思路。

项目开源地址

大体思路

在确定工作路径没有异常后存储参数和观测数据。根据样本参数生成样本，将每一个样本视为一个对象，完成拷贝，替换参数，执行，计算精度指标操作后汇总至Excel中。

模块

check()

check模块除了正常初始化外还要负责两个工作：测试运行模型以确保RawModel内的模型和Exes内部的文件可以正常运行，通过Parameters.txt和传递软件参数。

Parameters.txt具有一定的格式，名称行下一行为参数值，通过if判断所在行内容和参数，如遇多行内容则以-或=为中止符。如遇到土壤参数则查找对应土壤参数对应的csv。观测值文件在Parameters.txt内有定义，遵循原有格式进行扫描输入。这部分模块存在简易确定软件参数合理性功能。

测试运行方法为通过临时创建一个文件夹，并将原始模型和SWAT模型可执行文件复制到文件夹内。通过检测SWAT模拟完成文件`fin.fin`状态确定模拟是否成功。

sampleGeneration()

sampleGeneration为软件参数提供的参数进行模拟，根据上文2，可划归HRU参数，子流域参数，流域参数和土壤参数。土壤参数为一组，而其他三个参数则各有其尺度，方法可根据不同尺度和不同范围进行生成样本。

采用pyDOE2的内置方法进行拉丁超立方采样，使用两个参数，维度和每个维度内的样本数，得出应为m*n个样本。该方法返回值为0-1，因此还需要进行处理。将上下边界相减，乘以采样值加上下边界即可。

$Var_{tar}=Var_{sample}*(B_{upper}-B_{lower})+B_{lower}$

此处输出为一多维list，每条为一个样本，内部有多种参数。

ProcessSWAT()类

ProcessSWAT()是包括传递参数，模拟和计算的模块。其中传入参数主要为sampleGeneration()生成的参数。初始化模块过程包括创建实例文件夹，将参数输入至类中，同时创建变量记录各个尺度下参数数量以备使用。

修改参数模块为核心内容。出于模型物理意义考虑，将土壤参数单独剥离出来。在另一方面，土壤参数格式与其他格式也略有不同。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22

.Sol file Watershed HRU:2 Subbasin:1 HRU:2 Luse:FRST Soil: Haplic Luvisols2 Slope: 15-25 2022/5/7 0:00:00 ArcSWAT 2012.10_2.19
Soil Name: Haplic Luvisols2
Soil Hydrologic Group: C
Maximum rooting depth(mm) : 1000.00
Porosity fraction from which anions are excluded: 0.500
Crack volume potential of soil: 0.500
Texture 1                : L-GR--SL-GR--SL-FSL
Depth                [mm]:      300.00     1000.00
Bulk Density Moist [g/cc]:        1.51        1.52
Ave. AW Incl. Rock Frag  :        0.17        0.13
Ksat. (est.)      [mm/hr]:       15.62        4.04
Organic Carbon [weight %]:        0.74        0.36
Clay           [weight %]:       22.00       29.00
Silt           [weight %]:       37.00       34.00
Sand           [weight %]:       41.00       37.00
Rock Fragments   [vol. %]:        4.00        3.00
Soil Albedo (Moist)      :        0.01        0.01
Erosion K                :        0.17        0.17
Salinity (EC, Form 5)    :        0.10        0.10
Soil pH                  :        6.50        6.70
Soil CACO3               :        0.00        0.00
                             

土壤名称在第一行中即出现，可作为索引进行查找。土壤参数具有多层属性，由土壤数据库进行定义，需要根据现实情况进行修改。土壤参数本身不应使用v方法进行修改，因此在采用r方法修改时可以每一层都进行修改。参数名字与土壤数据库中略有不同，但可根据官方文档确定。

其他参数文件格式则相对统一，具有相似的格式。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17

.gw file Watershed HRU:1 Subbasin:1 HRU:1 Luse:FRST Soil: Dystric Podzoluvisols2 Slope: 5-15 2022/5/7 0:00:00 ArcSWAT 2012.10_2.19
      1000.0000     SHALLST : Initial depth of water in the shallow aquifer [mm]
      2000.0000     DEEPST : Initial depth of water in the deep aquifer [mm]
      61.500000     GW_DELAY : Groundwater delay [days]
       0.885000     ALPHA_BF : Baseflow alpha factor [days]
       0.763333     GWQMN : Threshold depth of water in the shallow aquifer required for return flow to occur [mm]
       0.025000     GW_REVAP : Groundwater "revap" coefficient
      45.916670     REVAPMN: Threshold depth of water in the shallow aquifer for "revap" to occur [mm]
       0.516281     RCHRG_DP : Deep aquifer percolation fraction
         1.0000     GWHT : Initial groundwater height [m]
         0.0030     GW_SPYLD : Specific yield of the shallow aquifer [m3/m3]
         0.0000     SHALLST_N : Initial concentration of nitrate in shallow aquifer [mg N/l]
         0.0000     GWSOLP : Concentration of soluble phosphorus in groundwater contribution to streamflow from subbasin [mg P/l]
         0.0000     HLIFE_NGW : Half-life of nitrate in the shallow aquifer [days]
         0.0000     LAT_ORGN : Organic N in the base flow [mg/L]
         0.0000     LAT_ORGP : Organic P in the base flow [mg/L]
         0.0100     ALPHA_BF_D : Baseflow alpha factor for deep aquifer [days]

gw文件，为HRU内地下水参数文件，可见其参数格式统一，各参数所占字符长度基本一致。对此可直接采用右侧说明进行查找。为排除说明中的参数引用导致的索引混乱，将：作为分隔符确保仅采用名称部分作为参数索引。

不同尺度下的修改方式也不同，采用初始化时的HRU结构建立目录结构。根据目录结构和修改尺度进行处理。最终得到修改后的模型。

执行模型方法和check()方法一致，不再赘述。

执行后提取结果文件，进行R^2和NSE计算，将结果输出。

printRes(results)方法

方法为输出结果到文件中，基本采用输出到文件的方法。

改进空间

1. 正则表达式的使用可有效提高输入输出操作；
2. 提高模型计算速度，参照SWAT-CUP，魔改了输入输出文件调用方法。

值得注意的事情

1. 不同SWAT版本计算结果偏差较大，且不同版本SWAT参数敏感性分析结果可能也不同；
2. 不同HRU大小下敏感性分析结果不同；
3. 讨论一种新的多区域调参方法，即：用子流域观测站控制调参该子流域后，在上一级子流域调参时该子流域参数不作改动。