Open Journal of Optimization
Vol.07 No.04(2018), Article ID:90154,34 pages

An Overview of Recently Developed Coupled Simulation Optimization Approaches for Reliability Based Minimum Cost Design of Water Retaining Structures

Muqdad Al-Juboori1,2 , Bithin Datta3

1College of Science and Engineering, James Cook University, Townsville, Australia

2College of Engineering, Wasit University, Wasit, Iraq

3Discipline of Civil Engineering, College of Science and Engineering, James Cook University, Townsville, Australia

Copyright © 2018 by author(s) and Scientific Research Publishing Inc.

This work is licensed under the Creative Commons Attribution International License (CC BY 4.0).

Received: December 13, 2018; Accepted: December 28, 2018; Published: December 31, 2018


This paper reviews several recently-developed techniques for the minimum-cost optimal design of water-retaining structures (WRSs), which integrate the effects of seepage. These include the incorporation of uncertainty in heterogeneous soil parameter estimates and quantification of reliability. This review is limited to methods based on coupled simulation-optimization (S-O) models. In this context, the design of WRSs is mainly affected by hydraulic design variables such as seepage quantities, which are difficult to determine from closed-form solutions or approximation theories. An S-O model is built by integrating numerical seepage modeling responses to an optimization algorithm based on efficient surrogate models. The surrogate models (meta-models) are trained on simulated data obtained from finite element numerical code solutions. The proposed methodology is applied using several machine learning techniques and optimization solvers to optimize the design of WRS by incorporating different design variables and boundary conditions. Additionally, the effects of several scenarios of flow domain hydraulic conductivity are integrated into the S-O model. Also, reliability based optimum design concepts are incorporated in the S-O model to quantify uncertainty in seepage quantities due to uncertainty in hydraulic conductivity estimates. We can conclude that the S-O model can efficiently optimize WRS designs. The ANN, SVM, and GPR machine learning technique-based surrogate models are efficiently and expeditiously incorporated into the S-O models to imitate the numerical responses of simulations of various problems.


Linked Simulation-Optimization, Water-Retaining Structures, Machine Learning Technique, Reliability Based Optimum Design, Multi-Realization Optimization Model, Heterogeneous Hydraulic Conductivity

1. Introduction

Construction of water-retaining structures (WRSs) [1] [2] , such as dams, barrages, regulators and weirs, is essential for stable and safe water management and to generate clean energy. However, significant considerations and hazards must be recognized in their design, such as economic factors and the risk of failure. Accordingly, the design and analysis of such structures require precise estimation and understanding of relevant variables and parameters, especially those related to seepage quantities, and their impacts on WRS safety. This study presents a coupled simulation-optimization (S-O) approach to identifying minimum-cost WRS designs while incorporating numerical seepage analysis, uncertainty in seepage quantities, and hydraulic design safety factors.

Hydraulic structures that impound a considerable amount of water (head) and are constructed on permeable soil foundations are associated with water seepage. Seepage forces threaten the hydraulic efficiency and structural stability of such structures. Another consequence of seeping water is pore-water pressure, which applies uplift (upward) pressure on the structure’s floor (apron) and may cause it to collapse.

The mathematical relationship between seepage design variables with flow domain characteristics is complex and nonlinear. An analytical solution may be obtained only for simple and symmetrical cases. Many approximation and empirical theories have been proposed for estimating seepage quantities (particularly uplift pressure and exit gradient). These theories include Bligh’s creep theory, Lane’s weighted creep theory, the flow-net method, the fragment method and Khosla’s theory [3] [4] . Their applications are associated with significant amounts of error compared to applications that use analytical solutions or experimental modelling. Additionally, the analytical methods and approximation theories apply to ideal general soil conditions (homogeneous and isotropic), which are rarely found in real-life cases [5] .

Recently, as a result of the development of numerical methods and computerized processes, many seepage problems related to WRS have been accurately simulated and solved by numerical methods such as the finite element method (FEM). Furthermore, hydraulic conductivity variation and other soil parameters can be integrated into a numerical model to study the consequences of varying soil parameters on WRS design. However, numerical techniques only provide solutions when certain parameters of the hydraulic structure are known, including the boundary conditions and geometry of the flow domain. This means that numerical models do not provide generalized performance equations similar to those provided by analytical solutions.

Considering the above-mentioned arguments, the contradicting goals of efficiency, accuracy, safety and cost must be simultaneously integrated into the design process to attain optimum, safe and economical WRSs. Hence, the optimization approach can be used to identify optimal WRS designs. As a result, the minimum cost and safest WRSs can be achieved. Directly coupling a numerical model to an optimization model to attain an optimum WRS design is computationally inefficient and time-consuming. Alternatively, the numerical model could be replaced by an approximate machine learning model (surrogate model) that accurately and expeditiously imitates the numerical model’s responses. Such surrogate models (meta-models) can be trained based on numerically simulated data (input and output) sets.

Linked simulation optimization (S-O) models based on surrogated models have been applied to solving groundwater and water resource problems and other applications related to groundwater management and identification of contaminant sources [6] - [18] . Specifically, in WRS design involving seepage effects, few studies have utilized S-O techniques. Singh [19] [20] formulated an optimization model to find the optimum dimensions of a barrage at minimum cost. The author used Khosla’s theory to obtain seepage characteristics that were processed by an optimization algorithm. Kholsa’s theory is only applicable to small hydraulic structures and the resulting solutions have a noticeable amount of error. Al-Suhaili and Karim [21] implemented an indirect S-O model based on an artificial neural network (ANN) model to find optimum solutions for hydraulic structures at minimum cost. In their study, WRS safety factors were only considered in terms of exit gradient and uplift pressure, while sliding, overturning and eccentric load effects were ignored. Also, the ranges of the implemented cases were only applicable to small WRSs (total head less than 10 m). Accordingly, there is limited potential to apply S-O models to the hydraulic design of WRSs. Additionally, no studies have yet considered the effects of uncertainty in some of the design parameters related to seepage analysis.

The current paper presents a few innovative applications of a newly developed S-O technique that can determine optimal WRS designs (safe, reliable and economic) based on adequately trained meta-models. Induced seepage forces, and many safety factors and design requirements related to WRSs such as overturning, sliding safety factors, and prevention of eccentric load conditions, are considered in the S-O approaches. Also, uncertainty in seepage characteristic estimates due to uncertainty in hydraulic conductivity is integrated into the reliability based optimum design using S-O model. For each S-O model, the type of machine learning technique and optimization solver are selected based on prediction accuracy and efficiency.

In the following sections, different scenarios and applications of linked S-O approaches were used to find the optimum WRS design. Section two encompasses the linked S-O models for simple WRS including two end cut-offs and apron. The coupled S-O models in section three is used for comprehensive (complex) WRS incorporating 10 cut-offs and stratified flow domain based on different hydraulic conductivity values. The fourth section incorporates the reliability concept in the S-O models using RBOD framework to address uncertainty in estimate the seepage characteristics due to uncertainty of hydraulic conductivity. The fifths section includes more advance application of RBOD using multi-objective multi-realization optimization (MOMRO) models.

2. Application of the S-O Approach to a Simple Conceptual Seepage Model Related to WRS Design

2.1. Conceptual Seepage Model and Data Generation

A simple conceptual seepage model was proposed to illustrate the WRS design problem [22] , as shown in Figure 1. This model uses fixed positions and variable depths for the cut-offs at each end of a WRS (upstream, downstream) and a variable apron length between the cut-offs. The proposed training ranges of the input design variables were: 1 - 40 for cut-off depths at the upstream (d1) and downstream side (d2); 1 - 60 for upstream water head (H) and half the WRS apron width (b).

Training of meta-models is based on datasets simulated by a numerical seepage modeling code (Geo-Studio/SEEP/W) [23] . Five hundred input data points are used as independent variables (d1, d2, 2b, H) and are randomly generated using the Latin hypercube sampling (LHS) method [24] . The output data is obtained

Figure 1. Conceptual seepage model.

from numerical seepage modeling for each input set. The most important output data are uplift pressure on the floor at the upstream cut-off (PC), uplift pressure on the floor at the downstream cut-off (PE) and the exit gradient (ie) at the toe of the WRS.

