Elastic forces drive nonequilibrium pattern formation in a model of nanocrystal ion exchange.

Chemical transformations, such as ion exchange, are commonly employed to modify nanocrystal compositions. Yet the mechanisms of these transformations, which often operate far from equilibrium and entail mixing diverse chemical species, remain poorly understood. Here we explore an idealized model for ion exchange in which a chemical potential drives compositional defects to accumulate at a crystal's surface. These impurities subsequently diffuse inward. We find that the nature of interactions between sites in a compositionally impure crystal strongly impacts exchange trajectories. In particular, elastic deformations which accompany lattice-mismatched species promote spatially modulated patterns in the composition. These same patterns can be produced at equilibrium in core/shell nanocrystals, whose structure mimics transient motifs observed in nonequilibrium trajectories. Moreover, the core of such nanocrystals undergoes a phase transition-from modulated to unstructured-as the thickness or stiffness of the shell is decreased. Our results help explain the varied patterns observed in heterostructured nanocrystals produced by ion exchange and suggest principles for the rational design of compositionally patterned nanomaterials.

Methods of chemical transformation are routinely used to alter the properties of nanocrystals postsynthesis (1). Among these, ion exchange—the replacement of one ion species by another in a crystal (2, 3)—is a common and effective technique for modifying nanoparticle composition. Cation exchange in particular has been used to produce a variety of heterostructured nanocrystals, chiefly metal chalcogenides, whose mixed compositions exhibit diverse morphologies (4, 5). Despite their utility, the mechanisms which govern the progress of cation exchange reactions remain unclear. Thus, cation exchange has likely not attained its full potential as a method for precisely tuning the composition and spatial organization of species within nanocrystals.

Experiments probing cation exchange reactions have yielded several key observations. Kinetic measurements of exchange among certain cation species suggest distinct behaviors on different time scales: rapid change in the composition at short times, presumably due to the fast introduction of guest ions at the nanocrystal surface, followed by much slower change at long times, thought to reflect internal diffusion of cations (5–8). Meanwhile, product morphologies vary widely, depending on the chemical identity of the species involved and on reaction conditions like stoichiometry (3, 4, 9). In cases where exchange is only partial, the resulting heterostructures can depend strongly on the relative sizes of the cations. Exchange of species with significant lattice mismatch tends to produce graded, alloy-like structures (as in Zn/CdSe) (5) or spatially modulated patterns with length scales of several lattice spacings (as in Ag2/CdS) (10–12), while similarly sized ions instead adopt core/shell or Janus morphologies (4). However, the limited spatial and temporal resolution of current experimental techniques has frustrated attempts to build a conceptual framework with which to understand these intriguing behaviors. For example, it is still not entirely clear whether ion exchange proceeds primarily via vacancies or interstitial defects (6, 13, 14). More broadly, the relative contributions of thermodynamics and kinetics in shaping the (potentially metastable) final products of cation exchange have yet to be elucidated in detail.

In contrast to experiments, computer simulations can fully resolve atomic motions but are limited to short time and length scales (especially for chemically detailed models). Two previous simulation studies, which rendered specific materials with a fairly high level of detail, have yielded some interesting mechanistic insights into cation exchange. Ott et al. (15) performed kinetic Monte Carlo (KMC) simulations of CdS →Ag2S exchange, parametrized based on density functional theory calculations of defect energies in the bulk crystal. They found evidence for cooperativity among charged defects and identified this as a key driver of exchange. Fan et al. (16) instead performed molecular dynamics simulations of PbS →CdS exchange, employing pseudoligands as a coarse-grained representation of solvent and ligands which typically decorate nanocrystals. They found that pseudoligands promote extraction of Pb and that interstitial defects mediate internal impurity transport, at least at the elevated temperatures necessary to observe significant exchange in the simulations. However, in both cases the authors were only able to monitor reactions on time scales far short of complete exchange. Moreover, neither study explicitly addressed lattice mismatch, which can produce long-ranged elastic interactions (11, 17, 18).

