Work hardening in micropillar compression: In situ experiments and modeling

Available online at www.sciencedirect.com

Acta Materialia 59 (2011) 3825–3840 www.elsevier.com/locate/actamat

Work hardening in micropillar compression: In situ experiments and modeling D. Kiener a,d, P.J. Guruprasad b,⇑, S.M. Keralavarma b, G. Dehm c,d, A.A. Benzerga b,e a

University of California at Berkeley, Lawrence Berkeley National Laboratory, National Center for Electron Microscopy, Berkeley, CA 94720, USA b Department of Aerospace Engineering, Texas A&M University, College Station, TX 77843-3141, USA c Erich Schmid Institute of Materials Science, Austrian Academy of Sciences, Leoben, Austria d Department of Materials Physics, University of Leoben, Leoben, Austria e Materials Science & Engineering Graduate Program, Texas A&M University, College Station, TX 77843, USA Received 9 December 2010; received in revised form 1 March 2011; accepted 4 March 2011 Available online 12 April 2011

Abstract Experimental measurements and simulation results for the evolution of plastic deformation and hardening in micropillars are compared. The stress–strain response of high-symmetry Cu single crystals is experimentally determined using in situ micropillar compression. Discrete dislocation simulations are conducted within a two-dimensional plane-strain framework with the dislocations modeled as line singularities in an isotropic elastic medium. Physics-based constitutive rules are employed for an adequate representation of hardening. The numerical parameters entering the simulations are directly identified from a subset of experimental data. The experimental measurements and simulation results for the flow stress at various strain levels and the hardening rates are in good quantitative agreement. Both flow strength and hardening rate are size-dependent and increase with decreasing pillar size. The size effect in hardening is mainly caused by the build-up of geometrically necessary dislocations. Their evolution is observed to be size-dependent and more localized for smaller sample volumes, which is also reflected in local crystal misorientations. Ó 2011 Acta Materialia Inc. Published by Elsevier Ltd. All rights reserved. Keywords: Discrete dislocation dynamics; Geometrically necessary dislocations (GNDs); Flow stress; Hardening; In situ pillar compression

1. Introduction Bulk materials harden when they are plastically deformed. In pure materials plastic resistance comes from dislocation interactions and intersections. The mechanisms for this resistance vary from one stage of deformation to another, with many details being material specific. However, the main features of hardening remain the same, as manifested, for example, by universal values of the workhardening rates when normalized by the material stiffness. This holds for both single crystals and polycrystals [1].

⇑ Corresponding author.

E-mail address: Guruprasad).

[email protected]

