Twinning-induced anisotropy of mechanical response of AZ31B extruded rods

Texture and twinning-induced anisotropy of the yield stress and hardening of AZ31B extruded rods is investigated. The multidirectional compression tests involving strain path changes are performed in order to: i. assess which slip and twinning systems are active in the polycrystalline sample with a strong texture, ii. analyze the influence of the preliminary deformation upon twin formation, iii. observe the resulting change of the mechanical response. In order to fulfil these goals mechanical testing is supplemented by microstructure analysis. Experimental observations are used to validate the proposed crystal plasticity framework when it is combined with the viscoplastic self-consistent scheme. On the other hand, the results of numerical simulations are used to confirm an advocated interpretation of experimental findings. Finally, the experimental and numerical results are discussed with respect to the theoretical study of slip and twinning activity on the basis of the generalized Schmid criterion. It is concluded that twinning activity influences the mechanical response predominantly by the texture change and to lesser extent by modification of strain hardening due to slip-twin interactions.


Introduction
In the case of hcp materials such as magnesium, zirconium or Ti alloys, twinning is an important mechanism of plastic deformation at room temperature, cf [1][2][3][4][5][6]. For these materials of high specific strength, twinning enables to compensate for the lack of five independent slip systems required for a general shape change of a single crystal. Activation of twinning is therefore, although to a limited extent, salutary for reduction of the excessive stress level and enhancement of ductility. On the other hand, twinning has a marked effect on texture evolution and hardening response of these metals and alloys, which is especially apparent when the strain path is changed. The limited number of easy slip systems results in the development of strong crystallographic textures upon mechanical processing, e.g. rolling, cf [1,[7][8][9]. or extrusion, cf [10]. It is particularly true for magnesium alloys such as AZ31B, for which the family of easy basal slip provides only two independent systems per grain. The pronounced texture is a source of strong anisotropy of mechanical properties of sheets or extruded rods.
As recognized in the literature [11] there are three main sources of twin-induced hardening in metals and alloys: i. appearance of difficult-to-deform orientations (textural hardening) due to the abrupt crystal reorientation [10,12]. This source of hardening was indicated as a most decisive one for magnesium [12,13]; ii. reduction of a mean free path for dislocations (Hall-Petch-like effect) due to the appearance of multiple matrix-twin boundaries [2,3,14]. This mechanism seems to be more important in fcc materials of low stacking fault energy, where twins cluster to form platelets a few micrometers thick with stable thin layers of matrix between them, than for magnesium in which twinning proceeds by the growth of initial platelets to encompass finally the whole parent grain [13,15]; iii. conversion of glissile dislocations into the sessile ones due to the twinning shear transformation of the lattice (Basinski's-type hardening [16]). This mechanism for AZ31B was investigated in Ref. [17].
Competing families of easy and hard to initiate slip systems and their mutual interactions with twinning as well as abruptly reorienting parts of crystallites constitute a challenge in developing successful crystal plasticity modelling frameworks. Particularly, in order to provide meaningful predictions of three hardening mechanisms indicated above the framework must involve three components: i. a reorientation scheme to account correctly for the appearance of twin-related orientations in the representative volume. This aspect is crucial for the description of textural hardening; ii. constitutive laws of coupling between slip and twinning mechanisms which is important for proper modelling of the Hall-Petch-type hardening; iii. modification of material parameters for reoriented crystallites to incorporate Basinski's-type hardening.
In relation to the first issue basically two types of modelling schemes can be distinguished in the literature: two-scale and three-scale approaches. In a two-scale approach only the volume effect of twins is accounted for and a key role is played by the reorientation condition [18][19][20][21]. The scheme should ensure the consistency between the twin volume fraction stemming from the constitutive model and the number of twin-reoriented grains in the representative aggregate (see discussion in Refs. [22,23]). In a three-scale approach a two-phase laminate-like substructure of twinned grain is directly accounted for [2,11,14,24,25]. Nevertheless, such models are usually more computationally demanding than two-scale ones. Note that terminology: two-scale and three-scale models is justified by the fact that in the first models only one scale transition rule is specified: this between the single grain micro-level and a macro-level of the polycrystalline aggregate, while for the second case two scale transition rules are needed: the first one between the matrix or twin lamella micro-level and a level of the twinned grain (called a composite grain in Ref. [26]), and the second one between this intermediate level and the polycrystalline aggregate macro-level. More details on this classification can be found in Ref. [23].
In the case of two-scale models, Basinski's hardening is taken into account in a simple way by modifying (usually increasing) slip and twin resistances after reorientation [3,20], also taking into account coplanarity or non-coplanarity of a given system with the primary twinning plane. In three-scale models where the laminate-like structure with a specified orientation of the matrix-twin boundary is present, modifications have usually directional character [26].
More discussion and critical review of all above aspects are presented in Ref. [23]. Recent contributions are reviewed for example in Ref. [11].
In Refs. [21,22] the two-scale framework of crystal plasticity with twinning was proposed. The formulation involved a novel probabilistic twin volume consistent (PTVC) scheme of reorientation, together with rules for the related modification of material parameters, as well as a systematic description of the critical shear stresses evolution in the form of the matrix of self and mutual slip and twin interactions. This formulation, combined with the viscoplastic self-consistent scale transition scheme, was first successfully verified for fcc materials in Ref. [21], as concerns predictions of both the texture development and the mechanical response. The model was also used for predicting texture evolution in hcp commercially pure titanium subjected to severe plastic deformation processes [6] as well as for simulation of the change in misorientation angle distribution for cold-rolled titanium in Ref. [37] where it was incorporated into the three-scale model of grain refinement [38]. Preliminary verification concerning the mechanical response of hcp magnesium alloy was conducted in Ref. [22] using the experimental data available in Ref. [14].
The goal of the present paper is twofold. Firstly, the combined effect of initial texture and twinning on the anisotropy of mechanical properties of AZ31B alloy is investigated experimentally. Observed differences in mechanical response upon the changes in loading directions are explained by the transfer of activity between two types of plastic deformation modes: slip and twinning. Concerning twinning activity, two effects are quantified: the combined textural and Basinski's effect of formation of twinned crystallites (lattice rotation towards the so-called hard orientation with respect to the applied loading) and the impact of twin platelets on slip activity in the parent crystallites. Secondly, identification and validation of the developed modelling framework for AZ31B rods is performed. It should be noted that as compared to the fcc materials the proposed crystal plasticity model involves relatively large number of parameters, so the evolutionary algorithm is developed to identify an optimal set of values for available experimental data.
The paper is organized as follows. After this introductory section, in Section 2 the experimental procedures are described and next results of the performed mechanical tests as well as texture and microstructure analysis for selected samples are presented. Section 3 summarizes the framework of crystal plasticity accounting for twinning [21] employed in the numerical calculations and outlines an algorithm of parameters identification. Next, the results of simulations and their comparison with experimental data are shown. Section 4 presents interpretation of experimental and numerical analyses in view of slip and twin systems activities that is based on the generalized Schmid criterion for initiation of these two modes of plastic deformation. The paper is closed by conclusions.

Experimental procedures
The experiments were performed on the commercial AZ31B magnesium alloy produced by Magnesium Electron Ltd. AZ31B extruded rods with diameter of 25.4 mm were annealed at 350 ∘ C for 2 h. From such prepared material 'large' samples with 25.4 mm diameter and 55 mm height were cut out and then divided in three groups. The first group was left in the initial (undeformed state) whereas the samples from the second and third groups were initially compressed using the MTS 810 testing machine in the extrusion direction up to 2% and 4% of plastic strain, respectively. Next, from each group of 'large' samples the two groups of 'small' samples (7.5 mm diameter and 11.25 mm height) with different orientations, i.e. 45 ∘ and 90 ∘ with respect to the compression axis were cut out. Additionally, the group of samples in the initial state with the compression axis parallel to the extrusion direction was prepared. Then, all 'small' samples were deformed in displacement controlled compression using the MTS 858 hydraulic testing machine until failure. The displacement rate was equal to 0.01125 mm/s what for a given sample geometry corresponded to mean strain rate 0.001 1=s. The scheme of all types of performed experimental tests is shown in Fig. 1.
To sum up, the following mechanical tests were carried out: During compression both force and displacement were recorded and tests were performed for four samples from each group. Moreover, in order to obtain the logarithmic strain-true stress curves an additional test was done in each case, where the strains were determined using the ThermoCorr software with built-in implementation of Digital Image Correlation (DIC) algorithm [39]. Due to the reproducibility of the force-displacement curves, additional tests with DIC were performed on a single sample from each group.
The multidirectional compression tests were performed for two reasons: i. in order to assess which slip and twinning systems are active in the polycrystal with the strong texture; ii. to analyze the influence of the preliminary deformation upon twin formation, and the change of the mechanical response following this. In order to fulfil these goals mechanical testing was supplemented by microstructure analysis.
The texture in both initial state and after deformation was determined using the Bruker D8 Discover X-ray diffractometer with filtered Co Kα radiation. The analyses were performed on the samples' cross sections perpendicular to the extrusion direction (for the non-deformed sample) or the subsequent compression direction (for deformed samples).
The orientation maps obtained using Electron Backscattered Diffraction (EBSD) were used in order to confirm twinning activity suggested by XRD measurements and to identify operating twin family. The analysis was performed using Hitachi S-3500 N Scanning Electron Microscope with the step 0.4 μm on the cross-section parallel to the extrusion axis and lying in the sample's symmetry axis. Fig. 2 shows the true stress -logarithmic plastic strain curves obtained from DIC measurements. The current cross section area was calculated assuming the material's incompressibility. The plastic part of the logarithmic strain was obtained by subtracting the true stress divided by Young's modulus from the logarithmic strain. Young's modulus was determined for the initial part of the measured stress-strain curve. The figure presents the absolute values of stress and strain. In order to facilitate analysis and interpretation of results the curves S2-90 and S2-45 were translated to start at the logarithmic plastic strain equal to 0.02. Similarly, the curves S4-90 and S4-45 were translated to begin at the logarithmic plastic strain equal to 0.04. Note that 0.02 and 0.04 are approximate values of accumulated plastic strain for prestrained samples at the beginning of secondary compression process.

Experimental results
The highest stress, and at the same time the lowest strain at failure, occurs in case of compression in ED (S-0). Three hardening stages can be seen on the curve: I decrease of the strain-hardening rate, II increase of the strain-hardening rate, III decrease of the strain-hardening rate.
According to Refs. [40,41], stage I is similar to the response of the material deforming by slip, and the decrease of the hardening rate results from dynamic recovery. The initiation of stage II is connected with twinning activity. The hardening due to twinning is explained by mechanisms mentioned in the introduction, e. g. the Hall-Petch mechanism or texture hardening (lattice rotation resulting in formation of hard orientations). These observations are in accordance with the texture measurements presented below (see Fig. 3) showing that in the sample compressed in ED many grains become reoriented due to twinning. The renewed decrease of strain-hardening rate observed in stage III is usually ascribed to the saturation of twinning. In this case, further deformation can proceed by slip only, which finally reaches the state of dynamic recovery.
In the case of compression at the angle 45 ∘ to ED (S-45) the occurrence of three hardening stages can be also seen, but not as evidently as during compression in ED direction. It appears that twinning takes place also in this case, but in fewer number of grains than in the previously discussed test (S-0). This results from the fact that during compression at the angle 45 ∘ to ED in the sample having strong basal texture only part of grains have the orientation favourable to T1 twinning. On the other hand, in case of the compression at the angle 90 ∘ to ED (S-90) even smaller number of grains can be easily twinned. This is reflected on the stress-strain curve where the stage II is nearly absent. These observations will be further elaborated and related to the theoretical predictions of crystal plasticity theory in the discussion section.
Initial 2% and 4% of pre-strain by compression results in the increase  of the yield stress in every pre-strained sample. It may be, to a limited extent, the result of the hardening due to dislocation slip, but it seems that the increase of the stresses necessary to plastically deform the sample should be rather attributed to twinning occurring during the initial compression. It can be seen in Fig. 2 that the yield stress is higher in case of the material subjected to 4% initial compression than in case of only 2% initial compression. Furthermore, during compression at an angle 90 ∘ to ED, the stresses in case of S4-90 test remain higher than in case of the S2-90, which is rather intuitive and agrees with a standard strain hardening phenomenon. On the other hand, in the case of compression at an angle 45 ∘ to ED, after 3% of strain in the less pre-strained S2-45 sample the stresses attain higher values than for the same strain level in the more pre-strained S4-45 sample. Moreover, not only stresses are higher in the S2-45 test, but also strain at failure is the largest among all the samples (about 26% strain).
In order to confirm the relation between twinning activity and the mechanical response the microstructure analysis has been performed.  There are two conventions used for the stereographic projection in the case of hexagonal close-packed (hcp) materials. As it was elaborated in Ref. [37] the choice of the convention is crucial for how the pole figures look like. Here, the y direction of the orthogonal coordinate system was specified to be parallel to the 〈2110〉 direction of the hexagonal lattice system.
The undeformed material has strong texture (see Fig. 3a), where the ED is nearly parallel to the most of the basal planes ð0001Þ. The texture only slightly deviates from the expected axial symmetry related to the extrusion process. The textures of the samples compressed up to 2% (cf. Figs. 3b) and 4% (cf. Fig. 3c) of axial plastic strain, are qualitatively similar to each other, but the ð0001Þ fiber occupies about 23% of the total volume in case of the sample compressed up to 2% strain, and about 40% in the sample compressed up to 4%. This texture component was not found in the undeformed material. Such abrupt texture changes can be ascribed to twinning, which reorients the lattice. The most common twinning mode in AZ31B deformed at room temperature is the T1 tensile twinning f1012g〈1011〉 [14,28], which results in lattice rotation by 86.4 ∘ [42]. Compression of the sample having initially the extrusion texture favours the activation of this twinning mode, because it results in the elongation of the crystal cell in the c direction.
On the left of Fig. 3 maps of the crystallographic orientations projected onto the extrusion direction (ED) of the rod are shown. Fig. 3a shows the map for the sample in the annealed state. The EBSD analysis revealed the bimodal grain structure in this case. Figure S2 presenting the grain size and corresponding histogram of grain size distribution for the larger area of the non-deformed material (1000 � 600 μm) is included in the supplementary material. At the histogram the maximum value corresponding to the small grain fraction is observed for approximately 15 μm and in the case of larger grains for approximately 130 μm. Fig. 3b and c shows the maps plotted for the samples compressed up to 2% and 4% in ED, respectively. The presence of twinning (shown in red and orange according to the orientation triangle) can be seen here even more readily than in the XRD pole figures presented on the right. The EBSD analysis qualitatively confirms also that the twin volume fraction increases with increasing strain. It seems that the increase stems both from the emergence of new twins as well as from the broadening of the existing ones. Moreover, the analysis confirms that the T1 tensile twinning is an active twinning mode. Detailed study of the point to point disorientation in the grains present on the orientation maps for the samples after 2% and 4% compression shows that solely a disorientation angle characteristic for T1 type tensile twins (i.e. 86.4 � ) are found. This result agrees with previous studies [36]. An example of the line profile of disorientation angle is presented in Fig. S3 in supplementary material.
It is observed that the change of compression direction seems to reduce twinning activity. Fig. 4a and b present pole figures projected onto the cross-section of the compressed sample after the strain path change. What is mainly seen in these pole figures is a rotated image of the texture resulting from the extrusion and initial prestrain by compression in ED direction.
Reported results of experiments have been used to identify and validate the proposed crystal plasticity modelling framework. Mutually, numerical simulations employing the developed theory have been applied to confirm provided interpretation of experimental observations.

Crystal plasticity framework
The rate-dependent velocity based formulation of crystal plasticity is used. Accordingly, the velocity gradient l is additively decomposed into elastic and plastic parts. The elastic part l e is equal to the lattice spin (the elastic stretches are not considered as they are small as compared to the plastic ones). Twinning is treated as a pseudo-slip [18,27,29]. The plastic part l p of the velocity gradient l is given as follows: In equation (1), M and N are the number of slip and twin systems, respectively, and _ γ kðlÞ � 0 is the rate of shearing on the kðlÞ-th slip or twinning system, _ f l is the rate of the volume fraction of twin mode no. l. γ TW;l is the characteristic twin shear for the twin mode l, which depends on the geometry of the lattice. It is equal to 0.169 for hcp magnesium (c= a ¼ 1:624) and the tensile twinning (T1).
For both, slip and twinning modes, the rate of shearing on system r is specified by the power law: where _ γ 0 is the reference shear rate. The non-negative resolved shear stress τ r is calculated using the projection of the Cauchy stress on the direction m r and plane n r of a given system r, namely The hardening laws governing the evolution of the critical resolved shear stresses (CRSSs) τ r c due to self and mutual slip-twin interactions are reminded in Table 1. The laws and their physical background are discussed in Ref. [22]. Let us only mention that the evolution of CRSSs due to twin is of sigmoidal character, contrary to its exponential evolution due to slip.
Salem et al. [2] and Wand and Agnew [17] pointed out that the mechanical properties of the twins are changed as compared to the parent matrix phase following the Basinski's hardening mechanism. This property is mimicked in a simplified way by multiplying the CRSSs of the twin-reoriented grain by a factor μ. This parameter, together with τ 0 c ; β; τ sat ; f sat and q ðαβÞ for each slip or twinning system constitute the set of material parameters that have to be identified. Concerning twinning activity in the reoriented grains it is assumed that only detwinning can occur.
As a polycrystal model the viscoplastic self-consistent (VPSC) scheme developed in Refs. [43,44] is applied. The code being used is the VPSC-7, with the probabilistic twin volume consistent (PTVC) reorientation scheme [21] and modified hardening laws [22] incorporated. Let us recall that the PTVC scheme, in a statistical way, ensures consistent predictions of the volume effect of twin-related orientations on the texture image. The latter modifications are crucial in the study of influence of twinning activity on the texture evolution and stress-strain response of AZ31B samples.

Evolutionary algorithm for model identification
For hcp magnesium up to seven families of slip and twin systems are relevant: basal ð0001Þ〈1120〉, prismatic f1100g〈1120〉, pyramidal 〈a〉 f1101g〈1120〉, pyramidal I 〈c þa〉 f1011g〈1213〉 and pyramidal II 〈c þa〉 f1122g〈1123〉 slip mode families as well as T1 tensile f1012g〈 1011〉 or compressive f1011g〈1012〉 twin systems [5,28]. In the present study, following e.g. Ref. [14], only three families of slip systems (basal, prismatic and pyramidal II 〈c þ a〉) and tensile twinning are considered for AZ31B alloy. 2 The single crystal model employs then 45 material parameters. Their identification directly by the user, based on his knowledge and expertise, is a valuable approach, because it gives the possibility to understand the physical basis of the acting deformation mechanisms and hardening laws relevant for the given material. This approach is, however, very laborious and time-consuming. Therefore, in the present study the evolutionary algorithm (EA) was applied for optimizing the parameters identification. Its main advantages are the easiness of implementation and the ability to overcome local minima.
EA was implemented in the Python programming language based on the algorithm described in Ref. [49]. Since the idea of EA is the biological evolution, such terms as generation, individual or mutation are used in its description. In order to perform the optimization in a finite time and to obtain results being in line with the physics, the maximum and minimum values of each material parameter have to be set. Moreover, in order to run the algorithm, the number of individuals N ind has to be specified. Each individual contains the set of parameters. In the first generation each of them is randomly chosen from the specified interval. The ranges of the parameters used for the optimization are presented in Tab. S1 included in the supplementary material. Then for each individual the VPSC simulation was performed and the results were evaluated based on their agreement with experimental data.
The evaluation was performed by numerically comparing the stressstrain curves. All 7 experimentally measured stress-strain curves (S-0, S-90, S-45, S2-90, S2-45, S4-90, S4-45) were used for the optimization. The experimental data were first discretized. The curves were compared by evaluating the sums of absolute values of differences in both the stress level and secants to the curve at each point. The secants were computed based on the stress values in the neighbouring points, cf. Fig. S4 in the supplementary material. The specified number of the highest ranked individuals (HRI) N hri , having the least values of those sums, were then Table 1 The hardening laws. See Ref. [22] for details and physical basis.  however, the volume fractions produced by these twin systems are much smaller than those resulting from tensile twins, therefore they do not contribute significantly to the texture image.
chosen. Next the crossover was carried out. It consisted in randomly choosing the pairs of individuals from the set of HRI. The number of pairs was equal to the number of individuals N ind . For each pair of parents, one child was generated. The parameters of the child were calculated using: a) inheritanceeach of the child's parameters was calculated as: where N p1 and N p2 denoted the values of parents' parameters and α was a drawn random real between 0 and 1.
b) mutationeach parameter of the child was mutated when a random real between 0 and 1 was grater than mutation probability P mut . In our implementation mutation consists in randomly generating the parameter from the interval specified at the beginning.
Next, the VPSC simulations, evaluation of results, selection of HRI and creation of the next generation were performed. After running the specified number of generations N gen , the algorithm was stopped and the individual with the minimal value of the difference with the experimental data was chosen.

