Table of Contents Author Guidelines Submit a Manuscript
Mathematical Problems in Engineering
Volume 2014, Article ID 587412, 14 pages
Research Article

An Elastic-Viscoplastic Model for Time-Dependent Behavior of Soft Clay and Its Application on Rheological Consolidation

1Ningbo Institute of Technology, Zhejiang University, Ningbo 315100, China
2Research Center of Coastal and Urban Geotechnical Engineering, Zhejiang University, Hangzhou 310058, China
3Zhejiang Provincial Institute of Communications Planning, Design & Research, Hangzhou 310006, China

Received 10 February 2014; Accepted 12 May 2014; Published 5 June 2014

Academic Editor: Gradimir Milovanović

Copyright © 2014 Jinzhu Li et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.


To describe the time-dependent behavior of soft clay, this paper extended one-dimensional Nishihara model to three-dimensional stress state based on the framework of Perzyna’s overstress theory and modified cam-clay model. The yield criterion of modified cam-clay model was used to describe the plastic properties of soft clay, and the overstress theory was used to describe the strain rate effect. Triaxial rheological tests were carried out on Ningbo soft clay and the rheological characteristics were studied. Based on laboratory results, the parameters of proposed model were determined by curve fitting, which show that this model is suitable for the rheological characteristics of Ningbo soft clay. The analysis of parameters shows that, the value of parameters changes slightly with different deviatoric stress when the confining pressure was constant, but changes notably with the increase of confining pressure. A user material subroutine of the proposed constitutive mode was coded on the platform of the FEM software ABAQUS and verified by triaxial compression of soil column. A plain strain problem was computed to analyze the rheological consolidation properties of soft clay, in which the rheological effect and the finite strain effect were considered.

1. Introduction

According to the classic consolidation theory presented by Terzaghi, when soft clay is consolidated under external loads, the excessive pore pressure dissipates and the effective stress increases, while the deformation grows until the excessive pore pressure dissipates completely. This theory has laid the foundation of soft soil mechanics and even became the main theoretical basis of the computation of soft soil foundation. However, in the former works, Buisman first discovered the long term deformation of soils after the excessive pore pressure was dissipated in experiments, which is also reported in the works of Zeevaert [1], Leonards and Ramiah [2], and also in lots of engineering practices. This long term deformation which cannot be explained by classic theories is usually called rheology, which contains creep, stress relaxation, long-term strength, elastic aftereffect, hysteresis effect, and so on. Because of strict requirement of post-construction settlement, researchers are paying increasing attention on the rheology of soft clay. To propose a suitable constitutive model to describe the rheology of soft clay and to determine the model parameters accurately by laboratory tests, these are the key points of the topic. Since 1960s, a number of constitutive models were proposed for time-dependent behavior of soft clay. Commonly used rheological constitutive models can be divided into three categories: empirical models, component models, and microscopic models, in which the first two categories were most widely used. Empirical models such as Singh-Mitchell model [3], Mesri et al. model [4, 5] are mostly one dimensional models, which were often used to calculate the rheological deformation in one dimensional problem but were limited in complex stress state. Component models are intuitive in physical meaning, which can describe numerous rheological models by different combinations of components and can be easily extended to multidimensional problems to apply in numeral computation for complex stress state, such as the work of Adachi and Oka [6], Fodil et al. [7], Yin et al. [8, 9], Hinchberger and Rowe [10], Hinchberger and Qu [11]. Based on the framework of Perzyna’s overstress theory [12] and modified cam-clay model, this paper proposed a three-dimensional rheological model. The model was verified by laboratory triaxial rheological tests of Ningbo soft clay, and the model parameters were determined by curve fitting. On the platform of ABAQUS, a user material subroutine of the model was coded and verified.

2. Model for Soft Clay

As one of commonly used component models, Nishihara model was connected by Hooke body, Kelvin body, and Bingham body, which can describe elastic deformation, creep, elastic aftereffect, and viscous flow comprehensively, such as Figure 1. When a constant external pressure was given, the stress and strain of the model can be written as follows:

Figure 1: Sketch of Nishihara model.

