A Comprehensive Mathematical Formulation of an Extended Taylor-Bishop-Hill Model Featuring Relaxed Constraints , the Renouard-Wintenberger Theory and a Strain Rate Sensitivity Model

The Taylor theory and the Bishop-Hill theory for the plastic deformation of polycrystals are expressed in a mathematical way which makes extensive use of vectors. These vectors represent either plastic strain rate tensors or deviatoric stress tensors, both in a unified five-dimensional stress-strain space. Such formulation permits a unified formulation of both theories, which can then easily be solved by means of linear programming. The computer implementation of this formalism (in Pascal) conserves this mathematical formalism to a high extent. Relaxed constraints (or "mixed boundary conditions") can very easily be incorporated in the method. The concept of a "relaxed constraint" is formulated in a much more general way than has ever been done before. It is not only shown why there are often multiple solutions for the slip rates, but also that such difficulties can arise for the stress state as well. A few methods for making an appropriate choice among these equivalent solutions are explained, one based on the Renouard-Wintenberger theory that proposes a secondary energy criterion, and another that takes strain rate sensitivity effects into account.


INTRODUCTION
The Taylor theory and the Bishop-Hill theory for the plastic deformation of polycrystals have been outlined several times (Kocks, 1970;Gil Sevillano, Van Houtte and Aernoudt, 1980;Van Houtte and Wagner, 1985). Some papers give details about solution methods such as linear programming (Van Houtte and Aernoudt, 1975;Gil Sevillano et al., 1980). Several extensions or variations have been proposed to this theory since the time that it has first been described and used for the prediction of deformations textures.
An important development has been the introduction of relaxed constraints. Several authors have nearly simultaneously felt that Taylor's first assumption of homogeneous strainstating that each crystallite must undergo exactly the same plastic deformationis too strict and have proposed less strict hypotheses concerning the distribution of plastic strain through the polycrystal Mecking, 1978, 1981;Kocks and Canova, 1981;Van Houtte, 1981; Kocks and Chandra, 1982;Van Houtte, 1982). These so-called "relaxed constraints theories" have been successfully used for the explanation of detailed features of the rolling textures of f.c.c, and b.c.c, metals. The "Hosford" model is a special variant to be used for materials that exhibit the "curling" phenomenon (Hosford, 1963) upon axisymmetric extension or compression (Van Houtte, 1984;Van Houtte, Wenk and Wagner, 1984).
A major problem of the Taylor-Bishop-Hill theory (especially of the "full constraints" version, with the original assumption of homogeneous strain) is, that several solutions are found for the set of slip systems that must be activated in order to achieve the increment of plastic deformation that at a given moment is imposed on a given crystallite. Many discussions have been held about how to choose among those energetically equivalent solutions. In recent times, a few new criteria have been proposed for this selection: to account for strain rate sensitivity: the various solutions may not be energetically equivalent when instead of adopting the same critical resolved shear stresses (CRSS) on all active slip systems, values are used that take the strain rate sensitivity of the CRSS into account (Canova and Kocks, 1984;Asaro and Needleman, 1985).
to account for a second-order term in the minimization of the plastic work: Taylor's second hypothesis is, that the slip system combination must be such, that the internal plastic work required to achieve the imposed strain is minimal. As said above, this still leads to several equivalent solutions. Renouard and Wintenberger (1981) have proposed to require in addition that the rate of variation of this plastic work as strain goes on should be minimized as well. This variation is due to the rotation of the crystal lattice of the grain towards crystal orientations that are more favourably oriented for plastic deformaton, or less. The model proposed by Renouard and Winterberger gives preference to those solutions that let the crystal rotate into more favourable orientations. Several authors have compared the results of this model to experimental data (Skalli, Fortunier, Fillit and Driver, 1985;Bacroix, Jonas, Montheillet and Skalli, 1986).
Several authors (amongst whom the present author) have made computer codes in which all the features mentioned above have been implemented. However, to the knowledge of the author, there is no paper yet that gives a comprehensive view on the mathematical details involved; until now, they had to be collected from several papers by different authors. So the purpose of the present paper is to present the theory underlying these complex computer codes in a modern, comprehensive way. The conventional assumptions and simplifying hypotheses of the Taylor-Bishop-Hill theory and its variants are adopted throughout this paper. They will not be discussed here. For discussions about such aspects, and in general for the physical background of these models, the reader is referred to several review papers (Kocks, Examples: A is a tensor, A its conjugate. A and A T are their matrix representations in a particular reference system (A T is the transposed of A). Without further specification, the crystal system is used as reference system. A superscript indicates that a different reference system is used; the superscript s refers to the sample system: A and AsT.
The elements of these matrices (or the components of the tensors) would then be: A,S. and AS. The transformation laws between A and A are: AS T-AT TTAT (la) and A TAT -1= TAST T (lb) in which T is the transformation matrix (see for example Eq. (1) in the paper by Van Houtte and Wagner (1985)).
The summation convention is used for repeated indices, except when they are put between rounded brackets. The ranges used for the indices are as follows: fori,]: lto3.
for k: 1 to rn + n; rn is the number of possible deformation mechanisms (slip or twinning systems); n is the number of relaxations that are considered. For full constraints models (see below), n is zero. forp, q: lto5. forh:lto2.
for o:: 1 to N; N is the number of equivalent optimal basis solutions of the Taylor theory (see below).
The conventional summation symbol E will however be used whenever other ranges are required than those mentioned above.