Here we report computer simulations of a lattice model for ion exchange in nanocrystals. We have in mind cation exchange reactions of metal chalcogenides, but the model does not specify a particular compound. Instead, we consider a minimal model which focuses on the key experimental features outlined above: differing time scales for surface exchange and internal diffusion, and lattice mismatch. Although our approach sacrifices a significant degree of microscopic chemical detail, the resulting reduction in computational expense carries benefits that are unprecedented in numerical studies of ion exchange. First, we are able to establish equilibrium properties of the microscopic model not just in reactant and product states but also across the entire range of intermediate compositions, providing a key reference for understanding reaction dynamics. Second, we can propagate exchange trajectories all the way from reactants to products, without heating the system unrealistically or artificially biasing its evolution. Finally, we can thoroughly sample ensembles of nonequilibrium exchange trajectories, enabling conclusions and comparisons that are statistically robust.

Model

Our model resolves compositional change on a scale that is microscopic but not atomic. We associate a binary variable with each site R of a finite lattice, roughly equivalent to dividing a nanocrystal into N unit cells that have either exchanged or not. thus indicates local majority of the incoming cation (labeled B and colored red in figures), while indicates majority of the outgoing species (labeled A and colored blue). In this way, we represent crystal structure at the scale of its simplest repeating units (a few ), slightly coarser than in previous simulation studies.

Exchange is driven by a difference between the chemical potentials of A and B cations, i.e., an effective energy . This thermodynamic bias could represent a difference between ion concentrations in a nanocrystal’s solution-phase environment, a difference in their solubilities (possibly mediated by dissolved ligands), and/or a difference in their lattice energies (19). In typical experiments, conditions for ion exchange are strongly favorable, so we set the magnitude of to be much larger than the thermal energy . We found that exchange trajectories were insensitive to the precise value of for .

Trajectories of this model advance through discrete, reversible events in which the state of one or more lattice sites changes. At a site on the nanocrystal’s boundary, one cation type can be exchanged for another (), with a rate that depends on the system’s compositional configuration . The nanocrystal’s interior evolves through diffusive steps that swap the states of adjacent lattice sites (, where and are nearest neighbors), with configuration-dependent rate . Such a swap may represent a series of microscopic barrier-crossing events, perhaps involving transient vacancies or interstitials (6, 13, 14); we resolve only the net transport of ion density. See Fig. 1 for illustrations of these events. The rates and are formulated to satisfy detailed balance with respect to the equilibrium probability distribution of the nanocrystal’s compositional state. Because they depend on configuration , their values generally change as the reaction proceeds. Similar to previous simulations by Ott et al. (15), these stochastic rate processes are numerically realized with a KMC algorithm ().

Fig. 1.

(A) Schematics illustrating exchange at the surface of a nanocrystal (Top) and diffusion within a nanocrystal (Bottom). Note that and depend on the compositional state of the nanocrystal before and after the corresponding event; the value of is therefore different for the two surface exchange events shown. (B) Representative configurations taken along KMC exchange trajectories for a statically uncorrelated model (noninteracting); a model with nearest-neighbor interactions (Ising); and a model with the nonlocal, mechanically mediated interactions of Eq. (elastic).

Ion Exchange Dynamics

If governed only by a difference in chemical potential (and not by configuration-dependent energies), the model exchange dynamics would amount to simple diffusion in a domain bounded by a sink. A trajectory exemplifying this noninteracting case is shown in Fig. 1 for a hexagonal lattice in two dimensions (2D). A corresponding plot of the net composition over time is shown in Fig. 2. With , A →B exchange at the perimeter is extremely facile and nearly irreversible. Mirroring experimental observations (5–8), the nanocrystal’s outermost layer is transformed immediately. Interior swap moves, which do not change the system’s net composition , leave the energy unchanged in this simple case. Neither facilitated nor hindered, diffusion proceeds slowly, transporting B inward and A outward with a high degree of randomness at this length scale.

Fig. 2.

Net composition versus time averaged over 100 KMC trajectories of (A) the noninteracting model with , (B) the Ising model with and , and (C) the elastic model with and . In each case, N = 271. Inset in A shows the same data with composition on a logarithmic scale, emphasizing the approximately exponential decay of the composition. Insets in B and C show the same data with time on a logarithmic scale. The steplike decay of composition apparent in B, Inset, reflects monolayer-by-monolayer advancement of the B-rich shell into the crystal interior. Adapted with permission from ref. 20.