The ANN technique is used to develop the meta-models in this section. ANN can explore complex, discontinuous and nonlinear relationships between datasets [25] . For seepage and groundwater problems related to hydraulic structures, ANN has been utilized to simulate and identify seepage characteristics [21] [26] - [32] . Based on the 500 generated data points related to the WRS seepage system, ANN was used to build three meta-models. These models provide accurate predictions of seepage characteristics without utilizing further numerical seepage simulation, i.e., Geo-Studio/SEEP/W code. To select the best parameters for the ANN models and provide accurate predictions, many designs of experimental models were implemented using the Taguchi method [33] to find the ideal ANN parameter [34] combination for each surrogate model of each seepage characteristic. The final PC model had 11 neurons, 60/40 training-to-validation ratio, logsig transformation function for the hidden layer, and purelin transformation function for the output layer. The final PE model had four neurons, 50/50 training-to-validation ratio, logsig transformation function for the hidden layer, and purelin transformation function for the output layer. Similarly, the final ie model had five neurons, 50/50 training-to-validation ratio, logsig transformation function for the hidden layer, and tansig transformation function for the output layer.

2.2. Formulation of the Optimization Model

The optimization model was formulated as follows:

Find X = { x 1 , x 2 , x 3 , x 4 } = { d 1 , d 2 , b , b * }


f ( X ) = c f b ( t 1 + t 2 ) 2 + t c s = 1 2 c s d s (1)

Subject to:

F S i e = ε ( H , d 1 , d 2 , b , i e c r t )

F S i e F S e x i t (2)

F s f l u s = ϵ ( H , d 1 , d 2 , b )

F s f l u s F S u p l i f t (3)

F s f l d s = γ ( H , d 1 , d 2 , b )

F s f l d s F S u p l i f t (4)

E c c b 3 (5)

E c c 2 b 3 (6)

M p a s = f ( H , b , b * , t 1 , t 2 , G c , G w , P C , P E )

M a c t = f ( H , b , b * , t 1 , t 2 , G c , G w , P C , P E )

V l = f ( H , b , b * , t 1 , t 2 , G c , G w , P C , P E )

H l = f ( H , G w )

E c c = M p a s M a c t V l o a d

F s o v e r 1.5 (7)

F s o v e r = M p a s M a c t

F s s l i d 1.5 (8)

F s s l i d = C × b + f × V l H l

where t 1 , t 2 represent the floor thicknesses at the upstream and downstream sides, respectively. c f is the construction cost of the floor per cubic meter ($400/m3), c c is the construction cost of the cut-offs per cubic meter, which is a function of the cut-off depth (d1, d2) as shown in Equation (9) and Equation (10), and t c is the cut-off thickness, which equals 1.0 m.

C 1 = x 1 3 + 20 x 1 2 + 200 x 1 + 400 (9)

C 2 = x 2 3 + 20 x 2 2 + 200 x 2 + 400 (10)

F s i e is exit gradient safety factor determined by the meta-models {ε( )}. F s f l u s , F s f l d s are safety factors imposing weight on the upstream and downstream floors of the WRS to safely counterbalance the uplift pressures (PC, PE) [35] [36] . Computing of F s f l u s , F F s f l d s is mainly based on the developed meta-models { ϵ ( ) } , { γ ( ) } , respectively. E c c is an eccentric distance (condition) to prevent eccentric loading on the WRS foundation. M p a s is the passive momentum obtained from all forces, which increases the stability of the WRS, M a c t is the active momentum obtained from all forces, which decreases the stability of the WRS, and V l o a d is the resultant of all vertical loads influencing the WRS. M p a s , M a c t , V l o a d are determined based on H , b , b * , t 1 , t 2 , G c , G w , P C and P E . F s o v e r is the overturning safety factor and F s s l i d is the sliding safety factor. C is the soil cohesion resistance, and f = tan , where is the internal friction angle [5] . The values of f and C are assumed to be f = tan = 0.7 and C = 20 kPa. Hl is the resultant of all horizontal loads affecting the WRS.

2.3. Results and Discussion

Simulation-Optimization Model

The S-O technique was implemented for different H values ranging from 2 m to 40 m. The design requirement of the WRS and all the constraints were satisfied for each optimum solution. The optimum solutions for the different H values show that d1 and d2 make considerable contributions to the hydraulic safety of a WRS (see Figure 2 and Figure 3). Nonetheless, the length of d2 is more important than that of d1, because d2 has a substantial impact on the ie value, which is the critical factor in the hydraulic design of WRSs [37] [38] .

The optimum width (2b) of a WRS effectively influences its optimum hydraulic design, because the total width is directly integrated into the specified safety factors. Reduction of the optimum WRS width (2b) could be attributed to increasing floor thicknesses at the upstream and downstream sides (t1, t2) with increases in head, which results in an expensive design. Therefore, the S-O model reduces the total width and simultaneously augments the depths of d1 and d2, which is an efficient and cost-effective solution to reducing the tremendous uplift pressure and exit gradient effects experienced at large H values. Moreover, b* considerably contributes to the optimum hydraulic design of the WRS. The weight of water head above b* counterbalances the uplift pressure and enhances the stability of the WRS. This means that the value of b* substantially contributes to WRS safety and provides a cheaper solution.

3. Coupled Simulation-Optimization Technique for Optimum WRS Hydraulic Design Incorporating Complex Seepage Scenarios

Often, the seepage prevention components of WRS projects are the end cut-offs

Figure 2. Optimum solutions (d1, d2) for different head values.

Figure 3. Optimum solutions (2b, b*) for different head values.

(upstream and downstream) and the apron between them. Limitations in the orientation, length, and number of cut-offs, and the width of the apron reduce the opportunity of finding a feasible optimum solution using the linked S-O technique. Including the effects of different scenarios of hydraulic conductivity and anisotropic ratios on WRSs is an important concept that must be considered in their optimum design. Moreover, studying soil stratification based on different values of hydraulic conductivity and quantifying its effects on optimum WRS design is another concept that needs to be considered.

Hence, a comprehensive conceptual model is proposed, as shown in Figure 4. This model includes ten cut-offs distributed along the apron of the WRS. The cut-off lengths, orientations, and distance between them (the apron) are considered as variables. These variables are used to build meta-models and, within the S-O model, their optimum values can be determined. Based on the optimum solutions for safe, minimum-cost WRS designs that consider seepage impacts, the most important and active sets of optimum values can be identified. This section investigates the effects of soil properties on the optimum solutions and identifies the most important and effective seepage control components for optimum WRS design.

The geometry of the assumed numerical model comprised ten cut-offs (sheet piles) with varied positions, lengths, and orientations. Additionally, three subsoil layers were assumed and the principle (horizontal) hydraulic conductivity (kx) and anisotropic ratio (ky/kx) were varied for each layer and each case. As a result, the contribution of each variable of the comprehensive model to the optimal design can be explored for various boundary conditions.

The prescribed range of each design variable and parameter is shown in Table 1.

Figure 4. Seepage conceptual model scheme.

Table 1. Summary of the variables used in the model.

The ranges of hydraulic conductivity and anisotropic ratios used cover a wide range of typical real-life values [39] [40] [41] [42] . The 41 input variables included in the conceptual model were: total upstream water head (H), ten cut-off depths ( d 1 , d 2 , , d 10 ), their angles ( β 1 , β 2 , , β 10 ), the distances (widths) between cut-offs ( b 1 , b 2 , , b 10 ), three subsoil layer depths (LD1, LD2, LD3), their hydraulic conductivities in the horizontal direction ( k x 1 , k x 2 , k x 3 ), and their anisotropic ratios (ky/kx)1, (ky/kx)2, (ky/kx)3, respectively. The most important seepage design characteristics for each numerical seepage model were the uplift pressures in front of (PEi) and behind (PCi) each single cut-off ( S 1 , S 2 , , S 10 ) in addition to the exit gradient (ie) at the toe of the hydraulic structure. Hence, 21 meta-models were required, one for each seepage characteristic.

3.1. Support Vector Machine Surrogate Model

The support vector machine (SVM) is one of the most popular machine learning techniques and has recently been used for various nonlinear and complex engineering problems [43] - [54] . The Matlab programing language was utilized to develop the meta-models used in the present study because it is a versatile tool that has many options that can be modified to build efficient SVM-based meta-models. Twenty-one models were built to determine the uplift pressures in front of and behind each cut-off, and the exit gradient near the toe of the WRS. These models were trained on 1500 scenarios of numerically simulated data. The input variables were randomly generated using LHS, then numerically simulated to determine seepage characteristics (output variables).

For each uplift pressure dependent variable, two different SVM models were built and trained on different training/testing datasets randomly selected from source data. A basic version of the ensemble surrogate model based on an average of the two models was developed. This procedure provides a more robust and accurate prediction. Also, any uncertainty arising from source data or meta-model prediction can be reduced by using the ensemble surrogate model. For exit gradient, three different models were developed. The coefficient of determination (RSQ) and mean square error (MSE) for the training and testing phases were used as measures of the predictive accuracy of the meta-models based on detached data. The parameters of each SVR model were carefully selected after performing trial-and-error iterations until the best RSQ and MSE values were obtained.