In the above equation, , , denote the stress of Hooke body, Kelvin body, and Bingham body, respectively; , , , denote the total strain, the strain of Hooke body, Kelvin body, and Bingham body, respectively.

Substituting the stress-strain relationship of Kelvin body and Bingham body into (1a) and (1b), the constitutive equation of one-dimensional Nishihara model can be obtained as follows: where denotes the elastic modulus of Hooke body, and denote the elastic modulus and viscosity coefficient of Kelvin body, and denote the viscosity coefficient and yield stress of Bingham body, denotes a switch function, which can be defined as follows:

Equation (2) contains the first two stages of creep, which are often called primary creep stage and steady creep stage. In fact, for complex stress state, (2) cannot be used to calculate the deformation of soil, it is necessary to extend (2) to three dimensional, so some hypotheses can be presented as follows:(1)the soil material is isotropic;(2)the volume deformation is just caused by spherical stress, and it is unrelated to shearing stress;(3)the Poisson’s ratio of soil is constant and does not change with the stress or time.

Based on the above hypotheses, the stress-strain relationship of Hooke body is given as where is the first stress invariant and and denote the bulk modulus and shear modulus of Hooke body.

For Kelvin body, the stress-strain relationship is given as where and denote the shear modulus and 3D viscosity coefficient of Kelvin body, in which , is the Poisson’s ratio of soil.

For Bingham body, the stress-strain relationship is given as

In (6), denotes 3D viscosity coefficient of Bingham body, and is a switch function to judge whether plastic yield is occurring and whether the amplitude of the plastic yield occurred, is the plastic potential function, when the associated flow rule is adopted, , where is the yield function. can be defined as follows: where is the current yield function, is the initial yield function, and is a constant derived from experiments; according to the work of Wei [13], the value of can be assigned to 1.0 approximately for soft clay.

According to the work of Adachi and Oka [6] and Yin et al. [8, 9], we define a static yield criterion , which represents a reference yield surface for the material, such as Figure 2. Its initial shape depends on the consolidation pressure . The expansion of the static yield surface, which describes the hardening of the material, is expressed by the variation of the consolidation pressure due to the inelastic volumetric strain as follows:

Figure 2: Yield criterion of model.

In (8), is the intercept of initial static yield surface in the axis, for soft clay, it can be assigned to the biggest consolidation stress in history; is the slope of normal consolidation curve in plane; is the slope of recompression curve in plane.

A dynamic yield criterion is defined to describe the current state of stress as follows:

From (8) and (9), the and in (7) can be written as

So, we can write (7) as follows:

Combining (6)~(11), (6) can be written as follows:

From (4), (5), and (12), when the stress is constant, the total strain of the presented model can be written as follows:

3. Experiments and Computation

3.1. Test Program

In order to verify the presented model, laboratory triaxial rheological tests were performed to observe the mechanical behavior of soft clay under long-term loading. The soil investigated in the experiments is Ningbo soft clay, which is a kind of problematic soil for low strength, high compressibility, and time-dependent behavior and deposits in Hangzhou Bay, China. The basic properties of the soil were summarized in Table 1. Both moisture content and moist unit weight of this material are significantly less than those of typical natural sedimentary deposits.

Table 1: Property of testing material.

The experimental equipment was refitted on the platform of strain controlled triaxial apparatus, the original strain loading method was changed to stress loading method, while the confining pressure system, the back pressure loading system, and the measurement system were reserved. The soil was cut into replicate specimens with a diameter of 39.8 mm and a height of 80.0 mm and placed in the triaxial test apparatus. Both soil specimens were consolidated 24 hours under a confining pressure of 100 kPa and 200 kPa, respectively; then, keep the confining pressure as constant and apply the deviatoric stress increment until the specimens were damaged or the total strain exceeds 15%. In the whole process, the free drainage conditions were kept. Table 2 showed the different loading rates applied in the experiments. In order to make sure that the creep deformation can be fully developed, each load increment lasted no less than 7 days until the deformation of the specimen is lower than 0.01 mm within 24 h.

Table 2: Loading scheme of triaxial rheological tests.
3.2. Experimental Results

