Regularized Yield Surfaces for Crystal Plasticity of Metals

: The rate-independent Schmid assumption for a metal crystal results in a yield surface that is faceted with sharp corners. Regularized yield surfaces round o ﬀ the corners and can be convenient in computational implementations. To assess the error by doing so, the coe ﬃ cients of regularized yield surfaces are calibrated to exactly interpolate certain points on the facets of the perfect Schmid yield surface, while the di ﬀ erent stress predictions in the corners are taken as the error estimate. Calibrations are discussed for slip systems commonly activated for bcc and fcc metals. It is found that the quality of calibrations of the ideal rate-independent behavior requires very large yield-surface exponents. However, the rounding of the corners of the yield surface can be regarded as an improved approximation accounting for the instant, thermal strain-rate sensitivity, which is directly related to the yield-surface exponent. Distortion of the crystal yield surface during latent hardening is also discussed, including Bauschinger behavior or pseudo slip systems for twinning, for which the forward and backward of the slip system are distinguished.


Introduction
A robust and general formulation with a yield surface for crystal plasticity can conveniently be applied in a framework similar to the continuum-plasticity framework for finite-element calculations, e.g., including elasticity and using a return-mapping algorithm. There are two major reasons for choosing this approach. Firstly, in perfectly rate-independent crystal plasticity finite element (CPFEM) formulations a nonunique stress versus strain-rate relation makes the calculation procedure break down, when it is applied into a finite element method for the solution of boundary value problems for polycrystals [1]. The reason is a nonunique slip-system selection. Uniqueness can be obtained with a positive definite hardening matrix in the hardening [2,3]. However, that kind of a limitation of the hardening law in order to avoid an ill-posed mathematical problem prevents the use of relevant constitutive laws. The rounded crystal yield surface does not suffer from this problem [4][5][6].
Secondly, the general metal plasticity framework enables crystal elasticity to be part of the model and allows arbitrary models for the critical resolved shear stress, as compared to the viscoplastic or rate-independent asymptotical cases. The activation of different slip systems with different glide planes is of interest, as well as complex models for latent hardening of the slip systems.
The nonuniqueness of the rate-independent approach disappears with the introduction of positive strain-rate sensitivities individually on each slip system. From a physical point of view, this is more realistic since the critical resolved shear stresses of the slip systems will have an instant strain-rate limit of the viscoplastic model. Until recently, return mapping algorithms for yield surfaces with high exponents have not been numerically stable. In a recent publication, Paux et al. [15] avoid this problem by combining the best part of the two available methods, using the regularized yield surface to obtain a smooth, unique tangent modulus, while using the highly efficient Schmid model for the integration of the single crystal constitutive equations. However, progress has recently been reported related to stable, efficient return-mapping algorithms for yield surfaces with high exponents. Thus far, this has been successfully tested for only for continuum plasticity yield surfaces, see [16]. However, the technology is generic and enables stable implementations also of high-exponent regularized yield surfaces for crystals. The details related to adapting this approach and implement return mapping algorithms for crystal yield surfaces will be reported elsewhere.
So far, crystal-plasticity implementations into the finite-element method or into Taylor type of models have mainly assumed either the strain-rate independent plasticity approach or simple visco-plastic behavior as described by power laws. The main reasons for choosing these approaches are the same as they were for choosing continuum-plasticity models some years ago, i.e., model simplicity and numerical convenience. However, the limitation to these specific models for the critical-resolved shear stress, nowadays become a hinder for dealing with more realistic models from a physical-metallurgical point of view, i.e., for dealing with real alloys, containing various types of particles and solid-solution elements applying models that distinguish the instant and the microstructurally induced parts of the strain-rate sensitivity. The framework discussed in this article, enables formulation of models with latent hardening and several sets of active slip systems. This allows yield-surface distortions as predicted by models dealing with strain-path changes related to Bauschinger and cross-hardening effects, e.g., [17][18][19][20]. Recently progress on similar distortional-hardening approaches in continuum plasticity models has been reported in [21][22][23][24][25][26]. The same framework can be applied to pseudo-slip systems applied in deformation-twinning models [27].
In the current work, the quality of the rounded crystal yield surface as an approximation to the rate-independent one will be explored for various set of slip systems activated. The theory is outlined in Section 2, relating the yield surface to iso-surfaces of the potential for the rate-dependent power-law theory and interpreting the regularized yield surface as an approximation accounting for the instant strain-rate sensitivity of the critical resolved shear stresses. Furthermore, calibration strategies are specified. In Section 3, calibrations are made to fcc slip systems, fcc with nonoctahedral slips and bcc with 24 slip systems. Based on the results and theory, the applicability of the crystal yield surfaces is discussed in Section 4 and conclusions are drawn in Section 5.

