# Microstructural Dependence and Reduction of the Size of the Representative Volume Element in Polycrystals: Case of Cubic and Hexagonal Elasticity

Romain Quey

► To cite this version:

Romain Quey. Microstructural Dependence and Reduction of the Size of the Representative Volume Element in Polycrystals: Case of Cubic and Hexagonal Elasticity. 2025. hal-05206608v2

HAL Id: hal-05206608

https://cnrs.hal.science/hal-05206608v2

Preprint submitted on 12 Nov 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Microstructural Dependence and Reduction of the Size of the Representative Volume Element in Polycrystals: Case of Cubic and Hexagonal Elasticity

Romain Quey¹

¹Mines Saint-Etienne, Univ Lyon, CNRS, UMR 5307 LGF, F - 42023 Saint-Etienne France

## Abstract

Closed-form expressions of the size of the representative volume element are determined for computing the elastic moduli of polycrystalline materials with cubic and hexagonal crystal symmetries. The size of the representative volume element corresponds to the number of grains needed in a polycrystal, *N*, to compute a property of the material (an elastic modulus) with a given error, *e*. The elastic moduli are computed using periodic polycrystals and finite element simulations. A typical, experimental grain-growth microstructure is considered as a reference case, and the microstructure properties are then modified. This allows us to determine the effect of the different microstructural properties on the error (*e*). The relationship between the error (*e*) and the size of the representative volume element (*N*) is shown to be *e* = *a*/√*N*, where *a* depends on both the crystal anisotropy of the material and the variability of the effective orientation distribution of the microstructure. This variability in turn arises from the cumulative effects of two main sources: the random sampling of the orientation distribution and the grain size distribution. Then, general expressions of *a* applying to all cubic and hexagonal-symmetry materials are determined. It is demonstrated that reducing the microstructural variability does not affect the computed elastic modulus values, thereby enabling the use of polycrystals with significantly fewer grains. For copper, the error in the computed shear modulus (*e*) can be reduced by a factor of 9.4 (for a given number of grains in a polycrystal, *N*), or equivalently, the number of grains needed in a polycrystal (*N*) can be reduced by a factor of about 90 (for a given error, *e*). Similar results are obtained for (hexagonal-symmetry) titanium and zinc, for which *N* can be reduced by factors of 140 and 45, respectively. The results also indicate that the size of the representative volume element (*N*) is smaller for hexagonal-symmetry materials than for cubic-symmetry materials, and smaller for the bulk modulus than for the shear modulus.

## 1 Introduction

The concept of “representative volume element (RVE)” plays a crucial role in the mechanics of random heterogeneous media such as polycrystalline materials. In the absence of crystallographic and morphological texture, the elastic response of a polycrystalline material is linear, isotropic, and can be characterized by two independent moduli: the shear modulus, *μ*, and the bulk modulus, *κ*. These values can then be used in upper-scale simulations if there is scale
separation, i.e. if the physical size of the RVE is smaller than the characteristic length of the structure. The concept of representative volume element is closely related to the notion of variability, which is an inherent property of the microstructures of polycrystalline materials. A typical polycrystalline microstructure is the grain-growth microstructure, which is found in a material that underwent grain boundary migration by capillarity effects. In addition to a uniform orientation distribution, a grain-growth microstructure is characterized by well-known morphological statistics. These include the grain size and shape distributions [1]. Actual polycrystals, which come from distinct locations inside a microstructure or from different (statistically equivalent) microstructures, have a finite numbers of grains, and vary in both grain morphology (and arrangement) and crystal orientations.<fnref n="1" /> Crystal orientations typically sample from the uniform orientation distribution through a random process. This is a first source of microstructural variability, since random sampling naturally leads to orientation clusters. Thus, some orientations (among all possible orientations) will be represented more than others. The varying grain volumes represent the second source of variability related to orientations, as the orientation of a larger grain will be overrepresented (in terms of volume fraction) relative to that of a smaller grain. In elastic polycrystals, the intrinsic mechanical response of a grain depends on its crystallographic orientation, while its contribution to the macroscopic response is proportional to its volume.<fnref n="2" /> Consequently, the variabilities in grain orientation and volume are the primary contributors to the overall variability of the effective elastic properties. This is particularly true when the number of grains in the polycrystal is small, as deviations from the theoretical uniform orientation distribution become more pronounced.

A quantitative definition of the representative volume element was provided by Kanit et al. in the context of heterogeneous (or “composite”) materials [2], but also previously discussed by, e.g., Huet [3]. Quoting Kanit et al. [2], “A RVE size can be associated with a given precision on the estimation of the wanted overall property [...]. It is shown to depend on the investigated morphological or physical property [...]”. Precision for a given property can be formally defined from the variability of measurements on different polycrystals, based on the theory of statistics [4]. Given a set of $n$ uncorrelated measurements, $\{z_i, i = 1, ..., n\}$, of mean $\bar{z}$ and standard deviation $\sigma_z$, and assuming normally distributed values (which is very generally the case of averaged properties, by virtue of the central limit theorem [4]), it is a common practice to define the (relative) error on a single measurement as $e = 2\sigma_z/\bar{z}$, which provides 95% confidence. Considering a specific property (shear modulus, bulk modulus, etc.), and aiming for a desired level of precision, the evolution of $\sigma_z$ as a function of the number of grains in a polycrystal can be determined through experiments or simulations. It is then possible to identify the number of grains that effectively represents the RVE size. Conversely, in a structural design approach, the error ($e$) can be interpreted as an uncertainty in the material’s response at the local scale. In Ref. [2], Kanit et al. also showed that periodic boundary conditions allow to reduce the RVE size by avoiding boundary effects. They then showed that the size of a representative volume element (RVE) associated to a given property is related to the integral range of the field associated with that property over the microstructure [2]. Further works discussed the existence of the representative volume element for specific composite microstructures [5], or concerned, more generally, the development of new numerical methods to compute

<footnote n="1">¹In this paper, “microstructure” is used in a general sense and is described by well-defined, “theoretical” statistical properties. Conversely, “polycrystal” represents a particular volume (with a finite number of grains) of a microstructure, whose statistical properties may deviate from the theoretical properties.</footnote>

<footnote n="2">²This is distinct from plasticity, where grain size has an intrinsic effect through the Hall–Petch relationship.</footnote>
the mechanical response [6–8]. The relationship between the RVE size of a mechanical property (e.g. an elastic modulus) and the underlying microstructural properties, as well as the potential to reduce the RVE size by tailoring these properties, has received comparatively little attention. Some works have focused on developing reduced or “compact” descriptions of (non random) crystallographic textures, by using a limited number of components [9, 10], especially in the context of machine learning [11, 12]. In previous works, we presented methods specifically focused on the representation of polycrystals, either for (i) fine control on the grain morphology (and its variability) [1] or (ii) (nearly) uniform sampling of orientation distributions [13]. The latter method eliminates the orientation clustering resulting from random sampling, ensuring that all possible orientations are represented equally. These methods enable the systematic variation of microstructural properties, ultimately allowing the effects on the predicted elastic moduli and RVE size to be assessed.

In this work, (single-phase) cubic and hexagonal crystal symmetry materials with different microstructures are considered. The microstructures exhibit no morphological or crystallographic texture, i.e. are (statistically) composed of equiaxed grains having a uniform orientation distribution. Their elastic moduli are computed using periodic polycrystals and finite element simulations. An expression of the RVE size as a function of the microstructural (variability) properties and the elastic anisotropy properties is first determined. It is then shown how, by reducing the variability of the microstructure (compared to a grain-growth microstructure), the RVE size can be significantly reduced while ensuring that no bias is introduced in the computed elastic modulus. The paper is organized as follows. Section 2 first outlines the methodology used to generate polycrystalline microstructures with controlled morphological and orientation properties. Section 3 then describes the simulation campaign and the different microstructures considered. Next, Section 4 presents the results and relates them to the microstructural properties. Section 5 finally closes the paper with conclusions.

## 2 Polycrystal representation

The grain structures and the grain orientations are generated independently, so that the morphological and topological properties of a grain, such as its volume and position, and its orientation are independent. Grain structures are obtained by optimization of Laguerre tessellations, providing complete control over the actual grain properties and particularly the variability in grain size, as detailed in Section 2.1. Grain orientations are obtained via (standard) random sampling or nearly uniform sampling, as described in Section 2.2. The methods are implemented in the free / open-source Neper software package [1, 13–15].

### 2.1 Optimization of Laguerre tessellations

Mathematically, a Laguerre tessellation is a partition of an $nD$ space in a collection of $nD$ entities that fill the space with no overlaps and no gaps. These entities are polyhedra and are formally defined as zones of influence of a particular set of weighted points, or seeds. Each seed, $S_i$, is defined by its coordinates, $x_i$, and a (non-negative) weight, $w_i$. Being given a spatial domain $D \in \mathbb{R}^n$ and a set of seeds $\{S_i(x_i, w_i), i = 1, \dots, N\}$, every seed $S_i$ is associated a Laguerre cell $C_i$ as follows,

$$C_i = \{ P(x) \in D \mid \|P, S_i\|^2 - w_i < \|P, S_j\|^2 - w_j \quad \forall j \neq i \}, \quad (1)$$
where “|||” denotes the Euclidean norm. In general, the larger the weight of a seed (relative to the weights of other seeds), the bigger its associated cell. Laguerre tessellations are normal tessellations, i.e. tessellations comprised of convex cells intersecting along planar faces, straight edges, and vertices [16].<fnref n="3" /> They can describe any normal tessellation, regardless of the size or arrangement of the cells. The Laguerre tessellation parameterization can therefore be considered as a general, continuous parameterization that allows any 3D tessellation composed of N (convex) cells to be described by only 4N (seed) parameters.<fnref n="4" /> Figure 1 provides an example of a (periodic) 3D Laguerre tessellation generated from seeds of random positions and weights. Note that, for a periodic tessellation, some fraction of a cell can be outside D, but it would fit inside D as a result of periodicity. The advantage of allowing a portion of a cell to extend beyond D, rather than dividing the cell into smaller parts that fit entirely inside D, is to avoid additional constraints for meshing. While Laguerre tessellations such as the one illustrated in Figure 1 properly reproduce the random character of polycrystalline microstructures, fine control of their morphological properties requires a specific method to set their seed coordinates and weights.

In this work, Laguerre tessellations are generated from prescribed cell properties, using an optimization approach. The coordinates and weights of the seeds represent the variables, and the objective function (to minimize) is defined as the difference between the prescribed properties and the current properties of the tessellation. The properties will be defined in terms of statistical distributions, but they can also be defined for individual cells in other applications [1]. The statistical description of polycrystalline morphologies first includes the (relative) grain size distribution. The (relative) size of a grain can be defined as the diameter of its sphere of equivalent volume, d, divided by the average diameter over all grains, $\bar{d}$. In the following, $d/\bar{d}$ will simply be referred to as “diameter” (or “size”). Generally, it is also important to describe the shapes of the grains, as they may influence the mechanical response of a polycrystal (and optimization could otherwise lead to unrealistic morphologies). The sphericity, s, is used, which is defined as the ratio between the surface area of the sphere of equivalent volume and the surface area of the grain, and takes a maximal value of 1 for a spherical grain [18]. For a grain-growth microstructure, $d/\bar{d}$ generally follows a lognormal distribution of average 1 (by definition) and standard deviation 0.35, and $1-s$ follows a lognormal distribution of average 0.145 and standard deviation 0.03 [18–21]. For finite element simulations, it is also desirable to avoid the presence of the very small edges found in standard random tessellations [14]. In this work, this is also achieved in the optimization (and not using “regularization” as introduced in Ref. [14]).

The size and sphericity metrics being defined, their contributions to the objective function (via the differences between their prescribed and current distributions) are computed separately. Let F denote a prescribed cumulative distribution function, and X denote the cumulative distribution function of the discrete cell values, {$x_i$}, $i \in [1, N]$. F is a continuous function while X is a step function which can be written as

$$X(x) = \frac{1}{N} \sum_{i=1}^{N} H(x - x_i) \quad \text{with} \quad H(\bullet) = \begin{cases} 0 & \text{if } \bullet < 0 \\ 1 & \text{if } \bullet \ge 0 \end{cases}. \qquad (2)$$

<footnote n="3">³This is topologically equivalent to grains of a polycrystal meeting along grain boundaries, triple lines and quadruple points.</footnote>

