Surface  energy is a fundamental material property that determines important functions such as catalytic, sensing, and imaging properties. Over the past century, various experimental studies and models including the broken bond theory and Wulff construction have been developed to analyze surface  free energies. However, it remains a challenge to measure or predict thermal fluctuation effects on surface  energies. In particular, crystals of functionalized building blocks, such as self-assembling proteins and DNA-functionalized nanoparticles, assembled via the specific surface interactions of the building blocks, are highly sensitive to thermal fluctuations. In the case of DNA-functionalized nanoparticles, it has been shown that the crystals are formed as a result of thermally active hybridizations. We show here that the surface  energy along different planes can be obtained from the ratio of hybridization events. The surface  energy fluctuations in these systems are shown to bear a nearly linear correlation with the fluctuations in DNA hybridization events in the bulk. We further demonstrate that short DNA chains and high DNA loading increase the volume density of the DNA sticky ends. The relationship between thermally active hybridizations and surface  energy found here can be used to aid the design of single crystals of functionalized colloids with active surface groups.

surface  energy has been a topic of extensive study since the seminal work by Gibbs.1 As one of the most fundamental properties of a solid, surface  energy determines the equilibrium shape of a crystal and plays an essential role in faceting and roughening, as well as in crystal growth phenomena.2 The fluctuations of surface  energy, while important for growing a stable and sharp single crystal,2,3 have seldom been reported or studied. This is because the accuracy and precision of most experimental methods fall short producing relative errors ranging from 5% to 30%.4 Unfortunately, this difficulty only increases when trying to measure thermal fluctuations at a surface. On the theory side, density functional theory has helped to establish a database of surface  energies,5 but energy fluctuations have not yet been reported using such an approach.

With the success of the first single crystal  growth composed of DNA-functionalized nanoparticles (DNA-NP), these man-made “programmable atom equivalents” provide a robust system for growing faceted single crystals in soft materials.6 Over the past two decades, this elegant strategy of using complementary DNA linkers to guide the assembly of nanoparticles has increased our ability to control crystal structure, resulting in a transition from disorganized amorphous products7,8 to organized crystalline symmetries,9–13 as well as micron-scale single crystals with fewer defects.6 Excitingly, the DNA-NP system provides us with a unique opportunity to explore the origin of surface  energy fluctuations, by probing the thermodynamics relevant to DNA chains via molecular dynamics (MD) simulations. Using DNA hybridizations to mimic atomic-scale chemical bonds, one can track the details of a single DNA attachment and detachment event.14–18 Using MD simulations, one is able to capture properties associated with thermal fluctuations, providing a deeper understanding of the atomic crystal systems as well.

It is widely accepted that surfaces bearing large thermal fluctuations do not benefit the stability of single crystal polyhedrons.2,19,20 Therefore, where do these surface  energy fluctuations come from and how can we design a system that minimizes their effect whilst optimizing crystal growth? In this work, we apply a molecular model with robust MD simulations to extract a mathematical relationship bridging design parameters (e.g., length of an individual DNA chain, DNA coverage on each NP, NP size) and surface  energy fluctuations. Based on this relationship, we discuss optimal design parameters that minimize these fluctuations for the design of single crystal polyhedra.

Herein we used a scale-accurate coarse-grained model with explicit DNA chains, which faithfully captures relevant sizes and stiffnesses for the double-stranded (dsDNA) and the single-stranded DNA parts in the design. The attraction between complementary DNA linkers (or sticky ends) is modeled by the Lennard-Jones potential. This model excellently represents DNA hybridization kinetics and the energy level per DNA hybridization, and it has been reported to be robust and reliable in accurately recovering the self-assembly processes of various symmetries with DNA-NPs.15,16