Introduction
In the Taylor-Bishop-Hill theory, the crystallites are treated one by one. The theory is normally used to identify the active slip systems of a particular crystallite (of which the crystal orientation is known) at a certain moment during a process of plastic deformation. The strain history of that process is assumed to be known in terms of macroscopic strains or strain rates. The identification of the active slip systems is accompanied by the calculation of the local plastic stress state, of the local rate of plastic work, of the slip rates on the individual slip systems and finally of the resulting rotation rate of the crystal orientation.
The basic assumptions of the Taylor theory are: 1) All crystallites are subject to the same plastic strain (= the prescribed strain).
2) The plastic strain is achieved by (multiple) slip.
3) Of all possible combinations of slips that achieve the prescribed strain, those that minimize the internal work must be chosen.
On the other hand, the basic assumptions of the Bishop-Hill theory are: 1) All crystallites are subject to the same plastic strain (= the prescribed strain).
2) The plastic strain is achieved by (multiple) slip.
3) The stress state must satisfy the yield locus of the crystal, i.e. for no slip systems, whether active or not, should the resolved shear stress exceed the momentaneaous critical shear stress of the system.
4) The Maximum Work Principle (Hill, 1950) is used to derive the stress state from the prescribed strain. This means that of all stress states permitted by the yield locus, the one that maximizes the external plastic work (as calculated from the stress tensor and the prescribed strain tensor) is the right one. Bishop and Hill (1951b) and several other authors (see the review papers by Gil Sevillano et al. (1980)) demonstrated that the two methods are strictly equivalent. This theory is now called the "full constraints" (FC) theory, in contrast to the so-called "relaxed constraints" (RC) theories. The latter somewhat "relax" the first hypothesis of both the Taylor and the Bishop-Hill theories. As a result, it is typically replaced by a hypothesis of the type (example for rolling): "All the crystallites are subject to the same plastic strain, apart from a shear component K3 which is free." In this example, g]3 stands for a shear in which the rolling plane is the shear plane and the rolling direction is the shear direction (lath model, see Van Houtte (1981)).
More details about these assumptions can be found in several review papers (Gil Sevillano et al. (1980), Van Houtte and Wagner (1985) for FC models, Van Houtte (1984) for RC models). In this paper, focus will be on the mathematical aspects. It is nevertheless worthwhile to mention that a rigid-plastic approach is used here: elastic strains are neglected; all strains or strain rates that are mentioned are plastic.

Full constraints models
In this section, n is zero.
The prescribed strain which is mentioned in the first assumption (see above) has often been described by a prescribed displacement gradient, since not only the strain tensors, but also the displacement gradient tensors that describe the deformation of each crystallite must be equal to each other in order to guarantee the geometrical integrity of the deforming polycrystalline aggregate.
Since it is not the intention here to use a finite strain theory for the description of stresses and strains, it is better to use momentaneous slip rates, strain rates etc. instead of finite slips or strains. For this reason, the prescribed plastic strain will be thought to be described by a displacement rate gradient: prescribed displacement rate gradient Ko/eq (2) in which/eq represents the macroscopic equivalent strain rate. Ko is a (not necessarily symmetric) tensor that characterizes the prescribed strain mode. In the case of rolling, its matrix representation could be (Xl is the rolling direction, x3 the sheet plane normal): In the case of simple shear (as is the case at the surface of the sample in a torsion test), it would be (Xl is the shear direction and x2 the shear plane normal): The expressions of these prescribed tensors K0 are known in the sample reference system. Equation (lb) shows how they can be converted to the crystal's system. Let /k be the slip rate on slip system k. Such slip rate causes a deformation of the crystal (a s!mple shear) expressed by the displacement rate gradient K(k)F(k). In a reference system g, associated to the glide system k, and for which x is the glide direction and x the slip plane normal, the components of the tensor K are given by: (K)ij 0 except (K)12 1 (5a) In the crystal system, the components of the tensor are given by (see Gil Sevillano et al. (1980) or Van Houtte and Wagner (1985)): (Kk)i r(kiv(k (5b) in which I' k is a unit vector in the glide direction of slip system k, and Vk is the unit vector normal to the slip plane. The same formalism can be used for a mechanical twinning system, for which we also refer to the review papers mentioned above. The condition that the strain caused by all slip systems together must be equal to the prescribed strain, now becomes: Kk/k q-RJeq-" Koeq (6) The tensor R represents the unknown rate of rotation of the crystal lattice per unit macroscopic strain rate.
With an eye on the com.puter implementation, it is advantageous to rep!ace the slip rates Fk (unknown) and the macroscopic strain rates Eeq (known) by the variables gk: gk -'k //eq (7) Equation ( (15) shows how the crystal's freedom (represented by the "free" variables gk, k ranging from rn + 1 to rn + n) is introduced in K. Moreover, in the case of relaxed constraints, K0 does not necessarily describe the average deformation of the polycrystal, as it does for the full constraints models. The average deformation of the polycrystal would be described by the weighted average of I, which depends on the sample symmetry of the texture at that moment. Let us illustrate some of these ideas for the example of the lath model for rolling. The shear that corresponds to (K))13 is relaxed. In practice, the component (K)13 of K0 is still given the value zero, as would be done in the full constraints case; however, it could be given any other value as well, since it will be "overruled" by the term Km+l gm+l anyway (see Eq. (14a) and Eq. (15)), in which gm+l is a "free" variable. The 13-shear component of the polycrystal as a whole will be zero on condition that at that moment, the texture has an orthorombic sample symmetry (i.e., the rolling direction, the rolling plane normal and the transverse direction are twofold axes of rotation). When present, this symmetry will be conserved, even when the texture evolution is simulated by a relaxed constraints model.
The two degrees of freedom that must be introduced in order to obtain the Hosford model (the "curling" model described by Van Houtte (1984) and Van Houtte, Wenk and Wagner (1984)) that gives the grains the kind of freedom that permits them to develop curly microstructures (Hosford, 1963) (11)(12)(13), neither to the way of solving them. The tensor K0 that describes the prescribed strain in the FC model can formally still be treated as a constant. The only difference of some practical importance is, that for the slip systems, the matrix representations of the K-tensors of the slip systems are known in the crystal's reference system (Eq. (5)), while on the other hand those of the relaxations are usually known in the sample's reference system. Equation (lb) can of course be used for the conversion.
Generally speaking, all that follows applies as well to FC models as to RC models. Each time g,, K,, A, etc. are mentioned, they may be associated as well to a true slip system, to a pseudo-slip system that simulates a mechanical twinning system (see e.g. Van Houtte and Wagner (1985)), or to pseudo-slip systems associated to a relaxed constraint.