<footnote n="4">⁴It is even possible to remove one parameter, since all weights are defined up to a constant [17].</footnote>
$H$ is the Heaviside step function. $F$ and $X$ are compared in terms of a metric, $A^s$, defined as,

$$A^s = \int_{-\infty}^{+\infty} \frac{[X^s(x) - F^s(x)]^2}{F^s(x) [1 - F^s(x)]} dx \quad \text{with} \quad \begin{cases} F^s = F \circ s \\ X^s = X \circ s \end{cases}, \qquad (3)$$

where “◦” is the convolution product and $s$ is a normal convolution kernel of standard deviation $\sigma_s$. The role of $s$ is to facilitate the convergence during optimization, and we use kernels of standard deviations equal to 0.05 for $d/d$ and 0.005 for $s$. Minimizing $A^s$ ensures that the distribution $X$ (and the values of $\{x_i\}$) converges toward the prescribed distribution $F$. Given a target minimal edge length, $l_r$, its contribution to the objective function, $B_r$, can be computed directly from the cell values:

$$B_r = \sqrt{\sum_{i=1}^{N} [l_i^{\min} - l_r]_+^2}, \qquad (4)$$

where $l_i^{\min}$ is the minimum length of the edges of cell $i$. A small value of $l_r = 0.02 \bar{d}$ is used, so that the overall grain morphology is almost unmodified.<fnref n="5" /> A last but important point is to avoid empty cells, which is done by means of a penalty term. The objective function, $O$, is finally computed from the $A^s$ values of the different variables, $A_d^s$ for the size and $A_s^s$ for the sphericity, the small edge component ($B_r$), and the penalty term, as

$$O = \sqrt{A_d^{s2} + A_s^{s2} + B_r^2} + p N_e, \qquad (5)$$

where $N_e$ is the number of empty cells and $p = 10^3$. The Subplex (gradient-free) algorithm [22], as available in the NLopt free (open-source) library [23], associated with a “tessellation update” algorithm were used to minimize $O$ [1]. Optimization terminates when

$$\Delta O_{\min} < 10^{-6}, \qquad (6)$$

where $\Delta O_{\min}$ is the decrease of the minimal value of $O$ evaluated over the last 40 $N$ iterations. The microstructures effectively used in this work will be illustrated in Section 3.

## 2.2 Crystal orientation distribution

The description and manipulation of crystal orientations require a parameterization, which is presented in Section 2.2.1. The method for random sampling is straightforward and is outlined within Section 2.2.1, while the method for nearly uniform sampling is detailed in Section 2.2.2. A measure of the variability of the resulting orientation distributions is defined in Section 2.2.3.

### 2.2.1 Orientations represented as unit quaternions

A crystal orientation can be described by a unit quaternion, $q$, defined as

$$q = \cos \frac{\theta}{2} + r \sin \frac{\theta}{2}, \qquad (7)$$

where $r$ ($\|r\| = 1$) and $\theta \in [0, \pi]$ are the (uniquely defined) rotation axis and angle from the reference (sample) coordinate system to the crystal coordinate system, respectively [24]. A unit quaternion is therefore composed of a scalar part and a vectorial part. A unit quaternion

<footnote n="5">⁵In practice, optimization is run first without the $B_r$ component, and then including it.</footnote>
<figure><img src="image_2.png" /><figcaption>Figure 1: Example of a 3D periodic Laguerre tessellation in (a) perspective view and (b) front view. The wireframe represents the cubic “domain” (D). Cells are arbitrarily colored.</figcaption></figure>

can equivalently be described as the 4D vector $(\cos(\theta/2), r_1 \sin(\theta/2), r_2 \sin(\theta/2), r_3 \sin(\theta/2))$, where $r_1, r_2$ and $r_3$ are the three components of $\mathbf{r}$. Unit quaternions lie on the unit 3-sphere, $\mathbb{S}^3$. It can be noted that opposite quaternions, $\mathbf{q}$ and $-\mathbf{q}$, which are antipodal points on $\mathbb{S}^3$, correspond to rotations of $(\theta, \mathbf{r})$ and $(\theta+2\pi, \mathbf{r})$, respectively, and so represent the same orientation. A key interest in unit quaternions lies in the proportionality between their intrinsic distance and the angle between their corresponding orientations [25]. It follows that a statistically uniform distribution of 3D orientations corresponds to a statistically uniform distribution of points on $\mathbb{S}^3$. These properties make random sampling of orientations particularly straightforward, since it reduces to randomly distributing points on $\mathbb{S}^3$. To distribute the points, a rejection algorithm is used. For each orientation, four numbers, $v_1, v_2, v_3$ and $v_4$, are drawn randomly in $[-1, 1]$. If $v = \sqrt{v_1^2 + v_2^2 + v_3^2 + v_4^2} \le 1$, the orientation is defined by the unit quaternion $(v_1, v_2, v_3, v_4)/v$. Otherwise, new numbers are drawn.

### 2.2.2 Nearly uniform orientation sampling

A nearly uniform set of orientations can be obtained by distributing points nearly uniformly on $\mathbb{S}^3$. To generate such a set starting from randomly distributed orientations, unit quaternions (or “points”) are considered as charged particles, and the energy of the system of particles is minimized, similarly to the Thomson problem for electrons on the 3D unit sphere [26].

Let us consider a set of two orientations, $\mathbf{q}_i$ and $\mathbf{q}_j$. The force at $\mathbf{q}_i$ due to the interaction with orientation $j$, $\mathbf{f}_i$, can be expressed in tangent space at orientation $\mathbf{q}_i$ as

$$ \mathbf{f}_i = - \left( \frac{1}{d_{ij}^2} - \frac{1}{(\pi - d_{ij})^2} \right) \mathbf{r}_{ij}, \quad \text{with } d_{ij} = 2 \theta_{ij}, \qquad (8) $$

where $d_{ij}$ represents the intrinsic distance between $\mathbf{q}_i$ and $\mathbf{q}_j$, and $\theta_{ij}$ and $\mathbf{r}_{ij}$ represent the misorientation angle and the misorientation axis between $\mathbf{q}_i$ and $\mathbf{q}_j$, respectively. The expression of $\mathbf{f}_i$ (Equation 8) includes two terms, which come from the interaction of $\mathbf{q}_i$ with the two quaternions of orientation $j$ ($\mathbf{q}_j$ and $-\mathbf{q}_j$). Additionally, $\mathbf{f}_i$ being aligned with $\mathbf{r}_{ij}$ indicates
that the two orientations would rotate about their misorientation axis ($r_{ij}$) as a result of the interaction. Equilibrium ($f_i = 0$) is reached for $d_{ij} = 2\pi$ or, equivalently, $\theta_{ij} = \pi$.

Equation 8 can easily be extended to account for crystal symmetry and interaction among a set of $N$ orientations ($N \ge 1$). First, in the presence of crystal symmetry, an orientation, represented by quaternion $q$, can be described by several, equivalent quaternions, $q^k = q u^k$, $k = 1, \dots, n_c$, where $u^k$ is an operator of the corresponding symmetry group [25]. Quaternion $q_i$ therefore interacts with all equivalent quaternions of orientation $j$. The generalization to $N$ orientations similarly involves a summation of the forces. We finally get

$$f_i = - \sum_{j=1, j \neq i}^{N} \sum_{k=1}^{n_c} \left( \frac{1}{d_{ij}^k} - \frac{1}{(\pi - d_{ij}^k)^2} \right) r_{ij}^k, \quad \text{with } d_{ij}^k = 2 \theta_{ij}^k. \qquad (9)$$

Knowing the forces ($f_i$), orientations are evolved iteratively until equilibrium is reached. During an iteration ($l$), each orientation evolves from its initial position, $q_i^{(l-1)}$, to its final position, $q_i^{(l)}$. The rotation during the iteration, $\Delta q_i^{(l-1)}$, is computed from the force at the beginning of the iteration, $f_i^{(l-1)}$, as

$$\Delta q_i^{(l-1)} = \alpha^{(l-1)} f_i^{(l-1)}, \qquad (10)$$

where $\alpha^{(l-1)}$ is the step size of the iteration ($l$) as computed by the Barzilai-Borwein method [27]. The orientation evolution during an iteration can then be written as

$$q_i^{(l)} = \overline{q_i^{(l-1)} + \Delta q_i^{(l-1)} q_i^{(l-1)}}. \qquad (11)$$

where the composition of $\Delta q_i^{(l-1)}$ with $q_i^{(l-1)}$ corresponds to its parallel transport from tangent space to orientation $q_i^{(l-1)}$, and “$^\wedge$” represents the normalization of the vector. Finally, optimization terminates based on a relative error on the forces:

$$\frac{\sqrt{\sum_{i=1}^{N} \|f_i^{(l)}\|^2}}{\sqrt{\sum_{i=1}^{N} \|f_i^{(0)}\|^2}} < 10^{-3}. \qquad (12)$$

### 2.2.3 Variability of the orientation distribution

The distribution of a discrete set of weighted orientations, $\{(q_i, p_i), i = 1, \dots, N\}$, where $p_i$ are the orientation weights ($\sum_{i=1}^n p_i = 1$), is formally described by the orientation distribution function (ODF), $f$. The ODF is computed using a discrete method, for which each orientation ($q_i$) (and its symmetry equivalents ($q_i^k$)) is represented by a continuous, Bingham distribution [25]. The ODF can therefore be written as

$$f(q) = \frac{1}{n_c} \sum_{i=1}^{N} p_i \left( \sum_{k=1}^{n_c} \left( {}_1F_1 \left( \frac{1}{2}, 2, K \right) \right)^{-1} \exp \left( w_i^{kt} U K U^t w_i^k \right) \right) \quad \text{with} \quad w_i^k = q_i^k q^{-1}, \qquad (13)$$

where $U$ is an orthonormal $4 \times 4$ matrix, $K$ is a diagonal matrix describing the shape of the distributions, and ${}_1F_1$ is the hypergeometric function with matrix argument normalizing the density of the Bingham distribution. In this work, we consider Bingham distributions with
isotropic spreads of (unidirectional) standard deviation equal to $d_r$, as obtained using $U = I_4$, $K_{11} = 2/d_r^2$, and $K_{22} = K_{33} = K_{44} = 0$. $d_r$ is taken as the “average orientation radius” and is defined as the radius of the hyperspherical cap of average equivalent surface area (on $S^3$) [13]. Using $d_r$ as standard deviation provides a level of smoothing that is independent of the number of orientations ($N n_c$). The variability of an ODF, $\sigma_f$, is simply defined as its standard deviation over orientation space. The lower the variability of an orientation distribution, the closer to 1 its ODF ($f$), and the lower $\sigma_f$. For the theoretical, perfectly uniform orientation distribution ($f = 1$), $\sigma_f = 0$. For randomly sampled orientations, $\sigma_f$ is equal to 0.3, while for nearly uniformly sampled orientations, it is equal to 0.007.<fnref n="6" /> In this section (and specifically in Section 2.2.2), all orientations have the same weight ($p_i$), in which case the distribution is simply referred to as “orientation distribution (function)” ($f$, of standard deviation $\sigma_f$). In contrast, in the case of actual polycrystals (and for future reference), orientations will effectively be weighted according to the volumes of their grains, and their distribution will be referred to as “effective orientation distribution (function)”, $f^e$, of standard deviation $\sigma_f^e$.

Finally, an ODF is visualized in the 3D Rodrigues orientation space, where an orientation is represented as

$$ \boldsymbol{x}_R = \boldsymbol{r} \tan \frac{\theta}{2}. \qquad (14) $$

The “fundamental region” or the space is used, which corresponds to the region closer to the origin in which each orientation is defined only once (except on the boundary). For cubic and hexagonal symmetry, the fundamental region is a truncated cube and a dodecagonal prism, respectively. The ODFs corresponding to random and nearly uniform samplings are illustrated in Figure 2.

# 3 Application to cubic and hexagonal elasticity

Given a microstructure and a material, the elastic moduli are computed by finite element simulation. The crystal behavior is detailed in Section 3.1. Six types of microstructures, differing both in their grain size distribution and grain orientation sampling, are used, as described in Section 3.2. The procedure for the numerical computation of their elastic moduli is finally detailed in Section 3.3.

## 3.1 Single crystal elastic anisotropy

