**3. Modification of Mesolevel Relations with an Explicit Separation of the Rigid Moving Coordinate System with Respect to Which Elastic Distortion Is Defined** *3.1. Relations in Terms of Lattice Unloaded Configuration (LU-Model)*

In [63,76], it was proposed to modify the decomposition of motion at a mesolevel, i.e., to realize a multiplicative representation of the deformation gradient *f* with an explicit separation of the rigid moving coordinate system (MCS) Op1p2p3:

$$f = f' \cdot f^p = \widetilde{f}' \cdot r \cdot f^p. \tag{14}$$

Here, *r* = *k*ˆ*<sup>i</sup>* ⊗ - *k*ˆ*i t*=0 ≡ *k*ˆ*<sup>i</sup>* ⊗ *k<sup>i</sup>* is the proper orthogonal tensor which converts the reference basis MCS) *<sup>k</sup><sup>i</sup>* into the current basis *<sup>k</sup>*ˆ*i*; in other words, *<sup>r</sup>* is the tensor of rotation of the MCS from reference to actual configuration; *f <sup>p</sup>* is the component of the deformation gradient that transforms the initial configuration into the plastically deformed configuration, *f e* is the component of the deformation gradient component that transforms the plastically deformed configuration into the actual configuration (this tensor also characterizes crystalline lattice distortion), and *f e* , *f <sup>p</sup>* are, in the general case, the non-symmetric second rank tensors.

The motion of the MCS is a rigid motion, and the motion relative to the MCS is a deformation motion: geometrically generalized nonlinear constitutive relations are formulated by an observer related to the MCS [65]. Relation (14) is a modification of the classical Kröner–Lee decomposition (1) via the explicit separation of the quasi-rigid motion *r*. Note that, in the general case, the tensor does not coincide with the orthogonal tensor *r<sup>e</sup>* in the polar expansion *f <sup>e</sup>* = *r<sup>e</sup>* · *u<sup>e</sup>* = *v<sup>e</sup>* · *r<sup>e</sup>* in the case when the spin of the corresponding MCS is defined (2).

For analytical calculations, we also introduce the crystallographic coordinate system (CCS) Oy1y2y3 with basis *qi*. Without loss of generality, we assume that, in the reference configuration, the CCS basis *q<sup>i</sup>* (0) is an orthonormal basis that coincides with *ki*. Note that the lengths of the CCS basis vectors and the angles between them vary during deformation.

The constitutive laws are formulated in the context of the observer related to the MCS. This MCS is assumed to be linked to the symmetry axes of the material. Probably, the idea of the necessity of determining corotational derivatives with reference to the symmetric elements of the material (directors) was initially put forward by J. Mandel [77], yet, unfortunately, his study contains no specific relations for describing a crystallite lattice spin.

Based on the introduced multiplicative decomposition (14) and assuming that the elastic properties of the crystallite in MCS are constant, the constitutive elastoviscoplastic

relation was formulated in terms of the unloaded configuration (unloading is carried out with (*f e* ) −1 , and the MCS remains fixed) [63,76]:

$$
\overline{\mathbb{K}} = \overline{\mathbb{K}} : \mathbb{E}^{\mathfrak{c}},
\tag{15}
$$

which includes the second Piola–Kirchhoff tensor *k*, defined in the unloaded configuration, *k* = *k ijk*ˆ*<sup>i</sup>* <sup>⊗</sup> *<sup>k</sup>*ˆ*<sup>j</sup>* <sup>=</sup> *<sup>J</sup>*(*<sup>f</sup> e* ) −1 · *σ* · (*f e* ) <sup>−</sup>*<sup>T</sup>* and the four-rank elastic tensor *<sup>π</sup>* <sup>=</sup> *<sup>π</sup>ijmnk*ˆ*<sup>i</sup>* <sup>⊗</sup> *k*ˆ*<sup>j</sup>* ⊗ *k*ˆ *<sup>m</sup>* ⊗ *k*ˆ *<sup>n</sup>*, *c<sup>e</sup>* is the elastic constituent of the right Cauchy–Green deformation tensor *c<sup>e</sup>* = *cijk*ˆ*<sup>i</sup>* ⊗ *k*ˆ*<sup>j</sup>* = 1/2((*f e* ) *<sup>T</sup>* · *<sup>f</sup> e* − *I*), and *I* is the unit tensor. The accepted definition of the elastic tensor (its components are constant in the MCS) conforms to the principle of independence of constitutive relations from reference choice [21].