3.2. Optimization Model

The optimization model described in this section includes 31 decision variables ( b 1 , b 2 , , b 11 , d 1 , d 2 , , β 1 , β 2 , , β 10 ) and several constraints. Also, the optimization solver evaluates the objective function and constraint values based on 21 ensemble surrogate model responses. This makes the optimization problem a complex task. Safety factors and other hydraulic design requirements represent the imposed constraints of the optimization model within the S-O model. The best value of each design/decision variable was selected by the optimization algorithm to provide a safe and economic design. For this complex optimization task, the hybrid genetic algorithm (HGA) was used. The HGA is a combination of two optimization algorithms: GA and the interior point algorithm (IPA) [55] [56] [57] . The HGA provides a global optimum solution and can deal with complex problems.

Formulation of the optimization model in this section is similar to that in Section 2. However, because many cutoffs were used, flotation constraints (safety factors) were applied at all key uplift pressure points associated with different cut-off locations. Moreover, the construction cost incorporates the effect of inclination angle. Hence, Equation (11) is used in determining the construction cost of the cut-offs.

C c = 0.05 d i 2 + 200 d i 2 + 0.0698 β i 2 12.558 β i + 565.93 i (11)

3.3. Results and Discussion

Head Variation Effects

The linked S-O was implemented for different head values ranging from 20 m to 100 m. Other parameters were kept constant, such as hydraulic conductivity for all layers (kx = 5 m/day), anisotropic ratio ((ky/kx)1 = 1) and the depth of soil layers (50 m). The obtained optimum solutions can determine the vital values of all the provided design (decision) variables. This means that the optimization solver selects the design variables that provide the safest and most cost-efficient WRS design for the optimum solution.

In general, the optimum solutions obtained demonstrate that the contributions of variables b1 to b8 and d2 to d8 to WRS safety are insignificant. The optimum values for these variables approached zero. In contrast, values b9, b10, b11, d9 and d10 had vital roles in optimum WRS design, as shown in Figure 5 and Figure 6. These variables, for most implemented cases, had considerable values which varied according to head values. The function of d9 is to reduce PC9 and PE10 uplift pressure and exit gradient value. More importantly, the function of d10 is to directly reduce the exit gradient value, which is the most critical seepage characteristic. The functions of b9 and b10 are to provide sufficient weight to stabilize the WRS, reduce uplift pressure, and provide sufficient width to counterbalance overturning and sliding forces.

Other important variables were β9 and β10, which are related to d9 and d10, respectively. The optimum value of β10 is 150˚. This is logical, as making the inclination angle of the last cut-off toward downstream (greater than 90˚) substantially decreases the exit gradient value. This can be attributed to the augmentation of

Figure 5. Optimum width between cut-offs for the implemented cases according to head values.

Figure 6. Optimum cut-off depths for the implemented cases according to head values.

the streamline length of seeping water, particularly when β10 reaches 150˚. Thus, the duration and distance of seeping water would increase, which can reduce the exit gradient value. The optimum value of β9 was 30˚ in all implemented cases. Such an inclination angle can reduce uplift pressure under b10. This helps decrease the construction costs of WRSs.

In general, it seems that effective and optimum WRS designs must include two upstream and downstream cut-offs, the width (b10) between them, and the width (b9) on the upstream side. These widths are necessary to provide sufficient weight for the WRS to resist external hydrostatic loads and uplift pressure and play a vital role in satisfying WRS design requirements (constraints) such as sliding, overturning and eccentric load conditions.

All the optimum solutions satisfied the safety factors and requirements of WRS design. For all implemented cases, the optimum solution attained the minimum allowable value of the exit gradient safety factor, which is five. This reflects the significance of the exit gradient safety factor in WRS design and how it affects the construction cost because it is mainly controlled by the depth and inclination angle of the last cut-off (d10, β10), which are indispensable and expensive components used to reduce the exit gradient impact.

3.4. Hydraulic Conductivity (kx1) and Anisotropic Ratio (ky/kx)1 Effects

The same procedure which was applied to study the effects of upstream water head was implemented to quantify the effects of the hydraulic conductivity (kx1) and anisotropic ratio (ky/kx)1 of the first layer. The first layer is the nearest layer to the foundation of the WRS (Figure 4). The effect of (ky/kx)1 was studied by assuming eight different values ranging from 0.1 to 1.5. Ten values of kx1 ranging from 0.01 m/day to 20 m/day were specified and processed using the S-O technique. The values of the other design variables and parameters were left constant.

Generally, the optimum solutions obtained revealed that increases in (ky)1 and the (ky/ky)1 ratio significantly decrease the total cost of the WRS, as shown in Figure 7 and Figure 8. The reason for this is that when kx1 increases with a constant anisotropic ratio ((ky/kx)1 = 1), seeping water can easily move from the high-pressure zone (upstream) to the low-pressure zone (downstream). Consequently, the pore-water pressure underneath the WRS and the exit gradient values decrease. Thus, great cut-off depths and significant widths between cut-offs are unnecessary. Similarly, when the anisotropic ratio (ky/kx)1 is large with a specified hydraulic conductivity (kx1 = 5), the seeping water motion in the vertical direction becomes faster and the exit gradient becomes smaller compared to that

Figure 7. Minimum-cost optimum WRS design for different values of kx1.

Figure 8. Minimum-cost optimum WRS design for different values of (ky/kx)1.

obtained for small (ky/kx)1 ratios. Hence, for high values of (kx)1 and (ky/kx)1, the optimum values of d9, d10, b9, and b10, which are the most effective variables, decrease and, consequently, the optimal cost declines.

For small values of (ky/kx)1 and kx1, the exit gradient safety factor and safe eccentric distance play crucial roles in the optimum solution, compared to the other safety factors. This is evident as these safety factors reached the maximum or minimum allowable limit to satisfy the design requirements, while the optimum design attained the minimum construction cost. Table 2 and Table 3 demonstrate the safety factor variations for different values of (ky/kx)1 and kx1, respectively.

The S-O results demonstrate that the contributions of b1 to b8 and d1 to d8 to WRS safety are insignificant because their optimum values approached zero. Optimum solutions for WRS were based on increasing the values of b9 and b10 to counterbalance the uplift pressures, and on augmenting d9, d10 and β10 to decrease the exit gradient, as shown in Figures 9-12. Also, there is a significant contribution of d9 associated with the minimum value of β9 required to decrease uplift pressure beneath b10, which represents a large portion of the WRS floor.

Table 2. Safety factors for the implemented cases for different values of kx1.

Table 3. Safety factors for the implemented cases for different values of (ky/kx)1.

Figure 9. Optimum widths between WRS cut-offs for various values of kx1.

Figure 10. Optimum WRS cut-off depths for various values of kx1.

Figure 11. Optimum width between WRS cut-offs for various values of (ky/kx)1.

Figure 12. Optimum cut-off depths for various values of (ky/kx)1.

In particular, the optimization solver increased the d10 and β10 values to achieve a safe exit gradient, even though it is a more expensive option (Equation (11)). These variables were more effective at reducing the exit gradient, which is the most critical safety factor. Hence, the optimum value of β10 was 150˚, which is the maximum specified limit for this variable. For the same reason, the inclination angle of the upstream cut-offs (β9) approached the minimum allowable limit (30˚) for all implemented cases.

Simultaneously, to corroborate WRS stability and satisfy related safety factors, the required optimum WRS width was provided by b9 and b10. Furthermore, the uplift pressure on the downstream side decreased with increases in total width, which contributed to reductions in exit gradient. Therefore, the values of b9 and b10 mainly provide an efficient cross-section and weight to resist external loads and uplift pressures, and partially reduce the uplift pressure and exit gradient.

4. Reliability Based Optimum Design of WRSs Constructed on Heterogeneous Porous Media

Soil properties, and especially hydraulic conductivity, have large covariance values of 200% - 300%, which means the level of uncertainty in hydraulic conductivity is high [58] . Therefore, many studies have been conducted to study the effects of uncertainty and variations in soil properties on the reliability of designs [59] - [66] . Specifically, for groundwater and seepage for hydraulic structures, most studies have concentrated on stochastic analyses of seepage characteristics based on different realizations of hydraulic conductivity generated from different probability distribution functions (PDF) or sets of means and standard deviations [61] [67] - [72] . The important conclusion of such studies is that the degree of uncertainty drastically influences the seepage characteristics, which may negatively affect a design’s performance and safety.