At every material point, the (Cauchy) stress tensor, $\boldsymbol{\sigma}$, is related to the (elastic) strain tensor, $\boldsymbol{\varepsilon}$, via Hooke's law:

$$ \boldsymbol{\sigma} = \boldsymbol{C} : \boldsymbol{\varepsilon}, \qquad (15) $$

where $\boldsymbol{C}$ is the (fourth-order) elastic stiffness tensor. Under the Voigt matrix notation, both the stress tensor ($\boldsymbol{\sigma}$) and the elastic strain tensor ($\boldsymbol{\varepsilon}$) can be expressed as $6 \times 1$ vectors, and $\boldsymbol{C}$ can be reduced to a $6 \times 6$ matrix.

<footnote n="6">⁶Since the values are derived from the actual (continuous) orientation distribution function ($f$), they can depend on the particular method used for its computation. The theoretical value for random sampling is actually $\sigma = \sqrt{1/12} \approx 0.289$, by analogy with a 1D uniformly distributed random variable on $[0, 1]$. The computed value (0.3) is close to this theoretical value.</footnote>
<figure><img src="image_3.png" /><figcaption>b</figcaption></figure>

<figure><img src="image_4.png" /><figcaption>Figure 2: Example orientation distribution functions represented over the fundamental region of Rodrigues orientation space, in the case of (a) cubic symmetry and (b) hexagonal symmetry. Random sampling appears on the left and nearly uniform sampling on the right.</figcaption></figure>
In the case of cubic symmetry, the C matrix contains only three independent elastic constants, $C_{11}$, $C_{12}$ and $C_{44}$, so that Equation 15 can be written as

$$ \begin{pmatrix} \sigma_{11} \\ \sigma_{22} \\ \sigma_{33} \\ \sigma_{23} \\ \sigma_{31} \\ \sigma_{12} \end{pmatrix} = \begin{pmatrix} C_{11} & C_{12} & C_{12} & 0 & 0 & 0 \\ C_{12} & C_{11} & C_{12} & 0 & 0 & 0 \\ C_{12} & C_{12} & C_{11} & 0 & 0 & 0 \\ 0 & 0 & 0 & C_{44} & 0 & 0 \\ 0 & 0 & 0 & 0 & C_{44} & 0 \\ 0 & 0 & 0 & 0 & 0 & C_{44} \end{pmatrix} \cdot \begin{pmatrix} \epsilon_{11} \\ \epsilon_{22} \\ \epsilon_{33} \\ 2 \epsilon_{23} \\ 2 \epsilon_{31} \\ 2 \epsilon_{12} \end{pmatrix}. \quad (16) $$

The degree of anisotropy is commonly described by the Zener ratio,

$$ A = \frac{2 C_{44}}{C_{11} - C_{12}}, \quad (17) $$

where a value of 1 denotes an isotropic behavior. It is also interesting to consider another metric, the universal anisotropy factor, $A_u$ [28]. $A_u$ can be computed from the Voigt and Reuss averages (see Appendix A) but is also related to A by

$$ A_u = \frac{6}{5} \left( \sqrt{A} - \frac{1}{\sqrt{A}} \right)^2. \quad (18) $$

By definition, $A_u \ge 0$, and a value of 0 denotes an isotropic behavior. Finally, and since A and $A_u$ are directly related through Equation 18 from the ratio $C_{44} / (C_{11} - C_{12})$, another metric is necessary to fully parameterize the ($C_{11}$, $C_{12}$, $C_{44}$) space. We use the ratio

$$ r = \frac{C_{44}}{C_{11}}. \quad (19) $$

So, given a value of A (or $A_u$) and r, relative values of $C_{11}$, $C_{12}$ and $C_{44}$ can be computed. The elastic constants and anisotropy ratios of different cubic-symmetry materials are provided in Table 1 [29]. From the values in Table 1 and for future reference, it can be seen that r shows a mild dependence on A among different materials, which can be approximated by

$$ r = 0.28 + 0.056 A. \quad (20) $$

In the case of hexagonal symmetry, the C matrix contains five independent elastic constants, $C_{11}$, $C_{12}$, $C_{13}$, $C_{33}$ and $C_{44}$, so that Equation 15 can be written as

$$ \begin{pmatrix} \sigma_{11} \\ \sigma_{22} \\ \sigma_{33} \\ \sigma_{23} \\ \sigma_{13} \\ \sigma_{12} \end{pmatrix} = \begin{pmatrix} C_{11} & C_{12} & C_{13} & 0 & 0 & 0 \\ C_{12} & C_{11} & C_{13} & 0 & 0 & 0 \\ C_{13} & C_{13} & C_{33} & 0 & 0 & 0 \\ 0 & 0 & 0 & C_{44} & 0 & 0 \\ 0 & 0 & 0 & 0 & C_{44} & 0 \\ 0 & 0 & 0 & 0 & 0 & C_{66} \end{pmatrix} \cdot \begin{pmatrix} \epsilon_{11} \\ \epsilon_{22} \\ \epsilon_{33} \\ 2 \epsilon_{23} \\ 2 \epsilon_{13} \\ 2 \epsilon_{12} \end{pmatrix}. \quad (21) $$

with $C_{66} = (C_{11} - C_{12})/2$. The degree of anisotropy can also be described by the universal anisotropy factor ($A_u$), as detailed in Appendix A. The elastic constants and anisotropy ratios of different hexagonal-symmetry materials are provided in Table 2. From the values, it can be seen that hexagonal-symmetry materials tend to have smaller elastic anisotropies compared to cubic-symmetry materials.
<table><caption>Table 1: Elastic constants ($C_{ij}$) [29,30] and anisotropy ratios of various cubic-symmetry materials, listed in increasing A values.</caption><thead><tr><th>Material</th><th>C<sub>11</sub> [GPa]</th><th>C<sub>12</sub> [GPa]</th><th>C<sub>44</sub> [GPa]</th><th>A [-]</th><th>A<sub>u</sub> [-]</th><th>r [-]</th></tr></thead><tbody><tr><td>Mo</td><td>460.0</td><td>176.0</td><td>100.0</td><td>0.8</td><td>0.1</td><td>0.22</td></tr><tr><td>W</td><td>501.0</td><td>198.0</td><td>151.4</td><td>1.0</td><td>0.0</td><td>0.30</td></tr><tr><td>Al</td><td>108.2</td><td>61.3</td><td>28.5</td><td>1.2</td><td>0.0</td><td>0.26</td></tr><tr><td>C (diamond)</td><td>1076.0</td><td>125.0</td><td>575.8</td><td>1.2</td><td>0.0</td><td>0.54</td></tr><tr><td>Si</td><td>165.7</td><td>63.9</td><td>79.6</td><td>1.6</td><td>0.2</td><td>0.48</td></tr><tr><td>Ge</td><td>128.9</td><td>48.3</td><td>67.1</td><td>1.7</td><td>0.3</td><td>0.52</td></tr><tr><td>Fe (α)</td><td>230.4</td><td>134.1</td><td>115.9</td><td>2.4</td><td>1.0</td><td>0.50</td></tr><tr><td>Ni</td><td>246.5</td><td>147.3</td><td>124.7</td><td>2.5</td><td>1.1</td><td>0.51</td></tr><tr><td>Au</td><td>186.0</td><td>157.0</td><td>42.0</td><td>2.9</td><td>1.5</td><td>0.23</td></tr><tr><td>Ag</td><td>124.0</td><td>93.4</td><td>46.1</td><td>3.0</td><td>1.6</td><td>0.37</td></tr><tr><td>Cu</td><td>168.4</td><td>121.4</td><td>75.4</td><td>3.2</td><td>1.8</td><td>0.45</td></tr><tr><td>Th</td><td>75.3</td><td>48.9</td><td>47.8</td><td>3.6</td><td>2.3</td><td>0.63</td></tr><tr><td>Pb</td><td>46.6</td><td>39.2</td><td>14.4</td><td>3.9</td><td>2.6</td><td>0.31</td></tr><tr><td>K</td><td>4.57</td><td>3.74</td><td>2.63</td><td>6.3</td><td>5.4</td><td>0.58</td></tr><tr><td>Na</td><td>6.03</td><td>4.59</td><td>5.86</td><td>8.1</td><td>7.5</td><td>0.97</td></tr><tr><td>Li</td><td>14.8</td><td>12.5</td><td>10.8</td><td>9.4</td><td>9.0</td><td>0.73</td></tr></tbody></table>

## 3.2 Microstructures and equivalent orientation distribution functions

To analyse the effect of the two factors contributing to the microstructure variability, namely the grain size distribution and the grain orientation sampling, several types of microstructures are used. First, three grain size distributions are considered, with standard deviations ($\sigma_d$) equal to 0, 0.155 and 0.35. The case for which $\sigma_d = 0$ corresponds to a delta Dirac grain size distribution, for which all grains have the same equivalent diameter (and volume), and will be referred to as “Dirac”. The case for which $\sigma_d = 0.35$ corresponds to a grain-growth microstructure and will be referred to as “grain-growth”. Finally, the intermediate case for which $\sigma_d = 0.155$ corresponds to the grain size distribution of a standard Voronoi tessellation, and will be referred to as “Voronoi”. This intermediate case is used to enrich the data set and because standard Voronoi tessellations, for which seeds are randomly distributed according to a Poisson point process and a uniform distribution over the domain, were standard in early studies [2, 32, 33]. As the grain sphericity distribution is not expected to influence the elastic moduli significantly, only the (experimentally representative) case for which $1-s$ follows a log-normal distribution of average 0.145 and standard deviation 0.03 is considered [14].<fnref n="7" /> Example polycrystals are provided in Figure 3. Then, for each grain size distribution, either random or nearly uniform grain orientation sampling is considered, which finally leads to a set of six microstructures. Among them, two are of particular interest. The association of $\sigma_d = 0.35$ to random orientation sampling corresponds to an experimentally representative grain-growth microstructure and will simply be referred to as “experimental” in the following. The association of $\sigma_d = 0$ to nearly uniform orientation sampling corresponds to the microstructure