Vector representation of stress and strain rate tensors
Suppose. that the symmetrical tensor I) is either a plastic strain rate tensor E, a tensor such as A in the equations above (and hence proportional to a plastic strain tensor), or a deviatoric stress tensor S. Then the following property holds for D: So i) has only five independent components. It is convenient then to represent such tensors as vectors in a five dimensional stressstrain space. Let us define the vector d that represents the tensor D as follows: Note that these definition make use of the components of D in the crystal's reference system. A different vector would be obtained when D were expressed in another reference system, for example the sample's system. The Eqs. (18)(19)(20)(21)(22) are chosen in such a way that the following expression holds: This for example means that the rate of plastic work per unit volume is equal to the scalar product of the vectors/ and s, which are the vector representations in stress-strain space of the tensors ! and S. Lequeu, Gilormini, Montheillet, Bacroix and Jonas (1987) have derived conversion formulas which are similar to Eqs. (18-22) (but different) and for which Eq. (26) also holds. It is also interesting to note that the equivalent stresses and strain rates according to von Mises are directly proportional to the lengths of the vectors/ and s: ]'VM /'2/" ,1,2 (, )1,2  (29) in which the right-hand side represents the prescribed strain.
Taylor's criterion of minimum internal work can be expressed as follows: in which the factor eq can be dropped since it is constant with respect to the minimization. For a full constraints model, W} is (after minimization) equal to W, the rate of plastic work per unit volume and per unit macroscopic strain rate Eeq (see Eq. (7)). ry, represents the critical resolved shear stress on slip system k. It is usually assumed to be zero for the pseudo-systems that correspond to the relaxed constraints. This expression has two disadvantages: it contains a non-linear element, and it does not allow different values for the critical resolved shear stress in positive direction and in negative directions. Both difficulties can be resolved by substituting the variables gk by ghe as follows: if gk -> 0 then gle ifgk<0 then g lk Note that and g2k 0 and g2k----gk (31) ghk >-0 (32) f W can also be called "the plastic work per unit volume and per unit macroscopic strain." Inversely, g g g2k (33) Now let the absolute values of the critical resolved shear stresses on slip system k be tlk in positive direction, and tzk in negative direction Equation (30) now becomes: W thkghk Min (34) A mechanical twinning system can be simulated by treating it as a "pseudo-slip system" (see e.g. Van Houtte and Wagner (1985)).
An undesired "pseudo-slip" in the antitwinning direction can be avoided by assigning a prohibitively high value to tzk.
The vectors ak are now replaced by vectors ahk: alk Ilk and azk --ak