Results
The identification of material parameters was carried out with the aid of the evolutionary algorithm described in the previous section. Number of individuals N ind was set to 40, and the number of generations was set to 90. After the simulation, for each generation N hri ¼ 8 highest ranked individuals were chosen. The probability of mutation was set to 0.02. The sum of the differences between simulated and experimental stress-strain curves was used as a value to be minimized.
The optimized parameters are collected in Table 2. Fig. 5 shows the stress-strain curves (absolute values) obtained in simulation and experiment (left column) and the current slip and twinning systems activities (right column). For each strain increment a relative current activity of a family of slip or twin systems is defined as, cf [14], The EA was able to give very good agreement between stress level for most of experimental loading programs (see also Fig. S5 in the supplementary material for S2-90 and S2-45 results). The changes in the slopes (or changes in the rate of hardening) were also predicted reasonably well, though in this case more differences can be seen, especially in the case of twinning activity. It indicates that the hardening model for twinning requires further development. However, it should be noted that the predicted twin volume fraction is consistent with the observations made on the basis of experimentally measured textures (see Fig. 5a). This verifies validity of identification procedure since this factor is not used in the identification procedure.
The simulated stress-strain curves can be also presented in one plot as in Fig. 2 to enable additional comparison with experimental data (see Fig. S6 in supplementary material). It is seen that some crucial features, like the highest hardening rate for S-0 curve and the increase of yield stress associated with preliminary compression are correctly reproduced. Unfortunately, the proper fitting of relative values of stress in the case of preliminarily compressed sample was not attained. According to the experimental data, for S-0, S2-90 and S4-90, an increase in preliminary compression stress was in line with the increase in maximal stress, while the opposite can be observed in simulations. Moreover, for S-45. S2-45 and S4-45 similar values of maximal stress were obtained, while the experimental data show significant differences. Fig. 6 shows the simulated ð0001Þ and f1010g pole figures corresponding to the measured ones in Figs. 3a and 4 (see also the results for 2% prestrain in Fig. S7 of the supplementary material). The optimized parameters enabled to achieve good agreement between predicted and measured textures, although the optimization of parameters for the crystal plasticity model was performed only on the basis of the stressstrain curves. Fig. 6 shows pole figures projected onto the cross section perpendicular to the current compression direction for the sample, in analogy to Fig. 4. Note that in case of S-0 test (Fig. 3a), this cross section is equivalent with the cross-section of the extruded rod, while in case S4-90 (Fig. 4a) and S4-45 ( Fig. 4b) tests, the cross section is rotated wrt. to the cross section of the extruded rod about an angle 90 ∘ and 45 ∘ , respectively. In order to properly assess the texture evolution, the PFs in the cases b and c were projected also on the cross-section perpendicular to the extrusion direction, see Fig. S8 in supplementary material. It has been observed that the dislocation slip, which dominates deformation up to 3.5% of strain for pre-strain samples, was not enough to result in a significant change of the texture resulting from the initial compression (see Fig. 6a). The twinning was acting only to a limited extent (see Fig. 5), thus the results of its activity are hardly noticeable in the pole figures. A slight asymmetry can be seen, which results from the change of loading direction.