This section concentrates on developing a reliability based optimum design (RBOD) framework to find optimum WRS designs at minimum cost and with a particular level of reliability while considering uncertainty in hydraulic conductivity and seepage. This can be achieved by formulating a constrained multi-realization optimization model linked with an S-O technique utilizing a GA optimization solver and incorporating many stochastic ensemble GPR meta-models. The minimum-cost objective function and stochastic constraints within the S-O model are based on the responses of ensemble meta-models. Reliability constraints are simultaneously integrated into the S-O model and are based on the ensemble surrogate responses to quantify the reliability of the design. Each meta-model in the ensemble model is trained and tested based on large datasets simulated by a numerical seepage modeling code (SEEP/W) [23] . Each meta-model imitates the numerical seepage modeling responses based on a particular field of heterogeneous hydraulic conductivity. The characteristics of each random field are based on means (μ) and standard deviations (σ) from a log-normal probability distribution function (PDF). Predictions of each surrogate model represent one of the seepage characteristics based on a particular random field involving different realizations of heterogeneous hydraulic conductivity. The process used to quantify the reliability of design within the RBOD framework is based on determining the number of stochastic responses and satisfying a particular constraint of the total number of meta-models (stochastic responses) in the ensemble. These meta-models were trained and tested based on various seepage datasets resulting from numerical simulation of various seepage models and scenarios of heterogeneous hydraulic conductivity.

Reliability level was formulated as an additional constraint. All stochastic constraints were continually controlled until the desired reliability level was achieved for each single iteration of the optimization model. Reliability constraints, stochastic constraints, and deterministic constraints were simultaneously evaluated with the objective function to attain the optimum solution. The majority of the constraints and the objective function were based on the ensemble surrogate model responses within the S-O model.

The optimization model is considered a complex one. Hence, the optimization solver and machine learning technique had to be efficient and accurate enough to provide reliable and accurate solutions. Therefore, GA was utilized as an optimization solver for this task. Based on the conceptual model shown in Figure 13, input data could be generated. The important design variables influencing seepage quantities were upstream cut-off (d1), downstream cut-off (d2), total width of WRS (b), upstream water head (H), and hydraulic conductivity.

The input data comprised 150 sets of seepage design variables (d1, d2, b, H), which were randomly generated utilizing the Halton sequences (HS) method

Figure 13. Conceptual model of a WRS.

[73] . Heterogeneous hydraulic conductivity was assumed to be a random field sampled from a log-normal distribution. Five standard deviations (0.85, 1.55, 2.25, 2.95 and 3.65) were assumed based on a constant mean of 2 m/day.

Based on each standard deviation value, a random field of hydraulic conductivity was generated and incorporated in the numerical seepage model. As unlimited realizations could be generated from a log-normal distribution for a certain standard deviation value, each input dataset (d1, d2, b, H) was simulated with four different random realizations (random fields) of the same standard deviation value. Some 600 simulated datasets were used for training each surrogate model for a particular seepage characteristic. This procedure ensures that the numerical responses for different hydraulic conductivity realizations are recorded and incorporated into the surrogate model’s training data. Figure 14 and Figure 15 represent different realizations of hydraulic conductivity for the same case and its effects on exit gradient (contours).

Figure 14. Different/random realizations of hydraulic conductivity for same standard deviation value.

Figure 15. Effect of different realizations (for the same σ value) of hydraulic conductivity variation on exit gradient contours.

Accordingly, the varied seepage quantities, such as uplift pressure on the upstream side (Pc1), downstream uplift pressure (Pe2) and exit gradient (ie), were determined by the numerical seepage modeling code four times for each input dataset (case). Furthermore, because exit gradient value is more critical than other quantities and hydraulic conductivity varies randomly, four points, shown in Figure 13, were selected at which exit gradient values were determined for each simulation. Determining four values of exit gradient and ensuring each value was within allowable limits ensured safety for WRS constructed on a heterogeneous flow domain. Hence, each training data set for a single meta-model included one set of input design variables (d1, d2, b, H) and four stochastically varied sets of output data (Pc1, Pe2, ie1, ie2, ie3, ie4). Therefore, the responses of meta-models reflect variation of seepage characteristics obtained from the four scenarios of random hydraulic conductivity. For each seepage design variable, five meta-models were trained to imitate different responses, reflecting the influence of five different hydraulic conductivity random fields drawn from the five log-normal distributions. As a result, thirty meta-models were built and grouped to develop six ensemble stochastic meta-models linked to the optimization model within the RBOD framework. Deterministic meta-models were developed separately to compare stochastic optimum solutions with deterministic solutions. Deterministic responses were used to train three meta-models (Pc1, Pe2, ie) based on expected hydraulic conductivity (σ = 0, μ = 2). Deterministic meta-models were incorporated in the deterministic S-O model to find the optimum solution of WRS for different head values.

4.1. Gaussian Process Regression (GPR) Surrogate Model

The meta-models were built based on Gaussian process regression (GPR) [74] [75] , a machine learning technique that imitates numerical model responses under various conditions. Many researchers dealing with geotechnical and civil engineering problems have demonstrated that GPR predicts certain responses more precisely than other machine learning techniques, such as SVM and backpropagation neural networks [76] - [84] . Building a surrogate model to use in the S-O approach is a delicate task. Although meta-models provide an expeditious alternative to numerical models, the training and testing phases need to be established carefully and accurately. The performance of meta-models must be precisely evaluated before being used in the S-O approach. The efficiency and accuracy of the developed meta-models increase the robustness of the linked S-O-based RBOD technique. The evaluation strategy is based on several statistical error measures (indices), such as the correlation coefficient (R), Nash-Sutcliffe efficiency (NSE) and mean square error (MSE), as shown in Table 4. These measures and others are briefly described by Gupta, Sorooshian, and Yapo [85] , Moriasi et al. [86] .

4.2. Formulating the Reliability Based Optimization Model

Formulation of the optimization model was based on the multi-realization

Table 4. Samples of surrogate model training testing error measures.

(“stacking”) optimization technique [18] [87] [88] . The reliability level is specified beforehand and the optimum design of the WRS satisfies the desired level of reliability at minimum cost. This can be achieved when the optimum solution satisfies a certain number of stochastic responses of all safety factors (constraints) of the total incorporated responses. This means that a particular reliability value (n/m) could be established within the S-O model by imposing a candidate design to satisfy n stochastic constraints of the total number (m) of constraints based on the WRS design’s safety factors. Each stochastic constraint is based on the responses of m meta-models within the stochastic ensemble surrogate model. For each safety factor, the reliability value n/m of the optimum design indicates that at least n stochastic constraints of all involved stochastic constraints (m) in the S-O model must be satisfied. Reliability is considered to be 100% when m/m of all constraints are satisfied and is considered to be 50% when 0.5 m/m of the stochastic constraints are satisfied, etc.

It is also important to note that some stochastic design variables, such as the thicknesses of upstream and downstream floors (t1, t2), involved in computation of the objective function are based on stochastic ensemble meta-models. Therefore, to provide a safe design, the maximum values of each thickness are considered in determining the objective function. The multiple-realization optimization-based RBOD using the stochastic S-O model is formulated as:

Find X = { x 1 , x 2 , x 3 , x 4 } = { d 1 , d 2 , b , b * }


f ( X ) = c f x 3 max ( t 1 m ) + max ( t 2 m ) 2 + t c s = 1 2 c s c x s (12)

The term { max ( t m ) } is a function returning the maximum of many stochastic values related to upstream and downstream floor thickness (t1, t2) of WRS. Similar constraints to those shown in Section 2 are utilized to formulate the optimization model described in this section. However, because we used stochastic responses due to uncertainty in seepage characteristics, there were 5 stochastic constraints for each safety factor. Based on a satisfaction percentage of these constraints, the desired reliability level can be achieved as shown below.

Z q l o g i c a l m = F s q m / F s q m a l l o w a b l e q , m (13)

g ( x ) q = m Z l o g i c a l m D R q (14)

where t 1 m , t 2 m represents the stochastic floor thicknesses at the upstream and downstream sides (Figure 13), respectively. These thicknesses were determined utilizing (m) stochastic meta-models related to uplift pressure. c f is the cost of constructing the floor per cubic meter ($400/m3); c c is the construction cost of the cut-offs per cubic meter, which is a function of cut-off depth, as shown in Equation (15), while t c is the cut-off thickness and equals 1.0 m.

c s c = x s 3 + 20 x s 2 + 200 x s + 400 s (15)

4.3. Results and Discussion

A RBOD framework based on stochastic S-O methodology was applied to hypothetical cases to study the effect of varying reliability on optimum WRS design. In these cases, the upstream head values (H) considered were 10 m, 20 m, 40 m, 60 m, 80 m and 100 m. The S-O models were implemented with different reliability levels (20%, 40%, 60%, 80% and 100%). The percentage of reliability only reflects the uncertainty of seepage quantities under WRS due to uncertainty associated with heterogeneous hydraulic conductivity.

