A generalized backward Euler scheme for the integration of a mixed isotropic-kinematic hardening model for clays

The mechanical response of natural clays is characterised by highly non-linear behaviour, memory of the past strain-history, evolving anisotropy, non-coaxiality and, when cemented, mechanically induced bond degradation phenomena. In recent years a number of constitutive models have been proposed to mathematically describe these features, often being characterised by complex formulations leading to non-trivial problems in their numerical integration. On the other hand, accuracy and stability are recognised as crucial requirements in the development of any integration algorithm for realistic material models, in order to ensure the necessary computational correctness and efficiency in their use within Finite Element codes. The paper describes a fully implicit stress-point algorithm for the numerical integration of a single surface mixed isotropic-kinematic hardening plasticity model for bonded clays. The model is characterised by an ellipsoid-shaped yield function, inside which a stress dependent reversible stiffness is accounted for by a hyperelastic formulation. The isotropic part of the hardening law extends the standard Cam-Clay one to include plastic strain-driven softening due to bond degradation, while the kinematic hardening part controls the evolution of the position of the yield surface in the stress space. This latter hardening term is formulated in such a way that along radial stress paths the centre of the yield surface initially moves to then achieve a stabilised position, corresponding to the imposed direction of the path. The Generalised Backward Euler algorithm is formulated in the space of elastic strain and internal variables, leading to a system of 14 unknowns and corresponding non-linear equations. The solution is obtained iteratively by means of the Newton's method. The proposed algorithm allows the consistent linearization of the constitutive equations guaranteeing the quadratic rate of asymptotic convergence in the global level Newton-Raphson iterative procedure. The accuracy and the convergence properties of the proposed algorithm are evaluated by numerical simulations of single element tests and represented by iso-error maps. A discussion on the theoretical implications of the kinematic hardening assumptions on the anisotropic and non-coaxial response of the model is also summarised in the final part of the paper.