Issuu

00_130_8 FHR reports

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 Implementation and testing of a framework for flexible hydrological modelling

www.flandershydraulicsresearch.be

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Tran, Q. Willems, P.; Pereira, F.; Nossent J.; Mostaert, F.

Cover figure © The Government of Flanders, Department of Mobility and Public Works, Flanders Hydraulics Research Legal notice Flanders Hydraulics Research is of the opinion that the information and positions in this report are substantiated by the available data and knowledge at the time of writing. The positions taken in this report are those of Flanders Hydraulics Research and do not reflect necessarily the opinion of the Government of Flanders or any of its institutions. Flanders Hydraulics Research nor any person or company acting on behalf of Flanders Hydraulics Research is responsible for any loss or damage arising from the use of the information in this report. Copyright and citation © The Government of Flanders, Department of Mobility and Public Works, Flanders Hydraulics Research 2021 D/2021/3241/023 This publication should be cited as follows: Tran, Q. Willems, P.; Pereira, F.; Nossent J.; Mostaert, F. (2021). Effect of climate change on the hydrological regime of navigable water courses in Belgium: Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling. Version 1.0. FHR Reports, 00_130_8. Flanders Hydraulics Research: Antwerp Reproduction of and reference to this publication is authorised provided the source is acknowledged correctly. Document identification Customer: Keywords (3-5): Knowledge domains: Text (p.): Confidential: Author(s):

Flanders Hydraulics research Ref.: WL2021R00_130_8 climate change; hydrological impact analysis; distributed & lumped models; Belgium; model uncertainty Watermanagement > hydrology > numerical modelling 94 ‫ ܈‬No

Appendices (p.): ‫ ܈‬Available online

Tran, Q. Willems, P.

Control Name Reviser(s):

Project leader:

Nossent, J.

Pereira, F.

Signature Getekend door:Jiri Nossent (Signature) Getekend op:2022-03-22 09:45:47 +01:0 Reden:Ik keur dit document goed

Getekend door:Fernando Pereira (Signat Getekend op:2022-03-22 08:48:17 +01:0 Reden:Ik keur dit document goed

Approval Mostaert, F. Head of Division:

F-WL-PP-16.03.01-02 Valid as from 6/11/2021

For the Head of Division, absent Patrik Peeters, Engineer, in charge with leading the Division Flanders Hydraulics Research

Digitaal ondertekend Patrik Peeters Patrik Peeters door (Signature) (Signature) Datum: 2021.09.08 08:15:16 +02'00'

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Abstract The innovative flexible hydrological modelling concept proposed by Tran et al. (2014) was further developed in this inventory. A generalized model structure which consists of the most basic and/or popular hydrological components was constructed. It includes the meteorological components, the storage components, the splitting components, the routing components and a component for closing water balance. Equations from the NAM, PDM and VHM models were used as demonstrations for each components. The additional components can be added by users depending on applications and on the available data. Since the concept allows flexibility in spatial resolution, the Grid-to-Grid approach for explicitly representation of the hydrological connectivity between gird cells with the catchment is needed. Next to that, the Source-to-Sink approach was also implemented but it was not discussed in details in this report due to its simplicity in implementation. Two schemes of Grid-to-Grid approach were developed: i) the normal routing scheme, in which the interaction between different subflow components are omitted and ii) the detailed routing scheme, which allows to describe more detail the infiltration and percolation processes. The step-wise methodology for calibration was developed for both lumped and spatial distributed models. For the distributed models, the calibration process starts from the well-calibrated model parameters for lumped case and then uses the catchment characteristics (i.e. the topography map, the land use map and the soil type map) and other additional data (e.g. the groundwater observation, the internal gauging station…) to disaggregate into spatial distributed parameters. Those spatial parameters can be fine-tuned to capture the distribution of the runoff in the whole catchment. The efficiency of the identified model structure is assessed using both traditional statistical indices and the graphical evaluation. Three conceptual model codes (NAM, PDM and VHM) spatially implemented in the proposed framework were applied for a case study of Grote Nete catchment in Belgium. They were calibrated using the predefined parameters from the lumped models. The comparison of the distributed outputs versus the lumped results and versus the observations at gauging stations indicated good matches, which show that the proposed framework has successfully implemented. Nevertheless, further research is still in need in order to complete the concept.

Final version

WL2021R00_130_8

III

F-WL-PP-16.03.01-02 Valid as from 6/11/2021

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Contents Abstract ............................................................................................................................................................ III Contents ............................................................................................................................................................ V List of tables..................................................................................................................................................... VII List of figures .................................................................................................................................................. VIII Abstract ............................................................................................................................................................. 1 1

Generalized hydrological model structure ................................................................................................ 3 1.1

1.1.1

Evapotranspiration component .................................................................................................... 5

1.1.2

Rainfall component (Rainfall contribution to surface storage) .................................................... 6

1.1.3

Net rainfall component ................................................................................................................. 6

1.2

Storage components ............................................................................................................................ 7

1.2.1

Surface storage component.......................................................................................................... 7

1.2.2

Soil storage component ................................................................................................................ 7

1.3

Meteorological components ............................................................................................................... 5

Splitting components ........................................................................................................................... 8

1.3.1

Rainfall constribution to slow flow (groundwater recharge) ....................................................... 8

1.3.2

Rainfall contribution to quick flow ............................................................................................... 8

1.3.3

Rainfall contribution to soil storage ............................................................................................. 9

1.4

Component to close water balance ................................................................................................... 10

1.5

Flow routing components .................................................................................................................. 10

1.6

Special component: Capillary rise component .................................................................................. 11

1.7

Sequence of solving equations .......................................................................................................... 11

Spatial linking........................................................................................................................................... 12 2.1

General considerations ...................................................................................................................... 12

2.2

The normal routing approach ............................................................................................................ 12

2.3

The detailed routing approach .......................................................................................................... 13

2.4

Modifications for detailed routing approach .................................................................................... 15

Model calibration method ....................................................................................................................... 17 3.1

Calibration concept for lumped models ............................................................................................ 17

3.1.1

Step 1: Prior time series processing............................................................................................ 18

3.1.2

Step 2: Step-wise model structure identification and calibration .............................................. 23

3.2

Calibration concept for distributed models....................................................................................... 23

3.2.1

Final version

Step-wise disaggregation concept .............................................................................................. 24

WL2021R00_130_8

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

3.2.2

Spatial disaggregation of model parameters based on spatial catchment characteristics ........ 24

3.2.3 Adjustment/calibration of the spatial disaggregation of model parameters based on spatial catchment characteristics ....................................................................................................................... 43 4

Model performance evaluation method ................................................................................................. 45 4.1

Graphical evaluation of runoff flows and response behavior at gauging stations ............................ 45

4.2

Statistical evaluation of runoff flows at gauging stations ................................................................. 49

4.3

Need for transformation of the flow variable: logarithmic or Box-Cox transformation ................... 51

4.4

Need to address serial dependence of flow ...................................................................................... 52

4.5

Evaluation of catchment response behavior ..................................................................................... 53

4.6

Effect of rainfall change ..................................................................................................................... 54

4.7

Evaluation of spatially variable model outputs ................................................................................. 55

Evaluation results of spatial models for Grote Nete catchment ............................................................. 62 5.1

Graphical evaluation of catchment behavior at the catchment outlet ............................................. 62

5.2

Spatial validation using internal gauging stations ............................................................................. 68

5.3

Spatial validation using groundwater observations .......................................................................... 74

5.4

Climate change effects on hydrological extremes given by distributed models ............................... 75

5.5

Spatial resolution effects on model performance ............................................................................. 78

Conclusions .............................................................................................................................................. 84

Bibliography ..................................................................................................................................................... 85

WL2021R00_130_8

Final version

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

List of tables Table 1: Soil hydraulic parameters characterized by soil textual classes (Liu and De Smedt, 2004) .............. 26 Table 2: Overview of the 47 land use classes in the NGI map and reclassification to the 9 main IGBP categories ......................................................................................................................................................................... 27 Table 3: Parameters characterizing land use classes in WetSpa ..................................................................... 29 Table 4: Potential runoff coefficient for different land use, soil type and topographical slope in WetSpa ... 30 Table 5: Runoff coefficient for different land use and vegetation types ........................................................ 31 Table 6: Runoff coefficient for some coarse topographical slope ranges, soil type classes and land use types ......................................................................................................................................................................... 32 Table 7: Runoff coefficient for different rainfall intensity classes, two hydrological soil type classes, and some land use types .................................................................................................................................................. 32 Table 8: Catchment averaged runoff coefficient for different classes of topographical slope, land use, soil type and geology ..................................................................................................................................................... 33 Table 9: CN-value for different land use types, agricultural crop types and hydrological soil type ............... 34 Table 10: Description of the hydrological soil types considered in the look-up table for CN ......................... 35 Table 11: Impervious percentage associated to specific land use types in WetSpa ....................................... 35 Table 12: Manning coefficient for different land uses .................................................................................... 37 Table 13: Combined lookup table for TAW_rootdepth as a function of soil type and root zone depth ........ 41 Table 14: NSE statistical performance values downstream the Grote Nete catchment for the calibration period (2003-2005) and validation period (2006-2008) using different models ............................................ 65 Table 15: Discharge stations in the Grote Nete catchment ............................................................................ 68 Table 16: NSE performance of different models for Tessenderlo and Vorst station ...................................... 73

Final version

WL2021R00_130_8

VII

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

List of figures Figure 1-1: Main components of the generalized hydrological model structure.............................................. 3 Figure 1-2: Generalized hydrological model structure where the QF and SF routing components are further detailed into two QF (OF and IF) and two SF (SF and deep SF) components .................................................... 5 Figure 2-1: Linking between two grid cells following the generalized hydrological model structure, and the normal routing approach ................................................................................................................................ 13 Figure 2-2: Linking between two grid cells following the generalized hydrological model structure, and the first attempt of detailed routing approach implementation .......................................................................... 14 Figure 2-3: Linking between two grid cells using the second attempt of detailed routing approach implementation ............................................................................................................................................... 14 Figure 2-4: Spatial linking between two grid cells using the modified detailed routing approach ................. 15 Figure 2-5: Spatial linking between two grid cells using the modified detailed routing approach with adaptation for the NAM model ....................................................................................................................... 16 Figure 2-6: Spatial linking between two grid cells using the modified detailed routing approach with adaptation for the PDM model ....................................................................................................................... 16 Figure 3-1: Steps in the lumped model structure identification and calibration procedure .......................... 17 Figure 3-2: Subflow filter results for the hourly river flows at Sint-Joris-Weert ............................................. 21 Figure 3-3: Parameters used in the criteria to select approximate independent peak flows from a river flow time series. ...................................................................................................................................................... 22 Figure 3-4: Different spatial modelling scales ................................................................................................. 24 Figure 3-5: Soil texture triangle chart .............................................................................................................. 26 Figure 3-6: Textural classes from Belgian system (A) and USDA system (B) ................................................... 27 Figure 3-7: Potential runoff coefficient for 74 gauged catchments in Flanders ............................................. 33 Figure 3-8: Simulation results (log-scale) of the Grid-to-Grid spatial linking model of the Dijle case for uniform parameter values (spatial disaggregation step 1) ........................................................................................... 38 Figure 3-9: Simulation results (log-scale) of the Grid-to-Grid spatial linking model of the Dijle case for uniform parameter values (spatial disaggregation step 1) ........................................................................................... 39 Figure 3-10: Spatial disaggregation of the quick runoff recession constants for overland flow (CKOF) and interflow (CKIF) (in hours) for the Dijle case ................................................................................................... 40 Figure 3-11: Spatial disaggregation of the maximum soil storage parameter (in mm) (without applying a controlling factor) for an example of the Dijle catchment: (left) based on the dominant land use type per 100m grid cell; (right) based on the weighted average of the parameter values obtained for the 5m land use map .................................................................................................................................................................. 42 Figure 3-12: Distribution of runoff coefficients over the Dijle catchment at different time steps of the simulation ........................................................................................................................................................ 43 Figure 3-13: Distribution of maximum soil storage parameter, before and after adjustments ..................... 44 Figure 4-1: Time series of total runoff flows, for the most extreme events in the simulated time series ..... 46

VIII

WL2021R00_130_8

Final version

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Figure 4-2: Time series of slow runoff flow (baseflow), for calibration and validation periods ..................... 46 Figure 4-3: Time series of cumulative total runoff volumes, for calibration and validation periods .............. 47 Figure 4-4: Scatterplot of peak flows .............................................................................................................. 47 Figure 4-5: Return period plot of peak flows, total simulation period............................................................ 48 Figure 4-6: Scatterplot of low flow extremes .................................................................................................. 48 Figure 4-7: Return period plot of low flow extremes, total simulation period ............................................... 49 Figure 4-8: Runoff coefficient based on overland flow (left) and combined interflow & overland flow (quick flow; right) versus relative soil moisture level downstream. .......................................................................... 53 Figure 4-9: Cumulative frequency distribution of the runoff coefficient for the NAM, PDM and VHM model based on overland flow (left) and combined interflow & overland flow (right) downstream. ...................... 53 Figure 4-10: Cumulative frequency distributions of the relative change in quick flow event volume per rainfall change class, for the Dijle catchment at Sint-Joris-Weert .............................................................................. 54 Figure 4-11: Quantile analysis for the relative change in overland peak flow versus change in rainfall, for the Dijle catchment at Sint-Joris-Weert ................................................................................................................ 55 Figure 4-12: Internal river flow gauging stations and calibration and validation groundwater wells considered by Vansteenkiste et al. (2012) for the Grote Nete catchment. ....................................................................... 55 Figure 4-13: Graphical evaluation of the temporal groundwater dynamics. Example of MIKE SHE head levels for calibration and validation wells in the Grote Nete catchment. ................................................................. 56 Figure 4-14: Graphical evaluation of the spatial groundwater dynamics. Example of MIKE SHE head levels along transects across the Grote Nete catchment. ......................................................................................... 57 Figure 4-15: Comparison of 60 cm depth averaged measured surface soil moisture content measured by 61 TDR sensors at weekly intervals, with saturation degree in the unsaturated zone simulated by a rainfall-runoff model, for the Weiherbach catchment in Germany ....................................................................................... 58 Figure 4-16: Correlations between measured runoff ratio and measured antecedent soil moisture for four catchments in the Walnut Gulch experimental watershed in Arizona, USA ................................................... 59 Figure 4-17: Soil moisture data from the ASCAT product rescaled to match the relative soil moisture simulated by the model, for Tiber river basin in Italy ...................................................................................................... 59 Figure 4-18: Example of runoff flow forecast results with observed flows, before and after model state updating based on data assimilation applied to the rescaled ASCAT soil moisture product, for Tiber river basin in Italy .............................................................................................................................................................. 60 Figure 4-19: Example of runoff flow forecast results with observed flows, before and after model state updating based on data assimilation applied to the rescaled ASCAT soil moisture product and change in NashSutcliffe efficiency, for Tiber river basin in Italy .............................................................................................. 60 Figure 4-20: Comparison of measured soil moisture based on DACOM sensors and model results based on the SPHY model for the Hunze catchment in Drente, The Netherlands ......................................................... 61 Figure 4-21: Runoff coefficient versus TDR measured 30 cm depth averaged soil moisture prior to the event, for the Rio Vauz basin in Italy .......................................................................................................................... 61 Figure 5-1: Total flow hydrograph simulated by lumped NAM model and spatial NAM model versus observation at the outlet station .................................................................................................................... 62 Figure 5-2: Total flow hydrograph simulated by lumped PDM model and spatial PDM model versus observation at the outlet station .................................................................................................................... 63

Final version

WL2021R00_130_8

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Figure 5-3: Total flow hydrograph simulated by lumped VHM model and spatial VHM model versus observation at the outlet station .................................................................................................................... 63 Figure 5-4: Observed and simulated discharges downstream the Grote Nete catchment for the extreme event in January 2003 for lumped and distributed NAM, PDM and VHM ................................................................ 64 Figure 5-5: Observed and simulated discharges downstream the Grote Nete catchment for the extreme event in November 2004 for lumped and distributed NAM, PDM and VHM ........................................................... 64 Figure 5-6: Observed and simulated discharges downstream the Grote Nete catchment for the extreme event in January 2007 for lumped and distributed NAM, PDM and VHM ................................................................ 64 Figure 5-7: Observed and simulated subflows at the outlet station using the spatial version of PDM, NAM and VHM models .................................................................................................................................................... 65 Figure 5-8: Observed and simulated cumulative discharge and baseflow at the catchment outlet using the spatial version of the PDM, NAM and VHM models ....................................................................................... 66 Figure 5-9: Peak flows versus return period at the outlet station for the spatial version of NAM, PDM and VHM models .................................................................................................................................................... 66 Figure 5-10: Simulated versus observed nearly independent peak flows at the outlet station after Box-Cox transformation (λ=0.25) for the spatial version of NAM, PDM and VHM models .......................................... 67 Figure 5-11: Low flow extremes versus return period at the outlet station for the spatial version of NAM, PDM and VHM models ............................................................................................................................................. 67 Figure 5-12: Simulated versus observed nearly independent low flow extremes at the outlet station after BoxCox transformation (λ=0.25) for the spatial version of NAM, PDM and VHM models ................................... 68 Figure 5-13: Total discharge simulated by NAM, PDM and VHM spatial models at Meerhout (top), Tessenderlo (middle) and Vorst (bottom) stations .............................................................................................................. 69 Figure 5-14: Peak flows versus return period at Meerhout (top), Tessenderlo (middle) and Vorst (bottom) stations for NAM, PDM and VHM spatial models ........................................................................................... 70 Figure 5-15: Simulated versus observed nearly independent hourly peak flows at Meerhout (top), Tessenderlo (middle) and Vorst (bottom) stations for NAM, PDM and VHM spatial models ........................ 71 Figure 5-16: Observed and simulated cumulative discharge at Meerhout (top), Tessenderlo (middle) and Vorst (bottom) stations for NAM, PDM and VHM spatial models .................................................................. 72 Figure 5-17: Low flow extremes versus return period at Tessenderlo (top) and Vorst (bottom) stations for NAM, PDM and VHM spatial models............................................................................................................... 73 Figure 5-18: Comparison between observed and calculated groundwater height at different piezometers 74 Figure 5-19: Hourly peak flows versus return period at the catchment outlet for the low, mean and high climate scenarios simulated with lumped and spatial versions of NAM, PDM and VHM models for 2050 ... 75 Figure 5-20: Hourly peak flows versus return period at the catchment outlet for the low, mean and high climate scenarios simulated with lumped and spatial versions of NAM, PDM and VHM models for 2100 ... 76 Figure 5-21: Hourly low flows versus return period at the catchment outlet for the low, mean and high climate scenarios simulated with lumped and spatial versions of NAM, PDM and VHM models for 2050 ................ 77 Figure 5-22: Hourly low flows versus return period at the catchment outlet for the low, mean and high climate scenarios simulated with lumped and spatial versions of NAM, PDM and VHM models for 2100 ................ 78 Figure 5-23: Total discharge simulated by NAM, PDM and VHM distributed models at outlet station using three different spatial resolutions................................................................................................................... 79

WL2021R00_130_8

Final version

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Figure 5-24: Observed and simulated hydrographs of the extreme events in January 2003 (top) and in January 2007 (bottom) by three models using 50m, 250m and 1km resolution maps................................................ 80 Figure 5-25: Observed (in back line) and simulated cumulative discharge given by three models using 50m (in red line), 250m (in blue line) and 1km (in green line) resolution maps .......................................................... 80 Figure 5-26: Observed and simulated peak distribution given by three models using 50m, 250m and 1km resolution maps ............................................................................................................................................... 81 Figure 5-27: Observed and simulated peak distributions for Meerhout (top), Tessenderlo (middle) and Vorst (bottom) stations for three distributed models using 50m, 250m and 1km resolution maps ....................... 82 Figure 5-28: Observed and simulated peak distributions for Tessenderlo (top) and Vorst (bottom) stations for three distributed models using 50m, 250m and 1km resolution maps .......................................................... 83

Final version

WL2021R00_130_8

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Abstract As concluded from the previous report on this project (Tran et al., 2021), the use of a single model for impact analysis in water management might not be sufficient to reach complete information on the impacts of given scenarios (e.g. climate change). An ensemble approach therefore was recommended to determine the most appropriate model structure and model details based on applications and available data. In this approach, the different plausible hydrological models and model assumptions are applied and inter-compare in order to quantify the impact of model uncertainty on climate change and other types of impact results (Vansteenkiste et al., 2012b; 2014). Therefore, it is a significant demand of a modelling system that allows to test different model structures and model assumptions and to assess their impact results. Such of that system has been proposed as an innovative flexible hydrological modelling concept and preliminarily tested by Tran et al. (2021). The proposed concept consists of several benefits including the flexibility not only in model structure (components and equations) but also in spatial resolution. It also couples the data based approach that helps to maximize the extracted information from the field measurements. Moreover, a model calibration methodology using catchment characteristics and a multi-criteria evaluation for model performance were integrated in this concept. The presented flexible hydrological modelling concept was then implemented in Python programming language using the PCRaster package. The first tests were done with three conceptual hydrological model codes, namely NAM, PDM and VHM. The preliminary results showed a very promising potential for a successful implementation of the concept. Furthermore, a further development of the concept of the framework for flexible hydrological modelling has been carried out and carefully tested by applying to Grote catchment, Belgium. The work is described in detail in this report as follows. The first chapter introduces the generalized hydrological model structure that was constructed and which consists of four basic components: meteorological components, storage components, splitting components and routing components. The generalized structure, moreover, allows users to add or remove other components if necessary. Other existing hydrological models (HBV or WetSpa for example) may also transform into the proposed structure. In the second chapter, two spatial linking approaches, representing the hydrological connectivity between the sources of runoff generated by the generalized structure and the downward locations are also described. The water is able to not only route from a grid cell to the neighboring downward cell and to the river cells but also move between different subflow components. Several modifications to adapt the spatial linking approach to specific model structures are also carried out. The calibration method for lumped and distributed models is described in detail in chapter 3. It gives the possibility to calibrate a spatial distributed model from a well calibrated lumped model while the consistency of the simulation results for different spatial resolution is maintained. Although extra calibration is hard to avoid when switching resolutions, the suggested method ensures an easier and less time consuming tasks for the users. The proposed framework was then applied for Grote Nete case study using three based model codes of NAM, VHM and PDM. Both lumped and distributed versions of these models were built within the framework. Their results were then validated by the model performance evaluation method, which can be found in chapter 4. This method makes use of both the classical method using statistical goodness of fit and the graphical (visual) method for a global assessment of model performance. The evaluation results for the Grote Nete case study can be accessed in Chapter 5.

