Bauschinger effect modelled by yield surface distortions

A model for distorting the yield surface by flattening the part in the reverse of the loading direction, is suggested. As the basis for the distortion, the model applies a pair of second-order back-stress tensors of similar type as in kinematic hardening models. The yield-surface formulation provides a flattening and shrinkage of a given first-order homogeneous yield surface in the reverse directions of the back-stress tensors. The mathematical formulation is based on similar ideas as the HAH (homogeneous yield function-based anisotropic hardening) model, for which the calibration of the equivalent stress-strain curve is independent of the Bauschinger part of the model. Severe mathematical and numerical challenges of the HAH model are pointed out, but are avoided in the new model. Furthermore, the yield surface doesn't have to contain the origin, and the r-value is conserved in stress reversals.


Kinematic hardening
Kinematic hardening occurs by that the yield surface, ( ), is shifted by the second-order back-stress tensor X. Here is the stress tensor. The adopted model is a combined isotropic and kinematic hardening formulation, firstly introduced by Armstrong and Frederick (Armstrong and Frederick, 1966;Frederick and Armstrong, 2007) and subsequently modified by Chaboche (1986).
Here¯(¯) is the equivalent isotropic hardening and¯is the equivalent accumulated plastic strain. The total stress in a tensile test will consist of contributions both from the isotropic yield surface expansion and from the kinematic back stress. Associated flow is assumed, hence Here denotes the stress gradient, i.e. / ij , and D p denotes the plastic rate-of-deformation tensor. Since the yield function (1) is homogeneous of order one, =¯p, where p is the equivalent plastic strain rate. In the original paper (Armstrong and Frederick, 1966;Frederick and Armstrong, 2007) the direction of the stress tensor, denoted X is following the direction of D p , i.e. at a given proportional strain path these two tensors will be aligned. This can be formulated as follows: For an undeformed material the initial condition is X 0. In this model the back stress will remain deviatoric. For application of an isotropic yield surfaces for a tensile test, and for steady back-stress solutions, = d dt X/ 0 , the tensile stress equals + X1 1 tensile . More general cases can be obtained with X ( ) p . The other parameter p controls the strain increment of the reversal transients. One interesting special case is when X remains proportional to the equivalent stress ( ) p . Then, only one proportionality constant is required to describe the back-stress attractor in Eq. (3). In a reversed strain path subsequent to a tensile test, the compression stress Note for the strainrate-space based formulations, where X is attracted by D p , that with a non-isotropic material with r-value different from unity, the r-value of a tensile test pre-strained in uniaxial compression, will have a different instant r-value than a tensile test without the Bauschinger prestrain.
An alternative that preserves the r-value in Bauschinger tests, is a stress-space based formulation of Eq. (3), where the direction of the back-stress tensor, X, is following the direction of the deviatoric stress tensor = S I tr( )/3, where I denotes the second order identity tensor: Also in this model formulation, the back stress remains deviatoric. Since the back-stress tensor and the deviatoric stress tensor are parallel during proportional load in the stress-space formulation, the interpretation of X (¯) as the stress contribution from the back stress to a solution for a tensile test in the equivalent direction with = d dt X/ 0, holds also for anisotropic yield surfaces in this case. When X is proportional to the equivalent stress,¯, the special case of the stress-space formulation (5) can be formulated Here the parameter g has the same interpretation ( / ) The strainrate-space and the stress-space based formulations, as given her, give similar results using the same parameters, but are strictly identical for all strain paths, only in combination with a von Mises yield surface. Both models can be generalized by applying multiple back-stress tensors, e.g. the stress-space formulation then can be written as Here X i are constants interpreted as saturation levels of the equivalent back stress contributions and¯i controls how large equivalent plastic strain each component needs to saturate.

