722
J. Phys. Chem. B 2008, 112, 722-726
Deformation of a Sheet with Coupling between Elasticity and Concentration Constantine Khripin,† Tian Tang,‡ and Anand Jagota*,† Department of Chemical Engineering, Lehigh UniVersity, 111 Research DriVe, Bethlehem, PennsylVania 18015, and Department of Mechanical Engineering, UniVersity of Alberta, Edmonton, Alberta, T6G 2G8 Canada ReceiVed: July 27, 2007; In Final Form: September 24, 2007
We investigate a simple model for equilibrium deformation of a sheet with concentration-dependent elasticity. The model is motivated by several physical situations where deformation of a sheet is modulated by concentration of a mobile species, for example, a quasi-2D array of carbon nanotubes. Elasticity of the sheet is modeled using a free energy functional that includes concentration-dependent potential energies, and a free energy of mixing. We show that the sheet responds to imposed distortion by rearranging its constituents and can do so in a smooth and stable manner only for sufficiently small distortion. For larger distortions, through phase plane and numerical analysis, we consider how the system can meet boundary conditions through nonsmooth solutions.
1. Introduction Imagine an elastic body subjected to some combination of external loads and kinematic constraints in contact with a source or sink of material. We are interested in how the body’s response is modulated by its ability to readjust the local concentration of its constituents in response to the applied external conditions. Concrete examples include crystalline solids into and out of which matter can diffuse,1-6 focal adhesions in a cell membrane embedded with a variable concentration of transmembrane protein molecules,7,8 and a proposed quasi-two-dimensional (Q2D) sheet of DNA-coated carbon nanotubes9 (DNA-CNTs) in equilibrium with a dilute suspension of DNA-CNT rods. Two consequences of stress on the chemical potential of the constituent can be discerned. The first arises if its insertion or removal is accompanied by volume change but there is little or no change in elastic modulus. Two very diverse examples are the phenomena of diffusional creep1-6 and force-driven evolution of focal adhesion structure in cells.7,8 In such cases, under applied load there is no equilibrium and the body responds by creep. Under applied kinematic constraints (e.g., applied displacement), matter diffuses until stress is relaxed. The second consequence, which is what we analyze in this work, occurs if the elastic moduli of the material increase due to insertion of constituents (but strains are not relaxed). Since regions with higher concentration become stiffer, there is a tendency for lower concentration in regions where there are higher imposed strains. However, strains change in response and concentration gradients set up reverse chemical potential gradients. We wish to know whether the system will attain an equilibrium distribution of strain and constituent concentration, and what it looks like. We are motivated specifically by recent experiments9 that have reported alignment of negatively charged DNA-CNTs near a hydrophobic silicon substrate in a drop of dilute DNA-CNT suspension. It was proposed9 that a Q2D nematic phase of DNACNT forms in the liquid near the substrate surface. The authors proposed a model for this phenomenon in which the elastic * Corresponding author. Tel. 610 758 4396. Fax 610 758 5057. E-mail
[email protected]. † Lehigh University. ‡ University of Alberta.
distortion of the sheet comprising DNA-CNT rods is coupled to their concentration. In this work, we study a simple model for coupled elasticity and concentration in detail and examine what kinds of stable equilibrium deformations (and concentration profiles) are possible. We develop the model specifically for a body comprising rods. However, the formulation can be readily applied to other physical situations with a reinterpretation of the variables. Consider a sheet of aligned rods that has a preferred orientation with respect to the substrate (Figure 1). The existence of a preferred orientation is a common observation for rods in a nematic liquid crystal near a surface.10 If the orientation is further fixed at some locations to be different from its preferred value, for example, due to electric fields at electrodes in its plane, the sheet will be distorted as a result. The sheet is simultaneously in chemical equilibrium with a source/sink so that rods are free preferentially to leave regions of high distortion and/or misorientation energy. The elastic distortion or misorientation is therefore coupled with local concentration. We derive a differential equation governing the alignment both from a microscopic moment balance and from a free energy formulation based on electrostatic repulsion between the rods. The continuity and stability requirements for the solutions of the governing equation are evaluated. Through analytical, phase-plane, and numerical analysis of this model, we show that the system responds with smooth solutions only for sufficiently low distortion. At a critical distortion, the governing equations show a singularity; for larger distortions, we show that smooth solutions are unstable, but nonsmooth or even discontinuous solutions may be possible. 2. Governing Equation Consider a two-dimensional sheet of aligned rods with onedimensional variation of orientation and concentration, as shown in Figure 1. The state of the sheet is characterized by two variables, the alignment angle, θ(x), and the rod concentration per unit volume, c(x). The rods in this plane are in contact with the bulk suspension with a much smaller concentration, and equilibrium between the two phases is assumed to be achieved rapidly. The rods are not adsorbed on the surface; they interact
10.1021/jp075965k CCC: $40.75 © 2008 American Chemical Society Published on Web 12/18/2007
Sheet Deformation
J. Phys. Chem. B, Vol. 112, No. 3, 2008 723
Figure 1. A quasi-2D sheet of rods deforms elastically in response to boundary conditions and has a preferred orientation θo with respect to the substrate. The state of the sheet is characterized by concentration c(x) and alignment θ(x). Moment balance on the distorted sheet relates body torque τ(x) and gradient of the bending moment M(x). Note that θ(x) and c(x) vary only in the x direction.
with the surface via electrostatic and van der Waals forces. The rods form an elastic sheet in which the concentration of rods can change. In the case of negatively charged DNA-CNT rods, this is because of mutual electrostatic repulsion that, at sufficiently high concentration, results in alignment of the rods. We now derive the differential equation describing the behavior of the sheet. We do this via two independent approaches, one by considering the torque balance in the sheet, and the other by proposing a free energy functional and finding its extrema. We show that the two approaches are equivalent and result in the same governing differential equation. 2.1. Torque Balance in the Sheet. When the local orientation, θ(x), is not aligned with the preferred direction, θo, the substrate exerts a restoring body torque τ(x), in the unit of moment per unit length, on the sheet. Moment balance of an element of the sheet with width dx, located at x (Figure 1), then yields
dM(x) ) τ(x) dx
(
)
f ) fsurface + fdeformation - TSmixing
(3)
Here fsurface ) c(x) µs where µs is the free energy of interaction between a single rod and the surface when θ ) θo, i.e., the standard chemical potential. The contribution due to deformation, fdeformation, has terms due to distortion and misorientation, equivalent to work done by the bending moment and body torque in eq 1. In other words, fdeformation is the sum of the elastic bending energy and the potential of the load by the substrate. Smixing is the entropy of mixing for rods in the sheet. A more detailed derivation of the distortional contribution to the free energy for the case of rods with mutual electrostatic repulsion is given in the Supporting Information. For small departures from the undeformed state (where θ,x ) θ - θo ) 0), the free energy density is
f ) µsc +
(ko + ck1) 2 k2 c θ,x + c sin 2(θ - θo) + kbTc ln 2 2 co (4)
(1)
where M(x) is the bending moment which is needed to balance the torque. If the body torque τ can be determined as a function of the orientation θ, and the constitutive relation between moment M and θ is specified, eq 1 becomes a differential equation for θ(x). Clearly, τ ) 0 when the rod is in its preferred orientation, i.e., θ ) θo + nπ (n ) 0, 1, 2, ...), and it is a periodic function in θ with a period of π. In addition, τ should be an odd function about θo. A simple choice for such a function is τ ) Kτ sin 2(θ - θo), where Kτ is a modulus representing stiffness against misorientation, which has the unit of moment per unit length and may depend on the local concentration. Moments in the sheet induce distortion in the form of gradient of orientation, θ,x. If one assumes a linear relation between the bending moment M and the distortion θ,x, i.e., M ) KMθ,x, where KM is a modulus representing the distortional stiffness of the sheet, then, eq 1 becomes
dθ d K ) Kτ sin 2(θ - θo) dx M dx
We now consider Kτ and KM as functions of the local concentration c. Kτ must incorporate the number of rods per unit volume, because the torque is exerted on each rod. Thus we take Kτ ) kτc, where kτ is a constant and c is the number of rods per unit volume. The distortional modulus KM also depends on concentration and in this work we take the relationship between the two to be linear: KM ) ko + k1c. This is discussed in more detail later in the paper and in the Supporting Information. Because we are interested in the situation where concentration adjusts in response to applied deformations, c(x) is unknown. It is specified in the next section by requiring uniform chemical potential for the rods. 2.2. Free Energy Functional. We propose a free energy functional for the DNA-CNT in the sheet. We take free energy density (energy per unit volume) in the sheet, f, to have three components
(2)
Note that KM has the dimension of force times square of length, the same as the moment of inertia for bending of an elastic beam.
where co is the maximum rod concentration; ko, k1, and k2 are material constants independent of concentration; kb is the Boltzmann constant; and T is temperature. Here, ko has the unit of energy per unit length, k1 has the unit of energy times square of length, and k2 has the unit of energy. In the case of rods repelling each other electrostatically, the term ko + ck1 is the distortional modulus, linearized about the concentration at zero distortion and misorientation. We can draw an analogy of this energy formulation with that of a strained composite sheet composed of an elastic matrix with a dilute concentration of rigid inclusions, if the kinematic variable θ is interpreted as a displacement.11 In this case, c represents the concentration of inclusions and ko represents the modulus of the matrix, while k1 represents the initial rate of increase of composite modulus with concentration of inclusions. In terms of examples cited in the Introduction, the composite sheet could be a cell membrane with protein molecules as inclusions.7,8 Specifically in the case of rods, the first term in eq 4 is the free energy density due to rod/surface interaction. The second and third terms are the internal strain energy (work of moment acting on distortion, θ,x) and misorientation energy (work of torque acting on misorientation θ - θo), respectively. Interestingly, the second and third terms have the same forms as the Frank-Oseen and the Landau-de Gennes free energies of a liquid
724 J. Phys. Chem. B, Vol. 112, No. 3, 2008
Khripin et al.
crystal, respectively. Specifically, the Frank-Oseen free energy is given by10,12
d dθ (k + ck1) ) k2c sin[2(θ - θo)]/2 dx o dx
1 1 2 2 fFO ) fFO o + K1(∇•n) + K2(n•∇ × n) + 2 2 1 K (n × ∇ × n)2 (5) 2 3
which is identical to eq 2 derived from moment balance. Writing eq 10 using normalized variables and substituting eq 9 eliminates c(x) as a variable resulting in the differential equation
fFO o
where n is the orientation vector of the liquid crystal, is a reference energy, and K1, K2, and K3 are elastic constants that may depend on concentration. For a two-dimensional liquid crystal, n is confined to a plane, i.e., n ) nxi + nyj ) cos θi + sin θj, where θ ) θ(x) is the orientation. Using this and 2 assuming K1 ) K3, eq 5 becomes fFO ) fFO o + K1θ,x /2, which has the same contribution as the second term in eq 4 if one identifies the Frank constants K1 ) K3 ) ko + ck1. The Landaude Gennes free energy for a liquid crystal can be written as10
fLD ) -
fLD o [3 cos2(θ - θo) - 1][3〈cos2(θ - θo)〉 - 1] 4 (6)
where 〈•〉 indicates the expectation value. Unlike the contribution from distortion, this depends on the local values of orientation, and so we may estimate it using the local value of θ(x). 2 Therefore, the orientational free energy becomes -fLD o [3 cos 2 (θ - θo) - 1] /4, which, neglecting the higher-order terms, contributes a sin2(θ - θo) term identical to the third term in eq 4 if fLD o ) k2c/6. The fourth term in eq 4 is the free energy of mixing. It should be noted that the layer of aligned rods exists in the liquid at some distance from the substrate. They are free to move between the bulk solution and the aligned phase. 2.3. Differential Equation from Free Energy Functional. We first impose the condition that, under equilibrium, the chemical potential for rods is uniform throughout the sheet and the bulk solution, i.e.,
∂f/∂c ) constant
(7)
Equation 4 describes f in the sheet; in the bulk, we take f ) kbTcb ln(cb/co) with cb the concentration in the bulk. Using these results and eq 7, we obtain
kbT ln
cb k1 k2 + kbT ) µs + θ,2x + sin2(θ - θo) + co 2 2 c kbT ln + kbT (8) co
( )
which gives
(
θ˜ ,2x˜ k˜2 2 c˜ ) exp - sin θ˜ 2 2
)
(9)
where we have introduced the normalization x˜ ) xxkbT/k1, θ˜ ) θ - θo, c˜ ) c/c∞, and k˜2 ) k2/(kbT). Here, c∞ ) cb exp(µs/kbT) is the concentration in the sheet far from any boundaries where the rods align with the preferred orientation and there is no distortion, i.e., θ˜ ,x˜ ) θ˜ ) 0. The differential equation governing alignment, θ, eq 2, can also be derived from the above free energy formulation. Specifically, we impose the condition that at equilibrium θ(x) should be the variational extremum of the total free energy F ) ∫L0 f dx. The Euler-Lagrange equation13 ∂f/∂θ - d(∂f/∂θ,x)/ dx ) 0 becomes
[
]
θ˜ ,x˜x˜ )
k˜2 1 + θ˜ ,2x˜ sin(2θ˜ ) 2 1 - θ˜ ,2x˜ - k˜o/c
(10)
(11)
where k˜o ≡ ko/[k1 ln(cb/co)]. In this work, we focus on the case where the constant term in the linear relationship between distortional modulus and concentration is small, i.e., ko , k1c. With this, eq 11 simplifies to
1 + θ˜ ,2x˜ k˜2 θ˜ ,x˜ x˜ ) sin(2θ˜ ) 2 1 - θ˜ ,2
(12)
x˜
When θ˜ ,x2˜ , 1, the equation yields smooth solutions, and if in addition θ˜ (x˜ ) is small, eq 12 reduces to θ˜ ,x˜ x˜ ) k˜2θ˜ , which yields solutions in terms of exponentials that decay away from the boundaries with a characteristic decay distance of xk1/k2. For θ˜ ,x2˜ . 1, there is another set of smooth solutions that, for small θ˜ (x˜ ), would be in terms of sinusoidal functions. We further study the behavior of the solution to eq 12 via the phase plane analysis presented below. 3. Phase Plane Analysis At this stage, it is helpful to examine the behavior of eq 12 on its phase plane (shown in Figure 2), where trajectories that solutions to eq 12 would follow are plotted, given some initial values of θ˜ (x˜ ) and θ˜ ,x˜ (x˜ ).14,15 Black arrows indicate the direction in which a solution will evolve with increasing x from any given initial condition. Blue lines are the contours of such solutions. Thus, if our boundary conditions are defined by two different values of θ˜ , the set of blue lines which connect the boundary conditions is the set of possible smooth solutions to the problem. Continuous but nonsmooth solutions will be represented by vertical jumps from one blue line to another. Suppose that two boundary conditions on θ˜ are imposed, distorting the sheet, as in Figure 1. Let the angle on the left, θ˜ L, be given by vertical green lines 2 or 3 in Figure 2 and the one on the right, θ˜ R ) 0, by vertical green line 1. Suppose that the system is large, |θ˜ R - θ˜ L|/L˜ , 1, where L˜ ) LxkbT/k1. Then, far from any boundary, θ˜ ,x˜ ) θ˜ ) 0; this is represented by the red filled circle. If the magnitude of |θ˜ L| is smaller than a critical value θ˜ c, then we have a smooth solution within the red lines. For example, if - θ˜ c < θ˜ L < 0 (line 2), the lower limit being represented by the red filled square, then the solution follows the purple contour a in Figure 2, and θ˜ ,x2˜ remains less than unity everywhere. For example, θ˜ c ) arccos(1 - 2(2 ln 2 - 1)/k˜2)/2 ≈ 0.67 for k˜2 ) 1. This represents a decaying solution that is exponential in x˜ for small θ˜ . Now, suppose we decrease θ˜ L to green line 3; we see that there is no longer a smooth solution that connects the boundary condition on the left with the condition far from the boundary. If the left boundary starts with |θ˜ ,x˜ | > 1, there may be a smooth solution connecting the two boundaries such as the one shown by purple contour c; such solutions are quasi-sinusoidal. In other cases, it will need to cross the singular boundary shown by red lines, e.g., purple line b, which represents a discontinuity in the gradient of the orientation. It is also possible that there will be solutions like d that contain a discontinuity in gradient
Sheet Deformation
J. Phys. Chem. B, Vol. 112, No. 3, 2008 725
Figure 2. Phase plane diagram for eq 12 with k˜2 ) 1. Black arrows indicate the direction of solution trajectories. The blue lines are possible solution paths. The red lines at θ˜ ,x2˜ ) 1 are the singularities that smooth solutions cannot cross. Smooth solutions that remain between the red lines represent variational energy minima, while those outside the lines are variational energy maxima. The red filled circle indicates the state far from any boundary. The red filled square indicates the maximum value of θ˜ that could be imposed on a boundary and still be consistent with a smooth decay to the conditions far away. This phase plane is periodic in θ˜ with a period of π.
but do not cross the singular boundary. Note also the possibility of oscillating solutions near |θ˜ | ) π/2; θ˜ ,x2˜ < 1, represented by the closed loops. If the rods are forced to be in a mean orientation perpendicular to the one they prefer, the orientation will oscillate in space. 3.1. Discontinuities in θ˜ and θ˜ ,x˜ . It can be shown that, if the solution contains a discontinuity in θ˜ , the total free energy F ) ∫L0 f dx is unbounded for finite c˜ (x˜ ). However, a discontinuity in θ˜ ,x˜ , or a “break point”, is allowed. Equating chemical potential on two sides A and B of the break point (see b and d in Figure 2) implies that θ˜ 2A,x˜ - θ˜ 2B,x˜ ) 2kbT ln(c˜ B/c˜ A), i.e., the jump in slope is balanced by a jump in concentration. 3.2. Stability of the Governing Euler-Lagrange Equation. To determine the stability of the solutions to the EulerLagrange eq 12, we consider whether they correspond to energy minima or maxima. Using the necessary condition for a relative minimum (maximum),13 ∂2f(θ˜ , θ˜ ,x˜ )/∂θ˜ ,x2˜ < 0(∂2f(θ˜ , θ˜ ,x˜ )/∂θ˜ ,x2˜ > 0), we find that a necessary condition for the stability of eq 12 is that θ˜ ,x2˜ < 1 everywhere. Quasi-sinusoidal solutions, like the trajectory c drawn in Figure 2, are therefore not stable. Once one allows breaks, the number of possible solutions increases dramatically. To examine whether such solutions can be stable, we restrict our attention in the next section to cases with a single break point, and search for the set of allowed solutions by numerical analysis of eq 12 in the two resulting domains. 4. Numerical Analysis We divide the domain (0, L˜ ) into two regions and look for smooth solutions in each, while allowing for a discontinuity in slope at the domain boundary. At the right boundary, we fix θ˜ R ) 0. This problem is fully specified by the location of the break point, x˜ m, and the value of angle, θ˜ m, at that point. In the numerical analysis, we vary x˜ m ∈ (0, L˜ ) and for each x˜ m, the lower and upper limits for θ˜ m are sufficiently large to capture a physically relevant range. On each side of the break point, a numerical solution of eq 12 was obtained using the finite
Figure 3. Energy map of numerical solutions with one break point (k˜2 ) 1). Each point corresponds to the location of and angle at a break point, (x˜ m, θ˜ m), and the color indicates the free energy of that solution. Equation 12 was solved numerically on each side of the break point, and eq 4 was used to calculate the free energy of the sheet. The areas of solid color (no gradients) are where no solution is possible with one break point. (a) θ˜ L ) -0.15: For any choice of x˜ m, the minimum energy solution corresponds to the case where θ˜ m is given by the smooth solution obtained by modeling (0, L˜ ) as a single domain. All the twodomain minimum energy solutions (the white line) are identical to each other and to the single-domain solution. (b) θ˜ L ) -0.8: Distortion is strong enough so that there is no smooth solution. The energy minima (blue regions) correspond to solutions with a discontinuity in θ˜ at the boundaries. (c) θ˜ L ) -2: Energy-minimizing solutions are as in case b. There is a smooth solution to eq 12, the white line, but it maximizes the free energy.
difference method with Newton-Raphson iterations on the resulting nonlinear algebraic equations. For each solution, we calculate the total free energy. Now, we can create a twodimensional contour plot of free energy in the space of x˜ m and θ˜ m (Figure 3, k˜2 ) 1). The color assigned to a point represents free energy for the solution associated with it. Thus, the contour plot studies the entire family of solutions with one break point for specified boundary conditions, and the color in the plot tells us which of these solutions are lower or higher in energy. For
726 J. Phys. Chem. B, Vol. 112, No. 3, 2008 some combinations (x˜ m, θ˜ m), there are no smooth solutions on one or both sides of x˜ m. Color in these regions appears to have no gradient. From among all the broken solutions with the same x˜ m but different θ˜ m, we plot the one that has the minimum energy (superimposed white lines in Figure 3a and black lines in Figure 3b,c). In Figure 3a, θ˜ L ) -0.15 and the distortion is small. One anticipates that the simple decaying solutions, like a in Figure 2, will be the energy minimum. This is confirmed by the fact that for any choice of x˜ m the minimum energy solution corresponds to the case where θ˜ m is given by the smooth solution obtained by modeling (0, L˜ ) as a single domain. As a result, all the two-domain minimum energy solutions (the white lines in Figure 3a) are identical to each other and to the single-domain solution. When θ˜ L ) -0.8 (Figure 3b), a smooth solution connecting the two sides is no longer possible. There are two local minima, one corresponding to a jump in θ˜ on the left boundary, and the second one corresponding to a jump at the right boundary. That is, if only one break point is allowed, the system minimizes free energy by making angles uniform everywhere except very close to the boundaries where it accommodates the abrupt jump in angle with vanishingly small concentration. Figure 3c represents the case of line c in Figure 2. In this case, distortion is very high. The solution that would appear to connect the boundary conditions in a smooth manner, highlighted by the white line in the figure, is actually an energy maximum, as evidenced by the contours in the plot. This is not a stable solution. Note that our analysis does not rule out the possibility of having stable solutions with multiple break points. Also, under other conditions, one can obtain stable solutions with break points removed from the boundaries like the one labeled d in Figure 2. 5. Conclusion and Discussion We have analyzed a simple model for the deformation of a sheet with concentration-dependent elasticity. The sheet is free to adjust the concentration of its constituent particles to maintain constant chemical potential. We have been motivated specifically by experimental observations of aligned DNA-CNT rods near a substrate. A free energy functional is formulated on the basis of the electrostatic repulsion between the rods, and it resembles the free energy for a continuum liquid crystal. If the elastic moduli are taken to be linear in concentration, smooth and stable solutions are possible only for sufficiently small distortions. A necessary condition is that the average slope, (θ˜ R - θ˜ L)/(x˜ R - x˜ L), be less than unity or that the gradient of normalized angle be less than unity everywhere. Nonsmooth “broken” solutions are allowed, and some combination of these can be stable. If only a single break is allowed, when boundary conditions exceed values where smooth, stable solutions result, we find that the minimum free-energy solutions have jumps in angle at the boundaries. Break points are possible because the energy formulation permits discontinuities in θ˜ ,x˜ (x˜ ) and c˜ (x˜ ). To explore the physical significance of the singularity at θ˜ ,x2˜ ) 1, consider first the limit |θ˜ | and θ˜ ,x2˜ , 1. In this limit, solutions decay exponentially from the boundaries and concentration is nearly uniform everywhere. In this case, where the Euler-Lagrange equation reads c∞k1θ,xx ) c∞k2(θ - θo), we
Khripin et al. can interpret c∞k1 and c∞k2 as the moduli that resist distortion and misorientation, respectively. The form of the EulerLagrange equation in the limit θ˜ ,x2˜ . 1 is (-ck1)θ,xx ) ck2 sin{2(θ - θo)]/2, from which it is clear that, effectively, one of the moduli is now negative and therefore the system must be unstable. It can be shown that the effect of nonuniform concentration is to reduce the distortional modulus by a factor (1 - k1θ,2x /kbT) so that the critical condition corresponds to the case where it vanishes in comparison with the misorientation modulus. The singularity we have found through analysis of our model should be interpreted as showing a tendency of the sheet as it is subjected to distortion. Our analysis is strictly valid either if the distortional modulus is linear in concentration or, in the more general case of nonlinear dependence on concentration, for small departures from the undistorted state. Inclusion of other terms, such as the effect of the constant contribution to the elastic modulus, ko, will very likely be needed to for more accurate analysis for large distortions. Acknowledgment. This work was supported in part by the National Science Foundation under grant CMS-0609050, and by NASA under award no. NNX06AD01A for the Lehigh University/Mid-Atlantic Partnership for NASA Nanomaterials. T. Tang acknowledges support from the Canada Research Chairs Program and National Science and Engineering Research Council (NSERC). The authors would like to acknowledge useful discussions with Prof. Philip Blythe, Prof. C-Y. Hui, Prof. M. K. Chaudhury, Dr. Ming Zheng, and Dr. N. J. Glassmaker. Supporting Information Available: Derivation of eq 4; an argument for how the electrostatic interaction between a pair of charged rods translates into the free energy functional, and the assumptions which are made in this derivation; derivation of the concentration eq 9. This material is available free of charge via the Internet at http://pubs.acs.org. References and Notes (1) Nabarro, F. R. N. Report of Conference on Strength of Solids; Physics Society: London, 1948. (2) Herring, C. J. Appl. Phys. 1950, 21, 437. (3) Coble, R. L. J. Appl. Phys. 1963, 34, 1679. (4) Nabarro, F. R. N. Metall. Mater. Trans. 2002, 33A, 213-218. (5) Asaro, R. J.; Tiller, W. A. Metall. Trans. 1972, 3, 1789-1796. (6) Grinfeld, M. A. J. Nonlin. Sci. 1993, 3, 35-83. (7) Shemesh, T.; Geiger, B.; Bershadsky, A. D.; Kozlov, M. M. Proc. Natl. Acad. Sci. U.S.A. 2005, 102, 12383-12388. (8) Nicolas, A.; Geiger, B.; Safran, S. A. Proc. Natl. Acad. Sci. U.S.A. 2004, 101, 12520-12525. (9) McLean, R. S.; Huang, X.; Khripin, C.; Jagota, Zheng, A. M. Nano Lett. 2006, 6, 55-60. (10) de Gennes, P. G.; Prost, J. The Physics of Liquid Crystals, 2nd ed.; Oxford University Press: Oxford, 1993. (11) Hashin, Z. J. Appl. Mech. 1983, 50, 481-505. (12) Oseen, C. W. Trans. Faraday Soc. 1933, 29, 883-899. (13) Sagan, H. Introduction to the Calculus of Variations; Dover Publications: Mineola, NY, 1969. (14) Arnold, D.; Polking, J. C. Ordinary Differential Equations Using Matlab; Prentice Hall: Upper Saddle River, NJ, 1999. (15) MATLAB V.7.0 by Mathworks; applet pplane7.m by Polking, J. C., Rice University, Department of Mathematics, http://math.rice.edu/∼dfield/ dfpp.html.
