# Graph Neural Network Modeling of Grain-scale Anisotropic Elastic Behavior
using Simulated and Measured Microscale Data

Darren C. Pagan1*††*Corresponding Author Calvin R. Pash1 Austin R. Benson2
Matthew P. Kasemer3 1Pennsylvania State University, University Park, PA 16802
USA  2Cornell University, Ithaca, NY 14853 USA  3University of Alabama,
Tuscaloosa, AL 35487 USA

###### Abstract

Here we assess the applicability of graph neural networks (GNNs) for
predicting the grain-scale elastic response of polycrystalline metallic
alloys. Using GNN surrogate models, grain-averaged stresses during uniaxial
elastic tension in Low Solvus High Refractory (LSHR) Ni Superalloy and Ti
7wt%Al (Ti-7Al), as example face centered cubic and hexagonal closed packed
alloys, are predicted. A transfer learning approach is taken in which GNN
surrogate models are trained using crystal elasticity finite element method
simulations and then the trained surrogate models are used to predict the
mechanical response of microstructures measured using high-energy X-ray
diffraction microscopy. The performance of using various microstructural and
micromechanical descriptors for input nodal features to the GNNs is explored.
The effects of elastic anisotropy on GNN model performance and outlooks for
extension of the framework are discussed.

††journal: arXiv

##  1 Introduction

The local micromechanical response of grains embedded within a polycrystal is
dictated not only by an isolated grain’s features (e.g., defect state or
crystallographic lattice orientation), but also by the features and mechanical
response of adjoining grains (i.e., the local grain ‘neighborhood’). Various
linking hypotheses or methods in which grains are embedded within homogeneous
matrices—most notably the Eshelby method [1]—attempt to capture the average
features of grain-scale response, but, by construction, do not consider
variation of behavior created by the interactions of specific grains within
their local neighborhoods. The relatively recent ability to explicitly model
grains and grain neighborhoods in three-dimensional (3D) polycrystals using
both finite element [2] and spectral [3] methods has allowed these
neighborhood effects to be more thoroughly explored. However, the use of these
methods comes at a sometimes-significant computational cost. Increasing
computational complexity, particularly with increases in microstructural
fidelity or the inclusion of various plastic deformation mechanisms in
modeling efforts, limits the number of microstructural configurations that can
be tested, which consequently limits the ability of these models to be
embedded within larger scale simulations (considering current computational
capabilities). To address these challenges, low-computational-cost surrogate
models which can rapidly evaluate micromechanical response and evaluate large
regions of microstructural parameter space are necessary. In this work we
propose and demonstrate that the fundamental network (or graph) structure of
polycrystals make them candidates for surrogate mechanical modeling through
graph neural networks (GNNs) [4, 5, 6]. We demonstrate the utility of this
surrogate modeling with GNNs trained to predict grain-scale elastic response
in two example alloy systems. The GNNs are trained with microscale crystal
elasticity finite element method (CEFEM) simulations, and then tested against
grain-scale elastic response measured using high-energy X-ray diffraction
microscopy (HEDM).

Supervised machine learning has become a primary choice for generating
surrogate models for predicting the mechanical properties and performance of
engineering alloys [7]. At the microscale, multiple efforts have utilized
convolutional neural networks (CNNs) to predict deformation fields at sub-
grain length scales in virtual polycrystals [8, 9, 10]. In these efforts,
polycrystals are represented by a grid of voxels containing microstructural
descriptors such as lattice orientation, and the CNN learns spatial
correlations between the voxels to generate full-field predictions. In
particular, CNNs take advantage of the grid structure of the data to learn
filter structures comprised of series of weights. These weights group (pool)
together neighboring grid values (in this case microstructural features) to
capture and predict the effects of neighborhood. The results of this approach
are very promising, but an issue is that (most) materials are not naturally
structured in a grid-like fashion (although due to data collection strategies,
materials are often represented as such). In the cases of polycrystalline
materials at the microscale, the grains themselves are more naturally
represented in an unstructured, connected format often referred to as a graph.

Broadly, a graph is a structure in which ‘vertices’ or ‘nodes’ are connected
through ‘edges’. Both nodes and edges can further be described via associated
attributes or features. In the case of a polycrystal represented via a graph,
grains are considered as nodes, while grain boundaries are considered as
edges. Features of nodes (grains) could include local micromechanical
properties such as stiffness or strength, as well as microstructural features
such as lattice orientation or dislocation density. Features of edges (grain
connections or boundaries) may be distance between grains or grain boundary
characteristics. GNNs, which are the focus of this work, adapt the data
pooling technique of CNNs on structured data to the pooling of neighboring
data on unstructured graph data. Instead of pooling features from neighboring
grid points, features are pooled from connected nodes [11]. A recent study has
demonstrated the utility of GNNs for predicting magnetorestriction response in
a polycrystalline material [12].

In this study, we apply GNNs to predict the elastic response in two alloy
systems: Low Solvus High-Refractory (LSHR) Ni superalloy and Ti 7wt%Al
(Ti-7Al), representing cubic and hexagonal elastic responses, respectively.
Two GNN surrogate models for predicting LSHR and Ti-7Al grain-averaged elastic
response (grain stress tensor components along the loading direction) that
utilize Gaussian Mixture convolutions are implemented, trained, and
preliminarily tested using CEFEM simulations that explicitly model grain
microstructure. The GNN models trained with CEFEM simulations are then
transferred to a new data domain to predict the grain-averaged elastic
response in polycrystals with microstructures experimentally measured with
HEDM. As part of the GNN model training effort, learning rates of the models
and accuracy of using various nodal features for stress predictions are
explored.

In this paper, vectors are generally lower-case bold characters ($\bm{a}$),
second order tensors are upper-case bold characters ($\bm{A}$), and fourth
order tensors are underlined, bold characters ($\underline{\bm{A}}$). Unless
otherwise noted, quantities are expressed in the sample frame. Prime
characters are generally used to indicate quantities in subsequent layers
within GNNs, i.e., $a^{\prime}$ for Hidden Layer n+1, while an overbar
$\bar{a}$ indicates an average.

##  2 Methods