Theory
The strain-rate independent crystal yield surface can mathematically be expressed as a smooth inner convex envelope of the linear facets from each slip system. The ideal Schmid assumption can mathematically be expressed as the following yield criteria: Here σ is the stress tensor at yielding, and τ s c the critical resolved shear stress for slip on slip system s. Furthermore, P s = (m s ⊗ n s + n s ⊗ m s )/2, where m s is the slip plane unity normal vector, n s the unity slip direction vector. Single slip system solutions are located on facets and obey the normality rule D p = P s . W p /τ s c , where D p is the symmetric part of the plastic velocity gradient and . W p = σ : D p is the rate of the external plastic work. Only the shape of the yield surface can be derived from Equation (2). The symmetric part of the plastic velocity gradient D p can be uniquely derived from the normality rule and Equation (2) only for cases of one single active slip system, which defines the facets (and hence indirectly also the corners where adjacent facets intersect). However, for solutions at corners or other facet interactions, the active slip systems are not provided by Equation (2), and then D p cannot be uniquely determined solely by Equation (2).

Rate-Independent Regularizations
The inner envelope of the yield surfaces by Equation (2) corresponds to a yield surface with sharp corners. By regularization it can be expressed as one yield criterion, i.e., by a smooth yield surface which corresponds to the yield criteria in Equation (2) in the limit when n → ∞ and with ξ s ≡ 1.
Here f is the yield function. In the absence of corners, the normality rule is now unique everywhere. It follows that the plastic-rate parameter can be found as Note that the yield surface in Equation (3) is centro-symmetric. For cases where the reverse and forward of a slip system have different critical-resolved shear stresses, a slightly more general yield surface can be adapted, counting both slip directions as two different slip systems: Here the McCaulay brackets are used, i.e., x = max(x, 0) = 1 2 (x + |x|).

Relation to Strain-Rate Dependent Formulations
It is interesting to compare Equation (3) to the stress potential for a rate-dependent crystal, obeying the power law [7,9].
In this special case, no true yield surface exists, but shape-invariant flow surfaces for constant values of . W p have the same shape and play a similar role as the regularized yield surface, Equation (3). Note the close correspondence between the strain-rate sensitivity exponent m and the yield surface exponent n.
Note that when ξ s and τ s c are the same for all the slip systems, the regularized yield surface with exponent n will give the same solution as the viscoplastic model with strain-rate sensitivity m at the same plastic work rate. It follows as a corollary that, in this case, the texture evolution predicted by the Taylor model will be the same for these two models.
The instant strain-rate sensitivity influences the yield stress in two distinct ways. Firstly, the magnitude of the critical resolved shear stress changes as a function of the strain rate, i.e., the volume of the yield surface changes. However, for cases with a small strain-rate sensitivity this effect is neglectable, as the critical resolved shear stress then depends only weakly on the logarithm of the strain rate. Secondly, the corners of the yield surface are rounded, locally. Corners are sharper in higher dimensions. Hence the local change of the stress at the rounded yield-surface corners in the five-dimensional stress space becomes significant, even in cases of a very small strain-rate sensitivity. Since most crystal-plasticity solutions are located near these corners, the magnitude of the resulting polycrystal stress tensor will decrease significantly as compared to as predicted by a perfectly rate-independent model with sharp corners. However, the stress itself still will be mildly sensitive to changes of the strain rate, i.e., it is caused by the changed strain-rate sensitivity, not the change of the strain rate. For a real metal at room temperature, the rate-independent model with rounded corners provides a more realistic approximation to the real behavior than a model with sharp corners.
The physical based mechanical threshold strength (MTS) model [28] provides estimates and discussion of the instant strain-rate sensitivity, from which the simplest power-law estimate m = g 0 µb 3 /kT, clearly reveals its dependence on the temperature T; see also [29] for a discussion of the formulation of the MTS model adequate for crystal plasticity implementations. Here µ is the elastic shear modulus, b the length of the Burgers vector, k the Boltzmann constant and g 0 nondimensional activation energy with magnitude of order unity. According to this estimate a realistic strain-rate sensitivity at room temperature for e.g., aluminum is m ≈ 10 −2 , hence a realistic yield surface exponent will be n ≈ 100, and the rounded corner solution will be 1-2% lower in stress than the ideal corner.
With increasing temperature, the strain-rate sensitivity of the critical resolved shear stress also needs to be accounted for, i.e., by a strain-rate dependent model, where the "volume" of the yield surface depends on the strain rate. Furthermore, with smaller yield-surface exponents the regularization itself will affect not only the corners locally, but also the magnitude of the entire yield surface will shrink. This can be accounted for by adjusting ξ s , so that the inscribed regularized yield surface touches the centres of the yield surface facets of the nonregularized yield surface. The distribution of the slip rates on the slip systems controls the shape of the yield surface, while the change of the magnitude of the yield surface mainly depends on the total amount of slip activity, i.e., τ s c depends on Γ, being the solution of s . Hence, efficient models can be formulated based on the regularized yield surface, also for cases at elevated temperatures. As discussed in [10], the rate-dependent plastic potential degenerates in its rate-independent limit m → 0 to one particular solution of the Taylor ambiguity, which corresponds exactly to the solution with the regularized yield surface in the limit n → ∞ . In practice, this limit can be approached with a stable numerical iterative algorithm that can handle a high exponent of the regularized yield surface or a low strain-rate sensitivity in the rate-dependent formulation. The use of line-search in addition to the Newton-iteration scheme seems very promising; see [16].

