- Email: [email protected]
Extracting single-crystal elastic constants from polycrystalline samples using spherical nanoindentation and orientation measurements
Available online at www.sciencedirect.com
ScienceDirect Acta Materialia 79 (2014) 108–116 www.elsevier.com/locate/actamat
Extracting single-crystal elastic constants from polycrystalline samples using spherical nanoindentation and orientation measurements Dipen K. Patel a, Hamad F. Al-Harbi a,b, Surya R. Kalidindi a,⇑ a
George W. Woodruff School of Mechanical Engineering, Georgia Institute of Technology, Atlanta, GA, USA b Mechanical Engineering Department, King Saud University, PO Box 800, Riyadh 11421, Saudi Arabia Received 12 March 2014; received in revised form 8 July 2014; accepted 8 July 2014
Abstract This paper describes a new approach for the extraction of single-crystal elastic stiffness parameters from polycrystalline samples using spherical nanoindentation and orientation measurements combined with finite-element (FE) simulations. The first task of this new approach involves capturing efficiently the functional dependence of the indentation modulus on the lattice orientation at the indentation site and the unknown single-crystal elastic constants. This step is accomplished by probing the function of interest using a suitably constructed FE model of spherical indentation, and establishing a compact spectral representation of the desired function using the discrete values obtained from the simulations. Note that this function needs to be established only once for a selected crystal lattice symmetry. In the second step of the approach presented here, the unknown single-crystal elastic constants for a selected phase are estimated through a regression technique that provides the best match between spherical nanoindentation measurements obtained on differently oriented grains of that phase in a polycrystalline sample (measured by orientation imaging) and the function established in the first step. The accuracy and viability of the proposed approach are demonstrated for an as-cast cubic polycrystalline Fe–3% Si sample. Ó 2014 Acta Materialia Inc. Published by Elsevier Ltd. All rights reserved.
Keywords: Spherical nanoindentation; Orientation imaging; Single-crystal elastic constants; Spectral representations; Finite-element models
1. Introduction Development of robust, physics-based, multiscale materials models is significantly hampered by the lack of validated tools and protocols for characterizing reliably the local (anisotropic) properties at length scales at or below the micron scale. Although numerical techniques such as the finite-element (FE) method have been shown to be successful in simulating complex interactions between microscale constituents of a composite material system [1–14], their predictive capabilities are strongly affected by assumptions made about the constitutive laws used to
⇑ Corresponding author.
E-mail address: [email protected] (S.R. Kalidindi).
describe the local response of the microscale constituents present in these systems. It is often very expensive, and sometimes impossible, to produce sufficiently large volumes of the microscale constituents of interest in their pure form to allow the application of traditional mechanical testing methods (e.g. compression or tensile testing). One approach explored in the literature involves the fabrication of micropillars [15– 17] using a focused ion-beam and testing these pillars in a scanning electron microscope. However, this approach requires access to highly sophisticated equipment and is not particularly well suited for extracting the elastic properties of the microscale constituents in composite material systems. In this paper, we present a new approach for estimating single-crystal elastic properties from polycrystalline
http://dx.doi.org/10.1016/j.actamat.2014.07.021 1359-6454/Ó 2014 Acta Materialia Inc. Published by Elsevier Ltd. All rights reserved.
D.K. Patel et al. / Acta Materialia 79 (2014) 108–116
samples using indentation methods, orientation measurements (by electron backscattered diffraction) and FE models. The approach presented here is formulated to solve the following inverse problem: find the values of the fundamental single-crystal elastic stiffness parameters (for a selected material phase) that are most consistent with a set of in-grain indentation measurements obtained on a range of lattice orientations in a polycrystalline sample. This work builds on prior work [18–22] from our research group. As noted above, the challenge posed in the present work is essentially an inverse problem. In the forward direction, it is relatively easy to build a FE model [18,23,24] that predicts the indentation modulus for a selected crystal lattice orientation and a selected combination of single-crystal elastic stiffness parameters. Such FE models are not ideally suited for addressing the inverse problem stated above. In this paper, we introduce a new approach for addressing this inverse problem that is built on the compact Fourier representations used extensively in our prior work to develop a microstructure-sensitive design framework [25–27]. This new approach is developed and presented first in a general framework that is applicable to any crystal lattice symmetry. Further details of the approach and its viability are then demonstrated with a specific case study involving previously published indentation measurements on a polycrystalline sample of Fe–3% Si [20]. 2. Spherical nanoindentation theory Spherical nanoindentation data analyses procedures are largely based on Hertz’s theory [28–30], which assumes frictionless, elastic contact between two isotropic quadratic surfaces. The main result of this theory can be expressed as: 1 4 3 P ¼ Eeff R2eff h2e ; 3
ð1Þ
where P is the indentation load at the elastic penetration depth, he . Reff and Eeff denote the effective radius and the effective indentation modulus of the sample and the indenter system, defined as: 1 1 t2s 1 t2i ¼ þ ; Eeff Es Ei
1 1 1 ¼ þ : Reff Rs Ri
ð2Þ
In Eq. (2), E and m denote the Young’s modulus and Poisson’s ratio of the indenter (subscript i) and the specimen (subscript s), and R denotes the radius. In the purely elastic indentation of a perfectly flat surface, Reff ¼ Ri . Although the theory for elastic indentation of isotropic materials is well established, it is not directly applicable to most in-grain nanoindentation studies on polycrystalline materials since the indentation zone size in such measurements is typically much smaller than the grain size in the sample. At this length scale, most crystalline materials exhibit significant anisotropy in their elastic response. A number of different approaches have been explored in the literature [31,32] to take the elastic anisotropy into account
109
in the analyses of the indentation measurements. Vlassak, Nix and co-workers [33–35] have developed a rigorous analytical framework based on Hertzian theory to address the elastic indentation of anisotropic samples. Their theory indicates that the inclusion of a crystal lattice orientationdependent parameter, b, into the definition of the effective indentation modulus will adequately capture the anisotropic elastic indentation response of cubic crystals for any arbitrary orientation of the crystal lattice in the indentation zone. More specifically, they suggest: 1 1 t2s 1 t2i ¼b þ ; ð3Þ Eeff Es Ei where Es and ms denote the effective Young’s modulus and Poisson’s ratio, respectively [34,35], for a randomly textured polycrystalline sample. In the present study, our goal is first to establish the main underlying features of the Vlassak–Nix theory [34,35] that are central to the indentation analyses protocols discussed in this paper. More specifically, Eq. (3) combined with Eq. (1) implies that the indentation load, P , 3 continues to be directly proportional to h2e , even when the sample exhibits an anisotropic elastic response. If, indeed, this relationship holds true for any arbitrary orientation of the crystal, then it allows us to extract an orientationdependent indentation modulus (the parameter b is expected to be orientation dependent) from the indentation measurement on any crystal in the sample. Note that this would be an extremely important step in addressing the goals of the present work. In order to establish whether or not an indentation modulus can be defined for any crystal lattice orientation (from 3 the expected linear relationship between P and h2e ), we developed and employed a FE model. Note that the investigation described above cannot be conducted directly with measurements as easily, because of the difficulties in isolating the elastic portions of the indentation load–displacement curves [19]. Furthermore, numerical approaches such as FE models circumvent many of the difficulties and uncertainties faced in experimental measurements (e.g. identification of the initial point of contact, accuracy in the description of the geometry of indenter, unavoidable plasticity, friction at the contact surface), and are ideally suited for critical validations of fundamental concepts [18]. 3. Finite-element model of spherical nanoindentation The 3-D FE model to simulate elastic spherical indentation of anisotropic crystals was produced using the commercial FE code ABAQUS [36]. The FE model developed for this study is comprised of two 3-D bodies: (i) an elastically deformable sample with an initially flat surface discretized into 57,000 eight-noded, 3-D, continuum (C3D8) elements, and (ii) a rigid hemispherical indenter of radius 13.5 lm (i.e. the same size indenter as in prior experimental studies [19,20,22]). The size of the sample was selected as 18 lm 18 lm 9 lm to ensure that it is much
110
D.K. Patel et al. / Acta Materialia 79 (2014) 108–116
larger than the typical indentation zone size (observed to be about 1.2 lm from our simulations). The sample has been discretized into seven regions as shown in Fig. 1a to permit the use of a progressively higher mesh density as we approach the indentation zone directly below the indenter–sample contact surface. Discretizing the mesh in this manner allowed us to capture the stress and strain fields with sufficient accuracy while keeping the total number of elements in the FE model (as well as the computational cost) relatively low. Significant computational resources were required and utilized to keep the runtime to just under 30 min per simulation. The FE simulations described in this work were executed on two nodes (32 cores per node) provided by Trestles (part of the XSEDE supercomputing facility). A hard surface-to-surface, frictionless, contact definition was used to model the contact behavior between the indenter (master) surface and the elastically deformable sample (slave) surface. Displacement boundary conditions in the z-direction were imposed on the reference node (placed at the center of the indenter), which is tied rigidly to the entire surface of the indenter. The top surface of the sample was free to move while the bottom surface of the sample was constrained along the z-direction (indentation direction). The inputs to the FE model are the single-crystal elastic stiffness constants of the material of interest in the crystal frame (i.e. C11, C12 and C44) and the crystal lattice orientation. The lattice orientation of a single crystal with respect to the sample reference frame (defined in Euclidean 3-D space) can be described by a set of three well-defined sequential rotations. The three ordered angles ðu1 ; U; u2 Þ associated with these rotations are commonly known as Bunge–Euler angles [37]. Since the only relevant crystallographic direction in the spherical indentation experiment is the one parallel to the indentation direction (i.e. the indentation response is independent of any in-plane rotation of the sample), the crystal lattice orientation with respect to the indentation direction can be described using only two of the three Bunge–Euler angles [20,37]. We shall use g ¼ ðU; u2 Þ to denote the lattice orientation in the indentation experiments and models discussed in this work.
4. New approach to estimation of single-crystal elastic constants In this section, we present the proposed approach as a general framework. In the next section, we demonstrate a specific application and provide many more details of the computations involved. 4.1. Spectral representation of indentation modulus As described earlier, the first step of the proposed protocol involves establishing the functional dependence of the effective indentation modulus on the single-crystal elastic constants and the crystal lattice orientation. In order to establish this functional dependence, a large number of data points can be accumulated from the FE model
(b) Load, P (mN)
(a)
Finite displacement in the z-direction was imposed on the indenter and the corresponding total reaction force exerted by the slave surface on the indenter is predicted and used to determine the effective indentation modulus of the sample. The FE model developed in this study was validated by comparing the simulated load–displacement curves against the analytical prediction from Hertz’s theory for purely isotropic elastic deformation (see Fig. 1b). Additionally, a detailed sensitivity study was conducted to ensure that the details of the mesh and the details of the far-field boundary conditions did not influence significantly the predicted load–displacement response. The FE model developed in this study was validated by comparing the predicted indentation modulus against the analytical expression reported by Vlassak and Nix [35] in Table 1. It is clear that the FE predictions are in excellent agreement with the theoretically expected values. Furthermore, the FE models predicted that a linear relationship 3 between P and h2e holds for the anisotropic elastic indentations, in full accord with the theories presented by Vlassak and Nix [34,35]. Some examples of the FE-predicted load– displacement relations are presented in Fig. 2. This is highly encouraging and provides a path forward for addressing the inverse problem described earlier.
Displacement, he (nm) Fig. 1. (a) The 3-D FE mesh used to simulate elastic anisotropic response under the indenter (isometric view of the region under the indenter tip). (b) Load–displacement response of FE simulation vs. Hertz theory for the isotropic response case.
D.K. Patel et al. / Acta Materialia 79 (2014) 108–116 Table 1 Indentation moduli values for the cubic crystals predicted from the FE model and theoretical values reported in Vlassak and Nix [35] for symmetric orientations. Material and corresponding elastic constants (GPa)
Orientation direction (h, k, l)
Theoretical Eind (GPa)
FE prediction Eind (GPa)
Aluminiumb C11 = 107.3 C12 = 60.90 C44 = 28.3
(1 0 0) (1 1 0) (1 1 1)
79 80 81
78.96 80.93 80.25
Copperb C11 = 170.2 C12 = 114.9 C44 = 61.0
(1 0 0) (1 1 0) (1 1 1)
129 138 141
129.65 137.19 141.05
b-Brass (CsCl)a C11 = 126.5 C12 = 107.7 C44 = 80.3
(1 0 0) (1 1 0) (1 1 1)
95 112 117
93.87 109.40 115.47
a b
Taken from Vlassak and Nix [35]. Taken from Simmon and Wang [38].
can be interpreted as results of probing single crystals of a selected phase through indentations in different crystallographic directions. The functional dependencies of the indentation modulus described above can therefore be succinctly expressed as Eind ðg; C 11 ; C 12 ; C 44 Þ. Our first goal is to capture this functional dependence of the indentation modulus using highly efficient spectral representations. In this paper, we have employed the spherical surface harmonics (SSH) [37] denoted as K ml ðgÞ to serve as a Fourier basis for compact representation of this function over the orientation space of interest. The advantages and compactness of using these basis functions on the orientation space have already been discussed in our prior publications [39,40]. For the variables denoting the elastic stiffness constants, we explore here the use of Legendre polynomials, P n ðÞ [41] as a Fourier basis. Note that the Legendre polynomials form a Fourier basis over the range [1, +1], and therefore the stiffness parameters need to be suitably rescaled. Each rescaled elastic stiffness constant can be defined as: ~ ¼ 2C C min C max ; C C max C min
ð4Þ
where C min and C max are the respective minimum and maximum values for the specific elastic stiffness constant. The ~ 11 ; C ~ 12 ; C ~ 44 Þ is then expressed function of interest, Eind ðg; C as ~ 11 ; C ~ 12 ; C ~ 44 Þ :¼ Eind ðg; C
MðlÞ X 1 X 1 X
~ 11 ÞP r ðC ~ 12 ÞP s ðC ~ 44 Þ: Amqrs K ml ðgÞP q ðC l
ð5Þ
l¼0 m¼1 q;r;s¼0
In Eq. (5), Amqrs represent the Fourier coefficients. The l function MðlÞ denotes the number of terms needed in the enumeration of index, m; these numbers are a function of the index l and are determined by crystal symmetries [37].
Load (mN)
predictions. Each FE simulation provides one data point as an estimate of the indentation modulus for a specific combination of one crystal orientation and one set of values for the elastic constants. A least-squares fit between the FE-predicted values of P 3 and h2e in each simulation yields an estimate of the indentation modulus (Eind ) (see Eq. (1)). Implicitly, the indentation modulus so obtained depends on the crystal lattice orientation (applied here uniformly to all elements in the sample) and the values of the single-crystal elastic constants used for the sample material (e.g. C11, C12 and C44 for cubic crystals) [38]. Therefore, the different measurements on the differently oriented grains in the polycrystalline sample
111
he3/2 (nm) Fig. 2. P vs. (he)(3/2) plots obtained for copper for different orientations spread over a cubic-transversely isotropic FZ (highlighted area in the subplot).
112
D.K. Patel et al. / Acta Materialia 79 (2014) 108–116
The most common approaches used in the literature to establish the Fourier coefficients generally exploit the orthonormal properties of the Fourier basis [25,42–44]. However, such approaches require evaluation of the function at an extremely large number of locations in the compounded space of the complete ranges of all the independent variables involved. Since the function evaluation here is being attempted through a computationally expensive FE model, the traditional approach of establishing Fourier coefficients is not very practical for the present problem. Consequently, we explore in this work an approach for establishing the values of the Fourier coefficients in Eq. (5) using ordinary least-squares (OLS) regression analysis. In order to use a regression technique, it is necessary to truncate the Fourier series in Eq. (5) to a finite number of terms in the expansion. Simplifying the notation a little, Eq. (5) can be expressed as: ~ 11 ; C ~ 12 ; C ~ 44 Þ Eind ðg; C
~ ~q L X X
~ ~ ~ Aqrs L K L ðgÞP q ðC 11 ÞP r ðC 12 ÞP s ðC 44 Þ:
dataset. Over-fitting needs to be avoided as it will likely produce unreliable estimates in subsequent application of the Fourier function (especially for new data points not included in the calibration dataset). In order to track the potential for over-fitting of the Fourier coefficients obtained at each level of truncation ~ and ~q), leave-one-out-cross-validation (i.e. selection of L (LOOCV) can be employed [47]. LOOCV allows an objective selection of the truncation levels by establishing the fit N times (this is size of the dataset), while leaving one data point out of the regression each time. In doing so, LOOCV will quantify the contribution of each data point to the Fourier coefficients established in each regression. Given a large N, for an over-fitted Fourier expansion, the exclusion of a single data point will cause significant change in the coefficients, whereas for a good fit this change will be negligible. In order to quantify the robustness of the fit, the following objective measures could be utilized:
ð6Þ
I. Mean absolute error of the fit, e, defined as:
L¼0 q;r;s¼0
where L enumerates each distinct combination of ðm; lÞ in ~ and ~ the SSH series expansion, and L q denote the truncation limits in the SSH and the Legendre polynomial expansions. Note also that a single limit was selected in Eq. (6) for all three Legendre polynomial bases used in this equation. This was merely done for simplicity here, and different limits can be selected for the different variables as needed in any specific application. As described earlier, each FE simulation of the indentation produces one discrete data point for a combination of one crystal orientation and one set of values for the cubic elastic stiffness constants. Consider a dataset with N such data points collected from multiple runs of the FE simulations for different combinations of crystal orientations and cubic elastic stiffness constants. n Such a FE-generated FEðnÞ ~ ðnÞ ; C ~ ðnÞ ; C ~ ðnÞ ; dataset can be denoted as Eind ; gðnÞ ; C 11 12 44 o n ¼ 1; 2; . . . ; N . Regression analysis is then formulated as a minimization of the difference (sum of the squares of the residual) between the FE-generated dataset and the truncated Fourier representation shown in Eq. (6). In other words, OLS [45,46] aims to solve the following minimization problem: ! ~ q ~ N L X FEðnÞ X X qrs ðnÞ ðnÞ ðnÞ ðnÞ ~ ~ ~ min AL K L g P q C 11 P r C 12 P s C 44 : Eind qrs AL n¼1 L¼0 q;r;s¼0
N 1X jeðnÞ j N n¼1 ~ q ~ N L 1X FEðnÞ X X qrs ðnÞ ðnÞ ðnÞ ðnÞ ~ ~ ~ AL K L g P q C 11 P r C 12 P s C 44 : ¼ E N n¼1 ind L¼0 q;r;s¼0
e ¼
ð8Þ
Median absolute deviation (MAD) of error of the fit as: MADe ¼ Median eðnÞ Median eð1Þ ; eð2Þ ; . . . ; eðN Þ ; n ¼ f1; 2; . . . ; N g:
ð9Þ
II. Mean absolute error of LOOCV, eCV , and MAD of LOOCV, MADCV , defined analogously as the mean and median of the error for the test data point in each of the N repeated fits. The two measures of error of fit defined above will show ~ and ~q, whereas improvement of fit with higher values of L the two measures of error of LOOCV are expected to show ~ and ~ a decline in robustness of fit with higher values of L q accounting for over-fit of data. Therefore, a compromise is often made in choosing the best fit based on the values of all four measures defined above. In the next step, Eq. (6) can be used to address the inverse problem of estimating the unknown single-crystal elastic stiffness constants for a selected phase in a polycrystalline sample.
ð7Þ
Note that the values of the Fourier coefficients established using the regression method described above can be sensitive to the truncation levels used, i.e. the values ~ and ~ of L q. Since we do not generally know a priori the ~ and ~ right values of L q in any specific application, this approach needs a few repeated trials. Although higher val~ and ~ ues of L q will generally produce a lower error (difference between the FE dataset and the fitted Fourier expansion), they may also result in over-fitting of the
4.2. Protocol to extract single-crystal elastic constants As described earlier, developing and validating an inverse solution methodology to extract grain-scale elastic properties for a selected phase in a polycrystalline material is the main objective of this work. The protocol described here aims to match the measurements of a set of indentation moduli on differently oriented grains (typically measured using orientation imaging microscopy [21]) of a
D.K. Patel et al. / Acta Materialia 79 (2014) 108–116
selected phase whose independent elastic stiffness constants are unknown (denoted as C ) with the function established in the previous step (Eq. (6)). For cubic materials, for sim~ , C ~ , C ~ ) correspond to rescaled plicity of notation, let (C 1 2 3 ~ , C ~ , C ~ ) using Eq. (4). unknown elastic constants (C 11 12 44 ðjÞ
Let g ; j ¼ 1; 2; . . . ; J denote the specific orientations ðjÞ
where the indentation measurements were made, and Eind denote the corresponding measurements of the indentation moduli. Let the corresponding theoretical value for each of the measurement, predicted from the Fourier representa ðjÞ ~ ; C ~ ; C ~ . The task tion (Eq. (6)), be denoted as, Eind gðjÞ ; C 1 2 3 of estimating the unknown single-crystal elastic constants then reduces to minimizing the difference between the val ðjÞ ~ ; C ~ ; C ~ and EðjÞ for all of the available ues of Eind gðjÞ ; C ind 1 2 3 measurements on the selected material phase. The corresponding minimization problem can be expressed as: J h X i2 ðjÞ ðjÞ ~ ;C ~ ;C ~ ¼ min ~ ; C ~ ; C ~ : min f C Eind Eind gðjÞ ; C 1 2 3 1 2 3
~ ;C ~ ;C ~ C 1 2 3
~ ;C ~ ;C ~ C 1 2 3 j¼1
ð10Þ
The above minimization problem can be addressed by evaluating and equating the relevant derivatives of the function to zero: ~ ;C ~ ;C ~ ¼ @f Fi C 1 2 3 ~ i @C J h i ðjÞ ðjÞ ~ ~ ~ X ðjÞ ðjÞ ~ ~ ~ @Eind g ; C 1 ; C 2 ; C 3 ¼ 2 EðjÞ ¼ 0: ind E ind g ; C 1 ; C 2 ; C 3 ~ i @C j¼1
ð11Þ
ðkÞ
1
~ i ½J ij F j ; ¼C @F i J ij ¼ : ~ j ðkÞ @C ~ ¼C ~ C Ci
n
ð12Þ
n
where the superscript (k) denotes the iteration number. 5. Case study: as-cast polycrystalline Fe–3% Si The new protocols presented in this paper to extract the single-crystal elastic stiffness parameters have been applied to an as-cast polycrystalline sample of Fe–3% Si. As described in the previous section, the first step involves ~ 11 ; C ~ 12 ; C ~ 44 in a spectral establishing the function Eind g; C form (see Eq. (6)) using regression methods on datasets assembled from results produced by the FE models described in Section 3. As described earlier, each FE simulation of the spherical indentation for a selected combination of crystal orientation and cubic elastic stiffness parameters yields one data point. In the present example, the desired dataset was accumulated by executing a total of 2700 FE simulations that employed 300 distinct sets of the cubic elastic stiffness constants [38] compounded with nine distinct crystal orientations. The elastic stiffness constants were selected in the ranges of 50 6 C 11 6 250, 40 6 C 12 6 150 and 15 6 C 44 6 125 (see Fig. 3), which were chosen so that they cover the typical ranges for most structural cubic metals and alloys of interest. This set covers a cubic anisotropy ratio, defined as A ¼ 2C 44 = ðC 11 C 12 Þ, in the range ð0 < A < 8Þ. The crystal orientations were selected such that they cover the relevant fundamental zone (FZ) of orientations of interest for the present application [25] (see also the inset of Fig. 2): " FZ ¼ ðU; u2 Þj cos1
cos u2 pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 1 þ cos2 u2
!
p 6 U2 6 ; 2
06u6
# p : 4
ð13Þ
C11 (GPa)
The system of nonlinear equations in Eq. (11) can be solved using various numerical approaches. In this paper, we have employed an iterative Newton–Raphson algorithm to solve this system of nonlinear equations. This approach can be summarized as [46]:
ðkþ1Þ
113
C44 (GPa) C12 (GPa) Fig. 3. 300 Distinct sets of independent elastic stiffness constants ðC 11 ; C 12 ; C 44 Þ used in the present study to establish the spectral representations.
D.K. Patel et al. / Acta Materialia 79 (2014) 108–116
The symmetrized spherical surface harmonics basis functions appropriate for this case study are those exhibit ing cubic-transversely isotropic symmetry denoted as K_ ml ðgÞ [37]. In this description of symmetry, the first term reflects the crystal symmetry and the second term reflects the sample symmetry (arises because the indentation direction is unaffected by any in-plane rotation of the sample). As described in the previous section, the Fourier coefficients of Eq. (6) can be determined using OLS regression (Eq. (7)). The values of e (Eq. (8)) and MADe (Eq. (9)) were com~ puted for the training set for different selections of ~ q and L, and are summarized in Fig. 4. The mean and median mea~ sures of the error decrease with increasing values of ~ q and L as expected (see Fig. 4a). The LOOCV analyses (Fig. 4b) indicates over-fitting of the data points at high values of ~ For instance, the eCV values for the fits both ~q and L. ~ > 12 and ~ obtained using L q > 3 start to increase significantly. Examination of the error measures shown in ~ ¼ 12 and ~ Fig. 4a and b suggests that L q ¼ 3 provide a good fit for the dataset acquired in this study using FE simulations. This level of truncation produces a total of 348 Fourier coefficients; these are utilized in this study to describe the functional dependence of indentation modulus on orientations and the independent elastic constants. The accuracy of the truncated spectral representation of the indentation moduli over the selected range of cubic elastic constants for several orientations was compared against the values predicted by the previously validated FE model and the analytical solution provided by Vlassak and Nix [35]. Representative results are shown in Fig. 5 for a selected range of elastic anisotropic ratios for two orientations that provide the extremal responses. The indentation modulus in this plot was normalized by the indentation modulus of a polycrystalline aggregate, following the approach outlined by Vlassak and Nix [34] and using the expressions provided by Hashin and Shtrikman [48]. The error between the FE prediction and the spectral representation for any given set of cubic elastic constants (within the selected bounds mentioned earlier) and any orientation was found to be within 4%. It is clear that the spectral representation developed in this work captures well the variations in the indentation moduli for different orientations and different degrees of elastic anisotropy in the sample. We now turn our attention to estimating the singlecrystal elastic stiffness parameters of Fe–3% Si, for which indentation measurements on differently oriented single crystals were already reported in prior work [20]. These results are summarized in Table 2, where the orientations are also plotted on the standard inverse pole figure. The effective indentation moduli reported in Table 2 for Fe– 3% Si were computed using Eqs. 1–3, where the Young’s modulus and Poisson ratio of the indenter were taken as 1000 GPa and 0.07, respectively [49]. Estimating the single-crystal constants of Fe–3% Si then reduces to minimizing the difference between the values predicted by the Fourier function established earlier and the discrete
Fig. 4. (a) The mean absolute error by varying numbers of SSH basis ~ and degree of Legendre polynomial ð~qÞ with MAD of functions ðLÞ absolute error as the error bars. (b) The mean absolute error of LOOCV ~ and degree of Legendre by varying numbers of SSH basis functions ðLÞ polynomial ð~qÞ with MAD of absolute error of LOOCV as the error bars.
1.1
Normalized Indentation Modulus
114
1.05 Spectral Function {100}
1
Spectral Function {111}
0.95
FE Prediction Analytical Function
0.9 0.85 0.8 0.75 1
3 5 Anisotropy Ratio (A)
7
Fig. 5. Comparison of the spectral representation of normalized indentation modulus as a function of anisotropic ratio vs. the predictions of the FE model and the analytical model of Vlassak and Nix [35].
D.K. Patel et al. / Acta Materialia 79 (2014) 108–116
115
Table 2 (a) Measured values of the effective indentation modulus, Eeff, for 11 different orientations, and (b) representation of the orientations in an inverse pole figure map.
(a)
(b)
measured indentation moduli reported in Table 2. This minimization was accomplished using a Newton–Raphson scheme described in the previous section (see Eqs. 10–12). The estimated values of the single-crystal elastic constants for the as-cast Fe–3% Si polycrystalline sample from this study are C 11 ¼ 215:7 GPa, C 12 ¼ 131:6 GPa and C 44 ¼ 122:2 GPa. The single-crystal elastic constants obtained using the inverse solution approach described in this paper are within 5% of the typical values reported in literature for Fe–3% Si ðC 11 ¼ 225:28 GPa; C 12 ¼ 135:1 GPa and C 44 ¼ 123:8 GPaÞ [50–52]. 6. Discussion The accuracy of the spectral framework mentioned above can be improved with additional measurements of indentation moduli on additional crystal orientations. Note that we used only 11 measured indentation moduli in the present case study. One might actually argue that we need only three measurements (on three distinct crystal orientations) to estimate the three unknown elastic stiffness constants in cubic crystals. However, because of the regression methods used in the approach described here, it is important to measure indentation moduli for orientations that are well spread out in the fundamental zone of orientations. For example, if we were to repeat the inverse solution presented here by considering only three orientations clustered in any one of the corners of the inverse pole figure (depicts a certain projection of the fundamental zone of orientations) shown in Table 2, it would produce highly
erroneous estimates of the single-crystal elastic stiffness constants. In fact, it is best to ensure that there at least some measurements in each of the three corners of the inverse pole figure shown in Table 2, as these orientations typically produce the maximum contrast in the indentation moduli. It is also likely that there are better representations of the indentation modulus as a function of the orientation and the elastic stiffness constants. Although we have had a lot of experience with the use of SSH as a Fourier basis for the orientation variable, we only have very limited prior experience in spectral representation of the elastic stiffness parameters. Therefore, it is entirely possible that one might accomplish a more accurate representation of the indentation modulus function (Eq. (5)) with a different Fourier basis, which might in turn improve the accuracy of the estimates of the elastic stiffness constants. 7. Conclusion A new two-step procedure has been developed to estimate the single-crystal elastic stiffness parameters from polycrystalline samples using spherical nanoindentation and orientation measurements combined with FE simulations. The first step involves a one-time capture of the functional dependence of the indentation modulus on the lattice orientation at the indentation site and the unknown single-crystal elastic constants for any presumed crystal symmetry. This step was accomplished in this study for cubic crystals by employing a FE model of spherical
116
D.K. Patel et al. / Acta Materialia 79 (2014) 108–116
indentation. In the second step of the approach presented here, the unknown single-crystal elastic constants for a selected phase are estimated through a regression technique. The accuracy and viability of the proposed approach were demonstrated for an as-cast cubic polycrystalline Fe–3% Si sample. Acknowledgements
[18] [19] [20] [21] [22] [23] [24] [25] [26]
This work was supported by ARO through Grant No. W911NF-10-1-0409. Computations described in this work used the Extreme Science and Engineering Discovery Environment (XSEDE), which is supported by National Science Foundation Grant No. OCI-1053575. The codes for fast computation of GSH coefficients were provided by Material Resources, LLC, Dayton, OH. References [1] Roters F, Eisenlohr P, Bieler TR, Raabe D. Crystal plasticity finite element methods in materials science and engineering. New York: Wiley; 2011. [2] Raabe D, Klose P, Engl B, Imlau K-P, Friedel F, Roters F. Adv Eng Mater 2002;4:169. [3] Haddadi H, Bouvier S, Banu M, Maier C, Teodosiu C. Int J Plast 2006;22:2226. [4] Zhao Z, Roters F, Mao W, Raabe D. Adv Eng Mater 2001;3:984. [5] Anand L, Kalidindi SR. Mech Mater 1994;17:223. [6] Bronkhorst CA, Kalidindi SR, Anand L. Philos Trans R Soc Lond Ser A Math Phys Eng Sci 1992;341:443. [7] Kalidindi SR, Bronkhorst CA, Anand L. J Mech Phys Solids 1992;40:537. [8] Kalidindi SR, Bhattacharyya A, Doherty RD. Proc R Soc Lond Ser A Math Phys Eng Sci 2004;460:1935. [9] Bhattacharyya A, El-Danaf E, Kalidindi SR, Doherty RD. Int J Plast 2001;17:861. [10] Kalidindi SR, Schoenfeld SE. Mater Sci Eng A Struct Mater Prop Microstruct Process 2000;293:120. [11] McDowell DL. Int J Plast 2010;26(9):1280. Special issue in honor of David L. McDowell. [12] Kalidindi SR, Bhattacharyya A, Doherty RD. Proc R Soc A 2004;460:1935. [13] Kalidindi SR, Anand L. J Mech Phys Solids 1994;42:459. [14] Bachu V, Kalidindi SR. Mater Sci Eng A Struct Mater Prop Microstruct Process 1998;257:108. [15] Shan ZW, Mishra RK, Syed Asif SA, Warren OL, Minor AM. Nat Mater 2008;7:115. [16] Uchic MD, Dimiduk DM, Florando JN, Nix WD. Science 2004;305:986. [17] Li N, Mara NA, Wang YQ, Nastasi M, Misra A. Scr Mater 2011;64:974.
[27] [28] [29] [30] [31] [32]
[33] [34] [35] [36] [37] [38] [39] [40] [41] [42] [43] [44] [45] [46] [47] [48] [49] [50] [51] [52]
Donohue BR, Ambrus A, Kalidindi SR. Acta Mater 2012;60:3943. Kalidindi SR, Pathak S. Acta Mater 2008;56:3523. Pathak S, Stojakovic D, Kalidindi SR. Acta Mater 2009;57:3020. Proust GKS. TMS Lett 2004;7:151. Pathak S. Development and validation of a novel data analysis procedure for spherical nanoindentation. Drexel University; 2009. Bhattacharya A, Nix W. Int J Solids Struct 1988;24:881. Suresh S, Giannakopoulos A. Acta Mater 1998;46:5755. Adams BL, Kalidindi SR, Fullwood DT. Microstructure sensitive design for performance optimization. Elsevier Science; 2012. Kalidindi SR, Houskamp JR, Lyons M, Adams BL. Int J Plast 2004;20:1561. Li DS, Garmestani H, Schoenfeld S. Scr Mater 2003;49:867. Hertz, H. On the contact of elastic solids. Hertz’s Miscellaneous Papers, Macmillan; 1986. Johnson KL, Johnson KKL. Contact mechanics. Cambridge: Cambridge University Press; 1987. Willis J. J Mech Phys Solids 1966;14:163. Cammarata R, Schlesinger T, Kim C, Qadri S, Edelstein A. Appl Phys Lett 1990;56:1862. Farthing L, Weihs T, Kisker D, Krajewski J, Tang M, Stevenson D. Hardness and elastic modulus measurements in CdTe and ZnTe thin film and bulk samples and ZnTe–CdTe superlattices. MRS proceedings, vol. 130. Cambridge: Cambridge University Press; 1988. Vlassak J, Ciavarella M, Barber J, Wang X. J Mech Phys Solids 2003;51:1701. Vlassak JJ, Nix W. Philos Mag A 1993;67:1045. Vlassak JJ, Nix W. J Mech Phys Solids 1994;42:1223. ABAQUS. Simulia Providence, RI: Dassault Systemes; 2009. Bunge HJ. Texture analysis in materials science: mathematical methods. Butterworths; 1982. Simmons G, Wang H. Single crystal elastic constants and calculated aggregate properties: a handbook. Cambridge, MA: MIT Press; 1971. Fast T, Knezevic M, Kalidindi SR. Comput Mater Sci 2008;43:374. Fullwood DT, Niezgoda SR, Adams BL, Kalidindi SR. Prog Mater Sci 2010;55:477. Copson ET. An introduction to the theory of functions of a complex variable. Oxford: Clarendon Press; 1950. Kalidindi SR, Duvvuru HK, Knezevic M. Acta Mater 2006;54:1795. Kalidindi SR, Landi G, Fullwood DT. Acta Mater 2008;56:3843. Knezevic M, Al-Harbi HF, Kalidindi SR. Acta Mater 2009;57:1777. Bishop CM. Pattern recognition and machine learning. Berlin: Springer Verlag; 2006. Press WH. Numerical recipes: the art of scientific computing. 3rd ed. Cambridge: Cambridge University Press; 2007. C ¸ ecßen A, Fast T, Kumbur EC, Kalidindi SR. J Power Sources 2014;245:144. Hashin Z, Shtrikman S. J Mech Phys Solids 1962;10:335. Field J, Swain M. J Mater Res (USA) 1993;8:297. Alberts H, Wedepohl P. Physica 1971;53:571. Machova A, Kadecˇkova´ S. Czech J Phys B 1977;27:555. Routbort J, Reid C, Fisher E, Dever D. Acta Metall 1971;19:1307.
ScienceDirect Acta Materialia 79 (2014) 108–116 www.elsevier.com/locate/actamat
Extracting single-crystal elastic constants from polycrystalline samples using spherical nanoindentation and orientation measurements Dipen K. Patel a, Hamad F. Al-Harbi a,b, Surya R. Kalidindi a,⇑ a
George W. Woodruff School of Mechanical Engineering, Georgia Institute of Technology, Atlanta, GA, USA b Mechanical Engineering Department, King Saud University, PO Box 800, Riyadh 11421, Saudi Arabia Received 12 March 2014; received in revised form 8 July 2014; accepted 8 July 2014
Abstract This paper describes a new approach for the extraction of single-crystal elastic stiffness parameters from polycrystalline samples using spherical nanoindentation and orientation measurements combined with finite-element (FE) simulations. The first task of this new approach involves capturing efficiently the functional dependence of the indentation modulus on the lattice orientation at the indentation site and the unknown single-crystal elastic constants. This step is accomplished by probing the function of interest using a suitably constructed FE model of spherical indentation, and establishing a compact spectral representation of the desired function using the discrete values obtained from the simulations. Note that this function needs to be established only once for a selected crystal lattice symmetry. In the second step of the approach presented here, the unknown single-crystal elastic constants for a selected phase are estimated through a regression technique that provides the best match between spherical nanoindentation measurements obtained on differently oriented grains of that phase in a polycrystalline sample (measured by orientation imaging) and the function established in the first step. The accuracy and viability of the proposed approach are demonstrated for an as-cast cubic polycrystalline Fe–3% Si sample. Ó 2014 Acta Materialia Inc. Published by Elsevier Ltd. All rights reserved.
Keywords: Spherical nanoindentation; Orientation imaging; Single-crystal elastic constants; Spectral representations; Finite-element models
1. Introduction Development of robust, physics-based, multiscale materials models is significantly hampered by the lack of validated tools and protocols for characterizing reliably the local (anisotropic) properties at length scales at or below the micron scale. Although numerical techniques such as the finite-element (FE) method have been shown to be successful in simulating complex interactions between microscale constituents of a composite material system [1–14], their predictive capabilities are strongly affected by assumptions made about the constitutive laws used to
⇑ Corresponding author.
E-mail address: [email protected] (S.R. Kalidindi).
describe the local response of the microscale constituents present in these systems. It is often very expensive, and sometimes impossible, to produce sufficiently large volumes of the microscale constituents of interest in their pure form to allow the application of traditional mechanical testing methods (e.g. compression or tensile testing). One approach explored in the literature involves the fabrication of micropillars [15– 17] using a focused ion-beam and testing these pillars in a scanning electron microscope. However, this approach requires access to highly sophisticated equipment and is not particularly well suited for extracting the elastic properties of the microscale constituents in composite material systems. In this paper, we present a new approach for estimating single-crystal elastic properties from polycrystalline
http://dx.doi.org/10.1016/j.actamat.2014.07.021 1359-6454/Ó 2014 Acta Materialia Inc. Published by Elsevier Ltd. All rights reserved.
D.K. Patel et al. / Acta Materialia 79 (2014) 108–116
samples using indentation methods, orientation measurements (by electron backscattered diffraction) and FE models. The approach presented here is formulated to solve the following inverse problem: find the values of the fundamental single-crystal elastic stiffness parameters (for a selected material phase) that are most consistent with a set of in-grain indentation measurements obtained on a range of lattice orientations in a polycrystalline sample. This work builds on prior work [18–22] from our research group. As noted above, the challenge posed in the present work is essentially an inverse problem. In the forward direction, it is relatively easy to build a FE model [18,23,24] that predicts the indentation modulus for a selected crystal lattice orientation and a selected combination of single-crystal elastic stiffness parameters. Such FE models are not ideally suited for addressing the inverse problem stated above. In this paper, we introduce a new approach for addressing this inverse problem that is built on the compact Fourier representations used extensively in our prior work to develop a microstructure-sensitive design framework [25–27]. This new approach is developed and presented first in a general framework that is applicable to any crystal lattice symmetry. Further details of the approach and its viability are then demonstrated with a specific case study involving previously published indentation measurements on a polycrystalline sample of Fe–3% Si [20]. 2. Spherical nanoindentation theory Spherical nanoindentation data analyses procedures are largely based on Hertz’s theory [28–30], which assumes frictionless, elastic contact between two isotropic quadratic surfaces. The main result of this theory can be expressed as: 1 4 3 P ¼ Eeff R2eff h2e ; 3
ð1Þ
where P is the indentation load at the elastic penetration depth, he . Reff and Eeff denote the effective radius and the effective indentation modulus of the sample and the indenter system, defined as: 1 1 t2s 1 t2i ¼ þ ; Eeff Es Ei
1 1 1 ¼ þ : Reff Rs Ri
ð2Þ
In Eq. (2), E and m denote the Young’s modulus and Poisson’s ratio of the indenter (subscript i) and the specimen (subscript s), and R denotes the radius. In the purely elastic indentation of a perfectly flat surface, Reff ¼ Ri . Although the theory for elastic indentation of isotropic materials is well established, it is not directly applicable to most in-grain nanoindentation studies on polycrystalline materials since the indentation zone size in such measurements is typically much smaller than the grain size in the sample. At this length scale, most crystalline materials exhibit significant anisotropy in their elastic response. A number of different approaches have been explored in the literature [31,32] to take the elastic anisotropy into account
109
in the analyses of the indentation measurements. Vlassak, Nix and co-workers [33–35] have developed a rigorous analytical framework based on Hertzian theory to address the elastic indentation of anisotropic samples. Their theory indicates that the inclusion of a crystal lattice orientationdependent parameter, b, into the definition of the effective indentation modulus will adequately capture the anisotropic elastic indentation response of cubic crystals for any arbitrary orientation of the crystal lattice in the indentation zone. More specifically, they suggest: 1 1 t2s 1 t2i ¼b þ ; ð3Þ Eeff Es Ei where Es and ms denote the effective Young’s modulus and Poisson’s ratio, respectively [34,35], for a randomly textured polycrystalline sample. In the present study, our goal is first to establish the main underlying features of the Vlassak–Nix theory [34,35] that are central to the indentation analyses protocols discussed in this paper. More specifically, Eq. (3) combined with Eq. (1) implies that the indentation load, P , 3 continues to be directly proportional to h2e , even when the sample exhibits an anisotropic elastic response. If, indeed, this relationship holds true for any arbitrary orientation of the crystal, then it allows us to extract an orientationdependent indentation modulus (the parameter b is expected to be orientation dependent) from the indentation measurement on any crystal in the sample. Note that this would be an extremely important step in addressing the goals of the present work. In order to establish whether or not an indentation modulus can be defined for any crystal lattice orientation (from 3 the expected linear relationship between P and h2e ), we developed and employed a FE model. Note that the investigation described above cannot be conducted directly with measurements as easily, because of the difficulties in isolating the elastic portions of the indentation load–displacement curves [19]. Furthermore, numerical approaches such as FE models circumvent many of the difficulties and uncertainties faced in experimental measurements (e.g. identification of the initial point of contact, accuracy in the description of the geometry of indenter, unavoidable plasticity, friction at the contact surface), and are ideally suited for critical validations of fundamental concepts [18]. 3. Finite-element model of spherical nanoindentation The 3-D FE model to simulate elastic spherical indentation of anisotropic crystals was produced using the commercial FE code ABAQUS [36]. The FE model developed for this study is comprised of two 3-D bodies: (i) an elastically deformable sample with an initially flat surface discretized into 57,000 eight-noded, 3-D, continuum (C3D8) elements, and (ii) a rigid hemispherical indenter of radius 13.5 lm (i.e. the same size indenter as in prior experimental studies [19,20,22]). The size of the sample was selected as 18 lm 18 lm 9 lm to ensure that it is much
110
D.K. Patel et al. / Acta Materialia 79 (2014) 108–116
larger than the typical indentation zone size (observed to be about 1.2 lm from our simulations). The sample has been discretized into seven regions as shown in Fig. 1a to permit the use of a progressively higher mesh density as we approach the indentation zone directly below the indenter–sample contact surface. Discretizing the mesh in this manner allowed us to capture the stress and strain fields with sufficient accuracy while keeping the total number of elements in the FE model (as well as the computational cost) relatively low. Significant computational resources were required and utilized to keep the runtime to just under 30 min per simulation. The FE simulations described in this work were executed on two nodes (32 cores per node) provided by Trestles (part of the XSEDE supercomputing facility). A hard surface-to-surface, frictionless, contact definition was used to model the contact behavior between the indenter (master) surface and the elastically deformable sample (slave) surface. Displacement boundary conditions in the z-direction were imposed on the reference node (placed at the center of the indenter), which is tied rigidly to the entire surface of the indenter. The top surface of the sample was free to move while the bottom surface of the sample was constrained along the z-direction (indentation direction). The inputs to the FE model are the single-crystal elastic stiffness constants of the material of interest in the crystal frame (i.e. C11, C12 and C44) and the crystal lattice orientation. The lattice orientation of a single crystal with respect to the sample reference frame (defined in Euclidean 3-D space) can be described by a set of three well-defined sequential rotations. The three ordered angles ðu1 ; U; u2 Þ associated with these rotations are commonly known as Bunge–Euler angles [37]. Since the only relevant crystallographic direction in the spherical indentation experiment is the one parallel to the indentation direction (i.e. the indentation response is independent of any in-plane rotation of the sample), the crystal lattice orientation with respect to the indentation direction can be described using only two of the three Bunge–Euler angles [20,37]. We shall use g ¼ ðU; u2 Þ to denote the lattice orientation in the indentation experiments and models discussed in this work.
4. New approach to estimation of single-crystal elastic constants In this section, we present the proposed approach as a general framework. In the next section, we demonstrate a specific application and provide many more details of the computations involved. 4.1. Spectral representation of indentation modulus As described earlier, the first step of the proposed protocol involves establishing the functional dependence of the effective indentation modulus on the single-crystal elastic constants and the crystal lattice orientation. In order to establish this functional dependence, a large number of data points can be accumulated from the FE model
(b) Load, P (mN)
(a)
Finite displacement in the z-direction was imposed on the indenter and the corresponding total reaction force exerted by the slave surface on the indenter is predicted and used to determine the effective indentation modulus of the sample. The FE model developed in this study was validated by comparing the simulated load–displacement curves against the analytical prediction from Hertz’s theory for purely isotropic elastic deformation (see Fig. 1b). Additionally, a detailed sensitivity study was conducted to ensure that the details of the mesh and the details of the far-field boundary conditions did not influence significantly the predicted load–displacement response. The FE model developed in this study was validated by comparing the predicted indentation modulus against the analytical expression reported by Vlassak and Nix [35] in Table 1. It is clear that the FE predictions are in excellent agreement with the theoretically expected values. Furthermore, the FE models predicted that a linear relationship 3 between P and h2e holds for the anisotropic elastic indentations, in full accord with the theories presented by Vlassak and Nix [34,35]. Some examples of the FE-predicted load– displacement relations are presented in Fig. 2. This is highly encouraging and provides a path forward for addressing the inverse problem described earlier.
Displacement, he (nm) Fig. 1. (a) The 3-D FE mesh used to simulate elastic anisotropic response under the indenter (isometric view of the region under the indenter tip). (b) Load–displacement response of FE simulation vs. Hertz theory for the isotropic response case.
D.K. Patel et al. / Acta Materialia 79 (2014) 108–116 Table 1 Indentation moduli values for the cubic crystals predicted from the FE model and theoretical values reported in Vlassak and Nix [35] for symmetric orientations. Material and corresponding elastic constants (GPa)
Orientation direction (h, k, l)
Theoretical Eind (GPa)
FE prediction Eind (GPa)
Aluminiumb C11 = 107.3 C12 = 60.90 C44 = 28.3
(1 0 0) (1 1 0) (1 1 1)
79 80 81
78.96 80.93 80.25
Copperb C11 = 170.2 C12 = 114.9 C44 = 61.0
(1 0 0) (1 1 0) (1 1 1)
129 138 141
129.65 137.19 141.05
b-Brass (CsCl)a C11 = 126.5 C12 = 107.7 C44 = 80.3
(1 0 0) (1 1 0) (1 1 1)
95 112 117
93.87 109.40 115.47
a b
Taken from Vlassak and Nix [35]. Taken from Simmon and Wang [38].
can be interpreted as results of probing single crystals of a selected phase through indentations in different crystallographic directions. The functional dependencies of the indentation modulus described above can therefore be succinctly expressed as Eind ðg; C 11 ; C 12 ; C 44 Þ. Our first goal is to capture this functional dependence of the indentation modulus using highly efficient spectral representations. In this paper, we have employed the spherical surface harmonics (SSH) [37] denoted as K ml ðgÞ to serve as a Fourier basis for compact representation of this function over the orientation space of interest. The advantages and compactness of using these basis functions on the orientation space have already been discussed in our prior publications [39,40]. For the variables denoting the elastic stiffness constants, we explore here the use of Legendre polynomials, P n ðÞ [41] as a Fourier basis. Note that the Legendre polynomials form a Fourier basis over the range [1, +1], and therefore the stiffness parameters need to be suitably rescaled. Each rescaled elastic stiffness constant can be defined as: ~ ¼ 2C C min C max ; C C max C min
ð4Þ
where C min and C max are the respective minimum and maximum values for the specific elastic stiffness constant. The ~ 11 ; C ~ 12 ; C ~ 44 Þ is then expressed function of interest, Eind ðg; C as ~ 11 ; C ~ 12 ; C ~ 44 Þ :¼ Eind ðg; C
MðlÞ X 1 X 1 X
~ 11 ÞP r ðC ~ 12 ÞP s ðC ~ 44 Þ: Amqrs K ml ðgÞP q ðC l
ð5Þ
l¼0 m¼1 q;r;s¼0
In Eq. (5), Amqrs represent the Fourier coefficients. The l function MðlÞ denotes the number of terms needed in the enumeration of index, m; these numbers are a function of the index l and are determined by crystal symmetries [37].
Load (mN)
predictions. Each FE simulation provides one data point as an estimate of the indentation modulus for a specific combination of one crystal orientation and one set of values for the elastic constants. A least-squares fit between the FE-predicted values of P 3 and h2e in each simulation yields an estimate of the indentation modulus (Eind ) (see Eq. (1)). Implicitly, the indentation modulus so obtained depends on the crystal lattice orientation (applied here uniformly to all elements in the sample) and the values of the single-crystal elastic constants used for the sample material (e.g. C11, C12 and C44 for cubic crystals) [38]. Therefore, the different measurements on the differently oriented grains in the polycrystalline sample
111
he3/2 (nm) Fig. 2. P vs. (he)(3/2) plots obtained for copper for different orientations spread over a cubic-transversely isotropic FZ (highlighted area in the subplot).
112
D.K. Patel et al. / Acta Materialia 79 (2014) 108–116
The most common approaches used in the literature to establish the Fourier coefficients generally exploit the orthonormal properties of the Fourier basis [25,42–44]. However, such approaches require evaluation of the function at an extremely large number of locations in the compounded space of the complete ranges of all the independent variables involved. Since the function evaluation here is being attempted through a computationally expensive FE model, the traditional approach of establishing Fourier coefficients is not very practical for the present problem. Consequently, we explore in this work an approach for establishing the values of the Fourier coefficients in Eq. (5) using ordinary least-squares (OLS) regression analysis. In order to use a regression technique, it is necessary to truncate the Fourier series in Eq. (5) to a finite number of terms in the expansion. Simplifying the notation a little, Eq. (5) can be expressed as: ~ 11 ; C ~ 12 ; C ~ 44 Þ Eind ðg; C
~ ~q L X X
~ ~ ~ Aqrs L K L ðgÞP q ðC 11 ÞP r ðC 12 ÞP s ðC 44 Þ:
dataset. Over-fitting needs to be avoided as it will likely produce unreliable estimates in subsequent application of the Fourier function (especially for new data points not included in the calibration dataset). In order to track the potential for over-fitting of the Fourier coefficients obtained at each level of truncation ~ and ~q), leave-one-out-cross-validation (i.e. selection of L (LOOCV) can be employed [47]. LOOCV allows an objective selection of the truncation levels by establishing the fit N times (this is size of the dataset), while leaving one data point out of the regression each time. In doing so, LOOCV will quantify the contribution of each data point to the Fourier coefficients established in each regression. Given a large N, for an over-fitted Fourier expansion, the exclusion of a single data point will cause significant change in the coefficients, whereas for a good fit this change will be negligible. In order to quantify the robustness of the fit, the following objective measures could be utilized:
ð6Þ
I. Mean absolute error of the fit, e, defined as:
L¼0 q;r;s¼0
where L enumerates each distinct combination of ðm; lÞ in ~ and ~ the SSH series expansion, and L q denote the truncation limits in the SSH and the Legendre polynomial expansions. Note also that a single limit was selected in Eq. (6) for all three Legendre polynomial bases used in this equation. This was merely done for simplicity here, and different limits can be selected for the different variables as needed in any specific application. As described earlier, each FE simulation of the indentation produces one discrete data point for a combination of one crystal orientation and one set of values for the cubic elastic stiffness constants. Consider a dataset with N such data points collected from multiple runs of the FE simulations for different combinations of crystal orientations and cubic elastic stiffness constants. n Such a FE-generated FEðnÞ ~ ðnÞ ; C ~ ðnÞ ; C ~ ðnÞ ; dataset can be denoted as Eind ; gðnÞ ; C 11 12 44 o n ¼ 1; 2; . . . ; N . Regression analysis is then formulated as a minimization of the difference (sum of the squares of the residual) between the FE-generated dataset and the truncated Fourier representation shown in Eq. (6). In other words, OLS [45,46] aims to solve the following minimization problem: ! ~ q ~ N L X FEðnÞ X X qrs ðnÞ ðnÞ ðnÞ ðnÞ ~ ~ ~ min AL K L g P q C 11 P r C 12 P s C 44 : Eind qrs AL n¼1 L¼0 q;r;s¼0
N 1X jeðnÞ j N n¼1 ~ q ~ N L 1X FEðnÞ X X qrs ðnÞ ðnÞ ðnÞ ðnÞ ~ ~ ~ AL K L g P q C 11 P r C 12 P s C 44 : ¼ E N n¼1 ind L¼0 q;r;s¼0
e ¼
ð8Þ
Median absolute deviation (MAD) of error of the fit as: MADe ¼ Median eðnÞ Median eð1Þ ; eð2Þ ; . . . ; eðN Þ ; n ¼ f1; 2; . . . ; N g:
ð9Þ
II. Mean absolute error of LOOCV, eCV , and MAD of LOOCV, MADCV , defined analogously as the mean and median of the error for the test data point in each of the N repeated fits. The two measures of error of fit defined above will show ~ and ~q, whereas improvement of fit with higher values of L the two measures of error of LOOCV are expected to show ~ and ~ a decline in robustness of fit with higher values of L q accounting for over-fit of data. Therefore, a compromise is often made in choosing the best fit based on the values of all four measures defined above. In the next step, Eq. (6) can be used to address the inverse problem of estimating the unknown single-crystal elastic stiffness constants for a selected phase in a polycrystalline sample.
ð7Þ
Note that the values of the Fourier coefficients established using the regression method described above can be sensitive to the truncation levels used, i.e. the values ~ and ~ of L q. Since we do not generally know a priori the ~ and ~ right values of L q in any specific application, this approach needs a few repeated trials. Although higher val~ and ~ ues of L q will generally produce a lower error (difference between the FE dataset and the fitted Fourier expansion), they may also result in over-fitting of the
4.2. Protocol to extract single-crystal elastic constants As described earlier, developing and validating an inverse solution methodology to extract grain-scale elastic properties for a selected phase in a polycrystalline material is the main objective of this work. The protocol described here aims to match the measurements of a set of indentation moduli on differently oriented grains (typically measured using orientation imaging microscopy [21]) of a
D.K. Patel et al. / Acta Materialia 79 (2014) 108–116
selected phase whose independent elastic stiffness constants are unknown (denoted as C ) with the function established in the previous step (Eq. (6)). For cubic materials, for sim~ , C ~ , C ~ ) correspond to rescaled plicity of notation, let (C 1 2 3 ~ , C ~ , C ~ ) using Eq. (4). unknown elastic constants (C 11 12 44 ðjÞ
Let g ; j ¼ 1; 2; . . . ; J denote the specific orientations ðjÞ
where the indentation measurements were made, and Eind denote the corresponding measurements of the indentation moduli. Let the corresponding theoretical value for each of the measurement, predicted from the Fourier representa ðjÞ ~ ; C ~ ; C ~ . The task tion (Eq. (6)), be denoted as, Eind gðjÞ ; C 1 2 3 of estimating the unknown single-crystal elastic constants then reduces to minimizing the difference between the val ðjÞ ~ ; C ~ ; C ~ and EðjÞ for all of the available ues of Eind gðjÞ ; C ind 1 2 3 measurements on the selected material phase. The corresponding minimization problem can be expressed as: J h X i2 ðjÞ ðjÞ ~ ;C ~ ;C ~ ¼ min ~ ; C ~ ; C ~ : min f C Eind Eind gðjÞ ; C 1 2 3 1 2 3
~ ;C ~ ;C ~ C 1 2 3
~ ;C ~ ;C ~ C 1 2 3 j¼1
ð10Þ
The above minimization problem can be addressed by evaluating and equating the relevant derivatives of the function to zero: ~ ;C ~ ;C ~ ¼ @f Fi C 1 2 3 ~ i @C J h i ðjÞ ðjÞ ~ ~ ~ X ðjÞ ðjÞ ~ ~ ~ @Eind g ; C 1 ; C 2 ; C 3 ¼ 2 EðjÞ ¼ 0: ind E ind g ; C 1 ; C 2 ; C 3 ~ i @C j¼1
ð11Þ
ðkÞ
1
~ i ½J ij F j ; ¼C @F i J ij ¼ : ~ j ðkÞ @C ~ ¼C ~ C Ci
n
ð12Þ
n
where the superscript (k) denotes the iteration number. 5. Case study: as-cast polycrystalline Fe–3% Si The new protocols presented in this paper to extract the single-crystal elastic stiffness parameters have been applied to an as-cast polycrystalline sample of Fe–3% Si. As described in the previous section, the first step involves ~ 11 ; C ~ 12 ; C ~ 44 in a spectral establishing the function Eind g; C form (see Eq. (6)) using regression methods on datasets assembled from results produced by the FE models described in Section 3. As described earlier, each FE simulation of the spherical indentation for a selected combination of crystal orientation and cubic elastic stiffness parameters yields one data point. In the present example, the desired dataset was accumulated by executing a total of 2700 FE simulations that employed 300 distinct sets of the cubic elastic stiffness constants [38] compounded with nine distinct crystal orientations. The elastic stiffness constants were selected in the ranges of 50 6 C 11 6 250, 40 6 C 12 6 150 and 15 6 C 44 6 125 (see Fig. 3), which were chosen so that they cover the typical ranges for most structural cubic metals and alloys of interest. This set covers a cubic anisotropy ratio, defined as A ¼ 2C 44 = ðC 11 C 12 Þ, in the range ð0 < A < 8Þ. The crystal orientations were selected such that they cover the relevant fundamental zone (FZ) of orientations of interest for the present application [25] (see also the inset of Fig. 2): " FZ ¼ ðU; u2 Þj cos1
cos u2 pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 1 þ cos2 u2
!
p 6 U2 6 ; 2
06u6
# p : 4
ð13Þ
C11 (GPa)
The system of nonlinear equations in Eq. (11) can be solved using various numerical approaches. In this paper, we have employed an iterative Newton–Raphson algorithm to solve this system of nonlinear equations. This approach can be summarized as [46]:
ðkþ1Þ
113
C44 (GPa) C12 (GPa) Fig. 3. 300 Distinct sets of independent elastic stiffness constants ðC 11 ; C 12 ; C 44 Þ used in the present study to establish the spectral representations.
D.K. Patel et al. / Acta Materialia 79 (2014) 108–116
The symmetrized spherical surface harmonics basis functions appropriate for this case study are those exhibit ing cubic-transversely isotropic symmetry denoted as K_ ml ðgÞ [37]. In this description of symmetry, the first term reflects the crystal symmetry and the second term reflects the sample symmetry (arises because the indentation direction is unaffected by any in-plane rotation of the sample). As described in the previous section, the Fourier coefficients of Eq. (6) can be determined using OLS regression (Eq. (7)). The values of e (Eq. (8)) and MADe (Eq. (9)) were com~ puted for the training set for different selections of ~ q and L, and are summarized in Fig. 4. The mean and median mea~ sures of the error decrease with increasing values of ~ q and L as expected (see Fig. 4a). The LOOCV analyses (Fig. 4b) indicates over-fitting of the data points at high values of ~ For instance, the eCV values for the fits both ~q and L. ~ > 12 and ~ obtained using L q > 3 start to increase significantly. Examination of the error measures shown in ~ ¼ 12 and ~ Fig. 4a and b suggests that L q ¼ 3 provide a good fit for the dataset acquired in this study using FE simulations. This level of truncation produces a total of 348 Fourier coefficients; these are utilized in this study to describe the functional dependence of indentation modulus on orientations and the independent elastic constants. The accuracy of the truncated spectral representation of the indentation moduli over the selected range of cubic elastic constants for several orientations was compared against the values predicted by the previously validated FE model and the analytical solution provided by Vlassak and Nix [35]. Representative results are shown in Fig. 5 for a selected range of elastic anisotropic ratios for two orientations that provide the extremal responses. The indentation modulus in this plot was normalized by the indentation modulus of a polycrystalline aggregate, following the approach outlined by Vlassak and Nix [34] and using the expressions provided by Hashin and Shtrikman [48]. The error between the FE prediction and the spectral representation for any given set of cubic elastic constants (within the selected bounds mentioned earlier) and any orientation was found to be within 4%. It is clear that the spectral representation developed in this work captures well the variations in the indentation moduli for different orientations and different degrees of elastic anisotropy in the sample. We now turn our attention to estimating the singlecrystal elastic stiffness parameters of Fe–3% Si, for which indentation measurements on differently oriented single crystals were already reported in prior work [20]. These results are summarized in Table 2, where the orientations are also plotted on the standard inverse pole figure. The effective indentation moduli reported in Table 2 for Fe– 3% Si were computed using Eqs. 1–3, where the Young’s modulus and Poisson ratio of the indenter were taken as 1000 GPa and 0.07, respectively [49]. Estimating the single-crystal constants of Fe–3% Si then reduces to minimizing the difference between the values predicted by the Fourier function established earlier and the discrete
Fig. 4. (a) The mean absolute error by varying numbers of SSH basis ~ and degree of Legendre polynomial ð~qÞ with MAD of functions ðLÞ absolute error as the error bars. (b) The mean absolute error of LOOCV ~ and degree of Legendre by varying numbers of SSH basis functions ðLÞ polynomial ð~qÞ with MAD of absolute error of LOOCV as the error bars.
1.1
Normalized Indentation Modulus
114
1.05 Spectral Function {100}
1
Spectral Function {111}
0.95
FE Prediction Analytical Function
0.9 0.85 0.8 0.75 1
3 5 Anisotropy Ratio (A)
7
Fig. 5. Comparison of the spectral representation of normalized indentation modulus as a function of anisotropic ratio vs. the predictions of the FE model and the analytical model of Vlassak and Nix [35].
D.K. Patel et al. / Acta Materialia 79 (2014) 108–116
115
Table 2 (a) Measured values of the effective indentation modulus, Eeff, for 11 different orientations, and (b) representation of the orientations in an inverse pole figure map.
(a)
(b)
measured indentation moduli reported in Table 2. This minimization was accomplished using a Newton–Raphson scheme described in the previous section (see Eqs. 10–12). The estimated values of the single-crystal elastic constants for the as-cast Fe–3% Si polycrystalline sample from this study are C 11 ¼ 215:7 GPa, C 12 ¼ 131:6 GPa and C 44 ¼ 122:2 GPa. The single-crystal elastic constants obtained using the inverse solution approach described in this paper are within 5% of the typical values reported in literature for Fe–3% Si ðC 11 ¼ 225:28 GPa; C 12 ¼ 135:1 GPa and C 44 ¼ 123:8 GPaÞ [50–52]. 6. Discussion The accuracy of the spectral framework mentioned above can be improved with additional measurements of indentation moduli on additional crystal orientations. Note that we used only 11 measured indentation moduli in the present case study. One might actually argue that we need only three measurements (on three distinct crystal orientations) to estimate the three unknown elastic stiffness constants in cubic crystals. However, because of the regression methods used in the approach described here, it is important to measure indentation moduli for orientations that are well spread out in the fundamental zone of orientations. For example, if we were to repeat the inverse solution presented here by considering only three orientations clustered in any one of the corners of the inverse pole figure (depicts a certain projection of the fundamental zone of orientations) shown in Table 2, it would produce highly
erroneous estimates of the single-crystal elastic stiffness constants. In fact, it is best to ensure that there at least some measurements in each of the three corners of the inverse pole figure shown in Table 2, as these orientations typically produce the maximum contrast in the indentation moduli. It is also likely that there are better representations of the indentation modulus as a function of the orientation and the elastic stiffness constants. Although we have had a lot of experience with the use of SSH as a Fourier basis for the orientation variable, we only have very limited prior experience in spectral representation of the elastic stiffness parameters. Therefore, it is entirely possible that one might accomplish a more accurate representation of the indentation modulus function (Eq. (5)) with a different Fourier basis, which might in turn improve the accuracy of the estimates of the elastic stiffness constants. 7. Conclusion A new two-step procedure has been developed to estimate the single-crystal elastic stiffness parameters from polycrystalline samples using spherical nanoindentation and orientation measurements combined with FE simulations. The first step involves a one-time capture of the functional dependence of the indentation modulus on the lattice orientation at the indentation site and the unknown single-crystal elastic constants for any presumed crystal symmetry. This step was accomplished in this study for cubic crystals by employing a FE model of spherical
116
D.K. Patel et al. / Acta Materialia 79 (2014) 108–116
indentation. In the second step of the approach presented here, the unknown single-crystal elastic constants for a selected phase are estimated through a regression technique. The accuracy and viability of the proposed approach were demonstrated for an as-cast cubic polycrystalline Fe–3% Si sample. Acknowledgements
[18] [19] [20] [21] [22] [23] [24] [25] [26]
This work was supported by ARO through Grant No. W911NF-10-1-0409. Computations described in this work used the Extreme Science and Engineering Discovery Environment (XSEDE), which is supported by National Science Foundation Grant No. OCI-1053575. The codes for fast computation of GSH coefficients were provided by Material Resources, LLC, Dayton, OH. References [1] Roters F, Eisenlohr P, Bieler TR, Raabe D. Crystal plasticity finite element methods in materials science and engineering. New York: Wiley; 2011. [2] Raabe D, Klose P, Engl B, Imlau K-P, Friedel F, Roters F. Adv Eng Mater 2002;4:169. [3] Haddadi H, Bouvier S, Banu M, Maier C, Teodosiu C. Int J Plast 2006;22:2226. [4] Zhao Z, Roters F, Mao W, Raabe D. Adv Eng Mater 2001;3:984. [5] Anand L, Kalidindi SR. Mech Mater 1994;17:223. [6] Bronkhorst CA, Kalidindi SR, Anand L. Philos Trans R Soc Lond Ser A Math Phys Eng Sci 1992;341:443. [7] Kalidindi SR, Bronkhorst CA, Anand L. J Mech Phys Solids 1992;40:537. [8] Kalidindi SR, Bhattacharyya A, Doherty RD. Proc R Soc Lond Ser A Math Phys Eng Sci 2004;460:1935. [9] Bhattacharyya A, El-Danaf E, Kalidindi SR, Doherty RD. Int J Plast 2001;17:861. [10] Kalidindi SR, Schoenfeld SE. Mater Sci Eng A Struct Mater Prop Microstruct Process 2000;293:120. [11] McDowell DL. Int J Plast 2010;26(9):1280. Special issue in honor of David L. McDowell. [12] Kalidindi SR, Bhattacharyya A, Doherty RD. Proc R Soc A 2004;460:1935. [13] Kalidindi SR, Anand L. J Mech Phys Solids 1994;42:459. [14] Bachu V, Kalidindi SR. Mater Sci Eng A Struct Mater Prop Microstruct Process 1998;257:108. [15] Shan ZW, Mishra RK, Syed Asif SA, Warren OL, Minor AM. Nat Mater 2008;7:115. [16] Uchic MD, Dimiduk DM, Florando JN, Nix WD. Science 2004;305:986. [17] Li N, Mara NA, Wang YQ, Nastasi M, Misra A. Scr Mater 2011;64:974.
[27] [28] [29] [30] [31] [32]
[33] [34] [35] [36] [37] [38] [39] [40] [41] [42] [43] [44] [45] [46] [47] [48] [49] [50] [51] [52]
Donohue BR, Ambrus A, Kalidindi SR. Acta Mater 2012;60:3943. Kalidindi SR, Pathak S. Acta Mater 2008;56:3523. Pathak S, Stojakovic D, Kalidindi SR. Acta Mater 2009;57:3020. Proust GKS. TMS Lett 2004;7:151. Pathak S. Development and validation of a novel data analysis procedure for spherical nanoindentation. Drexel University; 2009. Bhattacharya A, Nix W. Int J Solids Struct 1988;24:881. Suresh S, Giannakopoulos A. Acta Mater 1998;46:5755. Adams BL, Kalidindi SR, Fullwood DT. Microstructure sensitive design for performance optimization. Elsevier Science; 2012. Kalidindi SR, Houskamp JR, Lyons M, Adams BL. Int J Plast 2004;20:1561. Li DS, Garmestani H, Schoenfeld S. Scr Mater 2003;49:867. Hertz, H. On the contact of elastic solids. Hertz’s Miscellaneous Papers, Macmillan; 1986. Johnson KL, Johnson KKL. Contact mechanics. Cambridge: Cambridge University Press; 1987. Willis J. J Mech Phys Solids 1966;14:163. Cammarata R, Schlesinger T, Kim C, Qadri S, Edelstein A. Appl Phys Lett 1990;56:1862. Farthing L, Weihs T, Kisker D, Krajewski J, Tang M, Stevenson D. Hardness and elastic modulus measurements in CdTe and ZnTe thin film and bulk samples and ZnTe–CdTe superlattices. MRS proceedings, vol. 130. Cambridge: Cambridge University Press; 1988. Vlassak J, Ciavarella M, Barber J, Wang X. J Mech Phys Solids 2003;51:1701. Vlassak JJ, Nix W. Philos Mag A 1993;67:1045. Vlassak JJ, Nix W. J Mech Phys Solids 1994;42:1223. ABAQUS. Simulia Providence, RI: Dassault Systemes; 2009. Bunge HJ. Texture analysis in materials science: mathematical methods. Butterworths; 1982. Simmons G, Wang H. Single crystal elastic constants and calculated aggregate properties: a handbook. Cambridge, MA: MIT Press; 1971. Fast T, Knezevic M, Kalidindi SR. Comput Mater Sci 2008;43:374. Fullwood DT, Niezgoda SR, Adams BL, Kalidindi SR. Prog Mater Sci 2010;55:477. Copson ET. An introduction to the theory of functions of a complex variable. Oxford: Clarendon Press; 1950. Kalidindi SR, Duvvuru HK, Knezevic M. Acta Mater 2006;54:1795. Kalidindi SR, Landi G, Fullwood DT. Acta Mater 2008;56:3843. Knezevic M, Al-Harbi HF, Kalidindi SR. Acta Mater 2009;57:1777. Bishop CM. Pattern recognition and machine learning. Berlin: Springer Verlag; 2006. Press WH. Numerical recipes: the art of scientific computing. 3rd ed. Cambridge: Cambridge University Press; 2007. C ¸ ecßen A, Fast T, Kumbur EC, Kalidindi SR. J Power Sources 2014;245:144. Hashin Z, Shtrikman S. J Mech Phys Solids 1962;10:335. Field J, Swain M. J Mater Res (USA) 1993;8:297. Alberts H, Wedepohl P. Physica 1971;53:571. Machova A, Kadecˇkova´ S. Czech J Phys B 1977;27:555. Routbort J, Reid C, Fisher E, Dever D. Acta Metall 1971;19:1307.