J. Cent. South Univ. (2012) 19: 1764-1770

DOI: 10.1007/s11771-012-1203-2

Inverse procedure for determining model parameter of soils using real-coded genetic algorithm

LI Shou-ju(李守巨)¹, SHAO Long-tan(邵龙潭)¹, WANG Ji-zhe(王吉喆)², LIU Ying-xi(刘迎曦)¹

1. State Key Laboratory of Structural Analysis for Industrial Equipment (Dalian University of Technology),Dalian 116024, China;

2. The Second Affiliated Hospital of Dalian Medical University, Dalian 116027, China

? Central South University Press and Springer-Verlag Berlin Heidelberg 2012

Abstract: The hybrid genetic algorithm is utilized to facilitate model parameter estimation. The tri-dimensional compression tests of soil are performed to supply experimental data for identifying nonlinear constitutive model of soil. In order to save computing time during parameter inversion, a new procedure to compute the calculated strains is presented by multi-linear simplification approach instead of finite element method (FEM). The real-coded hybrid genetic algorithm is developed by combining normal genetic algorithm with gradient-based optimization algorithm. The numerical and experimental results for conditioned soil are compared. The forecast strains based on identified nonlinear constitutive model of soil agree well with observed ones. The effectiveness and accuracy of proposed parameter estimation approach are validated.

Key words: parameter estimation; real-coded genetic algorithm; tri-dimensional compression test; gradient-based optimization

1 Introduction

The identification of the material parameters, for a given constitutive model, can be seen as an inverse formulation. This approach consists of an optimization problem where the objective function is to minimize the gap between the experimental and the numerical results. The optimization variables are the material parameters that appear in the constitutive model [1]. This type of problem is usually called parameter estimation and is often solved by deterministic and gradient-based optimization methods. Unfortunately, the solution to the problem using these methods is usually around the local minima if there are more than one minimum available [2]. To overcome these local optimum points, a random search technique is developed. The derivative-free algorithms, also called direct search algorithms, are generally based on simple strategies and do not require the calculation of derivatives. The simplicity is the main attribute. However, direct search algorithms undergo the problem of converging to local minimums, and are also somehow user-dependent. The convergence of these algorithms is very time-consuming and involves the comparison of each trial solution with the best previous solution.

In order to perform the material parameter identification in a robust manner, the genetic algorithm was implemented [3-5]. The genetic algorithm (GA) is a selective random search algorithm designed to achieve a global optimum within a large space of solutions. Based on the type of modeling the natural evolution, GA can search for optimal or near-optimal solutions for an optimization problem over the search domain, and has superior performance over the local optimal techniques. This is due to searching for solution from only one single direction on the search space. However, GA can be regarded as a search method based on multiple  directions, and it inherently utilizes crossover and mutation operations for the searching process. This implies that it is easier to escape from a local minimum [6]. KHALIK et al [7] presented a parameter identification procedure based on real-coded genetic algorithm. GALDI et al [8] proposed a novel technique to identify the thermal parameters used for the estimation of the hot-spot temperature. TUTKUN [9] presented parameter estimation procedure in mathematical models using the real coded genetic algorithm approach. RECHEA et al [2] proposed two numerical procedures to quantitatively identify a set of constitutive parameters that best represent observed ground movement data associated with deep excavations in urban environments. LI and LIU [10] proposed a hybrid genetic algorithm for identifying vibration dynamic characteristics. PERRY  et al [11] proposed a modified GA strategy to improve the accuracy and reduce the computational time for parameter identification of multiple degree-of-freedom structural systems. HWANG et al [12] proposed a estimation method to determine the effective elastic constants of a woven composite plate and two personal computer boards (PCBs). HWANG and HE [13] developed a novel adaptive real-parameter simulated annealing genetic algorithm that maintains the merits of genetic algorithm and simulated annealing. Adaptive mechanisms were added to insure the solution quality and to improve the convergence speed.

The aim of this work is to develop a hybrid genetic algorithm to model parameter identification in a robust and fast manner, and to present a new approach to compute the calculated strains by multi-linear simplification approach instead of finite element method (FEM).

2 Parameter estimation procedures using real-coded genetic algorithms

