Mechanistic Principles of Colloidal Crystal Growth by Evaporation

Observations concerning colloidal crystal growth are interpreted in terms of the convective steering hypothesis, which posits that solvent flow into a...
0 downloads 0 Views 1MB Size
Langmuir 2008, 24, 13683-13693

13683

Mechanistic Principles of Colloidal Crystal Growth by Evaporation-Induced Convective Steering Damien D. Brewer,‡ Joshua Allen, Michael R. Miller, Juan M. de Santos, Satish Kumar, David J. Norris, Michael Tsapatsis,* and L. E. Scriven† Department of Chemical Engineering and Materials Science, UniVersity of Minnesota, Minneapolis, Minnesota 55455-0132 ReceiVed July 9, 2008. ReVised Manuscript ReceiVed September 21, 2008 We simulate evaporation-driven self-assembly of colloidal crystals using an equivalent network model. Relationships between a regular hexagonally close-packed array of hard, monodisperse spheres, the associated pore space, and selectivity mechanisms for face-centered cubic microstructure propagation are described. By accounting for contact line rearrangement and evaporation at a series of exposed menisci, the equivalent network model describes creeping flow of solvent into and through a rigid colloidal crystal. Observations concerning colloidal crystal growth are interpreted in terms of the convective steering hypothesis, which posits that solvent flow into and through the pore space of the crystal may play a major role in colloidal self-assembly. Aspects of the convective steering and deposition of highPeclet-number rigid spherical particles at a crystal boundary are inferred from spatially resolved solvent flow into the crystal. Gradients in local flow through boundary channels were predicted due to the channels’ spatial distribution relative to a pinned free surface contact line. On the basis of a uniform solvent and particle flux as the criterion for stability of a particular growth plane, these network simulations suggest the stability of a declining {311} crystal interface, a symmetry plane which exclusively propagates fcc microstructure. Network simulations of alternate crystal planes suggest preferential growth front evolution to the declining {311} interface, in consistent agreement with the proposed stability mechanism for preferential fcc microstructure propagation in convective assembly.

Introduction Colloidal self-assembly techniques are being actively investigated as routes to monolayer and multilayer crystals.1-3 Among the most promising applications for self-assembly processes are membrane separations4-8 and photonic crystals.9,10 Since early observations of phase behavior in colloidal suspensions,11 molecular and Brownian dynamics simulations have predicted a wide range of close-packed structures in dynamic and equilibrium states.12-15 Although equilibrium configurations are randomly close-packed in general, early experiments and simulations showed that face centered-cubic (fcc) packings are * Corresponding author. E-mail: [email protected]. † L. E. “Skip” Scriven conceived the model formulation and directed the early stages of this work. He should have been the corresponding author. This paper is dedicated to his memory. ‡ Author e-mail: [email protected]. (1) Jiang, P.; Bertone, J. F.; Hwang, K. S.; Colvin, V. L. Chem. Mater. 1999, 11, 2132–2140. (2) Dushkin, C. D.; Lazarov, G. S.; Kotsev, S. N.; Yoshimura, H.; Nagayama, K. Colloid Polym. Sci. 1999, 277, 914–930. (3) Xia, Y. N.; Gates, B.; Li, Z. Y. AdV. Mater. 2001, 13, 409–413. (4) Hedlund, J.; Sterte, J.; Anthonis, M.; Bons, A.; Carstensen, B.; Corcoran, N.; Cox, D.; Deckman, H.; De Gijnst, W.; de Moor, P.; Lai, F.; McHenry, J.; Mortier, W.; Reinoso, J.; Peters, J. Microporous Mesoporous Mater. 2002, 52, 179–189. (5) Lai, Z.; Bonilla, G.; Diaz, I.; Nery, J. G.; Sujaoti, K.; Amat, M. A.; Kokkoli, E.; Terasaki, O.; Thompson, R. W.; Tsapatsis, M.; Vlachos, D. G. Science 2003, 300, 456–460. (6) Lai, Z.; Tsapatsis, M.; Nicolich, J. P. AdV. Funct. Mater. 2004, 14, 716– 729. (7) Lee, J. S.; Ha, K.; Lee, Y.; Yoon, K. B. AdV. Mater. 2005, 17, 837–841. (8) Newton, M. R.; Bohaty, A. K.; Zhang, Y. H.; White, H. S.; Zharov, I. Langmuir 2006, 22, 4429–4432. (9) Arsenault, A.; Fournier-Bidoz, S. B.; Hatton, B.; Miguez, H.; Tetrault, N.; Vekris, E.; Wong, S.; Yang, S. M.; Kitaev, V.; Ozin, G. A. J. Mater. Chem. 2004, 14, 781–794. (10) Norris, D. J. Nat. Mater. 2007, 6, 177–178. (11) Perrin, J. Ann. Chim. Phys. 1909, 18, 1. (12) Alder, B. J.; Wainwright, T. E. J. Chem. Phys. 1957, 27, 1208. (13) Schmidt, M.; Loewen, H. Phys. ReV. Lett. 1996, 76, 4552–4555. (14) Zangi, R.; Rice, S. A. Phys. ReV. E 2000, 61, 660–670. (15) Pusey, P. N.; Van Megen, W. Nature 1986, 320, 340–342.

favored over hexagonally close-packed (hcp) structures.16,17 Because the free energy difference between hcp and fcc packings is a miniscule quantity (on the order of 0.005 RT/mol), hydrodynamic and colloidal interactions must direct crystal assembly in all but the most extraordinary phase-transition processes. Following the demonstration of thin, close-packed particle layers from spin coating,18 it was suggested that a highly ordered packing might be deposited on a substrate by contriving to make the meniscus contact line travel over the substrate at a constant speed.2,19 Various techniques aimed at realizing this evaporationdriven “convective assembly” have been investigated over the past decade.1,2,20-25 In these techniques, evaporation from quasistatic menisci engenders convective flow toward the nascent crystal. When the particles are 100 nm or greater, capillary and viscous forces appear to overcome Brownian motion and assemble these particles into ordered, predominantly fcc structures. “Synthetic” colloidal crystals fabricated by convective assembly and sedimented opals show similar tendencies toward fcc structure, although sedimented opals have time to reach equilibrium structure. (16) Pusey, P. N.; Van Megen, W.; Bartlett, P.; Ackerson, B. J.; Rarity, J. G.; Underwood, S. M. Phys. ReV. Lett. 1989, 63, 2753–2756. (17) Woodcock, L. V. Nature 1997, 385, 141–143. (18) Deckman, H. W.; Dunsmuir, J. H. Appl. Phys. Lett. 1982, 41, 377–379. (19) Dimitrov, A. S.; Nagayama, K. Langmuir 1996, 12, 1303–1311. (20) Lee, J. A.; Meng, L. L.; Norris, D. J.; Scriven, L. E.; Tsapatsis, M. Langmuir 2006, 22, 5217–5219. (21) Velikov, K. P.; Christova, C. G.; Dullens, R. P. A.; van Blaaderen, A. Science 2002, 296, 106–109. (22) Snyder, M. A.; Lee, J. A.; Davis, T. M.; Scriven, L. E.; Tsapatsis, M. Langmuir 2007, 23, 9924–9928. (23) Wong, S.; Kitaev, V.; Ozin, G. A. J. Am. Chem. Soc. 2003, 125, 15589– 15598. (24) Prevo, B. G.; Velev, O. D. Langmuir 2004, 20, 2099–2107. (25) Kim, M. H.; Im, S. H.; Park, O. O. AdV. Funct. Mater. 2005, 15, 1329– 1335.