The effect of reliability on optimum WRS design could be clearly seen by comparing the minimum construction costs obtained for different reliability levels, as shown in Figure 16. As logically expected, augmenting reliability significantly increases construction cost. For instance, the construction cost of a WRS which impounded 100 m water head and 100% reliability was around $143 million/m, whereas the cost was $102 million/m with 60% reliability. This infers that reliability substantially affects WRS design. Furthermore, ignoring uncertainty in hydraulic conductivity may result in unsafe designs. The deterministic optimum design, based on constant hydraulic conductivity (2 m/day), is presented

Figure 16. Optimum costs of WRSs with various reliability levels and head values.

in Figure 16. In general, the deterministic minimum cost curve was below the 60% reliability curve. This provides a general indication that the equivalent reliability of the deterministic design can be considered to be 50% - 60%, which is unsatisfactory for such an important structure.

The main role of d1 is to directly reduce the uplift pressure under the floor of the WRS and, indirectly, to reduce the exit gradient value. This is because the exit gradient is proportional to the uplift pressure located before the downstream cut-offs. In general, the optimum length of d1 decreases with reductions in head value, as shown in Table 5. In contrast, the optimum length of d1 is augmented by increasing the degree of reliability. However, for some values, especially with 100% reliability at H = 80, 40 m, the optimum length was less than for other reliability levels. This can be explained by considering that the objective function minimizes construction cost. Therefore, the optimization algorithm identifies the minimum construction cost solution for each case separately, as long as the design (decision) vector satisfies the constraints. On the other hand, because the meta-model responses are stochastic, it is difficult to expect the optimum value with different reliability ranks. Furthermore, if the optimization solver could provide an optimum solution that satisfies, for example, three out of five (60% reliability) of the stochastic constraints, this does not promise that the optimum solution with 80% reliability is in vicinity of the 60% solution. The justification is that additional stochastic constraints may require larger values of the variable, e.g. d1, which significantly increases the objective function value. Consequently, the optimization solver (GA) changes the direction of the search and continues in a more promising direction that provides a lower cost. Moreover, while the objective function is for minimum construction cost, the optimum solution for a certain reliability level does not necessarily follow the general trend of the other reliability levels. For instance, the optimum value of d1 at H = 80 m and 100% reliability was unexpectedly less than for other values. That may be logical if the values of d2, b and b* are considered simultaneously for this case. The value of d2, shown in Table 6, for the same case, was much larger than for other reliability levels because d2 is more important for reducing the crucial exit gradient value to the safe limit.

Table 5. Optimum length of upstream cut-off (d1) for various reliability (%) and head values.

Table 6. Optimum length of downstream cut-off (d2) for various reliability and head values.

On the other hand, the optimum value of d2, shown in Table 6, proportionally increases with increases in reliability. This design variable is the most important one as it controls the exit gradient value. In the reliability results, the majority of violated constraints were caused by the exit gradient safety factor. Therefore, the optimum solution for d2 with 100% reliability presents the highest value for different values of H to satisfy all stochastic seepage responses and relevant constraints.

Table 7 shows the optimum total WRS widths (b). The optimum width of b is the lowest for high reliability for different head values. Also, the optimization solver decreases the value of b and simultaneously increases the value of b*, which provided additional weight from upstream water (Table 7 and Table 8).

5. Optimum WRS Design Incorporating Uncertainty in Heterogeneous Hydraulic Conductivity Using MOMRO Models

To improve the efficiency and accuracy of the RBOD model and direct search optimization solver, a new approach was utilized. This approach is based on MOMRO model. The advantage of this approach is that some stochastic optimization constraints based on many ensemble meta-models are formulated as a second objective function to be minimized in the MOMRO model. The stochastic constraints used to achieve a safe exit gradient are formulated as a second stochastic objective function. The multi-objective optimization solver minimizes two stochastic objectives: exit gradient and construction cost. Desired reliability levels are implicitly incorporated into the objective functions and explicitly for the constraints. This significantly improves the solver’s search efficiency and aids in exploring more feasible candidate solutions in the search space. In this case, the multi-objective non-dominated sorting genetic algorithm-II (NSGA-II) was used.

The objective of this section is to find a safe, reliable and minimum-cost optimum WRS design that incorporates uncertainty in HHC estimates. The RBOD framework was implemented based on a more efficient and productive approach using the MOMRO technique. MOMRO integrates many stochastic responses

Table 7. Optimum WRS width (b) for various reliability and head values.

Table 8. Optimum width of b* for various reliability and head values.

from well-trained meta-models based on GPR machine learning techniques. These stochastic responses represent the uncertainties in estimates of particular seepage design variables, which are embedded in stochastic constraints and objective functions of MOMRO.

The steps used to generate data were the same as those mentioned in Section 4. The input design variables and seepage characteristics used were also the same. Some 150 random cases of input design variable data were generated. However, as random field HHC was used, the number of simulations for each input design variable was 20, including 20 different realizations of HHC for each case to cover a wide range of uncertainty in HHC. Each single realization represented a unique and randomly-varied distribution of hydraulic conductivity values from finite elements in the numerical model. Five log-normal distributions with standard deviations (σ) of 0.85 m/day, 1.55 m/day, 2.25 m/day, 2.95 m/day, and 3.65 m/day and constant a mean (μ = 2 m/day) were proposed to generate different HHCs. For each seepage quantity, a stochastic ensemble surrogate model was developed containing 20 meta-models. Therefore, for a single input dataset, 20 stochastic responses were obtained by the ensemble surrogate model for processing in the MOMRO model based on the RBOD technique. Similar to Section 4, GPR machine learning was used to develop 120 meta-models and build six stochastic meta-models. For each surrogate model, there were 150 sets of training and testing data (cases).

The training/testing performance of meta-models must be accurately evaluated before using them in the RBOD model. The developed GPR meta-models were evaluated using several error measures as described in Section 4. The majority of meta-models demonstrated perfect training and testing performance (Table 9).

5.1. Formulation of the Reliability Based MOMRO Model

Formulation of a multi-realization optimization model based on a single objective function with numerous stochastic constraints may lead to sub-optimum or infeasible solutions. The RBOD approach required imposition of a large number of explicit constraints, which needed to be satisfied as binding conditions for a feasible solution. Many attempts were made to formulate an RBOD model for this study with a large number (120) of stochastic meta-model-based stochastic constraints using a single objective function, but the majority of solutions were infeasible.

As multi-realization technique-based reliability requires a large number of stochastic constraints, an optimal solution search process based on evolutionary algorithms may produce an infeasible solution. Search efficiency decreases with increasing numbers of constraints and problem complexity [89] [90] . Therefore, a new formulation of the RBOD model was adopted in this study to improve the search process for such complex optimization tasks. The most important stochastic constraints are exit gradient constraints, as they are significantly influenced by HHC uncertainty and have critical impacts on WRS design and safety. These constraints were transformed as a second objective function to be minimized in addition to the first objective function related to WRS construction cost. Hence, a multi-objective optimization formulation was implemented to significantly decrease the number of constraints and improve searching efficiency. Reliability was included for exit gradient (objective function) and also implemented for stochastic constraints using a multi-realization technique.

The multi-realization optimization technique was based on formulating stochastic constraints based on the developed ensemble stochastic meta-models. For each safety factor or condition in the optimization model, there was one or more ensemble stochastic surrogate model/s encompassing 20 surrogate model responses for a specified seepage design variable. The desired reliability level was attained by allowing the optimum solution to satisfy any fraction (n) of the total number (m = 20) of constraints for all stochastic constraints, where n/m is equivalent to the required reliability level.

Table 9. Samples of surrogate model training testing error measure.

The second objective function, which minimizes the exit gradient value, integrates reliability by incorporating exit gradient stochastic responses in determining the objective function. As exit gradient is minimized, 20 stochastic exit gradient responses based on ensemble stochastic meta-models are determined and sorted in ascending order. The maximum value of all obtained exit gradients is then selected to be minimized. This is equivalent to 99.9% reliability because the resulting exit gradient is the safest of all other stochastic values. To attain 80% reliability, for example, the optimization solver is formulated to minimize the fifth maximum value (based on 20 responses) and allow up to four stochastic responses of exit gradient to be higher than the one selected for the objective function.

As there are four locations at which to determine the exit gradient values (ie1, ie2, ie3, ie4) the maximum value for each location is determined and their average is considered as the second objective function. The same technique is applied to the first objective function of minimizing construction cost. Construction cost is based on the upstream and downstream cut-off depths (d1, d2), and thicknesses (t1, t2), and floor width (b). Floor thicknesses are based on the stochastic responses of uplift pressure ensemble meta-models (pc1, pe2).