In this section, we give a broad overview of the various methods employed in
this work for GNN training and accuracy evaluation including: virtual sample
generation, CEFEM modeling, HEDM, and GNN surrogate modeling. In addition,
details regarding the various training and testing data used are provided. A
schematic of the various components of the effort, displaying LSHR data, are
given in Fig. 1. Briefly, a transfer learning approach [13] is taken in which
the GNN models for LSHR and Ti-7Al are trained using simulated data from
microscale CEFEM modeling (the Source Domain) and then transferred to predict
the mechanical response in a microstructure measured experimentally via HEDM
(the Target Domain).

Figure 1: Connectivity of the various efforts used to evaluate the accuracy of
GNN surrogate models for predicting the stress response in individual grains
during elastic loading.

###  2.1 Virtual Sample Generation

Virtual polycrystals upon which mechanical deformations are imposed and used
to create microstructure graphs are generated using the Neper polycrystal
generation and meshing software [14, 15, 16]. Generally, Neper utilizes
Laguerre tessellations to generate polycrystalline samples. Laguerre
tessellations produce convex, space-filling grains in a user-specified sample
domain. Specifically, Neper allows for user-specified target distributions of
grain size and grain shape. This affords the ability to create a wide range of
microstructures with various geometric features. An attendant finite element
mesh is then generated via Neper using Gmsh [17], in which the geometric
features of the microstructure (tessellation) are preserved.

Here, Neper is used to generate 50 virtual polycrystals (25 LSHR and 25
Ti-7Al) for elastic finite element simulations. The simulated domains are
$1\text{\,}\mathrm{mm}$ $\times$ $1\text{\,}\mathrm{mm}$ $\times$
$3\text{\,}\mathrm{mm}$, each containing 1,500 grains. Each polcyrystal had
approximately 16,500 shared grain boundaries which serve as edges. Grain size
and shape distributions are set to create nominally equiaxed grains with
diameters of approximately $150\text{\,}\mathrm{\SIUnitSymbolMicro m}$ and
minimal spread. Each polycrystal is meshed with ten node tetrahedral elements
and approximately 120,000 elements per sample. Orientations are assigned
randomly to grains from the cubic and hexagonal fundamental regions for the
LSHR and Ti-7Al virtual specimens respectively.

###  2.2 Polycrystal Anisotropic Elasticity Data

For generating training data for the GNN surrogate elasticity model
deformation response, we utilize the finite element solver within FEPX [18,
19]. Generally, FEPX considers the elasto-viscoplastic deformation response of
single crystals (grains) belonging to explicit representations of
polycrystalline aggregates. However, here, viscoplasticity is inhibited by
ensuring grain level stresses are significantly lower than the yield point. In
the model, grain-to-grain interactions are assumed to be rigid (no grain
boundary sliding or separation). As the implementation and use of these models
are well established [20, 2, 21, 22], and as this study considers elastic
deformation, only a truncated description of the model as implemented in FEPX
is included below. Please refer to [18] for a complete description of
kinematics, models, and finite element implementation. For the uniaxial
deformation studied in this work, loading is along the $\bm{z}$ direction in
the sample frame.

At each point in a finite element mesh, FEPX considers the elastic response to
be governed by the anisotropic form of Hooke’s law:

$\bm{\sigma}=\underline{\bm{C}}(\bm{r}):\bm{\varepsilon}\quad,$ (1)

where $\bm{\sigma}$ is the stress, $\underline{\bm{C}}$ is the lattice-
orientation-dependent elastic stiffness tensor, $\bm{r}$ is the lattice
orientation (a coordinate transformation from crystal frame to sample frame)
of a given crystal, and $\bm{\varepsilon}$ is the strain, assumed to be fully
elastic. Due to symmetry, crystals with cubic symmetry have three independent
constants ($C^{C}_{11}$ , $C^{C}_{12}$ , and $C^{C}_{44}$ in Voigt notation in
the crystal frame), while crystals with hexagonal symmetry have five
independent constants ($C^{C}_{11}$ , $C^{C}_{12}$ , $C^{C}_{13}$ ,
$C^{C}_{33}$ , and $C^{C}_{44}$ in Voigt notation in the crystal frame). We
note that due to the deformation decomposition formulation in FEPX, in
hexagonal crystals further dependence is required in the form of
$C^{C}_{33}=C^{C}_{11}+C^{C}_{12}-C^{C}_{13}$. The elastic moduli employed for
LSHR [23] in this work are: $C^{C}_{11}=247$, $C^{C}_{12}=147$, $C^{C}_{44}=$
$125\text{\,}\mathrm{GPa}$, while the Ti-7Al [24] moduli are: $C^{C}_{11}=162$
, $C^{C}_{12}=92$ , $C^{C}_{13}=69$ , $C^{C}_{33}=185$ , and $C^{C}_{44}=$
$45\text{\,}\mathrm{GPa}$.

In uniaxial deformation, the effective stiffness along the loading direction
can be approximated using a directional modulus $E(\bm{r})$. Generally, in
grains embedded in polycrystals, the stress along the loading direction is
correlated to the directional modulus [25]. Under uniaxial stress along the
$\bm{z}$ in the sample frame, the directional modulus is defined as:

$E(\bm{r})=\frac{\sigma_{zz}}{\varepsilon_{zz}(\bm{r})}\quad.$ (2)

To further describe the anisotropy of the grains embedded within a polycrystal
during uniaxial deformation, we define effective transverse contraction ratios
$\nu_{x}(\bm{r})$ and $\nu_{y}(\bm{r})$:

$\nu_{x}(\bm{r})=\frac{\varepsilon_{xx}(\bm{r})}{\varepsilon_{zz}(\bm{r})}\quad,$
(3)

and

$\nu_{y}(\bm{r})=\frac{\varepsilon_{yy}(\bm{r})}{\varepsilon_{zz}(\bm{r})}\quad.$
(4)

Each of the 50 virtual polycrystals generated (see §2.1) are deformed
elastically in FEPX. As all uniaxial loading linear elasticity solutions can
simply be scaled to account for increasing or decreasing deformation, each
polycrystal is only deformed with a single deformation increment to an applied
strain of 0.1% (chosen arbitrarily). Minimal displacement boundary conditions
are employed on the top and bottom surfaces of the specimen to prevent rigid
translation and rotation of the virtual specimens without impeding contraction
of the specimen perpendicular to the loading direction.