Cation exchange dynamics are made interesting by interactions, which correlate the composition of different sites, in space and in time, and can generate surprising collective behaviors. We consider two kinds of site–site interactions that have figured importantly in studies of alloy thermodynamics (21). One is spatially local,where the sum is restricted to nearest neighbors on the lattice and J is a constant. This Ising-like interaction describes short-ranged, material-specific chemical preferences for sites of the same type (J > 0) or of opposite type (J < 0), and its influence on the time dependence of phase change has been investigated extensively (21–24). The other energetic contribution is mechanical in origin, accounting for the elastic stress inherent to mixing solid materials with different lattice constants (25–27). We adopt a simple and general description of these elastic forces, in which neighboring sites of the lattice prefer a bond length that depends on the identities of their occupants (17). Integrating out fluctuations in lattice geometry yields an effective energy for the composition variables

where ϵ is a positive constant proportional to the Young’s modulus (17). The effective elastic interaction potential , whose form is detailed in , decays gradually and nonmonotonically with distance (Fig. 3). While is minimum for unmixed states (pure A or pure B), at intermediate compositions it can favor spatial patterns that are richly modulated (Fig. 4) (17). Its consequences for the dynamics of compositional phase change have, to our knowledge, not been thoroughly explored, especially for nanoscale systems.

Fig. 3.

Effective elastic interaction potential as a function of for two different values of R for an N = 271 hexagon-shaped nanocrystal with a triangular lattice structure. The site R is colored white and is outlined in black. The color of each other site indicates the value of . R is located (Left) at the center of the nanocrystal and (Right) on an edge. In both cases, the magnitude of the elastic potential is significant even when R and are separated by several lattice spacings. Adapted with permission from ref. 18.

Fig. 4.

(A) Representative configurations taken from bulk simulations of the elastic model at equilibrium (17). Across a range of net compositions c, unstructured (U, j = 1, 2) and superlattice (S, j = 1, 2) phases, as well as coexistence between these phases (U + S, S1 + S2), can be observed. (B) Seven-site plaquettes associated with different phases, together with corresponding values of the order parameter s and the color assigned to the central site of each plaquette. All other plaquettes are unclassified and assigned s = 0. (C) Representative configurations along the course of a KMC exchange trajectory of the elastic model shown in the original color scheme (Top) and in the plaquette color scheme (Bottom). Perimeter atoms were omitted in the plaquette color scheme renderings. Configurations were taken at times t 2.98 (Left), 25.6 (Middle), and 74.8 (Right) in units of (setting units of mass and distance to unity). (D) Fractional population of each phase α over time in KMC trajectories of the elastic model with N = 721. Solid-color lines represent individual trajectories; the thick lines with black outlines represent averages over 100 independent trajectories.

An effort to mimic a specific real material would include both Ising and elastic contributions, but here we consider their impact separately (either J = 0 or ϵ = 0). Our exploration of the combined case (both J and ϵ nonzero) suggests a similar range of behaviors. In both scenarios we set k to be small compared to characteristic energy scales, consistent with the high strength of covalent bonds in crystalline solids compared to thermal energy at room temperature. Ising interactions, for J > 0 and well below the critical temperature, greatly impede the overall kinetics of exchange, as illustrated in the trajectory shown in Figs. 1 and 2. In this case, rapid exchange at the perimeter is followed by much slower advance of a sharp A/B interface toward the center of the nanocrystal. Since A species can only be removed from the nanocrystal at its boundary, this process requires that A atoms enter and cross the B-rich shell as high-energy impurities. The costly creation of these defects is rate-limiting at the temperatures of interest and is especially hindered at smooth regions of the A/B interface where many favorable interactions must be disrupted simultaneously. Experiments on lead chalcogenide nanocrystals in which Pb is exchanged for Cd show strong kinetic trapping after exchange of the first few monolayers (14, 16, 28, 29), akin to our simulated trajectories and consistent with the bulk immiscibility of PbS and CdS (30). This exchange scenario is strongly foreshadowed by equilibrium thermodynamics of the Ising model at subcritical conditions (31). Whether in bulk or in a nanocrystal, nearly pure phases dominate here, and the line tension of the A/B interface is high, ensuring strong segregation and high barriers to impurity transport.