The inverse problem is often ill-posed. The ill-posedness is generally characterized by the nonuniqueness and instability of the identified parameters. The instability of the inverse solution stems from the fact that small errors in measurement data will cause serious errors in the identified parameters. The uniqueness problem has a great practical importance, because in the case of nonuniqueness, the identified parameters will differ according to the initial estimate of the parameters, and there will be no reason for the estimated parameters to be close to the true parameters [14]. Real-coded GA has been also introduced to a wide variety of applications in recent years [15-17]. All genes in a chromosome used in real-coded GA are real numbers. It is more suitable to directly represent genes as real values for most of real-optimization problems during genetic operations. Because the procedures of binary coding for a real number may suffer for the loss of precision depending on the number of the used bits, expectably, it will be very complicated and difficult to implement if the numerical values are large and have the decimal. For the real-coded GA, the length of chromosomes becomes much shorter than that using the traditional way [1].

The parameters estimation of mechanical model using a modified real-coded GA starts with a population with many chromosomes which are generated randomly. Each chromosome in the population represents a set of possible solution to the optimization problem of parameters estimation. The chromosomes are then evolved to generate better offspring according to the values of objective function by applying three genetic operations. The experimental database provides information for specific stress-strain paths, which allows parameter identification by minimizing the difference between the stress-strain values numerically simulated and the experimental curves. The material parameter identification is then a minimum problem of objective function, where the objective function can be defined as the difference between the results of the numerical simulation and the experimental data. The objective function can be expressed as

                    (1)

where ε_o is the measured strain vector, as shown in Fig. 1; ε_c is the computing strain vector, which is related to the identified parameter vector m. The optimization variables (here, model parameters) should be subjected to the bound constraints:

m_lu                                   (2)

where m_l and m_u are the lower and upper bounds of m, respectively.

Fig. 1 Experimental stress-strain curves for triaxial tests on soil

This objective function clearly depends on the measured data and the parameters of model. The objective function becomes complex, such as non-convex, or even multi-modal if errors contained in the model equation or/and errors in the measurement data are large. In such a case, the solution may vibrate or diverge when conventional gradient-based optimization methods are used, which gives rise to the necessary for a robust optimization method such that a stable convergence is always achieved. The solution of the inverse problem is contained in obtaining a minimum of an objective function which is defined taking into account the mathematical structure of the material model and asset of experimental data. The probability of survival of any individual is determined by its fitness when evaluating any individual, and this fitness function is computed for the individual. The fitting function for every individual is defined as

                                    (3)

where f_j is the fitting function of the j-th individual of population. There are two well-known selection mechanisms, i.e., roulette wheel and tournament selections used for reproduction operation. The roulette wheel selection can be visualized by imagining a wheel, where each chromosome occupies an area that is related to its value of objective function. When a spinning wheel stops, a fixed marker determines which chromosome will be selected to reproduce into the mating pool. When the roulette wheel is rotated, an individual has a chance of being selected corresponding to its share. One of the most commonly used is the roulette wheel selection, where individuals are extracted in probability following a Monte Carlo procedure. The extraction probability of each individual is proportional to its fitness as a ratio to the average fitness of all the individuals. In the selection process, the reproduction probabilities of individuals are given by their relative fitness:

                                 (4)

Another effective selection procedure is the rank- based selection. The problem of fitness-proportional selection is that it is directly based on fitness. In most cases, we cannot define an accurate measure of goodness of a solution, so the assigned fitness value does not express exactly the quality of a solution. Still, an individual with better fitness value is a better individual. In rank based selection, the individuals are ordered according to their fitness. The individuals are then selected with a probability based on some linear functions of their rank. Fitness-proportional selection is commonly used in the population reproduction [18]. When using this selection method, a solution has a probability of selection directly proportional to its fitness. The mechanism that allows fitness proportional selection is similar to a roulette wheel that is partitioned into slices. Each individual has a share directly proportional to its fitness:

p_k=q_max-(k-1)r                               (5)

where p_k is reproduction probability of the k-th individual, k is the rank of k-th the individual, q_max is reproduction probability of the best individual, in which it can be chosen as q_max=0.99.

The parameter r is defined as

                              (6)

where N is population size; q_min is reproduction probability of the worst individual, in which it can be chosen as

                                 (7)