Vector expression of the Bishop-Hill equations
Let S be the deviatoric stress tensor and s its vector representation. The corresponding resolved shear stress on slip system k is According to Eq. (5a) and Eq. (9), the tensor Ak has the following components in the glide system's reference system (see section (3.2)): (A)i 0 except (AD12 (A)21 1 / 2 (37) Hence (S12 + S,) sg2 (38) $" Ak S.(A)ij g which by definition is the resolved shear stress rk. Hence, according The classical way of expressing the Bishop-Hill equations is: together with the condition, that the external plastic work rate per unit volume must be maximal: /eqW3H Sij(Ao)ijeq Max (41) which is derived from the Maximum Work Principle of mathematical plasticity theory (Hill, 1950). In case of full constraints, WH is the rate of plastic work per unit volume and per unit macroscopic strain rate, or alternatively, the plastic work per unit volume and per unit macroscopic strain. Using Eq. (26), and dropping/eq since it is a constant with regard to the maximization, this becomes: Now using Eq. (39), and substituting r by thk and ak by ahk (Eq. (35)), Eq. (40) becomes: S ahk <--thk (43) Equations (42) and (43) (1969)) or used by the author in previous work (Van Houtte, 1982), these solutions can be subdivided in several classes" -general solutions; -valid solutions, for which the non-negativity conditions Eq. (32) hold; -so-called "basis solutions", for which no more than five ghk are different from zero. A basis solution can be valid or non-valid; -optimal basis solutions" these are valid basis solutions for which the minimum internal work condition (34) is satisfied.
The Taylor Eqs. (32), (34) and (36) form a set of equations that can be solved very efficiently by means of linear programming (see for example Chin and Mammel (1967), Lister (1974), Van Houtte and Aernoudt (1975)). The best known technique is the simplex algorithm. It consists of: -finding a first valid basis solution; -then performing an optimality test on it; when the solution is not optimal, a different valid basis solution is chosen using a special rule. The optimality test is then repeated, etc.
This algorithm makes extensive use of a special concept: the basis.