The relationship between thermodynamics and exchange kinetics is more subtle for the elastic model. In previous work (17, 18) we have broadly examined the equilibrium behavior of this model (reviewed in Fig. 4), which features microscopically patterned superlattice phases (denoted S1 and S2), in addition to unstructured phases (U1 and U2) with Ising-like symmetry breaking. Stabilities of these phases, and states of coexistence between them, are unusually sensitive to boundary conditions. In particular, the spatially modulated phases S1 and S2 dominate at intermediate composition and low temperature for periodic bulk systems but are unstable when boundaries are free to deform heterogeneously. At equilibrium the nanocrystals studied in this paper, due to their free boundaries, exhibit superlattice structure with a probability that is negligibly small, vanishing in the limit of large crystal size (18, 20). In the nanocrystal’s nonequilibrium exchange trajectories, however, superlattice patterns are pronounced. As shown in Fig. 1, motifs of S1 and S2 emerge shortly after the perimeter transforms and dominate the crystal’s interior as exchange continues. These trajectories appear not to involve U1/U2 interfaces (i.e., A/B interfaces), and they proceed to completion much more rapidly than for the case of Ising interactions (Fig. 2 ).

In order to quantify transient superlattice patterns in nanocrystal KMC simulations, we introduce an order parameter s that locally distinguishes among the four phases of the elastic model (). As illustrated in Fig. 4, the instantaneous value of s at each lattice site R is determined by the compositional arrangement of the seven-site plaquette comprising R and its nearest neighbors. The proportion of phase can then be judged by the fraction of sites with the corresponding value of s. This designation also allows a vivid rendering of KMC trajectories, with each site colored according to its current value of s (rather than the value of ). Fig. 4 shows an exchange trajectory visualized in this way, emphasizing the extent and shape of each phase’s domain as the reaction proceeds. The clearly evident decay of U1, through S1 and S2, to U2 is demonstrated quantitatively in Fig. 4 by plotting as a function of time for a large number of exchange trajectories. The proportion of U1 rapidly decays to nearly zero, and the populations of U2, S1, and S2 increase accordingly. A much slower decay of superlattice phases follows, accompanied by terminal growth of to 1. Although transient, modulated spatial order is thus long-lived during cation exchange of our simulated nanocrystals.

Why do superlattice phases characteristic of the bulk elastic model appear in our out-of-equilibrium nanocrystals, despite their thermodynamic instability under these boundary conditions? We argue that the rapid dynamics of exchange at the perimeter—the first step in all our exchange trajectories—produces an effective, transient change in boundary conditions experienced by the crystal’s interior. The stiffness of the exchanged shell mimics the influence of a bulk periodic environment, which favors modulated order. An adiabatic exchange dynamics, in which the nanocrystal adjusts completely after each small change in net composition, would begin in a similar fashion, accommodating impurities at the perimeter. However, in that reversible case, continued accumulation of B atoms would result in not a shell but instead compact domains of U1 and U2 whose interface spans the nanocrystal. Such shell disintegration is prohibitively slow on the time scale of nonequilibrium exchange at large , so that our irreversible trajectories follow a qualitatively different route.

Equilibrium Phases of Core/Shell Nanocrystals