Fitting Stress Points Projected Radially onto the Facets
With decreasing yield-surface exponents, e.g., at elevated temperatures, a calibration of ξ s becomes increasingly important. As an approximation, using equal values ξ for ξ s for all the slip systems, Gambin [5] proposed Here, N is the total number of slip systems. This estimate was obtained as an approximation for fitting one point on each yield surface facet, namely the one in the radial direction of the normal to the facet and being closest to the origin, i.e., S = 2P s τ s c , where S is the deviatoric part of the stress tensor. The approximately fit gives only one average value ξ for ξ s that is equal for all involved slip systems. For the case of only one family of slip systems, it can be shown that this corresponds to any of the solutions that exactly go through these points, and Equation (8) can be simplified to However, if a fit is made to several families of slip systems with each their constant critical resolved shear stress, there will be different values of ξ s for each of the involved slip systems. This idea can Crystals 2020, 10, 1076 6 of 17 be extended to cases where three or more families of slip systems are activated and is a more precise solution as compared to the estimate in Equation (8) made by Gambin. The solution, Equation (9), is an important special case of more general cases, for which one equation has to be solved for each involved slip system α.
An important special case is when cases involving two families of centro-symmetric slip systems are considered, and each family has one critical resolved shear stress. In these cases, only two different values of ξ s are required. These two values can be determined simply by fitting to one arbitrarily chosen facet of each type, i.e., fitted radially to S = 2P s τ s c for two selected slip systems, s 1 and s 2 , one from each family.
It follows that For general cases, the set of N linear Equation (10) can be solved to derive the coefficients ξ α . Note that the yield surface is convex only as far as all derived coefficients are positive.

Fitting the Midpoint of Each Involved Facet of the Crystal Yield Surface
Each facet of the crystal yield surface corresponds to one slip system and can be defined by the combination of one point on the facet and by that P s points in its normal direction. Equation (3) or Equation (5) are also recognized as a "Facet polynomial" [30]. Following the procedure proposed in that work, the yield surface can be fitted exactly to a certain number of stress points and corresponding symmetric part of plastic velocity gradient tensors. An improved choice, as compared to the approach by Gambin, is to match the center point of each facet, which corresponds to the average of the corner locations of the facet. However, the corner solutions must be found first.
A corner corresponds to a valid deviatoric stress-tensor solution S crn obtained by solving Equation (1) with only five of the slip systems. However, only the corners, for which for all the slip systems, are on the yield surface. All yield-surface corners belonging to a facet, i.e., to slip system s, can be found as the solutions of all possible combinations of five slip systems that obey Equation (2). A corner solution only involves some of the slip systems and is a common corner of all facets corresponding to these involved slip systems. The midpoint of a facet corresponds to the average of the stress tensors of all its corners, and the normal direction of the facet is parallel to P s . These midpoint stress points will be denoted S k , one for each slip system, i.e., 1 ≤ k ≤ N. The coefficients λ s of this calibration can be found by solving the following linear equation system N s=1 ξ s S k : P s τ s c n = 1 (14) Note that the yield surface is convex only if all coefficients ξ s ≥ 0, which is always the case with a sufficiently high exponent n. The system from Equation (12) can efficiently be solved iteratively. If the exponent n is too low, but one still wants to proceed, the points corresponding to the negative ξ s can simply be excluded to obtain a convex surface. When accounting for both slip directions of each slip system, one must include both directions as independent slip systems with each their ξ s .
For the important cases involving two families of slip systems, only two values of ξ s are required. Then, it is enough to fit the yield function to two facet center points, one to each type of facets, denoted S s 1 and S s 2 for the two families of slip systems, respectively.
It follows that Note that the corner-stress tensor is everywhere divided by the critical resolved shear stress, i.e., without loss of generality the coefficients ξ (1) and ξ (2) can be calculated based on corner solutions, for which one of the resolved shear stresses equals unity and depends only on the ratio τ

