Abstract
Understanding the complexity of the siltation process and sediment resuspension in shallow reservoirs is vital for maintaining the reservoir functionality and implementing sustainable sediment management strategies. The geometry of reservoirs plays an indispensable role in the appearance of various flow structures inside the basin and, consequently, the pattern of the morphological evolution. In this study, a three-dimensional numerical model, coupled with optimization algorithms, is used to investigate the morphological bed changes in two symmetric shallow reservoirs having hexagon and lozenge shapes. This work aims to evaluate the applicability, efficiency, and accuracy of the automatic calibration routine, which can be a suitable replacement for the time-consuming and subjective method of manual model calibration. In this regard, two sensitive parameters (i.e., roughness height and sediment active layer thickness) are assessed. The goodness-of-fit between the calculated bed levels and the measured topography from physical models are presented by different statistical metrics. From the results, it can be concluded that the automatically calibrated models are in reasonable agreement with the observations. Employing a suitable optimization algorithm, which finds the best possible combination of investigated parameters, can considerably reduce the model calibration time and user intervention.
HIGHLIGHTS
The flow structure and sedimentation pattern in symmetric expansions are numerically studied.
Local (GML) and a global (BB–BC) optimization algorithms are used to calibrate the numerical models of two symmetric shallow reservoirs.
GML outperforms BB–BC considering the convergence speed (efficiency) without trapping in local minima by having the same predicted parameter values (robustness).
INTRODUCTION
As general multi-purpose hydraulic structures, shallow reservoirs are used in different fields, serving as sediment and pollutant trapping tanks, retention ponds, water storage basins, and also for aquafarming (Dewals et al. 2020). Depending on several factors, such as geometrical aspects (i.e., expansion ratio: reservoir width/inlet width; and aspect ratio: reservoir length/reservoir width), hydraulic and boundary conditions, and sediment properties, the velocity field inside a reservoir emerges with various configurations, influencing the sedimentation along the basin. In turn, deposited sediments can further modify the flow pattern (Kantoush 2008; Dufresne et al. 2012). The flow behavior in shallow reservoirs, such that the horizontal dimensions are considerably larger than the water depth, can be considered a condition in which the flow passes through a narrow inlet channel into a sudden or moderate expanded basin. Thus, the flow field involves an entering jet associated with large-scale 2D horizontal coherent turbulent structures and recirculation zones, controlling the mass and momentum exchange. In addition, the flow field in symmetric expansions is specified by symmetry-breaking bifurcation under certain conditions. The importance of this behavior in engineering applications has been a subject of interest among many researchers in the field of fluid dynamics (Fani et al. 2012).
A comprehensive experimental study of laminar flow through a planar symmetric sudden expansion dates back to Durst et al. (1974), who found that the Reynolds number (Re) strictly governs the flow structure. They defined a critical Re value, beyond which the flow tends to become asymmetric. Cherdron et al. (1978) performed a similar study and confirmed that the reason for asymmetric flow formation lies in disturbances generated at the edge of the expansion, which amplify in the shear layers. Sobey & Drazin (1986) investigated the instabilities and bifurcation of flow in slight expansions by asymptotic analysis, numerical methods, and laboratory experiments. Shapira et al. (1990) carried out a linear stability analysis of symmetric flow considering different angles at the expansion section. A numerical and experimental investigation of the flow field in a symmetric expansion, indicating a symmetry-breaking (pitchfork) bifurcation point, was published by Fearn et al. (1990). They introduced a slight asymmetry in the inlet channel in their numerical model to reproduce the imperfections, which are inevitable in physical models. This small perturbation was also considered by Hawa & Rusak (2000). The effect of different expansion ratios on the occurrence of non-symmetrical flow regimes in symmetric sudden expansions was studied by Drikakis (1997). Many other researchers investigated the discussed phenomenon (i.e., the transition of the symmetric flow into the asymmetric state) in sudden expansions (e.g., Mizushima et al. 1996; Sarma et al. 2000; Quaini et al. 2016). The main focus of these studies is on determining a threshold for the critical Reynolds number as the bifurcation initiating factor and the effect of geometrical aspects on the flow structure.
In contrast to extensive investigations of laminar flow behavior downstream of sudden expansions, relatively limited studies exist about turbulent flow configurations (Escudier et al. 2002). As one of the earliest studies, Abbott & Kline (1962) experimentally investigated the turbulent flow through a symmetric plane sudden expansion. They noticed two same-length recirculation zones with a predominant central plane jet after the expansion. However, when the expansion ratio was increased, the asymmetric pattern was observed with two disproportionate recirculation regions by the deflection of the main jet toward one of the walls. They claimed that the reattachment length highly depends on the expansion ratio and is not sensitive to the Reynolds number value (turbulence intensity). Similar results were obtained by Mehta (1981), arguing the insensitivity of the asymmetric pattern for a range of high Reynolds numbers. Among other related studies on turbulent flow through symmetric expansions, the laboratory experiments of Casarsa & Giannattasio (2008) and the numerical simulations by De Zilwa et al. (2000) are worth mentioning. Most of these studies report the formation of the asymmetric flow pattern as a function of the expansion ratio.
Narrowing down the discussed phenomenon to shallow reservoirs, Kantoush (2008) and Kantoush et al. (2008a, 2008b) performed a comprehensive series of systematic investigations of the flow field with suspended sediments inside different geometries with suddenly expanded regions experimentally and numerically. They studied the effect of reservoir geometry alteration on the flow field, the silting process pattern, and the reservoir trap efficiency. The primary conclusion of these works is that an asymmetric flow pattern appears in symmetric geometries under certain conditions, affecting the process of suspended sediment transport and adjusting the spatial distribution of sediment deposits. Furthermore, the accumulation of deposited sediments due to the additional sediment supply can, in turn, alter the flow pattern. Regarding the related experimental studies in this field, Dufresne et al. (2010) classified flow patterns in shallow rectangular reservoirs into four stable categories. They found a symmetric flow pattern containing a straight jet from the inlet to the outlet with two equal-size recirculation zones at the sides for short-length reservoirs (called S0 type) and three asymmetric patterns with one, two, or three reattachment points (A1, A2, and A3) with unequal recirculation regions, depending on the length of the basin, for long-length reservoirs. In reservoirs with intermediate length, an unstable flow field was also identified (S0/A1), in which the flow randomly oscillates between the symmetric (S0) and asymmetric (A1) patterns. They investigated the effects of dimensionless length and flow depth, lateral expansion ratio, and the Froude number on the median reattachment lengths and defined a shape factor to predict the flow behavior based on the geometric aspects. Similar experimental results to those of Kantoush (2008) were obtained by Camnasio et al. (2011), who categorized the flow field inside shallow rectangular reservoirs into a channel-like (CH-L) type for reservoirs with very short width, two symmetric (S0 and S1), and two asymmetric (A1 and A2) stable types, based on the reservoir expansion and aspect ratios.
Regarding the numerical assessment of the flow behavior in rectangular shallow basins, 2D models using the shallow water equations were employed by Dewals et al. (2008) and Dufresne et al. (2011). On the other hand, the use of 3D numerical models, employing the Reynolds-averaged Navier–Stokes equations, can be found in Esmaeili et al. (2016) and Lakzian et al. (2020). Numerical models of symmetric expanding channels are complicated cases where the entering jet can randomly follow one side of the reservoir. It means that for such cases, multiple solutions for the Navier–Stokes equations can exist, where the presence of an obstacle (deposited sediments), the grid resolution, discretization schemes, and the turbulence model may cause the convergence of the equations with a different flow direction. Hence, in contrast to the extensively studied rectangular geometries, we investigated the performance of a fully 3D numerical model (SSIIM) coupled with an automatic calibration tool (PEST) to evaluate the morphological bed changes in two shallow reservoirs having hexagon and lozenge shapes. The reservoirs are selected based on available high-resolution measurements for calibration and to test the automatic model calibration for cases with complicated flow behavior. Since PEST uses a gradient-based local optimization algorithm, to confirm its ability in finding the global optimum point over the search space and avoid the local minima, the models are also calibrated by a global optimization algorithm.
Due to the fact that hydro-morphodynamic models involve a series of physically unmeasurable parameters, their accuracy and reliability highly depend on the calibration process, which means the adjustment of uncertain input parameters to minimize the misfit between simulation results and corresponding physical measurements. Generally, the calibration of such models has been manually carried out through trial-and-error parameter adjustments based on the user's understanding of the model structure and features of the environmental system until a satisfactory agreement between simulated and measured values can be achieved. By having multiple parameters for manual calibration, the typical way is to consider each parameter separately for tuning by keeping the others constant. The procedure can then be repeated for the rest of the parameters one at a time. However, in reality, the combination of the best values of every single parameter may not result in the overall best fit. Thus, having a complicated model with several uncertain input parameters, the manual calibration method becomes cost- and time-consuming, involving a high degree of subjectivity. Nevertheless, employing optimization algorithms is an ingenious approach to the model fitting process. The application of automatic model calibration in different fields of environmental studies, such as hydrologic or groundwater models, has gained popularity over the last four decades; however, there is a considerable gap in applying automatic calibration in hydro-morphodynamic studies. This work aims to evaluate the efficiency and accuracy of automatic model calibration based on mathematical optimization, which can be an innovative practice to overcome the time-consuming and subjective nature of manual model calibration.
MATERIALS AND METHODS
Experimental data
The experimental data obtained by Kantoush (2008) are used to set up the numerical models and calibrate them against measured bed levels. As an outline of the experimental work, a reference rectangular shallow reservoir with 4.0 m of width, 6.0 m of length, and 0.3 m of depth was constructed at the Laboratory of Hydraulic Constructions (LCH) of the Swiss Federal Institute of Technology (EPFL). There were up- and downstream rectangular channels centrally connected to the reservoir with 0.25 m of width and 1.0 m of length. The walls were made of movable PVC plates, which could be adjusted to create different geometrical shapes. The water-sediment mixture was supplied from a mixing tank into the basin by gravity. The thickness of deposited sediments was measured by a mini echo sounder. Fine ground non-uniform walnut shells were used as sediment particles with a median diameter of d50 = 50 μm, a geometric standard deviation of particle size gradation of , a density of
= 1.5 g/cm3, and an average entrance concentration of C = 3.0 g/l. Further details regarding the sediment characteristics can be found in Table 1. The discharge rate was constant during different experiments with a value of Q = 7 l/s and a corresponding water depth of h = 0.2 m, which was adjusted by a flap gate at the outlet.
Sediment characteristics used for numerical modeling
. | Sediment size classes . | ||||||
---|---|---|---|---|---|---|---|
Size (mm) | 0.025 | 0.03 | 0.05 | 0.06 | 0.09 | 0.125 | 0.28 |
Proportion (%) | 20 | 10 | 20 | 10 | 20 | 10 | 10 |
Cumulative proportion (%) | 20 | 30 | 50 | 60 | 80 | 90 | 100 |
Fall velocity (mm/s) | 0.2 | 0.25 | 0.7 | 1 | 2 | 4.2 | 20 |
. | Sediment size classes . | ||||||
---|---|---|---|---|---|---|---|
Size (mm) | 0.025 | 0.03 | 0.05 | 0.06 | 0.09 | 0.125 | 0.28 |
Proportion (%) | 20 | 10 | 20 | 10 | 20 | 10 | 10 |
Cumulative proportion (%) | 20 | 30 | 50 | 60 | 80 | 90 | 100 |
Fall velocity (mm/s) | 0.2 | 0.25 | 0.7 | 1 | 2 | 4.2 | 20 |
Photos of the physical models and their sketched plan views for (a and b) hexagonal and (c and d) lozenge-shaped reservoirs.
Photos of the physical models and their sketched plan views for (a and b) hexagonal and (c and d) lozenge-shaped reservoirs.
Numerical modeling

