The quantity of interest is the surface excess free energy per unit area of a specific crystal facet γ. By custom, γ is referred to as the “surface energy.” A direct measurement of γ requires to separate a bulk crystal into two parts along a particular plane. In the modeled system, periodic boundary condition along the z-axis was removed from the bulk crystal to expose the facet of interest on two sides (Figure 1). To calculate γ of the exposed facet, we measured the energy difference between the bulk system and the system with two exposed surfaces, then divided the value by twice the surface area. To keep within the capability of molecular modeling, several assumptions were made regarding entropic contributions. In particular, a perfect crystalline structure and constant lattice parameter are assumed.6 That is, we are neglecting the entropic contributions from the lattice vibrations. γ is computed taking the time average of the potential energy.

FIG. 1.

Simulation snapshots showing the setup for the surface energy calculation for a body-centered cubic DNA–NP superlattice. (a) Bulk system with periodic boundary conditions along the x, y, z directions. (b) System with two (100) surfaces exposed, top and bottom. γ ( 100 ) = E surface exposed E bulk / 2 × Area 100 . DNA hybridization events in these snapshots were circled in red. NP cores with different DNA linker types were colored differently for clarity purpose. All snapshots were generated with the visual molecular dynamics package21 and rendered with Tachyon ray-tracer.22 

FIG. 1.

Simulation snapshots showing the setup for the surface energy calculation for a body-centered cubic DNA–NP superlattice. (a) Bulk system with periodic boundary conditions along the x, y, z directions. (b) System with two (100) surfaces exposed, top and bottom. γ ( 100 ) = E surface exposed E bulk / 2 × Area 100 . DNA hybridization events in these snapshots were circled in red. NP cores with different DNA linker types were colored differently for clarity purpose. All snapshots were generated with the visual molecular dynamics package21 and rendered with Tachyon ray-tracer.22 

Close modal

The Wulff construction is an elegant method to determine the equilibrium crystal shape of fixed volume, known as the Wulff shape,23,24 if the orientational dependence of the surface  free energy is known. For each plane h k l , plotting a point in the h k l direction at a distance which is equal to the associated surface  energy γ(hkl) gives the Wulff plot. Planes perpendicular to the radial direction are drawn on each point of the Wulff plot. The inner envelope of all these planes gives the equilibrium shape, or the Wulff shape. With surface  energy values estimated from MD simulations, we apply here the WulffMaker software25 to calculate and display the corresponding Wulff polyhedron.

Wulff construction indicates that, instead of the absolute surface  energy values, it is the surface  energy ratios of different facets that determines the equilibrium shape of a crystal. For example, the body-centered cubic (BCC) metal should exhibit a ratio of γ(110):γ(100):γ(111) = 1:1.41:1.2226,27 corresponding to a rhombic dodecahedron; the face-centered cubic (FCC) metal should exhibit γ(111):γ(100):γ(110) = 1:1.15:1.22, which corresponds to a truncated octahedron. Naturally, surface  energy fluctuations propagate to the error bars of this ratio. In this section, we derive a correlation between the fluctuations on surface  energy ratios and fluctuations of DNA hybridization events with several assumptions made. We performed simulations on distinct systems to verify the model assumptions and to validate the obtained mathematical relationships. We then bridge the fluctuations of DNA hybridizations with design parameters directly controllable in experiments. As such, we found the relationship between design parameters and surface  energy ratios.

Due to the expandable nature of spherical nucleic acid,28 building blocks only interact with their first nearest neighbors12 thus low-index planes such as (100), (110), and (111) become cusp orientations in the Wulff plot,26 determining the equilibrium shape. In general, the surface  energy of any facet (hkl) in a DNA-NP superlattice is calculated by Eq. (1),

γ ( h k l ) = E h k l E bulk 2 × Area h k l = p H ( h k l ) p H ( bulk ) × N DNA pairs × E hybridization 2 × Area h k l = p H ( h k l ) ± Δ p H ( h k l ) p H ( bulk ) ± Δ p H ( bulk ) × η h k l

and

η h k l = N DNA pairs × E hybridization 2 × Area h k l ,
(1)