Recombination, also called crossover, is a process combined by information contained in two candidate solutions. In the recombination, each individual is first paired with an individual at random. New individuals are generally created as offspring of two parents (as such, crossover being a binary operator). One or more so-called crossover points are selected (usually at random) within the chromosome of each parent, at the same place in each. The parts delimited by the crossover points are then interchanged between the parents. The individuals resulting in this way are the offspring. In general, there are three different genetic operation modes of sexual recombination: one-point crossover, two-point crossover, and uniform crossover. The one-point crossover has only one pair of parameters to do the crossover process. Suppose that the parents could be described as

                     (8)

                    (9)

If the crossover point is k, then the new solutions are

          (10)

          (11)

Reproduction and crossover make the current population gradually stronger and individuals of the population resemble each other as the genetic process continues. This is expected trend, but if individuals become similar in early generations, the possible solutions will be trapped by local optima due to lack of genetic diversity [18]. Since the reproduction and the crossover operators have no straight forward mechanism to avoid local optima, the mutation operation is proceeded to increase the diversity. The mutation operation only occurs if a randomly generated number between 0 and 1 is less or equal to the mutation probability.

A new individual is created by making modifications to one selected individual. The modifications can consist of changing one or more values in the representation or in adding/deleting parts of the representation. In genetic algorithms, mutation is a source of variability, and is applied in addition to crossover and reproduction. Mutation is a process by which vectors resulting from selection and recombination are perturbed [19]. The mutation is conducted with only a small probability by definition:

                         (12)

where  is a random number, 0<<1;  and  are the lower bound and upper bound of the parameter x_k, respectively. The mutation operation follows the crossover and provides a possible mutation on some chosen chromosomes. Only smaller chromosomes in the current population will be randomly selected to mutate.

One feature that is currently missing in this selection procedure is that it does not guarantee that the best individual always survives into the next generation, particularly when many individuals have fitness close to that of the best individual. The elitist strategy, where the best individual always survives into the next generation on behalf of the worst individual, can compensate for some disadvantages of missing the best individual in selection operation or mutation operation. With the elitist strategy, the best individual always moves in a descent direction, thereby a stable convergence is obtained [20]. The gradient search algorithm adopted in genetic algorithm is the most popular quasi-Newton method with the Broyden-Fletcher-Goldfarb-Shanno (BFGS) algorithm. The gradient-based algorithm generates the following parameter sequence:

m^k+1=m^k-ρ^kd^k                               (13)

A^kd^k=g^k                                    (14)

A^k=[S(m^k)]^T[S(m^k)]                           (15)

g^k=[S(m^k)]^Te^k                               (16)

e^k=[ε_c(m^k)-ε_o]                               (17)

where e is error vector; S is Jacobian matrix of strain with respect to model parameters m; ρ^k is step size; d ^k is Gauss-Newton direction vector, k denotes number of iteration of gradient search. The main steps for parameter identification using genetic algorithm are shown as follows.

Step 1: Choose the size of population, the crossover probability, mutation probability, and stopping criterion.

Step 2: Determine the identified parameter domains according to prior information.

Step 3: Randomly generate an initial population of candidate solutions.

Step 4: Define a fitness function to measure the performance of an individual in the problem domain.

Step 5: Compute model responses with given model parameters by using numerical method.

Step 6: Calculate the fitness of each individual based on observed data and model responses computed by numerical method.

Step 7: Execute recombination operation by using real-coded genetic algorithm.

Step 8: Create new individuals by mutation operation.

Step 9: Execute selection operation according to the ranked selection.

Step 10: Perform elitist strategy in order to keep current best individual from missing and accelerate convergence speed of inverse problem.

Step 11: Replace the initial (parent) population with the new (offspring) population.

Step 12: Execute stopping criterion. If stopping criterion can not be reached, then, go to Step 5; otherwise, the inversion computation stops and best solution is recorded as the solution of the inverse problem.

3 Application of developed genetic algorithm to model parameter identification of soils

The stress-strain behavior of conditioned soil depends on a number of different factors including density, water content, structure, drainage conditions, strain conditions, duration of loading, stress history, confining pressure, and shear stress, as shown in Fig. 1. The hyperbolic equation of soils proposed by KONDNER [21] is expressed as

                            (18)

