1. Introduction

Appropriate constitutive models of soils are the key to a successful prediction of the behaviour of geotechnical structures. A large number of models based on various constitutive theories have been proposed in the last three decades. All of them assume a-priori, a mathematical framework of the model and the material parameters corresponding to the assumed framework have to be identified by physical material tests. Many material parameters in complex constitutive theories have no physical meaning, are difficult to determine and have to be identified by trial and error from numerical simulations. In spite of this, many features of soil behaviour such as stiffness at small strains, higher stiffness on reversal of stress path, influence of rotation of principal stress axes etc. have not been captured in a single model.

Thus, it is likely that models of greater complexity will have to be developed in the future. In recent years a number of applications of Artificial Neural Networks (ANNs) leading to Neural Network based Constitutive Models (NNCMs) have been proposed by a number of researchers as in [1], [2], [3]. There are three different ways in which NNCMs can possibly be used by engineers with considerable advantage. These are:

1. Firstly, NNCMs can be developed for any material from the raw test data without invoking any

constitutive theory as in [3]. This approach has many advantages, the most important being that one does not necessarily have to identify material parameters of the model. However, if one does need to identify them to get a feel of the conventional engineering prameters, they can be identified by carrying out what is known as ‘virtual tests’.

2. Secondly, NNCMs can be trained from incremental load and displacement data of ‘structures’. Here, by the term ‘structures’, it is meant that solids of arbitrary shape subjected to monotonically increasing loads having a non-uniform states of strain and stress. Thus, a cylindrical specimen of a geo-material having glued rigid platens, subjected to uniaxial load is a structure. This application of NNCMs leads to ‘intelligent finite elements’ as in [4] and is available for condition monitoring of real structures. It can also be used to identify material parameters for complex materials such as masonry from structural tests, as in [5].

3. Thirdly, since many constitutive models are very complex, incorporating them in a finite element code and using them for solving real life problems may not be a trivial task. Here, using ‘synthetic data’ generated from systematic exploration of strain space and corresponding stress response, a NNCM can be trained and plugged in a Finite Element (FE) code. This will certainly leads to computational efficiency since the computation of stress increment for a strain increment from a trained NNCM is almost instantaneous. Whilst the first two categories of applications have been reported by a number of researchers, there appear to have been no applications in this category.

This paper belongs to the third category of applications. The objective is to demonstrate how synthetic data from any constitutive model can be used to successfully train a NNCM. We choose two constitutive models for illustration. The first is a well-known model viz. the Hardening Soil Model (HSM) available in the commercial code PLAXIS. This model uses parameters obtained from established engineering practice and is suitable for situations of monotonic loading. The second model is a complex Two-surface Deviatoric Hardening Model (TDHM) in the multilaminate framework for cyclic loading developed as in [6] for the purpose of illustration. Section 2 presents briefly details of the architecture of NNCMs adopted here. Section 3 provides details of HSM and generation of ‘synthetic data’ along triaxial compression and extension paths under loading and unloading conditions. Comparison of the response of a trained NNCM with the synthetic data used for training is also given here. Section 5 gives details of TDHM and generation of ‘synthetic data’ for one-way cyclic triaxial compression loading under various specified volumetric strains including zero volumetric strain. The latter condition represents undrained conditions encountered in the liquefaction problems of loose sands under seismic loading. In this section the response of a trained NNCM with the synthetic data for a number of cycles of loading and unloading is also compared. Section 6 gives conclusions and indicates the implications for computational efficiency in three-dimensional practical problems.

2. Artificial Neural Network Based Constitutive Models

Artificial Neural Networks (ANNs) are pattern recognition algorithms using which relationship between a set of ‘causes’ and ‘effects’ can be captured. Any set of numeric data can be used to discover the pattern in it, if it exists. In the past two decades, a large number of applications of this methodology in almost all disciples of physical and biological sciences have been reported. Here we concentrate on developing nonlinear stress- strain relations for geo-materials.

2.1. Methodology and Architecture

ANNs can be used to simulate stress-strain response of any material by using appropriate data. In this case, the components of strain rates become the causes whilst the resulting stress rates are the effects. Description of stress-strain behaviour by a NNCM does not require checking for plastic flow, computation of flow vector, updating and reconstitution tangential stress integration matrix. Full details are contained in [2], [3], [4], [5]. Here a brief description for completeness and continuity is given.