The HAH model
The HAH model, as suggested by Barlat et al. (2011), modifies a given homogeneous yield surface ( ) of degree one. The resulting total yield function is written Here the McCaulay brackets are used, i.e.
. The second-order deviatoric tensor, ĥ , has to be initialized with length = ĥ | | 3/8. During monotonic loading it will align with the deviatoric stress-tensor direction. In a tensile-compression test with ĥ pointing in the deviatoric tensile-stress direction, the parameter g 1 corresponds to the ratio between the forward tensile yield stress with and without distortion of the stable yield surface, while the parameter g 2 corresponds to the ratio between the reverse compression stress with and without distortion. Associated flow is assumed, hence Since the yield function (9) is homogeneous of order one, = p , where p is the plastic strain rate. The second order deviatoric tensor ĥ describes a memory of the previous loading directions. An evolution equation that preserves h |ˆ| is required and the following one was given by Barlat et al. (2011) The attractor for the evolution of the tensor ĥ in Eq. (11) is parallel with the deviatoric stress tensor S. Since ĥ is normal to the yield surface facet (linear plane in the stress space), it means that the orientations of the facets during proportional load become radial in the deviatoric stress space. Because of this, the facet orientation will always correspond to that the r-value of a tensile test pre-strained in uniaxial compression, will have an instant r-value equal to unity, regardless the r-value prescribed by the steady yield surface. Hence the regularized yield surface, Eq. (9), will result in an r-value in-between these two, but closer to unity.
Additional evolution equations are formulated for the parameters g 1 and g 2 describing the magnitudes of the yield surface distortions in the forward and backward directions of ĥ .
The parameters g 3 and g 4 are internal variables that enables the model to capture permanent softening behavior. The initial conditions are = = = = g g g g 1 1 2 3 4 .

The new model
A new formulation is suggested as a mathematically and conceptually simpler alternative than the HAH model, combining best parts of the HAH model and the kinematic hardening models. The basic construction is similar to the HAH model, being based on a regularized combination of the stable, undistorted part of the yield surface, ( ), which is any homogeneous yield functions of degree one, and two additional facets. The distortion and its magnitude are controlled by two second order tensors, which conveniently can be modelled by similar models as in classical kinematic hardening models.