Regularized Crystal Yield Surfaces as Approximations to the Schmid Criteria
The application of regularized crystal yield surfaces in crystal-plasticity implementations can be convenient, and calibration of proper polynomial expressions is explained in Section 2. The simplest approach is to put all coefficients ξ s ≡ 1 in Equation (3) or Equation (5). This will be a good approximation for sufficiently large exponents n. The question is, however, how large n needs to be for a given set of slip systems. A more refined approach is the simplified solution proposed by Gambin,Equations (8) and (9), or the generalized version for several slip systems in Equations (10)-(12), which will be referred to as "radial" in the discussion below. Finally, the solution that intersects all the midpoints of yield-surface facets is uniquely defined but requires considerably more work to determine the coefficients, as the midpoints of the facets of the yield surface must be derived first. It is interesting and of practical relevance to look at the details for the most important examples for fcc and bcc slip systems.
An ideal regularization should match the facets of the crystal yield surface but slightly round off its sharp corners, i.e., the magnitude of the stress tensor becomes smaller at the corners. An analogy is an inscribed circle in a square, for which the difference is largest in the directions towards the corners. The corners of a cube are further away from an inscribed sphere. This illustrates the fact that corners become sharper from two to three dimensions, and even more in higher dimensions. Hence the corners in the five-dimensional stress space are sharp, and the error is larger close to the corners than elsewhere. By "error" it is here meant the difference between the regularized yield surface and the Schmid surface with sharp corners. Two estimates will be made to describe how well the regularized yield surface fits the Schmid criteria. The first one is the average error at the corners Here S crn,n is the stress in the same direction in the stress space like the corner of the Schmid yield surface S crn,∞ (i.e., with infinite exponent n), where N crn is the number of corners of the yield surface. The second measure is the average of the errors at the middle of each of the facets of the yield surface.
Here S fac,s is the stress tensor at the exact midpoint of the facet of the Schmid yield surface belonging to slip system s, and N s is the number of slip systems involved in the yield surface.

The fcc Case with {111} 110
The standard fcc case with only one set of slip systems is the simplest to deal with. The solution proposed by Gambin is identically equal to the "radial" solution, where ξ s takes the same value ξ for all slip systems. This solution is compared to ξ, being estimated by the midpoint fitting as functions of n in Figure 1. The average errors at the corners and at the midpoints of the facets are compared in Figure 2. Note that the midpoint fit has no error in the midpoints, i.e., e fac = 0 for all n. For an exponent n larger than about 15, both approaches coincide with the simple approach putting ξ s ≡ 1. As can be seen from Figure 2 this corresponds to about 15% corner error. Cases where n ≥ 15 correspond to a power-law strain rate sensitivity m ≤ 0.07, hence in such cases the simple asymptotic solution ξ s ≡ 1 applies.
Here crn, is the stress in the same direction in the stress space like the corner of the Schmid yield surface crn,∞ (i.e., with infinite exponent ), where crn is the number of corners of the yield surface. The second measure is the average of the errors at the middle of each of the facets of the yield surface.
Here fac, is the stress tensor at the exact midpoint of the facet of the Schmid yield surface belonging to slip system , and is the number of slip systems involved in the yield surface.

The fcc Case with 111 ⟨1 10⟩
The standard fcc case with only one set of slip systems is the simplest to deal with. The solution proposed by Gambin is identically equal to the "radial" solution, where takes the same value for all slip systems. This solution is compared to , being estimated by the midpoint fitting as functions of in Figure 1. The average errors at the corners and at the midpoints of the facets are compared in Figure 2. Note that the midpoint fit has no error in the midpoints, i.e., fac = 0 for all . For an exponent larger than about 15, both approaches coincide with the simple approach putting ≡ 1. As can be seen from Figure 2 this corresponds to about 15% corner error. Cases where ≥ 15 correspond to a power-law strain rate sensitivity ≤ 0.07, hence in such cases the simple asymptotic solution ≡ 1 applies.  Note that for the most rounded yield surfaces with the smallest exponents, the facet errors increase and only the midpoint fit provides a consistent regularization. The corner errors decrease very slowly, and even at = 90 the error remains about 2%.