We scrutinize this shell hypothesis by using constrained equilibrium Monte Carlo (MC) simulations (see ) to compute equilibrium averages (denoted as ) of when the properties of a nanocrystal’s shell are systematically varied. We first consider a multilayer shell consisting of B sites surrounding a core containing sites with fixed net composition c = 1/2 (i.e., half A and half B). A and B atoms may swap positions within the core, but the identity of the shell sites is held fixed. Fig. 5 show configurations representative of this constrained equilibrium for several shell sizes, along with in the nanocrystal core as a function of shell thickness. In the absence of a shell, the core is predominantly composed of U1 and U2 phases. As layers are added to the shell, regions of S1 and S2 gradually become larger, and those of U1 and U2 become smaller. At six layers, is roughly equal for all α. When an additional layer is added, U1 and U2 phases abruptly vanish, and S1 and S2 phases dominate, hinting at a phase transition. Demonstrating the existence of a phase transition in the thermodynamic limit is complicated by our inability to vary the shell thickness continuously in this model. As an alternative, we vary the stiffness of a shell with fixed size, which we choose here to be a single atomic layer. (In practice this is accomplished by changing the spring constant between bonds connecting shell sites; see for more details.) The resulting plots of versus (Fig. 6) exhibit the same general behavior, but the transition point can be located with much higher precision. The fraction of U1, U2 drops precipitously (and that of S1, S2 sharply rises) at for . An underlying phase transition would manifest as a peak in the related fluctuation quantity (where ) whose height grows with system size, diverging in the thermodynamic limit. Fig. 6 shows simulation results consistent with this scenario. The transient superlattices in our nonequilibrium simulations thus appear to be stabilized by a nearby equilibrium phase transition.

Fig. 5.

(A) Configurations taken from equilibrium MC simulations of a nanocrystal whose shell composition is fixed at c = 0, while the core fluctuates with net composition c = 1/2. Results are shown for and for three different shell thicknesses. Interior core sites are colored according to the plaquette scheme; the outermost core sites and the shell sites are translucent and colored according to the original scheme. Note the coexistence between U1 and U2 for zero shell layers; between U1, S1, S2, and U2 for six shell layers; and between S1 and S2 for eight shell layers. (B) Average fractional population of the different phases versus shell thickness, measured in atomic layers, for a core size of . Note the crossover between predominantly unstructured coexistence and predominantly superlattice coexistence that occurs around six shell layers.

Fig. 6.

(A) Average fractional population of the different phases as a function of shell stiffness for a single-layer shell and core size of . Note the sharp change in occurring at , consistent with a first-order phase transition. (B) Variance of the fraction f of unstructured phases (U1, U2) versus for different system sizes (here N refers to the size of the core). Peaks become narrower and higher with increasing system size, consistent with a first-order transition. Peak positions shift to higher as system size increases due to the decreasing ratio of shell size to core size.

The observed relationship between shell stiffness and modulated order can be explained from an understanding of the stability of superlattice phases in periodic bulk systems at equilibrium. In a fixed volume with regular shape, the coexistence of solid phases with different lattice parameters carries a cost that is extensive in system size (absent topological defects in lattice connectivity) (17, 18). The unstructured phases U1 and U2 of the elastic model in 2D have maximally different values of the area per lattice site A, so their coexistence is considerably disfavored. The superlattice phases S1 and S2 are more similar to one another in this respect and hence can dominate at intermediate composition. A nanocrystal with free boundaries can undergo collective elastic distortions which relax this energetic cost, stabilizing the coexistence of unstructured phases. Adding a large or stiff shell to the nanocrystal hinders such deformations and thus favors superlattice phases.

In Fig. 7 we view nanocrystal exchange trajectories through this mechanical lens, focusing on the local area per site at each position R (see for more details). Fig. 7 shows the deviation of this quantity away from its spatial average , for a few configurations taken from a single KMC trajectory; its instantaneous variance is plotted in Fig. 7. The spatial distribution is its most heterogeneous at short times: the shell is highly compressed compared to the core, reflecting the dominance of B sites at the perimeter of the crystal and A sites in the interior. As time progresses, this mechanical heterogeneity becomes less dramatic, as reflected in the decay of . Additionally, we observe distinct regions with characteristic values of . These regions align precisely with domains of the different compositional phases α. The growth of S1 and S2 phases (which have areas intermediate that of U1 and U2) relaxes strain at the core/shell interface, yielding the observed decrease in . The corresponding constrained equilibrium average decreases with increasing shell thickness (Fig. 7), reflecting how the shell hampers spatially varying elastic deformation of the core. Comparison of Figs. 5 and 7 shows that also closely tracks changes in . Local area distortions are thus strongly tied to local composition.

Fig. 7.

