Prediction of Anisotropic Elastoplastic Instability with Rice’s Criterion in Plane Compression

Objectives: The main objective of this work is the development of a mechanically coherent framework for the implementation of numerical simulation tools for sheet metal forming processes. The work at hand is part of a study of the mechanical behaviour of metallic materials subjected to large elastoplastic deformations during channel-die testing and a study of the mechanisms that limit their formability and cause their instability. We analyse the effect of induced rotation on the appearance of plastic instability predicted by Rice’s criterion. Methods: At the beginning, an adequate modelling of the mechanical behaviour of metals in large elasto-plastic deformations using the formalism in rotational referential is made. Then, a prediction of plastic instabilities is developed using Rice’s bifurcation criterion. Indeed, the integration of the law of anisotropic elastoplastic behaviour is presented in a three-dimensional kinematics framework in order to write the deformation tensor as a superior triangular matrix, with Hill’s anisotropic yielding function, using Runge Kutta’s method to integrate the obtained equations. Findings: we were able to show that Rice’s criterion can be used for any kind of material and that this criterion is effective in predicting what happens in metals. It predicts shear-band bifurcation in a continuous plastic medium. Hill’s criterion is suitable for what is weakly anisotropic. The model described above represents a good mathematical framework for analysing the behaviour of materials and plastic instability. Novelty: This study allows to apply Rice’s criterion in the case of a channel-die compression test, so as to locate plastic instability and see the appearance of shear bands for all types of materials; it demonstrates that the initial rotation of the material influences the appearance of shear-bands during a test simulating rolling. 
Keywords: Elastoplastic; Large deformations; Anisotropy; Plastic instability; Localisation


Introduction
The simulation of large deformations for elasto-plastic materials remains indispensable for forming and the calculation of structures. This is a problem that is always on the agenda due to its importance in applications, where industrial operations are limited by the appearance of deformation localization phenomena and plastic instabilities (1) (2) (3) (4) . Shear bands are a quite common phenomenon.in all sorts of solids. In geological materials, their role in fracture has been revealed by numerous experiences, not only in rocks (5) but also in granular media like sands (6) . Polymers, especially transverse anisotropic fiberreinforced polymers (7) , are sensitive to that form of instability. Even materials such as oxide layers that grow from a foam of iron in batteries offer analogous phenomena (8) .
In this paper, the focus is set on metallic materials. Within their grains many micro shear bands appear which eventually grow into macro shear bands, provoking failure. A distinctive phenomenon is that neither the micro bands or, all the more, the macro ones are crystallographic. They are not due to the coarsening of the sliding on slip planes. Now, the appearance of shear bands is often related to the geometrical features of the specimen. This is the case in (9) , which studies a notched bar submitted to tensile stresses in which shear banding takes place in the neck of the test-piece. This is not the approach taken in this paper. No boundary conditions favour the shear bands. They can initiate at any place in the material. In real solids, they often form in the vicinity of some material heterogeneity within the specimen, for example a metallic inclusion.
The appearance of shear bands is always related to the fact that the applied stresses do not have the same principal axes as the specimen. This happens in (9) where the tensile test is performed a few degrees away from the axis of symmetry of the testpiece. This is the view adopted in the present paper. An orthotropic material is considered (90 % of the metal produced in the world being rolled, orthotropy is quite common among metallic pieces). It is submitted to channel-die compression, but rotated 30 • from the direction of the flow. To predict instabilities, this work studies anisotropic elasto-plasticity in large deformations and rotations. It deal with this behaviour by a phenomenological approach, whose kinematics is described by a formulation in a rotating frame of reference, using the objective derivatives of Jaumann (10) and Green Nagdhi (11) . To facilitate numerical resolution, this study uses three dimensions kinematics in order to choose a rotating material frame of reference allowing the representation of the deformation gradient and the deformation velocity tensor by triangular upper matrices.
The proposed approach is validated by an application in the case of the channel-die compression test. Thus, the phenomenological elastoplastic behaviour with isotropic strain hardening has been numerically integrated and coupled with a criterion for predicting instabilities based on Rice's criterion. This criterion has the advantage to be three-dimensional and to apply in the elastoplastic case. Contrarily to other criteria used in early works, it is not limited in its field of application either by its two-dimensional character or by its limitation to the rigid plastic behaviour (12) . At the beginning of the research in this field, the work of band-based localization prediction was based on the study of stationarity conditions of acceleration waves in linear and non-linear elastic solids (13,14) . In (15)(16)(17) . Rice (4) adapted an approach equivalent to the first for elasto-plastic materials. In the case of incremental plasticity associated with a softening regime for some materials, the bifurcation is associated with the loss of ellipticity of the equations of the system. The use of non-associated plasticity laws makes it possible to predict early localization in the positive work-hardening regime (18) , for example: • Rudnicki et al. (15) were the first to propose a general criterion for the localization of the deformations of a thin sheet under a state of the plane stresses. Their analyses predict the localization of plastic deformations in regions of negative deformation. • in (19)   For this case of orthotropy we use the notion of the isoclinic configuration proposed by (21) whose orientation is parallel with the orthotropic axes ( Figure 1 ). https://www.indjst.org/ where: In the case of small elastic deformations So we can rewrite the relationship: In this isoclinic configuration, the tensors should be represented by Lagrangian quantities For that we will adopt here the revolving referential formulation (23) defining a matrix of rotation Q and the spin Ω =QQ T . According to this formulation, the equation of the rotation is written: where W is the spin of the global rotations. Writing the tensor w defines the choice of the rotating frame. Two cases will be considered in this work: -corotational Jaumann derivative such as: -spin W based on Green Naghdi's derivative.
where the tensor Γ i jkl is written in (44). In this work the orientation Q will be defined by means of three angles of Bunge(φ 1 , φ, φ 2 ) (24) : see Figure 2. These angles define the rotation Q through which the reference orthotropic , coincides with the Eulerian reference, such that: https://www.indjst.org/ The rotation rate is written:

Three-Dimensional Kinematics
We have adopted a three-dimensional kinematics writing which allows us to define a rotating material frame of reference whose deformation gradient and velocity gradient tensors are described in the form of an upper triangular(tr) matrix. The original idea of a material rotation frame was due to Mandel (21) . When researchers began to deal with plastic instabilities in finite element codes, its necessity appeared because of the large deflections and rotations they were simulating (25) . In codes such as Abaqus TM , which uses a Lagrangian formulation, the UMAT subroutine allows to program customized laws of behavior in such a rotational kinematic framework (7) . In order to avoid the problems encountered in the integration of the behavior laws in large deformations during non-triaxial tests, the main objective of this formulation is to write them in the form of an upper triangular matrix, as shown by (23) and (26) : where ( − → E i ) and ( − → e i ) are respectively two bases associated with the orthonormal coordinate system in the initial and current configuration.
The strain velocity gradient tensor L is then defined in the base ( − → e i ) in the form of an upper triangular matrix: For a gradient of deformation F defined by: we have the following writing in the deformed configuration: To configure the distorted configuration, we just use the orthonormal coordinate system − → e i . Consequently, let f i be: The transformation gradient F and the velocity gradient L can be written in the following triangular matrix form:

Constitutive Laws
The introduction of the rotating reference system allows us to define the tensors . . whose transformation laws are as follows: Where Q is the matrix defined in (8) by the three Bunge angles (φ 1 , φ, φ 2 ).
To describe the law of behavior, the criterion used is Hill's (15) defined by: where: and − H it is a (6x6) matrix given by: https://www.indjst.org/ so that: with Swift's law of hardening: where the coefficients F, G, H, L, M, N are Hill's coefficients. The formalism of the generalized standard material allows us to write.
These relations are supplemented by the elastic law which is written in the orthotropic reference system by: − C is the tensor of elasticity in the orthotropic system and − D e is the strain rate tensor. The plastic multiplier is written: and equation (23) becomes: and or: https://www.indjst.org/ In the Eulerian configuration this is written: with: Writing the decomposition of the strain rates: This can be written in the form: with: The objective derivative is written: or: with: and A ep is the elastoplastic tensor. Two cases can be distinguished: i) Corotational referential (Jaumann where w = 0, in which case the law is written: https://www.indjst.org/ and: ii) Referential based on Green Naghdi's derivative.with: Using the first Piola Kirchhoff tensor, we can write the incremental law in the form: with:

Application to the Channel-die Compression
Let us now consider the channel-die test in which the height of the sample changes from h to h 0 . We note ε 33 = ln ( h 0 h ) , (ε 33 > 0) . We assume that the deformation of the material is homogeneous. In particular we have ∂ x 33 33 . This is a test that fits well within the framework of three-dimensional kinematics. So we can then write: [F] =   e ε 33 γ 3 e ε 33 γ 2 e ε 33 0 We easily relate the quantities γ 2 and γ 3 to the angles χ 2 = α and χ 3 = β in Figure 3. These angles describe the inclination of the edges of the crystal, initially a rectangular parallele piped: Since F is a triangular tensor, the velocity gradient L is calculated using the expression F and it is equal to: The stress tensor is such that the faces perpendicular to axis 1 are free, i.e. T 11 = 0. We also assume the absence of friction T 12 = T 13 = 0. Hence:

Rice's Criterion
In this study we limit ourselves to the shear band phenomenon, for which the instability is characterized in terms of bifurcation of the equilibrium path (4) . From a homogeneous material deformation, the onset of a shear band will divide the material into two zones, the band and the matrix: see Figure 4. Kinematically, this is described by a discontinuity △Ḟ of the deformation gradient rate across the band from the matrix. These instabilities take place at time t = t c which is the moment of the passage from a homogeneous state to a non-homogeneous one and which corresponds to a discontinuity in the velocity gradient L. Now, L varies through a band inclined at an angle ψ with respect to the direction: see Rice (4) . The onset of shear bands has been linked to a number of material factors. For metallic alloys, a review has been proposed by Pineau et al. (28) . These authors examine phenomena such as material hardening (multiplication of dislocations), textural (or geometrical) hardening, presence of porosities, imperfections in the crystalline structures. Such issues certainly enhance the susceptibility of the alloys, but sole mechanical considerations suffice to predict shear banding without introducing other causes of instability. This is the view adopted in this paper with the implementation of Rice's criterion. It consider shear banding solely as a bifurcation in the response of a time-independent solid submitted to stress. When time t < t c the material is stable and the field L has no discontinuities.

Results and Discussion
This section presents two analyses: first the analysis of constraints based on the responses of stresses and shear strains as a function of the shortening deformation, then the analysis of instability followed by a discussion of the results. The numerical resolution of the system of differential equations is done using the 4 th order Runge Kutta method. The programming is written under the Matlab environment and the prediction of instability is done by a test of Rice's criterion.

Constraints Analysis
In this section, we focus on behavioural change and the effectiveness of the proposed criterion, in order to analyse and discuss the numerical results obtained and to simulate the anisotropic elastoplastic behaviour independent of rotation speed.
We illustrated the results in the case of a power type hardening law defined by: To demonstrate the feasibility of the computation of the mechanical behaviour, Hill's coefficients had to be specified. The example chosen in the literature was a mild steel sheet studied by Habracken et al. (29) . These authors adjusted the values of parameters F, G, H, L, M and N to the measurements of the Lankford coefficient. The values are typical of the mild orthotropy common among commercial products. Orthotropy is inevitable because sheets are rolled but manufacturers do their best to limit it.
To show the influence of anisotropy on the precocity of shear banding, the axes of orthotropy of the material must not coincide with those of the channel-die. Now, they are determined by three Bunge angles φ 1 , φ and φ 2 . As noted before, φ 1 was set to 30 • . There was no point in varying simultaneously the two other angles. So φ was set at 90 • and φ 2 ranged from 0 to 180 • .
The material parameters (29) were fitted on the Lankford coefficients Their values are: M=N=L=2.5091, F=0.5395, G=0.5526, H=1.4474 with: τ0=0.01, h0/E=0.007 and m=0.24 as hardening parameters. E=700 GPa is Young's modulus. On the figures, the values of the stresses are divided by Young's Modulus E. Figure 5 shows that the curves in the case of plastic rotation and in Jaumann's case are very close. The order of hardening observed for normal stress is qualitatively consistent with (27) experimental investigations, especially for the less hard materials. This remark shows that hardening plays an important role. Figure 6 shows the evolution of the lateral stress as a function of shortening in the two cases of Jaumann's and Green Naghdi's derivatives. At the beginning the curves are almost identical but from 1% strain on, the curves separate with a significant increase during the deformation. Nevertheless, the values of lateral stress remain lower than those of normal stress, a well expected result.
https://www.indjst.org/ As pointed out above, the use of an incremental constitutive law implies the choice of an objective derivative, for example Jaumann's and Green-Naghdi's. Both have been tested in the present work. Figures 5, 6, 7, 8 and 9 have shown that their influence on the prediction of the material behaviour is limited but real. On the contrary, they have practically no effect on the onset of shear bands. That is why, in the following section, Figures 10, 11, 12 and 13 will be presented without a reference to the choice of the objective derivative. Indeed, they have been drawn with Jaumann's scheme.
In Figure 7, we notice that the ratio σ 22/σ 33 obtained starts from 0.3 to stabilize afterwards around 0.52 in the case of Jaumann and 0.42 for Green Naghdi. This was observed in the case of channel-die compression applied to single crystals (27) (30) and explained by the phenomenon of relaxation due to lateral pressure at the edges of the channel-die as the deformation evolves. Figures 8 and 9 show the numerical curves of slip angles χ 2 and χ 3 as a function of shortening ε. They show a slight increase in the angles for the tangential deformation values in the Green Naghdi case and almost no change in the Jaumann case.