The fcc Case with 111 ⟨1 10⟩ and Nonoctahedral 001 ⟨1 10⟩
Nonoctahedral 001 ⟨1 10⟩ slips can be activated as long as the ratio of their critical resolved shear stress to the one for the 111 ⟨1 10⟩ slip systems is less than √3. The yield surface has 36 facets corresponding to the 18 slip systems and 186 corners. Figure 3 shows the fitted coefficients for two cases, when the critical resolved shear stress is equal for the two slip systems and when Note that for the most rounded yield surfaces with the smallest exponents, the facet errors increase and only the midpoint fit provides a consistent regularization. The corner errors decrease very slowly, and even at n = 90 the error remains about 2%.  3. The yield surface has 36 facets corresponding to the 18 slip systems and 186 corners. Figure 3 shows the fitted coefficients for two cases, when the critical resolved shear stress is equal for the two slip systems and when The corresponding errors are shown in Figure 4. The errors at a given exponent n are very similar for the two cases, and all the fitting procedures give similar errors for n larger than 10. When the critical resolved shear stresses are equal, the radial fit gives small facet errors and quite similar corner error as the midpoint fit.
When the ratio between the two critical resolved shear stresses equals 1.5, the facets due to the presence of nonoctahedral slips systems contribute to only a small portion of the yield surface, and their midpoints can only be fitted for n larger than about 7. However, their midpoints are still poorly fitted by the Gambin formula and by the radial approach up to values of n close to 30. are used, see e.g., [31]. Figure 5 shows the yield-surface parameter calibrations. With n > 50 the parameters are close to unity. Figure 6 shows that the error at the midpoints of the facets vanishes for n > 20, whereas the error in the corners is about 5% with n = 30, from which it slowly decays proportionally to n −1 .

The bcc
Nonoctahedral 001 ⟨1 10⟩ slips can be activated as long as the ratio of their critical resolved shear stress to the one for the 111 ⟨1 10⟩ slip systems is less than √3. The yield surface has 36 facets corresponding to the 18 slip systems and 186 corners. Figure 3 shows the fitted coefficients for two cases, when the critical resolved shear stress is equal for the two slip systems and when ⟨ ⟩ = 1.5 ⟨ ⟩ . The corresponding errors are shown in Figure 4. The errors at a given exponent are very similar for the two cases, and all the fitting procedures give similar errors for larger than 10. When the critical resolved shear stresses are equal, the radial fit gives small facet errors and quite similar corner error as the midpoint fit.     When the ratio between the two critical resolved shear stresses equals 1.5, the facets due to the presence of nonoctahedral slips systems contribute to only a small portion of the yield surface, and their midpoints can only be fitted for larger than about 7. However, their midpoints are still poorly fitted by the Gambin formula and by the radial approach up to values of close to 30.

The bcc Case with 2 11 ⟨111⟩ and 1 10 ⟨111⟩
It is commonly accepted that 48 slip systems consisting with 2 11 ⟨111⟩ and 3 11 ⟨111⟩ in addition to the basal 1 10 ⟨111⟩ are activated during plastic deformation of bcc steel at room temperature. However, for computational efficiency and convenience, calculations are commonly limited to two sets of slip systems, which provides many glide planes, i.e., close to the pencil glide assumption, hence only the 24 slip systems from 2 11 ⟨111⟩ and 1 10 ⟨111⟩ will be considered in the example here. When √3  [31]. Figure 5 shows the yield-surface parameter calibrations. With > 50 the parameters are close to unity. Figure 6 shows that the error at the

Latent Hardening or Twinning
In cases with strain-path changes, the latent hardening of the slip systems may be important. Models for the evolution of individual critical-resolved shear stresses for each slip system are required. Loss of symmetry will split the corners and introduce more corners of the yield surface. By the individual evolutions of the critical resolved shear stresses, the shape of the yield surface is updated at each integration point in time. In order to model complex Bauschinger behavior, as in [18], a distinction is made between forward and backward of the slip system. This requires the use of Equation (5). Similarly, this is important to distinguish twinning from de-twinning behavior described by pseudo slip systems, e.g., [27]. The stress-differential effect will also imply such a

Latent Hardening or Twinning
In cases with strain-path changes, the latent hardening of the slip systems may be important. Models for the evolution of individual critical-resolved shear stresses for each slip system are required. Loss of symmetry will split the corners and introduce more corners of the yield surface. By the individual evolutions of the critical resolved shear stresses, the shape of the yield surface is updated at each integration point in time. In order to model complex Bauschinger behavior, as in [18], a distinction is made between forward and backward of the slip system. This requires the use of Equation (5). Similarly, this is important to distinguish twinning from de-twinning behavior described by pseudo slip systems, e.g., [27]. The stress-differential effect will also imply such a