10.1021/la802180d CCC: $40.75  2008 American Chemical Society Published on Web 11/07/2008

13684 Langmuir, Vol. 24, No. 23, 2008

These observations26 prompt many questions regarding microstructure formation mechanisms; in particular, it is currently unclear how important solvent flow is to the convective assembly process.27 Because synthetic colloidal crystals form over timescales of hours and show fcc tendency similar to that observed in sedimented suspensions, solvent flow into and through the colloidal crystal may play a key role in microstructure evolution, propagation, or both. The associated length and time scales are many; spheres as they take their places at the crystal boundaries may experience capillary pressure and viscous forces that drive their motion, yet structural organization may follow at longer timescales.12 These hydrodynamic stresses could be computed with accuracy by full solution of the Navier-Stokes equations,28 yet the assembly mechanism must also be subject to mass transfer processes throughout the entire macroscopic domain. A useful compromise between full discretization of a solvent continuum and a coarse, volume-averaged model of the process is the network model of porous media.29-32 Among the predominant descriptions of fluid flow in both ordered and disordered porous media, the network model has provided insight into the link between microscopic structure and macroscopic transport properties. In this approach, an equivalent electrical network describes fluid motion through a connected graph of nodes. In his pioneering work on network models for petroleum applications, Fatt29 reduced an irregular, threedimensional structure to a regular, two-dimensional graph of nodes formed by the intersection of a regular array of channels. More recent studies have confirmed the validity of network models for the general case of unordered sphere packings.30,31 Additional network models of transport through porous media have been applied to a variety of problems, including drying,32 multiphase flow,33 percolation thresholds, and pressure drop in packed beds.34 While theoretical treatments of evaporation-induced assembly are available for drying sessile drops,35-37 we are not aware of any such studies in which the evaporation at menisci is the driving force for convective assembly. In this paper, a network model of evaporation-driven creeping flow in close-packed crystals of monodisperse spheres is presented. The relationships between a regular porous medium, its pore space, and the associated fluid motion are analyzed in hexagonally close-packed sphere arrays. Unlike previous studies, the present investigation focuses on meniscus configuration and mass transport at the network boundaries. We begin with the convective steering hypothesis,27 which posits that convective flow may preferentially guide spheres to regions of the crystal boundary;28 microstructure evolution, propagation, or both could be influenced by convective steering. In photonic crystal and monolayer crystal assembly, these mechanisms may be important to understand both in cases when the fcc structure is desired and when alternative crystal structures are required. Here, we focus (26) Meng, L.; Wei, H.; Nagel, A.; Wiley, B. J.; Scriven, L. E.; Norris, D. J. Nano Lett. 2006, 6, 2249–2253. (27) Norris, D. J.; Arlinghaus, E. G.; Meng, L.; Heiny, R.; Scriven, L. E. AdV. Mater. 2004, 16, 1393–1399. (28) Gasperino, D.; Meng, L. L.; Norris, D. J.; Derby, J. J. J. Cryst. Growth 2008, 310, 131–139. (29) Fatt, I. J. Pet. Technol. 1956, 8, 144–177. (30) Bryant, S. L.; King, P. R.; Mellor, D. W. Transp. Porous Media 1993, 11, 53–70. (31) Bryant, S. L.; Mellor, D. W.; Cade, C. A. AIChE J. 1993, 39, 387–396. (32) Nowicki, S. C.; Davis, H. T.; Scriven, L. E. Drying Technol. 1992, 10, 925–946. (33) Melli, T. R.; Scriven, L. E. Ind. Eng. Chem. Res. 1991, 30, 951–969. (34) Martins, A. A.; Laranjeira, P. E.; Lopes, J. C. B.; Dias, M. M. AIChE J. 2007, 53, 91–107. (35) Adachi, E.; Dimitrov, A. S.; Nagayama, K. Langmuir 1995, 11, 1057– 1060. (36) Deegan, R. D. Phys. ReV. E 2000, 61, 475–485. (37) Nguyen, V. X.; Stebe, K. J. Phys. ReV. Lett. 2002, 88, 164501.

Brewer et al.

on microstructure propagation by considering creeping solvent flow into rigid colloidal crystals; dynamic microstructure evolution is outside the present scope. Due to its formulation, however, the network model is quite general in treating ordered arrangements of spherical particles and could be extended to binary superlattice structures as well as irregular crystals. The main objectives of this article are to introduce a network model of flow in convective assembly processing and to highlight the potential role of the growing crystal interface in microstructure propagation. Specifically, a uniform solvent and particle flux near the crystal interface is proposed as the condition for “stable” interface structure. In order to simulate solvent flow into and through a colloidal crystal, we present an approximate representation of the pore space, referred to as the “equivalent network”. The associated network equations, following from the analogous matrix formulation for electrical networks, are enumerated along with a set of underlying assumptions. A set of crystal interfaces of varied crystal symmetry and relative inclination are simulated. We first consider solvent flow into and through the inclined {111} planar interface in detail, both as an extension to previous studies and as validation of the model. We then compare flow into the inclined {111} plane with flow into the {100} and {311} symmetry planes, and it is found that a slightly declining {311} interface imposes the most uniform solvent and particle flux.

Equivalent Network Development Packing of colloidal spheres into a regular array creates not only an ordered crystal but an accompanying porespace. The porespace can be reduced to an equiValent network by considering it as a series of nodes and branches as described below. farthermore, determination of flow into and through the porespace requires an equivalent representation of the porespace boundaries, including entrance channels, curved menisci, and drying interstices. Nodes and branches are each associated with specific sets of variables. The branch variables are currents that correspond to flow, while nodal variables are potentials that correspond to pressure. The network formulation also encompasses constitutive behavior of the curved menisci. In what follows, we describe each of these model components in turn. We then relate geometric aspects of a sphere packing and the accompanying porespace with the convective assembly process. Monodispersed spherical particles can be assembled into monolayer hexagonally coordinated arrays. If these, in turn, are stacked in the third dimension, then three mutually adjacent spheres of a given packing layer may support one sphere in the overlying layer. Arrays of spherical particles stacked in this manner belong to the hexagonally close-packed family of crystals, including fcc, hcp, and random hexagonal close-packed (rhcp) structures. The fcc structure is achieved with an ABCABC · · · sequence of stacking layers, while hcp is achieved with an ABAB · · · sequence. Deviations from either of these two sequences are referred to as stacking faults, also known as twinning defects. A random stacking of hexagonally packed planes is referred to as random hexagonal close-packed (rhcp). In twoand three-dimensional colloidal crystals, the solid matrix forms a regular array of close-packed spherical bodies, the interstices and constrictions between which fluid flows. These interstices and constrictions play a special role in network formulation of the equations of fluid motion: these respective regions are associated with large cross section, negligible velocity, or narrow cross section, finite velocity flow patterns. The interstitial spaces, or pore bodies, of a regular colloidal crystal form an equally regular network of nodes connected through converging-diverging constrictions, or pore throats.