where Ehkl and Ebulk are the energies of a system with two (hkl) surfaces exposed and the energy of the bulk system, respectively. pH is the percentage of hybridized DNA pairs, NDNA pairs is the maximum number of DNA hybridizations that could be formed, and Ehybridization is the hybridization energy per linker attachment (e.g., Ehybridization = − 42.3 k J/mol for DNA linker sequence −TTCCTT16). p H is the ensemble average of pH and ΔpH is the standard deviation of pH defined by Δ p H = p H 2 p H 2 1 / 2 .

By propagation of uncertainty (e.g., for f = c ( A + B ) , Δ f = c Δ A 2 + Δ B 2 , if A and B are uncorrelated.), Δ γ ( h k l ) = η h k l × Δ p H 2 ( h k l ) + Δ p H 2 ( bulk ) as the two measurements of pH are independent from each other.

As the bulk system and the (hkl) system have a similar system size, we can approximate that Δ p H 2 ( h k l ) Δ p H 2 ( bulk ) Δ p H 2 . We further define k h k l = p H ( h k l ) p H ( bulk ) , which is a mean field value consistent with the classical broken bond model.26,27 It is thus reasonable to consider khkl to be identical for spherical nanoparticles with different design parameters, as long as the crystalline systems bear the same dimension (i.e., the number of unit cells on each dimension of the system box keeps the same). As such, the relative surface  energy fluctuation of facet (hkl) can be approximated as

Δ γ ( h k l ) γ ( h k l ) 2 k h k l 1 × Δ p H ( bulk ) p H ( bulk ) .
(2)

Thus, surface  energy fluctuations originate from the fluctuations of DNA hybridization events. This is consistent with the thermally active hybridization crystallization process, which has been shown to be crucial for the successful crystallization of DNA-NPs.16 

surface  energy ratio is given respective to the reference facet (ref.) which has minimum surface  energy of all facets. For example, facet (110), the closest-packed plane in BCC, is used as the reference for BCC symmetry; similarly, facet (111) is the reference for FCC symmetry. surface  energy ratio regarding any facet (hkl) is expressed as γ ( h k l ) γ ( ref. ) , here after refer to as ghkl.

The uncertainty in ghkl is estimated by

Δ g h k l g h k l = Δ γ ( h k l ) γ ( h k l ) 2 + Δ γ ( ref. ) γ ( ref. ) 2 C h k l × Δ p H ( bulk ) p H ( bulk )

and

C h k l 2 × ( 1 k h k l 1 ) 2 + ( 1 k ref. 1 ) 2 .
(3)

(Propagation of uncertainty: for f = A B , Δ f f = Δ A 2 A 2 + Δ B 2 B 2 , if A and B are uncorrelated.)

With the above derivation, the relative error in ghkl is seen to bear an approximately linear correlation with the relative error in pH(bulk) with coefficient Chkl.

For the sake of comparability and simplicity, all of the systems tested below were constructed into the BCC symmetry. Nevertheless, the protocol we illustrate and the phenomena we discuss in this work is applicable to any crystalline symmetry in the DNA-NP system.

As shown in Figure 2, different BCC lattice were constructed to study different facets of interest, and all are composed of 96 building blocks. Since we have noted different values of k h k l = p H ( h k l ) p H ( bulk ) for different system geometries, we provide an estimation for k100, k110, and k111, as well as the coefficients C100 and C111 for the specific geometries study here. To better evaluate the model, these theoretically estimated parameters are validated in the section titled Validation by simulations by straightforward simulations.

FIG. 2.

Views of different BCC lattice constructions with 96 building blocks exposing surface (100), (110), and(111) on the z-direction. NPs with different DNA linker types were colored differently for clarity purpose.

FIG. 2.

Views of different BCC lattice constructions with 96 building blocks exposing surface (100), (110), and(111) on the z-direction. NPs with different DNA linker types were colored differently for clarity purpose.

Close modal

Consistent with the broken bond model,26,27 a number of building block connections are cut off when the bulk material is cleaved to expose a specific facet. As each building block’s nearest neighbors consist of eight building blocks with the complementary DNA linker type in the BCC lattice, we approximated that each neighbor contributes one eighth of its total connections. Note that building blocks do not connect with their second nearest neighbors in this DNA-NP design. We counted the number of neighbors cut off to estimate the coefficients as in Eq. (4). Parameters of 1/4, 1/2, and 3/8 are simply the portions of neighbors got disconnected by exposing the specific surface plane,