Basis in stress-strain space
The "basis" of the five dimensional stress-strain space is an important concept in this context. Any set of five independent vectors Ip can be used as a basis. It is then possible to express an arbitrary vector such as d as a linear combination of the basis: d dlp (44) in which the superscript b of the coefficients dp b refers to the particular basis Ip.
Let d be a column-matrix that consists of the components dp of the vector d, and d one that contains the coefficients dpb. The components of the basis vectors bp are bpq. The matrix U is called "the inverse of the basis"" Fbll b21 bs; - In what follows, it will prove to be very useful to choose five independent vectors among the ahk of the Taylor Eqs. (36) or the  in order to constitute a basis" bp ahp,p The indices hp and kp themselves depend on p (which ranges from 1 to 5), according to the choice that one wants to make. The analysis below will make it clear, that the Taylor theory has no solution if less than five independent vectors bp can be chosen among the ahk, except perhaps for some very special combinations of crystal orientations and prescribed strains.
In order to simplify the notations, the unknowns ghk that are associated with the abe vectors that have been chosen for the basis will be called yp: Ye ghe (48) Similarly, the coefficients th Of the internal work function that "belong" to the basis will be called u: up thke (49) 5.3 Finding a first valid basis solution The first thing to do is to find five independent vectors ahk that can constitute a basis, as has been explained in the previous section. This could for example be done by choosing at good luck five vectors bp among the ahk. They are independent if it proves possible to calculate the matrix U (Eq. (45)). If U does not exist, another choice must be made among the ahk etc. A more systematic approach is also possible, for example by using the first stage of the so-called simplex method (to be distinguished from the simplex algorithm, see for example Gass (1969) So a basis solution has been found. It is however not necessarily a valid solution, since a valid solution must satisfy Eq. (32), or: y _> 0 (52) Suppose that one of the yp is negative. A very easy remedy then is, to change the basis by replacing the corresponding ahk vector by the vector ah,k, with h'=3-h (53) This simply means a sign reversal of ahk. As a result, yp of the new basis solution will also undergo a sign reversal, and become positive.

Optimality criterion
A valid basis solution is not necessarily optimal, i.e. it does not necessarily minimize the internal work (Eq. (34)). Suppose now that a valid basis solution has been found. In order to be able to check optimality, a new concept is needed: the simplex multiplier.
Let the simplex multiplier associated to a particular basis be a vector n, which will be defined below. Whatever its value, it may be multiplied with the right-hand and the left-hand side of Eq. (36), after some manipulations leading to: -ahkghk q" , 10 0 (54) So the left-hand side of expression (54) may be added to the right-hand side of expression (34), which then becomes: W (thk n ahk)ghk qn" aO (55) Now for the present basis, the following "cost function coefficients" are defined: tk thk n ahk (56) Equation (55) now becomes: The tk that belong to the basis are given a shorter notation, similar as has been given to thk (Eq. (49)): " " Suppose now that one gradually moves on to another solution, respecting the Taylor equations. At the starting point, which is considered a valid basis solution, only the ghk that are associated to the basis can be different from zero; but they cannot contribute to W since their coefficients are zero (Eq. (60)). So their variation can exercise no influence on W.. The value of the latter will change because of the variation of the ghk outside the basis. Since these are initially zero, and can only evolve in positive direction because of Eq. (32), they will increase W when their coefficients are positive, and decrease it when their coefficients are negative. Because of Eq. (63), the coefficients are positive; so any departure from the considered valid basis solution increases W. Hence the minimum has been found.
Since all ghk outside the basis are zero, it follows from Eq. (57) that the value of this minimal internal plastic work rate is given by: W: n. a0 (64)

Changing the basis
Suppose that the optimality test led to the conclusion that the considered valid basis solution was not optimal. The following procedure then allows to go on to a new valid basis solution that has more chances to be optimal.
The principle is, that one of the vectors of the basis will be replaced by one that originally was outside the basis.
The hk indices of the vector ahk that must be brought in the basis are found from the following condition" t,k Min (65) Since the solution was not optimal and Eq. (63) was not satisfied, this minimum is negative.
Let us call this vector ahk and its associated variable gh, simply a and g. The more g can be increased from zero, the smaller W will become. Equation (50) gave the value of the non-optimal valid basis solution. In the next solution, the variable g will make its introduction, so that Eq. (50) must be adapted: bpy v + ag (ao)pbp (66) Let us express a as a function of the basis, using Eq. (44). Equation (66) now becomes: Ip(yp + abpg) (ao)bp (67) In the new solution, g will increase from zero to some positive value and the old basis variables yp will change accordingly. Their new values can be derived from Eq. (67): yp (ao)bp-abpg The fact that the yp should not become negative, sets an upper limit b is positive. Because of this, the new value of g will be" to g when ap g Min[ aP b j for ap Condition (68) also gives the index p of the variable yp that will become zero. It can be removed from the basis since it will become zero anyway. So the vector a will replace the vector bp in the basis.
When the new basis has been identified in this way, the new valid basis solution can be calculated, after which the optimality test can be repeated. After each change of the basis, a lot of coefficients must be updated" U, n, (ao)p b, (ahk)bp and tk. It is not necessary to carry out each time comparatively computer-time consuming calculations such as the evaluation of U by Eq. (45). Simple and efficient updating rules can be found in the specialized literature (e.g., Gass, 1969).

The Bishop-Hill stress
Once the simplex algorithm has found an optimal solution of the equations of the Taylor theory, the conditions (63) are satisfied. It is then seen from Eqs. (56) and (63) The external plastic work rate (per unit macroscopic strain rate and per unit of volume) is then given by s* a0, in analogy to Eq.  (49)) leads to the conclusion that the resolved shear stress is equal to the critical resolved shear stress on the slip systems that correspond to the basis. Hence in so far the stress state is concerned, these may be active. This corresponds to the solution that the Taylor theory finds for the ghk: these are always zero outside the basis, not necessarily inside the basis.
Relaxed constraints and stress conditions. In section 4.2 it was stated, that the th which are associated to pseudo-slip systems that represent relaxed constraints, are usually assumed to be zero. A more general assumption could be: This is an stress condition that is implicitely imposed when adopting a condition such as Eq. (74) in the Taylor theory. Indeed, this equation can be transformed into (using Eqs. (9) and (26), and realizing that the stress tensor S is symmetric): S:K, S c (76) When for example K, is defined by Eq. (14a), this leads to: S]3 S c (77a) i.e., the component 13 of the stress (in the sample's system) is imposed. The physical meaning of S c is obvious in this example; it is not always so clear. In general, a stress condition is imposed for every constraint that is relaxed" freedom of strain and stress are complementary. An interesting case is the "relaxed constraints" defined by Eq. (16a-b) ("Hosford model" for curling). Applications of Eq. (76) to Eq. (16a) leads to (with SC= 0)" S1 Sz2 (77b) On Eq. (16b) ( shear stress, in addition to those of the basis systems. As a result, it is quite probable that some of the corresponding ahk vectors may be brought in the basis (by exchanging them with some basis vectors bp) without affecting the optimality condition. Other solutions for the slips which are also optimal may be found in this way. Figure lb tries to give a graphical representation of this situation. The figure shows a part of stress-strain space, which for reasons of clarity is assumed to be three-dimensional. The stress conditions Eqs. (43) are represented by planes: the points that represent the plastic stress must stay on one side of them. Since the space of the figure is three-dimensional, three planes normally intersect in one point (Figure la). According to the Maximum Work criterion (Eq. (42)), the vector a0 which represents the imposed strain, must be contained within the cone of normals on the three planes (Bishop and Hill, 1951a). So the "corners" of the figure represent the plastic stress states (except in special cases when a0 happens to be normal to one of the planes or to an intersection of two planes). The three planes that meet in such corner indicate the active slip systems. The normals on them constitute the basis that leads to the optimal basis (al (b) Figttre 1 (a) In a three-dimensional stress-strain space, a "corner" of the single crystal yield surface would normally be the intersection of three planes. Each plane generates a "facet" and represents a particular slip system. The vector ao lies within the cone of normals on the three facets that intersect in the stress corner which satisfies the Maximum Work Principle. The active slip systems are identified by these three facets. (b) "By accident", the stress corner can also be contained by a fourth (or a fifth...) plane. Four slip systems may now be active, although three would be sufficient in a three-dimensional stress-strain space.
solution (in the space of Figure la, a basis consists of three vectors only). Suppose now that a fourth plane representing a slip system also passes through the point (Fig. lb). This would mean that the resolved shear stress on it also equals the critical shear stress. There are now four normals that can be used as basis vectors. This means that four different basises each consisting of three vectors can be constructed. Each of these basises leads to a valid optimal basis solution.
In true stress-strain space, a "corner" is the intersection of at least five hyperplanes representing slip systems. In the case of {111}(110) slip (f.c.c. metals), and for a FC model, up to eight hyperplanes intersect in most corners when the same critical resolved shear stress is assumed for all slip systems. It is still a subject of debate whether this is a good model assumption or not. It nevertheless remains necessary to find all optimal solutions when performing Taylor-Bishop-Hill calculations.
The procedure of finding all optimal solutions is based on the "change of basis" procedure of section 5.5. Each time an optimal basis solution is found, the thbg that do not belong to the basis but are nevertheless zero are sought, ahk vectors that belong to such zero tk can be brought into the basis without need to change , since condition (60) remains satisfied anyway. The procedure of section 5.5 is used in order to check which vectors could be removed from the basis in order to make room for the new vector.
The possibles changes of the basis that are detected in this way are not performed right away, but the identifications of these alternative basises (sets of 5 hk indices) are stored in a "linked list" in computer memory for later use. When all possible basises are stored in this way, an algorithm starts that examines if the linked list still contains basises for which the solution has not yet been computed and stored. If it finds one, the new optimal solution is computed, the resulting slips are stored, and the procedure described above is repeated. The algorithm stops if the linked list does not contain any more basises for which the solution has not yet been stored. For many applications, it is necessary to collect all those possible optimal solutions. Several methods based on physical arguments have been proposed to make a selection among them. But before some of these methods can be discussed, attention must be drawn to another type of non-uniqueness that may arise: degeneracy of the Taylor equations, leading to a non-unique Bishop-Hill stress. 6.2 Degeneracy of the Taylor Equation--Non-uniqueness of the Bishop-Hill stress Suppose that the vector a0 that describes the prescribed strain, happens to be perpendicular to a facet or an edge of the single crystal's yield locus (Fig. 2). Less than five slip system can in such case accommodate the prescribed strain. All the points of that facet or edge then represent stress states that are valid solutions: they do not violate the yield locus, and they satisfy the Maximum Work condition. In practice the set of possible solutions is fully described by the set of corner stresses that belong to the considered facet or edge.
In terms of the linear programming method for the solution of the Taylor theory, this situation is encountered when one or several of the components (a0)p (of the vector a0 when expressed in the basis associated to an optimal basis solution) happen to be zero. In such case the linear programming problem is said to be 'degenerated'. In the limiting case, when all but one of these components are zero, only one slip system will be activated. This can be directly Figure 2 In case that ao happens to be perpendicular to a facet of the single crystal yield locus, only one slip system must be active in order to accommodate the imposed strain. Several stress states may however achieve this: in the case of this figure, not only the corner stresses A, B, C, D, E or F but in fact all stress states that are represented by the points of the facet ABCDEF. concluded from Eq. (51); the vector a0 will be perpendicular to the facet that represents the active slip system. Any point on that facet, including the edges and corners that surround it, represents a valid solution for the plastic stress. The set of all possible solutions for the stress can be described by the set of corner stresses that surrounds the facet. Such corner stress is identified by a basis. In such case, it is possible to change the basis without changing the solution that the Taylor theory finds for the slips, i.e. the same slip system remains active and all others inactive, even when they belong to the basis. In a sense this is a inverse type of nonuniqueness as usually encountered in the Taylor theory (section 6.1): here the solution for the slips is unique but the solution for the stress is not.
In order to treat this problem, a procedure will be proposed that allows to find another equivalent basis to which a different stress state is associated. The procedure is intended to be used when an optimal basis solution is found for which at least one of the 5 slips Yv is zero. It is similar to the procedure described in section 5.5. Much in the same way as for the other non-uniqueness problem (section 6.1), its purpose is to collect all possible optimal basis solutions that correspond to different stress states.
Before it can be explained, the vectors c, (associated to a particular basis) must be introduced. The components Cpq of these vectors constitute by definition the rows of the matrix U defined by Eq. (45). Since UU -1= I (84) in which I is a 5 5 unit matrix, and since the components of the vectors b, constitute the columns of the matrix U-1, this means that lp" I}q q (85) in which (pq is the Kronecker symbol. The vectors Cp also constitute a "basis" that can be used to express vectors in stress-strain space, as explained in section 5.2. Let us now assume that an optimal solution of the Taylor theory (FC or RC) has been found, and that one of the slips is zero (Eq. (51)): Yb' (a0)p b, 0 Because of this, the decomposition of a0 according to the basis remains unchanged when the basis vector bp, is replaced by another one (which of course must also be independent of the remaining basis vectors). So the slips are not affected by this operation; neither is W, since it is given by Eq. (34). The Bishop-Hill stress s however changes. Assume that it is augmented by s', so the new stress is given by s + s'.
It will first be explained how acceptable values for s' can be calculated without changing the basis, since this is a costly operation in terms of computer time. The Bishop-Hill equations will be used.
Let us choose s' in the following way: Since both left-and right hand side of this condition are negative, the most negative value that can be chosen for s' is given by: All values of hk would typically be tried out in order to find this maximum. The hk-indices that are associated to it designate the vector ahk that will replace the vector b,, in the basis. A new basis that leads to the same solutions for the slips and for Ww and but to a different stress state, has now been identified. Once this basis is adopted, still other solutions may be detected. The set of all possible solutions can be obtained by a procedure very similar to the one described in section 6.1. Most generally, all optimal solutions of the Taylor theory (section 6.1) and of the Bishop-Hill theory (this section) must be collected before a selection based on physical arguments can be made among them.

Dealing with multiple optimal solutions
Suppose that N different basis solutions of the Taylor theory have been found. The particular optimal basis solution gives the solution (g) for the slip rates per unit macroscopic strain. It can easily be demonstrated, that any positive linear combination of these slip rates per unit macroscopic strain also leads to a valid, optimal but non-basis solution" g wo ( Eqs. (94-96) define a "cloud" of possible solutions.
One way of dealing with the problem is to use the point of gravity of this cloud, i.e. to choose the average of all possible solutions in order to calculate the rotations.
Another possibility is to make a random choice among the various optimal basis solutions. It is so avoided to introduce solutions that activate more slip systems than necessary. This method only makes sense when the calculation is a part of the simulation of the plastic deformation of a whole polycrystal, i.e. when many crystallites are involved. If different increments of strain have to be simulated, then the risk exists that for a given grain, the set of active slip systems changes after each strain increment. This can be avoided by giving to the slip systems that have been active in the previous step a slightly smaller critical shear stress than to those that were not active.
A variant of the previous method would be to avoid giving the same critical shear stresses to all slip systems. This can be done by using random numbers for the CRSS on the particular slip systems. These numbers must have a narrow Gaussian distribution around an average CRSS. As a result, only one optimal solution will be found in general. Chin (1969) has proposed to use the amount of cross slip or the amount of coplanar slip as selection criteria in order to choose among the optimal basis solutions of the Taylor theory. Such method would permit to take certain physical parameters (such as the stacking fault energy of the alloy) into account. Renouard and Winterberger (1981) have proposed to minimize not only W, but also the rate of evolution of it with macroscopic strain. This includes the hardening of the slip systems, but also the change of W due to the fact that the crystal rotates into orientations that are more favourable for plastic deformation, or less. Since the lattice rotation depends on the slips, it will be different for each optimal basis solution of the Taylor theory. One method of taking this effect into account will be explained in the next section.
-Several authors (Asaro and Needleman, 1985;Canova and Kocks, 1984) advocate to take the strain rate sensitivity of the critical shear stresses into account. It will be shown below how this can be done within the framework of the Taylor theory.
6.4 Implementation of the Renouard-Wintenberger method Renouard and Wintenberger (1981) have proposed to minimize not only W but also the rate of change of it with Eq, which then becomes a second order selection criterion.
It is proposed here to use such criterion only after a complete determination of all valid optimal basis solutions, using the principles explained in section 6.1 and section 6.2, has been carried out. Only such procedure avoids the generation of solutions that are not optimal with respect to the first order selection criterion (Eq. (34) or Eq. (42)).
Let us now consider a particular optimal basis solution of the Taylor theory. All slip rates that do not belong to the basis are zero; those that belong to the basis are the five yp; the corresponding critical resolved shear stresses are Up. From Eqs. (34), (48) (99) reduces to a set of 5 linear equations with 5 unknowns. Now even as the crystal rotates, the vectors Ip remain unchanged, as they are the vector representation of the strain caused by slip on a certain slip system as expressed in the crystal's reference system. The vector ao however is not constant as the crystal rotates, since it is the vector representation of the prescribed strain tensor (fixed in the sample's reference system) expressed in the crystal's system. So the derivation of Eq. (99) to Eeq leads to: dyp dao (100) Ip dEeq dEeq This is a system of 5 linear equations with 5 unknowns: the derivatives of the slip rates to the macroscopic strain. It can easily be inverted: the inverse of its coefficient matrix is the matrix U (Eq. (45)) which is directly available, since U is usually updated each time a new basis is chosen (see section (5.5)). The right hand side of the equation must of course be calculated before the inversion is possible. The vector a0 is obtained from A0, the matrix representation of A0 as expressed in the crystal system. Equation (lb) shows how A0 is obtained from A, which is constant when the crystal rotates. The columns of the transposed of the transformation matrix T can be regarded as vectors in physical space, representing the unit vectors on the reference axes of the crystal (Gil Sevillano et al., 1980). It follows from this that T T + dT T (I + R dEeq)T T (101) in which I is the unit matrix, and R the matrix representation in the sample's reference system of the antisymmetrical tensor R of the lattice rotation (defined by Eq. (6)). R can be calculated for each optimal basis solution by means of Eq. (8). It follows from Eq. can be calculated. The optimal basis solution for which this value is minimal is preferred to the other solutions. Only in rare cases, several optimal basis solutions are found that are also equivalent with respect to this second order criterion.
6.5 Strain rate sensitivity A method to use strain rate sensitivity in order to define a criterion for selection amongst the various optimal basis solutions of the Taylor theory will be proposed here. It should however not be forgotten that these optimal solutions have been obtained while neglecting strain rate sensitivity. So the model that is presented here is not a consistent strain rate sensitivity model as the one presented by Asaro and Needleman (1985) in which all slip systems are considered as active, not only five.
A variant of the classical power law for the strain rate dependency of the critical resolved shear stress would be: :' : oo (105) in which r and/o are constants (the first with the dimensions of a stress, the second with those of a slip rate). / represents the strain rate exponent. The conventional symbol rn has not been used because it already has an other meaning in this paper. A typical value for it is 0.01 for room temperature deformation.
In this paper, the values of thk will be used for r, which then is particularized for each slip system and for positive and negative slips. Let us furthermore define go =/0//eq (106) which is also a constant. Equation (105) can now be written in terms of ghk and thk (see Eq. (7)): "Chk t(hk) go Let us now define Ws*, "the rate of plastic work per unit volume and per unit macroscopic strain rate, according to the strain rate sensitivity model." It is associated to a particular optimal basis solution of the Taylor theory, and it is not equal to W or W: w: ( ao8) k=l using Eq. (107). This quantity is calculated for each optimal basis solution of the Taylor theory; the one that leads to the smallest value is preferred. Such procedure usually eliminates the nonuniqueness problem.
The Renouard-Wintenberger model (section 6.3) and the strain rate sensitivity model as presented here exclude each other, since they both provide in a criterion for choosing the "most appropriate" optimal solution for the slips.