Incremental (as distinct from total) stress vector can be computed from the corresponding incremental strain vector as follow:

(1)

Where:

(2)

(3)

Here, in the terminology of ANNs, NNCM stands for a neural network based constitutive model trained from appropriate incremental stress-strain data at various strain/stress levels. Thus, the input parameters must include current components of stress (σ) and increments of strain (dε) whilst the required outputs are increment of stress (dσ). In many geotechnical problems, soils are subjected to cyclic and transient loads. Even in quasi-static problems many elements are subjected to unloading and reloading. For constitutive models to be valid, in such situations, it is proposed to adopt ‘intrinsic time’ as in [7], [8] or the current length of strain trajectory, ξ, as an independent input parameter. Mathematically, intrinsic time, ξ, is defined as follows:

(4)

Where, dξ is an increment of deviatoric strain defined as:

(5)

ξ, is a monotonically increasing positive parameter. The above definition can be changed to include volumetric strain as well as real time as is the case in ‘endrochronic’ theories of constitutive behaviour.

Taking into consideration the above points in mind, the NNCM adopted in this paper has strain increments, stresses and ‘intrinsic time’ are the input variables whilst increments of stress are the output variables. For two-dimensional analyses, the optimal architecture of the NN has been obtained by trial and error and is presented in Figure 1. It is constituted by 9 input nodes two hidden layer of 18 and 8 nodes respectively and 4 output nodes. In order to train the neural network the resilient back-propagation (RPROP) algorithm is used. RPROP was first proposed by Reidmiller as in [9] and implemented by Shin as in [10]. It is a local adaptive learning scheme based on the standard back-propagation framework.

Figure 1. Architecture of NNCM [9-18-8-4] for two-dimensional analysis

2.2. Data Enrichment

The strain-stress pairs from the triaxial tests are actually principal stresses and strains since no shear stresses/strains are involved. If such data were used for training, NNCM will have to extrapolate wherever shear stress/strain components are involved. This will lead to large inaccuracies in stress-strain response of the NNCM.

To overcome this limitation, in [3] proposed a data enrichment strategy. They created additional data by transformation of stresses and strains in a co-ordinate other than the principal axes of soil sample being tested. This naturally leads to non-zero shear components of stresses and strains.

This method produces a large amount of training data depending on the number of co-ordinate transformations chosen to generate data. Among the expanded data there are many duplicated strain-stress pairs so an additional process of ‘data pruning’ is adopted. In this paper, the procedure described above has been adopted for data enrichment with an incremental angle, Δθ equal to 5º in order to rotate the strain-stress axes from -45^o to + 45^o.

3. Hardening Soil Model (HSM)

The Hardening Soil Model described in PLAXIS^TM Manual is a nonlinear elastic-plastic model with Mohr Coulomb failure criterion. It is an enhanced version of the nonlinear elastic hyperbolic model of Duncan & Chang as in [11] with deviatoric hardening operating on Mohr Coulomb yield surface. A non-associated flow rule defined by a dilatancy angle smaller than the peak friction angle is adopted. It is generally applicable to loose to medium dense sands and normally consolidated to lightly overconsolidated soils. The model captures apparently strong nonlinearity prior to failure, which is a drawback of conventional linear elastic-plastic models. A cap in deviatoric stress – mean effective stress space is also included.

3.1. Synthetic Data Generation

For heuristic purpose, we have chosen typical parameters for medium dense as given in Table 1. Though many parameters will be familiar to engineers, one should refer to the PLAXIS^TM software manuals for full details.

The data in Table 1 are used to generate stress-strain response of the sand under various experimental configurations, viz. Triaxial Loading in Compression (LC), Triaxial Loading in Extension (LE), Unloading in Compression (UC) and Unloading in Extension (UE). Data have been generated under stress controlled drained conditions. A single finite element subjected to uniform stress conditions was used in analysis with PLAXIS^TM software with HSM model with parameters given in Table 1. The stress paths, in deviatoric stress, q, and effective mean stress, p^’, space are shown in Fig. 2. Three different confining pressures of 50 kPa, 100 kPa and 150 kPa were used for each of the above stress paths. In engineering practice, if such test data were available for a soil, one would perhaps term them as ‘extensive’. These stress-strain data obtained from these simulations were then used for training the NNCM.

