Leon Bilton

Kinematic models of crystalline preferred orientation

Last modified on March 10, 2023

Kinematic models for simulating flow-induced CPO (crystalline preferred orientation) in minerals are attractive due to their relatively low computational cost. Rather than accounting for full stress compatibility between individual grains, aggregates are approximated as effective media. Results from kinematic CPO models have been shown to reproduce experimental observations for CPO in upper mantle materials to within acceptable error margins (e.g. Boneh et al., 2015). In this article, I provide a brief overview of kinematic CPO models, focusing on D-Rex in particular (Kaminski & Ribe, 2001), which I am currently (2022) using in my research.

For analytical derivations of single-grain rotation rates, see "Testing the D-Rex model of flow-induced preferred orientation". Mathematical notation should follow my usual conventions.

Contents:

  1. Flow-induced preferred orientation
  2. Polycrystal CPO models
  3. The D-Rex CPO model

Flow-induced preferred orientation

Geological materials experience extreme stresses and temperatures in the Earth. For example, in the Earth's mantle layer, physical conditions are sufficient to create convective currents that transport large volumes of material by the mechanism of plastic deformation. Plastic deformation occurs when stress surpasses a critical threshold, causing permanent changes in bulk material shape. Strain is a measure of such material deformations.

In crystalline minerals, significant strains (involving changes in length of more than ~1%) must be accommodated at the microscopic scale by appropriate restructuring of the crystal lattice. Various processes can occur in order for a mineral to change its internal geometry. The process of dislocation creep that can occur during large, prolonged strains, involves gradual motions within individual mineral grains along crystalline symmetry planes.

dislocation_creep_2d
Simple diagram of dislocation creep process using the ball-and-stick atom representation. The red ⊥ symbol shows a dislocation, which is moving through the crystal lattice. The shaded filament and finite strain ellipse (FSE) demonstrate the microscopic and macroscopic deformation, respectively. Harpoons indicate the macroscopic shear. The dashed red line lies in the slip plane, which rotates until it is coincident with the shear plane.

Dislocation creep is observed in sufficiently large grains, and manifests as slipping motion along an entire plane within the crystalline grain. The 2D slip diagram provides a visual resemblance of the process. As a consequence of planar slip, the crystalline lattice re-aligns in a predictable way. These changes are reflected in the bulk aggregate as a statistical crystalline preferred orientation (CPO). CPO is measurable either directly, in recovered mineral samples, or indirectly, using seismological methods. Interpretation of observed CPO is complicated by the fact that the speed of re-alignment in typical Earth materials is slow, and the bulk flow direction may change before re-alignment is complete. This means that full numerical simulations of the process need new information about the macroscopic flow at every timestep.

The direction of slip along a dislocation plane is limited by the molecular structure of the grain. Often, slip can only occur in a single direction. The combination of a slip plane and a slip direction is called a slip system. Slip on any single slip system can only induce rotations about one axis in a single direction, which is not sufficient to accommodate arbitrary strains in three dimensions. Denoting the normal and parallel unit vectors of a slip system as \(\bm{n}\) and \(\bm{l}\), the strain due to dislocation creep satisfies

\[ ε_{ij} = \sum_{k} \frac{1}{2}γ^{k} \left(n^{k}_{j}l^{k}_{i} + n^{k}_{i}l^{k}_{j}\right) \]

where the superscript \(k\) is a unique index for the slip system, and \(γ^{k}\) is an unknown slip activity on the \(k\)-th independent system. Because \(\bm{n}\) and \(\bm{l}\) are orthogonal, the sum \(n^{k}_{j}l^{k}_{i} + n^{k}_{i}l^{k}_{j}\) for slip system \(k\) in the grain-local frame will only be nonzero for one pair of \((i,j)\) indices. Volume conservation implies that

\[ \bm{ε}_{11} + \bm{ε}_{22} + \bm{ε}_{33} = 0 \]

which leaves five independent components of the symmetric strain tensor. This means that at least five independent slip systems are required to accommodate arbitrary strains exclusively by dislocation creep.

Other mechanisms may also contribute to strain accommodation. Crystal point-defects such as interstitial atoms may migrate through or around grains. Grains may also slide past each other to change the aggregate geometry. Furthermore, grain sizes can be altered by a variety of micro-physical processes. These other mechanisms are of less interest from a geodynamical point of view, because they do not encode flow patterns by signature lattice rotations. Fortunately, in the bulk of Earth's upper mantle layer, dislocation creep is the dominant deformation mechanism.

Polycrystal CPO models

Numerical methods of calculating CPO that track the deformation of individual grains are classified as "polycrystal" models. In contrast to the simpler finite strain ellipsoid (FSE) models, polycrystal models allow for the full orientation history of grains to be simulated. This is necessary in order to quantify parameters such as Grain Orientation Lag (GOL) (Kaminski & Ribe, 2002).

Polycrystal CPO models can be divided into three classes: full-field models, mean-field models and kinematic models (Fraters & Billen, 2021). Full-field models account for stress compatibility between grains, and provide the most realistic and detailed results. Due to the large number of computations involved, these models are not suitable for use within geodynamic models. Mean-field models are faster alternatives that retain the ability to describe interactions between mineral phases, such as the Visco-Plastic Self-Consistent (e.g. Lebensohn & Tomé, 1993) and Second-Order Self-Consistent (e.g. Ponte Castañeda, 2002) models. Finally, kinematic models are the most approximate, but also the fastest models. Mineral phases are assumed to respond independently to imposed stress, and CPO evolution is investigated by considering interactions between each grain and an effective medium with averaged properties for the chosen mineral phase.