CONCLUDING REMARKS
The concept of a "relaxed constraint" has been defined here in a more general way then has ever been done before.
The most general formulation until now has been the mixed boundary conditions as defined by Renouard and Wintenberger (1976) for the problem of the plastic deformation of a single crystal.
In short, they state that 5 boundary conditions are needed, of which there are n stress conditions and 5-n strain conditions. The stress conditions are of the following nature: the value of n components of the stress tensor as expressed in the external reference system (= the sample system in this context) are prescribed.-The nature of the strain conditions is as follows: the value of the 5n components of the strain tensor (in the external reference system) that correspond the non-prescribed stress components is prescribed.
In the terminology of the present paper, the n stress conditions are "relaxed constraints" since they correspond to non-prescribed strain components.
Although quite general, Renouard and Wintenberger's stress conditions (=relaxed constraints) are limited to the individual components of the stress and strain tensors as expressed in the sample's system. This type of conditions can be achieved by "relaxed constraints" of the concept defined in this paper. A stress condition for S3 can for example be achieved by introducing a relaxation for the component (K)13 Of the prescribed displacement rate gradient tensor (per unit macroscopic strain rate). See Eq. (14a) in section 3.3 for this. The new variable gk that is introduced in this way gives the component (A)3 of the prescribed strain rate (per unit macroscopic strain rate) its freedom back. This can be seen by taking the symmetric part of Eq. (15). The corresponding stress condition is obtained by applying Eq. (76) to the Kk tensor defined by Eq. (14a), leading to a Renouard and Wintenberger-type stress condition Eq. (77a). Note that the same result would have been achieved by a relaxation of the type (K))31; in the theory of Renouard and Wintenberger (1976), no distinction between these two cases is made.
The concept of a relaxed constraint as defined here (section 3.3) is by no means restricted to particular components of the stress or the strain tensors. Equation (16a) shows a more complex example; Eq. (77b) gives the corresponding stress condition. Nor is it restricted to tensors that are expressed in the sample's system. One can in principle introduce the relaxed constraints by defining a set of n tensors K with associated pseudo-slip systems (section 3.3) that " It would have been better to prescribe n independent components of the deviatoric stress tensor instead.
are subject to no restrictions at all, except that they must conserve the volume and that they better be independent of each other. They can for example represent a simple shear along the boundary between a cementite lammella and a ferrite crystal. Such boundary does in general not have a simple orientation with respect to either the crystal system or the sample system. The stress conditions that correspond to such general strain relaxations are given by Eq. (76).
An entirely new computer code has been developed in which the mathematical formalism presented in this paper has been implemented. The program has been written in Pascal, which has the advantage that the code is very readible" it reflects the mathematical formalism to a high extent by making extensive use of dedicated data types. Other advantages when compared with the old Fortran code are: greater flexibility when introducing for example strain hardening matrices or other special models; the use of pointer variables allowed an easier and more reliable processing of the "linked lists" mentioned in section 6.1 and 6.2 (non-uniqueness problem). For the rest of it, the new Pascal program approximately conserves the same speed as the previous Fortran version. Because of the use of the highly efficient linear programming technique, this speed is equivalent to that of programs based on the Bishop-Hill method (see e.g. Gil Sevillano et al. (1980) for details on this method). The latter have the disadvantage that they do not allow variations of the CRSS-ratios after each deformation step.
To the author's knowledge, the treatment of the non-uniqueness of the Bishop-Hill stress has been implemented for the first time in a linear-programming type of Taylor code. An interesting application would be the fully automatic and exhaustive generation all Bishop-Hill vertices for a material with any set of slip systems and any set of CRSS ratios. It would indeed be sufficient to use the program for finding the slips and the plastic stress for a single crystal, consecutively specifying a simple shear on each of the slip systems as imposed strain. For each case, the program would automatically generate the stress states surrounding the facet of the Bishop-Hill yield locus that correspond to the chosen slip system. After having done this for all slip systems, a list in which all vertices are present several times is obtained. This list must then finally be sorted by a sort program and the doubles removed.