Latent Hardening or Twinning
In cases with strain-path changes, the latent hardening of the slip systems may be important. Models for the evolution of individual critical-resolved shear stresses for each slip system are required. Loss of symmetry will split the corners and introduce more corners of the yield surface. By the individual evolutions of the critical resolved shear stresses, the shape of the yield surface is updated at each integration point in time. In order to model complex Bauschinger behavior, as in [18], a distinction is made between forward and backward of the slip system. This requires the use of Equation (5). Similarly, this is important to distinguish twinning from de-twinning behavior described by pseudo slip systems, e.g., [27]. The stress-differential effect will also imply such a distinction, but explanations involving non-Schmid effects will imply a nonassociated yield surface that must be accounted for, see [32].
The simple approach with all coefficients ξ s ≡ 1 is highly efficient when applicable, i.e., at room temperature applying a high exponent, n. Figure 7 shows corner and facet errors for fcc and bcc structures, where the critical resolved shear stress are given a ±5% random variation for cases, where forward and backward slip directions are not distinguished, i.e., with 12 slip systems and correspondingly 12 different ξ s to be calibrated. Here all ξ s are individually fitted. The result will be very similar if forward and backward were distinguished, except that twice as many ξ s would have to be calculated. The calibration quality varies similarly as without this variation. The asymptotic assumption ξ s ≡ 1 applies well for n > 30, and the corner errors slowly decay with increasing exponent, n.
Crystals 2020, 10, x 13 of 18 distinction, but explanations involving non-Schmid effects will imply a nonassociated yield surface that must be accounted for, see [32]. The simple approach with all coefficients ≡ 1 is highly efficient when applicable, i.e., at room temperature applying a high exponent, . Figure 7 shows corner and facet errors for fcc and bcc structures, where the critical resolved shear stress are given a ±5% random variation for cases, where forward and backward slip directions are not distinguished, i.e., with 12 slip systems and correspondingly 12 different to be calibrated. Here all are individually fitted. The result will be very similar if forward and backward were distinguished, except that twice as many would have to be calculated. The calibration quality varies similarly as without this variation. The asymptotic assumption ≡ 1 applies well for > 30, and the corner errors slowly decay with increasing exponent, .

Discussion
The results presented here provides an overview of how large the yield-surface exponent must be for the simplifying assumption = 1, to be valid. For the case of fcc with primary slip systems, this is > 10, i.e., for instant strain-rate sensitivities smaller than ≈ 0.1. Hence, the regularized yield surface can be applied at room temperature and at elevated temperatures up to this limit. In the interpretation of the yield surface exponent in terms of the instant strain-rate sensitivity, the regularized yield surface will provide a more accurate result than the perfectly rate-independent theory, and as a bonus, the Taylor ambiguity vanishes. However, when nonoctahedral slip systems are included, the number of corners increases, which gives the regularization a greater influence on the resulting yield surface. This must be compensated by decreasing , and this must also be done for larger , in the range between 10 and 50. For aluminum the nonoctahedral slip systems are activated at higher temperatures, for which < 10. The approach taken here is to adjust , so that the facet center is touched by the regularized yield surface from the inside. Note, that it is not necessarily the correct choice. Comparison to ratedependent models should be made, but the results here guides when such an adjustment must be made. Note that adjusting is mathematically equivalent to adjusting to /( ) ; hence, this adjustment is vital at elevated temperatures.

