In geotechnical engineering the optimization techniques can be categorized as: 1) mathematical programming; 2) optimality criteria methods; and 3) heuristic search algorithms. In mathematical programming, the major characteristic includes the linear and nonlinear programming. In the linear programming, the objective function and the associated constraints are represented in a linear combination of design variables. The linearization of objective function or constraints is not always easy and if the linearization techniques are used, the error in the representation of linearized problem is inevitable    . On the other hand, the nonlinear programming is introduced which was developed for unconstrained nonlinear problems. The incidence of the optimal solution is shown by the proof of Kuhan-Tucker conditions (KT)  which are the necessary conditions for the justification of the optimal solution. However, the application of KT conditions is enormously difficult for the most of engineering problems  .
The optimality criteria methods (OCM) are developed based on the combination of Kuhan-Tucker conditions from nonlinear mathematical programming and Lagrangian multipliers. In this approach, KT conditions support the necessary desires for the optimal solution while Lagrangian multiplier guarantees the satisfaction of constrains in the optimization problem. OCMs are used vastly in engineering problems including the continuous and discontinuous design variables   . OCMs basically describe continuous design variables. However, by considering two major steps, this procedure is compatible for discrete variables, as well. First, the optimization problem is defined based on the assumption of continuous solution. Second, a set of discrete values are estimated based on the solution from first step  . However, two problems are associated with this approach; for instance, if OCMs are used in structural design while a single cross sectional property is chosen as a design variable and the other properties are expressed as a function of the design variable. First, the relationship between the design variable and cross sectional properties is not unique. Second, the selected discrete variables may lead to a different structural response which does not fit the prescribed constraints  .
The techniques of the third group in the optimization problems are heuristic search algorithms which are relatively new optimization methods in geotechnical engineering. Evolutionary methods do not require an explicit relationship between the objective functions and constraints. Moreover, these methods are independent of objective function or constraint gradient. Instead, the set of random values for objective function is adjusted thorough an evolutionary procedure and subsequently, any violation of constraints is reflected in each iteration. Regarding the gradient based methods, considerable part of computational cost is served in sensitivity analysis phase whereas applications of evolutionary techniques which are based on probabilistic searching algorithm save computational time  ; this is because of evolutionary methods which do not require the gradient information and avoid sensitivity analysis. This property enables engineers to implement these techniques in the problems including discontinuity in the design variables space. In general, evolutionary methods work concurrently with a population of design points in the space of design variables instead of a single design point. These methods also can be easily implemented in discrete, continuous, and mixed optimization problems with minor adaptations. In addition, the open format for constraint statements and the possibility of defining multiple scenarios in the optimization process fascinated many researchers to implement heuristic search algorithms in engineering design optimization  .
The safety in designs is highly respected in geotechnical engineering. On the other hand, contractors and private sectors prefer economical design whit minimum labor effort and construction cost. Subsequently, implementation of robust optimization techniques which establish a trade-off between safety and total cost of geotechnical projects is necessary in practice.
In this study, we reviewed the two most implemented evolutionary methods in geotechnical engineering, Genetic Algorithm (GA) and Particle Swarm Optimization (PSO). The applications of GA and PSO in three common types of geotechnical problems including unconfined seepage analysis, slope stability analysis, and foundation designs are studied and subsequently, for each problem, the results from literature are regenerated, compared and discussed.
2. Genetic Algorithm (GA)
Genetic Algorithm which is based on the process of natural evolution is one the best-known heuristic search algorithms. GA has received substantial attention in engineering design optimization in the recent decades. The first time that GA was introduced in computer science goes back to sixties when a team of biologist tried to implement the process of evolution in nature in a computer code  . GA refers to any population-based algorithm which uses selection, crossover, and mutation over chromosomes to find the optimal solution. In general, a member of the population is called chromosome/genotype which is binary or real valued string. Following Barricelli  several types of GA have been introduced in optimization studies. However, one can simply define a given optimization problem in GA considering three main steps including  :
Step 0: Initialization.
Step 1: Selection.
Step 2: Generation.
Step 3: Stopping criteria.
Generation of initial chromosomes in step 0 (genotypes) is the first step in GA algorithm. In the most practical problems genotypes are generated randomly and goodness of each chromosome is evaluated via objective function and associated constraints.
Step 1: the selection operator is introduced and applied to the current population to create an intermediate one. The initial and intermediate populations are same in the first step (step 0). However, in the subsequent iterations, the imposition of selection operator forms the intermediate population.
Step 2: crossover-mutation operators are implemented to the results from Step 1 to create next population. The term crossover is used when the generator (operator) forms a chromosome by combining the properties of each of two parental chromosomes. However, the term mutation is applicable when a new chromosome is formed by introducing small alterations to the properties of single parent chromosome  .
There are many choices to design/encoding an optimization program based on the mentioned steps and the decision on design/encoding of GA depends on the nature of the problem  . The algorithm design depends on experience, the model specifications and associating the results from experiments with different heuristic search algorithms. However, a common design for GA with standard genetic operators can be defined as  :
1) Randomly generate an initial population including M chromosomes.
2) Calculate the fitness, , of all chromosomes in M.
3) Create the evolved population based on:
a) Using proportional fitness selection for selected two chromosomes, m1 and m2.
b) Apply crossover function to m1 and m2 to produce a new chromosome (m3).
c) Apply mutation function to m3 to produce.
d) Add to the next population.
4) Replace the new population with the old one.
5) If the stopping criteria have not been met, repeat the procedure from step 2.
Interested readers are encouraged to find more resources that are freely available online. For more details in GA, the Genetic Algorithm Archive by the US Navy Center for Applied Research in Artificial Intelligence (https://www.nrl.navy.mil/itd/aic/) or the available Toolboxes which provide instructive example for standard GA (https://sourceforge.net/projects/gatoolbox/) are introduced.
3. Particle Swarm Optimization (PSO)
Particle Swarm Optimization is a metaheuristic search algorithm introduced by Eberhart and Kennedy  to find optimal solution in engineering design optimization. PSO is based on the concept social models and swarm theories. The swarm consists of individual particles which mutually try to find the solution in the search space. In each iteration, individuals (particles) move toward the best solution which is experienced by them (Personal best) and concurrently to the best solution which is obtained by the other particles (Global best). PSO is established based on few or even no assumption on the search space; this feature enables PSO to search the optimum solution in a wide search space. In addition, PSO can be used in the optimization problems which are irregular, noisy, or dynamic  .
Regarding the PSO algorithm, the swarm of particles walks through the n-dimensional search space of a problem including n degree of freedoms. Each particle (i) in the mth iteration is identified by its position (), velocity (), and memory () vector. indicates the potential solution in the mth iteration, changes the position to investigate the solution in the search space, and stores the best personal solution for the ith particle. is updated when a better solution is found by the ith particle in the search space. Subsequently, each particle shares the best solution with the other particles in its vicinity and the Global best is updated (). Finally, the Personal and Global best solution is utilized simultaneously to update and. The updated velocity and position vector in iteration are defined as:
where and are constant values and defined as the acceleration factors, and are randomly generated weights with values between 0 and 1, is an inertia weight that controls the impact of velocity from the previous iteration on the newly computed velocity, and and are the Personal and Global best solutions in the mth iteration  . In general, these parameters are determined based on tuning and trial-error observations.
With respect to Equation  and aforementioned concepts, PSO algorithm covers the following steps to find the optimum solution:
1) Random generation the initial position and velocity vector including and, respectively.
2) Evaluate the objective function with respect to predefined constraints.
3) Comparing the fitness value of each particle to its local best. If the current value is better than the previous one, then it is replaced with the new value.
4) Comparing a particle’s fitness value to the population’s best. If the current value is better than the previous one, it is replaced as the global best.
5) Changing the velocity and position vectors using Equations (1a) and (1b).
6) Steps 2) - 5) are repeated until the stopping criteria are met or the number of iterations meets the predefined maximum number of iterations.
In each iteration, the predefined maximum velocity and position control the changes in velocity and position of each particle.
4. GA and PSO in Seepage Analysis
Determination of the seepage path through the earth dams is one the most difficult problems in geotechnical engineering. This phenomenon is recognized as unconfined seepage problem since the location of phreatic line is not determined. Subsequently, it represents a boundary value problem where at least one of the boundaries is not determined. In order to find the location of phreatic line (free water surface) an iterative procedure is needed  . Varity of optimization techniques are used in unconfined seepage problems to find the phreatic line in the earth dams    ).
In the first step, the governing equation of unconfined seepage problems is briefly introduced and followed by the definition of the objective function and related constraints. Based on the continuity of mass, the governing equation of laminar steady state flow in porous media is introduced as Laplace’s equation and Darcy’s law presents the relationship between the total head, , and velocity of fluid. Regarding the unconfined seepage problems, four types of boundary conditions are introduced (see Figure 1)  :
1) Impermeable boundaries: flow velocity is zero perpendicular to these
boundaries and streamlines are tangential to these boundaries (i.e., ,
is the normal vector on the impermeable boundary).
2) Equipotential lines: Hydrostatic pressure is applied on this boundaries and the total head is constant along these boundaries. and are the total
Figure 1. Schematic illustration of seepage through the earth dams.
head at the upstream and downstream boundaries respectively.
3) Phreatic line: this boundary represents the upper streamline inside the earth dam.
4) Seepage face: the pressure head is zero on this boundary and total head is equal to elevation.
The location where streamline intersects the dam body at the downstream is called exit point. Precise determination of the exit point is a critical step in the design of earth dams. Exit point is a critical location at the downstream where erosion is possible.
To summarize, the mathematical definition of governing equation and boundary conditions for a 2D seepage problem is defined:
where shows the hydraulic permeability in saturated porous media, and are the first and last points on the phreatic line, respectively,
The optimization procedure using evolutionary methods (GA or PSO) begins by generating randomly initial set of individuals which are chromosomes or particles regarding GA or PSO, respectively. The optimization scheme produces m individuals where each individual includes (). Next, the phreatic line is formed based on () for each individual. Then the governing equation (Equation  ) is solved using conventional numerical methods (i.e. Finite element method (FEM) or Finite Volume Method (FVM)). The pressure head is zero along the phreatic line and it represents that the difference between the total head and elevation of the points on the phreatic line is zero. Implementing the least square method and considering Equation (3c), the objective function can be defined as  :
Subsequently, the mathematical constraints are introduced in the optimization procedure as follow:
Fenton and Griffiths  analyzed the unconfined seepage problem using mesh deforming method based on nonlinear FEM. Their study represents an iterative procedure which adjusts the height of the nodes in the mesh geometry. However, Ouria and Toufigh  used linear FEM associated with Nelder-Mead simplex method. They approximated the phreatic line with a 4th degree polynomial and minimized the goal function based on the conditions on the phreatic line. Shahrokhabadi and Toufigh  used the same procedure but they used Natural Element Method (NEM) and GA in their solution procedure. In fact, they used a mesh deformation technique in NEM which is not sensitive to mesh geometry and GA to avoid the local optimum solution. Moreover, they showed that the average error in optimization and number of iterations are noticeably affected by the initial guess in gradient based optimization algorithms while GA is not sensitive to initial guess. Figure 2 represents the comparison among the results using gradient based optimization method and FEM  , nonlinear FEM analysis  , and NEM-GA  . It is shown that linear FEM with Nelder-Mead simplex and NEM-GA lead to similar results for the seepage path. However, the obtained exit point based on nonlinear FEM and NEM-GA is identical and it is reported as 3.63 m.
Precise estimation of the Exit point (see Figure 1) is a challenging problem in unconfined seepage problems and considered as a significant criterion that determines the privilege of a solution procedure to other solution algorithms. Determination of phreatic line in unconfined seepage problem in a rectangular domain is a benchmark problem which has been predominantly studied in literature (i.e.    ). The geometry of problem is defined by a rectangular shape with the length of is 16(m) and height of 24 (m). With respect to Figure 1,
Figure 2. Phreatic line obtained from three different methods (nonlinear FEM, linear FEM Nelder-Mead simplex method, and NEM-GA).
the upstream water level () is 24 (m) and downstream water level is () is 4 (m). In this review, three different methods are investigated, in which evolutionary optimization algorithms (GA and PSO) are used.
First, the results from
In order to present more accurate study, Table 1 shows the exit points reported by different numerical and analytical studies. It is observed that MFS- PSO-TCF represents closer results to analytical solution in comparison with other numerical-optimal solutions. Regarding the presented results in Table 1, MFS-PSO-TCF obtains the closest results to analytical solution introduced by Ozis  .
In the last part of study on application of evolutionary methods in unconfined seepage problems, the comparison between results from MFS-PSO-TCF and experimental results is studied  . Fu and Jin  made a tank of armor plates
Figure 3. Phreatic line in rectangular domain obtained from: NEM-GA, MFS-PSO, and MFS-PSO-TCF.
Table 1. Exit point (C) reposted from different studies.
and glass. The container included two armor plates, one base and two side panels. However, one side was made of glass to serve as observation window. In the container, an earth dam model was placed which was made of coarse sand and estimated porosity of 0.2. In order to protect the model from erosion, some small gravel was used to protect both side slopes. The dam was placed with nature sedimentation and natural density. Subsequently, the saturated hydraulic
conductivity is determined in the steady-state flow case as.
The comparison of experimental results and numerical-optimum results shows that the application of MFS-PSO-TCF is a successful algorithm to find phreatic line and exit point in unconfined seepage problems (See Figure 4). It is observed that NEM-GA, MFS-PSO, and MFS-PSO-TCF are three solution techniques in which a meshfree scheme is used for solving the governing equation
Figure 4. Comparison the results for phreatic line in earth dam based on experimental studies and MFS-PSO-TCF.
whereas an evolutionary optimization method is implemented to find the unknown boundary (phreatic line).
5. GA and PSO in Slope Stability Analysis
Slope stability is a major problem in geotechnical engineering that estimates the failure potential of geo-materials (soil/rock) covering the slopes. The stability is considered as the balance between shear stress and shear strength which are due to soil mass movement and material resistance, respectively. The stability of soil mass covering a slope depends on both slope geometry (slope angle, back slope angle, pore water pressure and etc.) and loading conditions (climatic events, loading and lateral pressure, and etc.).
Slope failure comes about in the failure zone (critical surface) in which the soil strength is less than shear stress. There are numerous techniques to determine the critical slope surface in the classical geotechnical engineering. In order to analyze and design slopes, there is a high demand to find either critical slip surface or optimum configuration    . Figure 5 shows a schematic slope in two dimensions (2D) with three different circular slip surfaces. A given circular slip surface (i) can be determined by its center and associate radius.
For a given slip surface (i), the factor of safety is defined based on the overloading definition in which the shear stress () is compared with the shear strength () along the slip surface  :
Figure 5. Schematic representation of failure surfaces in a given slope.
where is the differential length of a given slip surface. If the slip surface is divided into segments with discretized length, the evaluation of is described as:
The combination of NEM-GA has been used to calculate FOS in slope stability analysis. In the proposed framework, GA is used to find the minimum FOS for circular slip surfaces while NEM is used to analyze the shear stress and strength on the given slip surface  . The critical slip surface which is defined by and can be determined using the flowchart in Figure 6.
In the illustrative example using NEM-GA, a homogenous slope with height (H) 10 m and slopping 45 degrees is analyzed. The shearing strength is calculated based on Mohr Colum criteria:
where c is cohesion, is normal stress, and is friction angle.
The strength parameters for the illustrative example are: friction angle
degrees and, the unit weight, and the elastic parameters and the elastic parameters and Pois-
son’s ratio. Zheng et al.  analyzed the same slope by considering elastic behavior of material and Mohr-Coulomb failure criteria while they used
Figure 6. Solution from work to find critical slip surface.
Finite Element Method (FEM) to evaluate stress-strain behavior of the material. The minimum factor of safety reported from NEM-GA, FEM, and conventional methods based limit equilibrium (i.e. Bishop, Spencer etc.) is depicted in Table 2.
The results show FOS varies from 0.98 (Janbu simplified) to 1.09 (NEM-GA). In order to compare the results from classical limit equilibrium methods, FEM and NEM-GA, slip surface for above mentioned analysis is depicted in Figure 7. The proposed slip surface by FEM and NEM-GA are coarsely identical whereas the limit equilibrium methods find the critical slip surface in different coordinates. The difference in the results could be due to the failure criteria that have been introduced in different methods. Further details will be developed and explained in Discussions and recommendations.
The discussion in the previous methods is limited to circular slip surface while the heuristic algorithms can be expanded for finding the location of critical non-circular failure surfaces  . Zolfaghari et al.  introduced a simple genetic algorithm search to determine non-circular failure surface in slope stability analysis. They used Morgenstern-Price technique  to analysis the slopes. Despite of conventional studies, they considered the effect of internal-slice forces of adjacent slopes. This assumption significantly increases the computational efficiency in finding non-circular failure surfaces. This study suggested a robust algorithm to solve problems including finite and infinite slopes with heterogeneous materials.
In slope stability problems with heterogeneous geology the difference in the results between noncircular and circular slip surfaces is considerable. Figure 8
Table 2. Factor of safety obtained from different methods.
Figure 7. Critical slip surface proposed by classical limit equilibrium methods, FEM, and NEM-GA.
depicts the stability analysis in a finite slope including four different soil layers. Zolfaghari et al.  showed that the FOS in the circular slip surface analysis based circular slip surface using Bishop and Morgenstern methods are 1.475 and 1.5 respectively. However, the analysis with noncircular failure surface shows that FOS is 1.24. Moreover, it shows that failure surface mostly includes the weaker layers which are expected in real problems. They extended this idea in infinite slope stability and Figure 9 represents the results regarding an infinite slope stability analysis. The problem includes three infinite layers with different material properties. In this problem the seismic load is applied to the problem in pseudo-static conditions while the coefficient of quake is 0.1 which is applied horizontally. Regarding the proposed framework, the FOS is different in static and seismic conditions from 1.14 to 0.95, respectively. However, the location of slip surface is identical under both conditions.
In the abovementioned study the control variables were defined within static bounds. However, Cheng  introduced the dynamic bounds which are controlled by the kinematic requirement of failure. This idea later followed by
Figure 8. Comparison of failure surface with a noncircular and circular method.
Figure 9. Results of infinite slope analysis by Zolfaghari et al. (2005).
Cheng et al.  . They introduced a Modified PSO procedure in coupling with slip surface generator introduced by Greco  , Malkawi et al.  , and Cheng  . The MPSO procedure allows obtaining the optimized solution with fewer evaluations in comparison with standard PSO. Since particles with better objective function evaluation are permitted to search (fly) more within the given iteration step, subsequently MPSO suggested a faster approach in solution procedure. Moreover, in MPSO, the number of flies within the whole group of particles is limited. Figure 10 briefly depicts the MPSO algorithm by adding few simple steps to standard PSO algorithm.
In order to compare the results of presented methods, a natural slope with four soil layer that has been studied by Zolfaghari  is considered which implemented Spencer method within genetic algorithm to analyze the proposed slope under 4 loading conditions:
Case I: No water pressure and no earthquake loading in the system.
Case II: Water pressure and no earthquake loading.
Case III: No water pressure but earthquake loading.
Case IV: Water pressure and earthquake loading.
Cheng et al.  used the same number of slices in order to establish a fair comparison between GA, PSO, and MPSO. In all cases the minimum FOS is reported which is considerably lower than GA solution (reported by Zolfaghari et al.  ). Moreover, in all cases, the portion of slip surface which passes through
Figure 10. Proposed MPSO algorithm for non-linear optimization problems (Cheng, 2007).
the weakest layer is more than the solution by Zolfaghari et al.  . Table 3 shows the summary of the results from GA, PSO, MPSO under 4 loading cases.
For illustration intentions, the slope stability analysis results are shown for Case I in Figure 11. It is shown that MPSO suggests more critical conditions in comparison with standard PSO and GA. In order to compare the results for other cases (II, III, IV), interested readers are referred to Cheng et al.  .
6. GA and PSO in Foundation Design
The last part of this study introduces a brief overview on application of evolutionary optimization techniques (GA and PSO) in foundation design. Recently, in the concepts of foundation design a new idea which is called Robust Geotechnical Design (RGD) is the matter of interest    . RGD is interested for obtaining a certain level of design robustness in addition to safety satisfaction and cost requirements. Subsequently, a single best design does not exist anymore and a trade-off among multi objectives may be required. In these conditions, Juang et al.  implemented GA in a multi objective optimization program to
Table 3. Comparison of FOS, resulted from GA, PSO, MPSO.
Figure 11. Critical slip surfaces obtained from MPSO, PSO, and GA.
find a Pareto Frontier which describes a trade-off between cost and robustness at a given safety level. Following Juang et al.  , reliability-based RGD approach using GA is shown in six steps:
Step 1: Classify the design parameters and noise factors in the system (foundation). The dimension of system (B: width; L: Length; and D: Depth) are the design parameters. Subsequently, the design space is finite including M finite designs. On the other hand, uncertain soil parameters (: Friction angle;: cohesion) could be remarked as noise factors.
Step 2: This step includes an inner loop (Figure 12) which is used to compute variation of system response. In addition, coefficient of variation (COV) of noise
Figure 12. Flowchart representing reliability-based RGD approach using GA.
factors for each design is considered in this step. Note: The mean value of each noise is fixed while its COV is allocated based on point estimate method.
Step 3: It includes the reliability analysis for each design in the design space. The deterministic model for ultimate limit state (ULS) and serviceability limit state (SLS) for the system can be counted for probabilistic analysis. ULS and SLS probability may be calculated based on first order reliability method (FORM  ).
Step 4: This step accomplishes the reliability analysis for N times and uses the results to calculate the mean and standard deviation of probability of failure.
Step 5: The outer loop in Figure 12 that using point estimation method (PEM  ) for M times is accomplished in this step. Subsequently, the mean and standard deviation of the failure probability are obtained for all designs.
Step 6: The best design is determined in this step. Regarding the RGD, the multi-objective functions are defined as cost and robustness whereas safety is defined as a constraint. The cost of each design may be computed based on procedure introduced by Wang and Kulhawy  and the multi-objective optimization is done using non-dominated sorting genetic algorithm version II (NSGA- II) introduced by Pratp et al.  .
RGD is based on satisfying three major factors: Safety, cost, and robustness. Therefore, one can cast an optimization problem including robustness and cost as two objective functions which are subjected to safety as a constraint  . In order to find the Pareto Frontier, first we need to generate a random “parent population” in the design space which is similar to step 0 in standard GA (See Section 2). Here, we consider the number of initial “optimal” designs as n in the design space. Subsequently, a series of GA steps including selection, crossover, and mutation are executed to obtain the new “offspring population” (Step 1 in Section 2). Then the refinement of parent population is accomplished through an iterative process (Step 2 in Section 2). For a given generation, the new population is formed through the combination of parent and offspring population. Next, the non-dominant sorting is applied over new population and produces different sets of design points based on the non-dominant ranking. For instance, the best class is recognized as F1; the next is F2 and so on. Subsequently, the best n points are nominated for next generation. Crowding distance technique can be used in the generation on new population with respect to keep diversity in the selected design points when the number of points in exceed n  . This procedure continues until the convergence criterion is satisfied and stable solution obtained. Subsequently, the final solution in the proposed optimization process is determined as Pareto Frontier. Figure 13 shows the converged solution (Pareto Frontier) for a spread footing using NSGA-II optimization.
Pareto Frontier includes an optimal set of possible designs in RGD approach. This set plays a trade-off role in relationship between cost and robustness against uncertainty factors. The designer may opt for a greater robustness at a higher cost or vice versa. On the other hand, if a threshold is predefined for cost then the designer could decide on the most robust design from the Pareto Frontier. In
Figure 13. The Pareto frontier for the design of spread footing obtained using NSGA-II.
addition, based on feasibility robustness one can refine the decision for design objectives (Figure 14). For further information regarding feasibility robustness, interested readers are encouraged to study was done by Juang et al.  .
Application of evolutionary methods in foundation design is not limited just to spread footings. For instance, Lei  implemented PSO on the settlement fitting prediction of highway foundations. The settlement of highway foundations plays a major role in the deformation and stability of the roads. In general, there are three approaches to calculate settlement fittings:
1) Implementation of soil constitutive models to calculate the final settlement of road bases, this type of analysis is based on consolidation theories  .
2) Application of numerical methods (i.e. finite element method/finite difference method) using soil creep constitutive models  .
3) Employment of direct methods which are based on actual measured data (settlement) and presenting a mathematical model to predict settlement in road bases    .
Lei et al.  used the direct method which is based PSO and time series analysis to present a subgrade settlement model for Nanjing-Hangzhou highway. Since settlement in foundation follows a time dependent trend, one can re- present the settlement via a polynomial based time series. Time series analysis has been successfully applied in dynamic data processing phenomena (i.e. Prediction of dynamic deformation of landslides  :
Figure 14. Cost versus feasibility robustness for all designs on the Pareto Frontier.
If m and are determined with PSO, the above model turned to a predictive model which can be used to predict the future deformation prediction  .
In the proposed study, parameters in the model are generally determined by empirical methods. In the process of finding fitting parameters for the model, the measured data has a considerable effect on the predictive model. Lei et al.  used measured data from a part of Nanjing-Hangzhou highway settlement from October 2001 until November 2002 and the parameters of the algorithm are determined in the process of solving fitting model using PSO. The parameters for the optimal model are finally determined as:
Regarding the predictive model presented in Equation  , the learning value, predictive values are compared with the measured data in Figure 15. The comparison shows a good agreement in the analysis of measured data and predicted values.
7. Discussions and Recommendations
Seepage through earth dams or ground water flow represents a problem with an unknown boundary. In order to find a reasonable solution, this study suggests the simultaneous application of a robust technique for solving the governing
Figure 15. Settlement with time of the monitoring data.
equation while an optimization technique is required for finding the location of unknown boundary. Although the application of MFS-PSO-TCF in finding phreatic line/exit point in unconfined seepage problems shows privilege in comparison with presented methods in this study, it is limited only to saturated conditions. However, many studies show the location of phreatic line is a trade-off between unsaturated and saturated conditions  . Subsequently, implementation of evolutionary methods and adequate numerical techniques is suggested for future studies.
Slope stability analysis is a challenging problem in geotechnical engineering and it needs the application of constitutive models which include elasto-plastic- ity  . In NEM-GA framework only elastic behavior of material is considered which can result in error in slope stability analysis. However, it is encouraged to extend the proposed framework to elasto-plastic materials. In addition, it may be useful if the proposed studies could be implemented to retaining structures, MSE walls, and Concave Geosynthetic-Reinforced Soil Structures (CGRSS)  . The third part of this study briefly reviewed the implementation of GA/PSO in the foundation designs. The presented studies suggest the application GA/PSO in spread footings while there are numerous studies dealing with deep foundation and group piles    . Moreover, soil reinforcement and soil stabilization are considerably interested in foundation designs      . Subsequently, the application of evolutionary methods in the optimum design of deep foundations and soil reinforcements is strongly suggested for future studies.
In the presented review, we briefly introduced two well-known evolutionary methods (GA/PSO) and their application in geotechnical engineering. The implementation of GA/PSO in three major geotechnical problems is considered: 1) Unconfined seepage problems; 2) Slope stability analysis; and 3) Spread footing designs.
Regarding the unconfined seepage problems, at least one of the boundaries (phreatic line) is not determined and it is substantial to predict the location of phreatic line and exit point within the solution procedure. Three methods using GA/PSO in their solution procedure are studied in this review: NEM-GA, MFS- PSO, and MFS-PSO-TCF. In the comparison of results from above mentioned solution techniques, it is obtained those MFS-PSO-TCF results in more accurate solution, especially, in the determining exit point location.
In slope stability analysis, determination of the failure surface location is a challenging problem in geotechnical engineering. In this study, we found that GA, PSO, and Modified PSO are successful in finding critical slip surface and subsequently estimation of FOS. In real life scenario, usually critical slip surface is not a semi-circular, especially in problems including heterogeneous soil layers. In the presented study, MPSO performs adequately in the determination of critical slip surface with respect to number of iterations and computational cost.
In the last part of study, the application of GA/PSO in foundation design is investigated. It is shown that the combination of GA and RGD leads to find the trade-off between cost, safety, and design robustness. This criterion is depicted by a Pareto Frontier. Finally, the application of PSO in prediction of settlement in shallow foundations is presented by studying the spread footings in highways. It depicts that PSO is successful to predict the settlement trend based on observed data.
The safety in designs is highly respected in geotechnical engineering. On the other hand, contractors and private sectors prefer economical design with minimum labor effort and construction cost. Subsequently, implementation of robust optimization techniques which establish a trade-off between safety and total cost of geotechnical projects is necessary in practice.
 Rafsanjani, H.N., Shahrokhabadi, S. and Hadjahmadi, A. (2013) The Use of Linear Regression to Estimate the Actual Hourly Production of a Wheel-Type Loader in Construction Projects. International Conference on Sustainable Design, Engineering, and Construction 2012, 727-731.
 Camp, C., Pezeshk, S. and Cao, G. (1998) Optimized Design of Two-Dimensional Structures Using a Genetic Algorithm. Journal of Structural Engineering, 124, 551-559.
 Kuhn, H.W. (2014) Nonlinear Programming: A Historical View. In: Giorgi, G. and Kjeldsen, T.H., Eds., Traces and Emergence of Nonlinear Programming, Springer, Basel, 393-414.
 Keshanchi, B., Souri, A. and Navimipour, N.J. (2017) An Improved Genetic Algorithm for Task Scheduling in the Cloud Environments Using the Priority Queues: Formal Verification, Simulation, and Statistical Testing. Journal of Systems and Software, 124, 1-21.
 Holland, J.H. (1975) Adaptation in Natural and Artificial Systems: An Introductory Analysis with Applications to Biology, Control, and Artificial Intelligence. University of Michigan Press, Ann Arbor.
 Rada-Vilela, J., Johnston, M. and Zhang, M. (2014) Population Statistics for Particle Swarm Optimization: Resampling Methods in Noisy Optimization Problems. Swarm and Evolutionary Computation, 17, 37-59.
 Fernandez-Marquez, J.L. and Arcos, J.L. (2010) Adapting Particle Swarm Optimization in Dynamic and Noisy Environments. IEEE Congress on Evolutionary Computation, Barcelona, 18-23 July 2010, 1-8.
 Brunner, P., Cook, P.G. and Simmons, C.T. (2009) Hydrogeologic Controls on Disconnection between Surface Water and Groundwater. Water Resources Research, 45, W01422.
 Wang, Y., Hu, M., Zhou, Q. and Rutqvist, J. (2016) A New Second-Order Numerical Manifold Method Model with an Efficient Scheme for Analyzing Free Surface Flow with Inner Drains. Applied Mathematical Modelling, 40, 1427-1445.
 Fenton, G.A. and Griffiths, D.V. (1997) A Mesh Deformation Algorithm for Free Surface Problems. International Journal for Numerical and Analytical Methods in Geomechanics, 21, 817-824.
 Shahrokhabadi, S. and Toufigh, M.M. (2013) The Solution of Unconfined Seepage Problem Using Natural Element Method (NEM) Coupled with Genetic Algorithm (GA). Applied Mathematical Modelling, 37, 2775-2786.
 Westbrook, D.R. (1985) Analysis of Inequality and Residual Flow Procedures and an Iterative Scheme for Free Surface Seepage. International Journal for Numerical Methods in Engineering, 21, 1791-1802.
 Chen, J.T., Hsiao, C.C., Chiu, Y.P. and Lee, Y.T. (2007) Study of Free-Surface Seepage Problems Using Hypersingular Equations. Communications in Numerical Methods in Engineering, 23, 755-769.
 Chaiyo, K., Rattanadecho, P. and Chantasiriwan, S. (2011) The Method of Fundamental Solutions for Solving Free Boundary Saturated Seepage Problem. International Communications in Heat and Mass Transfer, 38, 249-254.
 Shahrokhabadi, S. and Ahmadi, A. (2013) Method of Fundamental Solution (MFS) coupled with Particle Swarm Optimization (PSO) to Determine Optimal Phreatic Line in Unconfined Seepage Problem. Scientia Iranica, 20, 1327.
 Shahrokhabadi, S., Vahedifard, F. and Yarahmadian, S. (2016) Integration of Thiele Continued Fractions and the Method of Fundamental Solutions for Solving Unconfined Seepage Problems. Computers & Mathematics with Applications, 71, 1479-1490.
 Hang, M. (2002) A 3D Slope Stability Analysis Method Assuming Parallel Lines of Intersection and Differential Straining of Block Contacts. Canadian Geotechnical Journal, 39, 799-811.
 Shahrokhabadi, S., Khoshfahm, V. and Rafsanjani, H.N. (2014) Hybrid of Natural Element Method (NEM) with Genetic Algorithm (GA) to Find Critical Slip Surface. Alexandria Engineering Journal, 53, 373-383.
 Zheng, H., Liu, D.F. and Li, C. (2005) Slope Stability Analysis Based on Elasto-Plastic Finite Element Method. International Journal for Numerical Methods in Engineering, 64, 1871-1888.
 Sun, J., Li, J. and Liu, Q. (2008) Search for Critical Slip Surface in Slope Stability Analysis by Spline-Based GA Method. Journal of Geotechnical and Geoenvironmental Engineering, 134, 252-256.
 Zolfaghari, A.R., Heath, A.C. and McCombie, P.F. (2005) Simple Genetic Algorithm Search for Critical Non-Circular Failure Surface in Slope Stability Analysis. Computers and Geotechnics, 32, 139-152.
 Cheng, Y.M., Li, L., Chi, S.C. and Wei, W.B. (2007) Particle Swarm Optimization Algorithm for the Location of the Critical Non-Circular Failure Surface in Two-Dimensional Slope Stability Analysis. Computers and Geotechnics, 34, 92-103.
 Malkawi, A.I.H., Hassan, W.F. and Sarma, S.K. (2001) Global Search Method for Locating General Slip Surface Using Monte Carlo Techniques. Journal of Geotechnical and Geoenvironmental Engineering, 127, 688-698.
 Juang, C.H., Wang, L., Liu, Z., Ravichandran, N., Huang, H. and Zhang, J. (2013) Robust Geotechnical Design of Drilled Shafts in Sand: New Design Perspective. Journal of Geotechnical and Geoenvironmental Engineering, 139, 2007-2019.
 Wang, Y. and Kulhawy, F.H. (2008) Economic Design Optimization of Foundations. Journal of Geotechnical and Geoenvironmental Engineering, 134, 1097-1105.
 Deb, K., Pratap, A., Agarwal, S. and Meyarivan, T.A.M.T. (2002) A Fast and Elitist Multiobjective Genetic Algorithm: NSGA-II. IEEE Transactions on Evolutionary Computation, 6, 182-197.
 Chan, C.M., Zhang, L.M. and Ng, J.T. (2009) Optimization of Pile Groups Using Hybrid Genetic Algorithms. Journal of Geotechnical and Geoenvironmental Engineering, 135, 497-505.
 Hekmat, A., Yusop, Z.B., Kashefizadeh, M. and Gholizadeh, R. (2013) Three Dimensional Dynamic Analysis of Interaction Soil-Pile with Coupling Finite Element Analysis (FEA) and Boundary Element Method (BEM). Caspian Journal of Applied Sciences Research, 2, 269-275.
 Armaghani, D.J., Raja, R.S.N.S.B., Faizi, K. and Rashid, A.S.A. (2015) Developing a Hybrid PSO-ANN Model for Estimating the Ultimate Bearing Capacity of Rock-Socketed Piles. Neural Computing and Applications, 28, 1-15.
 Shahrokhabadi, S. and Nazeryzadeh, N. (2013) Efficacy of Resinous Polypropylene (PP) Fibers on Strength Behavior of Reinforced Soils. In: Advanced Materials Research, Vol. 787, Trans Tech Publications, 75-80.
 Chaiyaput, S., Bergado, D.T. and Artidteang, S. (2014) Measured and Simulated Results of a Kenaf Limited Life Geosynthetics (LLGs) Reinforced Test Embankment on Soft Clay. Geotextiles and Geomembranes, 42, 39-47.
 Bazne, M.O., Vahedifard, F. and Shahrokhabadi, S. (2015) The Effect of Geonet Reinforcement on Bearing Capacity of Low-Compacted Soft Clay. Transportation Infrastructure Geotechnology, 2, 47-63.