Colloidal Crystal Growth

Langmuir, Vol. 24, No. 23, 2008 13685

Figure 2. (a) Illustration of two growth front entrance channels, which may be either clear or obstructed. Particle placement at the clear niches propagates fcc structure, while placement at obstructed niches propagates twinning defects. (b) Illustration of the convective assembly process.

Figure 1. Geometry of a (a) pore throat, (b) tetracoordinate pore body, and (c) octacoordinate pore body. (d) Spatial arrangement of tetracoordinate and octacoordinate pore bodies in hexagonally close-packed crystals. Tetracoordinate (sphere) and octacoordinate (cube) pore bodies and pore throats (rods) form a periodic void space through which solvent flows. Face-centered cubic lattice niches are those connected directly to octacoordinate pore bodies. These octacoordinate pore bodies engage two pore throats in the [111] direction and therefore provide minimal resistance to flow in fcc growth front niches.

Hence, an equivalent network representation for the porespace between spheres has two basic elements: pore bodies and pore throats. Illustrated in Figure 1a, pore throats between three mutually adjacent spheres may be characterized by an effective orifice dimension. In hexagonally close-packed arrays, all pore throats are defined by three adjacent spheres. Additionally, two types of pore bodies exist in close-packed hexagonal sphere arrays. Where three mutually adjacent spheres of a given packing layer also contact a single sphere of an adjacent layer, a tetracoordinate pore body forms as shown in Figure 1b. Where three mutually adjacent spheres of a given layer overlie three mutually adjacent spheres in the next, an octacoordinate pore body forms as shown in Figure 1c. Similarly, a crystal in the hexagonally close-packed family allows for distinction between boundary entrance or exit channels, hereafter referred to as niches. Once two or more rows of spheres appear in the sphere pack, two types of growth front niches become hydrodynamically distinct as illustrated in Figure 1d and Figure 2a.38 Particle placement at the first of these niches extends the fcc crystal (fcc niche), and particle placement at the (38) Arlinghaus, E. G. PhD. Thesis, University of Minnesota, 2004.

second creates a twin defect (twinning niche). The fcc niche is less obstructed than the twinning niche (Figure 2a). More specifically, there exists a clear path from the exterior of the crystal to an octacoordinate pore body and directly onward to additional layers in the crystal. The twinning niche leads to a tetracoordinate pore body, which then diverges in three directions defined by its tetracoordinate connectivity (Figure 1d). The transport of fluid (and with it, particles) to and across the sphere pack boundaries is of primary interest. In evaporationdriven self-assembly experiments, a liquid pool containing microspheres wets a substrate as it is evaporated; the contact line is typically pinned to the edge of a growing planar crystal (see Figure 2b).1 At the crystal boundary farthest from the substrate, the solvent forms menisci between exposed microspheres. Wetted crystal regions typically extend in the lateral and transverse directions on the order of millimeters or centimeters, several orders of magnitude greater than an individual crystal unit. Hence, we are concerned with three boundaries: the submerged growth front region at which microspheres are deposited, the exposed array of menisci from which liquid evaporates to drive the assembly process, and a rear drying front marking a transition from wetted to dry regions. The experimental system of interest has been described in detail elsewhere.26,27 Figure 3 illustrates the equivalent network for a small facecentered cubic array. At the interface between suspension and crystal, or growth front, solvent flows into and through the sphere packing shown schematically in Figure 3a. Hereafter, we shall refer to the periodic direction (x) as the lateral direction, the primary flow direction (y) as the transVerse direction, and the substrate normal (z) as the Vertical direction according to Figure 3. Pinned above the growth front interface is a free surface. We assume that the free surface is not highly curved, so that the growth front is at a uniform capillary pressure (Figure 2b). The mean curvature of this free surface is actually governed by capillary statics, H ≡ (d2 Z)/(d X2) ≈ (Fg/2σ)1/2, where 2H is the mean curvature of the pinned meniscus, Z is the liquid thickness

13686 Langmuir, Vol. 24, No. 23, 2008

Brewer et al.

Figure 3. Illustration of problem geometry: (a) close-packed bed of monodisperse spheres, the leading edge of which contacts suspension (left) and the rear end of which is the drying front region (right); (b) an equivalent network in which nodes are represented as small spheres inscribed within interstices defined by the sphere pack, including tetracoordinate pore bodies (green), octacoordinate pore bodies (blue), growth front niches (white), menisci (red), and drying front nodes (purple).

relative to the substrate, X is the distance along the substrate, F is density, g is gravity, and σ is surface tension.39,40 If the colloidal suspension is aqueous, then the capillary pressure is less than 100 Pa. Farthest from the substrate boundary, a curved liquid-vapor interface forms near the particle contacts. In the colloidal systems considered here, the pores between particles are small enough that capillary pressure dominates the process. Defined by the local region between three adjacent spheres, one section of the liquid-vapor interface is referred to simply as a meniscus. We consider processes for which the vapor is not saturated in solvent near the interface. At each meniscus, mass transfer draws solvent from the liquid to the vapor phase. In the absence of supplementary liquid flow, drainage into the pore space would occur.32,41 In the present study, we model steadystate mass transfer, in which the supplementary liquid flow matches the external mass transfer. At a macroscopic distance from the growth front interface is a drying front. Characterized by a transition from complete to partial liquid saturation in the pore space, the drying front represents another liquid-vapor interface embedded within the crystal pore space.32,41 Onset of the drying transition follows from meniscus instability at the crystal-air boundary, a condition that was described by Haines.42 Beyond the drying transition is a region of pore space partially saturated with solvent and partially invaded by air. Vapor diffusion mediates mass transport within the dry or partially wetted pore space, where partial vapor saturation and tortuosity hinder diffusion;32 hence, the transverse solvent flux is small compared with the solvent flux near the crystal-solution interface. Particle and solvent fluxes near the (39) Landau, L.; Levich, B. Acta Physicochim. 1942, 17, 42–54. (40) White, D. A.; Tallmadge, J. A. Chem. Eng. Sci. 1965, 20, 33–37. (41) Segura, L. A.; Toledo, P. G. Lat. Am. Appl. Res. 2005, 35, 43–50. (42) Haines, W. B. J. Agric. Soc. UniV. Coll. Wales 1927, 17, 264–289.