Final version

WL2021R00_130_8

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

NOTE The methodology outlined in this report is the result of joint research activities by KU Leuven under supervision of Prof. Patrick Willems for the Flemish Environment Agency (VMM) on “Next Generation hydrological modelling” in cooperation with IMDC and the present research project for Flanders Hydraulics Research (Waterbouwkundig Laboratorium) on “Effect of climate change on the hydrological regime of navigable water courses in Belgium”. See also our joint concept note in which next to the Grote Nete catchment results have also been reported for the Dijle catchment: Tran, Q. Q., Willems, P., 2014. Next Generation Tool for Flexible Hydrological Modelling - Concept note. KU Leuven – Afdeling Hydraulica i.s.m. IMDC, voor Vlaamse Milieumaatschappij – Afdeling Operationeel Waterbeheer & Waterbouwkundig Laboratorium, oktober 2014 / KU Leuven – section Hydraulics i.c.w. IMDC, for Flemish Environment Agency – division Operational Water Management & Flanders Hydraulics Research (WL), October 2014.

WL2021R00_130_8

Final version

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

1 Generalized hydrological model structure Based on the literature findings on the existing hydrological modelling summarized in the previous report (Tran et al., 2021), a generalized model structure is proposed in this section. The generalization is primarily based on the model structures of the existing lumped conceptual models PDM, NAM, VHM. Given that the different existing models have similar model components, their interlinks between these components are different (see Tran et al., 2021). Therefore, the interlinks in the generalized structure should also not be fixed. Only the components were generalized, but with equations that can be based on variables obtained as outputs from any of the components. This means that the equations selected by the users for each of the model components will determine the interlinks between the components. Figure 1-1 schematizes the generalized model structure and the different model components. They can be classified in four main types (which were given different colors in the figure): • • • •

Meteorological components (orange colour in the figure); Storage components (aqua colour in the figure); Splitting components (blue colour in the figure); Routing components (purple colour in the figure).

There is also an additional special type of components that describes how the water balance is closed in the entire model.

Figure 1-1: Main components of the generalized hydrological model structure

Final version

WL2021R00_130_8

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Each of the different components is more detailed described in the coming paragraphs, together with the equations that are used in PDM, NAM and VHM. The snow component has been left out for now. The descriptions make use of the following symbols for the different variables: x

catchment rainfall

net rainfall available for catchment runoff after surface storage

xQF

catchment rainfall amount that contributes to quick flow

xOF

catchment rainfall amount that contributes to overland flow

xIF

catchment rainfall amount that contributes to interflow

xSF

catchment rainfall amount that contributes to slow flow (groundwater recharge)

catchment potential evapotranspiration

catchment actual evapotranspiration

surface storage

soil storage

catchment outflow (total or for one of the subflows: yQF, yOF, yIF, ySF, …)

The following abbreviations are used: P

precipitation

ETp

potential evapotranspiration

ETa

actual evapotranspiration

quick flow

overland flow

interflow

slow flow

In principle, each component produces only one output variable (that is modelled by the component) but might require one or several input variables (mostly output variables of other components). If the time (t) variable is not specified in the equations, it means that the equations are applied to each time step. Note that the generalized model structure of Figure 1-1 can be further expanded, by replacing components of a given type to more components of the same type. An example is shown in Figure 2-2, where the QF routing component is further detailed into two QF routing components: one for OF and another one for IF. Also the SF routing component is further detailed in two SF routing components: a first SF routing component, and a deep SF routing component. The user has the flexibility to choose how many components he wants for the QF and SF routing components.

WL2021R00_130_8

Final version

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Figure 1-2: Generalized hydrological model structure where the QF and SF routing components are further detailed into two QF (OF and IF) and two SF (SF and deep SF) components

1.1 Meteorological components 1.1.1

Evapotranspiration component

The evapotranspiration model converts the potential evapotranspiration input into actual evapotranspiration, taking into account the soil water availability. In PDM, the evapotranspiration model depends on the soil storage u by a power relation with exponent be: be   u   max − u    ea = e p 1 −   umax     

(1)

In NAM, the following evapotranspiration model is applied: when ep>s

ea = e p ea =

u (e p − s ) u max

otherwise

(2) (3)

When be = 1 is applied to (1), model equation (3) is obtained. This means that the evapotranspiration model in NAM can be seen as a special case of the evapotranspiration model in PDM. Hence, only (1) needs to be implemented in the proposed framework to cover both the PDM and NAM evapotranspiration models.

Final version

WL2021R00_130_8

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

In VHM, also (3) is taken as the most parsimonious (linear) evapotranspiration model that is tested first. A more advanced version is: ea = f e

u u evap

e p when u < uevap

ea = fe e p otherwise

(4) (5)

Where fe is the linear fraction between actual and potential evapotranspiration. In the simplest VHM case, model (1) covers the evapotranspiration models of PDM, NAM and VHM. It has two parameters: the power exponent be and the maximum soil storage umax. The more advanced VHM version (4)-(5) (but in VHM with fe = 1) has parameter uevap, being the soil storage above which the actual evapotranspiration equals the potential evapotranspiration. 1.1.2

Rainfall component (Rainfall contribution to surface storage)

The rainfall component transfers the available rainfall observations, at the different rain gauges and the spatial radar images, to rainfall input, which for a lumped conceptual model is the catchment averaged rainfall. This is done by an interpolation method (commonly Thiessen polygon method), combined with a merging technique when radar data are used. The interpolation method can be included in the proposed framework, but this can be done at a later stage. For now, we assume that the interpolation and merging are done externally and that the rainfall component can be omitted from the software. This means the catchment averaged rainfall is provided as rainfall input (x) to the model. The rainfall component may also include a time delay, as considered in PDM. 1.1.3

Net rainfall component

The net rainfall component transfers the rainfall input (output from the rainfall component) to the net rainfall (xn). This is the rainfall left after losses at the surface (due to the interception storage, wetting losses, etc.). In NAM, the net rainfall is modelled as the overflow of the surface storage reservoir: x n (t ) = s(t − 1) + x (t ) − e p (t ) − x IF (t ) − s m ax when xn(t) > 0

(6)

xn(t) = 0 otherwise

(7)

Note that ep is used in (6) instead of ea for the evapotranspiration loss because in NAM the losses from the surface reservoir are taken at the potential rate; so this means that the actual evapotranspiration is taken equal to the potential one. This is different for the soil reservoir, where the actual evapotranspiration is lower than the potential one. That xIF is included as a term in (6) may sound as a surprise, but in NAM interflow is generated as outflow from the surface reservoir. So, following our generalized structure philosophy and when linking this to the NAM specific model structure, the interflow contributes to the surface losses. In PDM, the surface losses are not explicitly modelled. The net rainfall is simply calculated as the rainfall input minus the actual evapotranspiration. PDM moreover allows to increase or decrease the rainfall input with a factor fx, which indirectly may represent the effect of the surface storage: x n = f x x − ea

when ea(t) ≥ x

x n = 0 otherwise

WL2021R00_130_8

(8) (9)

Final version

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

In VHM, the net rainfall equals to the rainfall input (so the surface losses are neglected; this means that the submodels describing the different runoff fractions intrinsically should account for these losses instead):

(10)

xn = x

1.2 Storage components 1.2.1

Surface storage component

The surface storage component (s) describes the storage of water at the surface due to interception storage, wetting losses, etc. In NAM, surface storage is modelled as follows: s(t ) = s(t − 1) + x (t ) − e p (t ) − x IF (t )

when s(t) ≥ 0 and ≤ smax, (11)

s(t) = 0 when s(t) < 0

(12)

s(t) = smax when s(t) > smax

(13)

where smax is the parameter describing the maximum surface storage capacity, xIF is the rainfall amount that contributes to interflow. In PDM and VHM, surface storage is not modelled, which means that this component can be omitted because its output variable s is not used in the other components. Or, the component can be replace by a dummy equation: s(t) = 0 (14) That the surface storage is not modeled should not be a complication for the model flexibility. In the generalized model structure, a surface storage can be considered for the PDM and VHM models, but with storage capacity equal to zero. 1.2.2

Soil storage component

The soil storage component describes the soil storage u. In PDM, the soil storage model is simply closing the water balance taking the infiltration rate (xu in the case of PDM) as inflow and the groundwater recharge (xSF) as outflow, and considering a maximum soil storage capacity as parameter umax: u(t ) = u(t − 1) + xU (t ) − xSF (t ) when u(t) ≥ 0 and ≤ umax

(15)

u(t) = 0 when u(t) < 0 based on (15)

(16)

u(t) = smax when u(t) > umax based on (15)

(17)

In NAM and VHM, the actual evapotranspiration is additionally considered, which means that (15) is replaced by: u(t ) = u(t − 1) + xU (t ) − ea (t ) + ec

when u(t) ≥ 0 and ≤ umax

(17)

The term ec in (17) is only considered if capillary rise is modelled; otherwise ec=0. In HBV, the soil storage is modelled based on (17), but with ec=0.

Final version

WL2021R00_130_8

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

1.3 Splitting components 1.3.1

Rainfall constribution to slow flow (groundwater recharge)

The groundwater recharge component describes the rainfall amount that contributes to the groundwater or slow flow routing. In PDM, the groundwater recharge depends on the soil storage by a power relation: x SF =

(

1 u − u tr ,SF u ma x kg

(18)

The parameter utr,SF determines the threshold (considered as the relative u, hence defined as a fraction of umax) below which no groundwater recharge will occur. When utr,SF = 0, the equation (18) corresponds to the classical equation of a reservoir model, a linear reservoir model when bg=1 and a non-linear reservoir model when bg<>1. The parameter kg is the reservoir constant or recession constant. In NAM, the groundwater recharge also depends on the soil storage, this time in a linear way, but with a different relation, that requires the calculation of the rainfall contribution to quick flow first: (19)

x SF = fSF ( x n − xQF )

where the fraction fSF depends linearly on the soil storage: fSF =

u − u tr ,SF u max u max − u tr ,SF u max

(20)

The parameter utr,SF again determines the threshold (fraction of umax) below which no groundwater recharge will occur. In VHM, the groundwater recharge is not explicitly modelled; but considered as the component that closes the water balance. Note that PDM and the most parsimonious VHM model only consider one SF component, which means that only one groundwater recharge component is considered. In NAM, two SF components can be considered. In that case, there are two rainfall contributions to SF: xSF1 and xSF2. 1.3.2

Rainfall contribution to quick flow

This component describes the rainfall amount that contributes to the quick flow routing. In PDM, this rainfall amount is a fraction of the rainfall that is left after actual evapotranspiration and groundwater recharge: xQF = fQF ( x n − ea − x SF )

(21)

For the fraction fQF, there are different options depending on the distribution selected to describe the spatial variability in the soil storage capacity. For the Pareto distribution, which is most commonly selected, the fraction fQF depends as follows on the soil storage: fQF

 u = 1 − 1 −  u max

 b +1  

(22)

Note that the PDM model requires the groundwater recharge to be calculated first, before the rainfall contribution to quick flow is obtained. In NAM, the equations go as follows:

xQF = fQF x n

WL2021R00_130_8

(23)

Final version

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

fQF = ϕ

u − u tr ,QF u max

(24)

u max − u tr ,QF u max

In VHM, the same equation (23) is applied, but the fraction fQF can take different forms, depending on the data. The common forms are the linear and exponential equations: fQF = aQF ,1 + aQF ,2

u u max

fQF = aQF ,1 exp( aQF ,2

(25)

u ) u max

(26)

PDM only has one QF component, but NAM has two: OF and IF. For the OF component in NAM, equations (23) and (24) are applied. For the IF component, xIF is modeled in a special way based on both the surface storage and the soil storage: x IF = fIF s fIF =

(27)

1 u − u tr ,IF u max k IF u max − u tr ,IF u max

(28)

where kIF is the interflow recession constant and utr,IF the threshold (fraction of umax) below which no interflow will occur. In VHM, the user has the flexibility to choose between one QF component or two. In the latter case, the OF component is defined by (23) in combination with (25) or (26), and the IF component in a similar way: (29)

x IF = fIF x n

Where the fraction fIF again can have different forms, depending on the data. The common forms are the linear and exponential equations: fIF = aIF ,1 + aIF ,2

(30)

u max

fIF = aIF ,1 exp( aIF ,2

u u max

)

(31)

For the standard HBV, only one QF component is considered, but other versions, such as the HBV-light, consider both OF and IF. 1.3.3

Rainfall contribution to soil storage

This component describes the rainfall amount that contributes to soil storage. Both in PDM and NAM, the rainfall contribution to soil storage close the water balance. In VHM, the rainfall contribution to soil storage is explicitly modelled: (32)

xU = fU x n

Where the fraction fU can take different forms, depending on the data. The common forms are the linear and exponential equations: fU = aU 1 + aU 2

u u m ax

fU = aU 1 − exp( aU 2 (

Final version

(33)

u aU ,3 ) ) u m ax

WL2021R00_130_8

(34)

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

1.4 Component to close water balance As shown above, each model has one specific component for which no specific model equation is applied but the calculation is based on the water balance closing. This is for the rainfall contribution to soil storage in PDM and NAM, and for the groundwater recharge component in VHM. The water balance closing is done as follows: x n = xQF + x SF + xU

(35)

Or if there are two QF components: x n = xOF + x IF + x SF + xU

(36)

x n = xQF + x SF 1 + x SF 2 + xU

(37)

Or if there are two SF components:

1.5 Flow routing components The flow routing is in conceptual models most commonly done by means of the linear reservoir model: 1 1 c(t - 1) + c(t ) y (t ) = exp(- )y (t - 1) + (1 - exp(- ))( ) k k 2

(38)

where y(t) is the outflow at different time steps t, based on the inflow c (rainfall contribution to the subflow) and the recession constant k. An alternative solution for the linear reservoir is by means of the following two equations, hence explicitly considering the reservoir storage uY: y=

1 uY k

uY (t ) = uY (t − 1) + (

c (t − 1) + c (t ) ) − y (t ) − e c 2

(39) (40)

The term ec in (40) is only considered if capillary rise is modelled; otherwise ec=0. In case of a non-linear reservoir model, (39) is replaced by: 1 (uY ) m k

(41)

with exponent m. The linear reservoir model obviously is a special case of this for m=1. Some models use for some subflows two or more reservoirs in series. In that case, the reservoir equation is repeated two or more times. For instance, for two linear reservoirs in series this gives: 1 1 c (t - 1) + c (t ) )y 1(t - 1) + (1 - exp(- ))( ) 2 k1 k1

(42)

y (t - 1) + y 1(t ) 1 1 )y (t - 1) + (1 - exp(- ))( 1 ) 2 k2 k2

(43)