###  2.3 High-Energy X-ray Diffraction Microscopy Data

For final GNN evaluation, graph data derived from experimentally measured
microstructures and micromechanical response measured using HEDM are utilized.
HEDM is comprised of two variants, near-field and far-field, capable of non-
destructively characterizing the microstructure and micromechanical response
of polycrystalline materials [26, 27, 28]. The commonality between the two
techniques is the utilization of forward projection simulations of X-ray
diffraction peaks of individual grains to reconstruct information about the
local lattice state of crystalline materials in 3D. LSHR and Ti-7Al samples
were probed using the near-field variant to measure the grain structure,
orientation, and connectivity of the grains, and the far-field variant was
used to determined the stress in the same grains during elastic loading.
Detailed descriptions of the data collection for the LSHR sample can be found
in [29] and for the Ti-7Al sample in [30], but summaries of the methodology
and specimens are given below.

The near-field variant (nf-HEDM) is capable of reconstructing a 3D voxelized
distribution of lattice orientation with resolution on the order of
$\text{\,}\mathrm{\SIUnitSymbolMicro m}$ from series of diffraction images
collected as a specimen is rotated. In the near-field variant, the detector is
placed approximately $5\text{\,}\mathrm{mm}10\text{\,}\mathrm{mm}$ away from
the specimen, making the measurements sensitive to the spatial locations of
diffracting volumes [31]. For this work, the LSHR and Ti-7Al microstructures
were reconstructed with $5\text{\,}\mathrm{\SIUnitSymbolMicro m}$ voxel
spacing. The reconstructed volume for the LHSR specimen was
$1.0\text{\,}\mathrm{mm}$ $\times$ $1.0\text{\,}\mathrm{mm}$ $\times$
$0.5\text{\,}\mathrm{mm}$ with the short direction aligned with the loading
direction. Similarly, the reconstructed volume for the Ti-7Al specimen was
$1.0\text{\,}\mathrm{mm}$ $\times$ $1.0\text{\,}\mathrm{mm}$ $\times$
$0.5\text{\,}\mathrm{mm}$, again with the short direction along the loading
direction.

The far-field variant (ff-HEDM) can reconstruct the average grain orientation,
position, and elastic strain state of grains embedded in a polycrystal during
_in situ_ loading [32]. In this HEDM variant, a large area detector is placed
approximately $1\text{\,}\mathrm{m}$ away from the specimen. This positioning
provides more sensitivity to peak shifts due to changes in orientation and
strain state, as opposed to locations of diffraction events as utilized for
nf-HEDM reconstructions. Here ff-HEDM measurements collected prior to and
during _in situ_ elastic loading to applied strains near 0.001 — 0.0015 and
0.0009 for the LSHR and Ti-7Al respectively — are used for GNN prediction
evaluation. The full elastic strain tensors from individual grains are then
used to calculate the stresses in each of the grains, including the component
of stress along the loading direction which is of interest here. The
uncertainty per elastic strain component for HEDM measurements are generally
reported to be $10^{-4}$ [33].

###  2.4 Graph Neural Network Modeling

As previously described, the surrogate models employed in this work are graph
neural networks (GNNs). The GNNs are comprised of layered graphs in which the
features of each node in subsequent graph layers are weighted combinations
(convolutions) of $N_{F}$ nodal features from neighbors in the previous graph
layer. The weights may be independent or functions of nodal or edge features
and are learned throughout the training process. Similar to other neural
network formulations, nonlinear activation functions (functions which control
the flow of information between layers, i.e., nominally turning on and off)
control the passing of nodal feature information between layered graphs. The
input data is the graph itself — a set of nodes and a set of edges — along
with features associated with the nodes and edges expected to be useful for
predicting stress. Here graphs are generated from CEFEM or HEDM data using a
series of custom Python scripts.

Nodal features represent different microstructural descriptors or local
mechanical properties. Various nodal features, and combinations of nodal
features, were explored, including the directional modulus ($E(\bm{r})$),
contraction ratios ($\nu_{x}(\bm{r})$,$\nu_{y}(\bm{r})$), volume ($V$), and
lattice orientation ($\bm{r}$, Rodrigues paramaterization) of each grain. We
note that no appreciable accuracy difference was observed utilizing different
orientation parameterizations (e.g., Euler angles or other angle-axis
parameterizations). These nodal features alone could be useful for predicting
stress. However, rather than just treating all nodes independently, the GNN
combines these features with the same features of connected nodes in the
graph, along with additional edge features associated with the connections.
The edge features used here for a connected pair of grains $i$ and $l$ are the
coordinates of the vector between the grain centroids $\bm{p}_{i}$ and
$\bm{p}_{l}$:

$\bm{e}_{il}=\bm{p}_{i}-\bm{p}_{l}\quad.$ (5)

Rather than just looking at immediate connections, the GNN also works
recursively, learning a vector representation of node features at subsequent
layers. For this work, a common GNN architecture was utilized with two hidden
graph layers. A schematic of the architecture is given in Fig. 2. Each node
within the two hidden graph layers contains 16 features obtained from a
convolution of the previous layer’s features over neighboring nodes. For the
first hidden graph layer, the convolution is over the input nodal features.
Two hidden graph layers provide substantial flexibility for learning the non-
linear relationship between grain lattice orientation and elastic response.
After the two hidden graph layers, the second 16-dimensional vector
representation of the nodal features is passed through a final dense layer
that maps to a single nodal feature, which is the predicted stress tensor
component $\sigma_{zz}$. In total, there are four layers: one input graph with
the input nodal features, two hidden convolutional graphs with learned
features, and the final mapping layer. While the number of nodal features in
the two hidden graph layers could be slightly oversized for the problem at
hand, over-fitting was only observed when very small amounts of training data
were employed, as will be seen.

The graph convolution operator used in the elasticity surrogate models is a
Gaussian Mixture Model [34], implemented in PyTorch Geometric [35], with the
form (no sum over indices $i$ and $j$):