k 110 = p H ( 110 ) p H ( bulk ) = 96 24 × 1 4 × 2 96 = 84 96 , k 100 = p H ( 100 ) p H ( bulk ) = 96 16 × 1 2 × 2 96 = 80 96 , k 111 = p H ( 111 ) p H ( bulk ) = 96 16 × 3 8 × 2 96 = 84 96 .
(4)

Using Eq. (3), the coefficient for the BCC (100) system is calculated to be C100 = 14.14, and for (111) system it is C111 =16.00.

Until now the mathematical relationship (Eq. (3)) has been derived under several assumptions. In order to verify this formula, we simulate five distinct DNA-NP designs with different design parameters and directly calculate their surface  energy values. For each DNA-NP design, we carry out MD simulation runs on approximately five independently generated initial configurations with identical design parameters (DNA chains are randomly attached on each NP). Simulation statistics are summarized in Table I. The error bars in the table are calculated using the mean of the standard deviations of the five identical simulation runs.

TABLE I.

Surface energy values for BCC systems with DNA-NPs.

Design parametersa
No. NP size (nm) DNA coat dsDNA (bp) ρb (nm−3) pHc (%) ΔpH/pHc
10  60  18  0.017  37.84  0.0156 
10  60  58  0.0014  13.23  0.0368 
10  90  18  0.025  41.56  0.0110 
10  90  58  0.0021  16.84  0.0254 
20  150  18  0.023  32.90  0.0105 
Design parametersa
No. NP size (nm) DNA coat dsDNA (bp) ρb (nm−3) pHc (%) ΔpH/pHc
10  60  18  0.017  37.84  0.0156 
10  60  58  0.0014  13.23  0.0368 
10  90  18  0.025  41.56  0.0110 
10  90  58  0.0021  16.84  0.0254 
20  150  18  0.023  32.90  0.0105 
γ(100)d γ(110)d γ(111)d Ratios of surface energiesd
No. (mJ m−2) (mJ m−2) (mJ m−2) γ(110) : γ(100) : γ(111)
0.531 ± 0.077  0.354 ± 0.073  0.440 ± 0.089  1:1.50 ± 0.38:1.24 ± 0.36 
0.057 ± 0.018  0.039 ± 0.017  0.049 ± 0.021  1:1.45 ± 0.79:1.25 ± 0.77 
0.853 ± 0.082  0.573 ± 0.078  0.716 ± 0.095  1:1.49 ± 0.25:1.25 ± 0.24 
0.103 ± 0.022  0.070 ± 0.021  0.087 ± 0.026  1:1.47 ± 0.55:1.24 ± 0.53 
0.628 ± 0.058  0.434 ± 0.055  0.531 ± 0.068  1:1.45 ± 0.23:1.22 ± 0.22 
γ(100)d γ(110)d γ(111)d Ratios of surface energiesd
No. (mJ m−2) (mJ m−2) (mJ m−2) γ(110) : γ(100) : γ(111)
0.531 ± 0.077  0.354 ± 0.073  0.440 ± 0.089  1:1.50 ± 0.38:1.24 ± 0.36 
0.057 ± 0.018  0.039 ± 0.017  0.049 ± 0.021  1:1.45 ± 0.79:1.25 ± 0.77 
0.853 ± 0.082  0.573 ± 0.078  0.716 ± 0.095  1:1.49 ± 0.25:1.25 ± 0.24 
0.103 ± 0.022  0.070 ± 0.021  0.087 ± 0.026  1:1.47 ± 0.55:1.24 ± 0.53 
0.628 ± 0.058  0.434 ± 0.055  0.531 ± 0.068  1:1.45 ± 0.23:1.22 ± 0.22 
a