where σ₁ and σ₃ are the major and minor principal stresses, respectively; ε is the axial strain; a and b are constants whose values may be determined experimentally. The values of tangent modulus for any stress condition may be expressed as [22]

      (19)

where C and φ are the Mohr-Coulomb strength parameters; R_f is the failure ratio, which always has a value less than unity, and the value of R_f has been found to be between 0.75 and 1.00; p_a is the atmospheric pressure; K is a modulus number; n is the exponent determining the rate of variation of E_t with σ₃.

According to experimental stress-strain curves for triaxial tests on soil, the relationship between axial stress and strain can be simplified as multi-linear elastic model, as shown in Fig. 2. When σ₃ is constant, the relationship between the stress difference and axial strain can be approximately expressed as

                (20)

where σ_d is stress difference,  is increment of stress difference.

σ_d=σ₁-σ₃                                   (21)

                             (22)

Fig. 2 Simplified multi-linear elastic model when σ₃ is constant

The model solution, which depends on identified parameters because tangent modulus E_t includes model parameters, can be computed using Eq. (20) instead of using finite element method, which will save much computing time in parameter inversion.

To develop inverse procedure for evaluating the model parameters C, φ, K, n and R_f, a number of tests in laboratory have been conducted on a uniform conditioned sand. These tests are standard triaxial compression tests, which are used to determine the parameters representing the behaviors of the conditioned sand. The specimen tested is initially 39.1 mm in diameter and 80 mm in height, as shown in Figs. 3 and  4.

Fig. 3 Sketch map of laboratory test of conditioned soil for three-dimensional compression

Main characteristic parameters are listed as follows: population size 100, crossover probability 0.85, mutation probability 0.03, and maximum generation 100. Feasible domains of identified parameters are listed in Table 1 which are based on priori information. Based on genetic algorithm and objective function Eq. (1), the model parameters of soil are identified and listed in Table 2. Figure 5 shows the variation of fitness of best solution and average population solution versus generation.

Fig. 4 Dimension of test sample (Unit: mm)

Table 1 Feasible domains of identified parameters

Table 2 Identified model parameters of soil

Fig. 5 Variation of fitness of best solution and average population solution versus generation

In order to validate the effectiveness of proposed parameter identification procedure, the numerical and experimental results for conditioned soil are compared. The forecast strains based on identified nonlinear constitutive model of soil agree well with observed ones, as shown in Figs. 6-8. The identified model parameters are compared with these obtained by using artificial neural network and Gauss-Newton optimization [23-24]. The forecasting precision is higher and computing time is less than other parameter identification procedures.

Fig. 6 Comparison between forecast strains and observed ones at σ₃=100 kPa

Fig. 7 Comparison between forecast strains and observed ones at σ₃=200 kPa

Fig. 8 Comparison between forecasted strains and observed ones at σ₃=300 kPa

4 Conclusions

1) The multi-linear simplification approach instead of finite element method for computing model solutions can save much time in parameter inversion process because the forward problem should be solved 10 000 times. Lots of time will be spent in computing forward problem if finite element method is used.

2) The parameter inversion procedure based on genetic algorithm poses good robustness and fits to deal with inverse problem with observing noise because there are many local minimum values for objective function.

3) The investigation shows that the proposed parameter inversion procedure has high prediction precision. The forecast strains based on identified nonlinear constitutive model of soil agree well with observed ones.

References

[1] CHAPARRO B M, THUILLIER S, MENEZES L F, MANACH P Y. Material parameters identification: Gradient-based, genetic and hybrid optimization algorithms [J]. Computational Materials Science, 2008, 44(2): 339-346.

[2] RECHEA C, LEVASSEUR S, FINNO R. Inverse analysis techniques for parameter identification in simulation of excavation support systems [J]. Computers and Geotechnics, 2008, 35(3): 331-345.

[3] LIU G R, LEE J H, PATERA A T. Inverse identification of thermal parameters using reduced-basis method [J]. Comput Methods Appl Mech Engrg, 2005, 194(27/28/29): 3090-3107.

[4] HE R S, HWANG S F. Damage detection by an adaptive real-parameter simulated annealing genetic algorithm [J]. Computers and Structures, 2006, 84(31/32): 2231-2243.