$x^{\prime}_{ij}=\frac{1}{N_{N}}\sum_{k=1}^{N_{F}}\theta_{jk}\sum_{l\in
N(i)}w_{k}(\bm{e}_{il})x_{lk}\quad,$ (6)

where $x$ is a nodal feature (either the initial inputs, or the learned vector
representation in the hidden layers), $i$ indicates node (grain) number
ranging from 1 to $N_{G}$, $k$ indicates a nodal feature in a current layer
ranging from 1 to $N_{F}$, $j$ indicates a nodal feature in a subsequent layer
ranging from 1 to $N_{F}^{\prime}$, $w$ is a Gaussian weighting function, $l$
indexes a neighbor to grain $i$, $N(i)$ is the set of indices of nodes
connected to $i$, $\bm{e}$ is the vector of edge features between nodes $i$
and $l$ from Eq. 5, and $\theta$ is a (learned) weight for a given nodal
feature in the current layer.

The weighting functions are comprised of $N_{K}$ Gaussian kernels:

$w_{k}(\bm{e}_{il})=\sum_{M=1}^{N_{K}}\exp\left(-\frac{1}{2}(\bm{e}_{il}-\bm{\mu}_{m})^{T}\bm{\Sigma_{m}}(\bm{e}_{il}-\bm{\mu}_{m})\right)\quad,$
(7)

with $\bm{\mu}_{m}$ and $\bm{\Sigma_{m}}$ being Gaussian kernel and diagonal
covariance matrices. In total, $\theta$, $\bm{\mu}$, and $\bm{\Sigma}$ are
coefficients learned during the training process. Using the 3D vector between
grain centroid as edge features (Eq. 5), each layer in the neural network
requires $N^{\prime}_{F}\times(N^{\prime}_{F}\times(1+N_{K}\times(3+3)))$
coefficients to be learned. The choice of a Gaussian Mixture Model was
informed by our physical understanding of deformation compatibility and
mechanical equilibrium. In this model, GNN has sufficient degrees of freedom
to weight various positions around a grain differently (e.g., parallel and
transverse to applied uniaxial load).

For final fitting, seven Gaussian kernels ($N_{K}=7$) were used for the
Gaussian convolution operator in Eq. 7 with the goal of providing the
surrogate model the freedom to weight neighborhood grains arranged in various
positions differently (i.e., parallel and transverse to loading). Leaky ReLU
was chosen for the activation function, with a scaling value of -0.1 for input
values less than 0 as it provided better accuracy than a standard ReLU
activation function. Training of the surrogate models was performed by
minimizing the mean square error (MSE) between grain average stress components
in the training data and those predicted by the GNN surrogate model using Adam
stochastic gradient descent with a learning rate parameter of 0.01. Training
of a surrogate model for 10,000 epochs nominally takes 2-15 minutes depending
on the amount of training data (1,500 grains to 3,000 grains) using an NVIDIA
Quaddro GPU with 5 GB of memory.

Figure 2: Schematic of the GNN surrogate architecture used for this work. An
input graph with various sets of nodal features maps to a final output layer
for predicting stress along the loading direction $\sigma_{zz}$

##  3 Results

The results are divided into three subsections. The first subsection covers
the learning behavior of GNN surrogate models using CEFEM data. The second
subsection details how the trained GNN surrogate models perform predicting
micromechanical response in comparison to reserved CEFEM data is described.
The final subsection analyzes the performance of the trained surrogate GNN
models in predicting mechanical response of microstructures characterized
experimentally.

###  3.1 Model Training

After generating the various graph data from both CEFEM simulations and HEDM
results (Sec 2.2 and 2.3) for GNN surrogate model training, a study was
completed to examine the performance of the GNNs predicting the stress
response along the loading direction ($\sigma_{zz}$) in individual grains. For
GNN model training, graphs generated from various numbers of CEFEM simulations
(1, 5, 10, and 20 corresponding to 1,500, 7,500, 15,000, and 30,000 grains)
were used to train the LSHR and Ti-7Al GNN surrogates. Respectively for each
training scenario, 1, 2, 3, and 4 CEFEM simulations (1,500, 3,000, 4,500, and
6,000 grains) were reserved for evaluating the accuracy of the surrogate model
predictions. Accuracy of surrogacy models in predicting simulated stresses
were quantified using the mean of a 1-Norm error $e$:

$\bar{e}(\sigma_{zz})=\frac{1}{N_{G}}\sum\frac{|\sigma_{zz,i}^{\mathrm{SIM}}-\sigma_{zz,i}^{\mathrm{GNN}}|}{\sigma_{zz,i}^{\mathrm{SIM}}}\quad.$
(8)

Fig. 3 shows the learning rates (epoch vs error) for LSHR elastic response
using the various numbers of training data sets described above. Fig. 3a shows
the learning rate using the directional modulus $E(\bm{r})$ as a nodal
feature, while Fig. 3b shows the learning rates using components of the
Rodrigues vector $\bm{r}$ describing a grain’s lattice orientation as nodal
features. It is important to note that while these quantities are related (the
directional modulus is a function of the grain lattice orientation), the
directional modulus is a micromechanical property that should have a linear
relationship to the stress state, while the components of the lattice
orientation are a microstructural feature that should have a non-linear
relationship to the stress state. Solid lines show the mean error $\bar{e}$ of
grain stress predictions of the surrogate model on the training data while the
dashed lines show the mean error in comparison to the reserved testing data
sets. As a benchmark, mean errors as predicted by various mean-field linking
hypotheses, isostress (all grains have the same stress state as the
macroscopic stress, $\sigma_{zz}=\tilde{\sigma}_{zz}$) and isostrain (all
grains have the same strain state as macroscopic strain,
$\sigma_{zz}=E(\bm{r})\tilde{\varepsilon}_{zz}$), are provided. In most cases,
there is no evidence of over-fitting which would be seen as divergence between
mean errors on predictions of training and testing data, except for the case
when only one training data set is used with lattice orientation used as the
nodal feature input. The ‘spikes’ during the learning process are related to
the use of Adam stochastic gradient descent algorithm for fitting, which will
perturb the solution intermittently to try to ensure that a global minimum is
reached.