The yield function
The main idea of the new yield-surface construction is to assure that the reverse yield stress is similar as with a kinematichardening model with a back-stress tensor X, but by a distortion of the yield surface rather than by shifting it by kinematic hardening.
In the kinematic hardening model, the back-stress tensor X represents a memory of the recent deformation history and is part of the mathematical description of the current yield surface. The current magnitude and direction of X can be a result of any loading history, but only its current magnitude and direction matters. The further construction of the model is based on an equivalent history of proportional loading leading to X. At the end of this fictive proportional predeformation, the stress would correspond to the forward proportional yield stress point XF on the shifted kinematic hardening yield surface, = X ( )¯k in , in Fig. 1. This uniquely defines, XF for any given back stress X. Next, assume the same shape of the yield surface, ( ), but instead with isotropic hardening, and that the same equivalent preloading resulted in the same stress XF . In this case, the inscribed, smaller kinematic-hardening yield surface in Fig. 1 will touch the same forward stress point XF as the isotropic yield surface. Finally, a distortion of this larger isotropic yield surface is constructed, as shown in Fig. 1. Near the stress point XF (which is defined by X ) there is no yield surface distortion. In the reverse stress direction, one finds the backward stress point XB (which is also defined by X ) on the kinematic hardening yield surface. The new model is constructed to distort the yield surface, so that it exactly goes through XB . Furthermore, in this point the normal direction of the distorted yield surface is opposite of the normal direction in XF .
Similar as for the HAH model, the yield surface will be constructed as a regularized combination of the undistorted (envelope) yield surface = ( ) and two planes (facets) in the stress space. The undistorted yield surface is assumed to be homogeneous of degree one, i.e.
= k (k ) ( ), for any scalar k. The two planes are controlled by the two back-stress tensors X 1 and X 2 , which represent a memory of the previous loading. During monotonic loading, one of the back-stress tensors, say X 1 , will be active and try to catch up and align with the direction of the deviatoric stress tensor S, while X 2 then will be passive and shrink towards the origin. However, in general cases the yield surface is distorted by each tensors X 1 and X 2 independently. Firstly, the case of only one back stress tensor X 1 corresponding to an equivalent forward monotonic prestress XF , will be explained. The deviatoric part S XF of the forward yield stress XF is then found in the X 1 direction in the -plane. Since the undistorted yield surface is pressure independent, Fig. 1. The normal direction of the yield surface at XF is given by the gradient ( ) XF . For the pressure-independent, homogeneous yield surface of degree one, Knowing XF and ( ) XF , the tangent plane in Fig. 1 can be calculated. It will be identical for both the isotropic and the kinematic yield surface. The kinematic-hardening yield surface in Fig. 1 has the same shape but is shifted in the XF direction due to the back-stress tensor = X X 1 in Fig. 1b. Note that the center of the kinematic-hardening yield surface in the shown 11 22 section in Fig. 1a, then equals X X I 1 133 . In the reverse direction of XF , the yield stress corresponding to the kinematic-hardening yield surface, equals = X 2 XB 1 XF , and the normal direction of the tangent plane of the yield surface at this point is 1 . The corresponding tangent plane running through XB , can mathematically be expressed: The new yield surface will be a regularized combination of two tangent planes for X 1 and X 2 , as described by Eq. (16) for the case of X 1 , and the stable homogeneous yield surface, ( ). The resulting yield criterion is formally given similar as Eq (9) The regularization of the stable yield surface and the two tangent planes is made by the following mathematical formulation of the yield function ( ): It is assumed that the back-stress tensors X 1 and X 2 are deviatoric, and for a material in the undeformed condition, their initial conditions will be = X X 0 1 2 . The regularization of the stable yield surface and the two facets is mainly controlled by the regularization exponent p.
In the regularization of the facet described by Eq. (16) and the undistorted yield surface , the resulting yield surface will be the innermost one, except where the two regularized yield surfaces overlap or intersect, e.g. forming corners. The corners will be rounded and thereby the regularized yield surface will here locally be inside both the facet and the undistorted yield surface . In cases where X 1,2 vanishes, the regularized yield surface will be affected in the reverse direction of X 1,2 , something that is not desired by the model. (19) is included to assure that the yield surface distortion caused by X 1,2 completely vanishes when X 0 , utilizing that at this stress, As for the HAH model, associated flow is assumed as described by eq. (10). However, unlike the HAH yield function in Eq. (9), the total regularized yield function in Eq. (18) is not homogeneous of order one. Hence, in general cases,¯p. Similar as for the HAHmodel, the regularization of the two facets and the stable yield surface , as described by Eq. (17-19), ensures a convex yield surface . Note also, that alternatively, to account for corner effects, a non-associated yield criterion, e.g. by (Yoshida, 2017), can be applied instead of Eq. (10), but this will not be focus here.

Evolution of the microstructural tensors
In classical kinematic hardening models, the yield surface is shifted using only 1 s-order back-stress tensor. To obtain a similar behavior by the distortions of the yield surface by two facets, two back-stress tensors are applied, one in the forward and one in the backward direction. However, in monotonic strain paths, only one of them will be active and will then evolve in the same manner as in kinematic hardening models. The other one will simultaneously decay towards zero magnitude. The two kinematic back-stress tensors are suggested to follow two coupled differential equations in a stress-space formulation: The function X (¯) and the parameter p have the same interpretations and similar but not identical calibrations as in the kinematic hardening model, Eq. (5). In this model, the two back-stress tensors will remain deviatoric.
In Eqs. (20-21), the cases of monotonic deformation will develop a back stress by one of the back-stress tensors X 1 or X 2 . This one will be parallel with the deviatoric stress, while the other one will remain zero. When starting from an annealed material, i.e. with = X X 0 1 2 , the evolving tensor will be X 1 due to that initially then S X S X 1 2 . If the stress S is suddenly changed, the one amongst X 1 and X 2 that is closest will start evolving towards the direction of S, while the other one will start shrinking towards zero. The special case, for which the X is proportional to¯, is of interest for comparison with the HAH model.  (21): In applications, several back-stress tensors are sometimes adequate for a good description of the material behavior. The evolution Eqs. (20-21) for X 1 and X 2 can be generalized to a series as follows: Now, X i may not be functions but must be parameters, similar as i p .