crystal boundary will be proportional to the wetting length between growth and drying fronts.1,19 Convective steering of solvent and particle transport into the crystal may also be influenced by the crystal symmetry and inclination of the growth plane. On the basis of the hypothesis that the exposed crystal facet and structure may direct convective steering, the present model reduces the crystal pore space into a discretized network of nodes and solvent transport through the boundaries as a set of discrete currents. This equivalent network formulation implicitly defines an effective internal resistance based on the pore space between packed spheres. In what follows, we also propose a simple approximation for the external mass transfer resistance from each meniscus. On the basis of these approximations, the distributions of capillary pressure and internal and external solvent fluxes are calculated.

Mathematical Model Previous results in the science of porous media have established that interstitial spaces and the converging-diverging passages that connect them can be accurately reduced to pore bodies and pore throats.30 As such, pore bodies are modeled as pockets of uniform pressure distribution. Because pore throats feature a hyperbolic narrowing and widening contour about a minimum cross section, the resistance can be defined by the effective circular orifice dimension as in Figure 1a. In the creeping-flow regime, mechanical potential driven flow rate increases linearly with the pressure driving force as (r3/3µ)∆P, where r is the radius of a circular orifice, µ is the fluid viscosity, and ∆P is the pressure differential. In the case of three mutually adjacent and monodisperse spheres, a lower bound for an “effective” aperture radius is approximately the radius of an inscribed circle: rc ) 0.1547R. An upper bound is provided by the radius for which the circular

Colloidal Crystal Growth

Langmuir, Vol. 24, No. 23, 2008 13687

Specifically, we require that growth front boundary nodes remain at a uniform capillary pressure. Hence, the boundary condition for growth front niches is

Ph ) P0

Figure 4. (a) Schematic diagram of a meniscus defined by the region between three adjacent spheres. The mean curvature 2H ) 2/RM is approximated as the radius of an inscribed sphere (cross-hatched circle) of wetting height h above the particle centers; the effective mass transfer area is approximated by the area of a flat meniscus of the same wetting height h. (b) A cross-sectional view of the meniscus through the dashed plane in (a). Simple geometric equations describe the relationships between meniscus height h, meniscus curvature 2/RM, contact angle θ, and particle radius R. Quantities r, L, Ψ, and δ are geometric variables and parameters.

throat has equal area to the actual cross section minimum: re ) 0.2266R. A good estimate for the effective radius is simply the arithmetic average r ) (rc + re)/2 ) 0.1907R.30 All pore throats narrow into the same geometrical region defined by three spheres, regardless of orientation, and polydispersity is the principal contribution to any variations thereof. In the remainder of this section, we shall refer to pore bodies and pore throats as nodes and branches, respectively. The analysis is aided by a few simplifying assumptions. Most self-assembly experiments are run with nearly monodisperse 0.1 to 1 µm ((3%) diameter microspheres at concentrations of 0.1% to 4.0% by volume.1,26,43,44 Characteristic velocities and shear rates are typically on the order of 100 µm/s and 10 s-1, resulting in Reynolds and Peclet numbers on the order of 10-4 and 107 when R ) 0.5 µm. In this work, all particles are modeled as spherical and perfectly monodisperse. All crystals considered here belong to the hexagonal close-packed (fcc, hcp, or rhcp) crystal family, and spheres as they assemble on the growth front are fixed in place by local hydrodynamic and van der Waals forces. When the assembling particles are microspheres, the Peclet number is much greater than unity and Brownian motion is neglected. Spheres are considered to be rigid so that they pack closely without deforming,38 and the suspension is dilute so that particle motion does not perturb fluid streamlines in the vicinity of the growth front. As in many electrical network analysis problems, the majority of equations in the present system are linear. Mass conservation at each node of the network is represented in the linear matrix form

AP + BQ ) 0

(1)

where A is the analogue of a linear “admittance” matrix and BQ is the analogue of injected boundary currents, while P and Q are the vectors of nodal pressures and boundary node flow rates.45,46 As in the modified network analysis of Ho, Ruehli, and Brennan,46 we augment eq 1 with additional constitutive equations for each nonlinear element and expand the matrix formulation accordingly. (43) Wei, H.; Meng, L.; Jun, Y.; Norris, D. J. Appl. Phys. Lett. 2006, 89, 241913. (44) Malaquin, L.; Kraus, T.; Schmid, H.; Delamarche, E.; Wolf, H. Langmuir 2007, 23, 11513–11521. (45) Strang, G. Introduction to Applied Mathematics; Wellesley-Cambridge Press: Wellesley, MA, 1986. (46) Ho, C.; Ruehli, A. E.; Brennan, P. A. IEEE Trans. Circuits Syst. 1975, 22, 504–509.

(2)

where h is an index over all growth front niches, Ph is the associated pressure, and P0 is a reference pressure. In this work, the reference pressure could be interpreted as a uniform mechanical potential, including hydrostatic and stagnation pressures. farthermore, we require that pressure differentials across exposed menisci satisfy the Young-Laplace relationship and that the evaporation rate is consistent with the meniscus geometry. However, general analytical expressions for meniscus shape and capillary pressure are unavailable except in the simplest cases. In order to capture the meniscus behavior, we model the surface as a spherical interface.47 Then, a single radius describes the mean curvature, and capillary pressure is defined by

Pk ) Pamb - 2Hkσ

(3)

where k is an index over all menisci, Pk and 2Hk are the associated liquid pressure and mean meniscus curvature, and σ is the surface tension. A complete and accurate description of the geometry would pose a more expensive computational challenge, as the details of equilibrium or quasi-static meniscus behavior in confined geometries require numerical techniques.48 However, the proposed approximation, first described by Purcell,47 has been shown to yield good agreement with numerical simulation for wetting liquids (θ < 90°) as supposed in this study.48 Detailed measurement of effects due to roughness and surface functionalization has not yet been possible in this system, although such factors are known to control the wettability behavior.49 Hence, the present geometric argument provides a satisfactory description of the nonlinear influence of evaporation at exposed menisci. Geometric analysis of hexagonal and simple cubic sphere packs provides constitutive equations for the pressure and interfacial mass transfer rate at menisci. Schematic profiles of a trigonal meniscus are shown in Figure 4, including a cross-sectional view of the symmetry plane. In the special case of complete wetting (θ ) 0°), HkR ) cos(Ψk)/(1 - cos(Ψk) + r0/R) where sin(Ψk) ) hk/R, hk is the meniscus height, R is the particle radius, and r0 ) (δ - 1)R is the orifice radius. In the general case, HkR ) cos(θ + Ψk)/(1 - cos(Ψk) + r0/R) (see Supporting Information Figure S1). Rewriting this general equation in terms of the meniscus height yields