(A) Configurations taken from a KMC exchange trajectory of the elastic model with sites colored according to local area deviation. Configurations were taken at times t 2.98 (Left), 25.6 (Middle), and 74.8 (Right). Compressed regions at the exterior correspond to B-rich regions, while expanded areas in the core correspond to A-rich regions; different shades of color track regions populated by different phases. (B) rms deviations of (spatially averaged) versus time. Gray lines represent individual KMC trajectories (with N = 721), while the black line represents the average over 100 trajectories. (C) rms deviations (averaged over both space and sampled configurations) from equilibrium MC simulations of fixed-shell nanocrystals (with ) versus shell thickness. Note the drop between six and seven shell layers, corresponding to the transition from a state in which all phases coexist (at six layers) to a state of two-phase coexistence involving S1 and S2 only (at seven layers.)

Stability of Transient Modulated Structures

If , the driving force for ion exchange, is maintained indefinitely in our KMC simulations, then pure-A nanocrystals will inevitably transform into pure-B nanocrystals with only transient modulated structure. If, however, the driving force is removed before exchange is complete, such transient patterns can serve as an initial conditions for long-lived composition modulation. To demonstrate this, we have performed KMC simulations in which exchange moves are disallowed after a certain time, modeling relaxation of incompletely transformed nanocrystals. In Fig. 8, we show a series of snapshots from one such trajectory. The images shown are typical of relaxation trajectories initialized from configurations in which U1, S1, S2, and U2 domains are all present. As expected from equilibrium considerations, unstructured U1 and U2 domains do grow at the expense of modulated S1 and S2 ones. Nevertheless, regions of modulated order persist within the nanocrystal interior. Their occurrence primarily between flanking unstructured domains at opposing sides of the nanocrystal is consistent with the thermodynamic stability of S1 and S2 phases in the presence of a stiff nanocrystal shell, as discussed in Equilibrium Phases of Core/Shell Nanocrystals and demonstrated in Figs. 5 and 6.

Fig. 8.

Configurations along a relaxation trajectory of the elastic model, in which exchange moves were turned off after a time . Configurations are shown for times t –t 0 (first panel), 125 (second panel), 193 (third panel), and 221 (fourth panel).

Results in 3D

All the results we have presented thus far have been for 2D systems (which have the advantage of being easily visualized). Nanocrystals in the laboratory are of course mostly 3D entities. We therefore also performed simulations of a 3D version of our elastic model (Fig. 9). The results broadly agree with those in 2D: rapid surface exchange is followed by slower internal rearrangements (Fig. 9) characterized by spatially modulated patterns (Fig. 9). These patterns are more complicated than in 2D. Rather than characterizing them in detail, here we have simply identified unstructured regions (U1, U2) and labeled everything else as “unclassified.” As we show in Fig. 9, these unclassified regions dominate over unstructured phases at equilibrium in the presence of a sufficiently stiff shell. Thus, despite differences in the precise compositional morphologies that elastic forces produce, the impact of lattice mismatch on model ion exchange trajectories is similar in 2D and 3D.

Fig. 9.

(A) Configurations taken from an elastic model exchange trajectory of a 3D nanocrystal, whose atoms reside on a simple cubic lattice (with both nearest- and next-nearest-neighbor bonds). Here we set and . Configurations were taken at times t 2.58 × 10−26 (first panel), 1.12 × 10−4 (second panel), 0.0741 (third panel), and 0.579 (fourth panel). A atoms are colored blue; B atoms are colored red and are rendered as translucent for clarity. Note the spatially modulated patterns (columns of blue atoms) that develop over the course of time. (B) Fractional population of U1 and U2 phases versus time. The unclassified population here includes contributions from spatially modulated phases as well as the compositionally disordered phase. Thin, solid-color lines represent individual KMC trajectories with ; thick lines outlined in black represent averages over 100 trajectories. (Inset)The short-time behavior on a logarithmic scale. (C) Average fractional population of phases in equilibrium MC simulations of fixed-shell nanocrystals, as a function of (single-layer) shell stiffness for a core size and a fixed core composition c = 1/2. The predominantly A and predominantly B phases (U1 and U2) are present for small values of but are absent for large values.