Instability Analysis
In this section, we present an analysis of the results obtained when integrating Rice's criterion with anisotropic elastoplastic phenomenological behaviour and isotropic strain-hardening in the case of a planar channel-die compression test. The analytical approach proposed above led to the results below. Other authors have integrated only the forming limit diagrams flow https://www.indjst.org/ limit curves (FLC) in deformations by applying Rice's criterion to the elasto-plastic strain hardening in uniaxial and plane tension (31,32) .
Our approach consists in stopping when the determining function (58) passes through a minimum or else changes its sign. We trace the evolution of the stress and the force under control of the change of orientations. We start with a predictive calculation of instability. If Rice's condition of bifurcation is fulfilled, we pass to a new increment of the orientation by a jump of 15 • .
The results of the numerical simulations can be found in Figures 10 and 11 and Table 1, obtained for the case of a plane compression applied with a different initial orientation. It shows that the anisotropy has a significant influence on the behaviour and thus on the time of the appearance of plastic instabilities and that it is influenced by the form of the hardening evolution. Figures 11 and 12 show the decisive effect of anisotropy on the point of the deformation path at which shear banding occurs. It varies from ε = 0.25 to ε = 1.29. Since φ = 90 • , angle φ 2 gives the direction of the applied compression (vertical in the channel-die) in the (1,2) plane of the sheet, that is, between the rolling and transverse directions. For most applications, and especially deep-drawing, this is the key factor. There are definitively directions of this plane in which shear bands are more likely https://www.indjst.org/  The analysis of two initial orientations (φ 1 and φ are fixed and φ 2 varied between 0 • and 180 • with a step of 5 • ) studied and schematized in Figure 13 shows the effect of the evolution of the angle φ 2 on the appearance of plastic instabilities. The studies on single crystals confirm this evolution as well as the effect of microtextures. https://www.indjst.org/

Conclusion
The present paper shows how an orthotropic material submitted to channel-die compression can be modelled in the framework of large strains and how it is affected by shear banding. Several points came out when doing so.
Large strain formulation involves the use of numerous tensors but the calculations are greatly simplified by the introduction of a material rotation frame in which the deformation tensor takes the form of an upper triangular matrix.
Rice's criterion, which necessitates the large strain formulation, predicts efficiently shear banding as a bifurcation in the deformation path. It is linked only to considerations of continuum mechanics, regardless of material factors often presented as the source of these instabilities. In fact, they only accommodate it according to the specific microstructure of the considered solid.
Along with numerous previous works, the paper shows the sensitivity of shear banding to anisotropy. Its precocity relies heavily on the position of the applied stress with respect to the axes of the material. Here the latter was only mildly orthotropic, with a plastic behaviour represented by Hill's quadratic criterion. Nevertheless, their onset varied from ε = 0.25 to ε = 1.29 in a rotation around an axis. The next step in the research will be to introduce this formulation in a finite element code using a Lagrangian formulation, so as to retain the advantage of the formalism of large strains.