Isotropic hardening, permanent softening and reversed-stress transients
Transient hardening and permanent softening occurring after reverse strain paths may be modelled through modifying the work hardening of the equivalent stress. Following the models suggested by Manik et al. (2015), the equivalent stress is split into three contributions.
Here R is the path-independent isotropic hardening, commonly modelled as a function R (¯) p of the equivalent plastic strain. The following relation will be applied here: Here K, 0 and n are constants to be calibrated. The second term, R r p adds a negative contribution that leads to permanent softening of¯due to reduced work hardening during the strain-reversal transient. The last term, R r t , is a positive, transient hardening contribution to the isotropic stress¯and contributes only during strain-reversal transients. The following models are adapted from Manik et al. (2015).
Here k r and q r and r are constants. The mechanisms described by Eqns. (28) and (29) are believed to be caused by changes in the microstructure during the transient. Hence, in cases where the rapid variation is captured by X 1 (1) and X 2 (1) it sometimes may be better in Eqns. (28) and (29), corresponding to the specific microstructural mechanism described by this second expansion term.

Results and discussion
With prestrains by tensile tests in the x 1 direction, the yield surface in the 11 22 section will be shifted or distorted as shown in Fig. 2, where the isotropic Hosford/Hersey yield surface is applied as the steady yield surface.
Here i , i = 1,2,3 are the principal stresses. An exponent = a 6 is applied in Fig. 2. All yield surfaces have the same tensile stress point in common, and with increasing back stress the size of the kinematic-hardening yield surface (given by¯) decreases. For the HAH model and the new model, however, the part around the tensile direction of the yield surface, remains unaltered, while the reverse part is flattened and softened. Unless otherwise specified, the isotropic yield surface by Eq. (30) with = a 6 will be used throughout the examples.
The orientation of the flattening plane for the considered tensile pre-strain mode in Fig. 2, happens to be equal for the HAH model and for the new model. However, the orientation of the HAH facet will always be radial in the deviatoric stress space, hence for a general, proportional stress path, these planes will have equal orientation by the two models only with the von Mises yield surface. In the new model the plane will be prescribed by the yield surface tangent plane as explained in Fig. 1. Note, that a simplified version of the new model with radial orientation of the plane facets as in the HAH model, can easily be obtained by replacing all X ( ) i by X i in Eq. (18).
The distortion by the HAH model is limited by that the origin in the stress space must be inside the yield surface. Therefore, with increasing magnitude of the Bauschinger effect, the yield surface develops sharp corners due to the regularization. The roundness of these corners at the two ends of the flattened part of the yield surface, as seen in the 11 22 section in Fig. 2, decreases mainly with increasing magnitude of the Bauschinger effect and is only to a small extent controlled by the yield-surface regularization exponent q in Eq. (9).
On the other hand, the new yield-surface formulation is not limited to contain the origin in the elastic region. Furthermore, the roundness of the introduced corners of the flattened parts of the yield surfaces in Fig. 2 is controlled solely by the yield-surface regularization exponent p in Eq. (18), in a similar way as the roundness of the biaxial stress corner in the 11 22 section of the of the Hosford/Hersey yield surface in Eq. (30) is controlled by the yield-surface exponent a. It will not introduce new challenges in iterations and stability in numerical elasto-plastic implementations. Hence, the new model is well suited for implementations into simulation software.