[5] NYARKO E K, SCITOVSKI R. Solving the parameter identification problem of mathematical models using genetic algorithms [J]. Applied Mathematics and Computation, 2004, 153(3): 651-658.

[6] WEI D C. An improved real-coded genetic algorithm for parameters estimation of nonlinear systems [J]. Mechanical Systems and Signal Processing, 2006, 20(1): 236-246.

[7] KHALIK M A, SHERIF M, SARAYA S. Parameter identification problem: Real-coded GA approach [J]. Applied Mathematics and Computation, 2007, 187(2): 1495-1501.

[8] GALDI V, IPPOLITO L, PICCOLO A, VACCARO A. Parameter identification of power transformers thermal model via genetic algorithms [J]. Electric Power Systems Research, 2001, 60(2): 107-113.

[9] TUTKUN N. Parameter estimation in mathematical models using the real coded genetic algorithms [J]. Expert Systems with Applications, 2009, 36(2): 3342-3345.

[10] LI S J, LIU Y X. Identification of vibration loads on hydro generator by using hybrid genetic algorithm [J]. Acta Mechanica Sinica, 2006, 22(6): 603-610.

[11] PERRY M J, KOH C G, CHOO Y S. Modified genetic algorithm strategy for structural identification [J]. Computers and Structures, 2006, 84(8/9): 529-540.

[12] HWANG S F, WU J C, HE R S. Identification of effective elastic constants of composite plates based on a hybrid genetic algorithm [J]. Composite Structures, 2009, 90(2): 217-224.

[13] HWANG S F, HE R S. Improving real-parameter genetic algorithm with simulated annealing for engineering problems [J]. Advances in Engineering Software, 2006, 37(6): 406-418.

[14] YEH W W. Review of parameter identification procedures in groundwater hydrology: The inverse problem [J]. Water Resour Res, 1986, 22(2): 95-108.

[15] FRANULOVIC M, BASAN R, PREBI I. Genetic algorithm in material model parameters identification for low-cycle fatigue [J]. Computational Materials Science, 2009, 45(2): 505-510.

[16] CHANG W D. Nonlinear system identification and control using a real-coded genetic algorithm [J]. Applied Mathematical Modelling, 2007, 31(3): 541-550.

[17] CHANG W D. Coefficient estimation of IIR filter by a multiple crossover genetic algorithm [J]. Computers and Mathematics with Applications, 2006, 51(9/10): 1437-1444.

[18] XU T, ZUO W J, XU T S. An adaptive reanalysis method for genetic algorithm with application to fast truss optimization [J]. Acta Mechanica Sinica, 2010, 26(2): 225-234.

[19] REN Z W, SAN Y, CHEN J F. Hybrid simplex-improved genetic algorithm for global numerical optimization [J]. Acta Mechanica Sinica, 2007, 23(1): 91-95.

[20] LI S J, LIU Y X. Inverse determination of model parameters of nonlinear heat conduction problem using hybrid genetic algorithm [J]. Nonlinear Dynamics and Systems Theory, 2007, 7(1): 23-34.

[21] KONDNER R L. Hyperbolic stress-strain response: cohesive soil [J]. Journal of the Soil Mechanics and Foundations Division, 1963, 89(SM1): 115-143.

[22] DUNCAN J M, CHANG C. Nonlinear analysis of stress and strain in soils [J]. Journal of the Soil Mechanics and Foundations Division, 1970, 96(SM5): 1629-1653.

[23] SHANGGUAN Z C, LI S J, LUAN M T. Estimating model parameters of conditioned soils by using artificial network [J]. Journal of Software, 2010, 5(3): 296-302.

[24] YU H, LI S J, DUAN H X, LIU Y X. A procedure of parameter inversion for a nonlinear constitutive model of soils with shield tunneling [J]. Computers and Mathematics with Application, 2011, 61(8): 2005-2009.

(Edited by DENG Lü-xiang)

Foundation item: Project(2007CB714006) supported by the National Basic Research Program of China; Project(90815023) supported by the National Natural Science Foundation of China

Received date: 2011-03-28; Accepted date: 2011-05-23

Corresponding author: LI Shou-ju, Professor, PhD; Tel: +86-15040637078; E-mail: lishouju@dlut.edu.cn