For example, each building block consists of a nanoparticle sized at 10 nm, coated by 60 DNA chains, and each DNA chain contains 18 base pairs (bp) for the dsDNA part.

b

Volume density of DNA linkers, introduced in Eq. (5).

c

pH and ΔpH/pH refer to the value in bulk system, calculated from the DNA hybridization potential energy.

d

Surface energy values correspond to the 96 building block system discussed in the text. According to the central limit theorem, standard deviation decreases with sample size as σ 1 N , where σ is the standard deviation and N is the number of building blocks in the simulation system.

Note that the surface  energy values shown in Table I correspond to the 96 building block systems in Figure 2. According to the central limit theorem, standard deviation decreases with the number of building blocks: Δ g h k l g h k l 1 N , where N is the number of building blocks in the simulation system. Although N = 96 was the largest system we could simulate due to memory budget, the critical size for a single crystal to retain a stable structure is N ∼ 10 000.19 Structural changes were reported for clusters below a few thousand atoms. Importantly, we should note that in the range of N ≈ 10 000, the fluctuation levels of the surface  energy ratios are about one tenth of the values listed in the table.

Regression of Δ g h k l g h k l on Δ p H p H fit excellently into a linear model as shown in Figure 3, and the coefficients highly agree with the derived values. This statistical analysis has verified the linear correlation between the relative errors in the surface  energy ratios and the relative fluctuation of DNA hybridizations.

FIG. 3.

Linear regression of the relative error in the ratio of two surface energies Δghkl/ghkl on the fluctuation of DNA hybridizations ΔpH/pH. Surface (110) is used as the reference state as it is the plane with lowest surface energy for BCC system. Data points were obtained from five sets of design parameters in Table I. Red dots represent surface (100) and blue dots represent surface (111). The fitted model is Δ g 100 g 100 = 14 . 32 × Δ p H p H + 0 . 01 with the coefficient of determination R2 = 0.9974 for surface (100) and Δ g 111 g 111 = 16 . 30 × Δ p H p H + 0 . 02 with R2 = 0.9973 for surface (111).

FIG. 3.

Linear regression of the relative error in the ratio of two surface energies Δghkl/ghkl on the fluctuation of DNA hybridizations ΔpH/pH. Surface (110) is used as the reference state as it is the plane with lowest surface energy for BCC system. Data points were obtained from five sets of design parameters in Table I. Red dots represent surface (100) and blue dots represent surface (111). The fitted model is Δ g 100 g 100 = 14 . 32 × Δ p H p H + 0 . 01 with the coefficient of determination R2 = 0.9974 for surface (100) and Δ g 111 g 111 = 16 . 30 × Δ p H p H + 0 . 02 with R2 = 0.9973 for surface (111).

Close modal

We have confirmed that the minimization of DNA hybridization fluctuations is essential to minimize the errors in surface  energy ratios. We next predict DNA hybridization fluctuations on design parameters controllable in experiments. Accordingly, bulk systems with numerous sets of design parameters are simulated, including NP sizes at 5, 10, 15 nm, DNA chains composed of 18, 38, 58, 78 base pairs and various levels of DNA coats.

Based on tests with several variables, we find that the volume density of DNA linkers ρ (Eq. (5)) is an excellent indicator to predict ΔpH/pH,

ρ = r V layer ,
(5)

where Vlayer = VmaxVbb, V max = 4 3 π R max 3 , and V bb = 4 3 π R 3 .

In Eq. (5), R represents the hydrodynamic radius of a DNA-NP building block,12Rmax represents the maximum radius that a fully stretched building block can reach12 and r is the number of DNA chains coated on each nanoparticle. Correspondingly, Vbb is the hydrodynamic volume of a DNA-NP building block, Vmax is the maximum volume a building block can reach with DNA chains fully stretched out, and Vlayer is the estimated volume where sticky ends appear.