The sediment continuity equation (Exner equation) is used to calculate the bed elevation changes with respect to the continuity defect in cells close to the bed. The difference between the sediment inflow and outflow in a cell is multiplied by the time step and divided by the horizontal plane area of the cell. The bed is then raised or lowered after each time step, and the grid is regenerated.
Regarding the velocity distribution at the inlet, a minor linear perturbation is introduced to the velocity profile, described by Dewals et al. (2008), to provide an initiating factor and a condition for the genesis of the asymmetric flow in symmetric numerical models. This small perturbation can be considered as an initial boundary condition that is unavoidable in experimental set-ups, and can be further damped and become a straight jet regarding the reservoirs with stable-symmetric flow structure. In other words, a perfect symmetric numerical model needs a symmetry-breaking condition to reproduce an asymmetric flow field.
Model calibration
Model calibration is regarded as the initial stage of appraising the performance of a computational model that represents the complex behavior of a real-world system. This inverse process depends on how sensible the uncertain affecting parameters are adjusted according to the misfit minimization between calculated and corresponding measured values. The traditional manual trial-and-error method of model calibration has nowadays been overshadowed by the concept of mathematical optimization, which paves the way for automatic model calibration owing to its interesting features such as objective-based judgment of goodness-of-fit rather than relying on the user subjectivity and being less time-consuming (Evangelista et al. 2017). The principal elements of an automatic calibration routine are: an objective function to evaluate the model performance, an optimization algorithm to explore the parameter space by repetitive adjustment of the uncertain parameters, and a termination criterion to stop the search when the convergence of the objective function or the maximum allowable number of iterations is satisfied (Vidal et al. 2007).
In this study, the model-independent nonlinear Parameter ESTimation and predictive analysis package PEST (Doherty 2016), which employs the gradient-based Gauss–Marquardt–Levenberg (GML) local optimization algorithm, is used to calibrate the numerical models. This tool has demonstrated promising results for sensitivity analysis and automatic calibration of numerical models in different environmental studies (Shoarinezhad et al. 2020a, 2020b). Further, a global optimization algorithm, Big Bang Big Crunch (BB–BC) (Erol & Eksin 2006), is applied to validate PEST performance in finding the global optimum over the search space.
The numerical models are calibrated in accordance with the measured bed levels in different longitudinal and cross-sections, taking 8,600 and 16,500 points into account for the lozenge- and hexagon-shaped reservoirs, respectively. The residual sum of squares between measured and calculated bed levels, , is used as the optimization objective function for the automatic calibration of the models. The measured points were located in every 10 cm of cross-sections (61 cross-sections) with a lateral distance of about 1 cm from each other. Among the four phases of sediment feeding in the experiments, the measured bed levels from the first 1.5 hrs (considered a warm-up period) are used as the initial bed levels to simulate the following three phases. This is performed to provide an initial condition for the roughness calibration considering the surface friction development from the hydraulically smooth bed surface.
Gauss–Marquardt–Levenberg algorithm of PEST