The experimental results were presented to observe the rheological behavior of the soft clay. It is considered that deformation occurs during the whole process, but only the deformation that occurs after the dissipation of the excessive pore water pressure can be regarded as creep only. The full shearing process stress-strain-time curve of each specimen is illustrated in Figure 3, and the mutistage strain-time curves of each specimen are illustrated in Figure 4. The figures show that the creep process under triaxial loading exhibits attenuation characteristics when the deviatoric stress is small, which contains two stages of creep: primary stage and steady stage; while when the deviatoric stress is large, it may exhibit accelerated characteristics, which contains three stages of creep: primary stage, steady stage and accelerated stage, such as the last load increment of specimen 2.

Figure 3: Full-process press-strain-time curve.
Figure 4: Mutistage strain-time curves.
3.3. Computation and Curve Fitting

For triaxial rheological tests, , , so (13) can be written as follows: where, , , , .

Use (14) to fit the experimental data, the results were showed in Figures 5 and 6, in which the parameters SSE and denote the residual sum of squares and squares of correlation coefficient, respectively. From Figures 5 and 6, we can see that the model fits testing data of primary creep stage and steady creep stage very well, while for the accelerated creep stage, it needs further study. According to the value of , , , and , the value of the model parameters , , , and can be determined easily, such as in Tables 3 and 4.

Table 3: Model parameters of specimen 1 ( kPa).
Table 4: Model parameters of specimen 2 ( kPa).
Figure 5: Fitting curves of specimen 1 ( = 100 kPa).
Figure 6: Fitting curves of specimen 2 ( = 200 kPa).
3.4. Discussion

Table 3 shows that the model parameters of specimen 1 ( = 100 kPa) change slightly with the change of deviatoric stress except the first stage; Table 4 shows that the model parameters of specimen 2 ( = 200 kPa) change slightly with the change of deviatoric stress except . Compared with the model parameters of specimen 1 and specimen 2, it can be found that both and increase notably with the increase of , which means that the elastic strain reduces when increases, and increases notably while is reduced notably with the increase of , which means that the viscoelastic strain rate and viscoplastic strain rate change significantly when increases. For the lack of more experiment data, the further variation of model parameters was unable to find out, which needed further study.

Based on the parameters in Tables 3 and 4, the components of elastic strain, viscoelastic strain, and viscoplastic strain at any time can be calculated with (13). Tables 5 and 6 show , , and of all loading levels when = 7 d; it can be found that the elastic strain accounted for the major part of the total strain, but the viscoelastic strain and viscoplastic strain increase with the increase of , which takes about 15%~30% for each loading increment of the total and cannot be neglected.

Table 5: Percentage of strain components of specimen 1 ( kPa,  d).
Table 6: Percentage of strain components of specimen 2 ( kPa,  d).

Figure 7 shows several typical curves of the development of computed viscoelastic strain and viscoplastic strain of specimen 1 and specimen 2. It can be found that the viscoelastic strain increases quickly in the early stage and then reaches steady in a short time (about 30 hours), while the viscoplastic strain increases the whole process, which indicates that the long-term deformation after consolidation of soft clay foundation is mainly caused by the viscoplasticity rather than viscoelasticity

Figure 7: Development of computed viscous strain (VES: viscoelastic strain; VPS: viscoplastic strain).

4. Derivation and Verification of UMAT

4.1. Derivation of UMAT

ABAQUS is a FEM software which has powerful computing capabilities of strongly nonlinear problems. A user material subroutine (UMAT) is coded for the proposed model on the platform of ABAQUS to study the elastic-viscoplastic consolidation properties of soft clay. According to ABAQUS, the stiffness coefficient matrix (Jacobian matrix) and the equation to update the stress and strain should be defined in the UAMT. The Jacobian matrix is defined as follows: where and denote the stress increment matrix and strain increment matrix, respectively.

For the proposed model, (4)~(6) can be written as follows in 2D or 3D problems:

In (16a), (16b), and (16c) , , , and denote the Voigt matrix form of stress and strain component tensors; denotes the flexibility matrix, for plane strain problem and 3D problems, and it can be written as follows:

Suppose that the time step of the th increment is , so the viscoelastic strain increment in the increment step can be written as follows: where is differential coefficient, . Suppose that , ; (18) can be written as follows:

The viscoplastic strain increment in the increment step can be written as follows:

In order to calculate the viscoplastic strain rate at the end of the th increment, apply the Taylor series expansion and ignore the second order trace; then we can obtain the following equation:

Suppose that ; (21) can be written as

Substituting (22) into (21), we can obtain the following equation:

From (16c), can be written as where . According to (9), can be written as follows: where , , , .

From (25), the following equation can be derived: where , , . , , and are matrices written as follows:

Combining the above equations, the value of can be derived. Substituting into (22), the viscoplastic strain increment is obtained.

For current increment, is calculated by (28) as follows:

Substituting (19) and (23) into (28), we can obtain the following equation:

Equation (29) is the relationship between stress increment and strain increment, in which matrix is defined as follows:

Based on the above equations, a user material subroutine was coded for the proposed constitutive model on the platform of ABAQUS.

4.2. Verification of UMAT

In order to verify the validity of the UMAT, we study a soil column of triaxial compression such as Figure 8, the size of the column is 0.1 m × 0.1 m × 0.2 m. Supposing that = 100 kPa,  kPa, and the model parameters are given in Table 7. The column was cut into one element, the element type is C3D8, and the consolidation effect was not considered here. Embedding the UMAT coded in this paper on the platform of ABAQUS, we calculated the triaxial compression property of the column; the results were summarized in Figure 9~Figure 10.

Table 7: Material parameters of verification model.
Figure 8: Sketch of verification model.
Figure 9: Triaxial compression property of soil column.
Figure 10: Sensitivity of model parameters.

From Figure 9(a), we can see the deformation of the column exhibits apparent rheology, the rheological deformation takes about 37% of the whole deformation. Figure 9(b) shows that the static yield surface gradually expanded with the time until it reached steady state, when it coincided with the dynamic yield surface. From Figures 9(c) and 9(d), we can see that the viscoelastic strain and the viscoplastic strain both gradually increased with time until they reached their steady state; accordingly, the viscoelastic strain rate and the viscoplastic strain rate were gradually deduced to nearly zero. The above results fit the property of Nishihara model and the hardening law of modified Cambridge model very well, which proved the correctness of the UMAT.

In the proposed model, , , and are the most important parameters; Figure 10 showed the sensitivity of these parameters. From Figures 10(a) and 10(b), we can see the final viscous deformation is a constant regardless of the different values of , , but the time when it reaches steady varies greatly, so, the accuracy of and is directly related to the accuracy of computed displacement in the rheological process. Figure 10(c) shows that the value of has a great influence on the viscoplastic strain: the smaller is, the more rapidly the viscoplastic strain develops and the larger the final viscoplastic strain occurs, and so, the accuracy of is not only directly related to the accuracy of computed displacement in the rheological process, but also affected the final displacement, which needs to be paid more attention.

5. Analysis on Rheological Consolidation of Soft Clay

5.1. Computed Model

Consolidation problems can be classified as small strain theory and finite strain theory according to the difference of geometric assumption. Small strain theory assumed that the strain is little, but the consolidation deformation is usually very large for deep soft clay foundation. Take dredger fill for example, the deformation may exceed 80%, for which the small strain theory is no longer suitable [14, 15]. As the pioneer researchers studied on finite strain consolidation theory, Mikasa [16] and Gibson et al. [17] derived their one-dimensional finite strain consolidation equations, respectively, in 1960s, and many scholars developed further research on the topic in the following decades. Because of strong nonlinearity, most of existing researches on finite strain consolidation adopted some simplified assumptions without considering rheological characteristics. By means of ABAQUS, this paper studied the rheological consolidation of soft clay considering finite strain effect.