HkR )

εk(2 - εk2)1⁄2cos(θ) - (1 - εk2)sin(θ) δ - εk(2 - εk2)1⁄2

(4)

where εk ) (1 - hk/R)1/2 is a transformed meniscus height, introduced to avoid singular Jacobian matrices. The effective interfacial area is the projected cross section of wetted area from above, shown as the shaded region in Figure 4a, and naturally depends upon the meniscus height and particle radius: A ) RR2 + βh2 where R and β are geometric factors. When three mutually adjacent spheres form a meniscus, then R ) 3 - π/2, β ) π/2, and δ ) 23. Since the interfacial volume transfer rate is simply the product of evaporation flux JE and interfacial area, the evaporation flow rate is simply defined by (47) Purcell, W. J. Trans. AIME 1950, 189, 369–371. (48) Hilden, J. L.; Trumble, K. P. J. Colloid Interface Sci. 2003, 267, 463–474. (49) Wolansky, G.; Marmur, A. Langmuir 1998, 14, 5292–5297.

13688 Langmuir, Vol. 24, No. 23, 2008

Brewer et al.

Qk ) JE(RR2 + βhk2)

(5)

where Qk and hk are the (external) evaporation rate and meniscus height at meniscus k. On a liquid volume basis, note that the evaporation flux is JE ) kG(FG/FL), where kG is the gas-phase mass transfer coefficient and FG and FL are the respective gas and liquid densities. At the rear drying front nodes, we impose a no-flux boundary condition:

Qj ) 0

(6)

In the convective assembly process, the rear transverse boundary will end in a drying front at which capillary action and diffusion through the porespace determine the transport and pressure distribution.32 The mass transfer resistance is roughly an order of magnitude lower in the porespace than at the exterior, and partial solvent saturation reduces the driving force for vapor transport. On the lateral edges, we apply periodic boundary conditions to obtain symmetry about the central plane of the system, and the simulation approaches a two-dimensional domain. Isolated growth front defects and features can also be modeled using larger lateral dimensions as described in the next section. Constitutive eqs 3, 4, and 5 were combined to yield a set of 3m equations, where m is the number of menisci. An additional set of p Dirichlet conditions at the growth front nodes, q Dirichlet conditions at the rear transverse boundary nodes, and n mass conservation equations at each node yield a total of n + p + q + 3m equations in n + p + q + 3m unknowns. These unknowns are Pnx1, Hmx1, εmx1, and Q(m+p+q)x1 where P is pressure, H is meniscus curvature, ε is meniscus height, and Q is external flow rate. Data are summarized in Table 1; linear mass conservation at interior nodes is represented by eq 1 and boundary conditions are given in eqs 2-6. Equations 1, 2, 3, and 6 are linear in pressure and interfacial flow rates, while eqs 4 and 5 are nonlinear in pressure, interfacial flow rates, meniscus curvature, and meniscus height. These equations are solved iteratively with different evaporation flux until the meniscus height is sufficiently small at the drying front. In this study, a meniscus height below 3% of the particle radius is considered sufficiently small. The input parameters are the sphere positions and their radii, relevant solvent properties, contact angle, and ambient pressure. Although the evaporation rate is a parameter, only one value is consistent with meniscus instability at the rear transverse boundary of the simulation domain. All simulations reported here correspond to point-contact sphere packing of 1 µm diameter spheres. With respect to the substrate plane, spheres are stacked in the fcc (repeated “ABC”) sequence, and the {111}, {100}, and {311} growth planes are considered. Solvent viscosity, surface tension, and density are taken to be those of water at room temperature. Ambient pressure is taken to be atmospheric, the contact angle is taken to be 0° (complete wetting), and a typical evaporation flux is 1 cm/s. These evaporation rates correspond to mass transfer into a vacuum, which represents the highest physical rate of evaporation. We chose the evaporation flux to allow computationally accessible simulation sizes, and other properties were chosen to roughly match experimental conditions in previous studies.43,50

Any close-packed crystal geometry can be simulated using this model, including those with long-range disorder. However, we are primarily concerned with local flow behavior at the boundaries of fcc sphere packs. Although the network model is a simpler representation of local structure than Voronoi Tessellation,30,51 boundary integral,52 or meshed finite element53 solutions, a pore-network representation captures the microscopic aspects of solvent drainage32,41 and flow behavior in packed bed geometry while providing a link to macroscopic transport parameters.33 For instance, recent finite element simulations of the full Navier-Stokes equations also validated the network model for flow through close-packed arrays of monodisperse spheres.28 However, the same pore-network approach has also been highly successful in predicting such macroscopic parameters as Darcy permeability32,34,41 and liquid saturation curves.32 Hence, the pore-level network approach is expected to reproduce primary features of evaporation-induced assembly. Computational Methods. For each parameter set (including solvent properties and gas-phase mass transfer coefficient), we solve a large nonlinear set of equations given in Table 1 by the Newton algorithm. Within each Newton step, we solve the large linear sparse system using a package from the CHOLMOD library.54 Solutions for the desired gas-phase mass transfer coefficient were obtained by continuation over the dimensionless mass transfer coefficient kP0/Q0 ) (reff3/3µ)(Pamb/4JER2), beginning with a fully linear network as the initial guess. Simulations were performed on IBM LS21 nodes. Computations included between 10 000 and 20 000 degrees of freedom, and simulation times ranged from approximately 30 to 150 CPU hours.

Computational Results and Discussion Simulations were performed to study the capillary pressure and solvent flow distributions through assembling colloidal crystals over a range of thicknesses, lengths, evaporation fluxes, wetting behavior, and solvent properties. All of these parameters influence the capillary pressure and entrance flow profiles, although a detailed parametric study is outside the scope of this paper. Here, we first propose a mechanism for fcc selectivity in the context of a representative simulation. Then, we describe how isolated nucleation events may support fcc crystal structure propagation when twinning niches as well as fcc niches are present. In this case, the hydrodynamic distinction between growth front niches favors fcc propagation. Finally, we propose uniform particle flux as the main criterion for stability of crystal propagation. On the basis of this criterion, simulations of flow into and through alternative crystal facets suggest that a {311} interface may be the most stable to steady convective assembly. Selectivity of fcc over Twinning at the {111} Planar Interface. We find two main causes for selectivity of fcc over twinning structure propagation from {111} planar interfaces in convective assembly of colloidal crystals. The first cause arises from individual variations in resistance to solvent entry as depicted in Figure 1d and Figure 2a. More specifically, any local region of the growth front could feature two or more niche types as distinguished by their respective resistances to entrance flow. An intrinsic selectivity for convective steering toward the least

Table 1. Model Equations (1) (2) (3) (4) (5) (6)