The constraints for each safety factor, based on 20 developed meta-models, are formulated as shown in Section 1. However, the exit gradient safety factors (constraints) are formulated as a second objective function. Furthermore, the reliability constraints and all other variables and parameters are exactly the same as in Section 4. The modification in the formulation of the optimization model for MOMRO is as shown below

Find X = { x 1 , x 2 , x 3 , x 4 } = { d 1 , d 2 , b , b * }

Minimize, f 1 ( X ) = c f b max ( m ω ) ( t 1 m ) + max ( m ω ) ( t 2 m ) 2 + t c s = 1 2 c s c d s (16)


f 2 ( X ) = max ( m ω ) ( i e 1 m ) + max ( m ω ) ( i e 2 m ) + max ( m ω ) ( i e 3 m ) + max ( m ω ) ( i e 4 m ) 4 (17)

i e i m = ε i m ( H , d 1 , d 2 , b , k m ) i , m (18)

where max ( m ω ) is a function sorting stochastic responses in ascending order and returns the ( m ω ) th value of the sorted vector. m is the number of stochastic responses (20), and ω is based on the desired reliability level.

5.2. Results and Discussion

The MOMRO technique was applied to hypothetical design scenarios/cases to evaluate the performance of the RBOD-based MOMRO technique. These cases included five different upstream head values (100 m, 80 m, 60 m, 40 m, 20 m) each subject to four reliability levels. The obtained Pareto-optimum fronts for each head value, including different scenarios of reliability level, are presented in Figure 17 to Figure 18. To make an appropriate decision, minimum allowable

Figure 17. Optimum Pareto front for different reliability levels (H = 100 m).

Figure 18. Optimum Pareto front for different reliability levels (H = 80 m).

deterministic safe exit gradient [4] [91] values were used to locate safe and feasible optimum solutions, as shown in Figure 17 and Figure 18. There are two vertical lines showing the locations of safe exit gradient factors 5 and 3, considering a critical gradient value of 1.15. Based on these values, the minimum safe exit gradient can be allocated for different reliability levels. To provide greater safety related to exit gradient, many possible Pareto optimal solutions were available for consideration with ascending construction costs, and a WRS designer could select one according to their preference.

The effects of reliability on optimum WRS design were significant. Increasing reliability increased construction cost. For instance, the minimum construction costs for H = 100 m, reliability levels of 40%, 60%, 80% and 100%, and exit gradient safety factor = 5 were $112,191,378, $129,171,757, $162,166,799 and $268,206,048, respectively. Similarly, for the same reliability levels and an exit gradient = 3, the construction costs were $59,951,442, $79,158,696, $106,049,766 and $160,838,745. This means that considering reliability in WRS design significantly affects the optimum design attributes. Moreover, for high reliability levels, only a few applicable (feasible) scenarios can be obtained from the Pareto optimal front. For example, for H = 100 m, 99.9% reliability, and exit gradient safety factor = 5, only a few points were found at higher construction cost ($268,206,048.88).

The deterministic optimum Pareto front related to the expected hydraulic conductivity (2 m/day) was also considered in this study. In general, the deterministic Pareto optimal was located close to 60% reliability trade-offs. However, some deterministic optimum solutions approached 40% reliability solutions. The 60% or 40% reliability of the deterministic solutions mean that there is a high probability of the exit gradient approaching the critical value, which might lead to piping failure. Based on this, we can deduce that the deterministic safety factors of 3 and 5 are insufficient to provide adequate safety for such important projects (WRSs), and are inappropriate to measure the safety of seepage designs that incorporate a degree of uncertainty. This is true if we assume that the prescribed safety factor is used to quantify uncertainty in the HHC only.

One important benefit of using multi-objective optimization in RBOM is the diversity of provided optimum solutions. The multi-objective optimization solver provides many optimum solutions for the same objective function values (approximately). These solutions could not be obtained by a single objective optimization model. These solutions provide more flexible options because some optimum solutions are more applicable in terms of design requirements, such as field limitation and construction procedures, etc. Table 10 presents a few arbitrarily selected example solutions with the same objective function values including different optimum solution (X) scenarios.

6. Summary and Conclusions

This study presents several applications for a coupled simulation-optimization methodology that can optimize the design of WRSs. The method uses multiple meta-models of the seepage characteristics associated with such structures. These meta-models were trained on multiple datasets of simulated seepage scenarios to provide accurate seepage solutions. The methodology was extended to

Table 10. Optimum solution values for the same objective function values obtained by NSGA-II.

optimize the design of WRSs based on complex seepage models that consider non-homogenous and anisotropic hydraulic conductivity. Furthermore, the effect of uncertainty in seepage characteristics due to uncertain hydraulic conductivity was considered to determine reliable, minimum-cost designs.

The various meta-models are based on machine learning techniques such as ANN, SVM, and GRP. Many error measures were used to evaluate the predictive accuracy of the models to ascertain that their predictions were within allowable error ranges. The parameters of the optimization solver and machine learning-based meta-models were carefully selected according to multiple systematic experiments and trial-and-error iterations.

The results demonstrate that the S-O technique is a useful and efficient method of optimizing WRSs design with simple and complex seepage flow domains by integrating seepage modelling. The most important advantage of the S-O methodology is the incorporation the effect of uncertainty in seepage due to uncertainty in hydraulic conductivity on optimum design. This includes the incorporation of many levels of uncertainty or reliability based optimum design. Furthermore, several formulations based on multi-objective multi-realization optimization models were implemented to provide less restrictive formulation-based optimization models. This approach presents diverse optimum solutions, providing several design alternatives for the same objective functions (e.g. cost and exit gradient).

It is recommended for futures studies apply this S-O methodology to WRS design while incorporating different sources of uncertainty such as those of the upstream water table, critical exit gradient, and meta-model predictions. Furthermore, incorporating the requirement of structural design in WRS as additional constraints is a more advanced step to optimizing their designs.

Conflicts of Interest

The authors declare no conflicts of interest regarding the publication of this paper.

Cite this paper