Discussion
As indicated in Sec. 2.2 there are four families of slip and twin systems for magnesium alloy considered in the present study. At room temperature, among slip modes the prismatic and pyramidal 〈c þa〉 systems are the hard modes, while basal slip and tensile twinning are easier to initiate. Note that pyramidal 〈c þa〉 slip systems are necessary to fulfil the condition of five independent slip systems, and thus to realize some of the deviatoric deformation paths for single crystal. On the other hand a set of pyramidal slip systems itself is sufficient to perform arbitrary deviatoric deformation, which is not the case for other families considered. Thus, it is possible that for some assumed values of critical shear stress for a given family the respective slip systems are never activated. Such possibility is discussed below and next, possible configurations of compression loading with respect to the crystallite orientations applicable to the reported experimental data are analyzed as concerns the foreseen activity of deformation modes. Let us now, using the generalized Schmid law, calculate the yield stress τ k Y for a pure shear stress state τ k defined subsequently by the system k belonging to basal, prismatic and pyramidal 〈c þa〉 slip, and tensile twin families, namely (no summation over k): where P k ¼ 1 2 ðm k �n k þn k �m k Þ is the symmetric part of the Schmid tensor. From the Schmid law the resolved shear stresses on the system r  (Eq. (4)). The results were carried out for a) ED compression, b) compression perpendicularly to ED, c) compression at an angle 45 o to ED, d) compression perpendicularly to ED of the material preliminarily compressed in ED up to 4%, e) compression at an angle 45 o to ED of the material preliminarily compressed in ED up to 4% (results for 2% prestrain can be found in the supplementary material). are then given by (no summation over k in the second equation): Note that in this way we maximize the resolved shear stress for a slip or twin system r ¼ k. Taking into account that the plastic flow will start , substituting for k the slip or twin system of interest and evaluating maximum values of respective products 2P k ⋅P r for a given family we obtain: Fig. 6. Simulated ð0001Þ and f1010g pole figures projected on the plane of the sample's cross-section calculated for the texture: a) 4% ED compression, b) initial 4% ED compression and 3.5% compression perpendicular to ED and c) initial 4% ED compression and 3.5% compression at an angle 45 o to ED (additional results can be found in supplementary material). The figures were plotted using the ATEX software [50].
One concludes that basal slip is not initiated if α basal > 2:22, prismatic slip if α prism > 2:59 and tensile twinning if α twin > 1:24. As concerns the basal and prismatic slip the same bounds on their activity are obtained if a rigorous analysis of a deformation mode map is performed following the procedure formulated in Refs. [51,52]. According to the procedure the vertices of the Schmid envelope obtained for pyramidal slip are used for the verification of possible activity of other modes. As concerns activity of tensile twinning the rigorous limit value for α twin is equal to 1.86. More details concerning such analysis of the Schmid envelope are included in the supplementary material.
For the identified values of material parameters collected in Table 2 it is found that in the initial state the corresponding values of α calculated for τ c0 are One may also distinguish the state for which hardening for slip systems has already saturated, while twinning is still possible. If in order to simplify analysis it is assumed that τ twin c ¼ τ twin c0 then in such state � α basal ; α prism ; α twin � ¼ f0:233; 0:285; 0:399g : As it is seen in each case the respective value of α for prismatic and basal slip in Eqs. (9)-(11) is below the bounds specified above, so these slip families should be active for some deformation paths throughout the whole process. As concerns twinning it will be active at the initial and intermediate stage of proportional deformation process and then it will be no longer possible unless the strain path is changed. Now, let us analyze theoretically possible activity of four deformation modes in the uniaxial compression with a set loading direction n with respect to the hcp crystallite axes (see Fig. 7), namely σðnÞ ¼ À Y ϕ;φ nðϕ; φÞ � nðϕ; φÞ ; (12) where nðϕ; φÞ ¼ sin φ cos ϕk 1 þ sin φ sin ϕk 2 þ cos φk 3 ; and k 3 is coaxial with the ð1010Þ or ð0001Þ crystallographic direction for f1010Þ or ð0001Þ fiber texture components, respectively, k 1 with ð1210Þ and k 2 ¼ k 3 � k 1 . Three configurations considered in the following, namely: correspond to the set of compression experiments reported in previous sections. As concerns material parameters, for the initial state we have assumed the relation between the critical shear stresses τ r c given by specification (12), for the saturated hardening state by (13), and for intermediate state, when twinning is still possible but hardening of slip has saturated, by (14). It should be underlined that the analysis is performed assuming the ideal textures (without a scatter) and neglecting the continuous lattice rotation due to slip activity. Equivalently to the pure shear loading (5), the yield stress Y ϕ;φ of single crystal in uniaxial compression is specified as in Eq. (7) with replacement of 2P k ⋅P r by n⋅P r ⋅n.
First, let us consider the value of the yield stress when the compression axis is perpendicular to f1010g crystal plane (φ ¼ 0 o and k 3 jjð1010Þ in Eq. (12), cf. Fig. 7a). Such configuration corresponds to the one encountered in S-0 test for crystallites being the elements of ideal extrusion texture for which ED is perpendicular to f1010g crystal plane and parallel to the basal plane ð0001Þ. It can be verified that, using the Schmid law, the yield stress is then independent of angle ϕ, and at the initial stage activity of prismatic slip is foreseen, which with the advances of hardening for this mode is replaced by twinning. For the part of crystallites that are not reoriented by twinning, when the twinning saturates, then again prismatic slip will be active. The part of crystallites that were reoriented by twinning have the orientation with c-axis coaxial with compression direction (k 3 jjð0001Þ), cf. Fig. 7b. In such a case activity of pyramidal slip is predicted. Fig. 8a and b presents the value of the yield stress when the compression axis is parallel to f1010g crystal plane (φ ¼ 90 o in Eq. (12) and k 3 jjð1010Þ, c.f. Fig. 7a). It is a parametric plot where the yield stress is specified by the length of vector connecting the frame origin and the Schmid envelope curve (the orange line). The direction of vector represents the direction of compression axis in the f1010g crystal plane (note that n x ¼ n in Fig. 7a for this case). This configuration corresponds to the one encountered in S-90 test for crystallites being the elements of ideal extrusion texture. As it is seen initially (Fig. 8a) in the most of crystallites basal slip is active, however, small amount of crystallites require activation of prismatic slip (the compression axis aligned with a-axes of crystallites) or pyramidal slip (the compression axis aligned with c-axis of crystallites). In the intermediate regime, when hardening of slip mechanisms advances, in some crystallites, mainly those which required activation of prismatic slip, twinning is activated instead (Fig. 8b). Finally, when twinning is no more possible, those crystallites (non-reoriented part) again activate prismatic slip. Twinning reorients crystallite by the rotation close to 90 o in the plane of the presented parametric plots, so that the reoriented crystallites have their c-axis aligned with the loading direction n. As such they require activity of pyramidal slip.
In samples S2-90 and S4-90, the crystallites, which have not been reoriented due to twinning in the course of pre-strain, accommodate the plastic deformation in a similar fashion as grains in the S-90 sample (cf. Fig. 8a and b). Fig. 9a illustrates a situation for samples S2-90 or S4-90 and those crystallites which were reoriented due to twinning during prestrain process (φ ¼ 90 o , k 3 jjð0001Þ in Eq. (12), cf. Fig. 7b). All such grains have their c-axis nearly perpendicular to the compression direction. Initially (for the subsequent straining process) in most of those crystallites either prismatic or twinning mechanism will be active (depending on the advances in hardening). Those crystallites that again reorient due to twinning go back to the position in which c-axis is aligned with the loading direction (φ ¼ 0 o and k 3 jjð0001Þ in Eq. (12), cf. Fig. 7b). Fig. 8c and d illustrates compression of the sample S-45. As previously the value of the yield stress is specified by the length of vector connecting the frame origin and the Schmid envelope curve. The direction of vector represents in this case the direction of compression axis  12)). Initially in most crystallites basal slip is active (Fig. 8c), however, small amount of crystallites requires activation of prismatic slip (vector n x aligned with a axis of crystallites). As hardening advances, more and more crystallites activate prismatic slip (according to Table 1 basal slip hardens more than prismatic slip and at some stage of the process the critical shear stress for this mode becomes higher than for prismatic slip). In the intermediate regime, when hardening of this slip mechanism advances those crystallites which required activation of prismatic slip, activate twinning instead (Fig. 8d). Finally when twinning is no more possible those crystallites (non-reoriented part) again activate prismatic slip. The crystallites reoriented due to twinning in the course of S-45 test align their c-axis with the loading axis. They are therefore prone to pyramidal 〈c þa〉 slip (the case φ ¼ 0 o and k 3 jjð0001Þ in Eq. (12)). For samples S2-45 and S4-45 (the part of crystallites that have not been reoriented due to twinning in the course of pre-strain) the scenario is the same as described above for S-45 sample. Fig. 9b presents the case of those crystallites in the S2-45 and S4-45 samples that reoriented due to twinning during prestrain (φ ¼ 45 o and k 3 jjð0001Þ in Eq. (12), cf. Fig. 7b). It is seen that for this group of crystals basal slip is expected to be active during the whole process even though it hardens more than prismatic slip.
Above analyses were performed under simplifying assumption of equal stress state in all grains (Sachs hypothesis), however, it quite well explains the predicted evolution of current relative activity of modes obtained in the simulations presented on the right in Fig. 5. Indeed, as it can be seen in Fig. 5a, after initial activity of basal slip, which rotates the crystal towards the perfect f1010gkn orientation, the prismatic slip is activated, as predicted by the Schmid law. Then the activity of twinning rises more and more with the saturating CRSSs for slip, while after the saturation of twinning, the prismatic slip becomes again the dominant deformation mechanism, which also agrees with the analytical predictions presented above. In addition, in accordance with the Schmid law, in the reoriented grains the pyramidal slip becomes active and its activity increases with increasing strain.
In case of S-90 (Fig. 5b), most of the plastic deformation is accommodated by the basal slip and only limited twinning activity is predicted, consistently with Fig. 8a and b. The activity in S-45 (Fig. 5c) is similar to the one in case of S-90, although more prismatic slip and twinning is predicted, consistently with Fig. 8c and d.
An interesting case is also the plastic deformation of samples predeformed in compression up to 2 and 4% plastic strain. The samples predeformed to 2% plastic strain (S2-90 and S2-45) behave similarly to the samples without initial compression (respectively S-90 and S-45). Namely, most of the strain is carried out by the basal slip together with the prismatic one (see Fig. S5 in the supplementary material). This is because the twin volume fraction attained in the preliminary deformation is low and the untwinned grains carry most of the deformation. However, in case of S2-90 the activity of the prismatic slip is higher than in case of S-90. This seems to be a result of the prismatic slip activity in reoriented grains and is consistent with Fig. 9a where the prismatic slip and twinning were predicted to be the only active mechanisms both in the initial and the saturated states. In case of S2-45 the amount of prismatic slip is even slightly smaller than in case of S-45. On the other hand, the basal slip is visibly more active. We observe that this is the result of the basal slip activity in reoriented grains, consistently with Fig. 9b. Further twinning also plays some role in this case.
Finally, in samples prestrained up to 4% initial compression, the twin volume fraction is considerably larger than in the 2% case. An increase of the activity of the prismatic slip (in case of S4-90) and basal slip (in case of S4-45) is more pronounced than in case of S2-90 and S2-45, respectively, cf. Fig. 5d and e. The higher increase results from the larger twin volume fraction. The simulated slip activities are thus consistent with Fig. 9a and b also in case of 4% initial deformation.

