Environmental Modelling and Software 161 (2023) 105628

Contents lists available at ScienceDirect

Environmental Modelling and Software

journal homepage: [Link]/locate/envsoft

A CUDA Fortran GPU-parallelised hydrodynamic tool for high-resolution

and long-term eco-hydraulic modelling
Marcos Sanz-Ramos a, *, David López-Gómez b, Ernest Bladé a, Danial Dehghan-Souraki a
a
Flumen Institute, Universitat Politècnica de Catalunya (UPC), International Center for Numerical Methods in Engineering (CIMNE), 08034, Barcelona, Spain
b
Centro de Estudios Hidrográficos, Centro de Estudios y Experimentación de Obras Públicas (CEDEX), 28005, Madrid, Spain

A R T I C L E I N F O A B S T R A C T

Handling Editor: Daniel P Ames Eco-hydraulic models are wide extended tools to assess physical habitat suitability on aquatic environments.
Currently, the application of these tools is limited to short river stretches and steady flow simulations. However,
Keywords: this limitation can be overcome with the application of a high-performance computing technique: graphics
Fish habitat processing unit (GPU) computing. R-Iber is a GPU-based hydrodynamic code parallelised in CUDA Fortran that,
Weight useable area
with the integration of a biological module, performs as an eco-hydraulic numerical tool. R-Iber was validated
Numerical modelling
and applied to real cases by using an optimised instream flow incremental methodology in long river reaches and
High performance computing
GPGPU long-term simulations. R-Iber reduces the computation time considerably, reaching speed-ups of two orders of
magnitude compared to traditional computing. R-Iber allows for overcoming the current limitations of the eco-
hydraulic tools with the simulation of high-resolution numerical models calculated in a reasonable computation
timeframe, which provides a better representation of the hydrodynamics and the physical habitat.