Reversed stress direction
In the paper introducing the HAH model (Barlat et al., 2011), they tested the model for the "generic" case. The HAH-model parameters are listed in Table 1, along with calibrated parameters for the kinematic hardening model and the new model. The resulting compression-tensile curve is shown in Fig. 3. Note that in this special case, the result is independent of the yield surface regularization exponent, both for the HAH model (q) and for the new model p ( ). The shape of the stress-strain curve during the transient work hardening after the stress reversal in Fig. 3 is very similar for the HAH model and the kinematic hardening model, whereas the new model has a straighter hardening curve during this transient. Because the kinematic hardening term adds a stress contribution directly to 11 , and since X 11 increases fast at small strains, the prediction of a higher stress during the initial part of the prestrain by the kinematic hardening model could not be avoided in the calibration of the work hardening parameters of Eq. (27). The isotropic hardening coefficient n is in this case chosen as large as can be justified, i.e. = n 1, which gives linear isotropic hardening. The shape of the stress-strain curve in this initial part could be improved by introducing several back-stress tensors X i ( ) , but for the sake of simplicity, only one term is included in this calibration. This illustrates an important advantage of the HAH model and of the new model, where the calibration of the forward tensile-test work hardening is independent of the calibration of model parameters controlling the Bauschinger effect.
Note, that if the yield surface is anisotropic with an r-value different from unity, then the flattened reverse part of the HAH yield surface corresponds to an instant r-value close to unity, when the stress direction abruptly is reversed, from uniaxial compression to tensile mode. The total, regularized yield surface will then give r 1, where the flattened part is dominating the shape of the yield surface.
On the other hand, neither the new model nor the kinematic hardening as formulated in the stress space, will change the r-value. However, when formulated in the strainrate space Eqns. (3)-(4), the kinematic hardening model will change the r-value significantly, even during proportional loading. A calculation was made for the same compression-tension test as in Fig. 3 with the "generic" calibration but using the Hill48 yield surface with planar isotropy in the 11 22 plane and normal anisotropy. The Hill48 yield surface can then be written The calculations of the r-value evolution for the compression-tensile test are compared in Fig. 4a for the different models, using r = 4 in Eq. (31). Note that the r-value is here defined as = r D D | / | p 22 p 33 , so it is well defined also for the uniaxial compression mode. For the given model the strainrate-space based formulation of the kinematic hardening model saturates at = r 2.8 due to the nonradial direction of the back-stress tensor in the deviatoric stress space. The resulting shifts of the yield surfaces after the prestrain can be compared in Fig. 4b for the strainrate and stress-space based kinematic hardening models, amongst which the stress-space based one retains the r-value. The HAH-model keeps the r-value equal to 4 during the initial proportional loading, since the stress then remains at the undistorted part of the yield surface. However, right after the strain reversal, the flattened part of the yield surface dominates, and the instant r-value reaches almost unity, due to the radial orientation of the flattened facet. As seen in Fig. 4b, this direction is clearly different from the one by the new model, which doesn't alter the r-value. As the flattened part expands, the undistorted yield surface takes more and more over, and the r-value gradually changes back to the undistorted yield surface with = r 4. Fig. 5 shows the yield surface evolution in the π-plane for the HAH model subsequent to the same reverse strain-path change as in Fig. 3. Note, the asymmetry between the slower flattening in the positive 11 prestrain direction as compared to the faster expansion in Table 1 Calibrated model parameters for the cases considered: "generic" from Barlat et al. (2011), "TWIP" from Choi et al. (2018) and "AA1200" from Manik et al. (2015).  , applying the three models for the "generic" case based on the isotropic undistorted yield surface ( = a 6).
the opposite compression direction. The corresponding yield-surface evolution by the new model is shown in Fig. 6 for three different yield-surface regularization exponents = p 2, 4 and 6. Unlike the HAH surface, the flattening of the prestrain direction and the expansion of the opposite second loading direction happen at quite similar speeds, hence it appears as if the yield surface is shifted in this direction during the reversal transient, almost like in the kinematic hardening model, but undistorted in the orthogonal directions.
Right after the prestrain, the yield surface with regularization exponent = p 6, is close to the HAH model with exponent = q 2 for this particular case. As discussed above, the sharpness of the corners at the borders of the HAH flattened part, will depend on the magnitude of the yield surface distortion, as seen in Fig. 5. On the other hand, the corresponding corners of the distorted yield surfaces with the new model in Fig. 6 becomes rounder with decreasing exponent p.