Table 1. Chosen parameters for sand for HSM

Figure 2. Stress paths in q – p’ space used for generating synthetic data

Plots of q versus axial strain, ε_yy, as predicted by the trained NNCM and the original HSM data used for training, for various stress paths are shown in Figures 3 and 4. An excellent match is observed confirming that NNCM has been adequately trained. This, however, is not surprising since prediction of response is made for the confining pressure, which was also used in training data. We need to check if the response at confining pressures within the training range as well as outside that used in training is also reasonable. This issue is discussed in the next section of the paper.

Figure 3. Graph of q versus εyy under LC and LE conditions for a confining pressure of 100 kPa

Figure 4. Graph of q versus εyy under UC and UE conditions for a confining pressure of 100 kPa

3.2. Model Validation

In this part of the paper, we present the prediction of the trained NNCM for two additional triaxial tests with confining pressures of 10 kPa & 200 kPa which are out of the range of the training data. For these tests the predictions of the NNCM were poor. NNCM was then re-trained using additional data generated for the four stress paths described earlier and for confining pressures of 1 kPa and 200 kPa. The graphs of q versus ε_yy are presented in the figures 5 and 6 for the HSM, trained and re-trained NNCMs. This confirms that extrapolation by NNCM is always of poorer quality than interpolation.

Figure 5. Graph of q versus εyy under LC conditions for a confining pressure of 10 kPa

Figure 6. Graph of q versus εyy under LC conditions for a confining pressure of 200 kPa

4. Two-Surface Model in Multilaminate Framework (TDHM)

4.1. Background

A multilaminate framework for modelling the behaviour of soils was presented by Pande and Sharma as in [12] almost two decades ago. Yield and failure criterion for most elasto-plasticity based models are written terms of stress invariants. This type of formulation inhibits the development of plastic strains purely due to rotation of principal stresses. Such formulations also exclude plastic flow induced anisotropy. Multilaminate formulation of plasticity based models overcomes both the above mentioned short comings. Here, the yield and failure criteria as well as plastic potential and hardening softening functions are written in shear stress/ normal stress space for randomly oriented micro-planes. A framework similar to multilaminate framework was presented by Bazant as in [13] under the name ‘micro-plane model’. Extensive has been carried out on practical applications of the multilaminate framework by Scharinger & Schweiger, H. as in [14]. In the following a brief description of the formulation of a new two-surface model in multilaminate framework which is formulated for cyclic and transient loading, is presented as in [6]. Some details of multilaminate formulation are given for completeness and continuity. It is noted that the objective here is to develop an NNCM equivalent of this complex model.

4.2. Mathematical Formulation

Let us adopt a set of local co-ordinate axes (n, s, t) for each sampling plane. The n-axis is normal to the sampling plane whilst axes s and t are arbitrarily chosen on the sampling plane as shown in Figure 7 below.

Figure 7. Definition of a local system of axes on a typical sampling plane

In the local system of co-ordinates, the normal and shear stress components can be defined as follows:

(6)

(7)

(8)

(9)

Where

(10)

If the yield and bounding surfaces for frictional materials are plotted in resultant τ – σ_n’ (all stresses are effective and prime on the symbols will be omitted henceforth) space, they are represented by a pair of straight lines. The angle between the pair of lines representing the elastic domain is assumed as arbitrarily small while the angle between the lines representing bounding surface is related to the friction peak angle φ’. However, when two orthogonal components of shear stress on a sampling plane are considered, the yield and bounding surfaces are represented by two cones in τ_s – τ_t – σ_n space, the smaller one gyrating inside the larger one, as it is seen Fig 8.

Figure 8. Yield and bounding surfaces for a sampling plane in σ′n - τs - τt space

4.3. Equations of the Bounding, the Yield and the Plastic Potential Surfaces

The equation of the Bounding Surface is postulated below:

(11)

Where:

(12)

The equation of the Yield Surface:

(13)

Where:

(14)

The Plastic potential function is given by:

(15)