The plastic component of the deformation gradient *f <sup>p</sup>* is determined using the viscoplastic or elastoplastic relations applied to describe the intragranular dislocation slip. In the model, as in other compared models, we use Relations (4)–(6).

An important component of the constitutive model is the kinetic relation for the motion of the MCS. As noted above, it seems reasonable to relate the mesoscale MCSs and the symmetry elements of the crystal lattice [65,76,77]. Plastic deformations occur due to the movement of dislocations and they do not cause symmetry distortion, and therefore *r* can be completely determined by the tensor *f e* .

The elastic gradient *f <sup>e</sup>* = *f e* · *r* contains both the quasi-rigid motion (rotation) of the MCS and the distortion of the lattice relative to it. Generally speaking, there are different variants of this relation. In order to overcome this challenge, it is necessary to accept a hypothesis for the definition of *r*. We consider here one of the options. Let a relation between CCS and MCS be as follows: (1) the axes Oy<sup>1</sup> and Op<sup>1</sup> coincide at every instant of deformation (the vector *k*ˆ <sup>1</sup> is directed along the vector *q*1); (2) the vector *k*ˆ <sup>2</sup> is located in the plane Oy1y2 at each instant of deformation. In [65], the following relation for the MCS spin is presented:

$$
\varpi\_1 = \dot{\boldsymbol{r}} \cdot \boldsymbol{r}^T = \dot{\boldsymbol{k}}\_i \otimes \boldsymbol{\k}\_i = \begin{array}{c} -(\boldsymbol{\hat{k}}\_2 \cdot \mathbb{I}^\epsilon \cdot \boldsymbol{\hat{k}}\_1)\boldsymbol{\k}\_1 \otimes \boldsymbol{\k}\_2 + (\boldsymbol{\hat{k}}\_2 \cdot \mathbb{I}^\epsilon \cdot \boldsymbol{\hat{k}}\_1)\boldsymbol{\k}\_2 \otimes \boldsymbol{\k}\_1 - \\\ -(\boldsymbol{\hat{k}}\_3 \cdot \mathbb{I}^\epsilon \cdot \boldsymbol{\hat{k}}\_1)\boldsymbol{\k}\_1 \otimes \boldsymbol{\k}\_3 + (\boldsymbol{\hat{k}}\_3 \cdot \mathbb{I}^\epsilon \cdot \boldsymbol{\hat{k}}\_1)\boldsymbol{\k}\_3 \otimes \boldsymbol{\k}\_1 - \\\ -(\boldsymbol{\hat{k}}\_3 \cdot \mathbb{I}^\epsilon \cdot \boldsymbol{\hat{k}}\_2)\boldsymbol{\k}\_2 \otimes \boldsymbol{\k}\_3 + (\boldsymbol{\hat{k}}\_3 \cdot \mathbb{I}^\epsilon \cdot \boldsymbol{\hat{k}}\_2)\boldsymbol{\k}\_3 \otimes \boldsymbol{\k}\_2. \end{array} \tag{16}
$$

where *l <sup>e</sup>* <sup>=</sup> . *f e* · *f <sup>e</sup>*−<sup>1</sup> is the elastic component of the velocity gradient. In [78], the authors described the way to define *r* in finite form taking into account the above relation between the CCS and the MCS and showed the compliance of the obtained result with that obtained by integrating (16).

Expression (16) can take a more compact form [79]:

$$
\varpi\_1 = I \times (\pounds\_1 \otimes \pounds\_2 \otimes \pounds\_3 - \pounds\_2 \otimes \pounds\_1 \otimes \pounds\_3 + \pounds\_3 \otimes \pounds\_1 \otimes \pounds\_2) : I^c \tag{17}
$$

where «×» denotes the vector product.

A model that involves kinematic Relations (4), (14), and (17) and constitutive law (15) is termed here the "LU1-model" (lattice unloading configuration type of model with spin 1).

Besides, the MCS can be related to other symmetry elements of the crystal lattice. As an alternative, we consider the case when the axes Oy<sup>3</sup> and Op3 coincide at each instant of time (the vector *k*ˆ <sup>3</sup> is directed along the vector *q*3), the vector *k*<sup>1</sup> is positioned at each instant of deformation in the plane Oy1y3. Under these conditions, the spin can be calculated:

$$
\overline{\omega}\_2 = \mathbf{I} \times (-\hat{\mathbf{k}}\_1 \otimes \hat{\mathbf{k}}\_3 \otimes \hat{\mathbf{k}}\_2 + \hat{\mathbf{k}}\_2 \otimes \hat{\mathbf{k}}\_3 \otimes \hat{\mathbf{k}}\_1 + \hat{\mathbf{k}}\_3 \otimes \hat{\mathbf{k}}\_1 \otimes \hat{\mathbf{k}}\_2) : \mathbf{I}^\epsilon. \tag{18}
$$

We call the modified version of the LU-model with spin (18) the "LU2-model". Consider the issue regarding the difference between LU1-model and LU2-model.

We assume that, at equal kinematic effects *f*, the current positions of the MCS related to the symmetry axes of the lattice in different ways can be described by different tensors *r* and *r*∗ (due to the smallness of elastic distortions, these tensors are similar). So, the deformation gradient can be written as *f* = *f e* · *r* · *f <sup>p</sup>* = *f e* ∗ ·*r* ∗ ·*f p*∗.

Based on the constitutive Relation (15), in the case when the orientation of the MCS is defined via the application of the tensor *r*, the Cauchy stress are established as:

$$
\sigma = 1/\left| \overline{f}' \cdot \left( \overline{\pi} : 1/2 \left( \overline{f}'^T \cdot \overline{f}' - 1 \right) \right) \cdot \overline{f}'^T,\tag{19}
$$

where *π* = *πijmnk*ˆ*ik*ˆ*jk*ˆ *<sup>m</sup>k*ˆ *<sup>n</sup>* = *πijmn*(*r* · *ki*)(*r* · *kj*)(*r* · *km*)(*r* · *kn*). Relation (19) is easy to rearrange [80] into the form:

$$\sigma = 1/\left| f^{\epsilon} \cdot \left( \overline{\pi}^{ijmn} \mathbf{k}\_i \otimes \mathbf{k}\_j \otimes \mathbf{k}\_m \otimes \mathbf{k}\_n : 1/2 \left( f^{\epsilon T} \cdot f^{\epsilon} - I \right) \right) \cdot f^{\epsilon T} \right| \tag{20}$$

Analogously, when the position of the MCS position is defined via the application of the tensor *r*∗*,* the Cauchy stress can be obtained from the relation:

$$
\sigma\* = 1/\ulcorner \ulcorner \* \ulcorner \* \cdot \left(\overline{\pi}^{ijmn} \Bbbk\_i \otimes \Bbbk\_j \otimes \Bbbk\_m \otimes \Bbbk\_n : 1/2 \left(\overline{f}^T \ast \cdot \overline{f}^\epsilon \ast - I\right)\right) \cdot \overline{f}^{\epsilon^T \ast} \ast \tag{21}
$$