Big Bang Big Crunch algorithm

Here, is the position vector of the
th candidate in an n-dimensional search space, N is the population size, and
denotes the fitness value of the
th candidate. In the next Big Bang, new individuals are mainly generated around the former-calculated center of mass according to a normal distribution, where the standard deviation of the normal distribution decreases as the optimization proceeds. This step is followed by a contraction according to the recalculation of the new center of mass. In order to ensure the global convergence of the method, the algorithm always generates a number of new solutions far from the center of mass with a diminishing probability as iterations go forward. Eventually, with regard to a defined termination criterion, this successive two-phase scheme (explosion-contraction) converges to the optimum point (Kaveh & Bakhshpoori 2019).
Parameters selection
A series of sensitivity analyses are performed to find out the significance of uncertain input parameters (e.g., roughness coefficient, the porosity of compacted bed sediments, the coefficient for the bedform smoothing algorithm, the thickness of the upper active sediment layer, and the angle of repose for sediments) on the system behavior as the initial stage prior to the model calibration. SENSAN, a model-independent local SENSitivity ANalyzer as one of the PEST utilities, is used for the sensitivity analysis. SENSAN conducts several runs according to the pre-defined sets of parameter values and records the model output sensitivities to parameter changes. Among the assessed parameters, the two most sensitive ones are selected for auto-calibration as follows:
Nikuradse equivalent roughness height (
), which is generally considered to be proportionately related to the representative grain size, such that
. There are several values for
and various sediment sizes as the representative grain size
in the literature (e.g., see Dey (2014)). Here, the range of this parameter is set to be d50 ≤
≤ 10d90.
Active layer thickness (ALT), as the superficial exchange layer depth, engaged in the entrainment and deposition of sediment particles. The sorting mechanism occurs in the active layer, where the sediment continuity equation is computed separately for each size fraction inside a cell during each time step. Depending on the transport regime, ALT scales with the representative grain size of the sediment mixture or is defined as a fraction of bedforms height. The range for this parameter is selected to be d50 ≤ ALT ≤ 5dmax (Malcherek 2007).
Since the GML algorithm is a gradient-based approach, it might likely find the local optimum point on the objective function surface rather than the desired global one. Therefore, it is worthwhile to reassess the calibration procedure by using different random initial values within the parameter space in the case of using local optimization algorithms. Here, two different pre-defined starting values are considered for the investigated parameters:
GML#1:
= d90 = 0.013 cm ALT1 = dmax = 0.028 cm
GML#2:
= 3d90 = 0.039 cm ALT2 = 3dmax = 0.084 cm
RESULTS AND DISCUSSION
The geometry of a reservoir plays a vital role in the flow field development and, consequently, the sedimentation pattern. Regarding this fact, different relations in the literature can be found, the so-called shape factors, which correlate geometrical aspects with the appearance of various flow field categories (e.g., Kantoush (2008); Dufresne et al. (2011)). Above a critical value for the shape factor, the symmetric flow inside a symmetric reservoir evolves toward the asymmetric pattern due to the Coanda effect, where the flow deviates to one side of the reservoir according to a tiny imperfection in the physical symmetric configuration.
Contour maps of calculated velocity fields for the two reservoirs in different periods: 1.5 hrs of (a and d) phase 1 and (b and e) phase 2 for both reservoirs; (c) 4.5 hrs and (f) 3 hrs of phase 3 for the hexagonal and lozenge-shaped reservoirs, respectively.
Contour maps of calculated velocity fields for the two reservoirs in different periods: 1.5 hrs of (a and d) phase 1 and (b and e) phase 2 for both reservoirs; (c) 4.5 hrs and (f) 3 hrs of phase 3 for the hexagonal and lozenge-shaped reservoirs, respectively.
Central cross-sectional view of the final velocity field for (a) hexagonal and (b) lozenge-shaped reservoirs.
Central cross-sectional view of the final velocity field for (a) hexagonal and (b) lozenge-shaped reservoirs.
The auto-calibration process is based on the pairwise comparison of the calculated and measured bed levels. Table 2 shows the calibration results for both reservoirs using GML (with two different starting values for the selected parameters) and BB–BC (with randomly sampled initial values) algorithms. It can be seen that the calibrated values for the investigated parameters (calibrated and ALT calibrated) are very similar for GML#1 and GML#2 regarding both reservoirs, indicating the repeatability and robustness of the auto-calibration procedure based on the GML algorithm. Moreover, these values are almost the same as the results of the BB–BC algorithm. It means the gradient-based GML algorithm is not affected by local minima, and PEST can reliably calibrate the numerical models in a global manner. Furthermore, among innumerable parameter combinations, PEST calibrates the models with 30–40 runs depending on the starting values, which shows the method's efficiency compared to the BB–BC algorithm.
Initial values of the investigated parameters and calibration results (final calibrated values and the number of model runs)
Reservoir shape . | Algorithm . | Initial values . | Calibration results . | |||
---|---|---|---|---|---|---|
ksinitial . | ALT initial . | kscalibrated . | ALT calibrated . | Model runs . | ||
Lozenge | GML#1 | 0.013 | 0.028 | 0.0212 | 0.0308 | 32 |
GML#2 | 0.039 | 0.084 | 0.0214 | 0.0311 | 38 | |
BB–BC | Randomly sampled | 0.0212 | 0.0309 | 471 | ||
Hexagon | GML#1 | 0.013 | 0.028 | 0.0238 | 0.0382 | 31 |
GML#2 | 0.039 | 0.084 | 0.0235 | 0.0377 | 35 | |
BB–BC | Randomly sampled | 0.0237 | 0.0380 | 438 |
Reservoir shape . | Algorithm . | Initial values . | Calibration results . | |||
---|---|---|---|---|---|---|
ksinitial . | ALT initial . | kscalibrated . | ALT calibrated . | Model runs . | ||
Lozenge | GML#1 | 0.013 | 0.028 | 0.0212 | 0.0308 | 32 |
GML#2 | 0.039 | 0.084 | 0.0214 | 0.0311 | 38 | |
BB–BC | Randomly sampled | 0.0212 | 0.0309 | 471 | ||
Hexagon | GML#1 | 0.013 | 0.028 | 0.0238 | 0.0382 | 31 |
GML#2 | 0.039 | 0.084 | 0.0235 | 0.0377 | 35 | |
BB–BC | Randomly sampled | 0.0237 | 0.0380 | 438 |
Note: , roughness height at the bed; ALT, active layer thickness. Units are in centimeters.
Considering d90 = 0.013 cm as the representative grain size, the calibrated values can be rewritten as: ks ≈ 1.64d90 and ALT ≈ 2.4d90 for the lozenge-shaped reservoir; and ks ≈ 1.84d90 and ALT ≈ 2.94d90 regarding the hexagonal reservoir.
Comparison of bed level changes in the experiment (left) with the simulation results (right) in (a and d) phase 1, (b and e) phase 2, and (c and f) phase 3 for the hexagonal reservoir.
Comparison of bed level changes in the experiment (left) with the simulation results (right) in (a and d) phase 1, (b and e) phase 2, and (c and f) phase 3 for the hexagonal reservoir.
Comparison of bed level changes in the experiment (left) with the simulation results (right) in (a and d) phase 1, (b and e) phase 2, and (c and f) phase 3 for the lozenge-shaped reservoir.
Comparison of bed level changes in the experiment (left) with the simulation results (right) in (a and d) phase 1, (b and e) phase 2, and (c and f) phase 3 for the lozenge-shaped reservoir.
Despite the stable, straight, and symmetric nature of the flow in the hexagonal reservoir, which results in a central-longitudinal pattern of sediment depositions, the flow structure in the lozenge-shaped reservoir is unstable and fluctuating with main sediment deposits at the sides of the reservoir (Figure 5). During the first phase, the flow direction is clockwise, and sediments are mostly settled at the left part, resulting in an asymmetric bed topography (Figure 5(a) and 5(d)). The shift in the flow path during the second phase, due to the deposited sediments that can modify the unstable and sensitive flow structure in the lozenge-shaped reservoir, changes the bed topography to a semi-symmetric pattern (Figure 5(b) and 5(e)). Here, the numerical model underestimates the bed levels at the second half of the right side of the reservoir. Considering the final bed levels, sediments are almost symmetrically deposited along the left and right parts of the basin in the experiment, whereas the maximum deposition is along the first half of the right side in the numerical model (Figure 5(c) and 5(f)). Furthermore, in all phases, an erosion area can be seen just in front of the inlet in the experiment, which is deeper and larger in size compared to the results of the numerical model.
In order to quantitatively evaluate the overall performance of the automatically calibrated numerical models, various statistical metrics are used to compare the final measured and calculated bed levels, as shown in Table 3. The results of uncalibrated models using the initial parameter values of GML#1 are also presented to see the deviation between the results of the calibrated models and our initial guess. Mean Bias Error (MBE) is applied as a bias indicator, describing the degree of underprediction (negative values) or overprediction (positive values) of the model. Root Mean Squared Error (RMSE) and Mean Absolute Error (MAE) both reflect the average magnitude of the error. While MAE indicates a linear behavior of individual errors, RMSE gives more importance to big errors by giving higher weights to them. Pearson correlation coefficient (R) shows the linear correlation between estimated and measured bed levels. As a multi-component goodness-of-fit, Kling-Gupta efficiency (KG) is used, which combines Pearson correlation coefficient, bias, and variability within a single objective function.
Statistical performance of the automatic model calibration
Reservoir shape . | . | Goodness-of-fit . | ||||
---|---|---|---|---|---|---|
. | MBE (cm) . | RMSE (cm) . | MAE (cm) . | R (-) . | KG (-) . | |
Lozenge | Calibrated | −0.024 | 0.054 | 0.044 | 0.80 | 0.77 |
Initial guess | −0.061 | 0.288 | 0.212 | 0.56 | 0.52 | |
Hexagon | Calibrated | −0.013 | 0.055 | 0.039 | 0.86 | 0.85 |
Initial guess | −0.039 | 0.252 | 0.208 | 0.70 | 0.68 |
Reservoir shape . | . | Goodness-of-fit . | ||||
---|---|---|---|---|---|---|
. | MBE (cm) . | RMSE (cm) . | MAE (cm) . | R (-) . | KG (-) . | |
Lozenge | Calibrated | −0.024 | 0.054 | 0.044 | 0.80 | 0.77 |
Initial guess | −0.061 | 0.288 | 0.212 | 0.56 | 0.52 | |
Hexagon | Calibrated | −0.013 | 0.055 | 0.039 | 0.86 | 0.85 |
Initial guess | −0.039 | 0.252 | 0.208 | 0.70 | 0.68 |
Note: MBE, Mean Bias Error; RMSE, Root Mean Squared Error; MAE, Mean Absolute Error; R, Pearson correlation coefficient; KG, Kling-Gupta efficiency.
According to the negative values of MBE, the numerical models underestimate the final bed levels in both cases. Although having a similar calibrated RMSE, the lower MAE value for the hexagonal model represents its more precise calibration (the best value is 0). This can further be confirmed by higher values of correlation and Kling-Gupta efficiency metrics for the hexagonal model (the best value is 1). Nevertheless, the statistical performance of the calibrated lozenge-shaped model also agrees reasonably with the measured data.
SUMMARY AND CONCLUSIONS
In this study, we applied an automatic calibration tool (PEST), which uses the gradient-based Gauss–Marquardt–Levenberg optimization algorithm, to calibrate the 3D morphodynamic numerical models of two reservoirs (hexagonal and lozenge-shaped configurations) against the experimental data. The two most affecting input parameters (roughness height and active layer thickness) in the numerical models are selected for calibration through a sensitivity analysis. In order to verify the ability of the gradient-based optimization algorithm to find the global optimum values of the parameters over the search space rather than sticking to a local optimum point (evaluating the robustness and convergence of the algorithm), in addition to using the global optimization algorithm Big Bang Big Crunch, we tested the gradient-based method with two different pre-defined initial values with a reasonable range based on the literature. Achieving an almost identical set of calibrated values confirmed the accuracy and reliability of the optimization procedure. The parameters are calibrated with 30–40 model runs by PEST, which shows its efficiency and superiority over the trial-and-error manual calibration, suggesting its potential use for hydro-morphodynamic models calibration.
So far, most related research works have been on rectangular reservoirs with a 90° expansion angle. In this work, the effect of lower expansion angles on the flow field development and the sedimentation pattern is assessed numerically. Keeping the maximum inner dimensions constant, the lower expansion angle gives a higher chance of asymmetric flow appearance. Regarding the hexagonal reservoir, which has a stable flow structure, the main part of sediment particles is settled along the mid-longitudinal section. However, the unstable flow pattern inside the lozenge-shaped reservoir causes sediments to be distributed at the sides of the basin. The calculated bed levels are compared with the measured topography of the physical models in different time steps. Considering the primary features and specifications of both reservoirs regarding the flow field, jet direction, recirculation zones, and bed topography in the experimental work, the calibrated numerical models can reasonably reproduce similar patterns. The achieved results regarding the flow fields and erosion/deposition patterns are used as a criterion to judge the ability of the Gauss–Marquardt–Levenberg search algorithm to calibrate the numerical models (optimization performance). A difference between the bed levels in the physical set-up and simulation results is related to the scour hole in front of the inlet, likely due to boundary effects, which cannot be predicted by the numerical model.
The quality and performance of the calibrated numerical models are also investigated by different statistical metrics, comparing the predicted and measured bed levels. The mean bias error shows an underestimation of bed levels in both reservoirs. Nevertheless, according to the low values of the root mean squared error and the mean absolute error; and high values of the correlation coefficient and Kling-Gupta efficiency, the overall performance of the automatic calibration procedure is reasonable.
The manual trial-and-error model calibration approach with just one parameter subject to alteration may be reasonable and sufficient in most cases. However, if there are several uncertain parameters, their innumerable combinations cause the manual calibration method to become much more complex, time- and cost-consuming, and impractical. What is more, since the issue of subjectivity is involved in manual model calibration, there is no guarantee that the best possible combination of parameters can be achieved. The overall outcome of this study is that using suitable optimization algorithms for hydro-morphodynamic models calibration, which is not a common practice among researchers in this field, can considerably reduce the calibration time and user intervention/subjectivity and concurrently increase the precision of the process.
ACKNOWLEDGEMENTS
The first author received a Ph.D. scholarship from the German Federal Ministry of Education and Research (BMBF) through the Sustainable Water Management-NaWaM program of the German Academic Exchange Service (DAAD). The last author is indebted to the Baden-Württemberg Stiftung for the financial support from the Elite program for Postdocs. The authors also acknowledge Prof. Nils Reidar Bøe Olsen, the developer of SSIIM software, for his valuable suggestion regarding numerical modeling.
DATA AVAILABILITY STATEMENT
All relevant data are included in the paper or its Supplementary Information.
CONFLICT OF INTEREST
The authors declare there is no conflict.