Where:

are the sizes of the Bounding, the Yield surface and the Plastic Potential function respectively

5. Training an Equivalent NNCM of TDHM

TDHM is used in the situations of transient and dynamic loading. In this case, the computation of stress increment from stain increment is generally computationally intensive, especially for a complex model such as TDHM. In order to train the NNCM, strain-controlled triaxial test data conforming to seven different specified constant volumetric strains were generated using TDHM. The results of (TDHM) have been used to check the adequacy of the training of the equivalent NNCM. Additionally four new analyses which were not included in the training data have been used to compare and validate the NNCM equivalent to TDHM. The parameters of the model are given: Rf=0.52, R1=0.02, Rc=0.43, Coh=0, ge=40, β=70, γ =3.15, Acon=0.012. These parameters are described by Lee as in [15].

5.1. Synthetic Data Generation

Synthetic data of strain and stress increments in 7 hypothetical two-way strain-controlled cyclic triaxial configuration (dε_xx=dε_zz), have been generated using a point integration program of TDHM (DRIVER), described in the previous section. Various ratios of axial to radial strains have been chosen using the following strain paths trajectories:

Table 2. Strain paths trajectories

The NNCM was trained using the above data and the architecture shown in Fig. 1 and was used to predict the strain-stress response paths. These were very accurate. Figs.9-15 show stress paths in q-p^’ space together with stress paths obtained from the original DRIVER program of TDHM.

Figure 9. Stress path in q-p’ space for the test data SSTD1

Figure 10. Stress path in q-p’ space for the test data SSTD2

Figure 11. Stress path in q-p’ space for the test data SSTD3

Figure 12. Stress path in q-p’ space for the test data SSTD4

Figure 13. Stress path in q-p’ space for the test data SSTD5

Figure 14. Stress path in q-p’ space for the test data SSTD6

Figure 15. Stress path in q-p’ space for the test data SSTD7

5.2. Validation of the Trained NNCM

The following additional strain-stress data were generated using the DRIVER but were not used in training, were subsequently used for validation:

Table 3. Strain paths trajectories used for validation

The results of the equivalent NNCM prediction are compared with the validation data in Figs. 16 & 17, where stress paths are plotted in q-p^’ space. A very good agreement is observed.

Figure 16. Stress path in q-p’ space for the test data SSVD3

Figure 17. Stress path in q-p’ space for the test data SSVD4

6. Conclusions

In the authors’ opinion, the potential of applications of Artificial Neural Networks in engineering analysis and design is immense and is yet to be exploited. Appropriate constitutive models of materials are the key to the successful prediction of the behaviour of engineering structures. A large number of models based on various constitutive theories have been proposed in the last three decades for geomaterials, which show a large variation in their properties. All the models proposed assume a-priori, a mathematical framework for the model and the material parameters corresponding to the assumed framework have to be identified from physical material tests. Many material parameters in complex constitutive theories have no physical meaning, are difficult to determine and have to be identified by trial and error from numerical simulations. In spite of this, many features of soil behaviour such as stiffness at small strains, higher stiffness on reversal of stress path, influence of rotation of principal stress axes etc. have not been captured in a single model. It is expected that more complex models will be developed in the future.

In this paper a methodology for converting or recasting complex constitutive models for geomaterials developed in any mathematical framework into a fully trained neural network equivalent is proposed. This Neural Network based Constitutive Model can then be embedded in an appropriate finite element solution code. The length of strain trajectory traced by a material point, also called ‘intrinsic time’ is used as an additional input parameter in training. This is essential for situations of cyclic and transient loading. For the purpose of illustration, two constitutive models viz. Hardening Soil Model (HSM) available in the commercial software, PLAXIS and a Two-surface Deviatoric Hardening Model in the multilaminate framework (TDHM) developed by Lee and Pande (2004) have been cast in the form of an ANN. It is observed that equivalents for both models can be easily trained and can produce accurate results in all situations including a large number cycles.

It is perceived that real life problems in future will have to be solved with increasingly more complex constitutive models for geomaterials. This will lead to unacceptable computational processing times. The authors, therefore, provide a solution since response of a trained NNCM equivalent is instantaneous. Computational efficiency will be achieved even for simpler models.