In Fig. 3, we can see that using the GNN surrogate model with either the
directional modulus or lattice orientation as nodal features outperforms the
predictions of the mean-field theories reaching average errors of
approximately 0.05 (5%) for LSHR, with the lattice orientation performing
slightly better after 10,000 epochs of training. The fact that surrogate
models using directional modulus and lattice orientation have approximately
the same performance indicate that nominally the same information is encoded
into the two grain descriptors. This can be rationalized in that the
directional modulus is calculated using the lattice orientation and the single
crystal elastic moduli as previously mentioned. However, this relationship
between lattice orientation and directional stiffness is non-linear and
naturally the rate at which the relationship is learned is slower, but the GNN
surrogate model does in fact learn the relationship. We note that when using
only a single hidden layer, the surrogate model is incapable of learning this
relationship (not shown).

Figure 3: Learning rates of LSHR grain-scale elastic response using various
numbers of CEFEM simulations for training and testing. Insets show a magnified
view of the learning rate at larger epochs. a) Learning rates using the
directional modulus $E(\bm{r})$ as a nodal feature. b) Learning rates using
the lattice orientation (components of the Rodrigues vector $\bm{r}$) as nodal
features.

Similar to Fig. 3, Fig. 4 shows learning rates of the grain-scale elastic
response using different numbers of CEFEM simulations for training and using
directional modulus (Fig. 4a) and lattice orientation (Fig. 4b) as input nodal
features for Ti-7Al. In Fig. 4b, it can be seen that using 1 and 5 CEFEM
simulations for training along with lattice orientation as nodal features show
signs of over-fitting for Ti-7Al as opposed to only 1 FEM simulation for the
LSHR data. Like the LSHR, using directional modulus or lattice orientation as
input nodal features gives approximately the same performance (0.035 mean
error) at large epochs. The magnitude of the final mean error is lower than
that of the LSHR (0.035 for Ti-7Al versus 0.05 for LSHR) which is discussed
further in §4.1. In addition, the GNN surrogate model only performs slightly
better than using an isostress assumption to determine grain average stresses
in the Ti-7Al which will be discussed.

Figure 4: Learning rates of Ti-7Al grain-scale elastic response using various
numbers of CEFEM simulations for training and testing. Insets show a magnified
view of the learning rate at larger epochs. a) Learning rates using the
directional modulus $E(\bm{r})$ as a nodal feature. b) Learning rates using
the lattice orientation (components of the Rodrigues vector $\bm{r}$) as nodal
features.

###  3.2 GNN Surrogate Model Performance: Source Domain (CEFEM Data)

Besides directional modulus and lattice orientation, the performance of GNN
surrogate models in predicting stresses in the Source Domain (i.e., domain
where the surrogate models are trained using CEFEM data) using other sets of
nodal features were also tested. These include incorporating transverse
contraction ratios $\nu_{x}(\bm{r})$ and $\nu_{y}(\bm{r})$ in addition to the
volume $V$ of grains into the fitting, for a total of four different nodal
feature sets. A summary of the training results for the different nodal
feature sets using 20 CEFEM simulations (30,000 grains) for training and 4
CEFEM simulations for testing (6,000 grains) are shown in Fig. 5. Beyond what
can be seen in Fig. 3 and Fig. 4, the inclusion of transverse contraction
ratios generally do not improve accuracy. This is consistent with the
observation that the directional modulus and lattice orientation nominally
have nominally the same accuracy: if transverse contraction ratios improved
accuracy, it would be expected that lattice orientation would also perform
better (contraction ratios are calculated from orientation). Somewhat
surprisingly, inclusion of grain volume $V$ into the nodal features does not
improve the surrogate model prediction accuracy. For this work, however,
spread of grain size is minimal, and we note that this observation may not be
general. Volume may be a necessary descriptor for accuracy in other
microstructures with larger or bi-modal distributions of grain size.

Figure 5: Summary of learning rates for various nodal features included into
the GNN surrogate models for a) LSHR and b) Ti-7Al.

To examine if there were systematic errors in the GNN surrogate model
predictions, the stress predictions from the surrogate models using lattice
orientation (Rodrigues vectors) as nodal features were compared to stresses
from the CEFEM reserved testing data sets. Comparison of surrogate model
predictions using 20 CEFEM simulations for training and compared to stresses
from 4 testing CEFEM simulations are shown in 2D histograms in Fig. 6a and
Fig. 6b for LSHR and Ti-7Al respectively. The solid red diagonal lines
corresponds to perfect correspondence (100% accuracy) while the dashed and
dotted lines correspond to 5% and 10% error bounds respectively. For the LSHR,
64% of the GNN stress predictions fall within the 5% bound, and 93% of stress
predictions fall within the 10% bound. It should be noted that spread of
stresses across the LSHR grains is significantly larger than that seen in the
Ti-7Al due to the much larger elastic anisotropy. The LSHR fit shows no
evidence of systematic bias and is able to predict stresses withing the full
ranges of stresses generated by the significant elastic anisotropy in LSHR.
For the Ti-7Al, the GNN surrogate model captures most of the grain stress
spread, but appears to be truncating the bounds of the most extreme values of
the prediction. In Fig. 6b, 78% and 98% of the GNN stress predictions fall
within the 5% bound and 10% bound respectively. These observations will be
discussed further in §4.1.

Figure 6: Comparisons of stress predictions along the loading direction
($\sigma_{zz}$) from GNN surrogate models using lattice orientation as nodal
features compared to CEFEM simulations for a) LSHR and b) Ti-7Al.

###  3.3 GNN Surrogate Model Performance: Target Domain (HEDM Data)