Orthogonal change of the stress direction
A stress-path change can be measured as the angle S S S S arccos( : /(| | | |) between the deviatoric stress tensors S (0) before and S (1) after an abrupt change of loading, see Barlat et al. (2011). It is of interest to see the consequences of the model in the orthogonal°90 directions in the stress space. Some materials have pronounced transient behavior when the loading direction is abruptly changed into orthogonal directions, which requires dedicated model terms for e.g. cross hardening. Then, from a practical point of view, contributions from the terms in the model that are dedicated to take care of the Bauschinger type of behavior, are undesired, unless they describe the orthogonal transients without the need of additional terms. For a material without orthogonal transients, one does not want undesired effects to occur in this portion of the stress space. Without specific experimental observations, one rather wants as little as possible to change here, as compared to the isotropic yield surface. Fig. 7 compares the results for the considered models for an orthogonal stress-path change. Tensile loading in the 11 direction is followed by a biaxial loading, for which = 11 1 2 22 in the positive quadrant of the yield surface. The calibrations of the model parameters are the same as for the generic case in Fig. 3 for the reversed path.
Firstly, it is noted that the HAH model has two different solutions in this limit, obtained if the stress-path angle is infinitesimally larger or smaller than°90 , respectively. This type of mathematical issue for the HAH model was recently reported by Qin et al. (2018) for a more complex case, involving permanent softening. However, it occurs even in this simple "generic" case, taken from Barlat et al. (2011). Amongst the two solutions, the one approaching the°90 limit from below will in most cases be the desired one that does not change the yield surface in the orthogonal portions of the preloading direction. The other one introduces work softening during a strain transient, before approaching the isotropic yield surface solution at larger strains. Anyhow, it is conceptually undesired to have this type of discontinuity in the work-hardening behavior for neighboring stress paths, and it will cause numerical problems in finite element implementations.
The reason for the spurious ambiguity in the HAH model for the°90 stress-path change can be seen from Fig. 8 showing the HAH yield surfaces corresponding to the stress path in Fig. 7. After the abrupt stress-path change, the flattened part of the yield surface The von Mises strain of each yield surface is indicated. a) The "generic" case in Fig. 3 b) the "generic" case but with = p 4 and 6 (dotted lines). , for the three models for the generic case with the isotropic undistorted yield surface ( = a 6). In this limit two solutions exist for the HAH model. rotates as dictated by the rotation of the second order tensor ĥ . In the°90 limit, the part that is flattened in the reverse direction of the pre-deformation load starts vanishing by rapid expansion in this direction, while flattening of the yield surface starts in the opposite direction. This increasingly flattened part rotates clockwise during the stress-path change. In the°+ 90 limit, however, the flattening in the direction of the reverse of the preload stress direction does not vanish but rotates counter clockwise towards the reverse of the new loading direction during the subsequent straining. The temporary softening seen in the stress-strain curve in Fig. 7 occurs while the flattened part rotates through the°90 loading direction and vanishes when the expansion of the flattened part reaches the undistorted yield surface. In Fig. 7 the solution with the new model and = p 6, provides very similar stress strain curve as the°90 HAH curve in Fig. 8, and the influence of the yield-surface distortion is very small. With a lower exponent, = p 2, the solution is between the kinematic hardening and the°90 HAH curve. As seen in Fig. 9, the yield surface of the new model, corresponding to the stress path in Fig. 7, is not rotating much, but the flattened part developed during the prestrain, expands, while the yield surface is flattened in the reverse direction of the new stress path. This is due to the mutual exchange of the roles of the two back-stress tensors X 1 and X 2 . The sharper corner of the yield surface with = p 6 ensures that in this case, the orthogonal direction is little influenced by the flattened part. Note that the new model applies a somehow similar switch as the HAH model in the evolution of X 1 and X 2 in eqs. (20) and (21). However, two mechanisms are different than in the HAH model. Firstly, the discontinuity in the evolution of the back-stress tensors X 1 and X 2 doesn't occur at orthogonal stress paths, but at smaller stress-path changes, when = S X S X 1 2 . After monotonic  The von Mises strain strains are indicated. a) The "generic" case for the stress path in Fig. 7 with the isotropic undistorted yield surface ( = a 6), b) same, but = p 6.
preloading the discontinuity in the further work hardening will occur in directions closer to the preloading direction, safely in the undistorted part of the yield surface. Hence the solution is not influenced in single, abrupt stress-path changes. Secondly, the rotation of the yield surface is different. Unlike the HAH model, the flattened part of the yield surface will always rotate away from the new stress path. In theory, it may be rare occasions during complex loading histories, while both X 1 and X 2 are nonzero, that the transient work hardening rates at certain stress points at the yield surface (where X 1 and X 2 change roles), are slightly discontinuous for neighboring stress points. However, since the flattened part of the yield surface in the new loading mode always rotates out of this region, the transient difference between the two subsequent work hardening behaviors will be very small. This will occur very seldom and near the undistorted part of the yield surface, hence with a significantly smaller magnitude than with the HAH model. In practice, this seems unlikely to cause numerical issues. Further applications and testing of the model are required to tell.
In Fig. 10 the translation of the yield surface by the kinematic hardening model (stress-space based) is shown for the same case of orthogonal stress-path change as in Fig. 7. The stress direction corresponding to the new stress direction will be significantly affected by the shift, resulting in the transient stress-strain behavior in Fig. 7.

