Computational simulation of damage accumulation processes in cracked bodies by the UMAT procedure of SIMULIA Abaqus

The paper presents the experience of using the user subroutine UMAT for finite element package SIMULIA Abaqus/CAE for damage accumulation processes in the vicinity of the crack. A continuum damage mechanics model based on the constitutive relations of linear elastic isotropic materials with the incorporated damage tensor components is used to describe the material behavior. The material nonlinearity arising from the deformation process is modeled by introducing an anisotropic damage tensor of the second rank into the constitutive equation. The material model is described by means of user procedure UMAT of SIMULIA Abaqus. The finite element (FE) mechanical constitutive model is implemented in Abaqus/Standard via a UMAT routine. Numerical experiments for a large series of cracked specimens have been performed. Computed stress and damage tensor components were found. It is shown that they are not dependent on the FE mesh refinement. Distributions of the damage tensor components in the vicinity of the crack tip in cracked specimens of different configurations under mixed mode loading in a wide range of mixed mode loadings are found. The configurations of active damage accumulation process zone in the cracked specimens are obtained. It is shown that the damage accumulation process has substantial influence on the stress-strain state in the vicinity of the crack tip and leads to decrease of the stress concentration in cracked specimens.


Introduction
The accurate numerical modeling of fracture is a constantly pursued area in solid mechanics [1]. Damage initiation, accumulation and propagation processes in engineering materials are of particular importance for the design, production, certification, and monitoring of an increasingly large variety of structures [1]. Many techniques for fracture modeling have been developed over the time, which are broadly classified into two categories, namely discrete (finite element method (FEM), extended finite element method (XFEM), meshfree methods) and smeared approaches (phase field and gradient damage). Among discrete methods, the most popular one is the finite element method. Further developments in the discrete approaches include extended finite element method. Another class of models with continuum representation of crack are the gradient damage models which trace their roots to the continuum damage mechanics (CDM). Continuum damage mechanics models in provide a constitutive framework for representing damage induced stiffness reduction through effective damage parameters that represent overall material degradation. Both phenomenological and micromechanical CDM models, have been proposed in the literature for modeling composite failure. Phenomenological models in employ scalar, second order and fourth order damage tensors using mathematically and thermodynamically consistent formulations of damage mechanics. Damage parameters are identified through macroscopic experimental observations and data. A number of formulations based on continuum damage mechanics have been proposed for various applications in the literature [2][3][4][5]. Some of early contributions in this area were collected in three books by Kachanov [3], Lemaitre [4], as well as Lemaitre and Chaboche [5]. Stemming from the works of Lemaitre and Chaboche [5], the implementation of CDM using FEM for describing crack growth led to mesh dependent solutions. The dependency existed in terms of damage localization to a very thin band of vanishing volume and the direction of crack growth. Many remedies to overcome this pathological mesh dependence were proposed and discussed in Murakami and Liu [6]. Material libraries in FEM solvers need to be extended with modern, enhanced predictive fidelity models. As a result, the FEM user often undertakes the design of user-defined material subroutines (UMAT) to achieve this objective. They are optional interfaces that improve the usefulness of the FEM process in both academic and industrial sectors. The material response of a defined virtual domain within FEM is defined using carefully formulated constitutive models [7]. It is important to gain this theoretical, as well as practical, understanding in order for the user to interpret FEM outputs [7]. This is very important as not all known material behaviours have been modelled and validated as in-built material models within common FEM solvers [7]. Where new materials are discovered, or even modifications and improvements desired for existing in-built material models, the FEM user usually resorts to describing computationally their versions of the constitutive mathematics to reflect such changes. The UMAT subroutine within FEM solvers helps the researcher to do so. Material response is considered an important essential principle for any FEM framework. In fact, the material modules of FEM solvers continue to present the most challenges to FEM users. There is a need for the FEM user to understand how material models work and how they can also be created. The predictive fidelity of an FEM solution relies strongly on creating the most robust material model that captures experimentally observed mechanical behaviour. This starts with formulating the constitutive mathematics that define the experimentally observed behaviour before deploying the formulations into a numerical scheme, herein called the material model, which is then used within the FEM process to undertake simulations. Constitutive models are buzz words in computational mechanics that represent a numerical material model [7]. The keyword "constitutive" requires that the numerical models must capture the whole gamut of mechanical behaviour that constitutes known experimentally observed responses of the chosen material. It is therefore, a laudable objective for the constitutive modeller to create a model that represents comprehensively the behaviour of the test material. Although most material models do not necessarily meet this requirement, it remains the goal amongst model developers to continue to develop models that are truly 'constitutive.' All constitutive models are developed based on constitutive mathematics, which are a set of analytical formulations that describe the mechanics of the chosen material. UMAT is used in the ABAQUS/Standard analysis interface. This is a simulation engine designed specifically for principally implicit analysis and hence the included user-defined material subroutine is a UMAT. It is strictly used for describing material responses where the analysis procedure comprises of a mechanical behaviour.
The article is aimed at modeling the damage accumulation process near the stress concentrators and crack tips using the damage second rank tensor via UMAT procedure of SIMULIA Abaqus.