After evaluating the performance of the surrogate models in the Source Domain,
the surrogate models were then transferred to the Target Domain where stresses
during elastic loading was predicted from a graph generated from a 3D
microstructure measured using nf-HEDM. As the performance of various nodal
features used for prediction were fairly similar (Fig. 5), only GNN surrogate
models using lattice orientation as nodal features were tested as lattice
orientation is directly measured using HEDM. The stresses predicted by the GNN
surrogate models were then compared to the stresses measured from the same
microstructure during elastic loading using ff-HEDM. Comparison of the GNN
predicted versus ff-HEDM measured stresses for the LSHR and Ti-7Al specimens
are shown in Figs. 7a and 7b, respectively. For comparison, the stresses for
2,896 LSHR grains and 592 Ti-7Al grains were used. For prediction of stresses
using the GNN surrogates, the output stresses are scaled by the ratio of
macroscopic applied strain in the experiments (0.0015 for the LSHR and 0.0009
for the Ti-7Al) to the macroscopic applied strain for the CEFEM simulations in
the Source Domain (0.001). As in the Source Domain, there is no systematic
error in the prediction of stresses for LSHR seen in observed in Fig. 7a.
However, the mean of the errors are higher ($\bar{e}=13.3\%$), with 27% (783)
of grains having errors less than 5% and 49% (1,428) of grain having errors
less than 10%. There are multiple possible contributions to this elevated
error including measurement error, non-perfect application of uniaxial stress
(e.g., small bending or torsion stresses), or more complex grain geometries
and interactions than those created using Laguerre tessellations for the
simulations. In Fig. 7b, we can see the GNN stress predictions are fairly
tightly clustered around the average stress in the specimen (100 MPa). There
appears to be a systematic error in which the GNN predictions for Ti-7Al do
not capture the full spread of grain stresses in the simulation, similar to
observations in the Source Domain, but increased due to the various factors
likely contributing to errors with use of experimental data. Although the
systematic error is present, the Ti-7Al stress predictions still have a lower
mean error ($\bar{e}=9.2\%$) and higher percentage of grains within the 5%
bound (41%, 227 grains) and 10% bound (71%, 392 grains) in comparison to the
LSHR due to the lower elastic anisotropy of the Ti-7Al.

Figure 7: Comparisons of stress predictions along the loading direction
($\sigma_{zz}$) from GNN surrogate models using lattice orientation as nodal
features compared to those measured using ff-HEDM during in situ elastic
loading for a) LSHR and b) Ti-7Al.

##  4 Discussion

Here the efficacy of using GNN surrogate models to predict grain-scale elastic
response on both simulated and experimentally measured microstructures
(represented as graphs) was explored. The GNN models use a graph convolution
(Gaussian Mixture) to capture the inherent non-locality in both mechanical
equilibrium and deformation compatibility which govern mechanical response.
Elasticity was chosen as an initial test case for exploring prediction of
micromechanical response due to our generally strong understanding of the
constitutive relationship (Hooke’s Law) and primary material characteristics
(the single crystal elastic moduli and lattice orientation of a grain)
governing elasticity at the grain-scale. It was found that, once trained, the
accuracy of the GNN surrogate models exceeded that of both isostrain (Taylor)
and isostress (Sachs) mean-field theory assumptions in both a highly
anisotropic FCC alloy (LSHR) and a moderately isotropic HCP alloy (Ti-7Al).
Various sets of microstructural and micromechanical descriptors were tested
and found to have generally similar accuracy after full training (Fig. 5),
although it is noted that the choice of descriptors was informed by our
generally strong understanding of the elastic response. Importantly though,
the two hidden layer architecture employed was demonstrated to be capable of
learning the complex non-linear relationship between lattice orientation and
directional stiffness as evidenced by the similar final performances.

In a novel approach, transfer learning was successful in enabling the
prediction of the mechanical response of microstructures (as represented by
graphs) that were experimentally measured using HEDM, as opposed to only
predicting the response of synthetic microstructures. The ability to move
between various data domains opens up the possibility of the training schemes
that could further improve the accuracy of model. In particular, using both
simulated and experimental data to train GNN surrogate models is a promising
avenue for future efforts as both data types have strengths and weaknesses for
use in learning material response. Large amounts of simulated data is
generally more readily acquired than experimental data, however simulations
are generally not capable of capturing phenomena not included into the model.
Experimental data will naturally represent the physics of the phenomenon at
hand, but is susceptible to measurement error and bias. Creating training data
sets in the future containing both simulation and experimental data can help
to further to improve the accuracy of GNN surrogate model predictions.

###  4.1 Mechanical Anisotropy

At this point, it is worth discussing differences between the elastic
responses of LSHR and Ti-7Al and the utility of GNN surrogate models for
predicting elastic response in general. While the LSHR is crystallographic
more symmetric, it is actually significantly more elastically anisotropic. The
directional moduli in the LSHR varies by a factor of approximately 2.5 (ratio
of maximum to minimum), while the directional moduli of the Ti-7Al only vary
by a factor of approximately 1.5. The value of the use of the GNN surrogate
model in comparison to mean-field theory will increase as a material becomes
increasingly elastic (or plastically) anisotropic and the converse. At the
extreme, if a material is nearly elastically isotropic, assuming the
macroscopic stress state is equal to the grain stress state is highly
accurate. It was observed for the Ti-7Al, that there was a systematic error in
the GNN predictions in both the Source and Target domains. The predictions for
the Ti-7Al model are tightly bound around the average stress indicating the
model is only learning minor adjustments from an isostress mean-field theory,
as opposed to learning full neighborhood effects. However, this is simply a
consequence of the fact that any GNN model will converge to an isostress
prediction as the material approaches isotropy.

###  4.2 Extension of Modeling Efforts

With the ability of the GNNs to predict elastic response shown here, extension
of the approach to predicting plasticity should be considered. In this work,
prediction of mechanical response of each node is static and does not evolve
with time or deformation. To address material evolution that occurs during
plastic deformation, recurrent neural networks (RNNs) can be layered on top of
the GNN surrogate. In the RNN framework, the neural network has a ‘memory’
that transmits output nodal features through time, enabling dynamic processes
to be modeled. This approach has been successful in time-evolving graph
problems such as traffic forecasting [36], but has mostly been developed for
cases where the temporal data is largely stationary (e.g., regular traffic
data) and the nodal features are static (e.g., location information of roads
is static). However, existing approaches for modeling plasticity present a
path forward as there are significant commonalities between nodal features in
a graph and internal state variables characterizing the local microstructure
of material points [37]. One can imagine beyond predicting the mechanical
response at nodes after an loading increment (such as this work), nodal
features that are either implicit descriptors of microstructure (e.g., slip
system strengths) or explicit (e.g., dislocation density) can be predicted.