∑jn) 1 AijPj +∑jm) 1 BijQj ) 0 Ph ) P0 Pk ) Pamb - 2Hkσ HkR ) [εk(2 - εk2)1/2 cos (θ) - (1 - εk2) sin (θ)]/ [δ - εk(2 - εk2)1/2] Qk ) JER2[(R + β) - βεk2(2 - εk2)] Qj ) 0

nodes i ) 1, · · · , n niches h ) 1, · · · , p menisci k ) 1, · · · , m menisci k ) 1, · · · , m

mass conservation at nodes uniform pressure at growth front niches young-laplace equation at menisci meniscus configuration at menisci

menisci k ) 1, · · · , m rear drying front nodes j ) 1, · · · , q

interfacial mass transfer at mensici flow rate at rear drying front nodes

Colloidal Crystal Growth

obstructed niches was previously studied by Arlinghaus38 and Gasperino et al.28 When the exposed growth front is a planar {111} facet, network theory predicts selectivity of 4/7 to 3/7 in favor of the less obstructed fcc niches as in Figure 2a,38 and finite volume calculations predict the selectivity to be around 20% rather than 33%.28 The second cause is a more “coarse” hydrodynamic effect that arises due to a nonuniform internal resistance, or transverse pore space through which solvent must flow, along the entire vertical dimension of an inclined or declined growth front. In what follows, convective steering of solvent and particles due to internal pore structure, growth front inclination, and evaporative mass transfer are described. Figure 5a shows the volume-averaged capillary pressure distribution corresponding to the inclined crystal facet of an 18-layer crystal of transverse dimension 150 particle diameters, where JE ) 3.18 cm/s. The crystal is close-packed fcc, and the growth front is a {111} plane inclined at 70.5° relative to the substrate (Figure 3a). Several layers into the crystal, the capillary pressure distribution is nearly uniform in the vertical and lateral dimensions. At the entrance, however, solvent traverses a greater distance through porespace near the bottom (substrate) leading edge than at the top (free surface) leading edge; as a direct consequence, the (scalar) transverse pressure gradient increases with distance from the substrate. Because solvent flux is proportional to capillary pressure gradients in creeping flow, this result indicates that the associated solvent flux will also increase with normal distance from the substrate. The solid line in Figure 5b is a plot of the solvent flux corresponding to Figure 5a, averaged over the obstructed and clear niches alike (Figure 2a), versus the vertical distance from the substrate. Scatter points in the solvent flux reflect the relative resistance between clear and obstructed niches. Solvent and particle flux through clear niches is slightly higher than average, while the flux through obstructed niches is lower. Associated meniscus height and radius of curvature profiles show the greatest capillary pressure drops near the leading edge where the solvent flux is highest (Figure 4b and Figure 5c). Both the meniscus height and the radius of curvature decrease at a relatively large rate through the first 20% of the crystal, then decrease at nearly constant rates though the remaining wetted length. Figure 6a illustrates entrance flow into the inclined {111} interface corresponding to Figure 5. In this diagram, the front face is submerged beneath a free surface, which is pinned to a contact line at the top corner of the front face. The projected rods represent values of local influx or efflux rates at each of the growth front niches and menisci. This visual representation shows that entrance flow increases monotonically with vertical distance from the substrate. Figure 6b quantitatively shows how the selectivity for entrance flow into fcc over twinning niches varies with normal distance from the substrate. Selectivity is defined as the ratio of entrance flow into fcc and twinning niches. Each data point indicates a pairwise comparison between entrance flow into adjacent fcc and twinning niches with the normal distance reported as the average value of their vertical coordinates. Superposed onto the dichotomous entrance velocity,28,38 the velocity profile induced by the growth front inclination may enhance or diminish the inherent convective steering preference toward fcc niches. Note that the niche farther displaced from the substrate alternates between fcc and twinning, and for an inclined {111} plane solvent, velocity gradients favor flow into those displaced farthest from the substrate. As a result, selectivity for convective steering to fcc niches alternates between relatively high (150%) and low (105%) values when the fcc niche is above or below the adjacent

Langmuir, Vol. 24, No. 23, 2008 13689

Figure 5. Pressure, flow rate, and meniscus curvature distributions in an 18-layer crystal with JE ) 3.18 cm/s. (a) Volume-averaged pressure distribution in a cross section of the porespace. Given a uniform pressure distribution at its entrance, an inclined {111} growth front imposes a transverse pressure gradient profile. Pressure contours (in bars) are closely spaced near the top leading edge, indicating a proportional solvent flux. A uniform vertical pressure field develops within a few sphere diameters from the entrance. (b) Solvent flux profile at the entrance of the same inclined {111} growth plane. Solvent flux is normalized by the maximum flux. Scatter in the solvent flux data is due to the dichotomy between fcc and twinning entrance flow; the solid line represents the locally averaged solvent flux. (c) Meniscus height and radius of curvature profiles in the transverse direction, both scaled by the particle radius. All distances are scaled by particle diameter.

twinning niche, respectively. For comparison, the 4/3 (fcc versus twinning) selectivity into an unbounded {111} plane is shown, and the simulation results in Figure 6b approach this limiting theoretical value farthest away from the substrate and free surface boundaries. The entrance flow selectivity in Figure 6 also suggests that particle arrival to the hypothetical inclined {111} growth front

13690 Langmuir, Vol. 24, No. 23, 2008

Figure 6. Solvent flow into an 18-layer inclined {111} growth plane, where JE ) 3.18 cm/s. (a) Potential growth niches at the front (small colored spheres) direct flow, represented by protruding rods, preferentially into fcc niches. Protruding rods are scaled and shaded to emphasize the contrast in flow rates. All flow rates are normalized by the maximum flow rate. (b) Selectivity for entrance flow into fcc over twinning niches plotted against vertical distance from the substrate, where vertical distance is scaled by the sphere diameter. Each simulation point represents a comparison between two adjacent fcc and twinning niches whose vertical displacements from the substrate are unequal. In the limit of an infinite crystal facet, application of the network theory yields a 4/3 selectivity for flow into fcc niches over twinning niches.