The computed model is plane strain problem as shown in Figure 11. The soil foundation is homogeneous, anisotropic, saturated, and normally consolidated soft clay with a thickness of 10.0 m, and only the top boundary is permeable. In addition, the semilogarithmic empirical formula proposed by Tavenas et al. [18] was adopted to consider the variation of permeability in the process of consolidation as follows: where is the initial permeability coefficients, and is the permeability coefficient in the consolidation process; is the initial void ratio, and is the void ratio in the consolidation process; is the permeability index. According to the above work, the model parameters are summarized in Table 8 approximately.

Table 8: Material parameters of computed model.
Figure 11: Calculation model and mesh generation.
5.2. Consolidation Property of Soft Clay

From Figure 12 we can see the differences of the consolidation process between finite strain and small strain. Figure 12(a) shows that for any location the excess pore pressure of finite strain is less than that of small strain at the same time, and the difference is increasing with the increase of time in beginning period. Figure 12(b) shows that the displacement of finite strain is less than that of small strain and the difference increased significantly if the load is increasing. For example, when = 200 kPa, the difference may take about 20% while it is nearly 0 when = 50 kPa. So, it is essential to consider the error caused by small strain if the load is large and the compressibility of the soil is large.

Figure 12: Results comparison of large strain and small strain (SS: small train; FS: finite strain).

From Figure 13 we can see the influence of rheology on the consolidation process. Figure 13(a) shows that the dissipation of excess pore pressure is significantly slower if the rheological effect is considered; for example, when = 1000 d, the excess pore pressure is nearly 0 if the rheological effect is not considered, but it is still larger than 20 kPa if the rheological effect is considered. Figure 13(a) shows that the displacement is influenced by rheological effect too; for example, in the early stage ( ≤ 100 d), the displacement is nearly the same for two cases, but the final displacement may have a difference of more than 35%, which indicates that the displacement in the early stage is mainly caused by consolidation, while the displacement in the late stage is mainly caused by rheology.

Figure 13: Effects of rheology on calculated results. (NR: rheological effect not considered; R: rheological effect considered).

6. Conclusions

(1)An elastic-viscoplastic model developed to describe the time-dependent behavior of normally consolidated soft clay was presented in the paper on the framework of modified cam-clay elastic plastic model and Perzyna’s overstress viscoplastic theory.(2)A series of laboratory triaxial rheological tests of Ningbo soft clay with different confining pressure and deviatoric stress are performed; the results show that the presented constitutive model was suitable for the rheological characteristic of Ningbo soft clay, and it is achievable to determine the parameters of presented model by curve fitting.(3)The analysis on the model parameters shows that the value of all parameters is related to confining pressure and deviatoric stress. It may need further study on the relationship between parameters and the stress level so as to use the model conveniently.(4)In the proposed model, the viscoelastic viscosity coefficient , the viscoplastic viscosity coefficient , and the initial yield surface are the most important parameters. and have important influence on the time of viscous strain developed, and have important influence on both the viscoplastic strain rate and the final viscoplastic strain.(5)The existence of rheology will significantly slow down the dissipation of the excess pore pressure. Considering the rheological effect, the total displacement is larger, and the displacement in the early stage is mainly caused by consolidation, while in the late stage it is mainly caused by rheology.(6)For consolidation problems of high compressibility soft clay foundation, great error will be caused by small strain assumption, so finite strain effect should be considered. When finite strain effect is considered, the consolidation process is developed more rapidly, and the computed consolidation displacement will be less than that of small strain.

Conflict of Interests

The authors do not have any conflict of interests regarding the contents of the paper.