<footnote n="7">⁷The “Voronoi” microstructure therefore do not exactly correspond to a Voronoi tessellation, which generally has lower cell sphericities [1].</footnote>
<table><caption>Table 2: Elastic constants ($C_{ij}$) [31] and anisotropy ratios of various hexagonal-symmetry materials, listed in increasing $A_u$ values.</caption><thead><tr><th>Material</th><th>C<sub>11</sub> [GPa]</th><th>C<sub>12</sub> [GPa]</th><th>C<sub>13</sub> [GPa]</th><th>C<sub>33</sub> [GPa]</th><th>C<sub>44</sub> [GPa]</th><th>A<sub>u</sub> [-]</th></tr></thead><tbody><tr><td>Zr</td><td>144.0</td><td>74.0</td><td>67.0</td><td>166.0</td><td>33.4</td><td>0.09</td></tr><tr><td>Ti</td><td>160.0</td><td>90.0</td><td>66.0</td><td>181.0</td><td>46.5</td><td>0.19</td></tr><tr><td>Be</td><td>292.0</td><td>24.0</td><td>6.0</td><td>349.0</td><td>163.0</td><td>0.05</td></tr><tr><td>Mg</td><td>59.3</td><td>25.7</td><td>21.4</td><td>61.5</td><td>16.4</td><td>0.04</td></tr><tr><td>Sc</td><td>99.3</td><td>39.7</td><td>29.4</td><td>107.0</td><td>27.7</td><td>0.08</td></tr><tr><td>Co</td><td>295.0</td><td>159.0</td><td>111.0</td><td>335.0</td><td>71.0</td><td>0.20</td></tr><tr><td>Zn</td><td>165.0</td><td>31.1</td><td>50.0</td><td>61.8</td><td>39.6</td><td>1.82</td></tr><tr><td>Y</td><td>77.9</td><td>29.2</td><td>20.0</td><td>76.9</td><td>24.3</td><td>0.04</td></tr><tr><td>Ru</td><td>563.0</td><td>188.0</td><td>168.0</td><td>624.0</td><td>181.0</td><td>0.03</td></tr><tr><td>Cd</td><td>116.0</td><td>43.0</td><td>41.0</td><td>50.9</td><td>19.6</td><td>1.16</td></tr><tr><td>La</td><td>51.4</td><td>17.3</td><td>10.4</td><td>54.6</td><td>13.9</td><td>0.17</td></tr><tr><td>Pr</td><td>49.4</td><td>23.0</td><td>14.3</td><td>57.4</td><td>13.6</td><td>0.19</td></tr><tr><td>Nd</td><td>54.8</td><td>24.6</td><td>16.6</td><td>60.9</td><td>15.0</td><td>0.13</td></tr><tr><td>Gd</td><td>67.2</td><td>25.3</td><td>21.0</td><td>71.5</td><td>20.8</td><td>0.03</td></tr><tr><td>Tb</td><td>68.5</td><td>24.6</td><td>22.4</td><td>73.3</td><td>21.6</td><td>0.02</td></tr><tr><td>Dy</td><td>74.0</td><td>25.5</td><td>21.8</td><td>78.6</td><td>34.3</td><td>0.12</td></tr><tr><td>Ho</td><td>76.5</td><td>25.6</td><td>21.0</td><td>79.6</td><td>25.9</td><td>0.02</td></tr><tr><td>Er</td><td>84.1</td><td>29.4</td><td>22.6</td><td>84.7</td><td>27.4</td><td>0.02</td></tr><tr><td>Tm</td><td>92.5</td><td>33.5</td><td>21.0</td><td>81.5</td><td>28.2</td><td>0.05</td></tr><tr><td>Lu</td><td>86.2</td><td>32.0</td><td>28.0</td><td>80.9</td><td>26.8</td><td>0.01</td></tr><tr><td>Hf</td><td>181.0</td><td>77.0</td><td>66.0</td><td>197.0</td><td>55.7</td><td>0.03</td></tr><tr><td>Re</td><td>616.0</td><td>273.0</td><td>206.0</td><td>683.0</td><td>161.0</td><td>0.11</td></tr><tr><td>Os</td><td>765.0</td><td>229.0</td><td>219.0</td><td>846.0</td><td>270.0</td><td>0.01</td></tr><tr><td>Tl</td><td>41.4</td><td>36.0</td><td>29.4</td><td>53.9</td><td>7.2</td><td>1.88</td></tr></tbody></table>
with the least variability and will be referred to as “fully uniform” in the following. Other microstructures have intermediate variabilities and will be collectively referred to as “intermediate microstructures” in the following. They will allow to determine the relative effects of the grain size distribution and the grain orientation sampling. The variabilities of the different microstructures ($\sigma_d$, $\sigma_f$ and $\sigma_f^e$) are provided in Table 3. For example, for the experimental microstructure, $\sigma_d = 0.35$ and $\sigma_f = 0.304$ lead to $\sigma_f^e = 0.5$.

For each microstructure, polycrystals with different numbers of grains ($N$) are considered, from 3 to 1000. The special case of $N = 1$ is also considered for completeness and corresponds to single crystals, for which all microstructures are equivalent. Finally, for each value of $N$, $n$ different polycrystals are considered. This will allow to evaluate both the average ($\bar{z}$) and the standard deviation ($\sigma_z$) of the elastic moduli for a given microstructure and polycrystals with a given number of grains ($N$). For each value of $N$, the value of $n$ was chosen to ensure that the maximal standard error on both the mean and standard deviation of the elastic moduli computed from the polycrystals did not exceed 0.1%. This lead to values ranging from $n = 6500$ for $N = 3$ (or even $n = 9500$ for $N = 1$) to $n = 950$ for $N = 1000$.

<table><caption>Table 3: Variabilities of the microstructures. For each microstructure, the values represent averages over all polycrystals.</caption><thead><tr><th>Microstructure</th><th>σ<sub>d</sub> [-]</th><th>σ<sub>f</sub> [-]</th><th>σ<sub>f</sub><sup>e</sup> [-]</th></tr></thead><tbody><tr><td>Grain-growth – random</td><td>0.350</td><td>0.304</td><td>0.500</td></tr><tr><td>Grain-growth – uniform</td><td>0.350</td><td>0.007</td><td>0.397</td></tr><tr><td>Voronoi – random</td><td>0.155</td><td>0.304</td><td>0.334</td></tr><tr><td>Voronoi – uniform</td><td>0.155</td><td>0.007</td><td>0.137</td></tr><tr><td>Dirac – random</td><td>0</td><td>0.304</td><td>0.304</td></tr><tr><td>Dirac – uniform</td><td>0</td><td>0.007</td><td>0.007</td></tr></tbody></table>

### 3.3 Finite element simulations

The polycrystals are meshed into second-order tetrahedral elements, so that all elements have approximately the same size, with an average number of 160 elements and 230 nodes per grain. The shear ($\mu$) and bulk ($\kappa$) moduli are computed by subjecting the polycrystal to pure shear deformation and hydrostatic deformation, respectively. For cubic-symmetry materials, the bulk modulus can actually be analytically computed from the elastic constants $C_{11}$ and $C_{12}$, and does not depend on the orientation or microstructure (see Appendix A). It will therefore only be considered for hexagonal-symmetry materials. In principle, different shear moduli could be computed in the three shear planes, but only the 1–2 plane is considered to avoid any potential coupling effects in the results. An average strain of $10^{-3}$ was applied to mimic an experimental elastic loading, and the elastic moduli were determined from the computed average stress. The mesh of a 1000-grain grain-growth polycrystal and its shear deformation are provided in Figure 4. Finite element simulations are carried out with the Z-set software package [34].
<figure><img src="image_5.png" /><figcaption>Figure 3: Polycrystals of prescribed grain size distributions: (a) “Dirac” polycrystal ($\sigma_d$ = 0), (b) “Voronoi” polycrystal ($\sigma_d$ = 0.155) and (c) “Grain-growth” polycrystal ($\sigma_d$ = 0.35). Case of polycrystals containing 1000 grains. The wireframe represents the cubic “domain” (D). Grains are arbitrarily colored. Note the increase in grain size variability when going from (a) to (c).</figcaption></figure>
a

<figure><img src="image_6.png" /></figure>

b

<figure><img src="image_7.png" /><figcaption>Figure 4: Shear deformation of a grain-growth polycrystal containing 1000 grains. (a) Finite-element mesh. (b) Strain field, scaled 100 times. The wireframe represents the initial cubic “domain” (D). On (a), grains are arbitrarily colored.</figcaption></figure>
# 4 Results

Results are first obtained for copper, which is a cubic-symmetry material with a relatively high elastic anisotropy. In Section 4.1, the first results on the shear modulus ($\mu$) are presented, and the approach is validated. In Section 4.2, the evolution of the error on $\mu$ (e, as introduced in Section 1) as a function of the number of grains in a polycrystal ($N$) is analysed. In Section 4.3, the dependency of the error (and the RVE size) on the microstructural variabilities is determined. The analysis is then generalized to all cubic-symmetry materials in Section 4.4 and all hexagonal-symmetry materials in Section 4.5. In the following, “variability” ($\sigma_d$, $\sigma_f$, $\sigma_f^e$) is used to describe microstructures, while “error” (e) denotes the error on their elastic moduli.<fnref n="8" /> The error is considered relative to the mean value unless mentioned otherwise.

## 4.1 Validation

Figure 5 provides the evolution of the statistical properties of the $\mu$ values as a function of the number of grains ($N$), for the different microstructures. The average values of $\mu$, $\bar{\mu}$, are provided as solid lines, and the absolute errors are shown as ($\pm 2\sigma$) envelopes about the mean. It is evident from the graph that the different microstructures yield significantly different errors, varying by approximately one order of magnitude for a given value of $N$. However, before analyzing the errors, it is important to verify that all microstructures produce similar mean values ($\bar{\mu}$). To do so while avoiding potential bias due to the finite size of the polycrystal, the optimal (nearly asymptotic) case of $N = 1000$ is considered. The $\bar{\mu}$ values of the different microstructures are provided in Table 4, together with their standard errors, $2\sigma_\mu/\sqrt{n}$, at 95% confidence [4]. Given the large number of polycrystals ($n = 950$), the standard errors have very low values of 0.01–0.1% $\bar{\mu}$, for all microstructures. For the experimental microstructure, $\bar{\mu}$ has a value of $48\,541 \pm 47$ MPa, while for the fully uniform microstructure, $\bar{\mu}$ has a very similar value of $48\,495 \pm 5$ MPa. Among all microstructures, the Voronoi ones seem to have slightly different $\bar{\mu}$ values, but the difference remains very small (0.17%). It also appears from Table 4 that, for a given morphology (Dirac, Voronoi or grain-growth), nearly uniform orientation sampling provides slightly smaller values compared to random orientation sampling. However, any of these small differences on the asymptotic values can be reasonably neglected. In fact, and as will be shown in the following, they are significantly smaller than the errors caused by microstructure variability.

In Figure 5, it can be seen that the $\bar{\mu}$ values deviate from the asymptotic value when the number of grains becomes too small. This is better illustrated in Figure 6, which provides a close-up view in the range $N \in [1, 100]$. The extent of the deviation varies depending on the microstructure. In the literature, such a deviation is typically attributed to the effects of the boundary conditions, when applied statically or kinematically [2]. However, in the present case, periodic boundary conditions are used (for which there is no actual boundary), so that the deviation is rather to be attributed to the finite size of the polycrystals. It can be seen that limiting the deviation (or “systematic error”) to less than 1% requires 35 grains for the experimental microstructure, but only three grains for the fully uniform microstructure. In addition, only 10 to 15 grains are required to reach the asymptotic value (no deviation) for this last microstructure. This is consistent with the fact that, in a polycrystal, a grain typically has

<footnote n="8">⁸Of course, e is itself closely related to the variability of the moduli computed over several polycrystals, following the $e = 2\sigma_z/\bar{z}$ relation.</footnote>
<figure><img src="image_8.png" /><figcaption>Figure 5: Evolution of the average and absolute error (±2σ) as a function of the number of grains (N), for the different microstructures. “Grain-growth – random” corresponds to the experimental microstructure, and “Dirac – uniform” to the fully uniform microstructure. Case of the shear modulus (μ) of copper.</figcaption></figure>

<table><caption>Table 4: Microstructure dependency of the shear modulus for the asymptotic case of N = 1000 grains. For each microstructure, the averages (statistical expectations) and standard errors were calculated from a set of 950 different polycrystals. Case of the shear modulus (μ) of copper.</caption><thead><tr><td>Microstructure</td><td>Average, μ&#772; [MPa]</td><td>Standard error [MPa (%)]</td></tr></thead><tbody><tr><td>Grain-growth – random</td><td>48 541</td><td>47 (0.10)</td></tr><tr><td>Grain-growth – uniform</td><td>48 500</td><td>39 (0.08)</td></tr><tr><td>Voronoi – random</td><td>48 566</td><td>32 (0.07)</td></tr><tr><td>Voronoi – uniform</td><td>48 533</td><td>15 (0.03)</td></tr><tr><td>Dirac – random</td><td>48 529</td><td>29 (0.06)</td></tr><tr><td>Dirac – uniform</td><td>48 495</td><td>5 (0.01)</td></tr></tbody></table>
12 to 14 neighbors on average [14]. However, it is interesting to see that only minimal deviation is obtained for polycrystals containing as few as three grains in the case of a fully uniform microstructure. It can be concluded that the deviations observed in other microstructures are due to their greater variability, which influences their collective response when the number of grains in the polycrystal is small. No clear trend can be identified in Figure 6 regarding which microstructural parameter contributes most to the deviation.

<figure><img src="image_9.png" /><figcaption>Figure 6: Close-up view ($N \le 100$) of the average values computed for the different microstructures. The dashed line represents the asymptotic value, while the gray region represents a ±1% deviation from this value. “Grain-growth – random” corresponds to the experimental microstructure, and “Dirac – uniform” to the fully uniform microstructure. Case of the shear modulus ($\mu$) of copper.</figcaption></figure>

## 4.2 Influence of the number of grains

Figure 7 provides the evolution of the error ($e$) as a function of the number of grains ($N$), for the different microstructures. For each microstructure, the error decreases smoothly and can be approximated by a continuous function of the form $e = a/N^b$. The values of $a$ and $b$ were determined by fitting over the range $N \in [100, 1000]$ (i.e., avoiding the deviation regions). This provided similar $b$ values between 0.48 and 0.51, with no apparent dependence on the microstructure. A value of 0.5 actually corresponds to the standard law of statistical theory for averaging over a set of uncorrelated (identically distributed) measurements (similarly to the standard error introduced previously) [4]. This suggests that the shear modulus calculation adheres to the same law as a standard property averaging, regardless of the complexity of the actual stress distribution within the body. So, in the following, and to facilitate further analysis, we will consider the same value of 0.5 for all microstructures, which corresponds to a function of the form