Modulated compositional order has been realized experimentally in CdS nanorods via partial cation exchange, which can yield alternating stripes of Ag- and Cd-rich domains (10, 32). The equilibrium state of this mixed-cation system likely involves a single Ag-rich domain, separated from a Cd-rich domain by a single interface (11), yet stripes persist over times at least as long as experimental observation. Density functional theory calculations suggest that such patterns may be stabilized by lattice mismatch (11), whose elastic consequences are the essential ingredient of our model. We tested whether similar patterns could be obtained in our computational approach by changing the aspect ratio of 3D nanocrystals. In Fig. 10, we show the results of KMC simulations of an elongated nanocrystal. Fig. 10 show configurations taken from ion exchange and relaxation trajectories, respectively (the latter involving only swap moves, as described in Stability of Transient Modulated Structures). Fig. 10 show corresponding plots of c(z), the composition along the long (z) axis of the nanorod averaged over the x and y dimensions. Ion exchange dynamics in our elastic model is seen to generate transient concentration modulation along the long axis of the nanorod. Subsequent relaxation of these modulated configurations in the absence of a driving force for exchange can produce long-lived stripes of alternating composition, closely resembling those observed in experiments. In , we show that similar stripes arise at equilibrium in core/shell nanorods.

Fig. 10.

(A) Configurations taken from an elastic model exchange trajectory of an elongated 3D nanocrystal with dimensions . As in Fig. 9, we set and . Configurations are shown for times t 0.0724 (first panel), 7.23 (second panel), 11.8 (third panel), and 13.5 (fourth panel). (B) Composition profile along z, averaged over the x and y dimensions, for several times along the trajectory of A. The vertical line at signifies a uniform composition of c = 1 (all sites blue) at time t = 0. (C) Configurations taken from a relaxation trajectory of a partially ion-exchanged nanorod. The leftmost configuration, used to initialize the relaxation trajectory, was taken at time from the exchange trajectory of A and B. The relaxation portion of this trajectory consisted only of diffusion moves that do not change the net composition. Configurations were taken at times t –t (first panel), 1.93 (second panel), 9.54 (third panel), and 11.3 (fourth panel). (D) Composition profile along z, averaged over the x and y dimensions, for several times along the trajectory of C.

Discussion

Our model, while fairly simple, has proven able to reproduce some key features of cation exchange reactions. We have shown how a strongly favorable driving force for exchange creates an effective nonequilibrium boundary condition that shapes transport of material between the solution and the nanocrystal interior. Nanocrystal morphologies during reaction—observed to vary widely in experiments depending on cation identity—result in our models from the interplay between such nonequilibrium boundary conditions and the nature of interactions between the exchanging species. Elastic interactions, in particular, lead to rich spatial patterns in the composition, akin to modulated patterns observed in experiments (10, 11). Although our work was inspired by cation exchange within nanocrystals, our model could also apply to nanoparticle superlattices; in that context one could envision sites corresponding to individual nanocrystals, the length of whose DNA-mediated bonds could be tuned by pH, for example (33). Our results would then suggest altering the properties of the perimeter of such assemblies as a way of biasing spatial organization of the interior. Additionally, our finding that sufficiently stiff or large shells stabilize modulated phases in the nanocrystal core may explain the recent observation of labyrinth structures in a related elastic model of spin-crossover nanoparticles (34). Nanocrystal superlattice and spin-crossover systems, which are often extended in only 2D, may serve as ideal settings in which to test the detailed predictions of our 2D elastic model simulations.

As a model of cation exchange in particular, our model clearly has shortcomings. It does not resolve anion and cation sublattices and hence omits potentially important electrostatic interactions (15). Additionally, it has fixed lattice connectivity and hence cannot describe dislocations or account explicitly for interstitials or vacancies. One could imagine remedying this in part by introducing defect states (potentially carrying a net charge) in addition to the A and B states or by allowing bonds between adjacent sites to break or form. Importantly, our model entirely ignores the ligands which passivate semiconductor nanocrystals and whose density on different facets could affect the rate of surface exchange. A more elaborate model might account for their effect implicitly—perhaps by modifying the rates of surface exchange in the KMC scheme based on the facet and ligand type—or might represent them explicitly as fluctuating species which modify by their interactions with surface sites. Despite these limitations, the success of our model in capturing the basic phenomenology of nanocrystal compositional change can serve as a foundation for future, more microscopically detailed investigations.

Materials and Methods

Energy Function.