y 1(t ) = exp(y (t ) = exp(-

In VHM and HBV, the linear reservoir model is applied for each routing component (for each of the three subflows): OF, IF and SF. In NAM, also SF routing is based on a single linear reservoir model, but OF routing and IF routing are modelled by the same two linear reservoirs in series. Also in PDM, the QF routing goes with two linear reservoirs in series. Each reservoir may have a different recession constant. The recession constant are most often constants, except for one of the two OF or IF routing reservoirs in NAM, where the recession constant is depending on the rainfall amount contributing the overland flow as follows:

WL2021R00_130_8

Final version

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

x k =  OF  x*  OF

   

−0.33

(44)

* where k* is a parameter describing the recession constant when xOF equals a specific value xOF , which is taken equal to 0.4 mm/h in NAM.

1.6 Special component: Capillary rise component In some cases, some special additional components need to be added. One example is given below, for the component describing capillary rise in NAM. Capillary rise is modelled in NAM as follows, depending on the relative soil storage and the groundwater table height which (through the SF routing model) depends linearly on the SF outflow: y ec =  SF  ae1

  

ae 2

 u 1 − u max 

  

0 .5

(45)

To the parameters ae1 and ae2, a physical interpretation can be given that is related to the groundwater table height corresponding to 1 mm/h capillary rise.

1.7 Sequence of solving equations As explained above, the sequence in which the equations of the different components are solved depends on specific model. For example, the NAM model first solves the rainfall contribution to SF before computing the soil storage (because the soil storage is modelled as rest term, to close to overall water balance), whereas in VHM the opposite sequence is followed (the rainfall contribution to SF is considered as rest term). Question thus is whether first the storage models have to be solved, followed by the splitting components, or vice versa. Or, maybe first the QF related components have to be solved (splitting and storage components), followed by the SF related components. When the time step is taken very small or the equations solved by iteration per time step (till they converge to a solution), the sequence does not matter. However, in order to reach limited computational times, it would be best to avoid small time steps or iterations. Therefore, the following approach is proposed here (which might affect the results, but given that the model is calibrated at the end, this choice might not strongly affect the quality of the final model simulation results): First, solve the meteorological components (is very trivial; is done as such in all models). Second, solve the surface storage component. Third, solve the components that have the strongest changes in time. This means the QF splitting components. And in case of two QF components: OF first, IF next. Fourth, the SF splitting components. Fifth, the soil storage component. This component is taken last because it is assumed that the temporal variability of the soil storage is low. This means that the QF and SF splitting calculations can be based on the soil storage result at the previous time step. The soil storage results will indeed not vary much between the previous and the current time step. During the study, a sensitivity analysis can be conducted to investigate the influence of the model results on the selected sequence (and we may give flexibility to the tool; and allow the user to change the sequence). Final version

WL2021R00_130_8

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

2 Spatial linking 2.1 General considerations Two general routing approaches, i.e. source-to-sink and grid-to-grid routing, for implementing the spatial linking in the distributed hydrological modelling were briefly discussed in the previous report (Tran et al., 2013). The source-to-sink approach is simple and commonly used for lumped or semi-distributed models to route the flow indirectly from a certain source of the catchment to the outlet. It can be implemented simply by taking the runoff accumulation produced in the whole catchment or in every single cell or hydrological response unit and transferring to the outlet. The grid-to-grid method, however, is much more complicated and can be implemented in different ways as illustrated by Tran et al., (2013). In this method, the explicit spatial linking between grid cells, which is more physically based, is calculated using local drainage direction map. The main advantage of the explicit spatial linking scheme is that the hydrological connectivity is considered (the way grid cells are inter-linked or linked to the river system). The correct representation of the hydrological connectivity may be important when the model is used for studying the impact of spatial (e.g. land use) scenarios. Depending on the spatial pattern of land use changes (e.g. increase in pavement, which most often is spatially grouped), a change in land use of a given cell may produce additional surface runoff. But this extra surface runoff may not reach the river if there is no hydrological connectivity between that cell and the river network or if it gets infiltrated in the soil of neighboring cells. The same percentage of increment in pavement but with different spatial patterns or grouping may lead to same river impacts in the source-to-sink approach, but highly different impacts in the gridto-grid approach. However, when a coarse spatial resolution is considered for the distributed model, the role of the hydrological connectivity reduces. From this, it becomes obvious that the source-to-sink approach is more useful for coarser grid sizes. In the ultimate case of the lumped model, the approach becomes source-to-sink in any case. Nevertheless, the discussion on how to implement the source-to-sink approach will be skipped in this report due to its simplicity. The grid-to-grid approach, which is more complex, will be further developed. It can be constructed in the generalized structure in two ways: the normal routing and the detailed routing. Depend on the base model code (NAM, PDM or VHM) that is used in the framework, the detailed routing approach requires some modifications to represent the interaction between sub flows. The description of these approaches can be found below.

2.2 The normal routing approach The implementation of the normal routing approach is schematized in Figure 2-1. For each grid cell the generalized hydrological model structure of Figure 1-1 (so called rainfall runoff generator) will be implemented. Every grid cell does not only receive the local rainfall contributions but also the output from the different upward cells. In this approach, the quick runoff flow from a given upward cell will remain quick runoff flow also in the downward cell. This means that the amount of rain water contributing to QF in the downward cell will be based on the rainfall contributing to QF in the local cell itself plus the QFs received from the different upwards cells (based on the LDD map). The same applies to the other subflows (SF in Figure 2-1). For reasons of simplicity, Figure 2-1 is only based on two subflows (QF and SF), but can be extended to more subflows using the same principle.

WL2021R00_130_8

Final version

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Figure 2-1: Linking between two grid cells following the generalized hydrological model structure, and the normal routing approach

2.3 The detailed routing approach In the detailed routing approach, water originating from a given subflow (e.g. QF or SF) can move to other subflows or to soil storage depending on the properties of the neighbouring cells. Figure 2-2 proposes a change to the linking by the normal routing approach in Figure 2-1, and also a slight change/extension to the generalized model structure by introducing components where water is transferred/added (blue nodes). In the upper node, the net rainfall from the downward cell is added to QF from the different upward cells. The output from this node is the total amount of water available for runoff in the downward cell. How this water is divided over the different subflows and soil storage(s) is controlled by the different process equations, in the same way as proposed before for individual or lumped cells. For the non-OF subflows (SF in Figure 2-2), lower nodes are introduced (see lower blue node in Figure 2-2) to sum up the non-OF subflows from the different upward cells and add these to the total rainfall amount available for sub-surface or groundwater processes (storage or runoff). Slight modification of the net rainfall distribution process is then required to prevent a buildup of water in the SF subflow. The method can be extended even to include the effect of topography on the distribution of the different subflows, with a gradient-based alteration to the distribution of subflows from upward cells. In case of full soil saturation and depending on topographical differences between the upward and downward cells, the non-OF subflows may lead to QF more downstream. This is the “return flow” node link component shown in Figure 2-2.

Final version

WL2021R00_130_8

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Figure 2-2: Linking between two grid cells following the generalized hydrological model structure, and the first attempt of detailed routing approach implementation

However, when QF from upstream cells is added into the total water budget in the upper node, the amount of water infiltrating into the soil storage is changed. As a result, the soil moisture content, which is the core of the generalized structure, is not kept as it was after careful calibration of the lumped model. This consequently strongly affects the other subflow volumes, which are calculated based on the relative soil moisture content. To avoid additional calibration for the spatial distributed model, which is rather time consuming and complicated, it is necessary to modify the routing scheme to maintain the soil moisture storage as in the calibrated lumped model. The adjusted scheme is illustrated as in Figure 2-3. In this scheme, instead of introducing the new distributing node before the runoff generation, another nodes are established right after the routing processes.

Figure 2-3: Linking between two grid cells using the second attempt of detailed routing approach implementation

A part of the OF subflow generated in the upward cell is routed to participate in the infiltration process while the remaining keeps routing to OF in the downward cell. The infiltrated part together with the generated IF flow is again divided into two parts, one percolates to the groundwater reservoir and one is added to IF of the next downstream cell. The division ratios can be determined using the ratios between the rainfall contributions to each subflow. However, BF keeps transferring directly into the 14

WL2021R00_130_8

Final version

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

same component as in the normal routing approach. This leads to an extra water is coming to the BF volume, which was well calibrated in the lumped model, from upper layers. If there is no water that is taken away to surface or sub-surface zone, the BF volume will be overestimated. Hence, the OF and IF volumes will be underestimated. This problem can be solved if the “return flow” term, which can occur in reality due to slope or capillary rise, is introduced. However, it requires additional calibration, which might be complicated and time consuming. Another solution is to modify the structure of the runoff generator by applying two new parameters (FactorOF and FactorIF) for increasing the OF and IF volumes. By doing that, the volume of BF is reduced to compensate for the added amount. The FactorOF and FactorIF parameters can be easily calibrated by comparing the simulated and observed subflow volumes in a step wise approach, which is described in the next chapter.

2.4 Modifications for detailed routing approach In the detailed routing approach, the splitting of the subflows is done subsequently after receiving income from the upstream cell. To make this process more physically sound, a slight modification has been applied. In the modified scheme, the splitting processes for both OF and IF occur independently at the same time. The updated routing scheme is shown in Figure 2-4.

Figure 2-4: Spatial linking between two grid cells using the modified detailed routing approach

Nevertheless, the proposed spatial linking scheme still requires some specific adaptations when applying it for a specific hydrological model, due to the difference in structures. For example, there is no IF in this version of PDM or the IF term does not interact with other subflow components in NAM. Figure 2-5 and Figure 2-6 describe how the spatial linking scheme has been adapted for the NAM and PDM models. For the NAM model, the generated IF is simply the horizontal leakage from the surface storage and it has no interactions with other components. Hence, it does not take part in any equations related to OF and BF. It also participates in neither the infiltration, nor the percolation processes. In opposite, OF and BF originate from the same source (the excess water), and therefore they have significant interactions. The scheme in Figure 2-5 describes this relationship by allowing OF to infiltrate directly into the BF reservoir. In this case, only one additional parameter is needed (FactorOF) while the other (FactorIF) is taken equal to 1.0 (no effect on the generated IF).

Final version

WL2021R00_130_8

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Figure 2-5: Spatial linking between two grid cells using the modified detailed routing approach with adaptation for the NAM model

For the PDM model, the IF variable is omitted in the working version. Water in the sub-surface and surface zones (QF or OF in this case) is allowed to move to and interact directly with the groundwater component. The parameter FactorOF is also required to increase the OF volume in order to compensate the reduction after routing (see Figure 2-6).

Figure 2-6: Spatial linking between two grid cells using the modified detailed routing approach with adaptation for the PDM model

WL2021R00_130_8

Final version

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

3 Model calibration method The proposed model calibration method was briefly described in the previous report of Tran et al. (2021). In this inventory, the method is summarized and then further discussed for both lumped and distributed modelling.

3.1 Calibration concept for lumped models For the lumped model, a top-down methodology based on the VHM calibration approach (Willems, 2000, 2013 and 2014) is presented for the implementation and calibration of the model based on river flow data at a gauging station. This methodology aims to have better model evaluation, modelstructure selection; to avoid overparameterization; and to reduce the time needed for calibration. The methodology is an empirical and step-wise technique that includes a step by step examination of the various model components through a data-based analysis of response characteristics. Subresponses are separated and submodel components and related subsets of parameters are calibrated as independently as possible. At the same time, the model-structure identification process aims to reach parsimonious submodel-structures. Parsimony means that the best approach is the simplest that fits the purpose of the application (Harremoës, 2003). The flexibility in defining the detailed model structure is a virtue that fits with parsimony. More details about this VHM based model structure identification and calibration approach can be found in (Willems, 2000, 2013 and 2014). As explained above, this approach works by identifying and calibrating the different hydrological model components based on multiple and non-commensurable information derived from river flow series by means of a number of sequential time series processing tasks. These include separation of the river flow series in subflows, split of the series in approximate independent quick and slow flow hydrograph periods, and the extraction of independent peak and low flows. Prior time series processing tasks Step 1a: Reverse river routing

river flow series ⇒ rainfall-runoff series

Step 1b: Subflow separation

⇒ recession constants ⇒ series of overland flow, interflow and slow flow

Step 1c: Split in quick and slow flow events ⇒ event-based runoff volumes

Step-wise model-structure identification and calibration tasks Step 2: Identification and calibration routing models

⇒ event-based rainfall fractions to subflows

Step 3: Identification and calibration soil moisture storage model ⇒ closing the water balance

Step 4: Identification and calibration rainfall fractions to quick and slow flow events

Figure 3-1: Steps in the lumped model structure identification and calibration procedure

Final version

WL2021R00_130_8

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

3.1.1

Step 1: Prior time series processing

The VHM based approach to identify the lumped conceptual model structure starts from a series of river flow observations available in the catchment. In Step 1, a number of prior time series processing tasks are carried out: Transformation of the river flow series into a series of lumped (upstream catchment averaged) runoff as produced by rainfall, i.e. elimination of the river hydraulic aspects; Separation of the rainfall-runoff series in its subflow components, overland flow, interflow and baseflow using a numerical digital filter; Split of the rainfall-runoff series in approximated independent quick and slow flow events. These three prior time series processing steps are explained hereafter. Transformation of river flow series into rainfall-runoff series or elimination of the river hydraulic effects Catchment rainfall-runoff, i.e. the runoff as produced by rainfall, enters the catchment’s river network in spatially distributed manner. The flow observed at the river flow gauging station is the cumulated result of the runoff from the upstream catchment, and is moreover affected by flow attenuation along the river network. Due to spatial variability of the runoff, river flow routing, river flooding (floodplain storage) and hydraulic regulation, the flow observed at the river gauging station might differ from the rainfall-runoff of the upstream catchment. This effect is often neglected in current rainfall-runoff modelling practices, or included in the rainfall-runoff model (e.g. by including a river routing component or a time delay; although such component often will not address the effect of river flooding or hydraulic regulation). This project focuses on the modelling of rainfall-runoff which can also be used as input for river hydrodynamic models, thus with an approach which has the flexibility to exclude the river hydraulic effects. For the latter purpose, the observed river flow series need to be transformed in a series of “equivalent” upstream rainfall-runoff discharges (Willems, 2000, 2013, 2014). Given that these runoff discharges are spatially variable and distributed over the length of the river(s) and given the lumped nature of the rainfall-runoff model considered, a series representing the spatially averaged runoff over the river network upstream of the gauging station is required. Such runoff series can be derived from the observed river flow series by reverse modelling of the river flow hydrodynamics along the river network (reverse routing). When an accurate hydrodynamic model is available for the river under study, the river response behavior (response of river to upstream runoff) at the gauging station can be studied and parameterized. One option is to calibrate a conceptual routing model (e.g. linear reservoir model or cascade of linear or non-linear reservoirs) to the full hydrodynamic river model, which enables the runoff – river flow response to be reversed. In (Willems, 2013, 2014), the following conceptual model is proposed for describing the influence of river flooding, when the river flow is higher than the ‘flood threshold’ qF:

where:

q= y

y < qF

q = y - f ( y - z -1 q) + q d

y ≥ qF

f = 1 - exp(-( y - q F ))

y≥q

f = 1 - exp(-α hFP )

y<q

and where : q

: river flow series at the flow gauging station;

WL2021R00_130_8

Final version

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

: upstream averaged rainfall-runoff series (upstream of the flow gauging station);

hFP

: floodplain water depth in meter;

: river discharge capacity increase per time step (only for y < q; qd = 0 for y ≥ q);

α, q

: conceptual flood model parameters;

z-1

: backward shift operator ( z-1q(t) = q(t-1) ).

The conceptual model is based on the simplified concept that above the flood threshold, a fraction f of the difference between the rainfall-runoff and river discharge (assessed based on the river discharge during the previous time step) is stored in the floodplain. The remaining fraction contributes to the river flow q. In the rising limb of the flood hydrograph, the fraction f increases with the rainfall-runoff discharge in excess of the flood threshold. In the recession limb of the flood hydrograph, or more specifically after the rainfall-runoff discharge decreases below the river discharge, the stored floodplain volume is released to the river depending on the floodplain water depth. This release will go up to a maximum of the river discharge capacity, which is the current river discharge (based on previous time step) increased with a value qd per time step in periods of inundation. The flow difference between runoff and river discharge (y-q) is used to calculate flood volumes, hence the storage of water in the floodplain (vFP). From this storage, the water depth along the floodplain (hFP) can be assessed using a hypsographic curve, derived from digital terrain elevation data along the floodplain. It is clear that this river conceptual model structure might be strongly case dependent, and can be derived only when an accurate hydrodynamic river flow model is available. If no important river hydraulic effects are expected, this step of the methodology can be skipped. The same applies if no detailed river hydrodynamic model is available in a given case study, but the modeller has to take the limitation into account. Subflow separation The split of the observed flow series in quick and slow flow events can be done by means of subflow separation techniques (Chapman, 1999; Arnold & Allen, 1999; Eckhardt, 2005). A generalization of the recursive digital filter proposed by Chapman (1991) is recommended here, because of the direct linking between the parameters of this filter and the lumped hydrological characteristics of the catchment under study. The Chapman filter is derived from the general equation of a ‘low pass filter’, as applied in for instance signal processing, and makes additional use of the assumption of exponential recessions for the hydrological subflows, with nearly constant recession constants for different low flow periods. The filter aims to split the total flow time series q(t) in the subflow or slow flow component series b(t) and the quick flow series f(t). The original filter by Chapman (1991) has one parameter: the recession constant k of the subflow to be separated. It is shown by Willems (2000, 2009, 2013, 2014) that in the original form of that filter the intrinsic assumption is made that the total long-term volumes of the slow and quick flow series are identical (each 50% of the total runoff). The slow and quick runoff fractions, however, strongly vary between catchments depending on their (topographical, soil type…) characteristics. Therefore, a generalization of the original Chapman-filter was proposed by Willems (2000, 2009, 2013, 2014), where a new filter parameter w is introduced that represents the casespecific average fraction of the quick flow volumes over the total flow volumes. After this generalization, the filter equations become:

f (t ) = a1 f (t - 1) + a2 (q (t ) - αq (t - 1))

b(t ) = q (t ) - f (t ) = α b(t - 1) + a3 (1 - α ) ( f (t - 1) + f (t ))

Final version

WL2021R00_130_8

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

where:

a1 =

( 2 + v )α - v 2 + v - vα

a2 =

2 2 + v - vα

a3 = 0.5ν α = exp(−1 / k )

1− w w

The recession constant k equals the time in which the flow is reduced during dry weather flow periods to a fraction exp(-1) = 0.37 of its original discharge. Indeed, for (dry weather) periods with no quick flow (f(t)=0), the assumption is made that the subflow decreases in an exponential way, as is observed worldwide for many rivers (Nathan and McMahon, 1990):

1 b(t ) = exp(− )b(t − 1) = α b(t − 1) k When the outflow discharge at time t=0 is denoted by b(0), the outflow discharge at time t=k is a fraction exp(-1) of b(0):

b(k ) = exp(−1)b(0) The recession constant or recession time k of each subflow can be quantified as the average value of the inverse of the slope of the linear path in the subflow recession periods of a ln(q) - time graph. Willems (2000, 2099, 2013, 2014) explained how this can be done by visual inspection of this slope for a number of recession periods. The second parameter w can be calibrated by optimizing the height of the subflow during the recession periods (the subflow needs to be identical to the total flow during the subflow recession periods, and needs to be lower than the total flow during the other periods). Also this calibration can be done by trial-and-error through visual inspection in the time series (Willems, 2000, 2009, 2013, 2014). Given the strong difference in the order of magnitude of the recession constants for the three classes of subflows, the subflow separation can be carried out in a step-wise way. In a first step, the slow flow component is split from the total flow, and in a second step the interflow split from the remaining flow (total flow minus filtered slow flow). The rest fraction then represents the quickest flow component. Figure 3-2 shows the baseflow and interflow filter results based on the Dijle river flow results at Sint-Joris-Weert. The calibrated filter parameters (k and w) equal kSF = 3700 hours for the baseflow or slow flow (SF), kIF = 30 hours for the interflow (IF), kQF = 6 hours for the quickest flow (QF), wSF = 0.23 for the sum of the quickest flow and interflow fraction in the total flow (when obtaining the SF component in a first filter step), and wIF = 0.5 for the mean quickest flow fraction in the sum of the quickest flow and interflow runoff (when separating the IF component in the second filter step). Of course, the visual inspection based calibration of the k and w parameters involves some subjectivity and uncertainty in the estimation.

WL2021R00_130_8

Final version

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Figure 3-2: Subflow filter results for the hourly river flows at Sint-Joris-Weert (note that “filtered interflow” in the plot means interflow + baseflow)

Event separation In a given series of river flow discharges, independent peak flow discharges can be defined (or filtered from the flow series) using independence criteria. The dependence among subsequent peak flows depends on the autocorrelation in the flow series, which is determined primarily by the catchment recession. Subsequent rainfall-runoff peaks can be considered largely independent if the inter-event time exceeds the recession time and if the lowest inter-event discharge drops below a specific low flow level (see also USWRC (1976) and Lang et al. (1999)). Another approach, which is commonly used in statistical extreme value analysis (e.g. Madsen et al., 1997), is the Peak-Over-Threshold (POT) approach where peak flows (or peak flow extremes) are selected as the highest discharges in periods where the flow up-crosses a selected threshold. The threshold selection can be based either on physical or statistical considerations, as discussed extensively in Lang et al. (1999). The threshold is usually compared with the total discharge value of the peaks. It can also be applied to the discharge increment relative to the baseflow level preceding each peak event, as was done by Claps and Laio (2002). Another method, presented by Willems (2000, 2009, 2013, 2014), selects peak flows from the flow series based on criteria for the inter-event time, the inter-event low flow discharge and the peak height (see Figure 3-3). Two subsequent peak events are considered approximately independent when the following three conditions are fulfilled: (i)

the time length τ of the decreasing flank of the first event exceeds a time kp:

τ > kp (ii) flow:

the discharge drops down - in between the two events - to a fraction lower then f of the peak

qmin <f qmax

or close to the baseflow qbase:

(iii)

q min − qbase < f q max

the discharge increment qmax – qmin has a minimum height qlim:

q max − q min > qlim Final version

WL2021R00_130_8

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

qmax

τp

total flow series baseflow series

qmin qbase

Figure 3-3: Parameters used in the criteria to select approximate independent peak flows from a river flow time series.

This procedure for peak flow selection has three parameters: kp, f and qlim. It is based on the concept that a peak flow event can be considered largely independent from the next one, when the inter-event discharge drops down to a low flow condition or almost to the baseflow level (see criterion (ii)). Under this condition, the quick flow components attributed to the peak flow events are indeed approximately independent. It is clear that a near-baseflow condition can only be reached when the time length τ is longer than the recession constants of the quick flow components. The parameter kp thus can be taken equal to the recession constant of the quick flow, or higher (e.g. two or three times the recession constant) depending on how strong one wants to have the independence between subsequent quick flow events. When the dry inter-event period is longer than the recession constant, the quick flow indeed decreases down to a fraction lower than 1/e of the peak flow. Because inter-event periods are not always completely dry, also criterion (ii) needs to be checked. When the fraction f is based on

qmin < f , it has to be taken approximately as the upper limit of the baseflow fraction in the peak flow. qmax

When applying q min − qbase < f , small fraction values of 5%, 10% or 15% are considered, depending q max

on the accuracy of the baseflow filter results. Criterion (iii) is finally needed to avoid that small noise peaks are selected. Parameter qlim consequently is taken as an estimate of the upper limit of the highest noise peaks (high frequency perturbations which cannot be associated with quick runoff events). The latter is visually judged in the time series. Based on these three criteria, all significant peak events are selected from the flow series, and large physical and statistical independence reached between consecutive peak events. After this selection of peak flows qmax, the lowest flow values in between two consecutive peaks (qmin) can be defined and the time series divided in periods based on the time moments of these low flows. These periods can be considered as partial-duration series: PDS series. They can be considered as approximately independent quick flow hydrograph periods. Using a similar procedure, but using the recession constant for baseflow (or larger) for the parameter kp in the peak flow selection procedure, approximately independent slow flow periods (nearly independent baseflow events or wet periods) can be defined. More details on the method are given in Willems (2000, 2013, 2014).

WL2021R00_130_8

Final version

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

3.1.2

Step 2: Step-wise model structure identification and calibration

Using the results from the time series pre-processing, characteristics of individual rainfall - runoff response processes are derived from the series using a new data-based method. Lumped conceptual representations of individual rainfall-runoff process equations are identified in a separate and more or less independent way, and calibrated to related subsets of model parameters. This includes sequentially: -

identification and calibration of the routing submodels (Step 2 in Figure 3-1);

identification and calibration of the submodels describing the soil moisture storage: the splitting component describing the rainfall contribution to soil storage, the meteorological component transforming potential evapotranspiration in actual evapotranspiration (Step 3 in Figure 3-1);

identification and calibration of the splitting describing the rainfall contribution to the different subflows (quick and slow flow components) (Step 4 in Figure 3-1).

The identification and calibration of the routing models in Step 2 is done based on investigation of the hydrographs’ shape and the calibration of recession constants describing the catchment response time for specific subflows. This step comes first in the model-structure identification procedure because the identification of the hydrographs’ shape can be done independently of the hydrographs’ volumes or peak flow magnitudes, and because the catchment response times can be derived empirically and separately for the different subflows. In Step 3, the overall water balance is investigated and the main lumped processes affecting the soil moisture storage are identified and calibrated (e.g. lumped representation of the infiltration and percolation processes). Hereafter in Step 4 the results on the soil moisture storage are used to investigate the distribution of the total runoff over the different runoff subflows and to identify and calibrate the subflow separation processes. For fixed model structures, such as NAM and PDM, the basic principles of the approach still can be applied, as explained in (Willems et al., 2002) and (Willems et al. 2014a). Applications of the method can be found in (Willems, 2014) and (Willems et al., 2014a) for the Molenbeek at Erpe-Mere and the Grote Nete at Varendonk / Geel-Zammel, and in (Vansteenkiste et al., 2012) for the Grote Nete.

3.2 Calibration concept for distributed models Question is how to extend the VHM based calibration approach to the spatial model versions. As specified before, the grid size of the model can be chosen depending on the application. Figure 3-4 illustrates differences in spatial resolution for a catchment depending on the application. The user can change the spatial resolution depending on the needs and/or the evaluation of the accuracy of the model results, or use different spatial resolutions in parallel for different applications. However, one of the concerns here is to obtain consistent results at any spatial detail. It means there should not be a strong and sudden change in model parameter values and corresponding simulation results, when one moves from one level of spatial detail to another. Moreover, to ensure that the lumped models are able to produce results that are as accurate as distributed models, one may consider the option to check whether the structure of the lumped model can be adjusted in order to improve the results in light of the consistency aimed with our disaggregation/aggregation approach between the lumped and distributed models. These aspects need to be taken into account in the model calibration approach.

Final version

WL2021R00_130_8

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Figure 3-4: Different spatial modelling scales

Note that the requirement of consistency between spatial resolutions is a constraint that helps to reduce the over-parameterization problem. Additional model performance statistics based on comparison of model results with observations, e.g. flow at internal flow gauging stations, groundwater well levels, soil moisture products, etc., and would also be very useful in order to further decrease the equifinality. In order to reach such consistency, it would be useful if the calibrated parameters for the model at one spatial resolution can be used for the models at other spatial resolutions. 3.2.1

Step-wise disaggregation concept

Considering the principles of the top-down calibration approach, the parameters calibrated at the lumped scale are disaggregated to distributed parameters. This can be done in two steps. In the first step, a first disaggregation is conducted based on spatial calibration data such as the flow at internal flow gauging stations and the groundwater heads at wells. If later remote sensing would provide reliable soil moisture data, these can be considered as well. The second disaggregation step makes use of available spatial catchment characteristics such as land use, soil type and geology. Thanks to the disaggregation, one can make sure that the model implementation (and consequently also the simulation results) at the different spatial resolutions are consistent. Or, at each spatial resolution, model calibration can be done by optimization, conditioned on the model calibration and simulation results at the coarser spatial resolutions already conducted. It is clear that the second disaggregation steps will involve assumptions. Because many different assumptions can be made, at least a sensitivity study needs to be done to investigate the impact of these assumptions to the uncertainty in the model calibration and simulation results. Several plausible assumptions will exist, such that the second disaggregation step leads to equifinality. If the uncertainty in this “unknown” spatial disaggregation induced by the second disaggregation step can be assessed, this would allow the additional uncertainty in the spatially variable model simulation results due to the spatial disaggregation uncertainty to be quantified. For the model simulation results of the total runoff discharges, it is expected that these simulation results have the highest accuracy at the flow gauging stations. At other locations along the river, the uncertainty will be higher, and higher for locations at larger distance from the gauging stations. It would be useful to quantify this spatial variation in the runoff discharges. 3.2.2

Spatial disaggregation of model parameters based on spatial catchment characteristics

The spatial disaggregation of the calibrated model parameters at the lumped scale into higher spatial resolution (distributed or semi-distributed), can be done by using the basic spatial catchment characteristics on topography, land-use and soil type. While disaggregating the model parameter values, these values should not exceed certain limits, which represent the physical ranges of these parameters. These ranges can be estimated from previous modelling studies in the region (hence from empirical experiences) and from literature review. Several existing approaches may be used for giving 24

WL2021R00_130_8

Final version

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

parameter values spatial variations. Such approaches exist in relation to specific distributed and semidistributed hydrological models such as WetSpa (Liu & De Smedt, 2004), MIKE-SHE (DHI, 2012) and HEC-HMS (HEC, 2000). Also for the SWAT model, a common approach exists to link parameter values to hydrological response units (HRUs) based on their specific land uses and soil types. Given that in the framework proposed here, the distributed model parameters will be derived by disaggregation, starting from the conceptual lumped parameters, it has to be taken into account that there are no direct relationships existing to particular catchment physical properties. For lumped models, the parameters indeed characterize one particular watershed. They describe the integrated behavior of the various properties that cover the entire watershed. The disaggregation process aims to decompose this integrated description; and in this way, to get a more direct linking (will be more direct for the finer model resolutions) with the spatially variable catchment properties. The next sections first will discuss how parameters characterizing physical hydrological processes can be related to specific spatial catchment properties. Examples are provided from the existing distributed modelling systems WetSpa, MIKE-SHE and HEC-HMS. Afterwards, this knowledge is applied to characterize catchment characteristics by means of conceptual model parameters. Parameters of physical processes characterizing soil texture classes Soil-type maps are used to provide information on soil physical properties such as porosity, hydraulic conductivity, etc. These soil hydraulic properties can be estimated by using the information from soil texture. Or, in other words, the model parameters related to soil hydraulic properties can be represented as a function of soil texture classes. This estimation approach has been applied widely in many hydrological models using different soil texture classifications. For instance, the WetSpa model employs the soil texture classification of the US Department of Agriculture, which consists of 12 classes depending on the percentage of sand, silt and clay in the soil sample (Liu & De Smedt, 2004). Based on this classification, a look up table analyzed by Rawls et al. (1982) and Cosby et al. (1984) was established as shown in Table 1. A similar table (Rawls et al., 1982) is applied by the HEC-HMS hydrological model (HEC, 2000). Moreover, it includes an extra property describing the wetting front suction, which is calculated from the pore size distribution of soil texture. In MIKE-SHE, the model parameters determined by the soil texture are the retention curve, the hydraulic conductivity, the bulk density and some other basic properties. Also in its evapotranspiration component, soil texture is taken into account, because evaporation from soil is controlled by the soil wetness, soil hydraulic properties, whereas the transpiration is mainly linked to the plant physiology, but which again depends strongly on the water availability in the unsaturated zone (Grahams and Butts, 2005). MIKE-SHE enables users to freely define or modify the soil texture characteristics so that Table 1 or any other table also can be employed. However, instead of using lookup tables to estimate model parameters, MIKE SHE calculates them directly from the soil properties using equations. Additionally, it can be taken into account that the vertical soil profile may consist of several layers.

Final version

WL2021R00_130_8

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling Table 1: Soil hydraulic parameters characterized by soil textual classes (Liu and De Smedt, 2004)

The classification system given by Rawls et al. (1982) is quite similar to that commonly applied in Belgium, in which each soil type is specified by a notation composited by three letters describing the soil texture, natural drainage and vertical profile (“horizontenopeenvolging”) (Van Ranst and Sys, 2000). Nevertheless, only the soil texture (of the upper part of the soil), characterized by the first letter is of most concern. Based on the different fractions of the three main textures including sand, clay and silt, the soil classification consists of seven classes: Heavy clay (U), Clay (E), Silt (A), Sandy silt (L), Light sandy silt (P), Silty Sand (S) and Sand (Z). A triangle chart of soil texture, introduced by the Centre of Soil Mapping to identify the texture classes is commonly applied (Figure 3-5).

Figure 3-5: Soil texture triangle chart (Van Ranst and Sys, 2000)

By using this triangle chart, the classified soil types in Table 1, including the hydraulic characteristics, can be converted into the Belgian classification system. In Figure 3-6 the comparable triangle chart

WL2021R00_130_8

Final version

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

from USDA textural classes have been put next to the Belgian one with the same axis in order to relate the different systems to each other.

Figure 3-6: Textural classes from Belgian system (A) and USDA system (B) (Dondeyne et al., 2012).

Parameters of physical processes characterizing land use classes Similar to soil type maps, land use or land cover maps are also very important to many hydrological processes such as the evapotranspiration, the surface runoff formation, the infiltration, etc. Each land use type has a certain surface roughness, which directly controls the overland flow velocity and infiltration process. Different land use types result in different evapotranspiration rates due to their vegetation cover characteristics. Because of its importance, consideration of detailed land use information is preferred. However, various hydrological models apply different classification systems and approaches to estimate the effect of land use on simulated flows. The WetSpa model defines the vegetation types based on the classification by the IGBP (International Geosphere-Biosphere Program), which includes 17 classes. This IGBP system can be converted into other equivalent ones such as Corine or Anderson classification systems (Herold et al., 2009; Tomaselli et al., 2013). HEC-HMS employs a more complex land cover classification. When the land use map is used as obtained from the National Geographical Institute (NGI 2007) (1:10.000 vector-layers; spatial accuracy of 1 meter; based on aerial photographs from 1998 (1:21.000)), it contains 47 different land use classes. Table 2 shows the conversion from these 47 NGI classes to the 9 main IGBP classes: Evergreen Needle leaf Forest, Broad-leaved woodland, Mixed Forest, Open Scrublands, Grasslands, Permanent Wetlands, Croplands, Impervious area and Water Bodies. For each land cover type, several parameters are estimated by taking previous studies as reference. Table 3 provides an example based on the WetSpa model. Table 2: Overview of the 47 land use classes in the NGI map and reclassification to the 9 main IGBP categories (Vrebos et al., 2013).

Description Coniferous trees Orchard Tree nursery Deciduous trees Poplar plantation

Final version

IGBP vegetation Evergreen Needle leaf Forest Evergreen Needle leaf Forest Evergreen Needle leaf Forest Broad-leaved woodland Broad-leaved woodland

WL2021R00_130_8

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Mixed deciduous and coniferous trees without dominant Mixed deciduous and coniferous trees with dominance of deciduous trees Mixed deciduous and coniferous trees with dominance of coniferous trees Sand Bare ground Coppice Heath Heath with deciduous trees Heath with coniferous trees Scrubs Brushwood brushwood with scrubs Cemetery Beds Pasture Gardens Deep swamp Reedland Cropland Transformer station Railway Road Crossroad Industrial building (in use) Industrial building (abandoned) Warehouse Silo Greenhouse Cooling tower Non-university hospital Town hall Schoolhouse Firehouse Commercial building Religious building Sports hall Covered grandstand Not-covered grandstand Indoor swimming pool Building for drinking water supply Normal building Building for public use Watercourse Pond Sluice

Mixed Forest Mixed Forest Mixed Forest Open Scrublands Open Scrublands Open Scrublands Open Scrublands Open Scrublands Open Scrublands Open Scrublands Open Scrublands Open Scrublands Grasslands Grasslands Grasslands Grasslands Permanent Wetlands Permanent Wetlands Croplands Impervious area Impervious area Impervious area Impervious area Impervious area Impervious area Impervious area Impervious area Impervious area Impervious area Impervious area Impervious area Impervious area Impervious area Impervious area Impervious area Impervious area Impervious area Impervious area Impervious area Impervious area Impervious area Impervious area Water Bodies Water Bodies Water Bodies

WL2021R00_130_8

Final version

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Table 3: Parameters characterizing land use classes in WetSpa (Liu and De Smedt, 2004)

Obviously, the root depth and Manning’s coefficient might not be determined by only using land cover information but also by the soil type. This is a limitation of the present implementation in WetSpa. In contrast, HEC-HMS has introduced another lookup table, taking land use, soil type information and hydrologic conditions into account to define its parameters (e.g. Curve Number when applying the SCSmethod). In MIKE-SHE, the vegetation distribution is also used to calculate the spatial and temporal variation of actual evapotranspiration (Kristensen and Jensen, 1975). Parameters of physical processes characterizing integrated effect of different soil types, land use types and other catchment characteristics Several model parameters are controlled by a combination of soil type, land use type, topographical and other catchment characteristics. One of the most important ones is the runoff coefficient (defined as the ratio of runoff volume over the rainfall volume, and directly used in WetSpa) or the SCS Curve Number (as used in HEC-HMS), the maximum soil retention, surface retention, etc. In WetSpa, the runoff coefficient depends on land use, soil type, topographical slope, and on rainfall intensity and antecedent soil moisture, where the latter are a model input and a modelled state variable. Given a near saturated condition of soil moisture, a simple lookup table for the runoff coefficient has been created with reclassified land use types to reduce the size of the table (see Table 4). Other tables were produced by Meert and Willems (2013) based on Chow (1964), De Smedt (1999, 2006), Vandekerckhove et al. (2001), Viessman and Lewis (2003) providing indications on the dependency between the runoff coefficient and land use including agricultural crop type, soil type, topographical slope and rainfall intensity for upstream catchments: see Table 5, Table 6 and Table 7. Runoff coefficient have also been derived for Flanders by Cabus and De Jongh (2007) but for catchment averaged conditions. Dependencies were identified with geological formation, land use, topographical slope and soil type (Figure 3-7 and Table 7). These may be useful as well, but it should be taken into account that they are valid for larger spatial zones only. In HEC-HMS, lookup tables for CN or maximum retention are applied based on soil type, land cover, cumulative precipitation and antecedent soil moisture. A composite CN value is calculated per grid cell by a given equation in case of complex soil type and land use conditions. Based on the CN value, the soil storage capacity can be computed based on the SCS method:

S= Final version

1000 − 10 CN

when the soil storage capacity is computed in inches WL2021R00_130_8

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

25400 − 254 CN

when the soil storage capacity is computed in mm

For Flanders, look-up tables for CN based on soil type, land use and initial soil moisture level, derived from literature review, have been proposed by Vandekerckhove et al. (2001) and Meert and Willems (2013) (Table 5). Meert and Willems (2013) show that this SCS method, after some assumptions (e.g. that the initial surface loss is 20% of the soil storage capacity), corresponds to the following formula for the runoff coefficient: ϕ=

P − 0.2 S P + 0.8 S

where P is the daily rainfall intensity in mm. Table 4: Potential runoff coefficient for different land use, soil type and topographical slope in WetSpa (Liu and De Smedt, 2004)

WL2021R00_130_8

Final version

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling Table 5: Runoff coefficient for different land use and vegetation types (combining information from Chow, 1964; De Smedt, 1999; Vandekerckhove et al., 2001; Viessman and Lewis, 2003; De Smedt, 2006)

Landgebruik Verhard Wegen

Daken Woonzones

Runoff-coëfficiënt 0.8 – 1.0 0.7 – 0.95 0.8 – 0.95 0.7 – 0.95 0.7 – 0.85 0.75 – 0.95 0.3 – 0.7 0.3 – 0.5 0.4 – 0.75 0.25 – 0.4 0.5 – 0.7 0.5 – 0.9 0.5 – 0.8 0.6 – 0.9 0.5 – 0.95 0.7 – 0.95 0.5 – 0.7 0.2 – 0.4 0.2 – 0.4 0.1 – 0.25 0.2 – 0.35 0.1 – 0.6 0.2 – 0.8 0.1 – 0.3 0.4 – 0.8 0.1 – 0.4 0.1 – 0.4 0.0 – 0.1

Beton Geasfalteerd Straatstenen

Villawijk Vrijstaande woningen Voorstedelijke agglomeratie Appartementsgebouwen Industriezones Verspreid Geconcentreerd Handelscentra Stedelijk Gemeentelijk Spoorwegen Tuinen Parken en begraafplaatsen Speelvelden Weiden Akkerland Braakliggende gronden Boomgaarden-laagstam Boomgaarden-hoogstam Bos Waterpartijen

Final version

WL2021R00_130_8

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling Table 6: Runoff coefficient for some coarse topographical slope ranges, soil type classes and land use types (Vandekerckhove et al., 2001, based on Tourbier & Westmacott, 1974)

Landgebruik Weiden

Akkerland

Bos

Helling

Zand-leem

Bodemtype Klei-leem

Zware klei

Vlak (0.5% helling) Hellend (5-10%) Sterk hellend (10-30%)

0.1 0.15 0.2

0.3 0.35 0.4

0.4 0.55 0.6

Vlak (0.5% helling) Hellend (5-10%) Sterk hellend (10-30%)

0.3 0.4 0.5

0.5 0.6 0.7

0.6 0.7 0.8

Vlak (0.5% helling) Hellend (5-10%) Sterk hellend (10-30%)

0.1 0.25 0.3

0.3 0.35 0.5

0.4 0.5 0.6

Table 7: Runoff coefficient for different rainfall intensity classes, two hydrological soil type classes, and some land use types (derived from Vandekerckhove et al., 2001, based on Horn and Schwab, 1963, and Schwab et al., 1993). For the hydrological soil types, see Table 8

Landgebruik en hydrologische condities Akkers – gewassen Rijgewas, slechte condities Rijgewas, goede condities Graangewas, slechte condities Graangewas, goede condities Weiden Hooiweide (rotatie), goede condities Graasweide (permanent), goede condities Bos Bos (volgroeid), goede condities

Neerslagintensiteit, voor bodemtype B, in mm/h 25 100 200

Bodemtype, voor neerslagintensiteit 100 mm/h A C D

0.63

0.65

0.66

0.58

0.71

0.73

0.47

0.56

0.62

0.48

0.61

0.64

0.38

0.33

0.42

0.44

0.18

0.21

0.22

0.18

0.23

0.24

0.29

0.36

0.39

0.29

0.41

0.42

0.02

0.17

0.23

0.11

0.21

0.22

0.02

0.1

0.15

0.05

0.13

0.14

WL2021R00_130_8

Final version

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Figure 3-7: Potential runoff coefficient for 74 gauged catchments in Flanders (the size of the circles are proportional to the value of the runoff coefficient). The blue background colour shows the spatial variation in runoff coefficient as expected from Table 7 in relation to the topographical slope, land use, soil type and geology (Cabus & De Jongh, 2007)

Table 8: Catchment averaged runoff coefficient for different classes of topographical slope, land use, soil type and geology (Cabus & De Jongh, 2007)

Geologie

Ieperiaan_Maldegem Ieperiaan_Maldegem Ieperiaan_Maldegem Ieperiaan_Maldegem Ieperiaan_Maldegem Ieperiaan_Maldegem Rupeliaan/Landeniaan Rupeliaan/Landeniaan Rupeliaan/Landeniaan Rupeliaan/Landeniaan Rupeliaan/Landeniaan Rupeliaan/Landeniaan Diestiaan/Scaldisiaan Diestiaan/Scaldisiaan Diestiaan/Scaldisiaan Diestiaan/Scaldisiaan Diestiaan/Scaldisiaan Diestiaan/Scaldisiaan

Final version

Bodemgebruik

Helling [%]

Akkers Akkers Weiland/bos Weiland/bos Braakland/bebouwing Braakland/bebouwing Akkers Akkers Weiland/bos Weiland/bos Braakland/bebouwing Braakland/bebouwing Akkers Akkers Weiland/bos Weiland/bos Braakland/bebouwing Braakland/bebouwing

<5 >5 <5 >5 <5 >5 <5 >5 <5 >5 <5 >5 <5 >5 <5 >5 <5 >5

WL2021R00_130_8

Zand Lemig zand Lichte zandleem 0.34 0.51 0.73 0.74 0.54 0.71 0.12 0.17 0.14 0.18 0.15 0.17 0.27 0.39 0.31 0.41 0.35 0.39

Bodemtype Zandleem Leem

0.62 0.88 0.62 0.68 0.47 0.74 0.12 0.20 0.14 0.19 0.17 0.22 0.27 0.46 0.33 0.42 0.39 0.50

Klei Zware klei 0.73 0.91 0.82 0.88 0.81 0.85 0.16 0.22 0.19 0.20 0.20 0.23 0.38 0.51 0.44 0.46 0.46 0.54

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling Table 9: CN-value for different land use types, agricultural crop types and hydrological soil type (see Table 5 for a description of these hydrological soil types) (Vandekerckhove et al., 2001)

Landgebruik of Teelttechniek gewas Braak Rechte rij Rijgewas Rechte rij Rechte rij Contour Contour Terras Terras Graangewas Rechte rij Rechte rij Contour Contour Terras Terras Dichtgezaaide leguminosen Rechte rij of Rechte rij weiderotatie Contour Contour Terras Terras Begraasde weide

Permanente weide Bos

Erf Verharde wegen

Contour Contour Contour

Hydrologische toestand -

A 77

Slecht Goed Slecht Goed Slecht Goed

72 67 70 65 66 62

81 78 79 75 74 71

88 85 84 82 80 78

91 89 88 86 82 81

Slecht Goed Slecht Goed Slecht Goed

65 63 63 61 61 59

76 75 74 73 72 70

84 83 82 81 79 78

88 87 85 84 82 81

Slecht Goed Slecht Goed Slecht Goed

66 58 64 55 63 51

77 72 75 69 73 67

85 81 83 78 80 76

89 85 85 83 83 80

Slecht Matig Goed Slecht Matig Goed

68 49 39 47 25 6

79 69 61 67 59 35

86 79 74 81 75 70

89 84 80 88 83 79

Goed

Slecht Matig Goed

45 36 25

66 60 55

77 73 70

83 79 77

WL2021R00_130_8

Hydrologische bodemgroep B C D 86 91 94

Final version

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling Table 10: Description of the hydrological soil types considered in the look-up table for CN (Vandekerckhove et al., 2001)

Bodemgroep A B C

Omschrijving Laagste afvoerpotentiaal. Omvat diepe zanden met zeer weinig leem en klei, ook diepe, snel doorlatende löss Matige lage afvoerpotentiaal. Meestal zandige gronden minder diep dan groep A, maar met een meer-dan-gemiddelde infiltratie na grondig natmaken. Matig hoge afvoerpotentiaal. Omvat ondiepe gronden en gronden met aanzienlijke hoeveelheden klei en colloïden, maar minder dan deze van groep D. De bodems hebben een minder-dan-gemiddelde infiltratie na voorverzadiging. Hoogste afvoerpotentiaal. Meestal kleien met een hoog zwellingspercentage, maar ook ondiepe bodems met bijna ondoordringbare substraten nabij het oppervlak.

Both HEC-HMS and WetSpa models consider the impervious percentage according to various land use types, especially for urban areas since this percentage strongly effects the runoff formation and due to the size of the mapping grid. For several land-use categories, the % imperviousness is defined and estimated based on different approaches as can be seen in Table 11 below. For other types of landuse, this parameter is assumed to be zero.

Table 11: Impervious percentage associated to specific land use types in WetSpa (Liu and De Smedt, 2004)

For coarse grid cells, the catchment characteristics needs to be weighed according the relative fractions of the different characteristics (e.g. land use types). Next to the runoff coefficient (or parameters related to the splitting components in the conceptual model, see later), also the runoff travel time or concentration time (or parameters related to the flow routing components in the conceptual model) represents integrated effects of spatial differences in topographical characteristics, land use and soil type. Meert and Willems (2013) discuss different methods for estimating this travel time based on these characteristics. The methods are based on kinematic wave equation, surface runoff based on Manning’s equation or various types of existing empirical equations.

Final version

WL2021R00_130_8

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Manning’s equation (also called Manning-Strickler or Gauckler–Manning–Strickler equation) is an empirical formula for the resistance loss Sf: 4 2

Sf = n q

P3 A

where:

10 3

q: the surface runoff discharge [m3/s] A: the wet area perpendicular to the main flow direction [m2] P: the wet perimeter perpendicular to the main flow direction [m] n: the Manning coefficient of the surface [s/m1/3] Making use of the kinematic wave equation, the resistance loss Sf is replaced by the surface slope S. If the water depth (h) of the surface runoff water is neglected in comparison with the width of the surface runoff wet area, or the runoff section approximated by a broad rectangular section (perpendicular to the main flow direction), the following assumptions can be made: q=U.A, A=h.b, P=b. Based on these assumptions, the runoff velocity can be computed based on Manning’s equation as follows: U=

where:

1 2 / 3 1/ 2 h S n

U: mean runoff velocity [m/s] h: runoff depth [m] The surface runoff discharge per unit width perpendicular to the main flow direction is given by:

q OF = i L = U h =

where:

1 5 / 3 1/ 2 h S n

i: net rainfall [mm/h] L: runoff distance [m] The mean travel time T then can be derived from T= L / U: T =7

L0.6 n 0.6 i 0 .4 S 0 .3

(46)

where: T: mean travel time [min] The Manning coefficient is empirically assessed (Table 12). Table 12 shows two types of coefficients: the theoretical one and the effective one. The theoretical one is valid for a surface with one single type of land use, equal for the whole surface. Such surface does, however, never exist in practise. In that case, the effective coefficient is considered, which takes into the effect of other factors that influence the travel time of runoff water at the surface next to the surface roughness. Examples of these influencing factors are: disturbances at the surface, impact of rain drops, soil erosion and sediment transport, the influence of the drag force at the surface, etc. Due to these factors, the Manning coefficient most frequently reaches values higher than 0.1. The effective values were obtained through experiments by Engman (1986) applying the kinematic wave equation together with an analysis of the surface storage based on observed runoff discharges at field scale and by Crawford and Linsley (1966) after calibration of the Stanford Watershed Model.

WL2021R00_130_8

Final version

Type ruwheid Landoppervlak Gladde oppervlakken (beton, asfalt, grind of kale grond) Braakland Weiland Akker Struikgewas Bos Waterlopen Recht Meanderend Begroeid Overwoekerd

Manning-coëfficiënt [s/m1/3] Theoretisch Effectief 0.01

0.01 – 0.014

0.02 - 0.03 0.03 0.03 - 0.05 0.05 - 0.1 0.1 - 0.15

0.05 – 0.3 0.06 - 0.2 0.1 – 0.4 0.2 – 0.6 0.4 - 0.8

0.03 0.04 0.05 0.1

Conceptual model parameters Even being conceptualized or presented in some “macroscopic” way, the parameters in conceptual models are still representing physical properties that are related to actual catchment characteristics. When disaggregating, some of them show a quite apparent (or close) linkage to the physical information; some of them require further inference or calculations. The idea proposed below is one way to approach the former ones while the latter ones demand more investigation by means of trialand-error or calibration. Following the proposed step-wise spatial disaggregation concept, as outlined above, the spatial conceptual model parameter values are set in two steps: -

First, the lumped conceptual model parameters (as obtained from the lumped conceptual model calibration at catchment scale; e.g. following the VHM approach) are transferred to uniform spatial conceptual parameter values (every grid cell same parameter values). The conceptual model parameters that are related to storage or subflow splitting and if defined per unit area, can be taken identical in the lumped and spatial conceptual models. The conceptual model parameters that are related to the routing, such as the recession constants, need to be scaled down to account for the spatial difference between the catchment scale and the single grid scale Second, spatial variations are given to the model parameters. This is done in a relative way against the uniform parameter values. This means that for each parameter the spatial average is taken identical to the uniform parameter value, but that relative spatial variations are given to that averaged value, proportional to the geographical catchment characteristics discussed in previous sections.

Flow routing component related parameters In conceptual models, the most important parameters for the flow routing components (and often the only parameters considered) are the reservoir routing recession constants. The recession constant of the quick runoff is strongly dependent on the distance to the river cell or the size of the catchment cell (depending on whether a Source-to-Sink approach is applied or a Grid-to-Grid method), the topography (i.e. the topographical slope) and the land use type. The distance to the river can be assessed by using the LDD map and the river reach definition map (UPS). The average number of cells (distance following the water path) to the river cell can be calculated.

Final version

WL2021R00_130_8

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

For a model based on Grid-to-Grid spatial linking, the recession constant derived from the lumped model can be disaggregated by dividing it by the calculated average number of cells to the river, to obtain the uniform parameter value. Although the average number of cells to the river (in the Grid-to-Grid spatial linking) can be used to obtain an initial estimate of the uniform value for the quick flow recession constants, further fine-tuning of this number might be needed. This is done by comparing the recession slopes in the decreasing flanks of the hydrographs with these in the observed flow series or lumped model results. See Figure 3-8 for an example of such validation for different numbers applied to scale the lumped model overland flow recession constant in the Dijle case (Tran & Willems, 2014); idem Figure 3-9 for the baseflow recession constant. In Figure 3-8, 6h is the overland flow recession constant for the lumped case; 2.2h is the recession constant when 6h is divided by the mean number of cells to the river. The other recession constant values shown are for illustration purposes only, to demonstrate the sensitivity of the results to the selected recession constant. In Figure 3-9, 3700h is the baseflow recession constant for the lumped case; 1357h when 3700h is divided by the mean number of cells to the river. The tested results shown in Figure 3-8 and Figure 3-9 confirm that the method suggested above for assessing the uniform recession constant values (dividing the lumped recession constant by the average number of cells to the river) works well. The simulation results indeed confirm that this approach leads to similar hydrograph shape at the catchment outlet. That the peak differs is of no relevance here, because only the hydrograph shape is evaluated. When an average value is taken for the recession constant instead (6h for the overland flow recession constant in the Dijle case), the hydrograph slope becomes too flat as expected.

Figure 3-8: Simulation results (log-scale) of the Grid-to-Grid spatial linking model of the Dijle case for uniform parameter values (spatial disaggregation step 1): validation and fine-tuning of the recession constant for overland flow (Tran & Willems, 2014)

WL2021R00_130_8

Final version

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Figure 3-9: Simulation results (log-scale) of the Grid-to-Grid spatial linking model of the Dijle case for uniform parameter values (spatial disaggregation step 1): validation and fine-tuning of the recession constant for baseflow (Tran & Willems, 2014)

In the second step, the uniform values for the recession constants are spatially varied. For the overland flow recession constant, this is done based on topographical slope and land use. One way is to make use of the kinematic wave routing equation. Note that there is no need to apply the kinematic wave equation as such, because one can start from the lumped calibrated recession constant, but give it a spatial variation using a factor that is proportional to the influencing characteristics as in the equation. This means, based on equation (46), that the factor should be proportional to L0.6, and to n0.6, and to S-0.3, for instance by using a factor = ( L0.6 n0.6 S-0.3 ) / ( Lmean0.6 nmean0.6 Smean-0.3 ) where L, n and S are size, Manning coefficient and slope for each grid cell, and Lmean, nmean, Smean the uniform values (mean values over all grid cells). Using such approach, when coarser grid cells would be applied, the L, n, S values would become averages of the values of the underlying subgrid cells. This means that for the lumped case, L, n, S of the single grid cell would become the mean of L, n, S for all grid cells, hence leading to a factor of 1.0. This demonstrates how the proposed parameter disaggregation works. A similar approach can be applied to distribute the baseflow (slow flow) recession constant, but additionally making use of the soil layers’ information (hydraulic conductivity). Figure 3-10 illustrates the outcome of the spatial variation of the overland and interflow recession constants for the Dijle case (Tran & Willems, 2014).

Final version

WL2021R00_130_8

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

CKOF

CKIF

Figure 3-10: Spatial disaggregation of the quick runoff recession constants for overland flow (CKOF) and interflow (CKIF) (in hours) for the Dijle case (Tran & Willems, 2014)

Soil storage component related parameters The most important parameter of the soil storage component in conceptual models is the maximum soil moisture storage (soil water capacity). This maximum storage is physically linked to the difference between the field capacity and the soil moisture at the permanent wilting point (Total Available (soil) Water - TAW). It is also defined by the rooting depth of the vegetation, which is linked with land use. Therefore, the spatial variation in the lumped maximum soil storage parameter can be related to the spatial variation in maximum soil retention as expected from the integrated effect of soil type, land use and other catchment characteristics, as discussed in section 3.2.2.1. In the Dijle case, for example, the information in the 5m land use and soil texture maps is used to investigate the combined effect of rooting depth and soil water retention. The spatial variation reflected in the land use and soil texture maps is used to distribute the lumped maximum soil storage parameter to a 5m grid. This distributed map is resampled to a 100m grid, by either one of two methods: by selecting the dominant value per 100m grid cell or by calculating the weighted average value. The same approach is applied for the Grote Nete case, where the existing 250 m (and 100 m) spatial maps are used. The TAW map can be generated using the soil type map only. A general look up table which links the different Belgian soil classes to their TAW values was derived. Moreover, the spatial variation of root zone depth corresponding to different land use types was also mapped. For calculation, an intermediate variable, named “TAW_rootdepth” was created to combine these two lookup tables into one as shown in Table 13. This variable was applied to give the maximum soil storage capacity spatial variability.

WL2021R00_130_8

Final version

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling Table 13: Combined lookup table for TAW_rootdepth as a function of soil type and root zone depth S

(Silty

(Light sandy

(Sandy

(Heavy

U Non

Soil type

(Sand)

Sand)

silt)

(Silt)

(Clay)

clay)

class

TAW (%)

0.038

0.082

0.105

0.123

0.141

0.121

0.127

0.001

Root zone depth (m)

Land use types

1.0

Evergreen Needleleaf Forest

38.0

81.5

105.0

123.3

140.5

121.0

127.0

1.0

Evergreen Broadleaf Forest

38.0

81.5

105.0

123.3

140.5

121.0

127.0

1.0

Deciduous Needleleaf Forest

38.0

81.5

105.0

123.3

140.5

121.0

127.0

1.0

Deciduous Broadleaf Forest

38.0

81.5

105.0

123.3

140.5

121.0

127.0

1.0

Mixed Forest

38.0

81.5

105.0

123.3

140.5

121.0

127.0

1.0

0.8

Closed Shrublands

30.4

65.2

84.0

98.7

112.4

96.8

101.6

0.8

Open Shrublands

30.4

65.2

84.0

98.7

112.4

96.8

101.6

0.8

1.0

Woody Savannah

38.0

81.5

105.0

123.3

140.5

121.0

127.0

1.0

0.8

Savannah

30.4

65.2

84.0

98.7

112.4

96.8

101.6

0.8

Grasshlands

30.4

65.2

84.0

98.7

112.4

96.8

101.6

0.8

0.5

Permanent Wetlands

19.0

40.8

52.5

61.7

70.3

60.5

63.5

0.5

0.8

Croplands

30.4

65.2

84.0

98.7

112.4

96.8

101.6

0.8

0.5

Urban and Built-Up

19.0

40.8

52.5

61.7

70.3

60.5

63.5

0.5

0.8

Cropland/Natural vegetation

30.4

65.2

84.0

98.7

112.4

96.8

101.6

0.8

0.1

Snow and Ice

3.8

8.2

10.5

12.3

14.1

12.1

12.7

0.1

0.5

Barren or Sparsely vegetation

19.0

40.8

52.5

61.7

70.3

60.5

63.5

0.5

0.1

Water bodies

3.8

8.2

10.5

12.3

14.1

12.1

12.7

0.1

NAN

Non existing land use

The resulted parameter map for the Dijle catchment is shown in Figure 3-11 (Tran & Willems, 2014). The comparison was done between the map obtained based on the dominant land use type in each 100m grid cell and the map obtained using the weighting average values of the underlying 5m land use grid cells. It is clear in Figure 3-11 that the second map leads to more smooth spatial variations in parameter values, and the presence of urban spatial patterns such as roads and urban centers becomes more visible (such as the highway and the airport). This is because urban areas very often have a limited spatial size, hence 100m grid cells often do not show urban as the most dominant land use type. This leads to underestimated urban fractions in the coarser resolution maps when these are based on the dominant land use type. The problem, however, can be solved by calculating the parameter values of the coarser map based on a weighted average of the parameter values obtained for the highest resolution map available.

Final version

WL2021R00_130_8

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Figure 3-11: Spatial disaggregation of the maximum soil storage parameter (in mm) (without applying a controlling factor) for an example of the Dijle catchment: (left) based on the dominant land use type per 100m grid cell; (right) based on the weighted average of the parameter values obtained for the 5m land use map (Tran & Willems, 2014)

Surface storage related parameters The maximum water content in the surface storage mainly depends on topography and land use, e.g. crop type. It moreover is related to the crop growth stage, such that this model parameter in fact changes over time. This time variation is, however, most often neglected in models. This is also done for the Dijle and Grote Nete implementation (see Tran & Willems, 2014, for both cases). Splitting process related parameters The parameters controlling the splitting components are in most conceptual models the runoff coefficients or the equations describing the dependency between these runoff coefficients and model state variables such as the relative soil saturation level. The dependency of the (surface) runoff coefficient to land use type, agricultural crop type, hydrological soil type, etc. was extensively discussed in section 3.2.2.3. Through the spatial variations in maximum soil moisture storage and the corresponding variations in soil saturation levels, the runoff coefficients also vary in space. These can be computed from the model simulation results (e.g. mean long-term surface runoff coefficient is computed as the ratio of the cumulated surface runoff volume over the cumulated rainfall), and compared with the differences as expected in relation to e.g. land use and soil type (section 3.2.2.3). If these are insufficient, additional variations in the parameters of the equations describing the dependency between these runoff coefficients and the relative soil saturation level, such as the maximum runoff coefficient, can be applied. Some comparisons were carried out as following: Figure 3-12 shows the spatial variation of the runoff coefficient at different time steps of the simulation. A clear dependency on the land use map is visible: higher coefficient for the urban grid cells; lower for agriculture and forest. More detailed quantification and investigation is, however, required.

WL2021R00_130_8

Final version

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Original Corine land-use map

Runoff coefficient at time step 15000h

Runoff coefficient at time step 35000h

Runoff coefficient at time step 75000h

Figure 3-12: Distribution of runoff coefficients over the Dijle catchment at different time steps of the simulation (Tran & Willems, 2014)

3.2.3

Adjustment/calibration of the spatial disaggregation of model parameters based on spatial catchment characteristics

The spatial disaggregation described in section 3.2.2 is based on assumptions regarding the dependency between model parameters and geographical catchment characteristics. However, this dependency may vary from region to region or catchment to catchment, hence some adjustments (calibrations) may be needed. For instance, the factor = (L0.6 n0.6 S-0.3) / (Lmean0.6 nmean0.6 Smean-0.3) can be applied to give spatial variations to the recession constants (see 3.2.2.4). By doing that, even though the mean recession constant is remained equal to the lumped parameter value, the spatial model may lead to different results at the catchment outlet in comparison with the lumped model. This needs to be checked; if the spatial model results show improved results, this obviously would be good; if not, adjustments to the implemented spatial variations in model parameters are required. Two calibration parameters are introduced for that purpose: one parameter (F1) to increase or decrease the catchment mean value for each model parameter, and a second one (F2) to make the spatial catchment variations stronger or less strong (relative factor applied to these spatial variations). Also the spatial differences in runoff flows as observed at internal stations, and in groundwater and soil moisture results can be checked. Taking the recession constant again as an example, the two calibration parameters can be applied as described next.

Final version

WL2021R00_130_8

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

For example, for spatial disaggregation of the recession constant K (from its uniform model parameter Kuniform to the spatial one Kspatial), before the calibration, the following equation is applied: K spatial = K uniform

L0.6 n 0.6 S −0.3 Lmean nmean 0.6 Smean −0.3

(47)

0.6

The following modified equation is introduced including the adjustment/calibration parameters F1 and F2: K spatial = K uniform + F1 + F2 K uniform (

L0.6 n 0.6 S −0.3 − 1) Lmean nmean 0.6 Smean −0.3 0.6

Another example is the spatial disaggregation of the maximum soil moisture storage parameter (see Figure 3-13). The original equation is as below: SMmax_spatial= SMmax_lumped * RootTAW_ratio

(48)

However, this equation has been modified as follows: SMmax_spatial= SMmax_lumped + F1 + F2 * SMmax_lumped * (RootTAW_ratio-1) For F1=0 and F2=1, the original equation (47 or 48) is obtained (no adjustment of the original spatial disaggregation). For F1=0 and F2=0, the spatial variation is reduced to zero, hence the uniform parameter values are obtained. By taking F2>1, the spatial variations are increased; by taking F2<1, the spatial variations are reduced. F1 and F2 are optimized by maximizing the model performance. See Figure 3-13 for an example of the spatial distribution of the maximum soil storage parameter, before and after different F2 adjustments (F1=0).

Original map ( F2 = 1.0)

F2 = 0.8

F2 = 0.5

F2 = 0.2

Figure 3-13: Distribution of maximum soil storage parameter, before and after adjustments (Tran & Willems, 2014)

WL2021R00_130_8

Final version

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

4 Model performance evaluation method A model performance evaluation protocol that evaluates the hydrological model performance by different criteria is proposed for the framework. Both qualitative (based on visual, graphical techniques) and quantitative (based on statistical measures) methods will be applied. These methods aim to consider a multi-objective set of preferably unrelated or slightly correlated statistics and the consideration of supporting graphical criteria. The proposed protocol aims to support the serious analysis, testing and discussion of model performance, as recommended by Jakeman et al. (2006), Robson et al. (2008), Welsh (2008), among others, in order to increase the credibility and impact of results from the hydrological models. Firstly, the runoff results will be evaluated at the basin outlet or the flow gauging station(s). This evaluation is applicable also to the lumped models. Secondly, additional evaluation is done of the spatially variable model outputs for which also observations are available (aquifer heads, soil moisture).

4.1 Graphical evaluation of runoff flows and response behavior at gauging stations The runoff results at the flow gauging stations are evaluated by (i) comparison of the observed and simulated river flows at the gauging stations, (ii) evaluation of a set of statistical measures, (iii) control of the water balance and (iv) analysis of different time series simulation properties. This evaluation is applicable also to the lumped models. This model performance evaluation is hereafter summarized and results illustrated after application of the lumped models PDM, NAM and VHM to the Dijle case study (catchment zone upstream of the flow gauging station at Sint Joris Weert). The models were calibrated based on the calibration approach presented in section 3.1 for the lumped case. They are based on the results obtained in the TWOL research project on “Rainfall generator for flood risk assessment” for VMM-AOW (Willems et al., 2014b). The Dijle catchment upstream of Sint Joris Weert was one of the three case studies considered in that study. These results are considered as the reference to compare the results of the spatial model versions with. The following model evaluation plots are constructed: -

Time series of total runoff flows;

Time series of total runoff flows, separately evaluated for winter and summer events;

Cumulative volumes of total runoff flows;

Time series and cumulative volumes of runoff subflows: baseflow / slow flows and quick flow;

Scatterplot of simulated versus observed peak flow extremes (maximum flows during independent quick flow periods);

Scatterplot of simulated versus observed low flow extremes (minimum flows during independent slow flow periods);

Empirical extreme value distributions of peak flow extremes;

Empirical extreme value distributions of low flow extremes.

Final version

WL2021R00_130_8

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

The subflows (quick and slow flows) and peak and low flow extremes are extracted from the runoff time series by the numerical filter and the method to split the flow time series in nearly independent quick and slow flow events, as discussed in section 3.1.1. Examples of these graphical model evaluation plots are shown below for the PDM, NAM and VHM lumped model results for the Dijle catchment upstream of Sint Joris Weert. The models used an hourly time step. The period 2004-2006 was considered for model calibration, and the period 2010-2012 for model validation.

Figure 4-1: Time series of total runoff flows, for the most extreme events in the simulated time series

Figure 4-2: Time series of slow runoff flow (baseflow), for calibration and validation periods

WL2021R00_130_8

Final version

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Figure 4-3: Time series of cumulative total runoff volumes, for calibration and validation periods

Figure 4-4: Scatterplot of peak flows

Final version

WL2021R00_130_8

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Figure 4-5: Return period plot of peak flows, total simulation period

Figure 4-6: Scatterplot of low flow extremes

WL2021R00_130_8

Final version

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Figure 4-7: Return period plot of low flow extremes, total simulation period

In the scatterplots (examples Figure 4-4 and Figure 4-6), the model results qm (vertical axis) are compared against the flow observations qo (horizontal axis). The model residuals are shown in this plot as the horizontal or vertical differences between each point and the bisector. In the scatterplots, the model residual mean and standard deviation are represented by the solid and dotted lines. The axes of the scatterplot are transformed based on the Box-Cox transformation, in order to reach homoscedastic residuals (see section 4.3 for more details). Runoff model residuals are typically heteroscedastic. This means that the model residual standard deviation depends on the runoff flow magnitude. The residuals are homoscedastic when this standard deviation is constant. In the scatterplots of simulated versus observed peak flow or low flow extremes and in the empirical extreme value distributions of these extremes, the peak flows are obtained as the maximum flows during independent quick flow hydrograph periods, and the low flows as the minimum flows during independent slow flow periods. The full flow time series can be separated in these quick and slow flow periods using the methods presented in section 3.1.1.3. More on the need to use independent values here is discussed in section 4.4. It is the aim of the model calibration to minimize both the model residual mean and standard deviation, or combined by the MSE (see next in section 4.2), or to reduce the scattering of the points in the scatterplot. After a good calibration, the bias needs to be reduced to zero (to become an unbiased model). For an unbiased model, the overall model performance can then be represented by a single value: the model residual standard deviation (also the horizontal or vertical distance of the dotted lines in the scatterplot to the bisector.

4.2 Statistical evaluation of runoff flows at gauging stations Next to the graphical model performance evaluations, the following goodness-of-fit statistics are computed: the mean absolute error (ME), the root mean squared error (RMSE), the regression coefficient (R²) and the Nash-Sutcliffe (1970) coefficient of model efficiency (NSE). These statistics are based on the model residuals, being the difference between each flow observation qo(i) and the corresponding model result qm(i) based on a number of different (i=1,n) observations. The mean of the 2

model residuals (ME), the variance ( S E ) and the Mean Squared Error (MSE) can be calculated as standard statistics (see e.g. Neter et al. (1996) or any basic statistical handbook). The model residual Final version

WL2021R00_130_8

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

mean is commonly used to measure the average systematic difference between the simulated and the observed values. The MSE and the model residual variance or standard deviation measure the average random discrepancies between simulations and observations. When n is sufficiently large, MSE is practically equal to the sum of the squared model residual mean and the model residual variance. The MSE thus is built up of two terms, where the first one reflects the bias or systematic deviation in the model results and the second one the variance of the random model residual after bias correction. The model residual mean and variance therefore can be used for goodness-of-fit evaluation of the model, analyzing the bias and the random error in a separate way, or the two effects combined in the MSE. The above-mentioned statistics have the disadvantage that their magnitudes highly depend on the flow magnitude, and thus on the river under study. This makes it difficult to compare the relative performance of different models applied to different river cases. Therefore, other model performance indices have been proposed in literature, for which the values range between 0 and 1 (0 or 1 for a perfect model). The most popular one is the model efficiency coefficient EF of Nash & Sutcliffe (1970), which is a dimensionless and scaled version of the MSE. The EF index is depicted next, where q o and

S q2o are the average and variance of the observed values:

 n 2  ∑ (qm (i ) - qo (i ) )     = 1 - n MSE  EF = 1 - i =1n 2 2    (n - 1) S q o  q i q ( ) ∑ o o   i =1  

(

)

The MSE based term can be split in the separate contributions of the ME bias and S E variance terms as follows:

EF = 1 -

S E2

S q2o

n ME 2 (n - 1) S q2o

Other goodness-of-fit statistics used in the literature are the mean relative error, the mean cumulative error (as applied in e.g. Johnson et al., 2003; Jain & Srinivasulu, 2004), Young’s Information Criterion (Young, 1986), etc.; for an overview see Legates and McCabe (1999) and Dawson et al. (2007). Some modellers apply the statistics only to specific flow conditions (e.g. more focus on high or low flows; e.g. Wagener et al., 2001; Madsen, 2000; Fenicia et al., 2007). The use of these other statistics most often aims to overcome problems reported hereafter on the use of the basic goodness-of-fit statistics. If the recommendations of next paragraph are followed, there is however no need to use the other statistics, as they solve the problems only in an indirect way. The goodness-of-fit statistics, such as the model residual mean and variance, MSE and EF, are most widely used for evaluating model performance, but have the disadvantage that they largely summarize the goodness-of-fit information only in a few numbers and values. Therefore, the recommendation is made here to complement the statistics with graphical goodness-of-fit plots (See section 4.1). In general, these plots compare the qm and qo values, as used on the basis of the statistics in a graphical way, and provide the modeller with far more complete information about the goodness-of-fit. For instance, in some cases a good model can lead up to poor goodness-of-fit statistics due to one or a few outliers in the residuals. The plot with underlying values can make this clear. The plot might provide the modeller with information, in case of a real poor performance, on whether this poor performance is to be explained by few outliers, by the majority of points, by the high flow events, the low flow events, etc.

WL2021R00_130_8

Final version

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

4.3 Need for transformation of the flow variable: logarithmic or BoxCox transformation In rainfall-runoff models, model residuals typically increase with higher flow values. This means that the model residual variance or standard deviation typically increases with increasing flow. It also means that the higher flow values receive more weight in the goodness-of-fit statistics. Mainly for the MSE and the EF statistics, which have squared model residual terms in the equations, it is clear that they are oversensitive to peak values for the model residuals. To solve this problem, some authors have proposed variations to the mathematical expression of the statistics to overcome this problem (e.g. by applying weighting factors to the model residual terms; see also Legates and McCabe, 1999; Madsen, 2000, 2003; Vázquez, 2003). Another approach, which is recommended here, consists in the transformation of the variable of interest, without the need to modify the equations. Different transformation forms can be found in literature. A non-parametric normal quantile transformation can be applied to transform the model residuals into a normal distribution, see e.g. Kelly & Krzysztofowicz (1997) and Montanari & Brath (2004). Other transformations are the square root transformation (e.g. applied in Kokkonen & Jakeman, 2001), the logarithmic or ln-transformation (e.g. Oudin et al., 2006), the Box-Cox transformation (Box & Cox, 1964), etc. The Box-Cox (BC) transformation is a very flexible one because – depending on its parameter value – it can cover a wide range from weak to strong transformations. The BC-transformation, when applied to a variable q, is given by: qλ −1 BC ( q ) = λ where the parameter λ needs to be calibrated in an attempt to reach homoscedasticity in the residuals (constant model residual variance, independent on q). The value of λ can range from 0 to 1. Calibration of the parameter λ can be done graphically as illustrated in Willems (2009a, 2013). Evaluation on whether the standard deviation is increasing, constant or decreasing can be done by calculation of the standard deviation for different ranges of the flow (e.g. moving window technique), or simply in a visual way by checking the uniformity in the magnitude of scattering around the bisector. The visual check is most often sufficiently accurate, given that the parameter λ only needs to be assessed in an approximate way (differences in λ up to 0.2 only have minor influences on the goodnessof-fit evaluations). When visually checking the uniformity of the scattering magnitude, one has to take into account that the points may not be equally distributed along the flow range. Typically fewer points will occur for the higher flow range. For this reason, even when the standard deviation is uniform, the largest deviation from the bisector will typically occur at the low flow range. For runoff discharges the parameter λ usually adopts a value of around 0.25. The strongest transformation in the range [0, 1] for λ is the logarithmic transformation. The BC-transformation converges asymptotically to the ln-transformation for λ towards 0. This means that for the case λ=0 the standard deviation becomes uniform after ln-transformation of both axes in Figure 4-4 and Figure 4-6. This is the case when the standard deviation of the residuals is proportional to the flow value. The model then has a uniform relative error (the ratio of the standard deviation over the flow value is constant). When after ln-transformation, an increasing standard deviation with increasing flow value is turned into a decreasing standard deviation, the transformation is too strong and a λ value between 0 and 1 is required to reach homoscedastic errors. When λ equals 1, no transformation is needed, and the standard deviation of the absolute error is uniform. In most practical cases, both the standard deviation of the absolute error neither the standard deviation of the relative error is uniform and a BC transformation is needed. If no transformation would be applied in this case, higher weights would be given to the higher flow values in the goodness-of-fit statistics (as also discussed before), and if the logarithmic transformation would be applied, or relative

Final version

WL2021R00_130_8

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

errors used in the goodness-of-fit statistics, then higher weights would be given to the lower flow flows. From the above, it becomes clear that more specific transformations, as used by some modellers, such as the square root transformation, the logarithmic transformation, etc. only will be valid for datasets with specific error structures. The BC-transformation is a more general method to reach independence of the residuals on the flow value. Homoscedasticity in the residuals has the advantage that the model residuals can be represented by 1 single distribution (most often normal or Gaussian). In case the model is unbiased (zero model residual mean), this distribution has only 1 parameter (the model residual variance or standard deviation), and the total model uncertainty or goodness-of-fit can be represented by this single value. When this standard deviation is optimized during model calibration, equal weight is given to the high and the low flow values.

4.4 Need to address serial dependence of flow When series of observations are used, model residuals often have a serial dependence. This dependence will be higher for a series with a smaller time step, and will also depend on the flow magnitude. For peak flows, the serial dependence is strong only at small time scale (e.g. at hourly time steps, or smaller than the recession constant of the quickest subflow in the catchment rainfall-runoff). For low flows or dry periods, the serial dependence is higher, due to the longer recession constant of the baseflow component in comparison with the quick runoff components. This causes additional problems in case the flow values at all time-steps are used in the goodness-of-fit statistics. Higher weight would then indeed be given to the lower flow values. This problem is most often overlooked in current rainfall-runoff calibration applications. Some authors solve this problem by assuming an autocorrelation structure for the residuals, e.g. a lag-1 correlation (e.g. Sorooshian & Dracup, 1980). Often in this approach ignorance is made of differences in the autocorrelation for different flow magnitudes. The correlation structure moreover might be more complicated. For quick flow events, for instance, the correlation of the residuals is very high (almost 100%) within the same quick flow event (driven by the same rainfall event), and low (almost 0%) in between two successive events. Under these conditions, nearly independent observations can be selected from the flow series by means of hydrograph separation. Based on the techniques discussed in section 3.1.1.3, the flow series can be split in (nearly independent hydrograph) periods, and peak values selected from subsequent quick flow hydrographs, cumulative volumes during these hydrographs, low flow values during subsequent baseflow recession periods, etc. The sum (for I = 1,n) in the goodness-of-fit statistics then goes over the n independent hydrographs or flow periods in the flow series. For use in peak flow based goodness-of-fit statistics, the highest flow value can be selected for both the observed and the simulated series in each of the independent hydrograph periods, and compared. The peak flows for a given hydrograph period selected from the observed and simulated series do not need to occur exactly at the same time moment. This has an additional advantage for flow series with a small time step (e.g. hourly), where small time shifts might occur in the model simulation results due to rain gauge locations different from the central location of the catchment. In case small time shifts in modelled versus observed peak flows do not pose a problem, they should not influence the goodness-of-fit statistics. After selection of independent values and after BC transformation, it is feasible to derive homoscedastic and uncorrelated residuals. In comparison with other, more standard approaches in literature, where independence and BC transformation is not considered, this procedure allows for a more unbiased model performance assessment.

WL2021R00_130_8

Final version

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

4.5 Evaluation of catchment response behavior To further investigate the difference in peak flow extremes and related model structure in PDM in comparison with NAM and VHM, the overland and quick runoff coefficients, estimated as the ratio of overland or quick flow amount versus total rainfall during each overland flow event, are examined in function of the soil moisture content for the different models (Figure 4-8 and Figure 4-9). Figure 4-8 shows that the NAM model produces lower quick runoff coefficients than PDM and VHM. The PDM quick runoff coefficients increase in a more exponential way (increasing slope of the runoff coefficient for higher soil moisture level) but have less scatter. The exponential increase probably explains the higher simulated flow extremes by PDM. When the cumulative distribution of the quick runoff coefficients are compared (Figure 4-9), the distributions differ only slightly.

Figure 4-8: Runoff coefficient based on overland flow (left) and combined interflow & overland flow (quick flow; right) versus relative soil moisture level downstream. Example from the Grote Nete catchment for PDM, NAM and VHM.

Figure 4-9: Cumulative frequency distribution of the runoff coefficient for the NAM, PDM and VHM model based on overland flow (left) and combined interflow & overland flow (right) downstream. Example from the Grote Nete catchment for the NAM, VHM and PDM model.

The cumulative frequency of the VHM runoff coefficients matches the observations best, which is a direct consequence of the stepwise calibration procedure of VHM in which the overland flow and interflow are sequentially fit to their filtered components. The cumulative frequency of the NAM overland flow coefficient has a good match with the observed ones but when accounting for the interflow, the distribution differs from the observation. The lower runoff coefficients have a higher frequency, while the higher runoff coefficients occur less frequently in comparison with the observations. This indicates that NAM overestimates lower interflows and underestimates the higher interflows. Analysis of the PDM frequency distribution showed opposite results with a lower frequency for the low runoff coefficients and a higher frequency for the high runoff coefficients. These findings Final version

WL2021R00_130_8

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

are consistent with the previous observations/conclusions based on the extreme value distributions: PDM overestimates the highest peak flows.

4.6 Effect of rainfall change Next question after a thorough calibration and validation of the models based on the above-outlined model performance evaluation protocol, is whether the assumptions made by the model are limiting their applicability for scenario investigations. For this particular study, concern may arise in the justification of applying hydrological models under rainfall conditions that are more extreme than the ones considered during model calibration and validation, e.g. more extreme rainfall conditions generated by the stochastic rainfall generator or after simulation of changed climate conditions. The calibrated models might appear robust for the present-day variability of climate, but they might produce uncertain results regarding the extrapolation to more extreme rainfall conditions or changes in future runoff. In this context, a model evaluation technique has been developed by Van Steenbergen and Willems (2012) that analyses whether the models are suitable to simulate the observed changes in high and low flow extremes related to the observed changes in rainfall rates. For evaluating the accuracy of the model response during extreme events, the simulated high flows are firstly coupled to their underlying rainfall. Overland flow, quick flow or peak flow maxima are considered during the different nearly independent quick flow periods. These flow maxima will be linked to rainfall volumes during the recession period preceding the flow maxima. Next, relations between changes in rainfall amounts and the changes in overland flow are examined and compared to the observations. More precisely, different classes of rainfall increases (e.g. 5-15%, 15-25%,…) are considered and factors flow change computed and evaluated for each of these classes (Figure 4-10 and Figure 4-11).

Figure 4-10: Cumulative frequency distributions of the relative change in quick flow event volume per rainfall change class, for the Dijle catchment at Sint-Joris-Weert (Tran & Willems, 2014)

WL2021R00_130_8

Final version

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Figure 4-11: Quantile analysis for the relative change in overland peak flow versus change in rainfall, for the Dijle catchment at Sint-Joris-Weert (Tran & Willems, 2014)

4.7

Evaluation of spatially variable model outputs

Spatially variable model outputs will be evaluated based on the internal river flow gauging stations in the catchment, the groundwater well heads and the soil moisture levels. This will be done both graphically and statistically, using similar plots and statistics as discussed before for the evaluation of the runoff flows at the gauging stations. Examples of graphical evaluations of spatial groundwater head results are shown in Figure 4-12, Figure 4-13 and Figure 4-14.

Figure 4-12: Internal river flow gauging stations and calibration and validation groundwater wells considered by Vansteenkiste et al. (2012) for the Grote Nete catchment. Final version

WL2021R00_130_8

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Figure 4-13: Graphical evaluation of the temporal groundwater dynamics. Example of MIKE SHE head levels for calibration and validation wells in the Grote Nete catchment.

WL2021R00_130_8

Final version

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Figure 4-14: Graphical evaluation of the spatial groundwater dynamics. Example of MIKE SHE head levels along transects across the Grote Nete catchment.

Final version

WL2021R00_130_8

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

The usefulness of soil moisture products is highly questioned in literature, but there also exist very successful soil moisture validation examples in literature. The usefulness of the soil moisture products and the evaluation of the soil moisture results for catchments in Flanders is therefore very unclear to date. Some (both negative and positive results) are given hereafter. Lee et al. (2007) compared for a catchment in Germany 60 cm depth averaged TDR measured surface soil moisture content with saturation degree in the unsaturated zone, simulated by a rainfall-runoff model, but no useful comparison was found (Figure 4-15). Zhang et al. (2011b) studied the correlations between measured runoff ratio and measured antecedent soil moisture for four catchments but did not find clear correlations (Figure 4-16). As opposed to these negative results, very positive result are shown by Brocca et al. (2012b), supported by their results in Brocca et al. (2010, 2011, 2012a). They used the ASCAT soil moisture product for catchments in Italy. After rescaling these soil moisture product data, they matched the relative soil moisture simulated by their model very well (Figure 4-17). Use of these rescaled soil moisture data as input for real-time data assimilation to update the rainfall-runoff model soil moisture state variable revealed a significant increase in the river flood forecast results (Figure 4-18). The Nash-Sutcliffe efficiency of these results increased from 0.79 to 0.86 (Figure 4-19). They, however, concluded that a flexible rainfall-runoff modelling approach is required in order to match the available spatial data with the model state variables at the proper spatial and temporal scales. Another positive example is shown in Figure 4-20 for a catchment in Drente, The Netherlands, where Terink et al. (2012) found a good match between observed (DACOM sensors) and model simulated soil moisture. Good match between observed and simulated soil moisture has also been shown by Tavakoli (2012) and Tavakoli and De Smedt (2012, 2013) for the Wetspa model applied to the Fork river basin in Oklohoma, USA, the Vermilion catchment in Illinois, USA, and the Bibeschbach catchment in the Grand-duchy of Luxembourg. For the Fork river basin, the TDR soil moisture measurements were well reproduced by the model with errors that are within the range of the measurement uncertainty, although the model results show somewhat more abrupt temporal fluctuations as the whole root zone is considered to react to changes in catchment rainfall. Figure 4-21 shows an interesting example of how soil moisture data can be applied for evaluating the submodel describing the saturation excess dependent runoff coefficient (Penna et al., 2011). Using ground-penetrating radar (GPR) data for Belgium, recent research by Tran (2014) and Tran et al. (2014) has shown that the use of these data for data assimilation on the vertical soil moisture profiles in a hydrological model (in their case limited to field scale) may improve the accuracy of such model significantly. They also have shown that the vertical soil moisture profiles can be reconstructed in the most accurate way by combining the soil moisture remote sensing data with models. And, as for groundwater levels, study of the spatial variability of soil moisture would be very useful as well.

Figure 4-15: Comparison of 60 cm depth averaged measured surface soil moisture content measured by 61 TDR sensors at weekly intervals, with saturation degree in the unsaturated zone simulated by a rainfall-runoff model, for the Weiherbach catchment in Germany (Lee et al., 2007).

WL2021R00_130_8

Final version

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Figure 4-16: Correlations between measured runoff ratio and measured antecedent soil moisture for four catchments in the Walnut Gulch experimental watershed in Arizona, USA (Zhang et al., 2011b).

Figure 4-17: Soil moisture data from the ASCAT product rescaled to match the relative soil moisture simulated by the model, for Tiber river basin in Italy (Brocca et al., 2012b).

Final version

WL2021R00_130_8

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Figure 4-18: Example of runoff flow forecast results with observed flows, before and after model state updating based on data assimilation applied to the rescaled ASCAT soil moisture product, for Tiber river basin in Italy (Brocca et al., 2012b).

Figure 4-19: Example of runoff flow forecast results with observed flows, before and after model state updating based on data assimilation applied to the rescaled ASCAT soil moisture product and change in Nash-Sutcliffe efficiency, for Tiber river basin in Italy (Brocca et al., 2012b).

WL2021R00_130_8

Final version

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Figure 4-20: Comparison of measured soil moisture based on DACOM sensors and model results based on the SPHY model for the Hunze catchment in Drente, The Netherlands (WaterSense project; Terink et al., 2012).

Figure 4-21: Runoff coefficient versus TDR measured 30 cm depth averaged soil moisture prior to the event, for the Rio Vauz basin in Italy (Penna et al., 2011).

Final version

WL2021R00_130_8

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

5 Evaluation results of spatial models for Grote Nete catchment The proposed concept was implemented in PCRaster-Python with three different conceptual hydrological models (NAM, PDM and VHM) that were used as the based model. The models were disaggregated and calibrated into spatial distributed and applied to Grote Nete catchment in Belgium. The final results were evaluated by plotting versus the results of the lumped models results and the observations using the same approach in Chapter 4. The simulations were done in hourly time interval and making use of the 250m of catchment characteristic maps. The calibration period is from 13/08/2002 till the end of 2005 while three years of remaining period is used for validation purpose. After that, different resolutions (i.e. 50m and 1 km) were employed to check if the simulation results remain consistently. Finally, a climate change effects on hydrological extremes assessment was carried out using the spatial information of 250m resolution.

5.1 Graphical evaluation of catchment behavior at the catchment outlet These spatial distributed versions of three hydrological models shown hereafter were re-calibrated in order to derive the appropriate runoff at the outlet station. Their total runoff hydrographs were intercompared to each other and to the lumped model and observation as in Figure 5-1, Figure 5-2 and Figure 5-3.

Figure 5-1: Total flow hydrograph simulated by lumped NAM model and spatial NAM model versus observation at the outlet station

WL2021R00_130_8

Final version

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Figure 5-2: Total flow hydrograph simulated by lumped PDM model and spatial PDM model versus observation at the outlet station

Figure 5-3: Total flow hydrograph simulated by lumped VHM model and spatial VHM model versus observation at the outlet station

Some extremes events were selected to zoom in the peaks of the hydrographs as in Figure 5-4, Figure 5-5 and Figure 5-6. Final version

WL2021R00_130_8

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Figure 5-4: Observed and simulated discharges downstream the Grote Nete catchment for the extreme event in January 2003 for lumped and distributed NAM, PDM and VHM

Figure 5-5: Observed and simulated discharges downstream the Grote Nete catchment for the extreme event in November 2004 for lumped and distributed NAM, PDM and VHM

Figure 5-6: Observed and simulated discharges downstream the Grote Nete catchment for the extreme event in January 2007 for lumped and distributed NAM, PDM and VHM

The simulated total flows were separated using the Chapman filter (see the previous section 3.1.1.2) into three components (This version of lumped PDM model itself produces only two output components). They were used to compare with filtered subflows of the observation as shown in Figure 5-7.

WL2021R00_130_8

Final version

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Figure 5-7: Observed and simulated subflows at the outlet station using the spatial version of PDM, NAM and VHM models

The Nash-Sutcliffe coefficients were calculated for the spatial model results and are listed as in Table 14. Table 14: NSE statistical performance values downstream the Grote Nete catchment for the calibration period (2003-2005) and validation period (2006-2008) using different models

Period

Spatial NAM

Spatial PDM

Spatial VHM

2003-2008 (total)

0.63

0.74

0.52

2003–2005 (calibration)

0.62

0.77

0.56

2006-2008 (validation)

0.66

0.69

0.40

The cumulative runoff volumes of the rainfall runoffs were also plotted for both calibration and validation periods in Figure 5-8. It is seen that in the calibration period, the spatial PDM model slightly overestimated the cumulative runoffs and cumulative baseflow while the spatial NAM model in the opposite underestimated those volumes. Among three models, the spatial VHM model was produces the most accurate prediction.

Final version

WL2021R00_130_8

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Figure 5-8: Observed and simulated cumulative discharge and baseflow at the catchment outlet using the spatial version of the PDM, NAM and VHM models

For the peak flow extremes evaluation, the different models show close peak flow distributions (see Figure 5-9). Also the scatterplots of peak flows (Figure 5-10) after Box-Cox transformation show similar model performance.

Figure 5-9: Peak flows versus return period at the outlet station for the spatial version of NAM, PDM and VHM models

WL2021R00_130_8

Final version

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Figure 5-10: Simulated versus observed nearly independent peak flows at the outlet station after Box-Cox transformation (λ=0.25) for the spatial version of NAM, PDM and VHM models

For the analysis of the low flow extremes (see Figure 5-11 and Figure 5-12), the VHM model shows systematic underestimations while the NAM model has slight systematic overestimations and the PDM model gave the best fit to the empirical values. Because of the lower number of low flow extremes, the analysis is obviously less clear/more uncertain.

Figure 5-11: Low flow extremes versus return period at the outlet station for the spatial version of NAM, PDM and VHM models Final version

WL2021R00_130_8

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Figure 5-12: Simulated versus observed nearly independent low flow extremes at the outlet station after Box-Cox transformation (λ=0.25) for the spatial version of NAM, PDM and VHM models

5.2 Spatial validation using internal gauging stations The recorded series at three internal gauging stations in the Grote Nete catchment were collected in order to validate the spatial variation of the runoffs generated using the spatial version of the NAM, PDM and VHM models. The selected stations are shown in Table 15. Table 15: Discharge stations in the Grote Nete catchment

Location

Station ID

Varendonk (outlet station)

0761-1066

190430

198098

Meerhout

L10_078

199681

203396

Vorst

L10_086

194015

196486

Tessenderlo

L10_087

200840

197124

Period 29/12/1984 to 31/12/2010 25/06/1997 to 05/04/2013 14/05/1986 to 18/04/2013 26/11/1998 to 25/09/2013

Source HIC

Discharge

VMM

Discharge

VMM

Discharge

VMM

Discharge

These series were truncated into the period of simulation. However, they contained many gaps of missing data so that only the overlapping parts were used. The similar comparative plots were used as 68

WL2021R00_130_8

Final version

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

in the validation of the runoffs at the catchment outlet. Moreover, a “rescaled flow” from the outlet records was added to compare with the spatial results given by the distributed models. This flow was calculated from the total discharge at the outlet station (Varendonk) by using the proportion in area of the sub-watershed to the whole catchment and this might help to evaluate the internal model performance the distributed models.

Figure 5-13: Total discharge simulated by NAM, PDM and VHM spatial models at Meerhout (top), Tessenderlo (middle) and Vorst (bottom) stations

Final version

WL2021R00_130_8

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Figure 5-13 shows the total discharge hydrographs at three stations. Except a strange summer record noticed at Meerhout station, three spatial models in general captured the dynamic variation of the runoff at three stations. However, no model was able to capture the magnitude of the peak events although they predicted the quite similar response of the sub-catchments to the meteorological conditions. Only at Vorst station, the models produced a close match to the observed values. This is illustrated more clearly in the next figures (Figure 5-14 and Figure 5-15).

Figure 5-14: Peak flows versus return period at Meerhout (top), Tessenderlo (middle) and Vorst (bottom) stations for NAM, PDM and VHM spatial models

WL2021R00_130_8

Final version

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Figure 5-15: Simulated versus observed nearly independent hourly peak flows at Meerhout (top), Tessenderlo (middle) and Vorst (bottom) stations for NAM, PDM and VHM spatial models

Final version

WL2021R00_130_8

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Nevertheless, the cumulative volumes captured by three models were very close to that of observations (see Figure 5-16). At Meerhout station, the strange period of recording caused differences in the cumulative volume.

Figure 5-16: Observed and simulated cumulative discharge at Meerhout (top), Tessenderlo (middle) and Vorst (bottom) stations for NAM, PDM and VHM spatial models

WL2021R00_130_8

Final version

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Due to the limitation in observation, the Meerhout series was not used to extract the low flow extremes to compare to the observation. In general, the low flow extremes were underestimated by all three models, especially for the higher return periods. For the smaller period, the extremes distribution was well described. The VHM model gave the highest underestimation of low flow extremes, which might be caused by the calibration process (see Figure 5-17).

Figure 5-17: Low flow extremes versus return period at Tessenderlo (top) and Vorst (bottom) stations for NAM, PDM and VHM spatial models

The poor calibration resulted in the worst performance of VHM amongst the three models. This is also confirmed by the statistical goodness-of-fit (NSE), which is listed in Table 16 (The comparison for Meerhout station was skipped due to the limitation in recording at this station). Table 16: NSE performance of different models for Tessenderlo and Vorst station

Final version

Location

Spatial NAM

Spatial PDM

Spatial VHM

Area proportional rescaled flow

Tessenderlo

0.53

0.57

0.52

0.37

Vorst

0.71

0.70

0.57

0.63

WL2021R00_130_8

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

5.3 Spatial validation using groundwater observations Several groundwater piezometers were used to collect the data to validate the groundwater variation produced by three spatial versions of NAM, PDM and VHM models. Since these models originate from a conceptual model, they do not calculate groundwater head or groundwater depth explicitly but one can derive the approximated values thanks to the computed groundwater volumes using a linear reservoir transfer function. The inter-comparisons were done between the calculated groundwater depths and the observed ones at the field.

7-0209

1-0311

1-0154

1-0485

7-0195 Figure 5-18: Comparison between observed and calculated groundwater height at different piezometers

WL2021R00_130_8

Final version

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

The variations or seasonal fluctuations of the groundwater head given by the models matched well to those of the observations, except for piezometer number 1-0311, which overestimates the head during the years 2008 and 2009. However, note that the derived model results here were calculated from the lumped groundwater reservoir volumes without any knowledge about the aquifers and groundwater table, nor any geological information on the soil layers. Therefore, some further calibrations on the groundwater head (or groundwater table) or detailed implementation of the aquifers are necessary to carry out.

5.4 Climate change effects on hydrological extremes given by distributed models Three future climate scenarios were used in three spatial distributed models to investigate the impacts on hydrological extremes in the Grote Nete. The analysis was done only at the catchment outlet but the results were available at other locations in the catchment as well. Moreover, the impact results generated by the original conceptual lumped models (NAM, VHM and PDM) were also used to compare with the outputs from the spatial versions. It was seen that the spatial distributed models produced quite similar impacts by climate change to that given by the lumped models both for low flow and high flow extremes. However, the underestimations were still noticed in the analysis for low flow extremes given by VHM. It was again explained by the poor re-calibration of the distributed VHM model.

Figure 5-19: Hourly peak flows versus return period at the catchment outlet for the low, mean and high climate scenarios simulated with lumped and spatial versions of NAM, PDM and VHM models for 2050

Final version

WL2021R00_130_8

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Figure 5-20: Hourly peak flows versus return period at the catchment outlet for the low, mean and high climate scenarios simulated with lumped and spatial versions of NAM, PDM and VHM models for 2100

WL2021R00_130_8

Final version

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Figure 5-21: Hourly low flows versus return period at the catchment outlet for the low, mean and high climate scenarios simulated with lumped and spatial versions of NAM, PDM and VHM models for 2050

Final version

WL2021R00_130_8

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Figure 5-22: Hourly low flows versus return period at the catchment outlet for the low, mean and high climate scenarios simulated with lumped and spatial versions of NAM, PDM and VHM models for 2100

5.5 Spatial resolution effects on model performance In order to investigate the effect of spatial resolution of the performance of distributed models, two different resolutions of spatial maps were employed, i.e. 50m x 50m and 1km x 1km. The static maps required by the framework were prepared for three models without any change in configuration. No extra or further calibration were refined to skip the effect of recalibration out of consideration. The analysis was firstly done at the catchment outlet as illustrated in Figure 5-23, Figure 5-24, Figure 5-25 and Figure 5-26.

WL2021R00_130_8

Final version

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Figure 5-23: Total discharge simulated by NAM, PDM and VHM distributed models at outlet station using three different spatial resolutions

The extreme events in January 2003 and in January 2007 were selected to zoom in the peak of the hydrographs derived by three models using different resolutions.

Final version

WL2021R00_130_8

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Figure 5-24: Observed and simulated hydrographs of the extreme events in January 2003 (top) and in January 2007 (bottom) by three models using 50m, 250m and 1km resolution maps

Figure 5-23 and Figure 5-24 show that three resolutions resulted in relatively similar results in hydrographs at the outlet station. Also, the cumulative volume simulated by the distributed models were closely matched with that from the lumped models for all three based structures (Figure 5-25). Therefore, the change in using various resolutions (while keeps the original system configuration and model parameters) does not make significant differences in the dynamic fluctuation of the hydrographs and total discharge volumes.

NAM

PDM

VHM Figure 5-25: Observed (in back line) and simulated cumulative discharge given by three models using 50m (in red line), 250m (in blue line) and 1km (in green line) resolution maps

WL2021R00_130_8

Final version

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

However, for high flow extremes distribution, the finest resolution map of 50m produced overestimation in all three models (Figure 5-25). The 250m and 1km maps gave a close estimation to each other and to the observation. While the hydrographs and total cumulative volumes were closely preserved, this is not the case for the quick flows. This may be caused by the differences in stream system definition where the configuration was kept as in the original implementation. This will require further investigation.

Figure 5-26: Observed and simulated peak distribution given by three models using 50m, 250m and 1km resolution maps

Three internal stations (see section 5.2) were also used to validate the performance of the distributed models with the use of different spatial resolutions. Due to a large amount of comparisons, only the distributions of the high and low flow extremes, which are mainly able to present the differences caused by the map resolutions (as were seen previously in the validation of the gauging station at the outlet of the catchment) are plotted.

Final version

WL2021R00_130_8

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Figure 5-27: Observed and simulated peak distributions for Meerhout (top), Tessenderlo (middle) and Vorst (bottom) stations for three distributed models using 50m, 250m and 1km resolution maps

WL2021R00_130_8

Final version

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Figure 5-28: Observed and simulated peak distributions for Tessenderlo (top) and Vorst (bottom) stations for three distributed models using 50m, 250m and 1km resolution maps

The performance of three models at the Tessenderlo station indicates the similarities with that at the outlet station when the finest resolution causes overestimation in peak distribution for all three models. However, the Vorst validation shows the opposite situation when the coarser resolution produced higher peaks for the same return period. But it may be caused by the issue of difference in co-ordination of the Vorst station when switching from a resolution to another. The collecting location may receive more water than it should in reality (or than in the fine resolution maps). Nevertheless, the simulation results at internal stations for three spatial resolutions were quite close to each other for all distributed models, especially when no extra calibration was done.

Final version

WL2021R00_130_8

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

6 Conclusions Having an overall look at the results for the Grote Nete catchment simulated by the framework using three model based structures, the following general conclusions regarding the concept of the implemented framework can be drawn:

The generalized model framework enabled to work well with the three different model structures implemented. Other model structures may be converted to the generalized structure to offer more options for the users. Also, additional hydrological components may be added such as a snow melting component, a capillary rise component, etc.. Additionally, some components may be modified to increase the physically based detail, e.g. applying a complex groundwater storage (component) instead of a lumped groundwater reservoir.

The spatial linking is now modified to describe the interaction between neighbouring grid cells in the catchment and between other sub-flow components. But some hydrological processes need to be investigated such as the “return flow” or may be further advanced such as the infiltration process. Also the connection between river cells and between river cells and the outlet station requires further study.

The model calibration method with the step-wise and transparent approach was applied for several catchments in Belgium using lumped models. The disaggregation from lumped model parameters into distributed ones or from one spatial resolution to another shows very promising results, but may be further advanced. For example, other approaches to disaggregate the maximum soil moisture storage, the runoff coefficients, the baseflow recession… can be tested in order to conclude on the most accurate one.

The evaluation of model performance for the Grote Nete case study shows very promising results. The lumped versions of all models implemented in the framework produced identical results as the individual stand-alone models. The distributed versions, with a slight extra calibration, generated results with similar accuracy to that of the lumped models at the outlet of the catchment. The distributed models were able to simulate proper discharges at various internal stations when comparing with observations and with the results from the detailed physically based distributed model MIKE-SHE and the semi-distributed model WetSpa (see report of Vansteenkiste, 2012a). Moreover, the climate change impacts on the hydrological extremes in the Grote Nete catchment were similar for the different models and similar to the ones obtained before by Vansteenkiste (2012b).

However, as mentioned above, there is room for further improvements or modifications to the implemented framework. Other catchments need to be studied as well.

WL2021R00_130_8

Final version

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Bibliography Ajami, N.K., H. Gupta, T. Wagener, S. Sorooshian, 2004. Calibration of a semi-distributed hydrologic model for streamflow estimation along a river system. Journal of Hydrology, 298, 112–135. Andréassian, V., Le Moine, N., Perrin, C., Ramos, M.H., Oudin, L., Mathevet, T., Lerat, J., 2012. All that glitters is not gold: the case of hydrological models’ calibration. Hydrological Processes, 26(14), 2206-2210. Apip, Sayama, T., Tachikawa, Y., and Takara, K., 2012. Spatial lumping of a distributed rainfall sediment-runoff model and its effective lumping scale. Hydrological Processes, 26, 855–871. Arnold, J.G., Allen, P.M., 1999. Automated methods for estimating base flow and groundwater recharge from stream flow records. Journal of the American Water Resources Association, 35(2), 411424. Bae, D.-H., Jung, I.-W., Lettenmaier, D.P, 2011. Hydrologic uncertainties in climate change from IPCC AR4 GCM simulations of the Chungju Basin, Korea. Journal of Hydrology 401, 90–105. Beck, M.B., Halfon, E., 1991. Uncertainty, identifiability and the propagation of prediction errors: a case study of Lake Ontario. Journal of Forecasting 10 135-162. Bell, V.A., Moore, R.J., 1998a. A grid-based distributed flood forecasting model for use with weather radar data: Part 1. Formulation. Hydrology and Earth System Sciences, 2, 265–281. Bell, V.A., Moore, R.J., 1998b. A grid-based distributed flood forecasting model for use with weather radar data: Part 2. Case studies. Hydrology and Earth System Sciences, 2, 283-298. Bell, V.A., Kay, A.L., Jones, R.G., Moore, R.J., 2012. Use of a grid-based hydrological model and regional climate model outputs to assess changing flood risk. International Journal of Climatology, 27, 16571671. Bekele, E.G., Nicklow, J.W., 2007. Multi-objective automatic calibration of SWAT using NSGA-II. Journal of Hydrology, 341, 3–4, 165–176. Beven, K.J., 1989. Changing ideas in hydrology: the case of physically-based models, Journal of Hydrology, 105, 157-172. Beven, K.J., 1993. Prophesy, reality and uncertainty in distributed hydrological modelling, Advances in Water Resources, 16, 41-51. Bormann, H., 2011. Treating an artificial catchment as ungauged: Increasing the plausibility of an uncalibrated, process-based SVAT scheme by using additional soft and hard data. Physics and Chemistry of the Earth, 36, 615-629. Box, G.E.P., Cox, D.R., 1964. An analysis of transformations, Journal of the Royal Statistical Society, 26, 211-243, discussion 244-252. Boyle, D.P., Gupta, H.V., Sorooshian, S., 2000. Towards improved calibration of hydrological models: Combining the strengths of manual and automatic methods, Water Resources Research, 36(12), 36633674.

Final version

WL2021R00_130_8

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Breuer, L., Huisman, J.A., Willems, P., Bormann, H., Bronstert, A., Croke, B.F.W., Frede, H.-G., Gräff, T., Hubrechts, L., Jakeman, A.J., Kite, G., Lanini, J., Leavesley, G., Lettenmaier, D.P., Lindström, G., Seibert, J., Sivapalan, M., Viney, N.R., 2009. Assessing the impact of land use change on hydrology by ensemble modelling (LUCHEM) I: Model intercomparison with current land use. Advances in Water Resources, 32, 2, 129-146. Brocca, L., Melone, F., Moramarco, T., Wagner, W., Naeimi, V., Bartalis, Z., Hasenauer, S., 2010. Improving runoff prediction through the assimilation of the ASCAT soil moisture product. Hydrology and Earth System Sciences, 14, 1881-1893. Brocca, L., Melone, F., Moramarco, T., 2011. Distributed rainfall-runoff modelling for flood frequency estimation and flood forecasting. Hydrological Processes, 25(8), 2801-2813. Brocca, L., Moramarco, T., Melone, F., Wagner, W., Hasenauer, S., Hahn, S., 2012a. Assimilation of surface and root-zone ASCAT soil moisture products into rainfall-runoff modelling. IEEE Transactions on Geoscience and Remote Sensing, 50(7), 2542-2555. Brocca, L., Melone, F., Moramarco, T., Zucco, G., Wagner, W., 2012b. Soil moisture assimilation into rainfall-runoff modelling: which is the impact of the model structure on runoff predictions? Geophysical Research Abstracts, 14, EGU2012-2307. Cabus, P., De Jongh, I., 2007. Hydrologisch modelleren op onbemeten stroomgebieden: Een eenvoudig modelconcept op basis van PDM en waargenomen maximale afvoercoëfficiënten. Tijdschrift Water, Thema Congres Watersysteemkennis 2006/2007: Modellen voor integraal waterbeheer in Vlaanderen, vol. 32, 55-59. Chapman, T., 1991. Comment on “Evaluation of automated techniques for base flow and recession analyses” by R. J. Nathan and T. A. McMahon. Water Resources Research, 27(7), 1783-1784. Chapman, T, 1999. A comparison of algorithms for stream flow recession and baseflow separation. Hydrological Processes, 13, 701-714. Chapman, T, 1999. A comparison of algorithms for stream flow recession and baseflow separation. Hydrological Processes 13 701-714. Chow, V.T., 1964. Handbook of applied hydrology; a compendium of water-resources technology. McGraw-Hill, New York, 1494 p. Chow, V.T., Maidment, D.R., Mays, L.W., 1988. Applied hydrology. McGraw-Hill, New York. Claps, P., Laio, F., 2003. Can continuous streamflow data support flood frequency analysis? An alternative to the partial duration series approach. Water Resources Research, 39(8), 1216. Clark, M.P., Kavetski, D., Fenicia, F., 2011. Pursuing the method of multiple working hypotheses for hydrological modelling. Water Resources Research, 47, W09301. Cosby, B.J., Hornberger, G.M., Clapp, R.B., Ginn, T.R., 1984. A statistical exploration of the relationships of soil moisture characteristics to the physical properties of soils. Water Resources Research, 20, 682-690. Crawford, N.H., Linsley, R.K., 1966. Digital simulation in hydrology: Stanford watershed model IV. Technical report 39, Department of Civil Engineering, Stanford University, Stanford, USA. Das, T., Bárdossy, A., Zehe, E., He, Y., 2008. Comparison of conceptual model performance using different representations of spatial variability. Journal of Hydrology, 356, 106-118.

WL2021R00_130_8

Final version

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Dawson, C.W., Abrahart, R.J., See, L.M., 2007. HydroTest: A web-based toolbox of evaluation metrics for the standardised assessment of hydrological forecasts. Environmental Modelling and Software 22 (7) 1034-1052. De Smedt, F., 1999. Two- and three dimensional flow of groundwater. Chapter 3 in: The Handbook of Groundwater Engineering (Ed. J.W. Delleur), CRC Press, 3.1-3.27. De Smedt, F., 2006. Two- and three dimensional flow of groundwater”, Chapter 4 in: The Handbook of Groundwater Engineering (Ed. J.W. Delleur) - second edition, CRC Press, 4.1-4.36. DHI, 2012. MIKE SHE Reference Manual. DHI, Water and Environment, Hørsholm, Denmark. Dondeyne, S., Bouhon, A., Legrain, X., Baert, G., Vancampenhout, K., Ampe, C., Cools, N., Engels, P., Langohr, R., Van Ranst, E., Chapelle, J., Deckers, J., 2012. 4th International Eurosoil Congress, Bari, Italië. Eckhardt, K., 2005. How to construct recursive digital filters for baseflow separation. Hydrological Processes, 19, 507-515. Engman, E.T., 1986. Roughness Coefficients for Routing Surface Runoff. Journal of Irrigation and Drainage Engineering - ASCE, 112(1), 39-53. Fenicia, F., Savenije, H.H.G., Matgen, P., Pfister, L., 2006. Is the groundwater reservoir linear? Learning from data in hydrological modelling. Hydrology and Earth System Sciences, 10, 139-150. Fenicia, F., Solomatine, D.P., Savenije, H.H.G., Matgen, P., 2007. Soft combination of local models in a multi-objective framework. Hydrology and Earth System Sciences, 11, 1797-1809. Ferket, B.V.A., Samain, B., Pauwels, V.R.N., 2010. Internal validation of conceptual rainfall-runoff models using baseflow separation. Journal of Hydrology, 381, 158-173. Freer, J., Beven, K., Ambroise, B., 1996. Bayesian estimation of uncertainty in runoff prediction and the value of data: An application of the GLUE approach, Water Resources Research, 32(7), 2161-2173. Ghavidelfar, S., Alvankar, S., Razmkhah, A., 2011. Comparison of the lumped and quasi-distributed Clark runoff models in simulating flood hydrographs on a semi-arid watershed. Water Resources Management, 25(6), 1775-1790. Gong, L., Widén-Nilsson, E., Halldin, S., Xu, C-Y., 2009. Large-scale runoff routing with an aggregated network-response function. Journal of Hydrology, 368, 237-250. Gosling, S., Taylor, R.G., Arnell,N., Todd, M.C., 2011. A comparative analysis of projected impact of climate change on river runoff from global and catchment-scale hydrological models. Hydrology and Earth System Sciences, 15, 1, 279-294. Graham, D.N., Butts, M.B., 2005. Flexible, integrated watershed modelling with MIKE SHE. In: Watershed Models (Eds. V.P. Singh & D.K. Frevert), CRC Press, 245-272. Gullentops et al., 2013. Next Generation Tool for Flexible Hydrological Modelling - Literature review of hydrological models. Tussentijds rapport, IMDC & KU Leuven voor project L 2012 T 0001 X Perceel 2 Dijle / Vlaamse Milieumaatschappij – Afdeling Operationeel Waterbeheer, augustus 2013. Gupta, V.K., Sorooshian, S., 1983. Uniqueness and observability of conceptual rainfall-runoff model parameters: the percolation process examined, Water Resources Research, 19(1), 269-276.

Final version

WL2021R00_130_8

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Gupta, H.V., Kling, H., Yilmaz, K.K., Martinez, G.F., 2009. Decomposition of the mean squared error and NSE performance criteria: Implications for improving hydrological modelling. Journal of Hydrology, 377(1-2), 80-91. Habets, F., Boé, J., Déqué, M., Ducharne, A., Gascoin, S., Hachour, A., Martin, E., Pagé, C., Sauquet, E., Terray, L., Thiéry, D., Oudin, L., Viennot, P., 2013. Impact of climate change on the hydrogeology of two basins in Northern France. Climatic Change, 121(4), 771-785. Harremoës, P., 2003. The need to account for uncertainty in public decision making related to technological change, Integrated Assessment, 4(1), 18-25. HEC, 2000. Hydrologic Modelling System HEC-HMS. Technical reference manual, US Army Corps of Engineers. Horn, D.L., Schwab, G.O., 1963. Evaluation of rational runoff coefficients for small agricultural watersheds. ASAE Trans. 6(3), 195-198, 201. Huisman, J.A., Breuer, L., Bormann, H., Bronstert, A., Croke, B.F.W., Frede, H., Gräff, T., Hubrechts, L., Jakeman, A.J., Kite, G.W., Lanini, J., Leavesley, G., Lettenmaier, D.P., Lindström, G., Seibert, J., Sivapalan, M., Viney, N.R., Willems, P., 2009. Assessing the impact of land use change on hydrology by ensemble modelling (LUCHEM) III: Scenario analysis. Advances in Water Resources, 32(2), 159-170. Jain, A., Srinivasulu, S., 2004. Development of effective and efficient rainfall-runoff models using integration of deterministic, real-coded genetic algorithms and artificial neural network techniques. Water Resources Research 40 doi:10.1029/2003WR002355. Jakeman, A.J., Littlewood, I.G., Whitehead, P.G., 1990. Computation of the instantaneous unit hydrograph and identifiable component flows with application to two small upland catchments, Journal of Hydrology, 117, 275–300. Jakeman, A.J., Letcher, R.A., Norton, J.P., 2006. Ten iterative steps in development and evaluation of environmental models. Environmental Modelling & Software 21 602-614. Jakeman, A.J., Hornberger, G.M., 1993. How much complexity is warranted in a rainfall-runoff model? Water Resources Research, 29(8), 2637-2649. Johnson, M.S., Coon, W.F., Mehta, V.K., Steenhuis, T.S., Brooks, E.S., Boll, J., 2003. Application of two hydrological models with different runoff mechanisms to a hill slope dominated watershed in the northeastern US: a comparison of HSPF and SMR. Journal of Hydrology 284 57-76. Kavetski, D., Fenicia, F., Clark, M.P., 2011. Impact of temporal data resolution on parameter inference and model identification in conceptual hydrological modelling: Insights from an experimental catchment, Water Resources Research, 47, W05501. Kelly, K.S., Krzysztofowicz, R., 1997. A bivariate meta-Gaussian density for use in hydrology, Stochastic Hydrology and Hydraulics, 11, 17-31. Klemeš, V., 1983. Conceptualization and scale in hydrology. Journal of Hydrology, 65, 1-23. Kokkonen, T.S., Jakeman, A.J., 2001. A comparison of metric and conceptual approaches in rainfallrunoff modelling and its implications. Water Resources Research 37(9) 2345-2352. Kristensen, K.J., Jensen, S.E., 1975. A model for estimating actual evapotranspiration from potential transpiration. Nordic Hydrology, 6, 70-88.

WL2021R00_130_8

Final version

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Laloy, E., Bielders, C.L., 2009. Modelling intercrop management impact on runoff and erosion in a continuous maize cropping system: Part II. Model Pareto multi-objective calibration and long-term scenario analysis using disaggregated rainfall. European Journal of Soil Science, 60, 1022-1037. Lang, M., Ouarda, T.B.M.J., Bobée, B., 1999. Towards operational guidelines for over-threshold modelling. Journal of Hydrology, 225, 103-117. Lee, H., Zehe, E., Sivapalan, M., 2007. Predictions of rainfall-runoff response and soil moisture dynamics in a microscale catchment using the CREW model. Hydrology and Earth System Sciences, 11, 819-849. Lees, M.J., Price, N., Wheather, H.S., Peach, D., 1998. A rainfall-runoff simulation model for South Galway, Ireland, In: Hydrology in a changing environment, Ed. H. Wheater and C. Kirby, John Wiley & Sons, vol. III, 93-104. Legates, D.R., McCabe, G.J., 1999. Evaluating the use of “goodness-of-fit” measures in hydrologic and hydroclimate model validation. Water Resources Research 35 (1) 233-241. Li, H., Zhang, Q., 2009. Multiobjective optimization problems with complicated Pareto sets, MOEA/D and NSGA-II. IEEE Transactions on Evolutionary Computation, 13(2), 284-302. Li, H., Beldring, S., Xu, C-Y., 2013. Implementation and testing of routing algorithms in the distributed Hydrologiska Byråns Vattenbalansavdelning model for mountainous catchments. Hydrology Research, in press. Liu, Y.B., De Smedt, F., 2004. WetSpa extension, documentation and user manual, Department of Hydrology and Hydraulic Engineering, Vrije Universiteit Brussel. Liu, T., Willems, P., Pan, X.L., Bao, A.M., Chen, X., Veroustraete, F., Dong, Q.H., 2012. Climate change impact on water resource extremes in a headwater region of the Tarim basin in China. Hydrol. Earth Syst. Sci., 15, 3511-3527. Lobligeois, F., Andréassian, V., Perrin, C., Tabary, P., Loumagne, C., 2013. When does higher spatial resolution rainfall information improve streamflow simulation? An evaluation on 3620 flood events, Hydrology and Earth System Sciences Discussions, 10, 12485-12536. Ludwig, R., May, I., Turcotte, R., Vescovi, L., Braun, M., Cyr, J.-F., Fortin, L.-G., Chaumont, D., Biner, S., Chartier, I., Caya, D., Mauser, W., 2009. The role of hydrological model complexity and uncertainty in climate change impact assessment. Advances in Geosciences 21, 63–71. Madsen, H., Rasmussen, P.F., Rosbjerg, D., 1997. Comparison of annual maximum series and partial duration series methods for modelling extreme hydrologic events. 1. At-site modelling. Water Resources Research, 33(4), 747-757. Madsen, H., 2000. Automatic calibration of a conceptual rainfall-runoff model using multiple objectives, Journal of Hydrology, 235, 276-288. Madsen, H., 2003. Parameter estimation in distributed hydrological catchment modelling using automatic calibration with multiple objectives. Advances in Water Resources, 26(2), 205-216. Mantovan, P., Todini, E., 2006. Hydrological forecasting uncertainty assessment: incoherence of the GLUE methodology. Journal of Hydrology, 330, 368–381. Maurer, E.P., Brekke, L.D., Pruitt, T., 2010. Contrasting lumped and distributed hydrology models for estimating climate change impacts on California watersheds. Journal of the American Water Resources Association, 46 (5), 1024–1035. Final version

WL2021R00_130_8

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Meert, P., Willems, P., 2013. Evaluatie van berekeningsmethoden voor het bepalen van de benodigde buffercapaciteit van kleinschalige opvangsystemen in het kader van erosiebestrijding. KU Leuven Afdeling Hydraulica, voor Vlaamse Overheid – Departement LNE – Afdeling Land en Bodembescherming, Ondergrond, Natuurlijke Rijkdommen, augustus 2013. Mernild, S.H., Liston, G.E., 2012. Greenland freshwater runoff.Part II: Distribution and trends, 1960– 2010. Journal of Climate, 25(17), 6015–6035. Montanari, A., Brath, A., 2004. A stochastic approach for assessing the uncertainty of rainfall-runoff simulations, Water Resources Research, 40, W01106, doi:10.1029/2003WR002540. Moore, R. J., Bell, V. A., Cole, S. J., and Jones, D. A., 2007. Rainfall-runoff and other modelling for ungauged/low-benefit locations: Operational Guidelines, Environment Agency, R&D Technical Report SC030227/SR2 Muleta, M.K., Nicklow, J.W., 2005. Sensitivity and uncertainty analysis coupled with automatic calibration for a distributed watershed model. Journal of Hydrology, 306, 127–145. Najafi, M.R., Moradkhani, H., Jung, I.W., 2011. Assessing the uncertainties of hydrologic model selection in climate change impact studies. Hydrological Processes, 25 (18), 2814–2826. Nathan, R.J., McMahon, T.A., 1990. Evaluation of automated techniques for base flow and recession analyses. Water Resources Research, 26(7), 1465-1473. Neumann, M.B., Gujer, W., 2008. Underestimation of uncertainty in statistical regression of environmental models: influence of model structure uncertainty. Environmental Science and Technology, 42(11), 4037–4043. Oudin, L., Andréassian, V., Mathevet, T., Perrin, C., Michel, C., 2006. Dynamic averaging of rainfallrunoff model simulations from complementary model parameterizations. Water Resources Research 42 W07410 doi:10.1029/2005WR004636. Penna, D., Tromp-van Meerveld, H.J., Gobbi, A., Borga, M., Dalla Fontana, G., 2011. The influence of soil moisture on threshold runoff generation processes in an alpine headwater catchment. Hydrology and Earth System Sciences, 15, 689-702. Perrin, C., Michel, C., Andréassian, V., 2001. Does a large number of parameters enhance model performance? Comparative assessment of common catchment model structures on 429 catchments, Journal of Hydrology, 242, 275–301. Pokhrel, P., Gupta, H.V., 2011. On the ability to infer spatial catchment variability using streamflow hydrographs. Water Resources Research, 47(8), W08534. Rawls, W.J., Brakensiek, D.L., 1982. Estimating soil water retention from soil properties. Journal of Irrigation and Drainage Engineering (ASCE), 108, 166-171. Reed, S., Koren, V., Smith, M. B., Zhang, Z., Moreda, F., Seo, D., Dmipparticipants, A., 2004. Overall distributed model intercomparison project results, Journal of Hydrology, 298, 27–60. Robson, B.J., Hamilton, D.P., Webster, I.T., Chan, T., 2008. Ten steps applied to development and evaluation of process-based biogeochemical models of estuaries. Environmental Modelling and Software 23 369-384. Ruelland, D., Ardoin-Bardin, S., Billen, G., Servat, E., 2008. Sensitivity of a lumped and semidistributed hydrological model to several methods of rain interpolation on a large basin in West Africa. Journal of Hydrology, 361, 96-117. 90

WL2021R00_130_8

Final version

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Savenije, H.H.G., 2009. The art of hydrology. Hydrology and Earth System Sciences, 13, 157-161. Schwab, G.O., Fangmeier, D.D., Elliot, W.J., Frevert, R.K., 1993. Soil and water conservation engineering. Fourth Edition, John Wiley and Sons, Inc., New York, 507 p. Seibert, J., McDonnell, J., 2000. Towards a better process representation of catchment hydrology in conceptual runoff modelling. IAHS Freiburg conference on Meso-scale modelling, Proceedings of the international workshop on Runoff generation and implications for river basin modelling, 9-12 October 2000. Sivakumar, B., 2008. Dominant processes concept, model simplification and classification framework in catchment hydrology. Stochastic Environmental Research and Risk Assessment, 22(6), 737-748. Sivapalan, M., Blöschl, G., Zhang, L., Vertessy, R., 2003. Downward approach to hydrological prediction. Hydrological Processes, 17, 2101-2111. Smith, M.B., Koren, V., Zhang, Z., Zhang, Y., Reed, S. M., Cui, Z., Moreda, F., Cosgrove, B. A., Mizukami, N., Anderson, E.A., 2012. Results of the DMIP 2 Oklahoma experiments. Journal of Hydrology, 418-419, 17-48. Sorooshian, S., 1981. Parameter Estimation of Rainfall-Runoff Models with Heteroscedastic Streamflow Errors - Noninformative Data Case", Journal of Hydrology, 52(1/2), 127-138, 1981. Sorooshian, S., Dracup, J.A., 1980. Stochastic parameter estimation procedures for hydrologic rainfallrunoff models: correlated and heteroscedastic error cases. Water Resources Research, 29, 1185-1194. Tavakoli, M., 2012. Simulating soil moisture and climate change impacts in a watershed through application of the distributed hydrological WetSpa model. PhD dissertation, Faculty of Engineering, Free University of Brussels, December 2012. Tavakoli, M., De Smedt, F., 2012. Impact of climate change on streamflow and soil moisture in the Vermilion basin, Illinois. Journal of Hydrologic Engineering, 17(10), 1059-1070. Tavakoli, M., De Smedt, F., 2013. Validation of soil moisture simulation with a distributed hydrologic model (WetSpa). Environmental Earth Sciences, 69(3), 739-747. Terink, W., Droogers, P., Immerzeel, W., van den Eertwegh, G., 2012. SPHY – Een hydrologisch model gericht op de berekening van bodemvocht in de wortelzone en de actuele verdamping. Report FutureWater 115, Wageningen, August 2012. Tourbier, J., Westmacott, R., 1974. Water resources protection measures in land development – a handbook. National Technical Information Service, U.S. Department of Commerce, Delaware, Newak, Delaware, 237 p. Tran, Q., Pereira, F., Nossent, J., Willems, P., Mostaert, F., 2013. Effect of climate change on the hydrological regime of navigable water courses in Belgium: Subreport 7 – Development of a framework for flexible hydrological modelling. Version 1_0. WL Rapporten, 00_130. Flanders Hydraulics Research: Antwerp, Belgium. Tran, A.P., 2014. Full-wave inversion of near-field GPR data for hydrogeophysical characterization of soil. PhD manuscript, Université catholique de Louvain, Faculté d’Ingénierie Biologique, Agronomique et Environnementale. Tran, A.P., Vanclooster, M., Zupanski, M., Lambot, S., 2014. Joint estimation of soil moisture profile and hydraulic parameters by ground-penetrating radar data assimilation with maximum likelihood ensemble filter. Accepted for publication in Water Resources Research. Final version

WL2021R00_130_8

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Uhlenbrock, S., Seibert, J., Leibundgut, C., Rohde, A., 1999. Prediction uncertainty of conceptual rainfall-runoff models caused by problems in identifying model parameters and structures, Hydrological Sciences Bulletin, 44(5), 779-797. USWRC, 1976. Guidelines for determining flood flow frequency. United States Water Resources Council, Bull.17, Hydrol. Comm. Washington DC, 73 p. Vandekerckhove, L., Swerts, M., Leyman, N., Mennens, K., Neven, H., Desmet, J., De Vrieze, M., 2001. Code van goede praktijk voor het opmaken van een gemeentelijk erosiebestrijdingsplan. Rapport D/2001/3241/252, Vlaamse Overheid, AMINAL, Afdeling Land, 83 p. Van Hoey, S., Vansteenkiste T., Pereira F., Nopens, I., Seuntjens, P., Willems, P., Mostaert, F., 2012. Effect of climate change on the hydrological regime of navigable water courses in Belgium: Subreport 4 - Flexible model structures and ensemble evaluation. WL Rapporten, project 706_18, Waterbouwkundig Laboratorium, KU Leuven, VITO & U.Gent, September 2012, 29 p. Van Ranst, E., Sys, C., 2000. Eenduidige legende voor de digitale bodemkaart van Vlaanderen (schaal 1:20 000). Universiteit Gent, Laboratorium voor Bodemkunde. Vansteenkiste, T., Pereira, F., Willems, P., Vanneuville, W., Van Eerdenbrugh, K., Mostaert, F., 2011a. Effect van klimaatwijzigingen op afvoerdebieten in hoog- en laagwatersituaties en op de globale waterbeschikbaarheid: Deelrapport 4 - Impactanalyse met de ruimtelijke verdeelde modelcodes MIKE SHE en WetSpa - gevalstudie in Vlaanderen. WL Rapporten, 706_13a. Waterbouwkundig Laboratorium: Antwerpen, België. Vansteenkiste, T., Pereira, F., Willems, P., Vanneuville, W., Van Eerdenbrugh, K., Mostaert, F., 2011b. Effect of climate change on the hydrological regime of navigable water courses in Belgium: Subreport 1 - Literature review of the climate research in Belgium. WL Rapporten, project 706_18, Waterbouwkundig Laboratorium & KU Leuven, September 2011, 43 p. Vansteenkiste, T., Pereira, F., Willems, P., Mostaert, F., 2011c. Effect of climate change on the hydrological regime of navigable water courses in Belgium: Subreport 2 – Climate change impact analysis by conceptual models. WL Rapporten, project 706_18, Waterbouwkundig Laboratorium & KU Leuven, December 2011, 57 p. Vansteenkiste, T., Pereira, F., Willems, P., Mostaert, F., 2012a. Effect of climate change on the hydrological regime of navigable water courses in Belgium: Subreport 3 – Climate change impact analysis by distributed models. WL Rapporten, project 706_18, Waterbouwkundig Laboratorium & KU Leuven, February 2012, 71 p. Vansteenkiste, T., Pereira, F., Willems, P., Mostaert, F., 2012b. Effect of climate change on the hydrological regime of navigable water courses in Belgium: Subreport 5 – Comparison of lumped and distributed models in the climate change impact assessment. WL Rapporten, project 706_18, Waterbouwkundig Laboratorium & KU Leuven, August 2012, 33 p. Vansteenkiste, T., Pereira, F., Willems, P., Mostaert, F., 2012c. Effect of climate change on the hydrological regime of navigable water courses in Belgium: Subreport 6 - Urban expansion impact analysis by distributed models. WL Rapporten, project 706_18, Waterbouwkundig Laboratorium & KU Leuven, September 2012, 17 p. Vansteenkiste, T., Tavakoli, M., Van Steenbergen, N., De Smedt, F., Batelaan, O., Pereira, F., Willems, P., (2014), ‘Intercomparison of five lumped and distributed models for catchment runoff and extreme flow simulation’, Journal of Hydrology, 511C, 335-349.

WL2021R00_130_8

Final version

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Van Steenbergen, N., Willems, P., 2012. Method for testing the accuracy of rainfall–runoff models in predicting peak flow changes due to rainfall changes, in a climate changing context. Journal of Hydrology. 414-415, 425-434. Velázquez, J.A., Schmid, J., Ricard, S., Muerth, M.J., Gauvin St-Denis, B., Minville, M., Chaumont, D., Caya, D., Ludwig, R., Turcotte, R, 2012. An ensemble approach to assess hydrological models' contribution to uncertainties in the analysis of climate change impact on water resources. Hydrology and Earth System Sciences, 17, 565-578. Viessman, W., Lewis, G.L., 2003. Introduction to hydrology. Fifth Edition, Prentice-Hall. Viney, N.R., Bormann, H., Breuer, L., Bronstert, A., Croke, B.F.W., Frede, H., Gräff, T., Hubrechts, L., Huisman, J.A., Jakeman, A.J., Kite, G.W., Lanini, J., Leavesley, G., Lettenmaier, D.P., Lindström, G., Seibert, J., Sivapalan, M., Willems, P., 2009. Assessing the impact of land use change on hydrology by ensemble modelling (LUCHEM) II: Ensemble combinations and predictions. Advances in Water Resources, 32(2), 147-158. Vrebos, D., Vansteenkiste, T., Staes, J., Willems, P., Meire, P., 2013. Water displacement by sewer infrastructure in the Grote Nete catchment, Belgium, and its hydrological regime effects. Hydrology and Earth System Sciences Discussions, 10, 7425-7467. Vrugt, J., Gupta, H.V., Bastidas, L.A., Bouten, W., Sorooshian, S., 2003. Effective and efficient algorithm for multiobjective optimization of hydrologic models. Water Resources Research, 39, 1214, doi:1210.1029/2002WR001746. Vrugt, J.A., Diks, C.G.H., Gupta, H.V., Bouten, W., Verstraten, J.M., 2005. Improved treatment of uncertainty in hydrologic modelling: combining the strengths of global optimization and data assimilation. Water Resources Research, 41, W01017. Wagener, T., Wheater, H.S., 2002. A generic framework for the identification of parsimonious rainfallrunoff models, In: Integrated assessment and decision support, Proceedings of the 1st Biennial Meeting of the iEMSs (Ed. A.E. Rizzoli and A.J. Jakeman), 434-439. Wagener, T., Boyle, D.P., Lees, M.J., Wheater, H.S., Gupta, H.V., Sorooshian, S., 2001. A framework for the development and application of hydrological models, Hydrology and Earth System Sciences, 5(1), 13-26. Wagener, T., McIntyre, N., Lees, M.J., Wheather, H.S., Gupta, H.V., 2003. Towards reduced uncertainty in conceptual rainfall-runoff modelling: Dynamic identifiability analysis. Hydrological Processes, 17, 455-476. Welsh, W.D., 2008. Water balance modelling in Bowen, Queensland, and the ten iterative steps in model development and evaluation. Environmental Modelling and Software 23 (2) 195-205. Westerberg, I. K., Guerrero, J.-L., Younger, P. M., Beven, K. J.,. Seibert, J., Halldin, S., Freer, J. E., Xu, C.-Y., 2011. Calibration of hydrological models using flow-duration curves Hydrology and Earth System Sciences, 15, 2205–2227. Willems, P., 2000. Probabilistic immission modelling of receiving surface waters. Doctoral dissertation, Faculty of Engineering, KU Leuven. Willems, P., 2009a. A time series tool to support the multi-criteria performance evaluation of rainfallrunoff models. Environmental Modelling & Software, 24, 311–321.

Final version

WL2021R00_130_8

Effect of climate change on the hydrological regime of navigable water courses in Belgium Sub report 8 – Implementation and testing of a framework for flexible hydrological modelling

Willems, P., 2009b. Actualisatie en extrapolatie van hydrologische parameters in de nieuwe Code van Goede Praktijk voor het Ontwerp van Rioleringssystemen. KU Leuven – Afdeling Hydraulica voor VMM – Afdeling Operationeel Waterbeheer, september 2009, 79 p. Willems P., et al., 2002. Methodologie ter hydrologische modellering van Vlaamse hydrografische (sub)bekkens en de afleiding van maatgevende hydrogrammen. KU Leuven – Afdeling Hydraulica voor AMINAL – Afdeling Operationeel Waterbeheer, mei 2002, 309 p. Willems, P., 2013. Waterloopmodellering. Acco Leuven & Den Haag, 268 p., ISBN 978-90-334-9296-9, D/2013/0543/132, NUR 955. Willems, P., 2014. Parsimonious rainfall-runoff model construction supported by time series processing and validation of hydrological extremes – Part 1: Step-wise model-structure identification and calibration approach. Journal of Hydrology, 510, 578–590. Willems, P., Mora, D., Vansteenkiste, Th., Taye, M.T., Van Steenbergen, N., 2014a. Parsimonious rainfall-runoff model construction supported by time series processing and validation of hydrological extremes – Part 2: Intercomparison of models and calibration approaches. Journal of Hydrology, 510, 591–609. Willems, P., Ntegeka, V., Baguis, P., Roulin, E., 2010. Climate change impact on hydrological extremes along rivers and urban drainage systems. Final report for Belgian Science Policy Office, KU Leuven – Hydraulics Section & Royal Meteorological Institute of Belgium, December 2010, 110 p. Willems, P., Quan, T.Q., Van den Zegel, B., Akurut, M., Verhoest, N., Vanhaute, W-J., 2014b. Onderzoek naar het gebruik van neerslaggeneratoren voor het aanmaken van overstromingsrisicokaarten - 2. Deelproces hydrologie. Eindrapport, KU Leuven & U.Gent, voor TWOL onderzoeksproject L 2011 T 0005 / Vlaamse Milieumaatschappij – Afdeling Operationeel Waterbeheer, januari 2014, 160 p. Wheather, H.S., Jakeman, A.J., Beven, K.J., 1993. Progress and directions in rainfall-runoff modelling, In: Modelling change in environmental systems, Ed. A.J. Jakeman et al., John Wiley & Sons, 102-132. Xu, C-Y, 2001. Statistical analysis of a conceptual water balance model, methodology and case study, Water Resources Management, 15, 75-92. Young, P., 1986. Time series methods and recursive estimation in hydrological systems analysis. In: Craijenhoff, D.A., Moll, J.R., (Eds.), River flow modelling and forecasting 129-180. Reidel: Dordrecht. Young, P., Beven, K.J., 1994. Data-based mechanistic modelling and the rainfall-flow nonlinearity, Environmetrics, 5, 335-363. Young, P.C., Parkinson, S.D. Lees, M.J., 1996. Simplicity out of complexity in environmental systems: Occam’s Razor revisited, Journal of Applied Statistics, 23, 165-210. Yu, P.S., Yang, T.C., 2000. Fuzzy multi-objective function for rainfall-runoff model calibration, J. of Hydrol., 238, 1-14. Zhang, H., Huang, G.H, Wang, D., Zhang, X., 2011a. Multi-period calibration of a semi-distributed hydrological model based on hydroclimatic clustering. Advances in Water Resources, 34, 1292-1303. Zhang, Y., Wei, H., Nearing, M.A., 2011b. Effects of antecedent soil moisture on runoff modelling in small semiarid watersheds of southeastern Arizona. Hydrology and Earth System Sciences, 15, 31713179.

WL2021R00_130_8

Final version

DEPARTMENT MOBILITY & PUBLIC WORKS Flanders hydraulics Research Berchemlei 115, 2140 Antwerp T +32 (0)3 224 60 35 F +32 (0)3 224 60 36 waterbouwkundiglabo@vlaanderen.be www.flandershydraulicsresearch.be