$$e = \frac{a}{\sqrt{N}}. \qquad (22)$$

The corresponding values of $a$ for the different microstructures are listed in Table 5, and the continuous fits are shown in Figure 7. Note how they match the numerical results for all
$N$ values (including for $N < 100$, although these values were excluded from the fit). The values of $a$ decrease from 0.942 for the experimental microstructure to 0.1 for the fully uniform microstructure. Thus, for a given number of grains ($N$), using the fully uniform microstructure reduces the error ($e$) by a factor of 9.4. Equivalently, for a given error ($e$), using the fully uniform microstructure reduces the number of grains ($N$) by a factor of $9.4^2 \simeq 90$. This indicates, for instance, that an experimental polycrystal containing 1000 grains can be replaced by a fully uniform polycrystal containing only 11 grains (while still avoiding significant deviation, as seen in Section 4.1).

<figure><img src="image_10.png" /><figcaption>Figure 7: Evolution of the error ($e$) as a function of the number of grains ($N$). The continuous fits are of the form $\varepsilon = a/\sqrt{N}$, and the values of $a$ are provided in Table 5. “Grain-growth – random” corresponds to the experimental microstructure, and “Dirac – uniform” to the fully uniform microstructure. Case of the shear modulus ($\mu$) of copper.</figcaption></figure>

<table><caption>Table 5: Values of the $a$ parameter for the different microstructures. “Grain-growth – random” corresponds to the experimental microstructure, and “Dirac – uniform” to the fully uniform microstructure. Case of the shear modulus ($\mu$) of copper.</caption><thead><tr><td>Microstructure</td><td>a [-]</td></tr></thead><tbody><tr><td>Grain-growth – random</td><td>0.942</td></tr><tr><td>Grain-growth – uniform</td><td>0.746</td></tr><tr><td>Voronoi – random</td><td>0.646</td></tr><tr><td>Voronoi – uniform</td><td>0.281</td></tr><tr><td>Dirac – random</td><td>0.591</td></tr><tr><td>Dirac – uniform</td><td>0.100</td></tr></tbody></table>

### 4.3 Influence of the grain size and grain orientation distributions

Both the variabilities associated with the grain size distribution ($\sigma_d$) and the grain orientation distribution ($\sigma_f$) are expected to contribute to $a$. This is first investigated in Figure 8,
which represents the evolution of a as a function of $\sigma_d$ and $\sigma_f$. These plots are actually slices of what should be a 3D plot of $a = \varphi_1(\sigma_d, \sigma_f)$, but they are preferred for clarity, given the relatively small number of available data points. It can be seen from Figure 8 that a exhibits similar levels of dependency on both $\sigma_d$ and $\sigma_f$. Each graph also exhibits large variations of a for a given $x$ value ($\sigma_d$ or $\sigma_f$), which comes from the dependence on the other variable ($\sigma_f$ or $\sigma_d$, respectively). However, the expression of $\varphi_1(\sigma_d, \sigma_f)$ is not necessarily trivial. Since we know that the polycrystal response (and therefore $\mu$) directly depends on the effective orientation distribution ($f^e$), it is actually natural to relate a to its variability ($\sigma_f^e$). $\sigma_f^e$ cumulates the effects of both $\sigma_d$ and $\sigma_f$, and can be easily computed from $f^e$. Its values were provided in Table 3. Figure 9 shows the evolution of a as a function of $\sigma_f^e$. The correlation is evident, and the data can be closely approximated by a function of the form

$$a = a_1 + a_2 \sigma_f^e, \qquad (23)$$

with $a_1 = 0.069$ and $a_2 = 1.721$ (for copper). Interestingly, Equation 23 includes not only a dependence on $\sigma_f^e$ (via $a_2 \sigma_f^e$), but also a constant, non-zero term, $a_1$. This term represents the “residual” error that would exist even if there was no microstructural variability ($\sigma_f = \sigma_d = \sigma_f^e = 0$). Fundamentally, this error can be attributed to the discrete representation of the orientation distribution in the polycrystals. It is related to the fact that, even with perfectly uniform orientation sampling and no variability in grain size, two polycrystals would still be different (in their actual grain orientations and grain arrangement) and therefore provide different values of $\mu$. Only a homogeneous material (or an isotropic behavior) would provide the same value. $a_1$ is small compared to $a_2 \sigma_f^e$ for the experimental microstructure, but dominates for the fully uniform microstructure.

<figure><img src="image_11.png" /><figcaption>Figure 8: Influence of the variabilities of (a) the grain size distribution ($\sigma_d$) and (b) the grain orientation distribution ($\sigma_f$) on a. Case of the shear modulus ($\mu$) of copper.</figcaption></figure>

Finally, it is interesting to dissociate (and compare) the contributions of the three different sources of variability on a: the discrete representation of the orientation distribution (through $a_1$), the variability in grain size ($\sigma_d$) and the variability of the grain orientation distribution ($\sigma_f$).
<figure><img src="image_12.png" /><figcaption>Figure 9: Influence of the variability of the effective grain orientation distribution ($\sigma_f^e$) on a. Case of the shear modulus ($\mu$) of copper.</figcaption></figure>

The latter two contributions can be determined from the $a_2 \sigma_f^e$ term of Equation 23 by applying the standard composition law for standard deviations of uncorrelated variables:

$$ \sigma_f^{e2} = \sigma_f^2 + \varphi_2(\sigma_d)^2. \qquad (24) $$

Function $\varphi_2(\sigma_d)$ is unknown, but its values can simply be computed from those of $\sigma_f$ and $\sigma_f^e$ using Equation 24. Table 6 provides the contributions of $a_1$, $\sigma_f$ and $\sigma_d$ (through $\varphi_2(\sigma_d)$) to a, for all microstructures. For the experimental microstructure, it can be seen that the error ($e$) is mainly due to grain size variability (59%), while random orientation sampling accounts for 34%, and the discrete representation of the orientation distribution for 7%. Reducing the variability in grain size drastically changes the distribution, as the error ($e$) of the Voronoi – random microstructure is rather caused by the random orientation sampling (73%). In contrast, for the fully uniform microstructure, the error ($e$) is dominated by the discrete representation of the orientation distribution (69%).

<table><caption>Table 6: Relative influences of the microstructural variability on the error, $e$ (via $a$). $a$ is decomposed into its two contributions, $a_1$ and $a_2 \sigma_f^e$. The contributions of $\sigma_f$ and $\varphi_2(\sigma_d)$ to $a_2 \sigma_f^e$ correspond to the squared contributions, $a_2 \sigma_f^2/\sigma_f^{e2}$ and $a_2 \varphi_2(\sigma_d)^2/\sigma_f^{e2}$. Case of the shear modulus ($\mu$) of copper.</caption><thead><tr><th>Microstructure</th><th>a [-]</th><th>a<sub>1</sub> [- (%)]</th><th>a<sub>2</sub> σ<sub>f</sub><sup>e</sup> [- (% = σ<sub>f</sub>% + ϕ<sub>2</sub>(σ<sub>d</sub>)%)]</th></tr></thead><tbody><tr><td>Grain-growth – random</td><td>0.942</td><td>0.069 ( 7%)</td><td>0.862 (93% = 34% + 59%)</td></tr><tr><td>Grain-growth – uniform</td><td>0.746</td><td>0.069 ( 9%)</td><td>0.685 (91% = 0% + 91%)</td></tr><tr><td>Voronoi – random</td><td>0.646</td><td>0.069 (11%)</td><td>0.575 (89% = 73% + 16%)</td></tr><tr><td>Voronoi – uniform</td><td>0.281</td><td>0.069 (24%)</td><td>0.236 (76% = 0% + 75%)</td></tr><tr><td>Dirac – random</td><td>0.591</td><td>0.069 (12%)</td><td>0.525 (88% = 87% + 0%)</td></tr><tr><td>Dirac – uniform</td><td>0.100</td><td>0.069 (68%)</td><td>0.013 (32% = 32% + 0%)</td></tr></tbody></table>
## 4.4 Influence of elastic anisotropy

The $a_1$ and $a_2$ factors appearing in Equation 23 depend on the material. Indeed, their values were determined in Section 4.3 for copper, but we also know that, for an isotropic material ($A = 1$ or, equivalently, $A_u = 0$), we must have $e = 0$, and so $a = 0$ and $a_1 = a_2 = 0$ (all quantities being non-negative by definition). The values of $a_1$ and $a_2$ can be computed for materials other than copper by running simulations for a specific number of grains (e.g. 300). Indeed, the value of $a$ can be computed from the measured error ($e$) using Equation 22, and $a_1$ and $a_2$ can then be computed from $a$ as was done in Section 4.3. All the materials listed in Table 1 were considered, as well as “virtual materials” of Zener ratios ($A$) varying continuously in the range [0.5, 10], for which (relative) $C_{ij}$ values were obtained using Equations 17 to 20. The evolutions of $a_1$ and $a_2$ as a function of the Zener ratio ($A$) and the universal anisotropy factor ($A_u$) are provided in Figures 10 and 11, respectively. Both $a_1$ and $a_2$ show clear dependencies on $A$ (and $A_u$), with the expected minimal values of $a_1 = a_2 = 0$ for $A = 1$ (and $A_u = 0$). It should be noted that the data obtained for the different materials (copper, aluminium, etc.) deviate only slightly from the line corresponding to the virtual materials, despite their different values of the $r = C_{44}/C_{11}$ ratio. This clearly indicates that $a$ mainly depends on the anisotropy ratio ($A$ or $A_u$), and that the $r$ ratio can be ignored. Finally, materials with very large Zener ratio ($10^2-10^9$) were considered, for which $a_1$ and $a_2$ converged toward saturation values of $a_1^{\text{sat}} \simeq 1.34$ and $a_2^{\text{sat}} \simeq 4.67$, respectively. Such a saturation behavior, and the fact that all curves exhibit decreasing slopes as anisotropy increases, indicate that further increases in the anisotropy yield diminishing effects on the evolution of the error ($e$). Ultimately, even an infinitely anisotropic material produces finite values of both the error ($e$) and the RVE size ($N$). Even if the asymptotic behavior develops well beyond the practical anisotropy range, it will be considered when fitting the data (e.g. for application to materials with higher elastic anisotropy). Since the numerical data do not follow well-defined analytical laws, fits over both the practical range of $A \in [0.5, 10]$ and the entire range of $A \ge 0$ will be provided, yielding to different expressions. The same will be done for $A_u$.

### 4.4.1 Expressions of $a_1$ and $a_2$ as a function of $A$

Over the practical range of $A \in [0.5, 10]$, the $a_1$ values can be fitted by a function of the form

$$a_1(A) = \alpha_1 |\log A|^{\beta_1}, \qquad (25)$$

with $\alpha_1 = 0.044$ and $\beta_1 = 2$, as shown in Figure 10a. An expression of $a_1$ including saturation (valid for all $A$) is

$$a_1(A) = a_1^{\text{sat}} \frac{\alpha_1 |\log A|^{\beta_1}}{\alpha_1 |\log A|^{\beta_1} + a_1^{\text{sat}}}, \qquad (26)$$

with revised values of $\alpha_1 = 0.041$ and $\beta_1 = 2.3$. Over the practical range of $A \in [0.5, 10]$, the values of $a_2$ can be fitted by a function of the form

$$a_2(A) = \alpha_2 \left| 1 - \frac{1}{A^{\beta_2}} \right|, \qquad (27)$$

with parameters $\alpha_2 = 8.4$ and $\beta_2 = 0.2$, as shown in Figure 10b. Properly accounting for saturation (which is built into the expression) yields revised values of $\alpha_2 = a_2^{\text{sat}} = 4.67$ and $\beta_2 = 0.46$ (valid for all $A$). It is unclear why the data align with these types of functions.
However, the expression of $a_2$ (Equation 27) is a scaled Cox-Box transformation, which is related to the expression of $a_1$ (Equation 25) and becomes equivalent to a logarithmic transformation as the exponent ($\beta_2$) approaches zero.