Permanent softening
An example of a calibration of the parameters of the HAH-model with permanent softening for a TWIP steel is taken from Choi et al. (2018). A calibration of the new model was made for the same case, and the resulting stress-strain curves are compared in Fig. 11 for a compression-tensile test, for which the compression prestrain = 0.1 vM . The calibrated parameters of the two models are listed in Table 1. The permanent softening is enabled in the HAH model by modifying the amount of the distortions by evolving g 3 and g 4 in Eqs. (14-15). The new model, on the other hand, imposes an isotropic work softening term R r p to the isotropic stress (26) during the transient as prescribed by Eq. (28). In both cases the permanent softening can be well captured. An orthogonal stress-path change from a tensile test in the x 1 direction to biaxial stress with = 11 1 2 22 was also calculated with the TWIP calibrations for the two models. As explained by Qin et al. (2018), the HAH model has two solutions in this limit, and stress evolution for both cases are shown in Fig. 12. The HAH solution, for which the°90 limit is reached from below, the stress-strain curve is not affected by the local yield surface distortion but remains similar as with the isotropic yield surface. The new model predicts this same solution with the highest regularization coefficient = p 6, while the distortion is contributing with a small amount of the reverse  in a compressiontensile test of an AA1200 commercially pure aluminum alloy. The experimental data (grey lines) are taken from Manik et al. (2015). softening with = p 2. The other HAH solution, for which the°90 limit is reached from above, however, results in permanent softening that is rotated into the orthogonal stress path and causing transient softening similar as reported by Qin et al. (2018).
The commercially pure aluminum experiences a surprisingly complex Bauschinger behavior. Stress-strain curves for an AA1200 alloy by Manik et al. (2015) were calibrated and demonstrate the capability of the new model in Fig. 13. The parameters are given in Table 1. In this case the isotropic reverse stress term R r t contributes during the transient according to Eq. (29). In addition, the isotropic work softening term R r p evolves during the transient according to Eq. (28) and reduces the isotropic stress (26) permanently.
Note that in this calibration two terms are used for the back-stress tensors, for which < + > S X X : ( ) 1 2 is replaced by < + > S X X : ( ) 1 (2) 2 (2) in Eqns. (28) and (29), i.e. the slowest responding back-stress tensor is controlling the permanent softening. In a similar manner as by the model by Manik et al. (2015), the complex Bauschinger behavior can then be well reproduced by the new model.

Conclusions
A new model has been formulated as an alternative to the kinematic hardening models and to the HAH model, for dealing with Bauschinger behavior. The new model avoids the mathematical ambiguity in orthogonal stress-path changes and the sharp corners with large Bauschinger effects in the HAH model. The roundness of the corners in the new model is controlled solely by the yieldsurface regularization exponents. As in the kinematic hardening models, the yield surface can be located outside the origin of the stress space. The suggested model is well suited for implementations into simulations software. Permanent hardening and complex transient hardening during the Bauschinger transient are covered by modifying the isotropic hardening by the aid of the back-stress tensors, as have been demonstrated for selected examples of a TWIP steel and a commercial pure aluminum alloy.