Once trained, the GNN surrogate models are highly efficient in predicting
mechanical response as they do not require the inversion of large systems of
equations or any non-linear optimization. This provides an opportunity for
embedding into larger-scale models. Efforts to include grain-scale evolution
in the larger-scale deformation or forming operations generally requires
homogenization [38] or a linking hypothesis [39, 40], usually a Taylor
assumption, which will remove local neighborhood effects from consideration,
although inclusion of non-explicit neighborhood effects have been attempted
[41]. While successful for predicting a bulk property such as anisotropic
strength [42] or a microstructural feature such as phase fraction [43], this
approach will be unable to predict mechanical responses, such as fracture and
fatigue, which are dictated by rare, extreme events. Embedding GNN surrogate
models, which naturally incorporate neighborhood effects without a large
computational overhead, within a larger-scale simulation is a possible avenue
for predicting extreme events in full size components, opening a path for
truly microstructurally-sensitive predictions of material failure in
component-scale simulations.

##  5 Summary

A novel transfer learning approach was taken for training and evaluating the
performance of GNN surrogate models. These novel surrogate models provide a
new opportunity for implementation in multiscale mechanical modeling where
micromechanical response can be rapidly evaluated using a GNN surrogate.
Models were used to predict the elastic response of grains in samples deformed
under uniaxial tension. From this work we find that:

  * 1.

GNNs can exceed the accuracy traditional mean-field assumptions for predicting
anisotropic elastic micromechanical response through incorporation of
information regarding local neighborhood.

  * 2.

Micromechanical descriptors such as directional modulus and near equivalent
microstructural descriptors such as lattice orientation provide similar
accuracy for GNN surrogate model predictions of elasticity.

  * 3.

As mechanical response becomes increasingly isotropic, GNN surrogate model
predictions converge toward traditional mean-field predictions.

  * 4.

Using a transfer learning approach, prediction of micromechanical response
using GNNs is possible even in scenarios when little training data is
available (i.e., experimentally measured microstructures).

GNN surrogate modeling, beyond predicting grain-scale elastic response,
provides a framework for rapidly predicting more complex processes such as
plasticity which can in turn be embedded into larger-scale simulation for
prediction of complex material behaviors such as fracture and fatigue.

##  6 Acknowledgements

We would like to thank Dr. Mark Obstalecki and Dr. Paul Shade at the Air Force
Research Laboratory for providing the LSHR HEDM data used for GNN model
testing. ARB is supported in part by ARO Award MURI, NSF CAREER Award
IIS-2045555, and NSF DMS-EPSRC Award 2146079.

##  7 Data Availability

Data will be made available by the authors upon request.

##  8 Software Availability

The GNN model fitting code will be made available upon request. FEPX is
available at https://fepx.info/. Neper is available at https://neper.info/.
HEXRD is available at https://github.com/hexrd.

## References

  * [1] J. D. Eshelby, The determination of the elastic field of an ellipsoidal inclusion, and related problems, Proceedings of the royal society of London. Series A. Mathematical and physical sciences 241 (1226) (1957) 376–396.
  * [2] E. Marin, P. Dawson, On modelling the elasto-viscoplastic response of metals using polycrystal plasticity, Computer Methods in Applied Mechanics and Engineering 165 (1-4) (1998) 1–21.
  * [3] R. A. Lebensohn, A. K. Kanjarla, P. Eisenlohr, An elasto-viscoplastic formulation based on fast fourier transforms for the prediction of micromechanical fields in polycrystalline materials, International Journal of Plasticity 32 (2012) 59–69.
  * [4] W. L. Hamilton, Graph representation learning, Synthesis Lectures on Artifical Intelligence and Machine Learning 14 (3) (2020) 1–159.
  * [5] Z. Wu, S. Pan, F. Chen, G. Long, C. Zhang, S. Y. Philip, A comprehensive survey on graph neural networks, IEEE transactions on neural networks and learning systems 32 (1) (2020) 4–24.
  * [6] J. Zhou, G. Cui, S. Hu, Z. Zhang, C. Yang, Z. Liu, L. Wang, C. Li, M. Sun, Graph neural networks: A review of methods and applications, AI Open 1 (2020) 57–81.
  * [7] J. Fuhg, L. van Wees, M. Obstalecki, P. Shade, N. Bouklas, M. Kasemer, Machine-learning convex and texture-dependent macroscopic yield from crystal plasticity simulations, Materialia (2022) 101446.
  * [8] A. Frankel, K. Tachida, R. Jones, Prediction of the evolution of the stress field of polycrystals undergoing elastic-plastic deformation with a hybrid neural network model, Machine Learning: Science and Technology 1 (3) (2020) 035005\.
  * [9] J. R. Mianroodi, N. H Siboni, D. Raabe, Teaching solid mechanics to artificial intelligence—a fast solver for heterogeneous materials, npj Computational Materials 7 (1) (2021) 1–10.
  * [10] A. Pandey, R. Pokharel, Machine learning based surrogate modeling approach for mapping crystal deformation in three dimensions, Scripta Materialia 193 (2021) 1–5.
  * [11] W. Hamilton, Z. Ying, J. Leskovec, Inductive representation learning on large graphs, Advances in neural information processing systems 30 (2017).
  * [12] M. Dai, M. F. Demirel, Y. Liang, J.-M. Hu, Graph neural networks for an accurate and interpretable prediction of the properties of polycrystalline materials, npj Computational Materials 7 (1) (2021) 1–9.
  * [13] K. Weiss, T. M. Khoshgoftaar, D. Wang, A survey of transfer learning, Journal of Big data 3 (1) (2016) 1–40.
  * [14] R. Quey, P. Dawson, F. Barbe, Large-scale 3d random polycrystals for the finite element method: Generation, meshing and remeshing, Computer Methods in Applied Mechanics and Engineering 200 (17-20) (2011) 1729–1745.
  * [15] R. Quey, L. Renversade, Optimal polyhedral description of 3d polycrystals: Method and application to statistical and synchrotron X-ray diffraction data, Computer Methods in Applied Mechanics and Engineering 330 (2018) 308–333.
  * [16] Neper: Polycrystal Generation and Meshing.