Discussion
The results presented here provides an overview of how large the yield-surface exponent n must be for the simplifying assumption ξ s = 1, to be valid. For the case of fcc with primary slip systems, this is n > 10, i.e., for instant strain-rate sensitivities smaller than m ≈ 0.1. Hence, the regularized yield surface can be applied at room temperature and at elevated temperatures up to this limit. In the interpretation of the yield surface exponent in terms of the instant strain-rate sensitivity, the regularized yield surface will provide a more accurate result than the perfectly rate-independent theory, and as a bonus, the Taylor ambiguity vanishes.
However, when nonoctahedral slip systems are included, the number of corners increases, which gives the regularization a greater influence on the resulting yield surface. This must be compensated by decreasing ξ s , and this must also be done for larger n, in the range between 10 and 50. For aluminum the nonoctahedral slip systems are activated at higher temperatures, for which n < 10. The approach taken here is to adjust ξ s , so that the facet center is touched by the regularized yield surface from the inside. Note, that it is not necessarily the correct choice. Comparison to rate-dependent models should be made, but the results here guides when such an adjustment must be made. Note that adjusting ξ s is mathematically equivalent to adjusting τ s c to τ s c /(ξ s ) n ; hence, this adjustment is vital at elevated temperatures.
For the case of bcc, nonbasal slip systems should also be included at room temperature. The results show that the approximation ξ s = 1 holds for n larger than ≈ 50, which is the case for ferritic steels at room temperature and below. Similarly, in cases with latent hardening, either fcc or bcc, the critical resolved shear stresses will vary for each slip system, and there will be many corners; hence n larger than ≈ 50 is also required in these cases for the validity of the approximation ξ s = 1.
The 5-dimensional yield surface corners are sharp, hence the difference between the regularized yield surface and the rate-independent one is largest at the corner solutions of the Schmid yield surface. Furthermore, the number of corners increases faster than the number of facets in higher dimensions, e.g., the fcc {110} 111 crystal yield surface has 24 facets and 56 corners. Hence the corner solution will dominate a considerable fraction of the yield surface. Figure 8 shows the first threeand two-dimensional sections of the fcc {110} 111 crystal yield surface. Note that what appears as corners and facets in these sections are often projected edges of the surface in the 5-diemensional stress space. The section in Figure 8c is chosen because it contains many true corner solutions, which are marked. The section in Figure 8d contains no corner solutions, but one facet midpoint, as indicated. It can be seen, that at an exponent n = 100, which is realistic at room temperature, the regularized solution is not significantly affected by the regularization. At smaller exponents, the entire yield surface shrinks in Figure 8a,c, since these sections are strongly influenced by their many corner solutions. Still, with ξ s ≡ 1, the yield surface approximately goes through the midpoint of the facet in Figure 8b,d, except for the case with n = 5 in Figure 8b. At this low exponent, which is realistic at hot deformation, the corners significantly influence and shrink all parts of the yield surface, including the facet midpoints. According to the curve in Figure 1, ξ s ≡ 0.86 would be required to make the facet midpoint touch the Schmid yield surface when n = 5.
In the full-constraints Taylor model, for which the velocity-gradient tensor is prescribed, only the corner solutions of the yield surface are realized. Because of stress relaxations due to neighbor-grain interactions that are not accounted for in the Taylor model, some more noncorner solutions may be realized, as seen in recent relaxed-constrains Taylor models e.g., [33,34], similarly as CPFEM with a sufficiently fine mesh [11,35]. However, with any rate-independent model, most of the solutions still will be corner solutions. As pointed out in [10], the corner solutions that are chosen remain the same, as far as the strain rate sensitivity is increased up to less than m ≈ 0.1 (or corresponding lowering of the exponent as low as n ≈ 10). Hence the texture change is not sensitive to the choice of m in this range. However, the magnitude of the stress tensor will decrease rapidly with increasing m. In other words, the Taylor factor decreases sensitively with increasing m.
When calibrating the crystal regularized yield surface to the faceted rate-independent one, the error in the corners remains significant even for very large exponents, n. However, this is not really an error. Rather, interpreted in terms of the rate-dependent potential, it is due to the weak strain-rate sensitivity at room temperature, which still has a significant impact and should be included in the models, i.e., by an appropriate choice of the yield-surface exponent. The instant strain-rate sensitivity of a material can be measured by abrupt strain-rate changes; see [36] for a recent discussion. If assuming the viscoplastic power-law behavior, the application of a rounded yield surface is justified when the exponent n matches the instant strain rate sensitivity, i.e., n = 1 + 1/m.
In applications, large yield-surface exponents may cause severe numerical challenges, since the associated direction of the symmetric part of the plastic velocity gradient will alternate between the adjacent facets of the considered stress corner during return-mapping iterations. Hence it is unrealistic to apply the regularized yield surface as an approximation to truly rate-independent behavior, without carefully dealing with this numerical challenge first. The approach by Scherzinger [16] shows a way to achieve this by adding a line search to the Newton iterations, ongoing work will be reported soon. temperature, the regularized solution is not significantly affected by the regularization. At smaller exponents, the entire yield surface shrinks in Figure 8a,c, since these sections are strongly influenced by their many corner solutions. Still, with ≡ 1, the yield surface approximately goes through the midpoint of the facet in Figure 8b,d, except for the case with = 5 in Figure 8b. At this low exponent, which is realistic at hot deformation, the corners significantly influence and shrink all parts of the yield surface, including the facet midpoints. According to the curve in Figure 1, ≡ 0.86 would be required to make the facet midpoint touch the Schmid yield surface when = 5. In the full-constraints Taylor model, for which the velocity-gradient tensor is prescribed, only the corner solutions of the yield surface are realized. Because of stress relaxations due to neighborgrain interactions that are not accounted for in the Taylor model, some more noncorner solutions may be realized, as seen in recent relaxed-constrains Taylor models e.g., [33,34], similarly as CPFEM with a sufficiently fine mesh [11,35]. However, with any rate-independent model, most of the solutions still will be corner solutions. As pointed out in [10], the corner solutions that are chosen remain the same, as far as the strain rate sensitivity is increased up to less than ≈ 0.1 (or corresponding lowering of the exponent as low as ≈ 10). Hence the texture change is not sensitive to the choice of in this range. However, the magnitude of the stress tensor will decrease rapidly with increasing . In other words, the Taylor factor decreases sensitively with increasing .
When calibrating the crystal regularized yield surface to the faceted rate-independent one, the error in the corners remains significant even for very large exponents, . However, this is not really an error. Rather, interpreted in terms of the rate-dependent potential, it is due to the weak strain-rate sensitivity at room temperature, which still has a significant impact and should be included in the models, i.e., by an appropriate choice of the yield-surface exponent. The instant strain-rate sensitivity of a material can be measured by abrupt strain-rate changes; see [36] for a recent discussion. If assuming the viscoplastic power-law behavior, the application of a rounded yield surface is justified when the exponent n matches the instant strain rate sensitivity, i.e., = 1 + 1/ .
In applications, large yield-surface exponents may cause severe numerical challenges, since the associated direction of the symmetric part of the plastic velocity gradient will alternate between the adjacent facets of the considered stress corner during return-mapping iterations. Hence it is unrealistic to apply the regularized yield surface as an approximation to truly rate-independent behavior, without carefully dealing with this numerical challenge first. The approach by Scherzinger [16] shows a way to achieve this by adding a line search to the Newton iterations, ongoing work will be reported soon.
Deformation twinning by pseudo slip systems or the stress-differential effect are examples of mechanisms that distort the yield surface. In the continuum-plasticity theory, Bauschinger effects are Deformation twinning by pseudo slip systems or the stress-differential effect are examples of mechanisms that distort the yield surface. In the continuum-plasticity theory, Bauschinger effects are commonly handled by a back stress, even in complex models accounting for strain-path changes [37], but recently, significant progress has been reported on continuum models involving yield-surface distortions [21][22][23][24][25][26]. Similarly, crystal-plasticity models with latent hardening predict a shape change of the crystal yield surface during deformation [17,18,20]. The origin of the Bauschinger effect can be a combination of slip-system dependent mechanisms, e.g., dislocation pile-ups, that distorts the yield surface, and composite effects, e.g., strong, large elastic particles, that shifts the yield surface by kinematic hardening. In the proposed framework, both types can be handled simultaneously at the crystal scale, where the yield-surface distortions occur by the critical resolved shear stress changes during simulations. At room temperature, one may apply ξ s ≡ 1, as the exponent n will be larger than 50, i.e., the instant strain-rate sensitivity for most metals at room temperature is m > 0.02. The advanced continuum-plasticity models [21][22][23][24][25][26] can be calibrated based on virtual experiments, performed using models for distortions of the yield surface by virtual experiments calculated by physical-based crystal-plasticity models [17,18,20], which may complement rather demanding mechanical experiments.

Conclusions
Mathematical formulations of yield surfaces for metal crystals have been analyzed. Three calibration methods have been applied for fitting regularized yield surfaces to rate-independent faceted ones. Yield-surface calibrations for the fcc crystal, with or without nonoctahedral slip systems, and for the commonly considered 24 slip systems for bcc steels have been analyzed. An exact match of the midpoints of the facets can be obtained, provided the exponent n is sufficiently high. A good match for the Schmid stress corners, can only be obtained by applying a very high yield-surface exponent. However, the deviation at the corners at lower exponents can be justified by the interpretation of the yield surface exponent in terms of the instant strain-rate sensitivity. This makes the use of a regularized yield surface a precise, flexible method at room temperature without having to deal with the Taylor ambiguity. A return-mapping algorithm that can handle large yield-surface exponents is expected to be available soon, which will make the use of a regularized yield surface a flexible, robust and efficient candidate for CPFEM implementations, enabling the use of realistic physical-based models for the critical resolved shear stresses, allowing a distinction between thermal and a-thermal contributions, and optionally accounting for latent hardening and/or twinning.