The D-Rex kinematic model is an established kinematic model for olivine and olivine-enstatite aggregates. It is presented in detail in Kaminski & Ribe (2001) and Kaminski et al., (2004). At the time of writing (2022), the most recent implementation of the model, which includes an automatic fabric selector to handle changes in olivine fabric type, is discussed in Fraters & Billen (2021).

The D-Rex theory is valid for mineral phases with at most four independent slip systems. Because (1) does not completely describe the aggregate strain in this case, dynamic recrystallisation processes are assumed to be active as well. Dislocation densities are used to determine slip activation, so the model only applies for grains that are in a steady-state of dislocation energy. Therefore, the model should not be used in the case of very small strains.

The D-Rex CPO model

D-Rex is a kinematic model of CPO development in polycrystal aggregates developed by Kaminski & Ribe (2001). A reference implementation was originally provided in the Fortran language, and has been widely used or re-implemented by later researchers, e.g. Boneh et al. (2015), Hedjazian et al. (2017) and Fraters & Billen (2021). In this section I provide a brief overview of the basic theory.

Stresses in the bulk material induce dislocation creep at the microscopic scale. Combined with auxiliary recrystallisation mechanisms, this allows for microscopic strains to develop. Bulk strain is comprised of microscopic strain and passive rotation of entire grains (rigid-body rotation).

The model describes CPO evolution in a crystal aggregate located at some point in the macroscopic flow field. Deformation of the aggregate is described by a tensor of deformation gradients \(\bm{F}\) which can be defined as

\[ F_{ij} = δ_{ij} + ∂_{j}u_{i} \]

for the displacement vector \(\bm{u}\). Deformation evolves according to the differential equations

\[ D(F)^{ij} = L^{i}_{\phantom{i}k}F^{kj} \]

where \(D\bm{F}\) denotes the material derivative, describing changes within the Lagrangian reference frame. \(\bm{L}\) is a velocity gradient tensor that can be decomposed into symmetric and antisymmetric components as

\[ \dot{ε}_{0}\bm{L} = \dot{ε}_{0}(\dot{\bm{ε}} + \bm{Ω}) \]

where \(\dot{ε}_{0}\) is the strain rate scale, which is removed to achieve a convenient dimensionless representation. The tensors \(\dot{\bm{ε}}\) and \(\bm{Ω}\) are called strain rate and vorticity, respectively.

The aggregate consists of a fixed number of microscopic grains. The microscopic velocity gradient tensor \(\bm{d}\) for each grain is also composed of symmetric and antisymmetric components,

\[ \bm{d} = γ\bm{G} + \bm{ω} \]

where \(γ\) is a characteristic slip rate factor (for the softest slip system), \(\bm{G}\) is a relative slip rate tensor and \(\bm{ω}\) describes a rotation rate of the crystallographic axes. This includes a component due to passive, rigid-body rotation. Because of the assumption that less than five slip systems are active, it is not possible to achieve \(\bm{L} = \bm{d}\). However, dislocation creep in the grains will accommodate the aggregate strain rate as much as possible, which is expressed as the minimisation of

\[ \sum_{n=1}^{N} f_{n} \begin{Vmatrix} \bm{L} - \bm{d}_{n}\end{Vmatrix}^{2} \]

for an aggregate of \(N\) grains (Kaminski & Ribe, 2001). This is a sum of squared Frobenius norms of the residual velocity gradients. The prefactors \(\bm{f}\) are weighting factors which encode the relative grain sizes.

The remaining deformation is assumed to occur via grain-boundary migration and subgrain rotation. These mechanisms can modify the grain size (and thus \(\bm{f}\)), as they are instances of dynamic recrystallisation. Their combined effects are modelled by the formula

\[ d_{t^{\ast}}\bm{f} = - M^{\ast} X (E^{\ast} - \overline{\bm{E}}^{\ast}) \]

where \(d_{t^{\ast}}\) is a derivative in dimensionless time, i.e. \(t^{\star} = \dot{ε}_{0} t\). The asterisk is also used with other variables here to denote nondimensionalisation. The grain boundary mobility, \(M^{\ast}\), is an empirically determined constant. \(X\) is the volume fraction of the mineral phase, which is irrelevant (\(X=1\)) for single-phase aggregates. The term in parentheses denotes a residual between the strain energy of the grain, \(E^{\ast}\), and the average strain energy of the aggregate, \(\overline{\bm{E}}^{\ast}\). Strain energy is a potential energy that depends on the density and distribution of dislocations. Here, the relationships are once again empirical, and I refer the interested reader to the cited literature.

Computation of (7), which relies on calculation of the quantities \(γ\), \(\bm{G}\) and \(\bm{f}\), yields the rotation rate of the grain axes, \(\bm{ω}\). The updated grain axes orientations \(\bm{a}\) are given by

\[ d_{t^{\ast}}\bm{a} = \bm{ω}^{⊺} \bm{a} \]

I have provided a high-level overview of the theory as it is currently employed. Details are omitted but hopefully the conceptual basis is clear. I may expand on this article with another section about tracking finite strain of individual grains. This functionality is also present in the reference implementation, but the method is only discussed in an earlier publication (Ribe, 1992).

References

Bibliography

© Leon Bilton. Last modified: March 10, 2023. Website built with Franklin.jl and the Julia programming language.