Constitutive equations
Initially, in order to understand the creep damage evolution in the material Kachanov [3] and Rabotnov [8] chose to represent continuum damage as an effective loss in material cross-section due to internal voids, microcracks and cavities. According to this model the creep strain rates and the damage evolution law have the forms 1 3 , Eqs.
(1) are generalized to multi-axial stress states using J2-flow theory as follows [9] ( ) where and are the creep strain, the effective stress and time respectively. The material constant are the Norton law constants, which control secondary creep behaviour and they can be determined from a log-log of the creep strain rate vs the applied stress while the constant m controls the primary creep strain. The damage under a multiaxial stress state can be ascertained experimentally by identifying those combinations of stress which give the same rupture time [10]. These may be presented by surfaces in stress space which are called isochronous surfaces. Mathematically the isochronous rupture surfaces are presented by a scalar invariant of stress σ is the effective stress, kk σ is hydrostatic stress: α β are constants determined from the experimental surfaces. Since the pioneering works of Kachanov and Rabotnov continuum damage mechanics was extensively used to describe the stiffness degradation of a material undergoing damage [11] both for brittle materials such as concrete and rocks and ductile materials such as metals and alloys at high temperatures. Existing isotropic damage models were modified to take into account the damage-induced anisotropy in materials. A simple anisotropic damage model was developed in [12] to investigate the damage caused by an indenter on a soda-lime glass. The model predicted the cracking damage pattern and zone size adequately and was in good agreement with experimental results [12,13]. The stiffness degradation arising from the damage accumulation was simulated by adding a damage tensor into the glass constitutive equation [12] ( ) where ij σ and ij ε are the stress and strain tensor components, respectively, and where ij D are the damage parameters and ( ), 1, 2 i C T i = are two constants of the material. It is supposed that the terms of the damage tensor are function of the stress state and take values between zero and one meaning that the material is virgin and one meaning that the material is fully damaged: ≤ Though the original damage model [12,13] was meant to take into account the effect of the shear stresses on damage evolution, only the diagonal components of the damage tensor were accounted for in this study. Their values are assumed to follow a linear evolution law such that: The values of 1 C and 2 C are calculated [11] by stating that the stresses 11 σ and 22 σ drop to zero 33 Equations (7) and (8) must be equivalent. It allows to find 1 2 C µ = and 2 3 C µ = − . This choice of the constants is different from the original paper [12] where 1 , so that the stress 11 σ drops to zero when 11 1 D = for uniaxial tension in direction 1 x . Other approaches for damage modeling can be found, for instance, in [14 -19].

Damage in the plate with the central circular hole
A finite isotropic flat plate with a central circular hole under equi-biaxial and uniaxial loads has been considered. In this work, a second order symmetric tensor ij D is used as a state variable to account for the anisotropy of material damage. To estimate the damage distribution, first the plate with the central circular hole under equi-biaxial and uniaxial tension in the material with the constitutive equations (6) was considered. The distribution of the damage tensor component D 11 under equi-biaxial and uniaxial tension is shown in Figure 1 and

Anisotropic damage in the plate with the central crack
In this part the plate with the central crack under tension in the damaged material with constitutive equations (6) is considered. The typical mesh with singular finite elements around the crack tip is shown in Figure 3.

Damage in the plate with the central inclined crack: mixed-mode loading
In this part of the study numerical experiments are performed on the mixed-mode I/II brittle fracture by edge cracked semicircular bend specimens in the damaged material with the constitutive equations (6). The schematic presentation of the specimen geometry and load is shown in Figure 17.         Distribution of the anisotropic damage tensor components is shown in Figure 18.

Conclusions
The detailed finite element analysis was performed for analyzing the effect of damage accumulation process in the cracked specimens. The experience of using the user subroutine UMAT for FEM package SIMULIA Abaqus/CAE for damage accumulation processes in the vicinity of the crack is presented and discussed. A continuum damage mechanics model based on the constitutive relations of linear elastic isotropic materials with the incorporated damage second rank tensor components is used to describe the material behavior. The material nonlinearity arising from the deformation process is modeled by introducing an anisotropic damage tensor of the second rank into the constitutive equation. The material model is described by means of user procedure UMAT of SIMULIA Abaqus. The finite element (FE) mechanical constitutive model is implemented in Abaqus/Standard via a UMAT routine. Numerical experiments for a large series of cracked specimens have been performed. Computed stress and damage tensor components were found. They are not dependent on the FE mesh refinement. Distributions of the damage tensor components in the vicinity of the crack tip in cracked specimens of different configurations under mixed mode loading in a wide range of mixed mode loadings are found. The configurations of active damage accumulation process zone in the cracked specimens are obtained. The damage accumulation process has substantial influence on the stressstrain state in the vicinity of the crack tip and leads to decrease of the stress concentration in cracked specimens.

References
[1] Sarkar S, Singh I V, Mishra B K, Shedbale A S and Poh L H 2019 Finite Elements in Analysis