Comparative analysis of (20) and (21) shows that the difference in the definitions of stresses obtained using the LU1- and LU2-models can be attributed to the fact that *f <sup>e</sup>* and *f e* ∗ are determined in different ways in case of inelastic deformation. As noted above, the definition for shear stress can be expressed as *τ*(*k*) = *k* : *n*(*k*) ⊗ *b*(*k*) . The unit vectors of the slip direction and the normal to the plane are found in terms of the actual configuration as *b*(*k*) = *f e* · *o b* (*k*) /*f e* · *o b* (*k*) , *n*(*k*) = *f <sup>e</sup>*−*<sup>T</sup>* · *o n* (*k*) /*f e* · *o n* (*k*) . It is evident that, at the elastic deformation stage and in a transition to plasticity (as *f <sup>e</sup>* = *f e* ∗), the crystallographic characteristics will be the same, and therefore the reasons for the differences between *f <sup>e</sup>* and *f e* ∗ that may occur are absent.

Thus, the mesostresses obtained within the framework of both LU1- and LU2-models are equal. The final orientations of the lattices determined using the models are different, though due to the smallness of elastic distortions, they differ only slightly [80]. Since the LU1-and LU2-models give almost the same definitions of stresses, we will call them an "LU-model" (if the case in point concerns the definition of stresses).

It is easy to show that the constitutive relation of the LU-model is independent of the imposed rigid motion [21]. Actually, let's consider deformation gradient *<sup>f</sup>* <sup>=</sup> *<sup>O</sup><sup>T</sup>* · *<sup>f</sup>*, where *O* = *O*(*t*)(*O*(0) = *I*) is the orthogonal tensor describing the rigid rotation added to the initial motion. For pure elastic deformation, we get *f <sup>e</sup>* <sup>=</sup> *<sup>O</sup><sup>T</sup>* · *<sup>f</sup> e* , and the analysis of (20) yields *σ*(*f e* ) = *O<sup>T</sup>* · *σ*(*f e* ) · *<sup>O</sup>*. In view of *<sup>σ</sup>* <sup>=</sup> *<sup>O</sup><sup>T</sup>* · *<sup>σ</sup>* · *<sup>O</sup>*, we have *<sup>σ</sup>* <sup>=</sup> *<sup>σ</sup>*(*<sup>f</sup> e* ). Thus, measures of the stress-strain state with an imposed rigid movement are related by exactly the same function. It means that the principle of independence of the constitutive relation from the imposed rigid motion is fulfilled. Hence, it follows that at the onset of the plastic flow*f <sup>p</sup>* <sup>=</sup> *<sup>f</sup> <sup>p</sup>* is correct, and, as a result, the principle of independence of the constitutive law from the imposed rigid motion (i.e., *<sup>σ</sup>* <sup>=</sup> *<sup>σ</sup>*(*<sup>f</sup> e* ) will also be fulfilled in the case of inelastic deformation as well.

By virtue of the fact that the used stress and strain measures are energy conjugate, the requirement for the absence of stress hysteresis and energy dissipation in the arbitrary closed cycles of elastic deformations is fulfilled automatically in the LU-model [63], as in the U-model.

In the general case, other physically justified rotation models (in particular, those capable of taking into account the interaction of defects in neighboring crystallites [81] or the contribution of grain boundary sliding [13]) can be incorporated into the proposed approach to determine the MCS spin.

#### *3.2. Rate Form Relations in Terms of the Actual Configuration (LR-Model)*

From the formulation of the LU-model with the constitutive law in finite Form (15) through the corotational differentiation [64,76], we can pass to the rate analogue for the observer in the MCS:

$$
\overline{k}^{CR} = \overline{\pi} : \overline{\mathfrak{c}}^{CR},
\tag{22}
$$

where *k CR* = . *<sup>k</sup>* <sup>−</sup> *<sup>ω</sup>*<sup>1</sup> · *<sup>k</sup>* <sup>+</sup> *<sup>k</sup>* · *<sup>ω</sup>*1, *<sup>c</sup>eCR* <sup>=</sup> . *c e* − *ω*<sup>1</sup> · *c<sup>e</sup>* + *c<sup>e</sup>* · *ω*<sup>1</sup> (for implementation of this transition, we use the property of constancy of *π* for the observer in the MCS).

It was shown in [63,82] that, at *f <sup>e</sup>* ≈ *<sup>I</sup>* which is characteristic of metals, the following formulation in the rate form can be derived in terms of the current configuration.

Since, at *f <sup>e</sup>* ≈ *<sup>I</sup>*, the unloaded and actual lattice configurations are close, the following estimate is valid for the velocity gradient:

$$d \approx \dot{\tilde{f}}^{\varepsilon} \cdot \tilde{f}^{\varepsilon - 1} + \varpi\_1 + r \cdot \dot{\tilde{f}}^p \cdot f^{p - 1} \cdot r^T = \varpi\_1 + \dot{\tilde{f}}^{\varepsilon} \cdot \tilde{f}^{\varepsilon - 1} + \sum\_{k = 1}^{K} \dot{\gamma}^{(k)} \tilde{\mathsf{b}}^{(k)} \odot \mathfrak{R}^{(k)},\tag{23}$$

where *b* (*k*) , *n*(*k*) are the vectors of the slip direction and the normal to the slip plane in the unloaded lattice configuration, i.e., they are in the actual position of the MCS. If the unloaded and actual lattice configurations are in close proximity, then we can assume *b* (*k*) <sup>=</sup> *<sup>b</sup>*(*k*) , *<sup>n</sup>*(*k*) <sup>=</sup> *<sup>n</sup>*(*k*) for the latter relation and use *<sup>K</sup>* ∑ *k*=1 . *<sup>γ</sup>*(*k*) *b*(*k*) ⊗ *n*(*k*) as an approximation

for the inelastic strain rate.

According to (23), the clear motion decomposition is established in the rate form: the velocity gradient *l* is represented as a set of the MCS spin and the rate of elastic distortions and inelastic deformations determined by the observer related to this MCS.

Taking into account the described approximations, the constitutive law (22) is close to

$$\mathbf{k}^{cr} = \overline{\boldsymbol{\pi}} : \left( \dot{\overline{\boldsymbol{f}}}^{c} \cdot \overline{\boldsymbol{f}}^{r-1} \right) = \overline{\boldsymbol{\pi}} : \left( \boldsymbol{I} - \overline{\boldsymbol{\omega}}\_{1} - \sum\_{k=1}^{K} \dot{\boldsymbol{\gamma}}^{(k)} \mathbf{b}^{(k)} \odot \mathbf{n}^{(k)} \right) \tag{24}$$

where *<sup>k</sup>* is the weighted Kirchhoff stress tensor, and *<sup>k</sup>cr* <sup>=</sup> . *k* + *k* · *ω*<sup>1</sup> − *ω*<sup>1</sup> · *k* is its corotational derivative. Let us call the model involving the kinematic Relations (8), (17), and (23) and a constitutive law (24) the "LR1-model" (lattice rate type of model with spin 1). We also consider an analogous LR2-model in which Relation (18) is used for the spin. It was shown in [80] that the Cauchy stress determined using the formulation of this kind depends only slightly on the way they are related to specific material symmetry elements.

From the above discussion, we can draw a conclusion that the formulations of the Uand LU-models can provide the fulfillment of the conservative conditions for elastic deformation. On the other hand, for the study of technological processes of thermomechanical treatment with an a priori unknown configuration of the computational domain (including contact with the tool), the formulations in terms of the actual configuration in the rate form offer the most promise. Another advantage of such formulations is the possibility of an additive decomposition of the inelastic deformation rate into the contributions of different deformation mechanisms. This explains why the T-model is also popular.

The formulation based on the decomposition of motion with an explicit separation of the motion of the moving coordinate system (LU-model) makes possible a theoretically substantiated transition to a similar formulation in the rate form in terms of the actual configuration (LR-model). As will be illustrated below, this helps to easily show that all the models considered here give a similar response under equal influences.

#### **4. Results and Discussion**

In this section, the issue regarding a comparison of the formulations for crystal plasticity models described above is considered.

## *4.1. Analytical Comparison*

Based on the Equation of motion (14), constitutive law (3) for the U-model can be written as:

$$\left(\boldsymbol{\mathcal{I}}\boldsymbol{r}^{\mathrm{T}}\cdot\left(\boldsymbol{\widetilde{\boldsymbol{f}}}\right)^{-1}\cdot\boldsymbol{\sigma}\cdot\left(\boldsymbol{\widetilde{\boldsymbol{f}}}\right)^{-\mathrm{T}}\cdot\boldsymbol{r} = \boldsymbol{\overline{\boldsymbol{\pi}}}^{ijmn}\boldsymbol{k}\_{i}\otimes\boldsymbol{k}\_{j}\otimes\boldsymbol{k}\_{m}\otimes\boldsymbol{k}\_{n}:\boldsymbol{1}/2\left(\boldsymbol{\boldsymbol{r}}^{\mathrm{T}}\cdot\left(\boldsymbol{\widetilde{\boldsymbol{f}}}\right)^{\mathrm{T}}\cdot\boldsymbol{r}-\boldsymbol{I}\right),\tag{25}$$

where *k<sup>i</sup>* is the reference basis of the MCS related to the lattice symmetry.

The rearranged form of Relation (25) is:

$$J(\widetilde{f}')^{-1} \cdot \sigma \cdot (\widetilde{f}')^{-T} = \pi i^{jmn} (\boldsymbol{\tau} \cdot \boldsymbol{k}\_i) \otimes (\boldsymbol{\tau} \cdot \boldsymbol{k}\_j) \otimes (\boldsymbol{\tau} \cdot \boldsymbol{k}\_m) \otimes (\boldsymbol{\tau} \cdot \boldsymbol{k}\_n) : 1/2 (\left(\widetilde{f}'\right)^T \cdot \widetilde{f}' - I). \tag{26}$$

Relation (26) can be represented as:

$$\left(\left(\widetilde{\boldsymbol{f}}'\right)^{-1}\cdot\boldsymbol{\sigma}\cdot\left(\widetilde{\boldsymbol{f}}'\right)^{-T} = \mp^{ijmn}\hat{\boldsymbol{k}}\_i\otimes\hat{\boldsymbol{k}}\_j\otimes\hat{\boldsymbol{k}}\_m\otimes\hat{\boldsymbol{k}}\_n:\mathbf{1}/2\left(\widetilde{\boldsymbol{f}}'\right)^T\cdot\widetilde{\boldsymbol{f}}'-\boldsymbol{I}\right).\tag{27}$$

Relation (27) corresponds to the constitutive law of the LU-model (15). Thus, the stresses in the U- and LU-models are equally determined. As shown in Section 3.1, the determined stresses are independent of the portion of the MCS motion in the deformation gradient. When we use models of this kind, it is important to keep in mind this distinction and verify separately texture description.

If, by analogy with the transition from the LU-model to the rate LR-model described in Section 3.2, we take for the U-model an MCS that moves with spin (2), then we get the rate formulation of the LR-model type, but with a corotational derivative of the Green–Naghdi type (the GN-model used, for instance, in [75,83,84]). In line with this assumption, the proposed model will give stresses close to those determined in the U-model.

In [75], the following relation between spins (2) and (9) is given:

$$
\overline{\boldsymbol{\omega}}\_{\mathcal{T}} = \overline{\boldsymbol{\omega}}\_{\mathcal{T}} - (\mathcal{B} : \boldsymbol{\sigma}) \cdot \boldsymbol{d}^p + \boldsymbol{d}^p \cdot (\mathcal{B} : \boldsymbol{\sigma}), \tag{28}
$$

where *B* is the four rank elastic compliance tensor, *d<sup>p</sup>* is defined using Formula (8). According to (28), by virtue of the smallness of elastic deformations (*B* : *σ*), the use of spin Relations (2) and (9) in the mesolevel model will cause slightly different stresses to occur (in other words, the proximity of the T-model and the GN-model is their characteristic feature).

Hence, we can suggest that at small elastic distortions, all the models examined should give a close response (stresses). Table 1 presents information obtained when performing the analytical comparison of determination of stresses within the framework of different models.

In order to confirm the analytical results summarized in Table 1, we performed numerical calculations.

**Table 1.** Comparison of determination of stresses within the framework of different models (notation used here: «=» precise coincidence of stresse obtained in the indicated models, «≈»—stresses are close; by "described" is meant the description of the essence of calculations with references to publications).


#### *4.2. Illustrative Numerical Examples*

We now turn to the consideration of the results obtained via the analysis of the application of kinematic relations and constitutive laws in the statistical constitutive model for simulating some kinematic loads imposed on the fcc polycrystal.

The representative volume included a sample of 343 crystallites, the initial orientations of which were distributed randomly according to a uniform law. The nominal properties of the polycrystal corresponded to those of copper. The mesolevel elastic property tensor contained the following independent components (constant for the observer in a rigid moving coordinate system linked with the lattice): *π*<sup>1111</sup> = 168.4 GPa, *π*<sup>1122</sup> = 121.4 GPa, *<sup>π</sup>*<sup>1212</sup> <sup>=</sup> 75.4 GPa [85]; the viscoplastic relations included . *γ*<sup>0</sup> = 0.001 s–1, 1/*m* = 0.012; the hardening law parameters were *h*<sup>0</sup> = 180 MPa, *a* = 2.25, *τsat* = 148 MPa, and the initial values of critical stresses for all slip systems were *<sup>τ</sup>*(*k*) *<sup>c</sup>* (0) = *<sup>τ</sup>c*<sup>0</sup> <sup>=</sup> 16 MPa (*k* = 1, ... , *K*) [46,47]. Due to the significant nonlinearity of the systems of equations for all considered models, in particular, because of the presence of the Heaviside function in the viscoplastic law (5), which led to the necessity of discretization with a small time step to trace the activity of slip systems in crystallites, these equations were integrated using an explicit Euler method (in the calculations, the time step was 0.002 s). The results obtained with the aid of the LR1 model correspond satisfactorily to the experimental data for uniaxial compression and shear; the comparison of the results is given in [63].

We consider the affine deformations of the sample (corresponding to the representative volume of a polycrystal) subjected to uniformly distributed kinematic loads. Due to the uniform deformation, the radius vector of the material point at an arbitrary instant of time *t* is determined according to *r*(*t*) = *f*(*t*) · *q<sup>i</sup> pi* , where *f*(*t*) is the deformation gradient, *q<sup>i</sup>* stands for the Lagrangian coordinates of the considered point, and *p<sup>i</sup>* is the basis of the fixed laboratory coordinate system. The motion is determined by the deformation gradient (the choice of motion can be arbitrary):

$$\begin{array}{ll} f(t) = & \frac{1}{1 + 4h\sin(\omega t)} p\_1 p\_1 + (1 + 4h\sin(\omega t)) p\_2 p\_2 + p\_3 p\_3 + \\ & + 3h(1 - \cos(0.7\omega t)) p\_1 p\_2. \end{array} \tag{29}$$

where *ω* = 0.001π, *h* = 0.1 is the constant parameter (deformation is considered within the time interval *t* = [0, *T*], *T* = 1000 s). For illustration, the trajectory of the point with Lagrangian coordinates (0, L, 0) *m* is given in Figure 1.

**Figure 1.** Trajectory of the point with Lagrangian coordinates (0, L, 0) for Motion (29) (in the plane *OX1X2*); the position of the point at the initial instant of time is defined by the coordinates *X1/L* = 0, *X2/L* = 1, *X3* = 0.

Figure 2 presents the time dependence of the Cauchy stress tensor components at a macrolevel (in the laboratory coordinate system) for the LU-model (Figure 2a, precise coincidence with the results of the U-model) and for LR1-model (Figure 2b).

**Figure 2.** Time dependence of the Cauchy stress tensor components at a macrolevel (averaged mesostresses, **Σ** = *σ*) in the laboratory coordinate system for the fcc-crystal: (**a**) LU-model, U-model, (**b**) LR1-model.

We note the proximity of the results obtained (the curves plotted for the varying macrostress components are visually almost indistinguishable for all models examined in this study). For numerical assessment of the deviation in the results, we take the U-model as a base model and introduce the norm for comparing the results obtained using other models:

$$\Delta\_G = \left\| \Sigma\_{t \in [0, T]} \, ^{(II)} - \Sigma\_{t \in [0, T]} \, ^{(G)} \right\|\_{\mathcal{C}^{n}\_{L^2}} \, ^{\prime} \tag{30}$$

where **Σ***t*∈[0,*T*] (*U*) is the history of changes in stresses at *t* = [0, 1000 s], obtained via the use of the U-model, **Σ***t*∈[0,*T*] (*G*) is the history of changes in stresses at *t* = [0, 1000 s] in case of the G-model (one of the models from Table 1). The norm is given based on the Riemann integral *Yt*∈[0,*T*]*C<sup>n</sup> L*2 = 1 *T* (*T* 0 *<sup>Y</sup>*(*t*) : *<sup>Y</sup>T*(*t*)*dt*1/2 in the space *C<sup>n</sup> <sup>L</sup>*<sup>2</sup> continuous at

*t* ∈ [0, *T*] on the vector-function dimension *n*; in the calculations *n* = 9 (the number of the macrostress tensor component **Σ** that is equal to averaged mesostresses in crystallites).

Table 2 summarizes the results for deviations in the macrostresses obtained by the U-model from those determined by other models.

**Table 2.** Deviation of the history of changes in macrostresses obtained via the use of the U-model from the history of changes in macrostresses determined in the framework of other models.


(Zero deviations from Table 1 suggest that they do not exceed the order of the computer error when working with real numbers.).

Taking into account that the norms are integral over a sufficiently long period of time (T = 1000 s), the results summarized in Table 2 demonstrate the proximity of stresses obtained via the use of different formulations of mesolevel models, supporting thus the efficacy of theoretical statements given above (Table 1). Similar textures are obtained, which is indirectly confirmed by the proximity of macrostresses for the considered complex loading path. In [80], the issue of comparing the results of using some different spins (for other complex loading paths) was considered with an explicit analysis of the simulated lattice rotations.

Figure 3 shows the time dependence of the mesolevel stress components for a single crystallite (chosen arbitrarily). The initial orientations of the crystallographic and moving coordinate systems coincide and can be determined through the sequential rotation of the crystallographic system initially coincident with the coordinate system for the fcccrystal around the Ox1 axis with an angle *φ*1= 4.078, around the Ox<sup>2</sup> axis with an angle *φ*2= 0.0739, around the Ox<sup>3</sup> axis with an angle *φ*3= 2.149 for the U-model (Figure 3a) and the LR1-model (Figure 3b).

**Figure 3.** Time dependence of the Cauchy stress tensor components at a mesolevel in the laboratory coordinate system for crystallite No1: (**a**) LU-model, U-model, (**b**) LR2-model.

Analysis of the results shown in Figure 3 demonstrates that the mesostresses for the crystallite under study, which were obtained within the framework of the models mentioned above, are similar. Table 3 gives the norms for deviations in the mesostresses *δ<sup>G</sup>* = *σt*∈[0,*T*] (*U*) − *σt*∈[0,*T*] (*G*)*C<sup>n</sup> L*2 for the considered models, which are similar to those for deviations in macrostresses (30).


**Table 3.** Deviation of the history of changes in mesostresses obtained via the use of the U-model for crystallite No1 from the history of changes in macrostresses determined within in the framework of other models.

(Zero deviations from Table 1 suggest that they do not exceed the order of the computer error when working with real numbers.).

The results shown in Figure 3 and in Table 3 indicate the proximity of mesostresses. This is typical for the overwhelming number of crystallites (i.e., for almost all random orientations). In special cases, deviations of the mesostress components (difference modulus) may, at a certain instant of time, exceed 50 MPa. This can be attributed to some specific features in the motion of a representative point in the stress space near the yield surface [86] in transitions between the neighborhoods of the vertices. After transitions, in a great number of crystallites, the stresses become similar (there is a separate publication devoted to the study of these special cases). Due to the small fraction of such crystallites (37 crystallites out of 343) and averaging, the macroscale stresses obtained on the basis of the statistical constitutive model data using mesolevel models differ insignificantly (Table 2).

The differences in stresses obtained at the macrolevel are significantly lower than the experimental data scatter for samples from one production batch. Indeed, for each material, there is no complete identity of samples even from the same batch due to the presence of many random factors during manufacturing, which leads to differences in the microstructure (distributions of crystallite orientations, dislocation densities, etc.) and, as a consequence, to different experimental dependences of stresses on deformations. For example, such scatters for copper alloy samples are given in [87,88].

It should be noted here that the U-model and LU-model offer the most promise for modeling purely elastic behavior (e.g., elastic cyclic loading during exploitation) since the constitutive law is presented in the finite form, and there is no stress hysteresis in closed elastic cycles. Moreover, due to the use of energetically conjugate stress and strain measures, no energy is dissipated during the elastic deformation.

On the other hand, crystal plasticity models are primarily used to describe large inelastic deformations and material structure evolution. In this regard, formulations in terms of the actual configuration hold more promise because the physical meaning of the stress tensor is clearer in this case, and, consequently, the simulation of hardening using evolutionary equations is facilitated. These arguments lend additional support to using this formulation in case of the boundary value problem and in some advanced constitutive models at different deformation mechanisms. The theoretical predictions and the results of the numerical experiments indicate that the response produced by these models is similar to that obtained using the U-model. So, it turns out that the formulation in rate form would be most suitable for constructing advanced constitutive models (e.g., with complex hardening laws). The point is that the hypotheses adopted for constructing the evolutionary equations of the model via the use of the internal variables determined in the unloaded configuration can give more significant errors than the small deviations from the strict fulfillment of conservative conditions for the formulations written in terms of the actual configuration.

Geometrically nonlinear crystal plasticity models are applied to anisotropic crystallites, and therefore the decomposition of motion, used in rate formulations in terms of the current configuration, must be constructed by taking into account the symmetry of the material. In other words, the moving coordinate system determining the rigid motion in motion decomposition must be associated with the elements of material symmetry so that its symmetry properties can be taken into account correctly. Above, some options for constructing constitutive models with the introduction of such moving coordinate systems have been considered.