Software availability Alcázar, 2012; Wilkes et al., 2016) in all situations, not only during
extreme conditions.
• Name of tool: R-Iber Eco-hydraulics is a technique that analyses the effect of physical
• Developers: David López-Gómez, Marcos Sanz-Ramos, Ernest Bladé. environmental properties (flow depth, velocity, turbulence, tempera­
• Year first available: 2021. ture, substrate, etc.) on aquatic environments (Bovee, 1982), such as
• Hardware required: basic computer with a graphical power unit rivers. To characterise the physical habitat, several approaches have
(GPU) based on CUDA architecture (mainly, any NVIDIA GPU). been proposed based on hydrologic, geomorphologic or hydraulic
• Requirements: Windows OS x64. criteria (Tonina and Jorde, 2013). Particularly for hydraulics, a common
• Source Code Availability: the numerical tool is freely distributed method for assessing physical habitat is the use of normalized univariate
through [Link] curves, or suitability curves, that relate variables of the river environ­
• Data availability: the authors do not have permissions to share the ment with the inhabiting species’ use of the space.
data. Several numerical hydrodynamic tools have been developed or
• Cost: free enhanced by integrating the biological requirements of one or several
• Program languages: CUDA Fortran target species typically of fishes (Cassan et al., 2022; Hung et al., 2022;
Jowett, 2004; Meza Rodríguez et al., 2019; Nones, 2019; Sanz-Ramos
1. Introduction et al., 2019; Shim et al., 2020; Stamou et al., 2018; Steffler and Black­
burn, 2002), but not exclusively (Hamilton et al., 2015; Zohmann et al.,
Water resources and river management are no longer uniquely 2013). These tools are based on the solution of the one-dimensional
related to flood or drought scenarios. On the contrary, to properly study (1D), and most recently, on the two-dimensional (2D) Saint-Venant
the system from a holistic point of view, the relationship between hy­ equations or shallow water equations (SWE), which provide the evolu­
draulics and biota must be considered (Benjankar et al., 2018; Palau and tion of the water depth and the specific discharge or velocity of the

Abbreviations: SWE, shallow water equations; HPC, high-performance computing; GPGPU, general-purpose computing on graphics processing unit; GPU, graphics
processing unit; CPU, central processing unit; CUDA, compute unified device architecture; IFIM, instream flow incremental methodology; WUA, weight usable area.
* Corresponding author.
E-mail addresses: [Link]-ramos@[Link] (M. Sanz-Ramos), [Link]@[Link] (E. Bladé), [Link]@[Link] (D. Dehghan-Souraki).

[Link]
Received 11 August 2022; Received in revised form 9 January 2023; Accepted 10 January 2023
Available online 14 January 2023
1364-8152/© 2023 The Authors. Published by Elsevier Ltd. This is an open access article under the CC BY-NC-ND license ([Link]
nc-nd/4.0/).
M. Sanz-Ramos et al. Environmental Modelling and Software 161 (2023) 105628

water. Despite the current rapid increase in the use of three-dimensional

(3D) river flow modelling (Bermúdez et al., 2017, 2018; Meselhe et al., where h is the water depth, Ux and Uy are the x and y components of the
2012; Pisaturo et al., 2017), the computational cost of this method limits depth-averaged velocity, g is the gravitational acceleration, So,x = − ∂zb /
its application to very short river stretches and research purposes. ∂x and So,y = − ∂zb /∂y are the spatial variation of the bed elevation (zb )
√̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅
Apart from the simplifications applied to the generation of the suit­
on the x and y directions, Sf ,x = (n2 Ux U2x + U2y )/h4/3 and So,y =
ability curves (Boudreault et al., 2022; Palau et al., 2016), limitations of √̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅̅
2D-SWE-based eco-hydraulic models typically relate to computational (n2 Uy U2x + U2y )/h4/3 are the energy dissipation in the x and y di­
cost (Carlotto et al., 2021; Morales-Hernández et al., 2021). The time rections due to bed friction, which is computed using the Manning for­
required to simulate large river stretches, flow properties in natural mula, and νt is the turbulent viscosity.
conditions (e.g. long time series), and the use of a very fine domain The previous system of equations are the main part of the hydro­
discretisation with thousands or millions of calculation points (or mesh dynamics module of Iber (Bladé et al., 2014a, 2014b), a 2D numerical
elements) typically require high computational effort. tool that was initially developed for hydrodynamics and sediment
These limitations can be overcome with a good code design or/and transport modelling (Bladé et al., 2014a, 2014b, 2019b, 2019a; Sanz-­
with the application of high-performance computing (HPC) techniques. Ramos et al., 2020a). It solves the 2D-SWE on irregular geometries using
Currently, one of the most applied techniques is the parallelisation of the the conservative scheme based on the Finite Volume Method of Roe
computations by using general-purpose computing on graphics pro­ (1986), which consists of the Godunov method together with the Roe
cessing unit (GPGPU) instead of the traditional central processing unit Approximate Riemann Solver (Toro, 2009). The tool has been continu­
(CPU) computing (Buttinger-Kreuzhuber et al., 2022; Carlotto et al., ously enhanced, and the sequential version currently includes a series of
2021; Morales-Hernández et al., 2021). modules for different free surface flow processes, such as hydrological
A GPGPU technique commonly used in computational fluid dy­ processes (Cea and Bladé, 2015; Sanz-Ramos et al., 2020b, 2021),
namics (CFD) is that provided by the NVIDIA compute unified device pollutant propagation (Cea et al., 2016), large-wood transport (Ruiz-­
architecture (CUDA), which is currently available for the C++ and Villanueva et al., 2014), and physical habitat suitability assessment
Fortran programming languages (NVIDIA, 2022a). The benefit of using (Sanz-Ramos et al., 2019).
GPU computing is the significant reduction in computation time, with
speed-ups of one or two orders of magnitude compared to CPU
computing. Nevertheless, 2D-SWE-based eco-hydraulic numerical tools 2.2. Eco-hydraulics methods
that include an integrated GPU-parallelised hydrodynamic module and
physical habitat module are lacking. The combination of a hydrodynamic model and a biological model
This study aims to fill this gap with the presentation of a new GPU- can be oriented to develop eco-hydraulic tools (Nestler et al., 2016). One
based hydraulic code parallelised in CUDA Fortran, named R-Iber, and of the most extended biological models to assess physical habitat,
its performance as an eco-hydraulic numerical tool with the integration particular for freshwater fishes, is the use of suitability curves, a
of a physical habitat module. The code was first validated using two normalized relation between environmental variables (hydraulic, water
different benchmark tests for 2D-SWE-based models and a laboratory quality, geomorphological, etc.) and the use of the space of the target
experiment with a vertical slot fishway. Then, R-Iber was applied to two specie/stadium. The elemental habitat suitability of a physical variable
real cases to show not only the full integration of the hydrodynamic and ranges from 0 to 1, values being close to 1 those correspond to a high
physical habitat module for some fish species but also the benefits of degree of suitability. The multiplication of the elemental suitability
using a fully-integrated GPU-parallelised eco-hydraulic tool in high- (depth- and velocity-dependent) of each species/stadium was used as
resolution numerical models and long-term simulations. habitat criteria in the following sections, but other criteria can be used
without loss of applicability.
2. Materials and methods The Instream Flow Incremental Methodology (IFIM) is a technique
that uses the combination of hydrodynamic and biologic models to
2.1. Governing equations create a habitat or eco-hydraulic model. It consists of the assessment of
different constant discharges (Q) in a river stretch to obtain the Weight
Most existing eco-hydraulic simulation tools are based on the solu­ Useable Area (WUA) for each one (Bovee, 1982). With these Q-WUA
tion of the mass and momentum conservation equations that, when relations, it is possible to obtain the theoretical discharge that maximises
applied to a 2D framework, are named 2D Shallow Water Equations (2D- the WUA for the species and stadiums analysed, and which is usually
SWE). They are derived from Navier–Stokes equations through a time utilized to propose the environmental flows.
averaging to filtrate the turbulent fluctuations (Reynolds Averaged In 2D eco-hydraulic modelling, the use of the IFIM traditionally
Navier–Stokes equations [RANS]) and a depth averaging to obtain the implies the generation of as many models as the discharges to be eval­
final 2D equations. Two-dimensional SWE are a hyperbolic nonlinear uated. Being an analysis conducted under steady flow conditions, 2D-
system of three differential equations in partial derivatives (Toro, 2009), SWE-based models usually require simulating the model from the
which can be written as follows: beginning to an unknown maximum time of simulation. Additionally,
the time required to achieve each steady state of the series of discharges

∂h ∂hUx ∂hUy
+ + =0
∂t ∂x ∂y
( ) ( ) ( )
∂hUx ∂ h2 ∂( ) ( ) ∂ ∂U ∂ ∂U
+ hU 2x + g + hUx Uy = gh So,x − Sf ,x + νt h x + νt h y [1]
∂t ∂x 2 ∂y ∂x ∂x ∂y ∂y
( ) ( ) ( )
∂hUx ∂ ( ) ∂ h2 ( ) ∂ ∂U ∂ ∂U
+ hUx Uy + hU 2y + g = gh So,y − Sf ,y + νt h y + νt h y
∂t ∂y ∂y 2 ∂x ∂y ∂y ∂y

2
M. Sanz-Ramos et al. Environmental Modelling and Software 161 (2023) 105628

to evaluate is, in general, unknown. This time depends not only on the made in previous versions of Iber. Despite that, the simulation of models
geometry and length of the river, but also on the flow intensity. This of thousands or millions of elements, or long-time series, requires
issue is addressed with the definition of a series of constant discharges, several hours or even days of computational time with the sequential
or ‘stepped discharges’, as an inlet condition. The ‘stepped discharge’ version.
option allows for computing, in a unique model, a series of discharges R-Iber is a new code that implements a GPU parallelisation of the
that generates steady conditions considering a tolerance between the hydrodynamic (Bladé et al., 2014b) and habitat (Sanz-Ramos et al.,
inlet and the outlet. When the steady conditions of the target discharge 2019) modules of Iber based on CUDA Fortran architecture. An efficient
are achieved, the model automatically steps to the next discharge. The CUDA programming requires not only a hierarchical organisation of the
hydrodynamic conditions at the end of a step are used as initial condi­ computing units but also the programming of each thread block to
tions for the next step. Thus, this option uses the time strictly necessary efficiently process data in parallel (Zhang and Jia, 2013). To achieve
to reach the steady flow conditions of any of the constant discharges this, the original code of Iber was re-organised and re-written to follow
defined by the user, without having to know the time to reach the steady the main structure of CUDA programming, that is: transfer data from
flow conditions for any discharge and reducing the computational time CPU memory (host) to GPU memory (device); run parallel computing
to a minimum. through ‘kernels’; and transfer data back to the host from the device.
More details of the structure of the code are described in Appendix A.
2.3. CUDA Fortran code: R-Iber The main limitation of GPU-parallelised codes is the memory access
and, especially, the memory transfer (García-Feal et al., 2018; Vacondio
The sequential version of Iber is partially parallelised by means of the et al., 2014; Zhang and Jia, 2013). R-Iber reads, transfers and writes
open multi-processing (OpenMP) technique. However, the intrinsic memory data only when it is strictly necessary, with the aim of max­
limitations of this technique only allow for speed-ups corresponding to imising the code’s efficiency. Another issue, which also occurs in CPU
the number of CPU cores, at most. Despite the computing time of the computing, is that of the multi-core computing simulations: running
hydrodynamics being reduced notably, the simulation of long stretches several simulations at the same time in the same processor unit creates a
of rivers, even entire rivers, or long time series for fish habitat assess­ significant bottleneck, particularly for GPU computing, and reduces the
ment procedures remains unapproachable with techniques based solely global computing capacity of the device (Morales-Hernández et al.,
on CPU computing. 2021). To minimise heavy performance penalties, R-Iber allows for the
The integration of several calculation modules in the Fortran-based selection of the GPU device on which to carry out each calculation. Thus,
code of Iber was originally functionality-oriented, prioritising the full the option of a multi-processor computation on a single computer with
integration of all modules instead of multi-core processing efficiency. several GPUs or on a cluster of GPUs is available to practitioners. Fig. 1
Due to the growing application of more complex and detailed simula­ presents the main workflow of the R-Iber code.
tions in river hydrodynamics demands an increase in computational Additionally, to achieve maximum performance on GPU processors
capacity and, thus, a code re-ordering and optimization was already on arithmetic operations, it is necessary to use as many single-precision

Fig. 1. General flow chart of the R-Iber code for a particular GPU computation. The sketch also presents the flow for n simulations on n different devices.

3
M. Sanz-Ramos et al. Environmental Modelling and Software 161 (2023) 105628

arithmetic operations as possible (Hwu et al., 2009), as it is done in elevation composed of an increasing slope (≈0.02%), followed by a
R-Iber. Using single-precision operations implies round-off errors in the decreasing slope (≈0.25%), and ending in an increasing slope (≈0.3%).
calculation results, i.e. after a single floating-point operation the result is Thus, a depressed area was generated, where the water would be
rounded after the seventh digit. However, CPU-based simulations retained. As boundary condition, a water level was imposed, varying
perform double-precision arithmetic operations, as it is done in Iber, from 9.7 to 10.35 m in 1 h, keeping a constant value for 8 h to ensure the
with similar performance to single-precision (Narumi et al., 2011). depressed area was completely filled, then returning to 9.7 m in 1 h.
On the other hand, the habitat module of Iber, and also of R-Iber, is Two calculation scenarios were performed: one following the pro­
computed as a post-process because the currently implemented fish posed domain discretisation, with 1125 mesh elements, and the other
suitability curves only depend on the hydrodynamics. Thus, for opti­ using the full digital elevation model (DEM) resolution, which provides
mum memory access/transfer, the information exchanged between hy­ 26,341 calculation points. Simulations were also done with the CUDA
drodynamics and habitat must be done directly in the host and only C++ GPU code of Iber, called Iber+ (García-Feal et al., 2018).
when the results have to be written. This moment may occur when the Fig. 2 shows the water evolution at the two control points: P1
simulation time is greater than the user-defined time to write the results, (400,50) and P2 (600,50). The three numerical models show the same
or when the steady flow conditions are accomplished if the ‘stepped hydraulic behaviour: water accumulated in the depressed area at an
discharge’ option is used. Additionally, the process of writing the results elevation of 10.25 m, flush with the topography.
utilises an independent CPU core (Fig. 1). Therefore, once all data is The computational times of the reference model (Iber v3.1 sequen­
transferred back to the host, the habitat module is computed in the CPU tial) were 44 and 2064 s for the proposed discretisation and the full DEM
while the GPU continues calculating without any performance penalty. models, respectively. The speed-ups of R-Iber reached 3.4 and 39.4-
It is worth noting that 2D-SWE-based models do not solve the energy times, while Iber + reached 2.8 and 40.5-times, respectively, for the
equation, which in 1D only depends on the energy balance between two GTX 1660 Ti. Similar values were obtained for the GTX 980 Ti and the
consecutive sections. In 2D models, a steady state is reached when the RTX 3070 devices (Table 2).
inlet discharge is equal to the outlet discharge and, thus, the whole flood It is worth noting that for the proposed discretisation, in which only
front propagation process must be solved anyway. This implies an in­ 1125 elements were used, the speed-up was not significant because the
crease in the computational time because the simulation time needed to parallelised part of the code was faster than the time required to transfer
reach steady conditions in general is unknown. Therefore, longer information from the device to the host and to write the results file. This
simulation times are required to ensure a steady-state. bottleneck is practically neglected when the number of elements
This particular issue is solved with a specific option called ‘stepped increases.
discharge’, which allows for the simulation of steady flow conditions for
several discharges in the same model and directly provides flow-habitat 3.2. Test 2: Hydraulic jumps and flow obstructions
relationships such as the WUA. With this option enabled, GPU parallel
computation would potentially reach its maximum calculation capacity The second validation case is the ‘Test 6A’ of Defra’s benchmark
as the results would only be written a few times (at each step). (Néelz and Pender, 2013). This case aims to verify the capability of
2D-SWE-based models to simulate transcritical flow regimes, hydraulic
2.4. Hardware and characteristics jumps and wakes behind obstructions. For that purpose, the experiment
of Soares-Frazão and Zech (2007) was simulated in Iber and R-Iber using
Simulations were carried out with Iber v3.1, which contains the new two different mesh discretisation of rectangular elements: that proposed
code of R-Iber and are freely distributed through [Link]. by Defra’s benchmark (56,616 elements) and the full DEM resolution
Both algorithms use the Roe 1st order numerical scheme. Computations (226,352 elements).
with the sequential version (CPU) were done with a CPU Intel Core i7- An initial water elevation of 0.4 m was imposed from the beginning
9750H, while computations on the CUDA Fortran GPU version (R- of the flume to the narrowing of it, upstream of the obstacle. This
Iber) were launch with the same CPU but computed in different GPU configuration tries to reproduce a dam-break process in the downstream
devices (Table 1). The NVIDIA GPUs selected were the GeForce® GTX area, where the obstacle is placed non-symmetrically on a wet area in
and RTX series which, being quite common devices for standard desktop the flume. Results were extracted each 0.1 s, considering a fixed time
computers and mid-range to high-end gaming laptops, cover a wide step of 0.005 s. Six control points monitored the water depth and ve­
temporal range. locity during the experiment (G1 to G6).
All numerical approaches behaved similarly. Fig. 3 shows the evo­
lution of the water depth and velocity at points G1, G2 and G3. The
3. Validation
arrival time of the flood front was the same and suitably captured both
the supercritical and subcritical flow regimes. At point G2, the presence
3.1. Test 1: flooding of disconnected floodplains
of the obstacle generated a mobile hydraulic jump upstream at around
12 s, where water depth evolved from 0.03 to 0.1 m during the 10 s after
The first validation test was the so-called ‘Test 1’ proposed by the
the shock wave. In the case of the proposed mesh discretisation, the
United Kingdom Environment Agency within the Defra Flood and
hydraulic jump was produced at 9.4 s (yellow and grey dashed line),
Coastal Erosion Risk Management Research and Development Pro­
while with the full DEM the hydraulic jump was produced at 16.3 s for
gramme (Néelz and Pender, 2013). This test aims to assess the capability
Iber (black dashed line) and 11 s for R-Iber (red line). Iber generated a
of the 2D hydrodynamic models to reproduce the flooding of discon­
weak hydraulic jump at 11 s, which then transformed into a strong
nected floodplains and the wetting and drying processes.
hydraulic jump a few seconds later (16.3 s). This difference, in this
The domain consisted of a 700 × 100 m channel with a bottom
particular case of the simulation of a hydraulic jump with a very fine
discretisation of the calculation domain (grid of 0.05 × 0.05 m), can be
Table 1
attributed to the fact that R-Iber computes in single precision, while Iber
Graphical power units (GPU) utilized in the computations.
does it in double precision. Thus, after a single floating-point operation
Model Release date Architecture Cores Memory (type) the result is rounded after the seventh digit in R-Iber.
GeForce GTX 980 Ti 02/06/2015 Maxwell 2816 6 GB (GDDR5) In general, the reproduction of a hydraulic jump is quite complex to
GeForce GTX 1660 Tia 22/02/2019 Turing 1536 6 GB (GDDR6) reach with full precision, particularly in 2D-SWE-based modelling when
GeForce RTX 3070a 29/10/2020 Ampere 5120 8 GB (GDDR6)
the hydraulic jump evolves along time due to the discharge is not con­
a
This GPU is built for laptop. stant. A sharper solution is generally obtained (Néelz and Pender, 2013).

4
M. Sanz-Ramos et al. Environmental Modelling and Software 161 (2023) 105628

Fig. 2. Evolution of the water elevation at control points P1 (a) and P2 (b) for Iber (black dashed-line), Iber+ (blue dotted-line) and R-Iber (red line).

elements. In general, a mesh size reduction implies a better represen­

Table 2
tation of the hydrodynamics and, consequently, the eco-hydraulic
Computational time, in seconds, of the Test 1 for the sequential (Iber) and
assessment. However, higher computational times are required.
paralelized (R-Iber) versions. The speed-up is the quotien between parallelised
and sequential computational times.
The numerical discretisation of 2751 elements (M1) provided
smooth results, with only a weak recirculation area on the right side
Model Sequential [s] Parallelised [s] Speed-up [-]
downstream of the pool inlet and a considerably straight main jet
1125 26,341 1125 26,341 1125 26,341 (Fig. 5a). This recirculation area increased as the number of elements
els. els. els. els. els. els.
doubled (M2), and another two weak eddies appeared on the upstream
GeForce GTX 44 2064 13 61 3.4 33.6 and downstream corners on the left side, between the baffle and the wall
980 Ti (Fig. 5b). A computational mesh of eight times more elements (M3) than
GeForce GTX 13 52 3.4 39.4
the first mesh (M1) allowed for the proper definition of two separately
1660 Ti
GeForce RTX 10 57 4.4 41.1 recirculation areas on both sides of the main jet (Fig. 5c). This domain
3070 discretisation effectively defined the general flow pattern shown in
Fig. 4b, with a well-defined main jet. In the model with a mesh of 88,032
elements (M4), the flow patterns became more complex, representing
Despite that, the flow patterns were well captured by R-Iber, even using
small eddies near the bigger eddies to both sides of the main jet (Fig. 5d).
single precision and a Roe 1st order numerical scheme (Cea et al., 2010).
Numerical results were compared to the experimental results
The numerical simulation in the sequential code (Iber) spanned 172 s
throughout the mean depth of the middle traverse section (h0 ), the mean
for the proposed discretisation and 706 s for the full DEM model, while
flow depth of the pool (hm ), the maximum depth of the pool (hmax ), the
the parallelised R-Iber code computed the models in 35–45 s and from
minimum depth of the pool (hmin ), the depth of the slot (hs ), the absolute
130 to 197 s, respectively. This means an average speed-up of around
velocity of the slot (|vs |), the discharge coefficient of the slot (Cd =
3.6–5.4 times with respect of the sequential code (Table 3). It is worth
Q/(bhs vs ), where b is the slot width) and the energy dissipation rate in
noting that the necessity of writing results at very small-time steps (0.1
the pool (Ed = (ρgQS0 )/(Bh0 ), where ρ is the flow density, g is the
s) considerably reduces the global simulation time. If only computa­
gravitational acceleration, Q is the flow discharge, S0 is the bottom slope
tional time is considered, depending on the GPU used, the speed-up of R-
and B is the pool width). Table 4 summarises the experimental and
Iber reaches values from 65.4 to 75.3 and from 96.8 to 140.4-times for
numerical results obtained with each numerical mesh.
the proposed and full DEM models, respectively (Table 3, values in
In general, R-Iber performed suitably for all the simulations,
brackets).
particularly for the finer meshes (Table 5), which provided a good fit for
the observations reproducing the flow patterns, such as eddies and the
3.3. Test 3: fishway and recirculation areas shape and dimension of the main jet. The M1 mesh (2751 elements)
showed the poorest results, with absolute relative errors slightly above
Test 3 is the numerical reproduction of the laboratory experiment 10% in some hydraulic parameters (h0 , hm and Ed ). By contrast, the M3
presented by Puertas et al. (2004), which consisted of the hydraulic and M4 meshes (22,008 and 88,032 elements, respectively), provided
characterisation of a vertical slot fishway. The original flume device was suitable results with errors below 3.5%. Numerical results of the
a 12 m long channel divided into 11 pools (four for the design of Test 1 discharge coefficient (Cd ) and the energy dissipation rate (Ed ) also
and four for the Test 2 design, while the rest were the upper, interme­ adjusted suitably to the observations, with a relative error less than 3.5%
diate, and lower pools). Pool #7, which belonged to the Test 1 design for M3 and M4 meshes.
(Fig. 4a), was analysed considering half of the flume, a 5.7% slope and a The benefit of using the CUDA Fortran GPU parallelisation R-Iber
discharge of 0.0741 m3/s. A detailed description of the flume experi­ instead of Iber is the considerable reduction in computational time
ments can be found in Puertas et al. (2004). (Table 6). In this particular case, R-Iber reached a speed-up from 23.9
Previous numerical experiments confirm the capability of 2D-SWE- (GTX 980 Ti) to 45.1-times (RTX 3070) for mesh M4 in comparison with
based models for modelling the free surface flow in vertical slot fish­ the sequential version. Additionally, the efficiency of the code can be
ways (Bermúdez et al., 2010; Cea et al., 2007; Puertas et al., 2012). Thus, evaluated by accounting for the computational cost per number of ele­
this validation case aimed to reproduce main (jet) and particular ments. The computational time required by Iber increased from 24.7
(recirculation areas) flow patterns (Fig. 4b) by considering four mesh (M1) to 106.6 s (M4) per 1000 elements, while in the case of R-Iber, it
discretisation: 2751 (M1), 5502 (M2), 22,008 (M3) and 88,032 (M4) decreased from 5 to 2.9 s. Thus, it can be seen that the new parallelised

5
M. Sanz-Ramos et al. Environmental Modelling and Software 161 (2023) 105628

Fig. 3. Evolution of the water depth (up) and velocity (down) at control points G1, G2, and G3. Grey circles: experimental data; Dashed lines: Iber for the proposed
discretisation (grey) and the full DEM (black); Continuous lines: R-Iber for the proposed discretisation (yellow) and the full DEM (red).

code of Iber, R-Iber, provides not only a valuable computational cost

Table 3
reduction but also a notable improvement in the efficiency of the com­
Computational time, in seconds, of the Test 2 for the sequential (Iber) and
putations when the number of elements of the model increases.
paralelized (R-Iber) versions. The speed-up is the quotien between parallelised
and sequential computational times.
4. Study cases
Model Sequential [s] Paralleliseda [s] Speed-upa [-]

56,616 226,352 56,616 226,352 56,616 226,352 The performance of the model has been tested in two real study
els. els. els. els. els. els. cases: Eume River and Cinca River (Fig. 6). The Eume case study is
GeForce 172 706 45 (3) 197 s (7) 3.8 3.6 focused on showing the benefits of using not only the GPU-parallelised
GTX (65.4) (96.8) eco-hydraulic model of Iber, R-Iber, but also the noticeable improve­
980 Ti
ment of using the ‘stepped discharge’ inlet condition instead of creating
GeForce 42 (2) 158 s (6) 4.1 4.5
GTX (70.8) (110.1) a separate model for each discharge.
1660 An intrinsic benefit of the GPU parallelisation of 2D-SWE-based
Ti models is the capability to assess long-term hydrodynamics in long-
GeForce 35 (2) 130 s (5) 4.9 5.4 rivers in a reasonable timeframe. If a fish habitat module is integrated
RTX (75.3) (140.4)
3070
into the code, such as in R-Iber, the eco-hydraulic evaluation of rivers
will no longer be limited to a few metres or short river stretches. To
a
Net computational time, without the time required to write the results (in demonstrate this, the capacity of R-Iber to simulate long-term fish
brackets).
suitability variations in long-river stretches (even entire rivers) is

Fig. 4. (a) Geometry of the Test 1 design pool (distances in meters). (b) General flow pattern in Test 1 design pool (sketch adapted from Puertas et al. (2004)).

6
M. Sanz-Ramos et al. Environmental Modelling and Software 161 (2023) 105628

resulting WUA is reduced for all flows and stadiums of the target specie,
remaining below 14% of the total wetted area. The maximum suitability
conditions were achieved in the fry for a discharge of 0.560 m3/s
(12.78% of the wetted area) and in the juvenile and adult stadium for
2.176 m3/s (fry, 8.8%; adult, 6.2% of the wetted area).
Identical hydrodynamic and habitat results were obtained with Iber
and R-Iber, but with significant differences in the computational time
(Table 7). The simulation process of computing each discharge per
model required between 195 and 357 s in the sequential mode (CPU),
while the GPU version required less than 30 s for each simulation
(Fig. 8a). Accounting for the total computational time of all the models,
a global time of 2844 s was necessary for the CPU version. This time was
reduced to 207 s for the GTX 980 Ti, 189 s for the GTX 1600 Ti, and 124 s
for the RTX 3070 (Fig. 8a).
The benefits of using the ‘stepped discharge’ option instead of a
separate model for each discharge was evaluated as the quotient be­
tween the addition of the time of each model and the time of the model
that used the ‘stepped discharge’ option. The advantage of using the
‘stepped discharge’ option instead of a separate model for each
Fig. 5. Numerical results of the depth-averaged velocity vectors (white rows)
discharge is proved by a 16.6% acceleration in the computation time in
and contour map (colored) in the x direction (negative velocities are plotted in
CPU-based simulations (Fig. 8b). When GPU calculations were applied,
black). Mesh of 2751 (a), 5502 (b), 22,008 (c) and 88,032 (d) elements.
speed-ups of 32–52% were achieved (Fig. 8b). Using the RTX 3070 GPU
did not reduce the computation time as much as possible due to the
presented in the case study of the Cinca River. model’s low number of elements. However, when a higher number of
elements discretise the domain, the advantage of using a more powerful
4.1. Improved IFIM for 2D-SWE-based models: Eume River GPU clearly shows its potentialities (as shown in section 4.2).
It is worth noting that only the computation time was considered in
The analysed stretch of the Eume River, located in the north west of the previous comparison: the time needed for changing the boundary
Spain (Fig. 6a), is characterised by two well-differentiated geomorpho­ condition was not taken into consideration. If a constant time of 30 s for
logical facies: a discontinuous channel with jumps and pools upstream each simulated discharge is considered as the necessary time for
and a continuous channel downstream. All discharges were evaluated changing the boundary conditions and re-starting the simulation, a
with Iber and R-Iber considering a separate model for each discharge global speed-up of around 350% would be reached for the GPU com­
and also using the ‘stepped discharge’ option, which consisted of a putations of the ‘stepped discharge’ option.
unique model. In this case, the tolerance was fixed at 1%: a new ‘step’ of
the series of inlet discharges was automatically imposed when the ab­
solute value of the outlet discharge minus the inlet discharge, all divided 4.2. Long-term simulations in long-rivers stretches: Cinca River
by the inlet discharge, was less than the tolerance.
The study area of 3 km in length (Fig. 6b) was discretised into 72,426 Cinca is a 191 km long mountain river that begins in the central part
calculation elements. Twelve discharges were evaluated, ranging from of the south of the Pyrenees and drains into the Ebro River, in Spain
0.149 to 11.770 m3/s. Salmo trutta was the target specie, and the depth- (Fig. 6a). The analysed stretch of the Cinca River is 11 km long, starting
and velocity-dependent suitability curves described in CHE and ACA, 2.5 km upstream of the Bellós River junction and ending at the Ara River
(2008) were used as biological information. A detailed description of the junction (Fig. 6c). This stretch is characterised as a braided river with a
study area can be found in Sanz-Ramos et al. (2019). considerably active morphology (Béjar et al., 2018; Vericat et al., 2017).
As expected, the WUA depended on the hydrodynamics of each The flood event of 9 October 2014, which had two peak discharges
evaluated flowing discharge. Fig. 7 shows, in the transition zone be­ and spanned 10 days (CHE, 2021), was simulated with three different
tween two geomorphological facies (Fig. 6b, hollow black square), the meshes of 74k elements (M1), 299k elements (M2), and 1885k elements
elemental suitability evaluated for different discharges and the adult (M3). The river domain was discretised using an irregular mesh of
stadium (product of the depth- and velocity-suitability). Low discharges
produced a narrow-wet area, particularly at the facies characterised by Table 5
jumps and pools (Fig. 7a and b). The most fish suitability was obtained in Relative error (in absolute value) of the numerical results for T1 design, 5.7%
zones with high depths (pools). When the discharge increased (Fig. 7c), slope, and Q = 0.0741 m3/s.
the discontinuous channel became continuous, the velocity increased, h0 hm hmax hmin hb |v0 | Cd Ed
and the suitability of this area decreased. This is clearly seen in the lower [%] [%] [%] [%] [%] [%] [%] [%]
part of the river, where a moderate depth and low velocity area provided M1 11.7 10.2 7.4 8.3 8.4 9.7 2.3 10.5
a high suitability zone. However, high discharges (11.77 m3/s) modified M2 7.1 5.6 3.9 2.1 4.1 8.9 5.7 6.7
the hydrodynamics in this area, reducing the suitability of the central M3 3.5 2.4 2.3 1.0 2.6 4.4 2.1 3.4
M4 2.7 2.0 3.4 2.6 3.3 3.1 0.1 2.7
part of the channel and displacing it to the lateral sides (Fig. 7d). The

Table 4
Experimental and numerical results for the Test 1 design, 5.7% slope and Q = 0.0741 m3/s.
h0 [m] hm [m] hmax [m] hmin [m] hs [m] |vs | [m/s] Cd Ed [W/m2]

Exp. 0.604 0.608 0.652 0.562 0.604 0.880 0.870 69.350

M1 0.675 0.670 0.700 0.608 0.655 0.795 0.890 62.040
M2 0.647 0.642 0.677 0.574 0.629 0.801 0.919 64.675
M3 0.625 0.622 0.667 0.568 0.620 0.841 0.888 66.980
M4 0.620 0.620 0.674 0.577 0.624 0.852 0.871 67.490

7
M. Sanz-Ramos et al. Environmental Modelling and Software 161 (2023) 105628

Table 6
Computational time, in seconds, of the Test 3 for the sequential (Iber) and paralelized (R-Iber) versions. The speed-up is the quotien between parallelised and
sequential computational times.
Model Sequential [s] Parallelised [s] Speed-up [-]

M1 M2 M3 M4 M1 M2 M3 M4 M1 M2 M3 M4

GeForce GTX 980 Ti 68 136 1138 9387 22 26 94 392 3.1 5.2 12.1 23.9
GeForce GTX 1660 Ti 14 18 55 257 4.9 7.6 20.7 36.5
GeForce RTX 3070 13 15 47 208 5.2 9.1 24.2 45.1

Fig. 6. (a) General location of Eume and Cinca rivers (red lines). (b) Study area of the Eume River (yellow line) and the location of the results shown in Fig. 7 (hollow
black square). (c) Study area of the Cinca River (yellow line).

triangular elements with 220 (M1), 887 (M2) and 5577 (M3) elements
per hectare, the last one being one order of magnitude above the com­
mon values used in flood studies (Sanz-Ramos et al., 2020a). The mesh
was updated with the 2 × 2 m DEM provided by National Geographic
Institute of Spain (IGN, 2022). During these 10 days, where discharges
ranged from 10 (base flow) to 250 m3/s (maximum peak discharge),
variations in the fish habitat suitability were evaluated for salmonids,
Salmo trutta (CHE and ACA, 2008), and cyprinids, Chondrostoma poly­
lepis and Barbus bocagei (Martinez, 2000). These target species were
selected within the Cinca’s basin characteristics (CHE, 2002).
In contrast with the IFIM, the main result of this analysis is the
evolution of the WUA over time, besides the evolution of the fish suit­
ability distribution. Fig. 9a presents the evolution of the WUA (in %) for
the three abovementioned species and stadiums (fry, juvenile, and
adult) and the M2 domain discretisation. Changes in the hydrological
regime also modify the WUA. As expected, low discharges favoured a
more suitable area for fry and juvenile stadiums. Despite high discharges
generating high depths and velocities, in general, other areas flooded

Table 7
Computational time, in seconds, of the Study Case 1 for the sequential (Iber) and
paralelized (R-Iber) versions. The speed-up is the quotien between parallelised
and sequential computational times.
Fig. 7. Hydraulic suitability at the stream analysed for the adult stadium of Model Sequential [s] Parallelised [s] Speed-up [-]
Salmo trutta according to the suitability curves of Bovee (1982): (a) 0.149 m3/s,
72,426 els. 72,426 els. 72,426 els.
(b) 0.560 m3/s, (c) 1.567 m3/s, and (d) 11.77 m3/s (background image:
IGN (2021)). GeForce GTX 980 Ti 2844 s 207 s 13.7
GeForce GTX 1660 Ti 189 s 15.0
GeForce RTX 3070 124 s 23.0

8
M. Sanz-Ramos et al. Environmental Modelling and Software 161 (2023) 105628

Fig. 8. (a) Simulation time of each discharge computed independently (low discharges: dark green; high discharges: light green). (b) Time reducing obtained when
the ‘stepped discharge’ option is applied in comparison to the independent computing of each discharge.

higher resolution, which are the bases of WUA assessment. In general,

the M2 mesh provided a mean increment of 7.9% of the WUA in com­
parison to the M1 mesh, while the M3 mesh reached up to 18.1% more
WUA than the M1 mesh.
In terms of computational time, the sequential code (Iber) required
7.4 h for the mesh M1, 71.3 h (≈3 days) for the mesh M2, and 1081.3 h
(≈45 days) for the mesh M3. For the M1, M2 and M3 mesh configura­
tions, R-Iber reached speed-ups from 19.8 to 183.1-times, respectively
(Table 8).
The efficiency of the code can be evaluated by the computation time
required per 1000 elements. The time required by Iber increased when
the number of elements increased from 356 to 2062 s, while in the case
of R-Iber, the time remained almost constant, at around 10–20 s per
1000 elements. A similar trend was obtained in the validation tests.
Complementary to the traditionally used IFIM, the assessment of the
fish habitat suitability considering real discharge scenarios –even in real
time– can be a step forward in eco-hydraulics as it is not limited to
steady conditions (particular constant discharges). Furthermore, it also
encourages the study of the convenience, or inconvenience, of reintro­
Fig. 9. (a) Evolution of the weight useable area (WUA) during the event of 9 ducing native species and the potential prevalence of invasive species to
October 2014 for adult (continuous line), juvenile (dashed-line) and fry (dotted- the natural or induced hydrologic regime of a river. Additionally, since
line) stadiums of Chondrostoma polylepis (black), Barbus bocagei (green) and hydrodynamic conditions in rivers vary over time due to morphody­
Salmo trutta (red) (mesh M2). (b) Sensibility analysis of WUA for the stadiums namic changes, if sediment transport (Hung et al., 2022; Pisaturo et al.,
of Barbus bocagei (continuous line: adult; dashed line: juvenile; dotted line: fry) 2021) and water quality (Wang and Lin, 2013) processes are considered,
considering different mesh discretisation: M1 (black), M2 (green) and M3 (red). this kind of modelling strategy would allow for the analysis of such
The grey area represents the hydrograph. conditions from a holistic point of view.

provided more suitable zones for the Chondrostoma polylepis (black lines) 5. Conclusions
and Barbus bocagei (green lines), while the WUA of the Salmo trutta (red
lines) decreased. Currently, eco-hydraulic numerical tools are widely used for river
Fig. 9b exemplifies the variation in the WUA evolution for adult habitat rehabilitation, restoration and enhancement purposes. These
(continuous line), juvenile (dashed line) and fry (dotted line) of the tools undergo continuous development due to the advances in data-
species Barbus bocagei when the models of 220 (M1, black line), 887 acquisition and data-treatment and the use of high-performance tech­
(M2, green line) and 5577 (M3, red line) elements per hectare were niques to formulate more efficient, robust, and powerful numerical
considered. As has also been demonstrated in sections 3.2 and 3.3, a codes. Nevertheless, 2D-SWE-based eco-hydraulic numerical tools that
finer discretisation of the study area provides hydrodynamic results of a include an integrated GPU-parallelised hydrodynamic module and

Table 8
Computational time, in hours, of the Study Case 1 for the sequential (Iber) and paralelized (R-Iber) versions. The speed-up is the quotien between parallelised and
sequential computational times.
Model Sequential [h] Parallelised [h] Speed-up [-]

M1 M2 M3 M1 M2 M3 M1 M2 M3

GeForce GTX 980 Ti 7.4 71.3 1080 0.4 1.6 16.2 19.8 45.0 66.8
GeForce GTX 1660 Ti 0.3 1.1 11.6 26.6 64.3 92.8
GeForce RTX 3070 0.2 0.6 5.9 37.2 116.8 183.1

9
M. Sanz-Ramos et al. Environmental Modelling and Software 161 (2023) 105628

physical habitat module developed for overcoming the main computa­ calculation meshes providing high-resolution numerical models with a
tional limitations are lacking. better representation of the hydrodynamics and, consequently, the
R-Iber is a GPU parallelised hydrodynamic numerical tool that in­ habitat, with no performance penalties. Additionally, it has been proven
tegrates a physical habitat module; thus, it is a fully-integrated GPU- that the simulation of long river stretches, even entire rivers, or long-
parallelised eco-hydraulic tool. The code, based on Iber, was developed term habitat analysis can be computed with GPU-based numerical
in CUDA Fortran language for faster computations. The hydrodynamic codes within a convenient timeframe.
module of R-Iber was first validated using two benchmark test cases and
a laboratory experiment with a fishway. R-Iber was also applied to two Funding
real cases, one following the IFIM and the other simulating a 10-day real
flood event. The contract of the D.D.-S. is funded by the International Center for
The benefit of using a GPU-parallelisation eco-hydraulic tool instead Numerical Methods in Engineering (VAC-2021-1).
of a CPU-based tool is the significant reduction in computation time,
with speed-ups of one or two orders of magnitude (above 100-times with
respect CPU computations), but also a notable improvement in the ef­ Declaration of competing interest
ficiency of the computations when the model’s number of elements in­
creases is shown. Additionally, the option ‘stepped discharge’ The authors declare the following financial interests/personal re­
demonstrated to be a suitable solution to address Q-WUA computations lationships which may be considered as potential competing interests:
following the IFIM, being the computations 50% faster than an in- Danial Dehghan Souraki reports financial support was provided by
cascade simulation process. This option uses the time strictly neces­ Centre Internacional de Mètodes Numèrics en la Enginyeria.
sary to reach the steady flow conditions without having to know the time
to reach the steady flow conditions for any discharge and reducing the Data availability
computational time to a minimum.
R-Iber can simulate calculation domains discretised with finer The authors do not have permission to share data.

Appendix A. CUDA Fortran structure of R-Iber

This appendix shows the internal structure of the new code R-Iber, which is based on CUDA Fortran programming language (NVIDIA, 2022b). The
main structure of the code R-Iber is depicted in Fig. 1 of the main manuscript.
R-Iber code can be split in two main parts: the host and the device. The host relates to the central processing unit, or CPU, and controls the main
instructions of the code (variable reading, allocatation, saving, transfer, etc.) and directives (e.g. loops), even the execution of the graphics processing
unit, or GPU. This last action is done through Kernel Loop Directive (KLD).
The host part of code of R-Iber can be synthetized as follow.

Loops of the original code of Iber have been adapted to the new structure of the CUDA Fortran programming and both codes, Fortran and CUDA
Fortran, can be called from the host as follows.

10
M. Sanz-Ramos et al. Environmental Modelling and Software 161 (2023) 105628

A kernel may be invoked with many thread blocks, each with the same thread block size. The thread blocks are organized into a one-, two-, or three-
dimensional grid of blocks, so each thread has a thread index within the block, and a block index within the grid (Zhang and Jia, 2013). When invoking
a kernel, the first argument in the chevron <<<>>> syntax is the grid size, and the second argument is the thread block size. Thread blocks must be
able to execute independently; two thread blocks may be executed in parallel or one after the other, by the same core or by different cores.

References resolution simulations of rural and urban flash floods. Environ. Model. Software 156,
105480. [Link]
Carlotto, T., Borges Chaffe, P.L., Innocente dos Santos, C., Lee, S., 2021. SW2D-GPU: a
Béjar, M., Vericat, D., Nogales, I., Gallart, F., Batalla, R.J., 2018. Efectos de las
two-dimensional shallow water model accelerated by GPGPU. Environ. Model.
extracciones de áridos sobre el transporte de sedimentos en suspensión en ríos de
Software 145, 105205. [Link]
montaña (alto río Cinca, Pirineo Central). Cuadernos Invest. Geogr. 44, 641–658.
Cassan, L., Roux, H., Courret, D., Richard, S., 2022. Sensitivity of aquatic habitat
[Link]
modeling to hydrodynamic calibration. J. Ecohydraulics 1–9. [Link]
Benjankar, R., Tonina, D., McKean, J.A., Sohrabi, M.M., Chen, Q., Vidergar, D., 2018.
10.1080/24705357.2022.2049014.
Dam operations may improve aquatic habitat and offset negative effects of climate
Cea, L., Bermúdez, M., Puertas, J., Bladé, E., Corestein, G., Escolano, E., Conde, A.,
change. J. Environ. Manag. 213, 126–134. [Link]
Bockelmann-Evans, B., Ahmadian, R., 2016. IberWQ: new simulation tool for 2D
jenvman.2018.02.066.
water quality modelling in rivers and shallow estuaries. J. Hydroinf. 18, 816–830.
Bermúdez, M., Cea, L., Puertas, J., Conde, A., Martín, A., Baztán, J., 2017. Hydraulic
[Link]
model study of the intake-outlet of a pumped-storage hydropower plant. Eng. Appl.
Cea, L., Bladé, E., 2015. A simple and efficient unstructured finite volume scheme for
Comput. Fluid Mech. 11, 483–495. [Link]
solving the shallow water equations in overland flow applications. Water Resour.
19942060.2017.1314869.
Res. 51, 5464–5486. [Link]
Bermúdez, M., Cea, L., Puertas, J., Rodríguez, N., Baztán, J., 2018. Numerical modeling
Cea, L., Bladé, E., Sanz-Ramos, M., Fraga, I., Sañudo, E., García-Feal, O., Gómez-
of the impact of a pumped-storage hydroelectric power plant on the reservoirs’
Gesteira, M., González-Cao, J., 2010. Benchmarking of the Iber Capabilities for 2D
thermal stratification structure: a case study in NW Spain. Environ. Model. Assess.
Free Surface Flow Modelling. Universidade da Coruña. Servizo de Publicacións.
23, 71–85. [Link]
[Link]
Bermúdez, M., Puertas, J., Cea, L., Pena, L., Balairón, L., 2010. Influence of pool
Cea, L., Pena, L., Puertas, J., Vázquez-Cendón, M.E., Peña, E., 2007. Application of
geometry on the biological efficiency of vertical slot fishways. Ecol. Eng. 36,
several depth-averaged turbulence models to simulate flow in vertical slot fishways.
1355–1364. [Link]
J. Hydraul. Eng. 133, 160–172. [Link]
Bladé, E., Cea, L., Corestein, G., 2014a. Numerical modelling of river inundations [in
133:2(160).
Spanish]. Ing. del agua 18, 68. [Link]
CHE, 2021. SAIH Ebro [WWW Document]. Confed. Hidrográfica del Ebro. [Link]
Bladé, E., Cea, L., Corestein, G., Escolano, E., Puertas, J., Vázquez-Cendón, E., Dolz, J.,
[Link]. accessed 3.15.21.
Coll, A., 2014b. Iber: river flow numerical simulation tool [in Spanish]. Rev. Int.
CHE, 2002. Study of the Integral Ecological Quality of the Most Important Fluvial
Métodos Numéricos Cálculo Diseño Ing. 30, 1–10. [Link]
Sections of the Cinca River, vol. III (in Spanish).
rimni.2012.07.004.
CHE, ACA, 2008. Evaluation of Environmental Flows and Biological Validation in
Bladé, E., Sánchez-Juny, M., Arbat, M., Dolz, J., 2019a. Computational modeling of fine
Significant Sections of the River Network of Catalonia (in Spanish).
sediment relocation within a dam reservoir by means of artificial flood generation in
García-Feal, O., González-Cao, J., Gómez-Gesteira, M., Cea, L., Domínguez, J.,
a reservoir cascade. Water Resour. Res. 55, 3156–3170. [Link]
Formella, A., 2018. An accelerated tool for flood modelling based on iber. Water 10,
2018WR024434.
1459. [Link]
Bladé, E., Sanz-Ramos, M., Dolz, J., Expósito-Pérez, J., Sánchez-Juny, M., 2019b.
Hamilton, S.H., Pollino, C.A., Jakeman, A.J., 2015. Habitat suitability modelling of rare
Modelling flood propagation in the service galleries of a nuclear power plant. Nucl.
species using Bayesian networks: model evaluation under limited data. Ecol. Model.
Eng. Des. 352, 110180 [Link]
299, 64–78. [Link]
Boudreault, J., Bergeron, N.E., St-Hilaire, A., Chebana, F., 2022. A new look at habitat
Hung, H.-J., Lo, W.-C., Chen, C.-N., Tsai, C.-H., 2022. Fish’ habitat area and habitat
suitability curves through functional data analysis. Ecol. Model. 467, 109905
transition in a river under ordinary and flood flow. Ecol. Eng. 179, 106606 https://
[Link]
[Link]/10.1016/[Link].2022.106606.
Bovee, K.D., 1982. A guide to stream habitat analysis using the Instream Flow
Hwu, W.-M., Rodrigues, C., Ryoo, S., Stratton, J., 2009. Compute unified device
Incremental Methodology. Instream Flow Information, FWS/OBS¬82/26, 12, p. 248.
architecture application suitability. Comput. Sci. Eng. 11, 16–26. [Link]
Fort Collins, Color.
10.1109/MCSE.2009.48.
Buttinger-Kreuzhuber, A., Konev, A., Horváth, Z., Cornel, D., Schwerdorf, I., Blöschl, G.,
IGN, 2022. Digital Elevation Models [WWW Document]. Cent. Descargas. [Link]
Waser, J., 2022. An integrated GPU-accelerated modeling framework for high-
[Link]/CentroDescargas/. accessed 4.26.22.

11
M. Sanz-Ramos et al. Environmental Modelling and Software 161 (2023) 105628

IGN, 2021. Ortofotos e imágenes satélite [WWW Document]. Cent. Descargas. [Link] Ruiz-Villanueva, V., Bladé, E., Sánchez-Juny, M., Marti-Cardona, B., Díez-Herrero, A.,
[Link]/CentroDescargas/[Link]?Serie=PNOAH. accessed Bodoque, J.M., 2014. Two-dimensional numerical modeling of wood transport.
4.6.21. J. Hydroinf. 16, 1077. [Link]
Jowett, I.G., 2004. RHYHABSIM river hydraulics and habitat simulation computer Sanz-Ramos, M., Bladé, E., Escolano, E., 2020a. Optimización del cálculo de la Vía de
manual. Software Manual 3.2 version 77. Intenso Desagüe con criterios hidráulicos. Ing. del agua 24, 203. [Link]
Martinez, F., 2000. Régimen Ambiental De Caudales: Estimación De Las Condiciones De 10.4995/ia.2020.13364.
Habitabilidad Para La Ictiofauna, pp. 1–19. Sanz-Ramos, M., Bladé, E., González-Escalona, F., Olivares, G., Aragón-Hernández, J.L.,
Meselhe, E.A., Georgiou, I., Allison, M.A., McCorquodale, J.A., 2012. Numerical 2021. Interpreting the Manning roughness coefficient in overland flow simulations
modeling of hydrodynamics and sediment transport in lower Mississippi at a with coupled hydrological-hydraulic distributed models. Water 13, 3433. https://
proposed delta building diversion. J. Hydrol. 472–473, 340–354. [Link] [Link]/10.3390/w13233433.
10.1016/[Link].2012.09.043. Sanz-Ramos, M., Bladé, E., Palau, A., Vericatl, D., Ramos-Fuertes, A., 2019.
Meza Rodríguez, D., Martinez Rivera, L.M., Olguín López, J.L., Aguirre García, Á., 2019. IberHABITAT: assessment of physical habitat suitability and weighted useable area
Simulation of physical habitat in Ayuquila-Armeria river in the west of Mexico. for fishes. Application in the Eume River. Ribagua 6, 158–167. [Link]
Cienc. e Ing. Neogranadina 29, 53–68. [Link] 10.1080/23863781.2019.1664273.
Morales-Hernández, M., Sharif, M.B., Kalyanapu, A., Ghafoor, S.K., Dullo, T.T., Sanz-Ramos, M., Martí-Cardona, B., Bladé, E., Seco, I., Amengual, A., Roux, H.,
Gangrade, S., Kao, S.-C., Norman, M.R., Evans, K.J., 2021. TRITON: a Multi-GPU Romero, R., 2020b. NRCS-CN estimation from onsite and remote sensing data for
open source 2D hydrodynamic flood model. Environ. Model. Software 141, 105034. management of a reservoir in the Eastern Pyrenees. J. Hydrol. Eng. 25, 05020022
[Link] [Link]
Narumi, T., Hamada, T., Nitadori, K., Sakamaki, R., Yasuoka, K., 2011. Fast quasi double- Shim, T., Kim, Z., Seo, D., Kim, Y.-O., Hwang, S.-J., Jung, J., 2020. Integrating hydraulic
precision method with single-precision hardware to accelerate scientific and physiologic factors to develop an ecological habitat suitability model. Environ.
applications. Int. J. Comput. Methods 8, 561–581. [Link] Model. Software 131, 104760. [Link]
S0219876211002708. Soares-Frazão, S., Zech, Y., 2007. Experimental study of dam-break flow against an
Néelz, S., Pender, G., 2013. Benchmarking the Latest Generation of 2D Hydraulic isolated obstacle. J. Hydraul. Res. 45, 27–36. [Link]
Modelling Packages. Report - SC120002. Environment Agency, Horison House, 00221686.2007.9521830.
Deanery Road, Bristol, BS1 9AH. Stamou, A., Polydera, A., Papadonikolaki, G., Martínez-Capel, F., Muñoz-Mas, R.,
Nestler, J.M., Stewardson, M.J., Gilvear, D.J., Webb, J.A., Smith, D.L., 2016. Papadaki, C., Zogaris, S., Bui, M.D., Rutschmann, P., Dimitriou, E., 2018.
Ecohydraulics exemplifies the emerging “paradigm of the interdisciplines. Determination of environmental flows in rivers using an integrated hydrological-
J. Ecohydraulics 1, 5–15. [Link] hydrodynamic-habitat modelling approach. J. Environ. Manag. 209, 273–285.
Nones, M., 2019. Numerical Modelling as a Support Tool for River Habitat Studies: an [Link]
Italian Case Study, vol. 11. [Link] Water Steffler, P., Blackburn, J., 2002. Two-dimensional depth averaged model of river
(Switzerland). hydrodynamincs and fish habitat. In: User Manual River2D. University of Alberta,
NVIDIA, 2022a. NVIDIA HPC SDK Version 22.3 Documentation [WWW Document]. URL. Canada.
[Link] accessed 2.25.22. Tonina, D., Jorde, K., 2013. Approaches for ecohydraulic non-numerical models.
NVIDIA, 2022b. NVIDIA CUDA Fortran Programming Guide [WWW Document]. URL. Ecohydraulics An Integr. Approach 31–74.
[Link] accessed Toro, E.F., 2009. Riemann Solvers and Numerical Methods for Fluid Dynamics. Springer,
2.23.22. Berlin/Heidelberg, Germany. [Link]
Palau, A., Alcázar, J., 2012. The basic flow method for incorporating flow variability in Vacondio, R., Dal Palù, A., Mignosa, P., 2014. GPU-enhanced finite volume shallow
environmental flows. River Res. Appl. 28, 93–102. [Link] water solver for fast flood simulations. Environ. Model. Software 57, 60–75. https://
rra.1439. [Link]/10.1016/[Link].2014.02.003.
Palau, A., Alcázar, J., Rocaspana, R., Aparicio, E., Mariño, F., 2016. Misusing physical Vericat, D., Wheaton, J.M., Brasington, J., 2017. Revisiting the morphological approach.
habitat assessment techniques for environmental flows calculations. In: Proceedings In: Gravel-Bed Rivers. John Wiley & Sons, Ltd, Chichester, UK, pp. 121–158. https://
of 11th International Symposium on Ecohydraulics. Melbourne (Australia). [Link]/10.1002/9781118971437.ch5.
Pisaturo, G.R., Folegot, S., Menapace, A., Righetti, M., 2021. Modelling fish habitat Wang, F., Lin, B., 2013. Modelling habitat suitability for fish in the fluvial and lacustrine
influenced by sediment flushing operations from an Alpine reservoir. Ecol. Eng. 173, regions of a new Eco-City. Ecol. Model. 267, 115–126. [Link]
106439 [Link] ecolmodel.2013.07.024.
Pisaturo, G.R., Righetti, M., Dumbser, M., Noack, M., Schneider, M., Cavedon, V., 2017. Wilkes, M.A., Neverman, A.J., Casas-Mulet, R., Adeva-Bustos, A., McCluskey, A.H.,
The role of 3D-hydraulics in habitat modelling of hydropeaking events. Sci. Total Ouellet, V., Vanzo, D., Franklin, P.A., Silva, A.T., 2016. Early careers on
Environ. 575, 219–230. [Link] ecohydraulics: challenges, opportunities and future directions. J. Ecohydraulics 1,
Puertas, J., Cea, L., Bermúdez, M., Pena, L., Rodríguez, Á., Rabuñal, J.R., Balairón, L., 102–107. [Link]
Lara, Á., Aramburu, E., 2012. Computer application for the analysis and design of Zhang, Y., Jia, Y., 2013. Parallelized CCHE2D flow model with CUDA fortran on graphics
vertical slot fishways in accordance with the requirements of the target species. Ecol. processing units. Comput. Fluids 84, 359–368. [Link]
Eng. 48, 51–60. [Link] compfluid.2013.06.021.
Puertas, J., Pena, L., Teijeiro, T., 2004. Experimental approach to the hydraulics of Zohmann, M., Pennerstorfer, J., Nopp-Mayr, U., 2013. Modelling habitat suitability for
vertical slot fishways. J. Hydraul. Eng. 130, 10–23. [Link] alpine rock ptarmigan (Lagopus muta helvetica) combining object-based
0733-9429(2004)130:1(10). classification of IKONOS imagery and Habitat Suitability Index modelling. Ecol.
Roe, P.L., 1986. A basis for the upwind differencing of the two-dimensional unsteady Model. 254, 22–32. [Link]
Euler equations. Numer. Methods Fluid Dyn. 2, 55–80.

12