### 4.4.2 Expressions of $a_1$ and $a_2$ as a function of $A_u$

A similar fitting procedure can be applied to the evolution of $a_1$ and $a_2$ as functions of $A_u$. Since $A_u$ is related to A by Equation 18, one could use the expressions determined in Section 4.4.1 to calculate $a_1$ and $a_2$ as a function of $A_u$. However, it may be seen in Figure 11 that $a_1$ and $a_2$ show relatively simple dependencies on $A_u$, so that new expressions of $a_1$ and $a_2$ as a function of $A_u$ can be determined. Over the practical range of $A_u \in [0, 10]$, the $a_1$ values can be fitted by a function of the form

$$a_1(A_u) = \alpha'_1 A_u^{\beta'_1}, \qquad (28)$$

with $\alpha'_1 = 0.037$ and $\beta'_1 = 0.82$, as shown in Figure 11a. As done previously for $a_1(A)$ (Equations 25 and 26), saturation can then be incorporated, which leads to

$$a_1(A_u) = a_1^{\text{sat}} \frac{\alpha'_1 A_u^{\beta'_1}}{\alpha'_1 A_u^{\beta'_1} + a_1^{\text{sat}}}, \qquad (29)$$

with revised values of $\alpha'_1 = 0.054$ and $\beta'_1 = 0.69$. The values of $a_2$ can be fitted (over the full range of $A_u$) by a function of the form

$$a_2(A_u) = a_2^{\text{sat}} \frac{\alpha'_2 A_u^{\beta'_2}}{\alpha'_2 A_u^{\beta'_2} + a_2^{\text{sat}}}, \qquad (30)$$

with $\alpha'_2 = 2.14$ and $\beta'_2 = 0.62$, as shown in Figure 11b.

### 4.4.3 Expression of the error, e

Using Equations 22, 23, 25 and 27, we finally get the general expression of the error as a function of the number of grains ($N$), the variability of the effective orientation distribution ($\sigma_f^e$) and the Zener ratio ($A$). For $A \in [0.5, 10]$,

$$e(A, \sigma_f^e, N) = \frac{\alpha_1 |\log A|^{\beta_1} + \alpha_2 \left|1 - \frac{1}{A^{\beta_2}}\right| \sigma_f^e}{\sqrt{N}}, \qquad (31)$$

with $\alpha_1 = 0.044$, $\beta_1 = 2$, $\alpha_2 = 8.4$ and $\beta_2 = 0.2$. The expression including saturation is

$$e(A, \sigma_f^e, N) = \frac{a_1^{\text{sat}} \frac{|\log A|^{\beta_1}}{|\log A|^{\beta_1} + a_1^{\text{sat}}} + a_2^{\text{sat}} \left|1 - \frac{1}{A^{\beta_2}}\right| \sigma_f^e}{\sqrt{N}}, \qquad (32)$$

with $\alpha_1 = 0.041$, $\beta_1 = 2.3$, $\alpha_2 = 4.67$ and $\beta_2 = 0.46$.

Finally, a simplified expression of the error ($e$) specific to the experimental microstructure can be determined. As we know that $\sigma_f^e = 0.5$ and that the error ($e$) is dominated by the $a_2 \sigma_f^e$ term (as seen in Figure 9 and Table 6), a scaled Cox-Box transformation function can be used (as in Equation 27). This finally leads to the general expression of the error for the
experimental microstructure, $e_{\text{exp}}$, as a function of the number of grains ($N$) and the Zener ratio ($A$):

$$e_{\text{exp}}(A, N) = \frac{a_{\text{exp}}^{\text{sat}}}{\sqrt{N}} \left| 1 - \frac{1}{A^{\beta_{\text{exp}}}} \right|, \qquad (33)$$

with $a_{\text{exp}}^{\text{sat}} = a_{1}^{\text{sat}} + 0.5 a_{2}^{\text{sat}} = 3.67$ and $\beta_{\text{exp}} = 0.152$ (valid for all $A$). Following the same procedure, we can get the general expression of the error for the experimental microstructure ($e_{\text{exp}}$) as a function of the number of grains ($N$) and the universal anisotropy factor ($A_u$):

$$e_{\text{exp}}(A_u, N) = \frac{a_{\text{exp}}^{\text{sat}}}{\sqrt{N}} \frac{\alpha'_{\text{exp}} A_u^{\beta'_{\text{exp}}}}{\alpha'_{\text{exp}} A_u^{\beta'_{\text{exp}}} + a_{\text{exp}}^{\text{sat}}}, \qquad (34)$$

with $\alpha'_{\text{exp}} = 0.89$ and $\beta'_{\text{exp}} = 0.52$. The evolutions of $e_{\text{exp}}\sqrt{N}$ as a function of $A$ and $A_u$ are provided in Figures 12 and 13, respectively. For a specific material ($A$ or $A_u$) and a number of grains ($N$), it is possible to determine the error on $\mu$ ($e$). For example, considering copper ($A = 3.2$) and $N = 1000$, we get $e \simeq 3.6\%$, which is consistent with the results provided in Figure 7.

<figure><img src="image_13.png" /><figcaption>Figure 10: Influence of the Zener ratio ($A$) on (a) $a_1$ and (b) $a_2$. Case of the shear modulus ($\mu$) and cubic symmetry.</figcaption></figure>
a

<figure><img src="image_14.png" /></figure>

b

<figure><img src="image_15.png" /><figcaption>Figure 11: Influence of the universal anisotropy factor ($A_u$) on (a) $a_1$ and (b) $a_2$. Case of the shear modulus ($\mu$) and cubic symmetry.</figcaption></figure>

<figure><img src="image_16.png" /><figcaption>Figure 12: Influence of the Zener ratio ($A$) on $e_{\text{exp}}/\sqrt{N}$. Case of the shear modulus ($\mu$) and cubic symmetry.</figcaption></figure>
<figure><img src="image_17.png" /><figcaption>Figure 13: Influence of the universal anisotropy factor ($A_u$) on $e_{\text{exp}}\sqrt{N}$. Case of the shear modulus ($\mu$) and cubic symmetry.</figcaption></figure>
## 4.5 Case of hexagonal symmetry

The analysis can be repeated for the materials with hexagonal symmetry listed in Table 2, in order to compute their $a_1$ and $a_2$ values. The dependency of $a_1$ and $a_2$ on the microstructural properties can then be determined. The cases of the shear and bulk moduli are analysed successively.

### 4.5.1 Case of the shear modulus, $\mu$

The evolution of $a_1$ and $a_2$ as a function of $A_u$ is illustrated in Figure 14, for the different materials. As there is noticeable scatter of the $a_1$ and $a_2$ values for similar $A_u$ values, the corresponding numerical values are also provided in Table 7. A larger scatter compared to the case of cubic symmetry can be explained by the fact that hexagonal symmetry has five independent $C_{ij}$ values (four after normalization), and so $A_u$ does not represent the elastic anisotropy as well as in the case of cubic symmetry, which has only three independent $C_{ij}$ values (two after normalization). The numerical data can be fitted directly using the same analytical laws as in the case of cubic symmetry. Over the practical range of $A_u \in [0, 2]$, the $a_1$ values can be fitted by

$$a_1(A_u) = \alpha'_1 A_u^{\beta'_1}, \qquad (35)$$

with $\alpha'_1 = 0.042$ and $\beta'_1 = 0.53$. The saturation behavior was examined as in Section 4.4, using virtual materials with very high elastic anisotropies. Starting from the elastic parameters of cadmium ($C_{ij}$; see Table 2), the $C_{44}$ value was modified to obtain very large $A_u$ values ($10^2-10^9$), while all other $C_{ij}$ values remained unchanged. This lead to a value of $a_1^{\text{sat}} = 0.567$. As in the case of cubic symmetry, saturation develops well outside of the practical range of $A_u \in [0, 2]$, but can be represented by incorporating saturation (in the same way as in Equation 29). Over the practical range of $A_u \in [0, 2]$, the $a_2$ values could be fitted by

$$a_2(A_u) = a_2^{\text{sat}} \frac{\alpha'_2 A_u^{\beta'_2}}{\alpha'_2 A_u^{\beta'_2} + a_2^{\text{sat}}}, \qquad (36)$$

with $a_2^{\text{sat}} = 1.56$, $\alpha'_2 = 1$ and $\beta'_2 = 0.78$.

Interestingly, the dependence of $a_1$ and $a_2$ on $A_u$ follows the same form for both hexagonal and cubic symmetry, but the values of $\alpha'_1$, $\alpha'_2$, $\beta'_1$ and $\beta'_2$ differ. In the practical range of $A_u \in [0, 2]$, $a_1$ has values similar to those for cubic symmetry, but $a_2$ is smaller by approximately a factor of two. So, both the error (e) and the size of the representative element (N) are generally smaller for hexagonal symmetry than for cubic symmetry.

For titanium, the value of $a = a_1 + a_2 \sigma_f^e$ (and therefore the value of the error (e) for a given number of grains (N)) decreases from 0.183 for the experimental microstructure to 0.015 for the fully uniform microstructure, which corresponds to a reduction factor of about 12. For zinc, the value of a decreases from 0.394 for the experimental microstructure to 0.059 for the fully uniform microstructure, which corresponds to a reduction factor of 6.7. Equivalently, for a given error (e), using the fully uniform microstructure reduces the number of grains (N) by factors of about 140 and 45 for titanium and zinc, respectively.

Finally, the expression of the error (e) specific to a hexagonal-symmetry experimental microstructure, is

$$e_{\text{exp}}(A_u, N) = \frac{a_{\text{exp}}^{\text{sat}} \alpha'_{\text{exp}} A_u^{\beta'_{\text{exp}}}}{\sqrt{N} \alpha'_{\text{exp}} A_u^{\beta'_{\text{exp}}} + a_{\text{exp}}^{\text{sat}}}, \qquad (37)$$
with $a_{\text{exp}}^{\text{sat}} = 1.34$, $\alpha'_{\text{exp}} = 0.47$ and $\beta'_{\text{exp}} = 0.68$. The evolution of $e_{\text{exp}}/\sqrt{N}$ as a function of $A_u$ is provided in Figure 15. For example, for $N = 1000$, we get $e \approx 0.58\%$ and $1.25\%$ for titanium and zinc, respectively. These values are significantly smaller than those of copper, which is due to smaller values of $a_2$ (for a given value of $A_u$) and $A_u$ itself. For the fully uniform microstructure, these values are reduced to $e \approx 0.004\%$ and $0.009\%$ for titanium and zinc, respectively. Equivalently, for the experimental microstructure, an error ($e$) of $1\%$ can be obtained with $N = 335$ and $1552$ for titanium and zinc, respectively. Considering a fully uniform microstructure reduces these numbers to $N = 2-3$ and $34$ for titanium and zinc, respectively. For titanium, values of $N \ge 10$ should still be used to avoid systematic deviation.

<figure><img src="image_18.png" /><figcaption>Figure 14: Influence of the universal anisotropy factor ($A_u$) on $a_1$ and $a_2$. Case of the shear modulus ($\mu$) of hexagonal-symmetry materials. Labels are provided only for Mg, Zr, Ti, Cd, Zn and Tl for clarity.</figcaption></figure>

### 4.5.2 Case of the bulk modulus, $\kappa$