would proceed sequentially from the top of the sphere pack down. Because the highest influx rate is near the top leading edge, newly arriving particles are most likely to be steered toward growth front niches at the top leading edge. Since the uppermost growth front niches are also fcc niches, growth onto an inclined {111} facet tends toward fcc propagation. This result follows from the postulate that structure may induce flow near the growth front; more specifically, solvent flow due to local pore structure and the growth plane inclination complement one another in favor of fcc microstructure propagation over twinning defect formation. Crystal layering thickness at the leading edge also plays a role in convective steering. Because the liquid-vapor interface and the substrate break the growth front symmetry plane, the crystal thickness influences fcc selectivity. As shown in Figure 6, the dichotomy between local flux into clear and obstructed niches is enhanced near the liquid-vapor interface. Shallow crystals therefore show significant tendency toward fcc propagation on the basis of convective steering, while thick crystals induce flow fields similar to the case of an infinite planar facet. In a threelayer fcc simulation, for instance, the rate at which liquid enters adjacent fcc and twinning niches differs by a factor of about 50% in favor of fcc niches nearest the free surface. While this mechanism does not address dynamic microstructure evolution from crystal nucleation, the mechanism could help explain fcc lattice propagation once it has formed. In experiment, the fcc tendency of a colloidal crystal increases with its thickness,43 and this feature may also be related to the variability in niche selectivity. Growth Front Propagation of the {111} Planar Interface. As arriving particles nucleate at growth front niches, they may also perturb the local solvent velocity field. When an inclined {111} crystal face is disturbed by one or more nucleation events, the local solvent flow increases in the adjacent region. The flow field perturbation due to an outcropping sphere at a {111} fcc growth face is short-ranged (see Supporting Information Figure

Brewer et al.

S2). Regardless of the nucleation site (fcc or twinning niche), the adjacent lattice niches within one fcc position to the nucleated particle experience a moderate increase in flow rate (typically 20% increase) regardless of its distance from the substrate; submerged growth niches any farther experience a negligible change in entrance flow rate. Arranged doublets or triplets of nucleated particles show similar trends, where the nearest niches remain preferentially selected for convective steering. Since the immediately adjacent niches are sterically unavailable due to the nucleated spheres, the convective steering likely favors the nearest available nucleation niche. Two distinct possibilities arise here. First, the flow field can force incoming particles to approach and dislodge a settled particle, or the flow field can direct incoming particles to the nearest available niches. If the nucleated spheres are relatively fixed due to the asymmetric normal stress distribution at its surface, then the latter outcome prevails. In the event of particle rearrangement, a close-packed ordering still remains, but may introduce a twinning fault. Such faults are in fact observed43 with frequency depending on the crystal layer height. However, the more likely event is for the first set of particles on a growth plane to lodge, or to become lodged by rearrangement, upon the growth niches at the highest vertical positions available. These niches propagate the fcc lattice, and when one sphere rests in this location, the steady-state flux of liquid into the sphere pack strongly favors the two adjacent top niches over the two fcc niches below the lodged particle. Local flow perturbations could propagate either local fcc propagation or a twinning fault; however, a preference for nucleation at or near the top row of fcc niches, described in the previous section, makes the former more likely. Note that lateral capillary forces also draw particles into close contact at the free surface,55 farther aiding the proposed fcc propagation mechanism. In this case, convective steering favors the nearest available fcc niches, with preference first for those near the top face of the crystal. Hence, the model predicts the potential for an overall tendency to nucleate new rows of particles at or near the menisci, although it also allows for twinning faults. More importantly, the results described in the previous section imply that an inclined planar interface would induce nonuniform solvent and particle fluxes near the growing crystal front. Influence of Growth Front Inclination: {111}, {100}, and {311} Planar Interfaces. Indeed, simulations of other growth front interfaces such as the {100} and {311} symmetry planes show velocity profiles with similar nonuniform trends as in the inclined {111} case (see Figure 5b and Supporting Information Figure S3); however, shelf planes (inclination angle above 90° with respect to the substrate) of the same crystal symmetry such as {311} induce velocity profiles that are nearly uniform. Plots of the solvent velocity profile, averaged over local variations in individual entrance channel resistance, are shown in Figure 7 for JE ) 1.0 cm/s and θ ) 0°. All simulations shown are based on 12 layer crystals with transverse dimensions of 220 particle diameters. These normalized plots are very similar when JE is greater or smaller than 1.0 cm/s, although the solvent velocity magnitude varies with the evaporation flux as discussed below. As the inclination angle relative to the substrate decreases, the solvent velocity profile becomes more nonuniform (Figure 7a). Conversely, as the growth front inclination increases relative to (50) Meng, L. Ph.D. Thesis, University of Minnesota, 2008. (51) Bryant, S.; Blunt, M. Phys. ReV. A 1992, 46, 2004–2011. (52) Zick, A. A.; Homsy, G. M. J. Fluid Mech. 1982, 115, 13–26. (53) Saeger, R. B.; Scriven, L. E.; Davis, H. T. J. Fluid Mech. 1995, 299, 1–15. (54) Davis, T. A.; Duff, I. S. ACM Trans. Math. Software 1999, 25, 1–20. (55) Kralchevsky, P. A.; Nagayama, K. Langmuir 1994, 10, 23–36.

Colloidal Crystal Growth

Langmuir, Vol. 24, No. 23, 2008 13691

Figure 7. Variation in entrance flux profile with growth front inclination relative to the substrate. Normalized solvent flux profiles are shown for entrance flow into faces of (a) {100} crystal symmetry and 54.7° inclination; (b) {111} crystal symmetry and 70.5° inclination; (c) {311} crystal symmetry and 80.0° inclination; (d) {311} crystal symmetry and 100.0° inclination where JE ) 1.0 cm/s. Vertical distance from the substrate is scaled by particle diameter and the local average solvent flux is normalized by the maximum solvent flux into discrete entrance channels (omitted for clarity).

the inclined {111} case (Figure 7b), the solvent velocity profile becomes more uniform (Figure 7c). As the inclination angle increases above 90°, a competition arises between mass transfer at the top edge and the nonuniform transverse pressure gradients imposed by the crystal face (Figure 7d). Mass transfer due to evaporation favors convective steering farthest from the substrate, by virtue of its proximity to that region of the crystal interface. However, the pore space traversed in a declining (or “shelf”) interface is largest in that region. As a direct consequence, evaporation competes with the imposed nonuniform transverse flow resistance. Local evaporative mass transfer is important when the ratio of transverse solvent flux to vertical solvent flux is not much greater than unity; when the evaporation rate is very small, a perpendicular growth front (inclination angle 90°) would yield the most uniform solvent and particle velocity profile. In either case, a nearly uniform velocity profile could develop near the exposed crystal facet. This result is important because a uniform velocity profile is expected to mediate steady growth front propagation. According to this criterion, the declining (“shelf”) front induces the most stable crystal growth. These simulations of {111}, {100}, and {311} symmetry planes also afford us a qualitative framework for describing the convective assembly process, although experimental evidence to date is not conclusive regarding the growth front structure. In

a previous study by Meng et al., the inclining {100} and {311} symmetry planes were observed by scanning electron microscopy of ex situ dry crystal interfaces, suggesting that these planes could be more stable than the {111} interface.26,50 It should be noted that the {100} and {311} planes support only one type of growth front niche, which means that they inherently propagate fcc structure as long as they are sustained as stable growth faces. However, the drying process could have introduced structural transitions from the in situ crystal interface. Additional visual measurements of the growth front interface as a function of distance from the substrate demonstrated a “shelf” shape when the substrate was immersed vertically in the suspension,50 but this result does not conclusively show that a {311} crystal plane was present in situ. This study provides one possible rationale for the preference toward the “shelf” crystal interface. If the {311} symmetry plane were indeed favored, as suggested by a highly simplified network analysis presented here, this conceptual framework could also help explain the strong tendency of convective assembly processing toward synthetic fcc crystals.