The total energy E(C) of a configuration is given by

In simulations of the Ising model, we set the parameter ϵ in the elastic energy to zero; similarly, in simulations of the elastic model, we set the parameter J in the Ising energy to zero.

In simulating the noninteracting model, we set . For our Ising model simulations, we set and . This temperature is well below the triangular-lattice Ising model critical temperature (35). Simulations of the 2D elastic model used and . For simulations of the 3D elastic model, we used and . The choice of temperatures for elastic model simulations would place our systems in the low-temperature regions of the corresponding bulk phase diagrams, where modulated phases are thermodynamically stable (17, 20).

Elastic Interaction Potential.

The effective elastic potential is determined by integrating out the mechanical fluctuations of an elastic model which resolves atomic motions. To be concrete, we take the small-mismatch limit of that model’s Hamiltonian (17, 18):where are bond vectors which may vary by site (surface sites have fewer bonds than interior sites), , a is the lattice parameter, and (where and l and l are the rest lengths of A–A and B–B bonds, respectively). The displacement variables measure deviations of atoms from their perfect lattice sites R. The energy scale depends on the spring constant K of harmonic bonds between sites, as well as the lattice mismatch . Integrating over the harmonic mechanical fluctuations is equivalent to minimizing with respect to and . The former simply contributes an additive constant to the free energy, which we ignore. The latter yields minimum-energy displacements for a given configuration of composition variables,

The effective energy function that results from evaluating at those displacements iswhere the matrices D, C, and S are given by

We identify an effective pair potential as

, which is sensitive to the nanocrystal’s boundaries, was evaluated and tabulated for all site pairs in advance of each simulation.

Form of Rate Constants.

The rate of an event depends on the initial and final configurations C and through the energy function E. We chose rate constantsthat are consistent with the Boltzmann distribution, satisfying as a condition of detailed balance (and similarly for ). The bare rate constants and set the fundamental time scales for exchange and diffusion. The mechanism of ion exchange depends only on their ratio, which we set as for simplicity. Surface dynamics are nonetheless much faster than diffusion for this case, as suggested by experiments (5–8), due to the large chemical potential difference driving exchange. Accordingly, our results are qualitatively insensitive to the value of this ratio so long as .

KMC Algorithm.

KMC trajectories were generated using the Gillespie algorithm (36). Specifically, given an initial configuration, we determine all possible exchange and diffusion moves and compute their rate constants. Denoting the rate constant of event i (out of n total) as k, we compute . We then generate a number r uniformly at random from the interval and execute the event j for which . Finally, we update the time , where is sampled from the exponential distribution . Repeating this process many times yields a single trajectory.

Order Parameter.

Consider a site R with nearest neighbors , and denote the mutual nearest neighbors between R and N as . The order parameter s associated with R is defined as

Calculation of Local Area Deviations.

In the small-mismatch limit to which we restrict ourselves, the minimum-energy displacements are uniquely determined by the compositions according to Eq. . From the displacements, we can measure local deviations from the average crystal area (in 2D) due to heterogeneous composition configurations. Specifically, we compute the local area deviation associated with site R asexpressed in units of a2. Here Z is the coordination number of the lattice, which is 6 for all calculations presented.

Equilibrium MC Simulations of Fixed-Shell Nanoparticles.

Equilibrium configurations were generated at fixed net composition using Kawasaki dynamics (37). Proposed moves consist of attempts to swap the identities of two randomly selected sites R, within the nanocrystal core. Such moves are accepted with probability (38)

A single MC sweep consists of attempted moves.

In each equilibrium fixed-shell simulation, core sites were initialized with random identities consistent with a net core composition of c = 1/2. Systems were then equilibrated by running 105 MC sweeps without collecting data. After the equilibration period, configurations were recorded every sweep for 105 sweeps. Reported observables were averaged over these configurations.

Stiff-shell nanocrystals were modeled by changing the spring constant of the bonds connecting shell atoms; explicitly, terms in Eq. corresponding to shell–shell bonds were multiplied by a factor before the energy minimization described previously to obtain .

To mitigate potentially large surface contributions to quantities like , equilibrium observables were spatially averaged only over interior core sites—the outermost core atoms were excluded.