In the case of the bulk modulus ($\kappa$), it is useful to examine the Voigt, Reuss, and Hill averages ($\kappa_V$, $\kappa_R$, and $\kappa_H$, respectively) before analysing the $a_1$ and $a_2$ values. The computation of these averages from the elastic constants $C_{ij}$ is described in Appendix A. The corresponding values for the materials listed in Table 2 are provided in Table 8. The error on the Hill average is computed from the difference between the Voigt and Reuss averages. For most materials, the error on the Hill average is very small, indicating that computing $\kappa$ using polycrystalline simulations is not necessarily required. Among all materials, only cadmium and zinc exhibit large errors on their Hill averages (about 10%) and actually benefit from a numerical computation of $\kappa$ and an RVE size analysis. These materials are also characterized by particularly high $A_u$ values (thallium, however, shows relatively close Voigt and Reuss averages despite its high $A_u$ value). The values of $a_1$ and $a_2$ computed for all materials are provided in Table 8. Materials with close Voigt and Reuss averages exhibit very small $a_1$ and $a_2$ values, on the order of $10^{-6}$ to $10^{-3}$. According to the relationship $e = a/\sqrt{N}$ (see Equation 22), values as small as $10^{-4}$ indicate that $\kappa$ can be estimated with a relative error of $10^{-4}$ (0.01%) using only a single grain. This is consistent with
<table><caption>Table 7: Values of $a_1$ and $a_2$ in the case of the shear modulus ($\mu$) of hexagonal-symmetry materials.</caption><thead><tr><th>Material</th><th>A<sub>u</sub> [-]</th><th>a<sub>1</sub> [-]</th><th>a<sub>2</sub> [-]</th></tr></thead><tbody><tr><td>Lu</td><td>0.01</td><td>1.7 × 10<sup>-3</sup></td><td>1.1 × 10<sup>-2</sup></td></tr><tr><td>Os</td><td>0.01</td><td>4.8 × 10<sup>-3</sup></td><td>4.7 × 10<sup>-2</sup></td></tr><tr><td>Er</td><td>0.02</td><td>6.8 × 10<sup>-3</sup></td><td>6.1 × 10<sup>-2</sup></td></tr><tr><td>Ho</td><td>0.02</td><td>5.8 × 10<sup>-3</sup></td><td>6.3 × 10<sup>-2</sup></td></tr><tr><td>Tb</td><td>0.02</td><td>6.3 × 10<sup>-3</sup></td><td>4.6 × 10<sup>-2</sup></td></tr><tr><td>Gd</td><td>0.03</td><td>8.5 × 10<sup>-3</sup></td><td>6.8 × 10<sup>-2</sup></td></tr><tr><td>Hf</td><td>0.03</td><td>6.8 × 10<sup>-3</sup></td><td>1.1 × 10<sup>-1</sup></td></tr><tr><td>Ru</td><td>0.03</td><td>8.4 × 10<sup>-3</sup></td><td>5.7 × 10<sup>-2</sup></td></tr><tr><td>Mg</td><td>0.04</td><td>9.3 × 10<sup>-3</sup></td><td>6.9 × 10<sup>-2</sup></td></tr><tr><td>Be</td><td>0.05</td><td>3.1 × 10<sup>-3</sup></td><td>1.9 × 10<sup>-1</sup></td></tr><tr><td>Tm</td><td>0.05</td><td>7.7 × 10<sup>-3</sup></td><td>4.9 × 10<sup>-2</sup></td></tr><tr><td>Sc</td><td>0.08</td><td>1.5 × 10<sup>-2</sup></td><td>9.5 × 10<sup>-2</sup></td></tr><tr><td>Zr</td><td>0.09</td><td>1.5 × 10<sup>-2</sup></td><td>1.0 × 10<sup>-1</sup></td></tr><tr><td>Re</td><td>0.11</td><td>1.7 × 10<sup>-2</sup></td><td>1.1 × 10<sup>-1</sup></td></tr><tr><td>Dy</td><td>0.12</td><td>7.2 × 10<sup>-3</sup></td><td>2.7 × 10<sup>-1</sup></td></tr><tr><td>Nd</td><td>0.13</td><td>1.8 × 10<sup>-2</sup></td><td>1.5 × 10<sup>-1</sup></td></tr><tr><td>La</td><td>0.17</td><td>1.9 × 10<sup>-2</sup></td><td>1.3 × 10<sup>-1</sup></td></tr><tr><td>Pr</td><td>0.19</td><td>2.1 × 10<sup>-2</sup></td><td>2.0 × 10<sup>-1</sup></td></tr><tr><td>Ti</td><td>0.19</td><td>1.3 × 10<sup>-2</sup></td><td>3.4 × 10<sup>-1</sup></td></tr><tr><td>Co</td><td>0.20</td><td>2.1 × 10<sup>-2</sup></td><td>2.1 × 10<sup>-1</sup></td></tr><tr><td>Cd</td><td>1.16</td><td>3.7 × 10<sup>-2</sup></td><td>6.3 × 10<sup>-1</sup></td></tr><tr><td>Zn</td><td>1.82</td><td>5.4 × 10<sup>-2</sup></td><td>6.8 × 10<sup>-1</sup></td></tr><tr><td>Tl</td><td>1.91</td><td>7.3 × 10<sup>-2</sup></td><td>1.0 × 10<sup>+0</sup></td></tr></tbody></table>
<figure><img src="image_19.png" /><figcaption>Figure 15: Influence of the universal anisotropy factor ($A_u$) on $e\sqrt{N}$, for the experimental microstructure. Case of the shear modulus ($\mu$) and hexagonal symmetry. Labels are provided only for Mg, Zr, Ti, Cd, Zn and Tl for clarity.</figcaption></figure>

the close values of the Voigt and Reuss averages, and with the Hill average itself being defined to within about 0.01% error. In contrast, cadmium and zinc show large $a_1$ and $a_2$ values, on the order of $10^{-2}$ to $10^{-1}$. For zinc in particular, the value of $a = a_1 + a_2 \sigma_f^e$ (and therefore the error (e) for a given number of grains ($N$)) decreases from 0.406 for the experimental microstructure to 0.036 for the fully uniform microstructure, corresponding to a reduction factor of about 11. Equivalently, for a given error (e), the number of grains $N$ can be reduced by a factor of approximately 120. For example, for $N = 1000$, one obtains $e \simeq 1.2\%$ for the experimental microstructure and $e \simeq 0.2\%$ for the fully uniform microstructure. Equivalently, an error (e) of 1% can be obtained with $N = 1552$ grains for the experimental microstructure, but with only $N = 35$ grains for the fully uniform microstructure.
<table><caption>Table 8: Values of $a_1$ and $a_2$ in the case of the bulk modulus ($\kappa$) of hexagonal-symmetry materials. $\kappa_V$, $\kappa_R$ and $\kappa_H$ are the Voigt, Reuss and Hill averages, respectively.</caption><thead><tr><th>Material</th><th>A<sub>u</sub> [-]</th><th>κ<sub>R</sub> [GPa]</th><th>κ<sub>V</sub> [GPa]</th><th>κ<sub>H</sub> [GPa]</th><th>a<sub>1</sub> [-]</th><th>a<sub>2</sub> [-]</th></tr></thead><tbody><tr><td>Lu</td><td>0.01</td><td>47.6</td><td>47.7</td><td>47.6 ± 0.12%</td><td>5.3 × 10<sup>-4</sup></td><td>1.5 × 10<sup>-3</sup></td></tr><tr><td>Os</td><td>0.01</td><td>411.6</td><td>412.2</td><td>411.9 ± 0.08%</td><td>3.4 × 10<sup>-4</sup></td><td>9.5 × 10<sup>-4</sup></td></tr><tr><td>Er</td><td>0.02</td><td>44.6</td><td>44.7</td><td>44.7 ± 0.05%</td><td>2.3 × 10<sup>-4</sup></td><td>6.3 × 10<sup>-4</sup></td></tr><tr><td>Ho</td><td>0.02</td><td>40.9</td><td>40.9</td><td>40.9 ± 0.00%</td><td>1.5 × 10<sup>-5</sup></td><td>4.4 × 10<sup>-5</sup></td></tr><tr><td>Tb</td><td>0.02</td><td>38.8</td><td>38.8</td><td>38.8 ± 0.01%</td><td>5.4 × 10<sup>-5</sup></td><td>1.5 × 10<sup>-4</sup></td></tr><tr><td>Gd</td><td>0.03</td><td>37.9</td><td>37.9</td><td>37.9 ± 0.00%</td><td>0.0 × 10<sup>+0</sup></td><td>0.0 × 10<sup>+0</sup></td></tr><tr><td>Hf</td><td>0.03</td><td>108.5</td><td>108.6</td><td>108.5 ± 0.01%</td><td>3.0 × 10<sup>-5</sup></td><td>8.4 × 10<sup>-5</sup></td></tr><tr><td>Ru</td><td>0.03</td><td>310.6</td><td>310.9</td><td>310.7 ± 0.05%</td><td>2.1 × 10<sup>-4</sup></td><td>5.7 × 10<sup>-4</sup></td></tr><tr><td>Mg</td><td>0.04</td><td>35.2</td><td>35.2</td><td>35.2 ± 0.01%</td><td>5.3 × 10<sup>-5</sup></td><td>1.4 × 10<sup>-4</sup></td></tr><tr><td>Be</td><td>0.05</td><td>111.3</td><td>111.7</td><td>111.5 ± 0.15%</td><td>6.9 × 10<sup>-4</sup></td><td>2.1 × 10<sup>-3</sup></td></tr><tr><td>Tm</td><td>0.05</td><td>45.8</td><td>46.4</td><td>46.1 ± 0.65%</td><td>3.0 × 10<sup>-3</sup></td><td>8.1 × 10<sup>-3</sup></td></tr><tr><td>Sc</td><td>0.08</td><td>55.8</td><td>55.8</td><td>55.8 ± 0.01%</td><td>2.7 × 10<sup>-5</sup></td><td>7.1 × 10<sup>-5</sup></td></tr><tr><td>Zr</td><td>0.09</td><td>96.5</td><td>96.7</td><td>96.6 ± 0.09%</td><td>4.4 × 10<sup>-4</sup></td><td>1.1 × 10<sup>-3</sup></td></tr><tr><td>Re</td><td>0.11</td><td>365.0</td><td>365.0</td><td>365.0 ± 0.00%</td><td>0.0 × 10<sup>+0</sup></td><td>0.0 × 10<sup>+0</sup></td></tr><tr><td>Dy</td><td>0.12</td><td>40.5</td><td>40.5</td><td>40.5 ± 0.00%</td><td>5.6 × 10<sup>-6</sup></td><td>1.7 × 10<sup>-5</sup></td></tr><tr><td>Nd</td><td>0.13</td><td>31.8</td><td>31.8</td><td>31.8 ± 0.01%</td><td>4.5 × 10<sup>-5</sup></td><td>1.2 × 10<sup>-4</sup></td></tr><tr><td>La</td><td>0.17</td><td>25.9</td><td>26.0</td><td>26.0 ± 0.04%</td><td>2.1 × 10<sup>-4</sup></td><td>5.4 × 10<sup>-4</sup></td></tr><tr><td>Pr</td><td>0.19</td><td>28.8</td><td>28.8</td><td>28.8 ± 0.00%</td><td>7.2 × 10<sup>-6</sup></td><td>1.9 × 10<sup>-5</sup></td></tr><tr><td>Ti</td><td>0.19</td><td>105.0</td><td>105.0</td><td>105.0 ± 0.00%</td><td>1.3 × 10<sup>-5</sup></td><td>3.6 × 10<sup>-5</sup></td></tr><tr><td>Co</td><td>0.20</td><td>187.4</td><td>187.4</td><td>187.4 ± 0.01%</td><td>2.7 × 10<sup>-5</sup></td><td>7.1 × 10<sup>-5</sup></td></tr><tr><td>Cd</td><td>1.16</td><td>48.9</td><td>59.2</td><td>54.0 ± 9.56%</td><td>3.6 × 10<sup>-2</sup></td><td>8.4 × 10<sup>-2</sup></td></tr><tr><td>Zn</td><td>1.82</td><td>59.5</td><td>72.7</td><td>66.1 ± 9.98%</td><td>3.1 × 10<sup>-2</sup></td><td>7.5 × 10<sup>-2</sup></td></tr><tr><td>Tl</td><td>1.91</td><td>36.1</td><td>36.3</td><td>36.2 ± 0.16%</td><td>8.3 × 10<sup>-4</sup></td><td>2.3 × 10<sup>-3</sup></td></tr></tbody></table>
## 5 Concluding remarks

The dependence of the RVE size on the properties of polycrystalline materials, as well as the potential to reduce it through microstructural modification, was analyzed using dedicated polycrystal generation methods. The RVE size corresponds to the number of grains $N$ required in a polycrystal to compute an average property with a specified level of precision (or error, $e$). The shear modulus ($\mu$) was of particular interest, as it exhibits greater variability compared to the bulk modulus ($\kappa$). A relationship of the following form was observed:

$$e = \frac{a}{\sqrt{N}}, \qquad (38)$$