Al-Juboori, M. and Datta, B. (2018) An Overview of Recently Developed Coupled Simulation Optimization Approaches for Reliability Based Minimum Cost Design of Water Retaining Structures. Open Journal of Optimization, 7, 79-112.


  1. 1. Griffiths, D. and Fenton, G.A. (1993) Seepage Beneath Water Retaining Structures Founded on Spatially Random Soil. Géotechnique, 43, 577-587.

  2. 2. Harr, M.E. (2012) Groundwater and Seepage. McGraw Hill, New York.

  3. 3. Garg, S.K. (1987) Irrigation Engineering and Hydraulic Structures. Khanna Publishers, Delhi.

  4. 4. Khosla, A.N., Bose, N.K. and Taylor, E.M. (1936) Design of Weirs on Permeable Foundations. Central Board of Irrigation, New Delhi.

  5. 5. Lambe, T.W. and Whitman, R.V. (1969) Soil Mechanics. John Wiley & Sons, New York.

  6. 6. Das, A. and Datta, B. (1999) Development of Management Models for Sustainable Use of Coastal Aquifers. Journal of Irrigation and Drainage Engineering, 125, 112-121.

  7. 7. Wagner, B.J. and Gorelick, S.M. (1986) A Statistical Methodology for Estimating Transport Parameters: Theory and Applications to One-Dimensional Advectivec-Dispersive Systems. Water Resources Research, 22, 1303-1315.

  8. 8. Gorelick, S.M. (1983) A Review of Distributed Parameter Groundwater Management Modeling Methods. Water Resources Research, 19, 305-319.

  9. 9. Willis, R. and Finney, B.A. (1988) Planning Model for Optimal Control of Saltwater Intrusion. Journal of Water Resources Planning and Management, 114, 163-178.

  10. 10. Ayvaz, M.T. (2016) A Hybrid Simulation-Optimization Approach for Solving the Areal Groundwater Pollution Source Identification Problems. Journal of Hydrology, 538, 161-176.

  11. 11. Bhattacharjya, R.K. and Datta, B. (2009) ANN-GA-Based Model for Multiple Objective Management of Coastal Aquifers. Journal of Water Resources Planning and Management, 135, 314-322.

  12. 12. Bhattacharjya, R.K., Datta, B. and Satish, M.G. (2007) Artificial Neural Networks Approximation of Density Dependent Saltwater Intrusion Process in Coastal Aquifers. Journal of Hydrologic Engineering, 12, 273-282.

  13. 13. Datta, B., Chakrabarty, D. and Dhar, A. (2011) Identification of Unknown Groundwater Pollution Sources Using Classical Optimization with Linked Simulation. Journal of Hydro-Environment Research, 5, 25-36.

  14. 14. Dhar, A. and Datta, B. (2009) Saltwater Intrusion Management of Coastal Aquifers. I: Linked Simulation-Optimization. Journal of Hydrologic Engineering, 14, 1263-1272.

  15. 15. Heydari, F., Saghafian, B. and Delavar, M. (2016) Coupled Quantity-Quality Simulation-Optimization Model for Conjunctive Surface-Groundwater Use. Water Resources Management, 30, 4381-4397.

  16. 16. Sreekanth, J. and Datta, B. (2015) Review: Simulation-Optimization Models for the Management and Monitoring of Coastal Aquifers. Hydrogeology Journal, 23, 1155-1166.

  17. 17. Sreekanth, J. and Datta, B. (2010) Multi-Objective Management of Saltwater Intrusion in Coastal Aquifers Using Genetic Programming and Modular Neural Network Based Surrogate Models. Journal of Hydrology, 393, 245-256.

  18. 18. Sreekanth, J. and Datta, B. (2011) Coupled Simulation-Optimization Model for Coastal Aquifer Management Using Genetic Programming-Based Ensemble Surrogate Models and Multiple-Realization Optimization. Water Resources Research, 47, 1-17.

  19. 19. Singh, R.M. (2010) Design of Barrages with Genetic Algorithm Based Embedded Simulation Optimization Approach. Water Resources Management, 25, 409-429.

  20. 20. Singh, R.M. (2011) Genetic Algorithm Based Optimal Design of Hydraulic Structures with Uncertainty Characterization. International Conference on Swarm, Evolutionary, and Memetic Computing, 19-21 December 2011, 742-749.

  21. 21. Al-Suhaili, R.H. and Karim, R.A. (2014) Optimal Dimensions of Small Hydraulic Structure Cutoffs Using Coupled Genetic Algorithm and ANN Model. Journal of Engineering, 20, 1-19.

  22. 22. Al-Juboori, M. (2018) Seepage Criteria Based Optimal Design of Water Retaining Structures with Reliability Quantification Utilizing Surrogate Model Linked Simulation Optimization Approach. James Cook University, Australia.

  23. 23. Krahn, J. (2012) Seepage Modeling with SEEP/W: An Engineering Methodology. GEO-SLOPE International Ltd., Calgary.

  24. 24. Lin, C.D. and Tang, B. (2015) Latin Hypercubes and Space-Filling Designs. In: Dean, A., Morris, M., Stufken, J. and Bingham, D., Eds., Handbook of Design and Analysis of Experiments, CRC Press, Lowa.

  25. 25. Jain, A. and Kumar, A. (2006) An Evaluation of Artificial Neural Network Technique for the Determination of Infiltration Model Parameters. Applied Soft Computing, 6, 272-282.

  26. 26. Garcia, L.A. and Shigidi, A. (2006) Using Neural Networks for Parameter Estimation in Ground Water. Journal of Hydrology, 318, 215-231.

  27. 27. Ersayin, D. (2006) Studying Seepage in a Body of Earth-Fill Dam by (Artificial Neural Networks) ANNs. Institute of Technology, Izmir University, Turkey.

  28. 28. Szidarovszky, F., Coppola Jr., E.A., Long, J., Hall, A.D. and Poulton, M.M. (2007) A Hybrid Artificial Neural Network-Numerical Model for Ground Water Problems. Groundwater, 45, 590-600.

  29. 29. Kim, Y.-S. and Kim, B.-T. (2008) Prediction of Relative Crest Settlement of Concrete-Faced Rockfill Dams Analyzed Using an Artificial Neural Network Model. Computers and Geotechnics, 35, 313-322.

  30. 30. Joorabchi, A., Zhang, H. and Blumenstein, M. (2009) Application of Artificial Neural Networks to Groundwater Dynamics in Coastal Aquifers. Journal of Coastal Research, SI56, 966-970.

  31. 31. Nourani, V., Sharghi, E. and Aminfar, M.H. (2012) Integrated ANN Model for Earthfill Dams Seepage Analysis: Sattarkhan Dam in Iran. Artificial Intelligence Research, 1, 22-37.

  32. 32. Santillán, D., Fraile-Ardanuy, J. and Toledo, M. (2013) Dam Seepage Analysis Based on Artificial Neural Networks: The Hysteresis Phenomenon. The 2013 International Joint Conference on Neural Networks (IJCNN), Dallas, 4-9 August 2013, 1-8.

  33. 33. Cavazzuti, M. (2012) Optimization Methods: From Theory to Design Scientific and Technological Aspects in Mechanics. Springer Science & Business Media, London.

  34. 34. MathWorks (2018) Matlab Statistics and Machine Learning Toolbox. Natick.

  35. 35. Bligh, W.G. (1915) Dams and Weirs: An Analytical and Practical Treatise on Gravity Dams and Weirs; Arch and Buttress Dams; Submerged Weirs; and Barrages. American Technical Society, Chicago.

  36. 36. US Army Corps of Engineers (1987) Engineering and Design Flotation Stability Criteria for Concrete Hydraulic Structures.

  37. 37. Al-Juboori, M. and Datta, B. (2018) Performance Evaluation of a Genetic Algorithm-Based Linked Simulation-Optimization Model for Optimal Hydraulic Seepage-Related Design of Concrete Gravity Dams. Journal of Applied Water Engineering and Research, 1-25.

  38. 38. Al-Juboori, M. and Datta, B. (2018) Linked Simulation-Optimization Model for Optimum Hydraulic Design of Water Retaining Structures Constructed on Permeable Soils. International Journal of GEOMATE, 14, 39-46.

  39. 39. Beckwith, C.W., Baird, A.J. and Heathwaite, A.L. (2003) Anisotropy and Depth-Related Heterogeneity of Hydraulic Conductivity in a Bog Peat. I: Laboratory Measurements. Hydrological Processes, 17, 89-101.

  40. 40. Burger, R.L. and Belitz, K. (1997) Measurement of Anisotropic Hydraulic Conductivity in Unconsolidated Sands: A Case Study from a Shoreface Deposit, Oyster, Virginia. Water Resources Research, 33, 1515-1522.

  41. 41. Greenkorn, R., Johnson, C. and Shallenberger, L. (1964) Directional Permeability of Heterogeneous Anisotropic Porous Media. Society of Petroleum Engineers Journal, 4, 124-132.

  42. 42. Terzaghi, K., Peck, R.B. and Mesri, G. (1996) Soil Mechanics in Engineering Practice. John Wiley & Sons, New York.

  43. 43. Liu, J., Chang, J.-X. and Zhang, W.-G. (Year) Groundwater Level Dynamic Prediction Based on Chaos Optimization and Support Vector Machine. 3rd International Conference on Genetic and Evolutionary Computing, Guilin, 14-17 October 2009.

  44. 44. Bashi-Azghadi, S.N., Kerachian, R., Rez, M. and Soloukia, B. (2010) Characterizing an Unknown Pollution Source in Groundwater Resources Systems Using PSVM and PNN. Expert Systems with Applications, 37, 7154-7161.

  45. 45. Fisher, W.D., Camp, T.K. and Krzhizhanovskaya, V.V. (2016) Anomaly Detection in Earth Dam and Levee Passive Seismic Data Using Support Vector Machines and Automatic Feature Selection. Journal of Computational Science, 20, 143-153.

  46. 46. Mahani, A.S., Shojaee, S., Salajegheh, E. and Mohsen, K. (2015) Hybridizing Two-Stage Meta-Heuristic Optimization Model with Weighted Least Squares Support Vector Machine for Optimal Shape of Double-Arch Dams. Applied Soft Computing, 27, 205-218.

  47. 47. Parsaie, A., Yonesi, H.A. and Najafian, S. (2015) Predictive Modeling of Discharge in Compound Open Channel by Support Vector Machine Technique. Modeling Earth Systems and Environment, 1, 1-6.

  48. 48. Ranković, V., Grujović, N. and Milivojević, N. (2014) Development of Support Vector Regression Identification Model for Prediction of Dam Structural Behaviour. Structural Safety, 48, 33-39.

  49. 49. Su, H., Chen, Z. and Wen, Z. (2016) Performance Improvement Method of Support Vector Machine-Based Model Monitoring Dam Safety. Structural Control and Health Monitoring, 23, 252-266.

  50. 50. Azamathulla, H.M., Ghani, A., Chang, C., Hasan, Z. and Zakaria, N. (2010) Machine Learning Approach to Predict Sediment Load—A Case Study. CLEAN—Soil, Air, Water, 38, 969-976.

  51. 51. Bhagwat, P.P. and Maity, R. (2012) Multistep-Ahead River Flow Prediction Using LS-SVR at Daily Scale. Journal of Water Resource and Protection, 4, 528-539.

  52. 52. Cimen, M. (2008) Estimation of Daily Suspended Sediments Using Support Vector Machines. Hydrological Sciences Journal, 53, 656-666.

  53. 53. Goel, A. and Pal, M. (2009) Application of Support Vector Machines in Scour Prediction on Grade-Control Structures. Engineering Applications of Artificial Intelligence, 22, 216-223.

  54. 54. Han, D., Chan, L. and Zhu, N. (2007) Flood Forecasting Using Support Vector Machines. Journal of Hydroinformatics, 9, 267-276.

  55. 55. Lesaja, G. (2009) Introducing Interior-Point Methods for Introductory Operations Research Courses and/or Linear Programming Courses. The Open Operational Research Journal, 3, 1-12.

  56. 56. Mulligan, A.E. and Ahlfeld, D.P. (2002) A New Interior-Point Boundary Projection Method for Solving Nonlinear Groundwater Pollution Control Problems. Operations Research, 50, 636-644.

  57. 57. Liu, M., Tso, S.K. and Cheng, Y. (2002) An Extended Nonlinear Primal-Dual Interior-Point Algorithm for Reactive-Power Optimization of Large-Scale Power Systems with Discrete Control Variables. IEEE Transactions on Power Systems, 17, 982-991.

  58. 58. Baecher, G.B. and Christian, J.T. (2005) Reliability and Statistics in Geotechnical Engineering. John Wiley & Sons, West Sussex.

  59. 59. Baroni, G., Zink, M., Kumar, R. and Attinger, S. (2017) Effects of Uncertainty in Soil Properties on Simulated Hydrological States and Fluxes at Different Spatio-Temporal Scales. Hydrology and Earth System Sciences, 21, 2301-2320.

  60. 60. Christian, J.T., Ladd, C.C. and Baecher, G.B. (1994) Reliability Applied to Slope Stability Analysis. Journal of Geotechnical Engineering, 120, 2180-2207.

  61. 61. Deng, Z.-P., Li, D., Cao, Z. and Phoon, K (2017) Reliability Evaluation of Slope Considering Geological Uncertainty and Inherent Variability of Soil Parameters. Computers and Geotechnics, 92, 121-131.

  62. 62. Duncan, J.M. (2000) Factors of Safety and Reliability in Geotechnical Engineering. Journal of Geotechnical and Geoenvironmental Engineering, 126, 307-316.

  63. 63. Hicks, M.A., Nuttall, J.D. and Chen, J. (2014) Influence of Heterogeneity on 3D Slope Reliability and Failure Consequence. Computers and Geotechnics, 61, 198-208.

  64. 64. Hicks, M.A. and Spencer, W.A. (2010) Influence of Heterogeneity on the Reliability and Failure of a Long 3D Slope. Computers and Geotechnics, 37, 948-955.

  65. 65. Popescu, R., Deodatis, G. and Nobahar, A. (2005) Effects of Random Heterogeneity of Soil Properties on Bearing Capacity. Probabilistic Engineering Mechanics, 20, 324-341.

  66. 66. Banerjee, S. and Datta, B. (1991) Reliability Analysis of Thaw-Induced Pore Pressures. Journal of Cold Regions Engineering, 5, 125-141.

  67. 67. Ahmed, A.A. (2012) Stochastic Analysis of Seepage under Hydraulic Structures Resting on Anisotropic Heterogeneous Soils. Journal of Geotechnical and Geoenvironmental Engineering, 139, 1001-1004.

  68. 68. Griffiths, D. and Fenton, G.A. (2004) Probabilistic Slope Stability Analysis by Finite Elements. Journal of Geotechnical and Geoenvironmental Engineering, 130, 507-518.

  69. 69. Griffiths, D. and Fenton, G.A. (1997) Three-Dimensional Seepage through Spatially Random Soil. Journal of Geotechnical and Geoenvironmental Engineering, 123, 153-160.

  70. 70. Le, T.M.H., Gallipoli, D., Sanchez, M. and Wheeler, S. (2012) Stochastic Analysis of Unsaturated Seepage through Randomly Heterogeneous Earth Embankments. International Journal for Numerical and Analytical Methods in Geomechanics, 36, 1056-1076.

  71. 71. Zhu, X., Wang, X., Li, X., Liu, M. and Cheng, Z. (2017) A New Dam Reliability Analysis Considering Fluid Structure Interaction. Rock Mechanics and Rock Engineering, 51, 2505-2516.

  72. 72. Jiang, S.-H., Li, D., Zhang, L. and Zhou, C. (2014) Slope Reliability Analysis Considering Spatially Variable Shear Strength Parameters Using a Non-Intrusive Stochastic Finite Element Method. Engineering Geology, 168, 120-128.

  73. 73. Loyola, D., Pedergnana, M. and García, S.G. (2016) Smart Sampling and Incremental Function Learning for Very Large High Dimensional Data. Neural Networks, 78, 75-87.

  74. 74. Rasmussen, C.E. (2004) Gaussian Processes in Machine Learning: An Advanced Lectures on Machine Learning. Springer, Berlin, Heidelberg, 63-71.

  75. 75. Shi, J.Q. and Choi, T. (2011) Gaussian Process Regression Analysis for Functional Data. CRC Press, Chapman and Hall, New York.

  76. 76. Chen, T. and Ren, J. (2009) Bagging for Gaussian Process Regression. Neurocomputing, 72, 1605-1610.

  77. 77. He, P., Xiao, J., Zhang, Q., Xu, F. and Zhang, J. (2017) Shallow Sliding Failure Prediction Model of Expansive Soil Slope Based on Gaussian Process Theory and Its Engineering Application. KSCE Journal of Civil Engineering, 22, 1-11.

  78. 78. Kang, F., Han, S., Salgado, R. and Li, J. (2015) System Probabilistic Stability Analysis of Soil Slopes Using Gaussian Process Regression with Latin Hypercube Sampling. Computers and Geotechnics, 63, 13-25.

  79. 79. Kang, F., et al. (2017) Slope Stability Evaluation Using Gaussian Processes with Various Covariance Functions. Applied Soft Computing, 60, 387-396.

  80. 80. Kim, K., Lee, D. and Essa, I. (Year) Gaussian Process Regression Flow for Analysis of Motion Trajectories. 2011 International Conference on Computer Vision, Barcelona, 6-13 November 2011.

  81. 81. Li, S.-C., et al. (2017) Gaussian Process Model of Water Inflow Prediction in Tunnel Construction and Its Engineering Applications. Tunnelling and Underground Space Technology, 69, 155-161.

  82. 82. Pal, M. and Deswal, S. (2010) Modelling Pile Capacity Using Gaussian Process Regression. Computers and Geotechnics, 37, 942-947.

  83. 83. Samui, P. and Jagan, J. (2013) Determination of Effective Stress Parameter of Unsaturated Soils: A Gaussian Process Regression Approach. Frontiers of Structural and Civil Engineering, 7, 133-136.

  84. 84. Xu, J.W. and Suzuki, K. (2011) Massive-Training Support Vector Regression and Gaussian Process for False-Positive Reduction in Computer-Aided Detection of Polyps in CT Colonography. Medical Physics, 38, 1888-1902.

  85. 85. Gupta, H.V., Sorooshian, S. and Yapo, P.O. (1999) Status of Automatic Calibration for Hydrologic Models: Comparison with Multilevel Expert Calibration. Journal of Hydrologic Engineering, 4, 135-143.

  86. 86. Moriasi, D.N., Allah, O. and Najafian, Y. (2007) Model Evaluation Guidelines for Systematic Quantification of Accuracy in Watershed Simulations. Transactions of the ASABE, 50, 885-900.

  87. 87. Chan, N. (1993) Robustness of the Multiple Realization Method for Stochastic Hydraulic Aquifer Management. Water Resources Research, 29, 3159-3167.

  88. 88. Feyen, L. and Gorelick, S.M. (2005) Framework to Evaluate the Worth of Hydraulic Conductivity Data for Optimal Groundwater Resources Management in Ecologically Sensitive areas. Water Resources Research, 41, W03019.

  89. 89. Dorsey, R. and Mayer, W. (1995) Genetic Algorithms for Estimation Problems with Multiple Optima, Nondifferentiability, and Other Irregular Features. Journal of Business & Economic Statistics, 13, 53-66

  90. 90. Kolda, T.G., Lewis, R.M. and Torczon, V. (2003) Optimization by Direct Search: New Perspectives on Some Classical and Modern Methods. SIAM Review, 45, 385-482.

  91. 91. Tanchev, L. (2014) Dams and Appurtenant Hydraulic Structures. CRC Press, London.