Conclusions
In the paper texture and twinning-induced anisotropy of the yield stress and hardening of AZ31B extruded rods have been analyzed. First, response of samples under the multidirectional compression tests, involving strain path changes, has been investigated experimentally. The set of tests has been designed with a goal to activate different set of slip and twinning systems in the polycrystalline sample with the strong texture. This enables us to analyze the influence of the preliminary deformation upon twin formation and their impact on the hardening phenomenon. Mechanical testing has been supplemented by microstructure analysis. Next, the two-scale crystal plasticity model accounting for twinning and employing the viscoplastic self-consistent scheme has been validated using the acquired experimental data. In the same time the results of simulations verified proposed interpretation of experimental observations. It should be underlined that the analysis of the AZ31B rods subjected to the 7 different multidirectional deformation paths has been performed using a single set of parameters for the crystal plasticity model. These parameters have been identified using the developed identification procedure based on the evolutionary algorithm. The set of parameters, established using stress-strain curves, made it possible to correctly predict the experimentally observed textures and twin volume fractions. To the authors' best knowledge presented set of experimental data accompanied by the crystal plasticity analysis is not yet found in the literature. Finally, all results are discussed in view of the theoretical study of slip and twinning activity on the basis of the generalized Schmid criterion.
The performed studies enable us to draw the following conclusions.
� The anisotropy of the mechanical response of strongly textured AZ31B extruded rods beyond the yield point cannot be neglected and its proper identification requires extended set of tests involving the strain path changes. � Twinning activity plays a crucial role in the observed anisotropy. It influences the mechanical response predominantly by the texture change that invokes the orientational hardening, and to lesser extent by modification of hardening moduli due to slip-twin interactions. � Both mean field crystal plasticity model and the generalized Schmid law for slip and twinning are useful tools to analyze the plastic behaviour of AZ31B. While the Schmid law is beneficial for overall qualitative understanding of a deformation scenario, the mean field crystal plasticity model enables analysis of the slip-twin interactions, texture evolution and their effect upon the mechanical behaviour in a more thorough quantitative fashion. � The simple evolutionary algorithm is an effective means for the identification of crystal plasticity model parameters. � Although the proposed computationally efficient two-scale approach based on [21] provides the predictions of satisfactory agreement with experimental outcomes, the observed discrepancies for the tests involving the strain path changes have indicated that the component of the model related to the modification of material parameters of the re-oriented grains requires further elaboration and/or modification.

Declaration of competing interest
The authors declare that they have no known competing financial