where $a$ depends on the material and microstructural properties. Computing an elastic modulus therefore follows the same law as averaging of uncorrelated (identically distributed) measurements. Knowing the value of $a$ allows one not only to determine the number of grains ($N$) required to estimate an elastic modulus within a given error ($e$), but also to predict the error ($e$) expected in a polycrystalline volume composed of a known number of grains ($N$). Equation 38 therefore represents a structure-property relationship that can be used before (or even without) computing the elastic modulus (or the corresponding mechanical field over the microstructure). When the number of grains is too small ($N < 3$ to $N < 35$, depending on the microstructure), a systematic bias develops in addition to the statistical error $e$. This deviation is related to the degree of microstructural variability, which influences the collective response of polycrystals when the number of grains is too small.

The expression of $a$ involves both the elastic anisotropy of the crystal ($A$ or $A_u$) and the variability of the effective orientation distribution of the microstructure, $\sigma_f^e$:

$$a = a_1(A) + a_2(A) \sigma_f^e \quad \text{and} \quad a = a_1(A_u) + a_2(A_u) \sigma_f^e, \qquad (39)$$

where $A$ is the Zener anisotropy ratio and $A_u$ is the universal anisotropy factor. Interestingly, factor $a$ does not directly involve information on the microstructure morphology, such as the grain size distribution. This information is instead directly contained in the effective orientation distribution ($f^e$). This tends to demonstrate an equivalency between the over-representation of an orientation due to statistical sampling and its increased weighting resulting from a larger grain size. This equivalency holds by definition for $f^e$, but was not necessarily expected for the error ($e$).

The expressions of $a_1$ and $a_2$ as a function of elastic anisotropy ($A$ or $A_u$) were determined by fitting of the numerical values. The same $a_1(A_u)$ and $a_2(A_u)$ laws apply to both cubic and hexagonal symmetries, albeit with different values of the parameters. For hexagonal symmetry, the dependence on $A_u$ is less pronounced than for cubic symmetry, and so larger deviations from the analytical law should be expected. For the same degree of anisotropy ($A_u$), hexagonal symmetry generally leads to smaller values of $a$ compared to cubic symmetry. Moreover, hexagonal-symmetry materials generally exhibit lower elastic anisotropy ($A_u$) compared to cubic-symmetry materials, which further reduces the error ($e$) and, ultimately, the RVE size ($N$). In the case of the bulk modulus ($\kappa$) of hexagonal-symmetry materials, even smaller values of $a_1$ and $a_2$ were obtained, which can be explained by the close values of the Voigt and Reuss averages.

It was finally demonstrated that the error ($e$) can be significantly reduced by modifying the microstructure, in particular by artificially reducing the variability of its effective orientation
distribution ($\sigma_f^e$). Importantly, this approach does not introduce any significant bias on the predicted value of the elastic modulus. A “fully uniform” microstructure combines a “Dirac” grain size distribution (i.e., all grains having the same size) with nearly uniform orientation sampling and reduces $\sigma_f^e$ to nearly zero. This enables, for instance, a decrease of the error (e) by factors of 9.4 for copper, 12 for titanium and 6.7 for zinc. Equivalently, the number of grains in a polycrystal (N) can be reduced by factors of 90 for copper, 140 for titanium and 45 for zinc. A typical 1000-grain polycrystal can therefore be replaced by polycrystals containing 11, 8 and 22 grains for these three materials, respectively. Consequently, although realistic microstructures are typically employed in numerical simulations, this work demonstrates that, for certain applications such as the prediction of elastic moduli in bulk, single-phase polycrystals, reduced microstructural variability can provide equally accurate results with a significantly smaller error and an even greater reduction in RVE size.

# 6 Acknowledgments

This work was funded in part by the ANR AFOTI project (ANR-20-CE08-0020). This concerned Darcey Britton’s internship, which was co-supervised by François Villette during his postdoc. Both are acknowledged for valuable discussions.

## A Definition of the Voigt, Reuss and Hill averages, and computation of $A_u$

For cubic symmetry, the Voigt averages of the bulk and shear modulus, $\kappa_V$ and $\mu_V$, respectively, are

$$ \kappa_V = \frac{C_{11} + 2 C_{12}}{3} \qquad (40) $$

and

$$ \mu_V = \frac{C_{11} - C_{12} + 3 C_{44}}{5}. \qquad (41) $$

The Reuss averages on the bulk and shear modulus, $\kappa_R$ and $\mu_R$, respectively, are

$$ \kappa_R = \kappa_V \qquad (42) $$

and

$$ \mu_R = \frac{5 (C_{11} - C_{12}) C_{44}}{4 C_{44} + 3 (C_{11} - C_{12})}. \qquad (43) $$

For hexagonal symmetry, the Voigt averages of the bulk and shear modulus, $\kappa_V$ and $\mu_V$, respectively, are

$$ \kappa_V = \frac{2(C_{11} + C_{12}) + 4C_{13} + C_{33}}{9} \qquad (44) $$

and

$$ \mu_V = \frac{1}{5} (G_{\text{eff}}^V + 2 C_{44} + 2 C_{66}) \qquad (45) $$

with

$$ G_{\text{eff}}^V = \frac{C_{11} + C_{33} - 2 C_{13} - C_{66}}{3}. \qquad (46) $$

The Reuss averages of the bulk and shear modulus, $\kappa_V$ and $\mu_V$, respectively, are

$$ \frac{1}{\kappa_R - C_{13}} = \frac{1}{C_{11} - C_{66} - C_{13}} + \frac{1}{C_{33} - C_{13}} \qquad (47) $$
and

$$ \mu_R = \left[ \frac{1}{5} \left( \frac{1}{G_{\text{eff}}^R} + \frac{2}{C_{44}} + \frac{2}{C_{66}} \right) \right]^{-1} \qquad (48) $$

with

$$ G_{\text{eff}}^R = \frac{\kappa_R G_{\text{eff}}^V}{\kappa_V}. \qquad (49) $$

For both cubic and hexagonal symmetries, the Hill averages, $\kappa_H$ and $\mu_H$, are defined as the arithmetic averages of the Voigt and Reuss averages:

$$ \kappa_H = \frac{\kappa_V + \kappa_R}{2} \qquad (50) $$

and

$$ \mu_H = \frac{\mu_V + \mu_R}{2}. \qquad (51) $$

The universal anisotropy factor is then defined as [28]

$$ A_u = 5 \frac{\mu_V}{\mu_R} + \frac{\kappa_V}{\kappa_R} - 6. \qquad (52) $$

By definition, $A_u \ge 0$.
## References

[1] R. Quey and L. Renversade. Optimal polyhedral description of 3D polycrystals: method and application to statistical and synchotron X-ray diffraction data. *Comp. Methods Appl. Mech. Eng.*, 330:308-33, 2018.

[2] T. Kanit, S. Forest, I. Galliet, V. Mounoury, and D. Jeulin. Determination of the size of the representative volume element for random composites: statistical and numerical approach. *Int. J. Solids Struct.*, 40:3647-3679, 2003.

[3] C. Huet. Application of variational concepts to size effects in elastic heterogeneous bodies. *J. Phys. Mech. Solids*, 38:813-841, 1990.

[4] David Lane. *Introduction to Statistics*. David Lane, 2003.

[5] J. Dirrenberger, S. Forest, and D. Jeulin. Towards gigantic RVE sizes for 3D stochastic fibrous networks. *Int. J. Solids Structures*, 51:359-76, 2014.

[6] H. Moulinec and P. Suquet. A computational method based on augmented Lagrangians and fast Fourier transforms for composites with high contrast. *Comput. Methods Appl. Mech. Engrg.*, 157:69-94, 1998.

[7] A. Marano, L. Gélébart, and S. Forest. Intragranular localization induced by softening crystal plasticity: Analysis of slip and kink bands localization modes from high resolution FFT-simulations results. *Acta Mater.*, 175:262-75, 2019.

[8] M. Schneider. A review of nonlinear FFT-based computational homogenization methods. *Rev. Perspect. Mech.*, 232:2051-2100, 2021.

[9] G. Rousselier, F. Barlat, and J.W. Yoon. A novel approach for anisotropic hardening modeling. Part I: Theory and its application to finite element analysis of deep drawing. *Int. J. Plasticity*, 25:2383-2409, 2009.

[10] R.E. Marki and M. Knezevic. Nonlinear optimization for compact representation of orientation distributions based on generalized spherical harmonics. *J. Mech. Phys. Solids*, 187:105609, 2024.

[11] L. van Wees, K. Shankar, J.N. Fuhg, N. Bouklas, P. Shade, M. Obstalecki, and M. Kasemer. Establishing the relationship between generalized crystallographic texture and macroscopic yield surfaces using partial input convex neural networks. *Materialia*, 36:102151, 2024.

[12] F. Villette, L. Van Wees, M. Kasemer, and R. Quey. Compact representation of texture for the prediction of macroscopic yield surfaces using machine learning methods. In *ICOTOM20*, 2024.

[13] R. Quey, A. Villani, and C. Maurice. Nearly uniform sampling of crystal orientations. *J Appl Crys*, 51:1162-1173, 2018.

[14] R. Quey, P.R. Dawson, and F. Barbe. Large-scale 3D random polycrystals for the finite element method: generation, meshing and remeshing. *Comput. Methods Appl. Mech. Eng.*, 200:1729-45, 2011.
[15] Neper. https://neper.info.

[16] C. Lautensack and S. Zuyev. Random Laguerre tessellations. *Adv. Appl. Prob.*, 40:630–650, 2008.

[17] Q. Duan, D.P. Kroese, T. Brereton, A. Spettl, and V. Schmidt. Inverting Laguerre tessel- lations. *Comp. J.*, 57:1431–1440, 2014.

[18] A. Spettl, T. Werz, C.E. Krill III, and V. Schmidt. Parametric Representation of 3D Grain Ensembles in Polycrystalline Microstructures. *J. Stat. Phys.*, 154:913–928, 2014.

[19] C. Zhang, A. Suzuki, T. Ishimaru, and M. Enomoto. Characterization of three-dimensional grain structure in polycrystalline iron by serial sectioning. *Metal. Mater. Trans. A*, 35A:1927–1953, 2004.

[20] K.M. Döbrich, C. Rau, and C.E. Krill III. Quantitative Characterization of the Three-Dimensional Microstructure of Polycrystalline Al-Sn using X-Ray Microtomography. *Metal. Mater. Trans. A*, 35A:1953–1961, 2004.

[21] D.J. Rowenhorst, A.C. Lewis, and G. Spanos. Three-dimensional analysis of grain topology and interface curvature in a β-titanium alloy. *Acta Mater.*, 58:5511–5519, 2010.

[22] T.H. Rowan. *Functional Stability Analysis of Numerical Algorithms*. PhD thesis, 1990.

[23] NLopt. http://ab-initio.mit.edu/wiki/index.php/NLopt.

[24] W.R. Hamilton. On quaternions or a new system of imaginaries in algebra. *Philos. Mag.*, 1844–1850.

[25] A. Morawiec. *Orientations and Rotations*. Springer, 2004.

[26] J.J. Thomson. On the structure of the atom: an investigation of the stability and periods of oscillation of a number of corpuscles arranged at equal intervals around the circumference of a circle; with application of the results to the theory of atomic structure. *Philos. Mag.*, 7:237–265, 1904.

[27] J. Barzilai and J Borwein. Two-point step size gradient method. *IMA. J. Numer. Anal.*, 8:141–148, 1988.

[28] S.I. Ranganathan and M. Ostoja Starzewski. Universal elastic anisotropy index. *Phys. Rev. Lett.*, 101:055504, 2008.

[29] H.B. Huntington. The elastic constants of crystals. *Solid State Physics*, pages 213–351, 1958.

[30] JJ Adams, DS Agosta, and RG Leisure. Elastic constants of monocrystal iron from 3 to 500 K. *J. Appl. Phys.*, 100:113530, 2006.

[31] D. Tromans. Elastic anisotropy of hcp metal crystals and polycrystals. *Int. J. Recent Res. Appl. Stud.*, 6:462–83, 2011.

[32] F. Barbe, L. Decker, D. Jeulin, and G. Cailletaud. Intergranular and intragranular behavior of polycrystalline aggregates. Part 1: F.E. model. *Int. J. Plast.*, 17:513–536, 2001.
[33] F. Barbe, S. Forest, and G. Cailletaud. Intergranular and intragranular behavior of poly-crystalline aggregates. Part 2: Results. *Int. J. Plast.*, 17:537–563, 2001.

[34] Z-set. http://www.zset-software.com.