Outlook We have presented a network model of evaporation-driven creeping flow in close-packed crystals of monodisperse spheres. We considered the situation in which a free meniscus configures

13692 Langmuir, Vol. 24, No. 23, 2008

itself so as to pin the contact line at or near the leading crystal edge. In the “convective assembly” processes modeled here, particle motion is dominated by convective forces, and Brownian motion plays a minor role in colloidal crystal formation. In this work, we assert that capillary pressure and solvent transport in the crystal porespace follow from conservation of mass and interfacial mass transfer at exposed menisci. These assumptions amount to a restatement of the convective steering hypothesis, which posits that solvent flow into and through a colloidal crystal govern the assembly process. Central to this hypothesis is the competition between timescales for structural reorganization and convection; if structural organization is a slow process, then a quasi-rigid microstructure steers solvent flow near the crystal interface. Under these limiting conditions, a simple pseudosteady (creeping) flow model allowed us to simulate microstructure propagation in assembling colloidal crystals, and the hypothetical process was found to favor fcc structure propagation. Additional opportunities exist for study of defect formation and propagation or repair of point and line defects. These could be partly addressed by accounting for the geometry of, and constitutive relations for, flow in pore bodies and throats that arise near point or line defects. In any case, these defects are of high interest for photonic crystal applications.3,9,56 A more detailed approach would employ force balances about each of the spherical particles, thereby tracking the rearrangement and packing. One of the primary challenges would be to incorporate van der Waals, elastic, and lubrication forces into the model as particles come into contact. Second, this model would dynamically track a series of pseudosteady states and would require a time integration scheme to propagate meniscus rearrangement. Such a model would introduce additional complexities that we do not address here, but may explain such phenomena as sphere pack restructuring14 that are probably at least as important to the final crystal structure as is niche selectivity at the growth front. Any time-integrated model might in fact reveal interesting features of the growth dynamics. One might consider modeling the drying and air-liquid flow transport through the colloidal crystal assembly. This would require imposing special, “averaged” time step integrations over Haines jump events42 and meniscus rearrangement processes. Since these events are governed by minute configurational details that are not captured at this level, they must be dealt with as stochastic in nature. After the rearrangement process, the regular time integration scheme would be resumed. Other processes of interest, such as structural reorganization at the growth front region, could be simulated using Brownian dynamics57 or Stokesian dynamics58 algorithms. It should be stressed that this work represents an initial approach to modeling the complex network flow behavior in sphere packs of nearly monodisperse colloids. Additional refinements to account for the effects of polydispersity on packing, the description of Haines jumps and air invasion within the porespace, lubrication flow between spheres as they settle into crevices on the crystal front, and electrostatic and elastic forces between particles are outside the scope of the present work. We have shown, however, that a highly simplified model grounded in matrix analysis and topology and with close ties to modern circuit analysis can qualitatively describe aspects of the convective assembly process.

Conclusions A uniform solvent and particle flux near the growing crystal interface was proposed as the primary criterion for stability of (56) Taton, T. A.; Norris, D. J. Nature 2002, 416, 685–686. (57) Ermak, D. L.; Mccammon, J. A. J. Chem. Phys. 1978, 69, 1352–1360. (58) Brady, J. F.; Bossis, G. Annu. ReV. Fluid Mech. 1988, 20, 111–157.

Brewer et al.

the interface. Pore network simulations of solvent flow into and through crystals in which the growth plane assumed {111}, {100}, and {311} crystal symmetries of varied inclination relative to the substrate showed that a {311} declining “shelf” interface yielded a nearly uniform solvent flux profile. On the basis of this observation, we proposed that selective fcc structure propagation in the “convective assembly” process may be due to the stability of the {311} symmetry plane as the crystal interface, a growth facet that propagates only fcc microstructure. Additionally, it was found that solvent and particle flux profiles near alternative growth planes tended to direct crystal growth so as to generate the declining {311} interface. Hence, the mechanism proposed here to explain preferential fcc microstructure propagation was found to be consistent with the pore network simulations based on a stability argument. While visualization experiments provide clues as to the growing crystal interface, these simulations cannot be compared with experimental measurements due to the inherent difficulties in precise in situ structural characterization. This work represents a contribution to the understanding of a special class of porous media: regular porous networks in which the evolution of boundary conditions by addition of structural units to the solid matrix influences the flow into and through the porespace. Previous work has focused on describing transient processes with multiple-phase flow in such structures, but has not considered the special circumstance that arises when the boundary topology changes. In particular, the effect of mass transfer from curved liquid-vapor interfaces on convective steering mechanisms has not been studied systematically before. The present work is a first approximation to modeling colloidal crystal growth via convective assembly. Due to its formulation, however, the network model is quite general in treating ordered arrangements of spherical particles and could be extended to binary superlattice structures as well as irregular geometries and dynamic structure evolution. Acknowledgment. Support for this work was provided by NSF (CMMI-NIRT 0707610) and the Industrial Partnership for Research in Interfacial and Materials Engineering (IPRIME). Computational support from the Minnesota Supercomputing Institute (MSI) is gratefully acknowledged. D.D.B. gratefully acknowledges that this material is based upon work supported by a National Science Foundation Graduate Research Fellowship. Supporting Information Available: Additional information regarding the approximate relationship between meniscus curvature and shape in Purcell’s toroidal pore approximation; illustrations of {100}, {111}, and {311} crystal facets and growth front niches; illustrations of solvent flow near isolated nucleation sites and incipient point defects at the {111} planar interface. This material is available free of charge via the Internet at http://pubs.acs.org.

Glossary

Q r µ P A B Q P σ θ r R

Nomenclature Symbols external flow rate on liquid basis, m3/s effective orifice radius, m solvent viscosity, Pa-s absolute pressure, Pa linear admittance matrix, Pa-s/m3 boundary topology matrix vector of external flow rates, m3/s vector of absolute pressure, Pa surface tension, N/m contact angle geometric parameter, m radius, m

Colloidal Crystal Growth

R β δ RM H h L Ψ ε A

geometric parameter geometric parameter geometric parameter radius of curvature, m meniscus curvature, 1/m meniscus height, m radius of curvature at complete wetting, m meniscus filling angle transformed meniscus height mass transfer area of a meniscus, m2

Langmuir, Vol. 24, No. 23, 2008 13693

JE JT kG FG FL

evaporation flux on liquid volume basis, m/s transverse solvent flux on liquid volume basis, m/s gas-phase mass transfer coefficient, m/s gas-phase density, kg/m3 solvent density, kg/m3

0 h, k, j

reference value node indexes

Subscripts

LA802180D