URL https://neper.info/

  * [17] C. Geuzaine, J.-F. Remacle, Gmsh: A 3-d finite element mesh generator with built-in pre-and post-processing facilities, International journal for numerical methods in engineering 79 (11) (2009) 1309–1331.
  * [18] P. Dawson, D. Boyce, FEpX – finite element polycrystals: Theory, finite element formulation, numerical implementation and illustrative examples (2015).
  * [19] FEPX: Finite Element Polycrystal Plasticity.
URL https://fepx.info/

  * [20] R. Asaro, A. Needleman, Texture development and strain hardening in rate dependent polycrystals, Acta Metallurgica 33 (6) (1985) 923–953.
  * [21] E. Marin, P. Dawson, Elastoplastic finite element analyses of metal deformations using polycrystal constitutive models, Computer Methods in Applied Mechanics and Engineering 165 (1-4) (1998) 23–41.
  * [22] F. Roters, P. Eisenlohr, L. Hantcherli, D. Tjahjanto, T. Bieler, D. Raabe, Overview of constitutive laws, kinematics, homogenization and multiscale methods in crystal plasticity finite-element modeling: Theory, experiments, applications, Acta Materialia 58 (4) (2010) 1152–1211.
  * [23] T. Turner, P. Shade, M. Groeber, M. Miller, M. Uchic, Two integrated experimental and modeling approaches to study strain distributions in nickel and nickel-base superalloy polycrystals, Models for Processing and Properties, The Minerals, Metals and Materials Society, Seven Springs, PA (2012) 643–652.
  * [24] E. S. Fisher, C. J. Renken, Single-crystal elastic moduli and the hcp → bcc transformation in Ti, Zr, and Hf, Physical Review 135 (2A) (1964).
  * [25] U. F. Kocks, C. N. Tomé, H.-R. Wenk, Texture and anisotropy: preferred orientations in polycrystals and their effect on materials properties, Cambridge university press, 1998.
  * [26] H. F. Poulsen, Three-dimensional X-ray diffraction microscopy: mapping polycrystals and their dynamics, Vol. 205, Springer Science & Business Media, 2004.
  * [27] J. V. Bernier, N. R. Barton, U. Lienert, M. P. Miller, Far-field high-energy diffraction microscopy: a tool for intergranular orientation and strain analysis, The Journal of Strain Analysis for Engineering Design 46 (7) (2011) 527–547.
  * [28] K. E. Nygren, D. C. Pagan, J. V. Bernier, M. P. Miller, An algorithm for resolving intragranular orientation fields using coupled far-field and near-field high energy x-ray diffraction microscopy, Materials Characterization 165 (2020) 110366.
  * [29] W. D. Musinski, P. A. Shade, D. C. Pagan, J. V. Bernier, Statistical aspects of grain-level strain evolution and reorientation during the heating and elastic-plastic loading of a ni-base superalloy at elevated temperature, Materialia 16 (2021) 101063.
  * [30] D. C. Pagan, K. E. Nygren, M. P. Miller, Analysis of a three-dimensional slip field in a hexagonal ti alloy from in-situ high-energy x-ray diffraction microscopy data, Acta Materialia 221 (2021) 117372.
  * [31] D. C. Pagan, M. P. Miller, Connecting heterogeneous single slip to diffraction peak evolution in high-energy monochromatic x-ray experiments, Journal of applied crystallography 47 (3) (2014) 887–898.
  * [32] M. P. Miller, D. C. Pagan, A. J. Beaudoin, K. E. Nygren, D. J. Shadle, Understanding micromechanical material behavior using synchrotron x-rays and in situ loading, Metallurgical and Materials Transactions A 51 (9) (2020) 4360–4376.
  * [33] R. C. Hurley, E. B. Herbold, D. C. Pagan, Characterization of the crystal structure, kinematics, stresses and rotations in angular granular quartz during compaction, Journal of Applied Crystallography 51 (4) (2018) 1021–1034.
  * [34] F. Monti, D. Boscaini, J. Masci, E. Rodola, J. Svoboda, M. M. Bronstein, Geometric deep learning on graphs and manifolds using mixture model cnns, in: Proceedings of the IEEE conference on computer vision and pattern recognition, 2017, pp. 5115–5124.
  * [35] M. Fey, J. E. Lenssen, Fast graph representation learning with pytorch geometric, arXiv preprint arXiv:1903.02428 (2019).
  * [36] Y. Li, R. Yu, C. Shahabi, Y. Liu, Diffusion convolutional recurrent neural network: Data-driven traffic forecasting, in: International Conference on Learning Representations, 2018.
  * [37] U. Kocks, A. Argon, M. Ashby, Thermodynamics and kinetics of slip, Progress in materials science 19 (1975) 1–281.
  * [38] M. Zecevic, I. J. Beyerlein, M. Knezevic, Coupling elasto-plastic self-consistent crystal plasticity and implicit finite elements: Applications to compression, cyclic tension-compression, and bending to large strains, International Journal of Plasticity 93 (2017) 187–211.
  * [39] A. Beaudoin, P. Dawson, K. Mathur, U. Kocks, A hybrid finite element formulation for polycrystal plasticity with consideration of macrostructural and microstructural linking, International Journal of Plasticity 11 (5) (1995) 501–521.
  * [40] P. Dawson, S. MacEwen, P. Wu, Advances in sheet metal forming analyses: dealing with mechanical anisotropy from crystallographic texture, International Materials Reviews 48 (2) (2003) 86–122.
  * [41] G. B. Sarma, P. R. Dawson, Texture predictions using a polycrystal plasticity model incorporating neighbor interactions, International journal of plasticity 12 (8) (1996) 1023–1054.
  * [42] N. R. Barton, J. V. Bernier, R. A. Lebensohn, D. E. Boyce, The use of discrete harmonics in direct multi-scale embedding of polycrystal plasticity, Computer Methods in Applied Mechanics and Engineering 283 (2015) 224–242.
  * [43] M. Zecevic, M. V. Upadhyay, E. Polatidis, T. Panzner, H. Van Swygenhoven, M. Knezevic, A crystallographic extension to the olson-cohen model for predicting strain path dependence of martensitic transformation, Acta Materialia 166 (2019) 386–401.