This research was supported by the National Natural Science Foundation of China (Grant no. 51378469) and the Ningbo Municipal Natural Science Foundation (Grant no. 2013A610196). The authors wish to express their gratitude to the above for financial support. (1) Project (51378469) is supported by the National Natural Science Foundation of China; (2) Project (2013A610196) is supported by the Natural Science Foundation of Ningbo City, China.


  1. L. Zeevaert, Consolidation of Mexico City Volcanic Clay Conference on Soils For Engineering Purpose, ASTM, Philadelphia, Pa, USA, 1958.
  2. G. A. Leonards and B. K. Ramiah, Time Effects in the Consolidation of Clays Symposium on Time Rate of Loading in Testing Soils, ASTM, Philadelphia, Pa, USA, 1960.
  3. A. Singh and J. K. Mitchell, “General stress-strain-time function for clay,” Journal of the Clay Mechanics and Foundation Division, ASCE, vol. 94, no. 1, pp. 21–46, 1968. View at Google Scholar
  4. G. Mesri and P. M. Godlewski, “Time and stress-compressibility interrelationship,” Journal of Geotechical Engineering Division, ASCE, vol. 103, no. 5, pp. 417–430, 1977. View at Google Scholar · View at Scopus
  5. G. Mesri, E. Febres-Cordero, D. R. Shields, and A. Castro, “Shear stress strain time behaviour of clays,” Geotechnique, vol. 31, no. 4, pp. 537–552, 1981. View at Google Scholar · View at Scopus
  6. T. Adachi and F. Oka, “Constitutive equations for normally consolidation clay based on elastic-viscoplastic,” Soils and Foundations, vol. 22, no. 4, pp. 57–70, 1982. View at Google Scholar · View at Scopus
  7. A. Fodil, W. Aloulou, and P. Y. Hicher, “Viscoplastic behaviour of soft clay,” Geotechnique, vol. 47, no. 3, pp. 581–591, 1997. View at Google Scholar · View at Scopus
  8. Z. Y. Yin, P. Y. Hicher, Y. Riou, and H. W. Huang, “An elasto-viscoplastic model for soft clay,” in Proceedings of the Geoshanghai Conference of Soil and Rock Behavior and Modeling, pp. 312–319, Geotechnical Special Publication, Shanghai, China, June 2006. View at Publisher · View at Google Scholar · View at Scopus
  9. Z. Y. Yin and P. Y. Hicher, “Identifying parameters controlling soil delayed behaviour from laboratory and in situ pressure meter testing,” International Journal for Numerical and Analytical Methods in Geomechanics, vol. 32, no. 12, pp. 1515–1535, 2008. View at Publisher · View at Google Scholar
  10. S. D. Hinchberger and R. K. Rowe, “Evaluation of the predictive ability of two elastic-viscoplastic constitutive models,” Canadian Geotechnical Journal, vol. 42, no. 6, pp. 1675–1694, 2005. View at Publisher · View at Google Scholar · View at Scopus
  11. S. D. Hinchberger and G. Qu, “Viscoplastic constitutive approach for ratesensitive structured clays,” Canadian Geotechnical Journal, vol. 46, no. 6, pp. 609–626, 2009. View at Publisher · View at Google Scholar · View at Scopus
  12. P. Perzyna, “Fundamental problems in viscoplasticity,” Advances in Applied Mechanics, vol. 9, pp. 243–377, 1966. View at Publisher · View at Google Scholar · View at Scopus
  13. L. M. Wei, Double non-linear fluid-solid coupled simulating analysis and settlement prediction for soft soil subgrade [Ph.D. thesis], Central South University, Changsha, China. View at Publisher · View at Google Scholar
  14. W. G. Weber, “Performance of embankments constructed over peat,” Journal of Soil Mechanics and Foundations, ASCE, vol. 95, no. 1, pp. 53–76, 1969. View at Google Scholar
  15. K. W. Cargill, “Prediction of consolidation of very soft soil,” Journal of Geotechnical Engineering, ASCE, vol. 110, no. 6, pp. 775–795, 1984. View at Google Scholar · View at Scopus
  16. M. Mikasa, “The consolidation of soft clay—a new consolidation theory and its application,” Civil Engineering in Japana, JSCE, vol. 1, no. 1, pp. 21–26, 1965. View at Google Scholar
  17. R. E. Gibson, G. L. England, and M. J. L. Hussey, “The theory of one-dimensional consolidation of saturated clays II: finite nonlinear consolidation of thick homogeneous layers,” Geotechnique, vol. 17, no. 2, pp. 261–237, 1967. View at Google Scholar
  18. F. Tavenas, P. Jean, P. Leblond, and S. Leroueil, “The permeability of natural soft clays. Part II: permeability characteristics,” Canadian Geotechnical Journal, vol. 20, no. 4, pp. 645–660, 1983. View at Google Scholar · View at Scopus