As shown in Figure 4, ρ and ΔpH/pH fits into a power law model log ( r × Δ p H / p H ) = 0 . 35 × log ( ρ ) 1 . 54 , with the coefficient of determination R2 = 0.983. Taking the cue from this analysis, high volume density of sticky ends is preferred for small errors in surface  energy ratios. From Eq. (5), high DNA coat, short DNA chains as well as large NP if DNA loading density retain constant will benefit a high volume density of sticky ends. Interestingly, experiments fail to grow single crystal Wulff shapes with 20-nm NPs coated by long DNA chains (e.g., 58 or 98 dsDNA pairs).6 

FIG. 4.

Log-log graph showing the power law fitting y(x) = kxn of y = r × Δ p H / p H on x = ρ, where r is the DNA coat on each nanoparticle, ΔpH/pH is DNA hybridization fluctuation, and ρ is the volume density of sticky ends. The fitted model is log ( r × Δ p H / p H ) = 0 . 35 × log ( ρ ) 1 . 54 , with the coefficient of determination R2 = 0.983.

FIG. 4.

Log-log graph showing the power law fitting y(x) = kxn of y = r × Δ p H / p H on x = ρ, where r is the DNA coat on each nanoparticle, ΔpH/pH is DNA hybridization fluctuation, and ρ is the volume density of sticky ends. The fitted model is log ( r × Δ p H / p H ) = 0 . 35 × log ( ρ ) 1 . 54 , with the coefficient of determination R2 = 0.983.

Close modal

We applied the WulffMaker25 to display Wulff polyhedra for surface  energy ratios with different fluctuation levels, as shown in Figure 5. Ideally, the Wulff shape for a BCC metal is the rhombic dodecahedron. We note that in Figure 5, we are not constructing Wulff shapes with fluctuations in surface  energy since Wulff shapes are constructed from ratios of mean surface  energies. In the hypothetical case that the surface  energies lower bounds could be accessible experimentally by growing Wulff shapes while suppressing the fluctuations (say with different solvent conditions and/or grown at low temperatures), and then the Wulff shapes are allowed to relax say by annealing them, the relaxed Wulff shape will become more globular.

FIG. 5.

Left: Wulff polyhedra constructed for surface energy ratios with different fluctuation levels for BCC system, which we denote here as degree of error to quantify the uncertainty on the surface energy ratio since fluctuations are not included in Wulff polyhedra construction. Only polyhedra constructed with the lower bounds of the ratios are displayed. Polyhedra constructed with the upper bounds give the ideal rhombic dodecahedra. This figure is simply to show the importance of the uncertainty on the surface energy ratios if and only if the Wulff shapes could be constructed when the fluctuations are suppressed (say by growing them via a rapid quenching to lower temperatures in different solvent conditions) and then the resulting Wulff shapes are heated up. The Wulff shapes for DNA functionalized NPs are today constructed by slowly cooling; therefore, these shapes are not realizable in the laboratory today.

FIG. 5.

Left: Wulff polyhedra constructed for surface energy ratios with different fluctuation levels for BCC system, which we denote here as degree of error to quantify the uncertainty on the surface energy ratio since fluctuations are not included in Wulff polyhedra construction. Only polyhedra constructed with the lower bounds of the ratios are displayed. Polyhedra constructed with the upper bounds give the ideal rhombic dodecahedra. This figure is simply to show the importance of the uncertainty on the surface energy ratios if and only if the Wulff shapes could be constructed when the fluctuations are suppressed (say by growing them via a rapid quenching to lower temperatures in different solvent conditions) and then the resulting Wulff shapes are heated up. The Wulff shapes for DNA functionalized NPs are today constructed by slowly cooling; therefore, these shapes are not realizable in the laboratory today.

Close modal