(P.J.

Over the past few years, new experimental methods have been developed to probe the mechanical response of materials at the scale of their microstructures. These techniques thus permit fundamental issues in crystal plasticity to be addressed and the limits of current models to be examined. Among such methods, micropillar compression has been extensively used [2–8] (see Ref. [9] for a recent review). In general, a common trend emerges from pillar compression experiments: smaller is stronger. However, the strength of the scaling of flow stress with pillar diameter varies from one experimental data set to another. In addition, very few attempts have been made to analyze the hardening behavior. Therefore, there is need for further experimental investigation coupled with analysis of plasticity in micronand submicron-sized objects, especially in the absence of imposed strain or stress gradients. In particular, design of

1359-6454/$36.00 Ó 2011 Acta Materialia Inc. Published by Elsevier Ltd. All rights reserved. doi:10.1016/j.actamat.2011.03.003

3826

D. Kiener et al. / Acta Materialia 59 (2011) 3825–3840

experiments that allow an investigation of hardening at the micron and submicron scales has far-reaching implications on physics-based plasticity modeling and simulation efforts. Phenomenological models of plasticity do not include adequate representation of microstructural effects at the dislocation scales. In addition, current continuum models are incapable of providing a rationale for micropillar plasticity and size effects. Under such circumstances, recourse to lower-scale, higher-resolution analyses is necessary. Fully discrete atomic-level methods, such as molecular dynamics (MD), have been used for understanding plasticity in nanoscale domains (e.g. [10,11]). However, MD is incapable of resolving sample sizes ranging from 100 nm to over 10 lm, which is the typical range of pillar diameters considered in the experiments thus far. Alternatively, semidiscrete analyses may be used which are based on dislocation theory, i.e. linear elasticity for long-range dislocation interactions as well as suitably specified atomic-level input. Progress on discrete dislocation dynamics (DD) simulations of micropillar plasticity has recently been reviewed by Uchic et al. [9]. Subtleties aside, three-dimensional (3D) simulations have essentially confirmed two strengthening mechanisms: (i) the role of source strength distribution when sources are available [12–15]; and (ii) the imbalance between rates of dislocation generation and dislocation annihilation/immobilization when there is paucity of sources (e.g. [12,16]). As noted by Uchic et al. [9], it is remarkable that some two-dimensional (2-D) DD simulations [17] have identified such strengthening mechanisms; see Ref. [18,19] for a more recent discussion. None of the DD simulations above have discussed the transition from the features observed in micropillar compression to bulk-like behavior, where forest-hardening processes are generally expected to result in a size-independent response. However, the experimental results of Kiener et al. [8] as well as the DD simulations of Guruprasad and Benzerga [20] have independently revealed that size-dependent, steady hardening can be obtained up to very large strains in micron-scale samples. Such behavior cannot be rationalized in terms of previously established strengthening mechanisms. The lack of investigations centered on size-affected hardening is not commensurate with the critical need for improved hardening models in continuum descriptions, and is in part due to the absence of clear trends in most previously published experiments. Driven by previous investigations [8,20], the objective of this work is to combine experiments and DD modeling to investigate (i) the propensity of micropillars to harden; and (ii) the size dependence of hardening. In situ microcompression experiments were carried out on high-symmetry Cu single-crystalline micropillars made by focused ion beam (FIB) milling, with either circular or square cross-sections and diameters or side lengths from 8.2 down to 0.9 lm [8]. The choice of Cu is motivated by its technological use in micro- and nanoelectronics applications, the availability of tensile data for micron-thick thin films [21–

24] as well as microtension specimens [25–28], and by the vast literature on the hardening behavior of macroscopic samples [29–32,1]. In an additional set of experiments, pillars coated with a thin TiN film and pillars on a stiff MgO substrate were also tested to investigate the effect of boundary conditions [33]. Strong effects of size on flow strength and hardening and a weak effect of cross-section shape were evidenced. The discrete DD formulation follows that of Guruprasad and Benzerga [20]. Pillars with a square cross-section are modeled using a plane-strain approximation and the paradigm of 2.5-D DD [34], where atomic-level input is incorporated through a set of constitutive rules for close-range interactions. Chief among these are rules that lead to dynamic multiplication at junction-anchored Frank–Read (FR) sources and to effective dislocation storage at dynamically formed junctions. This 2.5-D DD framework, which enhances the standard 2-D model [35], has predicted a range of features observed experimentally in bulk plasticity including Taylor hardening, stage I and stage II hardening with rates in keeping with experimental measurements, and refinement of the dislocation structure upon hardening [34]. In carrying out the quantitative comparison between experiments and modeling, a scheme was developed whereby few simulation parameters were calibrated so as to obtain a good representation of the behavior of the largest specimens. The calibrated parameters are related to the initial dislocation-source and dynamic dislocation-junction populations. Then, with all constitutive parameters fixed, the mechanical response beyond yielding is predicted for all other sizes. The computations are compared quantitatively with the experimental observations and provide insight into the mechanisms leading to the observed size dependence in hardening by the build-up of a density of geometrically necessary dislocations (GNDs) in a sizedependent manner, a result which could not be predicted from classical theories. 2. Experiments 2.1. Methods The mechanical response of four sets of Cu micropillars and a set of Cu microtensile specimens (Table 1) was investigated using flat punch indentation [36] and a custom-made microtension machine [26,27], respectively. In both techniques, FIB machining is used to cut out specimens of desired shape and dimension. The nominal stress r and strain e are determined from the measured load, P, and tip displacement, u, using the simple expressions: r ¼ P =A;

e ¼ juj=H ;

ð1Þ

where H is the initial height of the pillar measured from the base and A is the initial cross-sectional area measured at half the pillar height after FIB machining.

D. Kiener et al. / Acta Materialia 59 (2011) 3825–3840

3827

Table 1 Geometry and cross-section shape, crystal orientation, and the minimum and maxmimum equivalent diameter, D, of the specimens tested. Test Compression

Tension

Geometry/ cross-section

Orientation/ constraint

Dmin (lm)

Dmax (lm)

Standard deviation (%)

h1 0 0i Cu h1 0 0i Cu on MgO

0.8 0.98

6.7 1.08

11 6

h1 1 1i Cu h1 0 0i Cu on MgO h1 1 1i Cu/TiNcoated

0.95 0.41 0.72

5.71 1.07 7.27

4 4 2

h2 3 4i Cu

3.0

–

–

Both cylindrical and “square” micropillars were milled using 30 keV Ga+ ions and a final milling current of 100 pA in a dual-beam FIB/scanning electron microscope (SEM) workstation. The pillars with square cross-section are taper-free, while cylindrical ones are tapered. Each sample is represented by a dimensional parameter, D, corresponding to the diameter of a circle of area A. Table 1 lists ranges of D for the specimens tested. Most pillars originated from high-purity (99.999%) melt-grown single-crystal rods. The first two sets of specimens were h1 0 0i cylinders (minimum diameter Dmin = 0.8 lm) and h1 1 1i pillars with a square cross-section (Dmin = 0.95 lm). The third and fourth sets consisted respectively of TiN-coated h1 1 1i pillars and h1 0 0i pillars on an MgO substrate. To prepare the latter, 1 lm thick single-crystal films grown on MgOh1 0 0i [37] were FIB structured. The samples had an aspect ratio of about 2:1 as suggested in Ref. [38], except for pillars from the fourth set where aspect ratios of 1–1.3 were achievable. Finally, taper-free 3  3  6 lm3 microtensile h2 3 4i Cu specimens were fabricated and tested following the procedure detailed in Ref. [26,27]. The important feature of the tensile test is the low lateral stiffness, allowing the sample to deform in single slip. Stress and strain are evaluated using the counterpart of Eq. (1) and H = 6 lm as the gauge length. This data will be used solely for the purpose of model parameter calibration, as will be explained in Section 3.2. Immediately after fabrication, the samples were transferred from the FIB into an SEM to minimize exposure to air. Sample testing was performed in situ in this tungsten filament SEM using a microindenter equipped with a flat conical diamond tip. Here, most experiments were displacement-controlled at an applied strain rate of

3  103 s1. Further details may be found in Refs. [8,33,36] for pillar fabrication and testing, respectively, and in Ref. [26,27] for microtensile experiments. Two compression samples were further analyzed using electron backscatter diffraction (EBSD) to study local crystal orientation changes. One sample was deformed in load control up to a strain of 0.28, the other in displacement control to a strain of 0.18. The high compressive strains yielded large slip steps on the sample surface that complicated EBSD investigation. Therefore, the sample surfaces were FIB polished using a 100 pA ion current. EBSD scans were performed with a step size of 25 nm. 2.2. Experimental results Representative nominal stress–strain curves for the first three sets of micropillars (Table 1) are shown in Fig. 1. The 18 tested h1 0 0i Cu samples on MgO (set 4) cover only a limited size range. Therefore, their responses are not shown here but are reported in Ref. [39]. It should be mentioned that the curves in Fig. 1 differ in the used data output rate. The data for h1 0 0i Cu (Fig. 1a) was recorded with 32 data points per second, whereas the data for h1 1 1i Cu (Fig. 1b) and TiN-coated h1 1 1i Cu (Fig. 1c) had 16 and 4 points per second, respectively. Noise tends to obscure the fine details of the curves with 32 data points per second. This could be changed by meaningful binning of the data. In the case of only 4 data points per second (Fig. 1c), the resulting stress– strain curves are devoid of noise, but lose some of the fine details. For a discussion of size-dependent hardening, the stress values at various amounts of strain were extracted and converted from normal stress to shear stress using the nominal Schmid factors of fs = 0.408 for h1 0 0i Cu and fs = 0.278

3828

D. Kiener et al. / Acta Materialia 59 (2011) 3825–3840

Fig. 2. Size dependent engineering shear stress for strain values ranging from 0.05 to 0.20 for: (a) tapered circular shaped h1 0 0i Cu; (b) straight square shaped h1 1 1i Cu and straight square shaped h1 1 1iCu with TiN top coating. The straight lines represent a best fit to the data sets.

Fig. 1. Representative stress–strain curves showing a sample size effect for tested samples of: (a) tapered circular shaped h1 0 0i Cu; (b) straight square shaped h1 1 1i Cu; and (c) straight square shaped h1 1 1i Cu with a TiN top coating.

for h1 1 1i Cu. Changes of these values during straining were not taken into account. The nominal shear stresses for strain values of 0.05, 0.10, 0.15 and 0.20 are plotted

in Fig. 2a for h1 0 0i Cu and in Fig. 2b for h1 1 1i Cu and TiN-coated h1 1 1i Cu. Straight lines representing a best fit to the data were obtained if at least three data points were available for the crystal orientation and strain level. The results in Figs. 1 and 2 clearly show the trend of increasing flow stress with decreasing specimen size and with increasing strain at fixed D. The size-scaling exponent deduced from Fig. 2 is shown in Fig. 3 for all material systems and strain ranges with sufficient data points. No discrimination between h1 0 0i Cu and h1 0 0i Cu on MgO was made because of the limited size range of the h1 0 0i Cu on MgO data. In all cases, the scaling exponent increases with strain. This observation is in good agreement with recent micro-Laue experiments on Au and Ni micropillars [40]. For strain levels up to 0.10 the h1 1 1i Cu samples show a steeper increase in scaling compared to the h1 0 0i Cu samples, indicating that the h1 1 1i direction is the more hardenable crystal direction,

D. Kiener et al. / Acta Materialia 59 (2011) 3825–3840

3829

Fig. 3. The power exponent deduced from the best fit to the data in Fig. 2 as a function of strain for the different samples investigated.

consistent with bulk measurements [29,30]. Moreover, an effect of the TiN top coat is observed for low strain values, where higher hardening compared to uncoated samples is found. This was not seen in the previously reported flow stress values at a strain of 0.10 [6,27,33], where the data for both samples closely merges. 3. Modeling 3.1. Formulation and simulation methods The mechanical response of the Cu specimens is modeled assuming that plane-strain conditions prevail. This approximation is reasonable in pillars with a square cross-section, as it is generally in macroscopic specimens (e.g. [41]). Since the response of rectangular pillars is not fundamentally different from that of cylindrical ones, we conclude that plane-strain analyses may suffice to shed light on the fundamental aspects of micropillar deformation, as far as the states of stress and strain are concerned. The compression simulation setup is sketched in Fig. 4. Deformation takes place in the x1–x2 plane with x1 the direction of compression. Self-similar pillars having width D, height H and a fixed aspect ratio of H/D = 2 are considered. Each pillar has two potential slip systems oriented at ±35.25 from the x1-axis. The boundary conditions are such that surfaces at x2 = ±D/2 are traction free and the shear stress vanishes at x1 = ±H/2. A constant displacement rate u_ 1 ¼ U_ is uniformly prescribed on the top surface x1 = H/2 while the vertical displacement of all nodes of the bottom surface is precluded. To eliminate rigid-body motion, one node is fixed on the bottom surface, as shown in Fig. 4. All other nodes of that surface are thus free to move laterally. This model allows for the rotation of the crystal axis. In tension, the same set-up is used with u_ 1 ¼ U_ .

Fig. 4. Schematic showing the plane strain model of a micropillar of width D and height H, being oriented for symmetric double slip with two slip systems oriented at ±u0 from the x1 axis. The displacement boundary conditions applied allow the rotation of the crystal axis.

The DD formulation follows that of Ref. [20], which specializes to the pillar problem the superposition method in Ref. [35] for solving general small-strain boundary-value problems and the 2.5-D paradigm in Ref. [34] for representing 3-D dislocation interactions in a 2-D setting. In this framework, the dislocations are modeled as Volterra line singularities embedded in a linear elastic isotropic medium. Although not considered here, the effect of elastic inhomogeneities such as the TiN top coat or MgO substrate could be represented following the formulation in Ref. [42]. Consistent with the plane-strain approximation, dislocation loops are modeled as dipoles of edge dislocations of opposite sign and restricted to gliding in their slip plane. Long-range interactions are directly accounted for through elasticity with the Peach–Koehler glide force on dislocation i given by " # X i i j ^þ f ¼m  r r  bi ð2Þ j–i

i

where m is the slip plane normal, bi the Burgers vector with signed length bi and b = jbij, and rj the infinite medium ^ is the image stress obstress field of dislocation j. Also, r tained by solving the complementary problem to impose the boundary conditions specified above. As in previous works [20,43] the finite-element method is used to obtain the image fields at each time increment. The superposition method does not apply to dislocation cores where nonlinear effects are prominent. The latter are represented through a set of rules [34] describing the closerange dislocation interactions. Dislocation mobility, which

3830

D. Kiener et al. / Acta Materialia 59 (2011) 3825–3840

is governed by the dislocation–phonon interaction, is described through: Bvi ¼ f i  a

lb i b Sid

ð3Þ

where vi is the glide velocity, B the drag factor, a the line tension parameter, Sid a signed dipole separation and l the shear modulus. In addition, rules for dislocation annihilation at a critical distance Le and glide out of the crystal are used [34] along with “2.5-D” rules for dislocation junctions resulting from short-range reactions. Junction formation is taken to occur when non-coplanar mobile dislocations fall within a distance d. Based on arguments regarding the importance of cross-slip in stabilizing junctions (see Ref. [34] and references therein), we distinguish two populations. A junction can act as an anchoring point for a dynamic FR source, with probability p. It is postulated that p represents an effective cross-slip probability. Alternatively, the junction acts as a breakable obstacle. The breaking stress for junction I is given by: I fbrk =b ¼ bbrk

lb SI

ð4Þ

where SI is the nearest junction spacing and bbrk is a scaling factor for the junction strength. A dislocation pinned at a junction is released only if the junction is destroyed. The activation of a dynamic FR source follows the model further developed by Benzerga [18]. The critical stress is given by a formula similar to Eq. (4) but with bbrk replaced with bnuc. The critical time for nucleation is approximately given by: tInuc ¼ c

SI jf I j

ð5Þ

where c is a constant having units of a drag factor. As discussed by Benzerga and co-workers [17–19] there are two ways to initialize the dynamics. Initial dislocation sources and obstacles can be represented as sessile dislocations. In the above-cited works, forest dislocations were used as “centers” of multiplication or as static obstacles. However, in the large source-density case, a point representation of sources and obstacles is more adequate and will be followed here. A random distribution of point FR sources and obstacles is considered. A dislocation dipole is nucleated from a source when the Peach–Koehler force acting on it exceeds a critical value snucb for a prescribed time t0n. The source strengths are randomly assigned from a Gaussian distribution with average snuc . A glide dislocation may get pinned at a static obstacle and is released when the Peach–Koehler force at the location of the obstacle attains the obstacle strength sobsbi. This representation has been extensively used in the literature (e.g. [20,43–47]). Within this framework, plastic flow arises due to the collective motion of dislocations. In particular, the stress– strain response, including any emerging strain hardening, is an output of the simulations. Microscopically, junction formation results in DD source and obstacle evolution dur-

ing the deformation process, and this eventually leads to hardening. The average stress and strain are computed from: Z 1 D=2 U ð6Þ r11 ðH =2; x2 Þdx2 ; e ¼  r¼ D D=2 H where r11 is the normal traction. These quantities are directly comparable with experimental measurements, Eq. (1). Hence, the same notation is used. The value of the applied U_ is scaled with specimen size so that the resulting nominal strain rate j_ej ¼ U_ =H is the same for all specimens. In the simulations the dislocation population evolves: the number of dislocations increases, because of initial and dynamic sources, and some level of self-organization takes place. To detect the potential development of a dislocation substructure, the spatial distribution of GNDs is monitored following a procedure first outlined in Ref. [20] and later developed in Ref. [43,48]. In general, both the dislocation density, q, and the GND density, qGND, viewed as continuous fields, are non-uniformly distributed, but the latter is in addition resolution-dependent. To illustrate this, the dislocation density in any volume V, say Z 1 ðV Þ ¼ q qdV V V with q the local density, takes values that only depend on V. In particular, it does not depend on the way in which the density is distributed within V. Thus, it is resolutionindependent. On the other hand, the GND density in V: Z 1 GND ðV Þ ¼ q q dV ð7Þ V V GND depends on both V and how it is locally distributed. This distinction between q and qGND is rooted in the polarity of dislocations and takes its simplest form in single slip for which q = q+ + q and qGND = jq+  qj, q+ and q being the positive and negative dislocation densities, respectively; see Fig. 5. For multiple slip conditions, one has [20]:

(a)

(b)

Fig. 5. Resolution dependence of the GND density, illustrated in single slip in a unit volume. (a) Case of statistical storage: q = 6, qG = 0 and GND ¼ 0. (b) Case of an emergent substructure: q = 6, qG = 0 but q GND –0ð¼ 2Þ. q

D. Kiener et al. / Acta Materialia 59 (2011) 3825–3840

rhffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi i2 hX i2 X ðnÞ ðnÞ ðnÞ Þ cos uðnÞ ðnÞ Þ sin uðnÞ qGND ðxÞ ¼ ðq  q þ ðq  q þ þ   n n for any subdomain x  V. Here, u(n) refers to the orientation of slip system n with respect to the loading axis. In practice, dislocation densities are not continuous fields and the (non-Riemann) integral in Eq. (7) is evaluated as follows. Taking V to be the domain occupied by the specimen, then discretizing V into a grid of joint equi-sized elements Ve, Eq. (7) becomes: GND ¼ q

X Ve q ðV e Þ V GND e

ð9Þ

where qGND(Ve) is calculated from Eq. (8). The resolution GND manifests through an inherent dependependence of q GND ðV ; V e Þ dence upon the size of Ve. Thus we may view q e as a function of both V and V . The two extreme limits: GND ðV ; V e Þ; qG lim q e V !V

GND ðV ; V e Þ q ¼ lim q e V !0

ð10Þ

have clear physical meaning. qG is the net GND density within the specimen, whereas q is the total density in V. For any intermediate resolution, i.e. value of the size of Ve, we choose the term “effective GND density” for GND . This terminology is a parallel to the notion of effecq tive plastic strain in continuum plasticity theories. With respect to our investigation of the size effect, it is expected that qG = 0 under uniaxial loading, irrespective of specimen size V. However, the effective GND density GND does not necessarily vanish. If a dislocation substrucq GND ture emerges, then a strong resolution dependence of q is expected. What is of particular importance is that, at GND during straining fixed resolution, the monitoring of q will allow examination of its dependence upon the specimen size V. In addition, the spatial distribution of GND density at desired resolution can be analyzed by using a grid of fixed element size.

3831

ð8Þ

strength parameter, bbrk, and the anchoring point formation probability p. Parameter calibration proceeds in two steps (Fig. 6a). In the first step, the group fqobs ; sobs ; snuc g is identified so as to obtain a target overall “yield” stress and hardening rate in single slip (dashed line in Fig. 6a). The microtension data on a D = 3 lm specimen was used for that purpose. Further details are provided in Appendix A. In a second step, parameters bbrk and p are chosen so as to obtain a good fit, in the average sense, to the stress–strain response of the largest pillar tested, i.e. for D = Dmax (solid line in Fig. 6a). The rationale behind this procedure is that parameters of the second group do not affect the single-slip response, while parameters of the first group weakly affect the hardening rate under multislip. Thus, the two calibration steps are decoupled to a large extent. One set of parameter values that satisfies the above criteria is listed in Table 2.

(a)

3.2. Choice of parameters Simulation parameters representative of Cu were used with l = 47 GPa, m = 0.34 (Poisson’s ratio), b = 0.255 nm and B = 104 Pa s. An initial density of FR sources qnuc = 20 lm2 was estimated based on transmission electron microscopy (TEM) data [39,49]. Potential slip planes of a given system were equally spaced at 40b, irrespective of specimen size. Key parameters related to the initial source/obstacle population and to the dynamic junction population are inferred from a calibration process. One subset of the experimental data, including the tension data, is used to that end. Then, with all material parameters fixed, the effect of pillar width D is analyzed allowing the model to be assessed against a broader set of experimental data. The parameters to be identified fall under two groups: (i) the initial obstacle density, qobs, the obstacle strength, sobs, the mean nucleation strength, snuc ; and (ii) the junction

(b) Fig. 6. Schematic highlighting the steps followed in the calibration of the parameters entering the simulations: (a) The initial source and obstacle properties ðsnuc ; qobs ; sobs Þ were calibrated to match the apparent yield and hardening observed in the micro-tension experiment on a crystal oriented for single slip. The parameters governing the junction/dynamic obstacle strength (bbrk) and the probability of junctions stabilizing (p) were calibrated to achieve hardening observed in the largest micro-compression experiment on a multiple slip oriented crystal having Dmax. (b) The parameters obtained from the calibration were used to investigate size effects in crystals by reducing the width from Dmax down to Dmin.

3832

D. Kiener et al. / Acta Materialia 59 (2011) 3825–3840

Table 2 Material parameters used in the simulations. Calibrated

qobs (lm2) 20

sobs (MPa) 150

snuc (MPa) 20

bbrk (–) 5

p (–) 0.01

Other

d = Le (nm) 1.5

t0n (ns) 10

a (–) 0.3

bnuc (–) 1

c (Pa s) 101

Additional parameters of atomistic character are assigned values based on estimates from 3-D DD analyses, atomistic calculations or theory [19]. In view of the universality of scaling laws in bulk plasticity, the key trends are not sensitive to particular choices of many such parameters. Their values are also included in Table 2. The applied nominal strain rate is j_ej ¼ 104 s1 in all calculations. A time step of Dt = 0.5 ns was used to resolve the dynamics of dislocation nucleation and motion. When GND densities are discussed, two values of the resolution are used: 50 nm for density maps and 200 nm for effective GND density calculations. The coarser resolution used in the latter aims at ensuring a minimum number of dislocations in the majority of grid elements. For a more systematic analysis of resolution effects the reader is referred to Ref. [20].

in the simulations for the D = 2.08 lm pillar is 112.9 ± 6.8 MPa at a strain of e = 0.10. The corresponding value in the experiment is 137 MPa. In the crystal with width D = 1.08 lm the shear stress averaged over three realizations is 128.2 ± 6.6 MPa and the corresponding value in the experiment is 164.4 MPa at a strain of e = 0.10. Representative stress vs. strain responses for pillar widths in the range D = 0.4–3.2 lm are shown in Fig. 9. A general trend of increasing flow stress and hardening rate with decreasing size is observed. Subsequent to yielding, all specimens exhibit some hardening, with the exception of the D = 0.4 lm pillar, which does not harden much until

45

Micro-tension, D = 3.0 μm

Simulation – Experiment:<234> Cu

40 35

3.3. Simulation results

25 20 15 10 5 0 0.00

(a) 0.02

0.04

0.06

0.08

0.10

ε (-) 140

Micro-compression, D = 5.7 μm

Simulation Experiment: <111> Cu

120 100

τ (MPa)

The results of simulations using the material parameters in Table 2 were compared with experimental data of the calibration set, Fig. 7. This task was undertaken to check that the calibration procedure outlined above delivers the expected results. In single-slip tension, the expected yield stress is sY = 23.6 MPa; see Eq. (A.1). This is in keeping with the results in Fig. 7a. In addition, the value used for sobs delivers a good representation of hardening, which is very weak. Notice the difference in elastic stiffness between experiments and simulations, as the experimental stiffness is affected by the system compliance. In double-slip compression with Dmax = 5.7 lm, the computed and measured hardening rates are in good agreement (Fig. 7b). The difference in the hardening rate in the simulations, averaged over three different realizations of initial source and obstacle distributions, and the experiment is found to be less than 16% of the experimental value. Next, a systematic investigation of size effects is carried out with D values ranging from 0.4 to 9.6 lm. This range contains the values used in the experiments. For each pillar size, at least three realizations, corresponding to a fixed source and obstacle density but different initial distributions of sources and obstacles, were simulated. In Fig. 8a and b the shear stress vs. strain response is shown for pillars with width D = 2.08 lm and D = 1.08 lm, respectively. Good agreement between the simulation and experimental shear stress vs. strain response is observed. Immediately after the onset of yield, the simulations capture the hardening noticeable in the experiments. The shear stress averaged from three realizations

τ (MPa)

30

80 60 40 20 0

(b) 0

0.02

0.04

0.06

0.08

0.1

ε (-) Fig. 7. Comparison of shear stress (s) versus strain (e) response between simulations and experiments in the calibration step: (a) Micro-tension simulation response in comparison to micro-tension experiment on single slip oriented h2 3 4i Cu for a crystal with width D = 3.0 lm. (b) Microcompression simulation response in comparison to micro-compression experiment on h1 1 1i Cu for a crystal with width D = 5.7 lm.

D. Kiener et al. / Acta Materialia 59 (2011) 3825–3840

160 140

180 Experiment: round <100> Cu Simulation

D = 5.7 μm

160

D = 2.08 μm

120

140 Slip-system 1

120

100

-2

ρ (μm )

τ (MPa)

3833

80 60

D = 0.4 μm

100 80

Slip-system 2

60 40

40

20 0 0.00

(a) 0.02

0.04

0.06

0.08

0.10

ε (-) 180 160 140

0 0.00

0.02

0.04

0.06

0.08

0.10

0.12

0.14

ε (-) Fig. 10. Evolution of the dislocation density (q) in the simulations with strain (e) for two slip systems is shown for crystals with width D = 0.4 lm and D = 5.7 lm. A symmetric evolution in the two slip systems is observed for the larger sample, while for the smaller sample distinct differences between the slip systems and notably lower densities compared to the larger sample are observed.

Experiment: round <100> Cu Simulation D = 1.08 μm

120

τ (MPa)

20

100 80 60 40

(b)

20 0 0.00

0.02

0.04

0.06

0.08

0.10

ε (-) Fig. 8. Prediction of shear stress (s) versus strain (e) response from the simulations as compared to experiments for selected samples: (a) Comparison between the circular h1 0 0i Cu and simulation for crystal with width D = 2.08 lm; (b) Comparison between the circular h1 0 0i Cu and simulation for crystal with width D = 1.08 lm.

500

D = 0.4 μm

400

σ (MPa)

0.8 μm

300

200 3.2 μm

100

0 0.00

0.02

0.04

0.06

0.08

1.6 μm

0.10

ε (-) Fig. 9. Representative stress (r) versus imposed strain (e) responses from simulations are shown for crystals with width varying from D = 0.4 lm to 3.2 lm.

a strain of e = 0.02. However, beyond this strain the crystal hardens at a much faster rate than the larger specimens.

This feature is similar to the trend observed in the experiments. The evolution of dislocation density is shown in Fig. 10 for a small pillar (D = 0.4 lm) and a large pillar (D = 5.7 lm). In the D = 5.7 lm pillar the dislocation density builds-up at a rapid rate in both slip systems from the onset of plastic deformation. Within the strain range shown, the density continues to increase and the rate of increase is similar in both slip systems. The early activation of both slip systems precludes the observation of a two stage stress–strain response in the simulation. On the other hand, in the D = 0.4 lm pillar, the dislocation density does not pick up in any of the two slip systems until a strain of

0.02. Beyond this point the rate of increase in the dislocation density is different in the two slip systems. This shows that in smaller pillars the discreteness of source distribution can lead to localization of slip, even in samples that are potentially oriented for multiple slip, in agreement with in situ TEM [50] and 3-D DD [14] observations. In addition, this explains the lack of classical Taylor hardening in the smaller specimen at small strains. Moreover, from Fig. 10 it is observed that the dislocation density is greater in the D = 5.7 lm pillar. Yet, the flow stress is smaller than that reached in the D = 0.4 lm pillar, in the ensemble average sense. This indicates the breakdown of the Taylor hardening law, where the flow stress scales with the square root of the dislocation density. Deformed configurations at a strain of e = 0.10 are shown for pillars of width 0.4, 1.6 and 5.7 lm in Fig. 11. All the specimens show evidence of double slip. However, in the D = 0.4 lm pillar we observe localization of slip on one system, in fact on very few slip planes. This correlates with the dislocation density evolution shown in Fig. 10. In the larger D = 1.6 lm and D = 5.7 lm specimens, slip is

3834

D. Kiener et al. / Acta Materialia 59 (2011) 3825–3840

Fig. 11. Deformed configurations of crystals with width D = 0.4 lm, 1.6 lm and 5.7 lm at a strain level of e = 0.1.

more evenly distributed along the height of the specimens. Additional microstructural information provided by the discrete DD simulations will be discussed in the following section. 4. Comparison of experimental and simulation results 4.1. Flow stress First, comparison is made on the size-scaling of flow strength. Since the simulations were all performed using displacement control, only experiments in displacementcontrol mode have been analyzed. In fact, as shown in Ref. [33], the stress–strain curves of specimens tested under load and displacement control are in agreement. In addition, while the range of pillar diameters in the experiments is 0.9–7 lm, smaller values of D were investigated in the DD analyses. The simulation results and experimental data for flow stress vs. pillar diameter are summarized in Fig. 12. Overall, the simulations capture the tendency for the flow

τ f (MPa)

1000

round <100> Cu, ε=0.05 round <100> Cu, ε=0.1 square <111> Cu, ε=0.05 square <111> Cu, ε=0.1 Simulation: ε=0.05 Simulation: ε=0.1

4.2. Hardening rate

100

10

1

stress to increase with decreasing pillar size, as well as the stochastic character of plastic flow in micropillars, as reflected in the amount of scatter (also see Fig. 2). An inverse size effect is noteworthy in the D = 0.4 lm pillar up to a strain of 0.05. This behavior is related to the scarcity of point dislocation sources. It indicates a limitation of the point-source model in source-limited plasticity. This limitation can be remedied as in Refs. [17–19]. However, as soon as a dislocation structure develops, as is the case at e = 0.1, the experimental trend is retrieved down to D = 0.4 lm. The flow stress power-scaling exponent in the simulations is 0.09 ± 0.04 at e = 0.05 and increases to 0.17 ± 0.03 at e = 0.10. These values are lower than the exponents measured experimentally; but the trend of increasing exponent with increasing strain is consistent with the experiments (compare with Fig. 3). This trend is due to the effect of pillar size on hardening. The large standard deviation observed in the values of the power exponents determined from the simulations reflects the increase in the scatter of the flow stress values in smaller pillars. For example, at a strain of e = 0.10 the shear flow stress in the D = 0.4 lm pillars is found to vary between 238 and 118 MPa.

10

D (μm) Fig. 12. The shear stresses (sf) from the simulations for various crystal widths (D) are shown in comparison to experimental data from circular h1 0 0i Cu and square h1 1 1i Cu at strains of e = 0.05 and 0.1.

The calculation of the hardening rate is based on true stress–strain response. To calculate the hardening rates between two reference true strains (etrue), the values of true stress (rtrue) were calculated under the assumption of homogeneous pillar deformation and volume conservation, as is commonly done [5,7]. The hardening rate is then given true by: H ¼ Dr . In the present case, this procedure is justified Detrue by the investigated multiple slip orientation and the confirming in situ observations. A similar approach was followed in the simulations. The hardening rate data is determined for a true strain range of 0.02–0.1. In those simulations and experiments where the true strain reached is below 0.1, the true strain

D. Kiener et al. / Acta Materialia 59 (2011) 3825–3840 round <100> Cu square <111> Cu square <111> Cu - coated Simulation

Θ (GPa)

10

1

0.1

1

10

D (μm) Fig. 13. Hardening rate (H) data determined from true stress–strain response for crystals of various sizes are shown from simulations and experiments on circular h1 0 0i Cu, square h1 1 1i Cu, and square h1 1 1i Cu coated with TiN. The hardening rate in the simulations and experiments are determined between the true strain range of etrue = 0.02  0.1.

3835

range considered for the determination of hardening is between 0.02 to the maximum true strain reached. Fig. 13 shows the H vs. D data in a logarithmic plot. Data from the pillar with size D = 5.7 lm was used for calibration. For all other pillar sizes, measured and computed hardening rates compare well and the size effect on the hardening is well captured by the simulations. The simulations reveal that with decrease in pillar size there is an increase in the scatter of the hardening value. For the pillar with size D = 0.4 lm hardening is found to be as high as 4.6 GPa and as low as 1.7 GPa. Larger pillars demonstrated less scatter and evidence of forest mechanisms expected under multislip deformation. This is supported by examination of deformed pillars in SEM; see explanation to Fig. 14 below and in situ TEM observations [50]. 4.3. Microstructural features Junction formation results in dynamic dislocation source and obstacle evolution during the deformation

Fig. 14. (a,b) and (f,g) show in situ SEM images of h1 0 0i Cu samples with size D = 6 lm compressed by a flat conical diamond punch. The specimen in the upper row was tested under displacement control (d.c.) to a strain of 0.28, while the sample in the lower row was loaded under load control (l.c.) to a strain of 0.18, respectively. Multiple slip on several slip planes is observed in both cases, where fewer slip planes were activated in the l.c. tested sample. (c,h) Nearest neighbor misorientation maps and (d,j) global misorientation maps with respect to the underformed sample base determined by electron backscatter diffraction (EBSD) with 25 nm step size. Note that different color coding was applied. (e) and (k) show (0 0 1) pole figure maps of the samples with the same color coding as in (d) and (j). The compression axis corresponds to the (0 0 1) projection direction. Fragmentation in the misorientation maps and much stronger peak broadening in the pole figure maps is observed for the upper row sample deformed to a higher compressive strain.

3836

D. Kiener et al. / Acta Materialia 59 (2011) 3825–3840

process. This was also observed in quantitative in situ TEM experiments [50] and 3-D DD simulations [14] and is key to adequate simulation of small-scale hardening. Fig. 14a, b, f and g show images captured during in situ SEM testing of a D = 6 lm pillar in displacement-

controlled (a, b) and load-controlled (f, g) mode. There is clear evidence of the crystals deforming by double slip. Fig. 14c and h show the corresponding nearest-neighbor misorientation maps determined from electron backscatter diffraction (EBSD) scans with 25 nm step size, with the

(a)

(b)

(c)

Fig. 15. Contour plots at a nominal strain e = 0.1 in crystals with width (left to right) D = 0.4 lm, 1.6 lm and 5.7 lm: (a) Axial stress, r11, with superposed dislocation structures; (b) GND density, qGND (at 50 nm resolution); and (c) Lattice rotations, j.

D. Kiener et al. / Acta Materialia 59 (2011) 3825–3840

500

(a)

-2

ρ (μm )

400

300

200

100

0

0

1

2

3

4

5

6

7

D (μm) 1

(b)

0.8 – /ρ ρ GND

color code ranging from 0 to 8. Several highly misoriented regions are observed, especially for the sample in the upper row deformed to a strain of 0.28. Furthermore, the global misorientations with respect to the undeformed sample base increase with increasing maximum strain, as depicted in Fig. 14d and j, respectively. Note that in these images the color code ranges from 0 to 20. The increasing crystal fragmentation observed in the previous images is also reflected in the (1 0 0) pole figure shown in Fig. 14e and k, respectively, where increased peak broadening in multiple directions is seen with increasing compressive strain. The same color code as for the global misorientation maps was applied. Simulated dislocation structures are shown in Fig. 15 along with contours of internal stress fields, GND density and in-plane lattice rotation. At fixed pillar size, the inner region is harder than the surface. This is due to dislocation intersections being more likely in the central region. The formation of stable junctions pins glide dislocations and prevents them from escaping at the free surfaces. In the D = 0.4 lm pillar, dislocation activity is limited to a few slip planes, consistent with the deformed configuration in Fig. 11. This explains why it can take time for hardening to pick up in the smallest pillar (see Fig. 9 and 10). To explain the size effect, we first focus on a large pillar for which microstructural information is available from both experiments (Fig. 14) and simulations (corresponding to D = 5.7 lm in Fig. 15). Contours of GND density, determined from Eq. (8), are shown in Fig. 15b (right). These contours can be thought of as representing the nearest-neighbor misorientations (Fig. 14c and h). The resolution, i.e. the size of x = Ve in Eq. (8) is 50 nm, which is twice the step size of the EBSD maps. Pockets of high GND density are found in the central region. This feature is in keeping with the large number of local misorientation boundaries observed in Fig. 14c at the center of the crystal. A consideration of the lattice rotation field, j, in Fig. 15c reveals domains of large rotations (as high as 3). These domains are not oriented along the slip systems but are the signature of organized dislocation structures. In addition, there is a mismatch in the sign of the rotation fields, which is reminiscent of the fragmentation processes observed in the experiments, as highlighted in Fig. 14d and j. Experimental observations of the type shown in Fig. 14 were not available for smaller pillars. However, based on the good qualitative comparison between observed and simulated microstructural features, we have sought an interpretation of the size effect by examining, in the simulations, the changes that occur to the GND density and j fields when D decreases. It is worth mentioning that the lattice rotation and GND density fields are generally self-similar in different pillars; also see the results in Refs. [48,51] for crystals as large as D 13 lm. However, below a critical D, the bulk-like lattice rotation and GND fields are disturbed by the presence of free surfaces. In Fig. 15 a clear

3837

0.6 0.4 0.2 0 0

1

2

3

4

5

6

7

D (μm) Fig. 16. (a) Total dislocation density (q) in crystals at a strain of e = 0.1 as a function of crystal width D; (b) Effective GND density ( qGND ) normalized by the total dislocation density (q) in crystals at a strain of e = 0.1 as a function of crystal width D. A resolution size of 200 nm is used in the calculation of the effective GND density ( qGND ) in all the samples.

deviation from bulk-like contours of qGND and j is visible when D decreases (right to left in the figure). As explained in Section 3.1, the net GND density qG is expected to vanish in all specimens, irrespective of size. It is worth noting that this has been computationally verified for pillars of all sizes, up to some statistical variations associated with random dislocation escape events. In order to quantify the emergent dislocation substructures, the variations with pillar size of the total dislocation GND , defined in density, q, and the effective GND density, q Eq. (9) are reported in Fig. 16 at e = 0.10. The resolution GND is 200 nm. This choice is motivated used to determine q by a detailed study of the effects of resolution, grid topology and the average number of dislocations per grid element [52]. For the range of dislocation densities reached in the simulations, the 200 nm value is also well within the range of appropriate characteristic lengths suggested by Kysar et al. [53] to ensure the calculation of a meaningful value of GND density. The general trend in Fig. 16 is a slight decrease in q with decreasing D and a more significant drop for some D = 0.4 lm specimens. By way

3838

D. Kiener et al. / Acta Materialia 59 (2011) 3825–3840

of contrast, the mean and standard deviation of the ratio GND =q are both found to increase with decreasing D. Conq GND =q, corsistent with previous work [48], this measure, q relates best with the qualitative trend of smaller being harder. These observations suggest that the emergence of GND density plays an important role in the size effect in micropillars. 5. Discussion The combined micro-compression experiments and DD simulations show a propensity for micropillars to harden and that a size effect is induced by work-hardening under multislip conditions. Experimental and theoretical analytical tools were used to uncover the origin of the observed size effect. Our findings strongly suggest that the size effect is mainly caused by the build-up of GNDs. Relatively large strains and dislocation densities are needed for such behavior to become detectable. Multislip conditions favor this scenario. In the simulations, the emergent dislocation substructure is a result of self-organization at some level. As the specimen size decreases, the process of self-organization is disturbed by the presence of free surfaces, hence the size effect. What is of particular importance is that both longrange and short-range dislocation interactions contribute to the size effect. In particular, the atomic-scale events of junction formation result in dynamic dislocation-source and obstacle evolution during the deformation process, which are the fundamental processes of hardening in pure single crystals. While comparison was made in Section 4 between the experimental and computational results, it is instructive to point out some differences between the methods used. An effort was undertaken in order to identify the DD model parameters based on a subset of experimental data. The predictions from the DD model were then directly validated by comparison to other experimental data. However, some details are worth discussing.

There is a difference in the initial dislocation structure between experiment and simulation. While the model assumes statistically distributed dislocation sources in the volume, the FIB-prepared samples also contain near-surface dislocations from the fabrication process [54]. This FIB damage was shown to drop the theoretical strength of defect-free Mo down to what is observed for prestrained samples [55]. On the other hand, in situ TEM studies on (1 0 0) oriented submicron Cu pillars showed that a significant amount of this FIB damage escapes the sample during early loading [50]. Combined experimental and analytical [54] as well as computational [56] approaches concluded that the strengthening contribution of the FIB damage to the observed size effect is negligible for micron-sized samples. This is especially true in the present study focused on the hardening response and related size effects that develop over large strains. Furthermore, there are differences in slip orientation and the plane-strain condition assumed in the simulation. While for the calibration step the experimental situation was closely mimicked (see Appendix A), two orientations were considered in the experiments but only one was investigated in the simulation. As all orientations are highly symmetrical and the classical work by Diehl [29] and Suzuki et al. [30] showed that athermal hardening rates are to first order not affected by the crystal orientation, the simulation results are representative, and thus comparable to experimental ones. Moreover, while out-of-plane shear components are observed experimentally (Fig. 14), the simulations assume plane strain. This could in principle be improved by 3-D DD simulations. However, those reach strains that are smaller than needed for hardening-induced size effects to emerge [12–16]. Referring again to the classical works [29,30], we assume that variations in the out-of-plane shear strain do not significantly influence the determined hardening rates. Hardening due to dislocation–dislocation interactions, in the classical macroscopic sense of strain hardening, is

Fig. 17. Successive in situ SEM images of a h1 0 0i Cu specimen during compression. Note in (c) the local loss of contact at the center. The inclination of the top surface of the diamond punch is a result of the electron beam scanning from left to right during image acquisition while the sample was compressed.

D. Kiener et al. / Acta Materialia 59 (2011) 3825–3840

discussed based on true stress vs. true strain data. On a (sub)micron scale, an experimental determination of true stress vs. strain curves is desirable and seemingly within reach when using in situ SEM testing approaches. However, as depicted for a h1 0 0i Cu sample in Fig. 17, there are situations where, even during continuous observation of sample deformation, it remains a challenge to determine the actual contact area or the smallest cross-section. These strain localizations at the sample/punch interface should be minimized in microtensile testing [26], allowing a more accurate determination of true stress values. The significance of the above findings is now put in perspective of other recent work. First, the decrease of dislocation density with decreasing pillar size (Fig. 16a) is in contrast with previous 2.5-D DD simulations [20]. Despite that, smaller is harder in both series of investigations. This only strengthens the conclusion that the simulated, and most likely the measured, hardening is not of the Taylor kind. Experimentally, TEM studies can be carried out to observe the dislocation substructure [57], but are demanding. The DD simulations readily allow the quantification of the dislocation substructure evolution during deformation. The decrease in q with decreasing D is in contrast with findings based on TEM measurements in Ni microcrystals [57]. In addition, in the present investigation a continuous increase in q with strain is observed (Fig. 10) as in Ref. [20]. By way of contrast, in Ref. [57] the dislocation density was not found to vary much with the applied strain, which may be due to the single-slip conditions prevailing in their experiments. In addition, the presence of a significant density of dislocations in pillars with D as small as 0.4 lm suggests that the crystals are not “starved” of dislocations. This observation is in agreement with in situ TEM investigations on h1 0 0i Al [49] and h1 0 0i Cu [50]. Both the experimental and simulation results indicate that dislocation starvation is not necessary for a size effect to emerge. The interpretation of the size effect on hardening in the GND micron regime based on an effective GND density q bears some resemblance to the classical interpretation of the indentation size effect by Nix and Gao [58], where the GND density scales with the inverse indentation depth and the square root of the flow stress. While this analogy offers a simple interpretation of the pillar size effect in the investigated size range (>0.4 lm), a difference is noteworGND is not the net GND density thy. We emphasize that q in the pillar, which was denoted qG in Section 3. Under nominally uniform compression, qG 0. We check that qG/q = 0.3% in the D = 5.7 lm sample but larger values are found in some D = 0.4 lm specimens. Because of this GND to some macroessential difference, one cannot relate q GND is scopic strain gradient as in the Nix–Gao picture. q thought to relate to some evolving microdeformation field. Finally, the origin of the size effect suggested by the present findings may not be the only one. Below a critical pillar size, and at a given initial dislocation density, the

3839

behavior may become nucleation-controlled [50]. The micropillars of the current study are still in a size range where nucleation is not the only dominant mechanism controlling the flow strength. A further reduction in sample size may shift this situation to a dislocation nucleation-controlled state and accordingly source-controlled hardening [50], since dislocation interaction on different glide planes becomes less probable. For example, recent results indicate a breakdown of the Nix–Gao model in nanoindentation of very shallow indents, where dislocation nucleation effects become the controlling mechanism [59,60]. 6. Conclusions A quantitative comparison has been carried out between experimental measurements and DD simulations for the size dependence of plastic flow in Cu micropillars. Key parameters entering the DD model were calibrated based on a subset of the experimental data. Then with all DD parameters fixed, the DD simulations capture the experimental trends and are in good quantitative agreement with the experimental measurements. Specifically, the main findings are:  A size effect on the flow stress is confirmed and the sizescaling exponent increases with further straining. The size effect is not consistent with conventional Taylor hardening.  Although more challenging to quantify, the hardening rate increases with decreasing pillar size and thus explains the variation of the size-scaling exponent.  Evidence from both experiments and simulations strongly suggests that the size effect on hardening is caused by the increasing amount of stored GNDs in smaller pillars. Based on these observations, it can be concluded that the size effect observed in work hardening of micropillars does not result from either dislocation starvation or strain gradient plasticity effects. It originates from the deformation induced development and evolution of dislocation substructures that give rise to a size-dependent build-up of GNDs in the pillars. Acknowledgments Support from the National Science Foundation through the Faculty Early Career Development Program (CMMI0748187) and from the Texas A&M University Supercomputing Facility are greatly acknowledged. Appendix A. Calibration Step 1 In the DD simulations, the flow stress at incipient plasticity and small-strain hardening are governed by the average initial dislocation source strength ðsnuc Þ, the standard deviation of the source strength distribution (Rnuc), the

3840

D. Kiener et al. / Acta Materialia 59 (2011) 3825–3840

obstacle density (qobs) and obstacle strength (sobs). These four parameters are calibrated using the h 2 3 4i microtension data [27]. Since the standard deviation only affects the scatter, we set Rnuc ¼ 0:25snuc . In addition, the value sobs = 150 MPa is chosen so as to match the hardening rate in the above experiment. The number of parameters to be calibrated is thus reduced from four to two, snuc and qobs. Chakravarthy and Curtin [47] developed an analytical formula relating the yield stress sY to the above parameters through: sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi Lobs lb sobs ðA:1Þ þ s2nuc sY ¼ m  Lobs pð1  mÞ Lobs The yield stress sY in Eq. (A.1) is defined as the stress at which dislocations nucleated from sources achieve flow past obstacles in their path. The average obstacle spacing Lobs is related to the obstacle density as Lobs = 1/(qobsd) with d the slip plane spacing. It is the weakest of the sources which nucleates first and hence the active sources are typically from the lower set of the source strength distribution; thus, snuc ¼ snuc  2Rnuc . The ratio Lobs =Lobs is a material-independent parameter based on the statistical considerations of the obstacle spacing, having a value of 6.7. The value of the numerical factor m for obstacles randomly distributed around the sources is m = 4.5 [47]. Using Eq. (A.1), trial estimates of parameters snuc and qobs are obtained. The parameters are then used in the simulation of a D = 3.0 lm crystal oriented for single slip with u0 = 28.8. This slip configuration has a Schmid factor of fs = 0.422, which is the same as for the h 2 3 4i specimen. From the trial estimates the values of the parameters are iteratively modified until a good fit is obtained between simulated and measured yield stresses. References

[15] [16] [17] [18] [19] [20] [21] [22] [23] [24] [25] [26] [27] [28] [29] [30] [31] [32] [33] [34] [35] [36] [37] [38] [39] [40] [41] [42] [43] [44] [45] [46]

[1] Argon AS. Strengthening mechanisms in crystal plasticity. Oxford: Oxford University Press; 2008. [2] Uchic MD, Dimiduk DM, Florando JN, Nix WD. Science 2004;305:986. [3] Dimiduk DM, Uchic MD, Parthasarathy TA. Acta Mater 2005;53:4065. [4] Greer JR, Oliver WC, Nix WD. Acta Mater 2005;53:1821. [5] Volkert CA, Lilleodden ET. Phil Mag 2006;86:5567. [6] Kiener D, Motz C, Scho¨berl T, Jenko M, Dehm G. Adv Eng Mater 2006;8:1119. [7] Frick CP, Clark BG, Orso S, Schneider AS, Arzt E. Mater Sci Eng: A 2008;489:319. [8] Kiener D, Motz C, Dehm G. J Mater Sci 2008;43:2503. [9] Uchic MD, Shade PA, Dimiduk DM. Annu Rev Mater Res 2009;39:361. [10] Horstemeyer MF, Baskes MI, Plimpton SJ. Acta Mater 2001;49:4363. [11] Deng C, Sansoz F. Nano Letters 2009;9:1517. [12] Rao SI, Dimiduk DM, Parthasarathy TA, Uchic MD, Tang M, Woodward C. Acta Mater 2008;56:3245. [13] El-Awady JA, Wen M, Ghoniem NM. J Mech Phys Solids 2009;57:32. [14] Motz C, Weygand D, Senger J, Gumbsch P. Acta Mater 2009;57:1744.

[47] [48] [49] [50] [51] [52] [53] [54] [55] [56] [57] [58] [59] [60]

Akarapu S, Zbib HM, Bahr DF. Int J Plasticity 2010;26:239. Tang H, Schwarz KW, Espinosa HD. Acta Mater 2007;55:1607. Benzerga AA, Shaver NF. Scripta Mater 2006;54:1937. Benzerga AA. Int J Plast 2008;24:1128. Benzerga AA. J Mech Phys Solids 2009;57:1459. Guruprasad PJ, Benzerga AA. J Mech Phys Solids 2008;56:132. Hommel M, Kraft O, Arzt E. J Mater Res 1999;14:2373. Huang HB, Spaepen F. Acta Mater 2000;48:3261. Xiang Y, Vlassak JJ. Acta Mater 2006;54:5449. Gruber PA, Bohm J, Onuseit F, Wanner A, Spolenak R, Arzt E. Acta Mater 2008;56:2318. Brenner SS. J Appl Phys 1957;28:1023. Kiener D, Grosinger W, Dehm G, Pippan R. Acta Mater 2008;56:580. Kiener D, Grosinger W, Dehm G. Scripta Mater 2009;60:148. Richter G, Hillerich K, Gianola DS, Monig R, Kraft O, Volkert CA. Nano Letters 2009;9:3048. Diehl J. Z Metallkd 1956;47:331. Suzuki H, Ikeda S, Takeuchi S. J Phys Soc Jpn 1956;11:382. Argon AS, Brydges WT. Phil Mag 1968;18:817. Prinz F, Argon AS. Phys Status Solidi 1980;57:741. Kiener D, Motz C, Dehm G. Mater Sci Eng: A 2009;505:79. Benzerga AA, Bre´chet Y, Needleman A, Van der Giessen E. Model Simul Mater Sci Eng 2004;12:159. Van der Giessen E, Needleman A. Model Simul Mater Sci Eng 1995;3:689. Kiener D, Motz C, Dehm G, Pippan R. Int J Mater Res 2009;100:1074. Purswani JM, Spila TDG. Thin Solid Films 2006;515:1166. Zhang H, Schuster BE, Wie Q, Ramesh KT. Scripta Mater 2006;54:181. Dehm G, Worgotter HP, Cazottes S, Purswani JM, Gall D, Mitterer C, et al. Thin Solid Films 2009;518:1517. Maaß R, Van Petegem S, Ma DC, Zimmermann J, Grolimund D, Roters F, et al. Acta Mater 2009;57:5996. Benzerga AA, Tvergaard T, Needleman A. Int J Fract 2002;116:275. Keralavarma SM, Benzerga AA. Model Simul Mater Sci Eng 2007;15:239. Guruprasad PJ, Benzerga AA. Acta Mater 2008;56:5477. Cleveringa HHM, Van der Giessen E, Needleman A. Acta Mater 1997;45:3163. Deshpande VS, Needleman A, Van der Giessen E. J Mech Phys Solids 2003;51:2057. Nicola L, Bower AF, Kim KS, Needleman A, Van der Giessen E. J Mech Phys Solids 2007;55:1120. Chakravarthy SS, Curtin WA. J Mech Phys Solids 2010;58:625. Guruprasad PJ, Benzerga AA. Phil Mag 2008;88:3585. Oh SH, Legros M, Kiener D, Dehm G. Nature Mater 2009;8:95. Kiener D, Minor AM. Acta Mater 2011;59:1328. Guruprasad PJ, Benzerga AA. IOP conference series: materials science and engineering, vol. 3; 2009. p. 012008. Guruprasad PJ. Scale effects in crystal plasticity. Ph.D. thesis, Texas A&M University, USA; 2010. Kysar JW, Saito Y, Oztop MS, Lee D, Huh WT. Int J Plast 2010;26:1097. Kiener D, Motz C, Rester M, Jenko M, Dehm D. Mater Sci Eng: A 2007;459:262. Shim S, Bei H, Miller MK, Pharr GM, George EP. Acta Mater 2009;57:3895. El-Awady JA, Woodward C, Dimiduk DM, Ghoniem NM. Phys Rev B 2009;80:110410. Norfleet DM, Dimiduk DM, Polasik SJ, Uchic MD, Mills MJ. Acta Mater 2008;56:2988. Nix WD, Gao H. J Mech Phys Solids 1998;46:411. Rester M, Motz C, Pippan R. Phil Mag Lett 2008;88:879. McLaughlin KK, Clegg WJ. J Phys D 2008;41:074007.