The MD simulations were performed with the HOOMD-Blue package29,30 (available at http://codeblue.umich.edu/hoomd-blue) on graphics cards. Constant volume and temperature ensemble (NVT) was applied with the temperature controlled via a Nosé-Hoover thermostat.31,32 All the simulations were performed at the same temperature, providing a universal reference state. Additional details of the model and the simulation protocol can be found elsewhere.6,15,16 The molecular model and the simulation protocol have been shown to capture the observed crystal symmetries (BCC, FCC, CsCl, AlB2, Cr3Si, Cs6C60).15,16 In the present study, since a crystal composed of ∼100 nanoparticles may not stay stable, the positions of the nanoparticle centers are fixed while performing MD simulations yet the particles are free to rotate.6 To keep the system within the memory limitation of a Tesla K40 graphics card, 96 independently generated DNA-NPs were simulated for each design with DNA chains randomly attached on each NP. Nonetheless, the system is larger than the critical nuclei size which is reported to be ∼40 NPs for the DNA-NP system.33 

This work provides a method for understanding and predicting surface  energy fluctuations in DNA-NP systems that originate from thermally active hybridization events, while the latter is found to be easily predicted using the volume density of DNA sticky ends. Indeed, these simulation results are the first quantitative measurements that demonstrate the relationship between thermal fluctuations at surfaces and the thermal fluctuations of DNA hybridization events. We found that shorter DNA chains and higher loading of DNA chains on NPs are the main factors to reduce the errors in surface  energy ratios. Our studies show that the larger the number of particles in a single crystal, the less surface  energy fluctuation errors. Therefore, finding ways to grow larger crystals is essential to develop single crystals of nanoparticles with specific optical and electronic properties. Furthermore, we provide a general method to understand surface  energy fluctuations. This method can be adapted to describe the directed-assembly of anisotropic particle shapes11,34 and colloids covered by multiple types of ligands.35 It would be interesting to determine how these fluctuations could be introduced in other level of coarse-grained models such as those in Refs. 36 and 37.

M.O.d.l.C. and T.I.N.G.L. acknowledge support from the Air Force Office of Scientific Research (AFOSR) Award Nos. FA9550-11-1-0275 and FA9550-10-1-0167. T.I.N.G.L. gratefully acknowledges the Ryan Fellowship and the Northwestern University International Institute for Nanotechnology.

1.
J. W.
Gibbs
,
Trans. Conn. Acad. Arts Sci.
2
,
382
(
1873
).
3.
D. C.
Schlößer
,
L. K.
Verheij
,
G.
Rosenfeld
, and
G.
Comsa
,
Phys. Rev. Lett.
82
,
3843
(
1999
).
4.
V. K.
Kumikov
and
K. B.
Khokonov
,
J. Appl. Phys.
54
,
1346
(
1983
).
5.
L.
Vitos
,
A.
Ruban
,
H.
Skriver
, and
J.
Kollar
,
Surf. Sci.
411
,
186
(
1998
).
6.
E.
Auyeung
,
T. I. N. G.
Li
,
A. J.
Senesi
,
A. L.
Schmucker
,
B. C.
Pals
,
M.
Olvera de La Cruz
, and
C. A.
Mirkin
,
Nature
505
,
73
(
2014
).
7.
C. A.
Mirkin
,
R. L.
Letsinger
, and
R. C.
Mucic
,
Nature
382
,
607
(
1996
).
8.
A. P.
Alivisatos
,
K. P.
Johnsson
,
X.
Peng
, and
T. E.
Wilson
,
Nature
382
,
609
(
1996
).
9.
S. Y.
Park
,
A. K. R.
Lytton-Jean
,
B.
Lee
,
S.
Weigand
,
G. C.
Schatz
, and
C. A.
Mirkin
,
Nature
451
,
553
(
2008
).
10.
D.
Nykypanchuk
,
M. M.
Maye
,
D.
van der Lelie
, and
O.
Gang
,
Nature
451
,
549
(
2008
).
11.
M. R.
Jones
,
R. J.
Macfarlane
,
B.
Lee
,
J.
Zhang
,
K. L.
Young
,
A. J.
Senesi
, and
C. A.
Mirkin
,
Nat. Mater.
9
,
913
(
2010
).
12.
R. J.
Macfarlane
,
B.
Lee
,
M. R.
Jones
,
N.
Harris
,
G. C.
Schatz
, and
C. A.
Mirkin
,
Science
334
,
204
(
2011
).
13.
E.
Auyeung
,
R. J.
Macfarlane
,
C. H. J.
Choi
,
J. I.
Cutler
, and
C. A.
Mirkin
,
Adv. Mater.
24
,
5181
(
2012
).
14.
D. M.
Hinckley
,
J. P.
Lequieu
, and
J. J.
de Pablo
,
J. Chem. Phys.
141
,
035102
(
2014
).
15.
T. I. N. G.
Li
,
R.
Sknepnek
,
R. J.
Macfarlane
,
C. A.
Mirkin
, and
M.
Olvera de la Cruz
,
Nano Lett.
12
,
2509
(
2012
).
16.
T. I. N. G.
Li
,
R.
Sknepnek
, and
M.
Olvera de la Cruz
,
J. Am. Chem. Soc.
135
,
8535
(
2013
).
17.
C.
Knorowski
and
A.
Travesset
,
Curr. Opin. Solid State Mater. Sci.
15
,
262
(
2011
).
18.
C.
Chi
,
F.
Vargas-Lara
,
A. V.
Tkachenko
,
F. W.
Starr
, and
O.
Gang
,
ACS Nano
6
,
6793
(
2012
).
19.
P. M.
Ajayan
and
L. D.
Marks
,
Phase Transitions
24
,
229
(
1990
).
20.
L. D.
Marks
and
P. M.
Ajayan
,
J. Mater. Res.
5
,
1496
(
1990
).
21.
W.
Humphrey
,
A.
Dalke
, and
K.
Schulten
,
J. Mol. Graphics
14
,
33
(
1996
).
22.
J. E.
Stone
, “
An efficient library for parallel ray tracing and animation
,” M.Sc. thesis, Computer Science Department,
University of Missouri-Rolla
,
1998
.
23.
G.
Wulff
,
Z. Kristallogr. - Cryst. Mater.
34
,
449
(
1901
).
24.
25.
R. V.
Zucker
,
D.
Chatain
,
U.
Dahmen
,
S.
Hagège
, and
W. C.
Carter
,
J. Mater. Sci.
47
,
8290
(
2012
).
26.
27.
B.
Sonderegger
and
E.
Kozeschnik
,
Metall. Mater. Trans. A
40
,
499
(
2009
).
28.
J. I.
Cutler
,
E.
Auyeung
, and
C. A.
Mirkin
,
J. Am. Chem. Soc.
134
,
1376
(
2012
).
29.
J. A.
Anderson
,
C. D.
Lorenz
, and
A.
Travesset
,
J. Comput. Phys.
227
,
5342
(
2008
).
30.
T. D.
Nguyen
,
C. L.
Phillips
,
J. A.
Anderson
, and
S. C.
Glotzer
,
Comput. Phys. Commun.
182
,
2307
(
2011
).
32.
W. G.
Hoover
,
Phys. Rev. A
31
,
1695
(
1985
).
33.
C.
Knorowski
and
A.
Travesset
,
Soft Matter
8
,
12053
(
2012
).
34.
P.
Akcora
,
H.
Liu
,
S. K.
Kumar
,
J.
Moll
,
Y.
Li
,
B. C.
Benicewicz
,
L. S.
Schadler
,
D.
Acehan
,
A. Z.
Panagiotopoulos
,
V.
Pryamitsyn
,
V.
Ganesan
,
J.
Ilavsky
,
P.
Thiyagarajan
,
R. H.
Colby
, and
J. F.
Douglas
,
Nat. Mater.
8
,
354
(
2009
).
35.
S. K.
Kumar
and
R.
Krishnamoorti
,
Annu. Rev. Chem. Biomol. Eng.
1
,
37
(
2010
).
36.
S.
Angioletti-Uberti
,
P.
Varilly
,
B. M.
Mognetti
,
A. V.
Tkachenko
, and
D.
Frenkel
,
J. Chem. Phys.
138
,
021102
(
2013
).
37.
W. B.
Rogers
and
J. C.
Crocker
,
Proc. Natl. Acad. Sci. U. S. A.
108
,
15687
(
2011
).