Conformations of a Dipolar Solute in a Stockmayer Solvent Channel

Department of Mechanical Engineering, Northwestern University, Evanston, Illinois 60208, United States. Langmuir , 0, (),. DOI: 10.1021/la303252c@proo...
0 downloads 0 Views 3MB Size
Article pubs.acs.org/Langmuir

Conformations of a Dipolar Solute in a Stockmayer Solvent Channel Taeil Yi,* Qian Wang, and Seth Lichter* Department of Mechanical Engineering, Northwestern University, Evanston, Illinois 60208, United States S Supporting Information *

ABSTRACT: A wide range of molecules, from inorganic to biological, self-assemble on surfaces. Previous studies have elucidated many features of solute reorganization on surfaces using coarse-grained modeling, implicit solvents, and constraints such as chemically bonding the solute to the surface. Using molecular dynamics simulations under various combinations of interaction parameters, solute fractions, and solute dipole moment, we study the redistribution of freely-rotating dipolar solute molecules solvated by a water-like Stockmayer solvent initially adhered to a face-centered cubic substrate. The balance of attractive and repulsive forces is essential for acquiring a particular stable conformation. Here we show that the adsorbed molecules redistribute into different conformationswetting film, nonwetting, partial wetting, and pseudopartial wetting dropsdepending on the parameter values. We observe that the pseudopartial wetting drop is transient and its rate of spreading fluctuates, slowing to nearly zero as it passes through particular conformations before reaching an equilibrium thin film. Strong attraction between solute molecules yields a droplet with a net dipole moment. A high solute fraction leads to a pancake-like conformation arising from a balance of surface tension and van der Waals forces. This study augments our understanding of the evolution of aggregates in biological systems and also the design of polymers for self-assembled monolayers for industrial applications.



film, a partial wetting drop, and a pseudopartial wetting drop.10 The first two terms are in common usage, and the third refers to the case that is a hybrid of the first two in which a drop coexists with a film.11 Nonwetting drops are spherically shaped under zero gravity with a contact angle of 180°. One approach to determine the drop profile and test its stability is to solve the Young−Laplace (Y−L) equation augmented with disjoining pressure and/or minimize the free energy of the drop.12−15 The singularity of the moving contact line of the augmented Y−L equation is avoided by using the disjoining pressure to incorporate intermolecular interactions.16,17 These continuum models, based on mean-field theories and a sharp Gibbs dividing surface, may not be suitable for providing detailed information about nano- and atomic-scale phenomena. Molecular simulation is a technique to support experimental and theoretical studies by providing detailed atomic-level behavior and structure.18 There is a vast literature on assemblages of molecules on a substrate, such as water droplets on either hydrophobic or hydrophilic substrates19−22 and polymers on a substrate.23 SAM simulations assume chemiosorption between a solute and substrate.24,25 Other studies investigated single protein adsorption onto a hydrophobic surface in water with no prebuilt chemical bond between the protein and the surface.26 Xue et al. discussed conformational

INTRODUCTION The distribution and conformational changes of molecules such as proteins, DNA, and inorganic polymers as they adhere onto, detach from, and reorganize on a substrate have been studied using computational simulations,1,2 theoretical analysis,3 and experiments.4,5 Advances in experimental techniques to construct and probe nano- and atomic-scale systems have allowed theoretical and simulation expectations to be compared with observed properties.6 The self-assembled monolayer (SAM) has been shown to be an experimentally useful and application-prone realization using modern technology to achieve well-ordered molecular and atomic-scale structures.7 However, it is still not easy to physically manipulate specimen samples for systematic studies. On the theoretical side, coarsegraining methods such as mean-field theory, while admitting their short-comings such as an ignorance of fluctuation offsets,8 have been successful. Starting in the early 19th century, macroscopic studies of drop conformation used concepts of surface energy.9 Many equilibrium drop shapes are possible, depending on the contact angle, defined as the angle between the substrate surface and the drop interface at the contact point. de Gennes classified drop shapes as total wetting for zero contact angle and partial wetting for positive contact angle.3 Brochard-Wyart et al. used the Hamaker constant, which represents the van der Waals attraction, and the spreading coefficient, defined by the surface energy difference between the wet and dry substrates, to introduce three types of droplet shapes: a completely wetting © 2012 American Chemical Society

Received: August 10, 2012 Revised: September 23, 2012 Published: October 11, 2012 15286

dx.doi.org/10.1021/la303252c | Langmuir 2012, 28, 15286−15293

Langmuir

Article

changes of a chain polymer depending on the number of monomers.27 Investigations typically focus on a specific molecule or molecule−substrate pair. Often, especially for simulations of water droplets, only vacuum is used as the environment. For studies of polymers on surfaces, the most frequent choice has the solute chemically bonded to the substrate. Recently, studies of physisorbed monolayers at liquid−solid interfaces using modern microscopy have been introduced to investigate how factors such as chirality affect solute packing.28 Much remains to be learned about how flexible molecules organize, adhere, and detach from a substrate to which they are physisorbed. Our simulations consider a dipolar solvent containing a dipolar freely-rotating chain solute initially physiosorbed on a substrate. Inter- and intra-interaction energy parameters of the solute, solvent, and substrate are varied, along with the number of solute molecules and solute dipole strength. Molecular geometries and interactions are kept simple, so the system is described by a small number of parameters and thus a systematic parameter exploration is feasible. Furthermore, the simplicity facilitates finding the relationship between the observed conformations and the underlying interactions. Despite the choice of simple molecular geometries and interactions, the conformations are similar to those seen in real experiments. Therefore, the model systems capture fundamental aspects of different conformations and can be used to identify critical interactions needed to achieve particular arrangements of solute and solvent molecules.

Figure 1. (a) The simulation cell is bounded on top (purple) and bottom (gray) by an FCC crystal of area 1600. Periodic boundary conditions are applied to the unit cell in the x and z directions, whereas in the y direction the periodicity is enforced on the unit cell plus a buffer region (bounded in blue). The long buffer region prevents the electrostatic field above the top (bottom) boundary from propagating through the periodic y boundary condition to the bottom (top) boundary. The unit cell is shown with Stockmayer solvent (mauve) and 400 solute molecules (blue) from a snapshot of one of our simulations. (b) Initial position of 400 solute molecules (blue) on a substrate, from a snapshot before starting the simulation.



SIMULATION SETUP, METHODS, AND PROCEDURES System Configuration. Quantities in this article are normalized by the fundamental units, energy εll, length σll, and mass ml shown in Table 1, where the subscript l represents

z component of the solute dipole moment and N is the number of solute molecules), is extremely weak. When nondimension−3 −2 alized as shown in Table 1, Vcore C ≈ 6 (10 −10 ). Hence, the dipole charge induces a field whose effects are much weaker than those due to the Lennard-Jones interactions. Number density ρs = 4.0 for the substrate, and the number densities of the solute and solvent vary with the number of solute molecules used; for example, for 400 solute molecules, ρl = 0.73 and ρp = 1.18. Each solute chain is initialized with its long axis normal to the substrate with a square lattice spacing of unity (Figure 1b). Solvent. The Stockmayer solvent is a Lennard-Jones sphere combined with a point dipole at its center of mass.29−31 This simple solvent model can be used to study the characteristics of polar molecules such as water and ammonia. For this study, the Stockmayer solvent can be considered to be a coarse-grained water molecule that retains both the pairwise interaction and dipolar characteristics but loses information about hydrogenbond geometry. The Stockmayer solvent’s Lennard-Jones parameters are taken to match the TIP4P water model. The phase diagram for a Stockmayer solvent with this dipole moment strength is used to select the thermostat temperature of 2.5, which gives a liquid phase at a solvent density of ρl ≃ 0.7.32,33 The initial orientations of the Stockmayer dipole moments are randomly assigned from a uniform distribution of direction cosines in the x, y, and z directions. Typically, there are approximately 19 000 solvent molecules. Solute. The solute is coarse grained using a united-atom model. A finite extensible nonlinear elastic (FENE) potential is used for the chemical bonding between the four monomers of the freely-rotating chain solute using Kremer and Grest parameters, (KFENE, RFENE) = (30, 1.5).34 The head monomer is positively charged +|q|, the tail is negatively charged −|q|, and the two intervening monomers are uncharged. For a water

Table 1. Variables Presented in the Article Are Nondimensionalized with Respect to Characteristic Energy, Mass, and Size of the Solvent, {εll, σll, ml}, Respectivelya quantity

normalization factor

energy E (or U) length L mass M charge q temperature T time τ |dipole moment| |p⃗| translational entropy Syθ

εll σll ml (4πε0εrσllεll)1/2 εll/kB (mσll2/εll)1/2 (4πε0εrσll3εll)1/2 kB

kB, ε0, and εr are the Boltzmann constant, vacuum permittivity, and dielectric constant, respectively. a

the solvent (liquid) and subscripts p and s will subsequently be used for (polymer) solute and substrate components, respectively, and εij and σij denote the Lennard-Jones energy and length scale between components i and j, where i, j = l, p, s. Figure 1a shows the simulation geometry of the Stockmayer solvent plus the dipolar linear-chain solute confined between parallel substrates. The size of the periodic unit cell is (Lx, Ly, Lz) = (40, 250, 40). A vacuum slab is used to reduce the interference across the periodic boundary conditions in the y direction. Its length is chosen to be long enough such that the average y component of the total dipole moment of the N i i simulation cell, Vcorr C ≡ −2π/(4NLxLyLz) ∑i=1pz (where pz is the 15287

dx.doi.org/10.1021/la303252c | Langmuir 2012, 28, 15286−15293

Langmuir

Article

where p⃗i is the dipole moment vector of solvent molecule i. Parameters values are based on those used for similar polymer molecules.35−37 Simulation Protocol. Simulations were implemented using the large atomic/molecular massively parallel simulator (LAMMPS), which solves the equations of motion of Newtonian mechanics using the Verlet algorithm.38 Visualizations in this section were produced by the visual molecular dynamics (VMD) package.39 We consider an NVT ensemble in which the system is in contact with a heat reservoir at constant temperature, the system volume is constant, and the number of each constituent (solute and solvent) is fixed. The Berendsen thermostat40 is applied to regulate the solvent temperature by solving the equations of motion as

solvent at room temperature (300 K) and a dielectric constant of εr ≈ 80, |q| = 1 is equivalent to 0.45 times the elementary charge, which is approximately the charge assigned to the hydrogen of the TIP4P water model. Hence, the values used, |q| = 0.2, 1, 2, include those relevant to water and other polar solvents. Additionally, each of the monomers is given a Lennard-Jones interaction with other monomers, the solvent, the solute, and the substrate. Substrate. Atoms of the upper and lower substrates confining the mixture are assigned to fixed positions on a face-centered cubic (FCC) lattice in the xz plane. The distance between the centers of mass of the two substrates is 35.



INTERACTIONS Interactions among the solvent, solute, and substrate are modeled by a combination of nonbonding pairwise potentials. Superimposed on the Lennard-Jones interaction, we add a dipole−dipole potential for solvent−solvent pairs resulting from their permanent dipole moment; a Coulomb potential between paris of partial charge on the head and tail of the solute; and a charge-dipole potential between solute partial charges and the solvent permanent dipole moment. The Lennard-Jones potential VLJ describes steric repulsion and longer-range attraction by ⎧⎛ σij ⎞12 VLJ(r ) = 4ϵij⎨⎜ ij ⎟ − ⎩⎝ r ⎠ 12 ⎧ ⎪⎛ σij ⎞ ⎨ − 4ϵij ⎜ ⎟ ⎪ ⎩⎝ rc ⎠ ⎪

ij



VLJ(r ij) = 0,

⎛ σij ⎞6 ⎫ ⎜ ⎟ ⎬ ⎝ r ij ⎠ ⎭





dvi⃗ F⃗ T −T =− i − B vi⃗ dt mi ξT

where ν⃗i is the velocity vector, F⃗i is the force vector given by −∇Vtot,i, where Vtot,i is the total potential energy of the ith atom, T is the instantaneous temperature of the system, and the damping coefficient ξ = 0.5 determines the relaxation time. The solute temperature is regulated by a Langevin thermostat, dvi⃗ F⃗ = − i − γvi⃗ + R i dt mi

ij

r < rc

(1)

where rij is the distance between monomers i and j, rc = 2.5 is the cutoff distance, εij is the strength of the Lennard-Jones interaction, and σij is the separation distance between monomers i and j for which VLJ(σij) = 0. Subscripts continue to represent the component {l, p, s}, and each superscript for the distance rij (and subsequently for qi and p⃗i) identifies each solvent atom, solute monomer, or substrate atom. For example, rij might be the distance between head monomer i of one of the N solute molecules and any one j of the thousands of Stockmayer atoms. Interaction parameters between the solvent and solute are computed using the Lorentz−Berthelot combining rules, εij = (εiiεjj)1/2, and σij = (σii + σjj)/2, where i ≠ j. The Coulomb potential VC between charges qi and qj is given by VC(r ij) =

qi · q j r ij

(2)

The particle−particle particle-mesh (PPPM) algorithm is used to calculate the electrostatic potential. Charge−dipole (VCD) and dipole−dipole (VDD) interactions for the interactions between the solute charge and solvent dipole and between the solvent dipoles, respectively, are calculated by18

νk =

(p ⃗ j · r ⃗ ij) (r ij)3

VDD(r ij) =

1 3 p′ i ·p ⃗ j − ij 5 (p ⃗ i · r ⃗ ij)(p ⃗ j · r ⃗ ij) ij 3 ⃗ (r ) (r )

1 Nj 1 Nk

N

∑ j j Var(t j , k)

⎤ N ⎡1 N ∑k k ⎢ N ∑ j j Var(t j , k)⎥ ⎦ ⎣ j

1/Ni∑iNi(xi⃗ (tj,k)

(6)

where Var(tj,k) = − the center of mass of solute at the kth time window xCM = 1/(NiNj)∑jNj∑iNixi⃗ (tj,k), k⃗ Ni is the total number of solute molecules, Nj is the number of snapshots in each time window, Nk is the total number of time windows, and xi⃗ (tj,k) is the position vector of the ith chain molecule of the jth snapshot in the kth time window.

qi

VCD(r ij) =

(5)

where γ = 0.5 is the damping coefficient representing the coupling strength to the heat bath and Ri is a time-dependent stochastic noise force per unit mass with zero mean and with intensity ⟨Ri(t) Rj(t + Δt)⟩ = 2γikBTBδ(Δt)δij.36,41 In eq 5, subscript i identifies each solute monomer. We use the following protocol to relax the system to thermal equilibrium at the target temperature T = 2.5 (300 K).33 The Stockmayer solvent is first relaxed for a short time of approximately 30 000 time steps with the solute fixed in its initial conformation. Next, the solute is permitted to diffuse and flex for between 105 and 2 × 105 time steps. The time step size begins at 0.0007τ for the initial relaxation and shifts up to 0.003τ for the remainder of the simulation where τ = 2.155 ps.33 Postprocessing of the data uses a variety of software, notably VMD, as noted above, and MATLAB for image processing and general data handling. Steady State. We introduce parameter ν, which is the variance normalized by the average variance within sequential time windows, to determine whether conformations are at steady state. This variance is based on the distance squared of all solute molecules from the position of the solute center of mass. A constant value of ν suggests no spreading, and ν increasing with time indicates that the solute is spreading away from the center of mass. Normalized variance ν in the kth time window is given by

⎛ σij ⎞6 ⎫ ⎪ − ⎜ ⎟ ⎬, ⎪ r ⎝ c⎠ ⎭ r ij > rc

(4)

(3) 15288

2 x⃗CM k ) ,

dx.doi.org/10.1021/la303252c | Langmuir 2012, 28, 15286−15293

Langmuir



Article

RESULTS Effect of Intermolecular Parameters. Because interaction strengths εij for dissimilar molecules i ≠ j are calculated using the Berthelot mixing rules, they are not all independent; only two independent interaction strengths are required to parametrize all cases. Figure 2 shows the long-time distributions of

Figure 3. Evolution over 1200τ (starting after an arbitrary time) of normalized variance ν (defined in eq 6) about the center of mass for the nonwetting drop (red triangle), partial wetting drop (blue circle), and pseudopartial wettinglike drop (green square) conformations. Whereas the other two conformations have reached steady state, the pseudopartial wetting drop continues to spread as shown by the three side views of its conformation and by its increasing variance. Error bars show the standard deviation of variance ν. Figure 2. Simulation snapshots (viewed along the z direction) of the steady-state distributions of solute molecules from an initial state as shown in Figure 1b at the values of interaction strengths indicated by the neighboring colored symbols and specified in Table S1. In all simulations shown, there are 400 solute molecules and |q| = 1. To reveal the distribution of solute molecules better, the Stockmayer atoms are not shown. The colored symbols defined for each set of parameter values are used in subsequent figures to help identify each case.

shows the transient evolution between these two states, with the pseudopartial wetting drop as a transient intermediate state. The thick disklike drops of the eight conformations across the top of Figure 2, εpp/εpl = 3, 4, are achieved for strong solute−solute interaction. For these conformations, we consider the angle between the dipole moment of the solute p⃗p and the y axis θ = cos−1(j ̂ ·pp⃗ )

(7)

where j ̂ is the unit vector along the y axis. The upper row in Figure 4 shows the orientational distribution of solute dipole moments as a function of the height above the substrate. Colors represents the local orientational distribution density from low (blue) to high (red). For εpp/εpl ≤ 2.23, the orientational distributions are symmetric with respect to θ = 90°. This symmetry breaks for εpp/εpl = 3, for which the orientation of

the solute at the various combinations of interaction strengths given in Table S1. Parameters εps/εpl and εpp/εpl were found to lead to the simplest clustering of the types of distributions. We follow the usual naming of the distributions as originally developed for macroscopic droplets, namely, nonwetting and partial wetting drops and perfectly wetting films. To these three types, we add the dispersed distribution in which the solute does not aggregate. For example, in Figure 2, nonwetting drops are seen near (εps/εpl, εpp/εpl) ≃ (0.25, 2), partial wetting drops are seen near (εps/εpl, εpp/εpl) ≃ {(0.5, 2), (1.5, 2.23)}, and nonwetting films are formed near (εps/εpl, εpp/εpl) ≃ (1,1.72). The wettability of the solute depends on the interactions among the solute, solvent, and substrate as described by the ratio εps/εpl. Proceeding from low to high values of εps/εpl leads to an increasing tendency to spread over the solute, as is most easily seen across the middle of the figure (i.e., εpp/εpl ≃ 2) from left to right. For large enough εps/εpl and small εpp/εpl (i.e., across the bottom of the figure) even as it disperses into the bulk, the solute maintains a monolayer on the substrate, as is most easily seen for εpp/εpl ≃ 1.4 from left to right. Once the solute−solute interaction εpp/εpl is large enough, the solute aggregates with no dispersed phase for εpp/εpl ≳ 1.7 and with thick disklike conformations for εpp/εpl ≳ 3. As εps/εpl increases, the aggregate tends to maintain contact with the substrate. As shown in Figure 3, variance ν, defined in eq 6, of both the nonwetting drop and the partial wetting drop converge in time to unity, which implies that both conformations have achieved steady-state conformations. However, ν increases in time for the transient pseudopartial wetting case. Because it is not a steady state, the pseudopartial wetting drop is not shown in Figure 2. It would be placed as a transient state for (εps/εpl, εpp/ εpl) ≃ (1.0, 2.0). That is, it would occupy a place between the partial wetting drop and the perfectly wetting film. Figure 3 also

Figure 4. (a) Solute orientation θ as a function of height y, where the angles are mapped onto the interval 0 ≤ θ ≤ 180°. The number density of solute molecules in a small slab of size 0.1 centered at y with orientation θ ± 0.5° is indicated by color from zero (blue) to three (red). The cases shown here, εps/εpl = 2.0 and |q| = 1.0, correspond to the following distributions: dispersed (left column), nonwetting drop (center column), and highly flattened nonwetting drop (right column). The orientational order within the solute changes with the type of distribution, as described in the text. (b) The blue line represents the number ratio N± of positively to negatively charged monomers. The green and red lines show the fraction of positive and negative charges, respectively, relative to the total number of charged monomers as a function of height. For εpp/εpl = 1.0 and 2.23, the total numbers are equal at each height; for εpp/εpl = 3.0, negative charges are displaced closer to the substrate relative to the positive charges. 15289

dx.doi.org/10.1021/la303252c | Langmuir 2012, 28, 15286−15293

Langmuir

Article

the dipole moment is predominantly in the 90° < θ < 180° quadrant. Consequently, conformations at εpp/εpl = 3 have a net dipole moment in the y direction. The lower row in Figure 4 shows the polarization of the three conformations. Recall that the solute head and tail monomers carry a positive and negative charge, respectively. For εpp/εpl = 1.0 and 2.23, there are equal numbers of positive and negative charges at each height above the substrate. For εpp/εpl ≥ 3, the net positive and negative charges are slightly displaced relative to one another. Overall, then, the solute aggregate carries a net dipole moment of approximately 60, which increases slightly as εps/εpl increases and is approximately an order of magnitude greater than the net dipole moment for εpp/εpl ≤ 2.23. Solute Concentration. The effect of solute concentration is shown in Figure 5. Across each row, as the number of solute

increases, the growth of the aggregate size on the substrate can be interpreted using concepts from percolation theory. The percolation threshold is defined as the coverage at which the network extends from one edge of the area to the opposite side. For a 2D lattice, the percolation threshold is at a fractional coverage of approximately 0.41, so the available area is 0.59.42,43 We find that a chain of solute molecules stretches from one side of the unit cell to the other between N = 160 and 240, for which the coverage is between 0.39 and 0.61, in agreement with the percolation threshold. Percolation theory also predicts that the largest contiguous area drops precipitously near the percolation threshold. In our work, the network can be considered to be formed by the solute molecules adhered onto the substrate, so 3 measures the available portion of the substrate. We quantify the largest contiguous area 3 on the substrate that remains unpopulated by any solute molecules as a function of the surface coverage Φ in Figure 6f. The fraction of unpopulated contiguous area decreases linearly for 0 < Φ < 0.4 and then resumes a slow decrease for 0.7 < Φ < 1.0. Our results show a large decrease in the largest contiguous area within the interval 0.4 < Φ < 0.7, as predicted from percolation theory.42 Row b in Figure 5 is for partial wetting conditions (εps/εpl, εpp/εpl) = (1.08, 1.58)). The shape of the aggregated solute and its connectivity to the substrate change with the number of solute molecules. For the smallest N, the solute forms a sphere detached from the substrate. For such a small N, the large surface area/volume ratio allows the solute−solvent interaction to hinder adsorption onto the substrate. For N ≥ 49, the solute is adsorbed onto the substrate. The shape variations of the adsorbed solute are similar in appearance to macroscopic droplets under the combined influence of surface tension and gravity. For macroscopic dipoles, the capillary length κ−1 ≡ (γ/ ρg)1/2 indicates the relative effects of surface tension γ and gravity, where ρ is density and g is gravitational acceleration. If the system size is much smaller than κ−1, then gravity is negligible. For molecular-sized systems, a κ−1-like characteristic length Δ arises from the relative contributions of surface tension and the attraction to the substrate. To formulate an expression for Δ, the role played by γ is taken on by εpp, and g is replaced by a function of εps,

Figure 5. Effect of the number of solute molecules, shown at the top of each column, on the distribution of solute for (a) perfect wetting, (b) partial wetting, and (c) nonwetting. All cases have |q| = 1.0.

molecules N is varied, all other parameter values are kept constant. The top row shows the results for perfectly wetting conditions, (εps/εpl,εpp/εpl) = (1.08,1.58)). For N = 25, 49, and 100, sufficient solute is not present to cover the substrate. Most of the molecules adsorb onto the substrate, forming networks of solute molecules (Figure 5a−c). The substrate is completely covered at N ≃ 400. Figure 6a−e shows the view looking down on the adsorbed solute for perfectly wetting conditions. As N

Eγ EvdW

=

2πδ 2γ Δ⎡

ps

∫1 ⎢⎣ A6H (δ

2

− h2) ⎤ ⎥⎦ h3

dh

(8)

where we have approximated the volume of N solute molecules composed of n monomers as nNσpp3 and where the Hamaker constant AHps is given by 4π2ρpρsεpsσps, with ρp being the number density of the solute and ρs being the number density of the substrate.44 For example, we find Δ = 2.6 when the surface tension γ is calculated using the surface energy per area, App H /24πσ0, where σ0 is the smallest distance between two solute molecules and the interaction parameters in App H are for the partial wetting drop. The length scale, defined by balancing forces due to van der Waals interaction from a semi-infinite substrate and the surface tension force on a circular drop of radius Δ, can be compared to the observed drop radius R ≃ (nN/4)1/3. It is expected that small enough drops Δ/R ≫ 1 will be nearly spherical because solute−solute forces dominate. In the limit of Δ/R ≪ 1, the solute−substrate interaction causes the drop to flatten.

Figure 6. Nonhomogeneous substrate coverage for a perfectly wetting film, shown as a function of the number of solute molecules for N = (a) 25, (b) 100, (c) 160, (d) 240, and (e) 400. (f) Fraction of the total substrate area of the largest contiguous unpopulated area 3 as a function of surface coverage Φ for N = 25, 49, 100, 120, 160, 240, 400, 784, and 1296. The slow linear decrease as Φ increases is interrupted by a precipitous drop in 3 near the percolation threshold of 0.41 for 2D lattices. Dots are data points, and the dashed lines are aids to the eye. 15290

dx.doi.org/10.1021/la303252c | Langmuir 2012, 28, 15286−15293

Langmuir

Article

molecules at these angles increases for the partial wetting drop case. These angles are the limiting geometries for the freelyrotating solute molecule. Figure 7b shows a few conformations at these limiting angles. Although increasing charge from |q| = 0.2 to 2.0 has no noticeable effect on the perfectly wetting film, small changes can be seen in the other three conformations. In general, the orientation is broadened for larger |q|. For example, for the partial wetting drop, the peaks at 60 and 120° are reduced in amplitude while the surrounding background of orientations is increased. Similarly, for larger |q|, the distribution of orientations for the nonwetting drop is less peaked and has a broader region on either side of the maximum. The partial wetting cases have hybrid characteristics combining the peak at 90° of the perfect wetting cases with a broad background of orientations as in the nonwetting cases. This combination of attributes can be understood by looking at the distribution of orientations as a function of height (Figure 8). For the partial

The bottom row of Figure 5 is for nonwetting conditions εpp/εpl = 2 and εps/εpl = 0.2. A substrate with small εps/εpl weakly attracts the solute, so solute molecules form nonwetting drop aggregates in the bulk solvent. The greater the number of solute molecules, the larger the size of the aggregates. For the greatest number of solute molecules, N = 1296, the large aggregation experiences interactions through the periodic boundary condition in the z direction. In contrast to small N aggregates in Figure 5c, the profile of the aggregate for N = 1296 is cylindrical, linking one side of the periodic cell to the other. This case, then, is contaminated by interactions through the boundary condition. Alternatively, we can imagine a substrate patterned with a periodic array of initial conditions, the unit cell of which is as shown in Figure 1b. Then, this case can be considered to represent the evolution of such a periodically repeated pattern. Charge |q|. For each of four pairs of Lennard-Jones parameter values, the charge +q and −q on the head and tail of the solute dipole was varied over the values of q = 0.2, 1.0, and 2.0. The physical shapes of the conformations are not visibly affected over the range of |q|. However, the organization of the solute molecules within the conformations is affected. Results are not shown for q = 1.0 because they are similar to the results for q = 2.0. We again use the angle defined in eq 7 between the dipole moment of the solute p⃗p and the y axis. In Figure 7a, the number of solute molecules within ±0.5° of angle θ is plotted versus θ. In going left to right across the figure, from a perfectly wetting film to a nonwetting drop, the distribution changes from sharply peaked near 90° (with nearly all solute molecules parallel to the substrate surface) to bell-shaped. There are two small peaks near 60 and 120°, and the proportion of solute

Figure 8. Color bar indicates the number of solute molecules at location y ± 0.5 with dipolar angle θ ± 0.5°. The perfectly wetting film is highly ordered, with nearly all solute molecules in the layer adjacent to the substrate having their axes parallel to the plane of the substrate, θ = 90°. The partial wetting drop also shows layering, where it appears that the layers are somewhat thicker than the thickness of the perfectly wetting film. The greater prevalence of θ = 60 and 120° suggests that the increase in thickness of the layers arises from the greater y extent of the solute molecules when they are no longer lying in the plane of the substrate. The layering and the concentration of solute orientation near the discrete angles of 60, 90, and 120° decrease with height, and near the top of the drop, the distribution of angles is uniform across nearly the entire range of 0° ≤ θ ≤ 180°. The angular distribution of the dipoles is highly localized near the substrate. The nonwetting drop shows the most disorder in its interior, where distinct layers appear, with little preference for any particular orientation.

wetting cases in the second column, the orientation is highly ordered near the substrate (2 ≲ y) with a peaked distribution similar to that of the perfectly wetting film shown in the first column; further from the substrate, the partial wetting cases have distributions similar to that of the nonwetting drop shown in the third column.

Figure 7. (a) Distribution of solute dipole orientation for dipole charges |q| = 0.2 (top row) and 2.0 (bottom row) for three typical configurations. (See the corresponding colored symbol in Figure 2). The fraction > is the ratio of the number of solute molecules within ±1/2° of the angles shown divided by the greatest number of solute molecules at any angle. (b) Fully folded solute conformations. Blue dots are the solute monomer, and red arrows are along the solute dipole moment. The solute with two or three of its monomers on the substrate is at orientations of 60 and 120° as observed under partially wetting conditions.



CONCLUSIONS This research has explored the spreading behavior and final equilibrium configurations of short, flexible dipolar polymers physiosorbed on a surface. Our simulations use a Stockmayer liquid solvent, which provides a water-like environment to a dipolar flexible-chain solute initially adhered on an FCC substrate. The parameters for a Stockmayer solvent are selected 15291

dx.doi.org/10.1021/la303252c | Langmuir 2012, 28, 15286−15293

Langmuir

Article

Notes

to represent water at room temperature. The redistribution of solute molecules is investigated by the systematic variation of the interaction strengths εpp/εpl and εps/εpl, solute number, dipolar charge |q|, and solute fraction. The normalized variance ν is introduced to determine the stability of solute formations. The pseudopartial wetting drop shows an increasing variance in time. We observe this unstable conformation spreading to the perfectly wetting film at steady state.45−50 We observe intervals, during which the rate of spreading slows to nearly zero, as shown in Figure 3 for 258 < time < 774. During periods of reduced spreading rate, the central portion of the droplet reorganizes its shape by coarsening, while the precursor film remains at a constant area. Once the smooth drop extends over the entire precursor film, spreading, as measured by ν, resumes. The distribution of solute on the substrate under film conditions can be interpreted using 2D percolation theory. The largest contiguous area normalized by the total surface area drops beyond a critical surface coverage given by the percolation threshold. For the partial wetting drop case, we introduce a length scale Δ to interpret the plateau region generated when RN/Δ ≳ 1: the plateau is due to a large enough substrate interaction compared to surface tension forces. The effect of dipolar charge |q| does not appear in the visual appearance of the conformation but in the distribution of solute dipole-moment orientations. Peaks in the orientation are dominated by the tendency for the solute to be aligned with the substrate and to fold into fully folded configurations. For large enough solute−solute interaction, solute aggregates as a whole possess a dipole moment.51 As shown in the last column of Figure 6, the solute dipole moment distributes asymmetrically and the strength of the net dipole moment is at least an order of magnitude higher than for other conformations. Gregory et al. showed how water clusters enhance the strength of their net dipole moment by relative arrangements of water orientations and their asymmetry.52 In general, the large change in the dipole moment reflects a change in the solvent structure.51 This structural change is important in understanding phase transitions.52 Furthermore, near proteins, the structure of vicinal water is intrinsic to protein function, affecting folding, molecular recognition, enzymatic action, and protein−protein interactions.53 The relationship between the dipole moment and structure remains an area of active research.54 From a practical point of view, this work can contribute to the design of periodic patterns on a substrate. Controlling the solute fraction can yield isolated cells of solute versus solute connected from one cell to its periodically repeated neighbors. For example, in nanolens fabrication on a substrate, results such as those given here can be used to find the fraction of solute molecules needed for a given pattern separation and lens curvature.55



The authors declare no competing financial interest.



ACKNOWLEDGMENTS Simulations were performed on QUEST, Northwestern University’s high-performance computing facility operated by Northwestern University Information Technology. Q.W. and T.Y. gratefully acknowledge support from the National Science Foundation and the Army Tank-Automotive and Armaments Command (TACOM). S.L. thanks the Lillian Sidney Foundation for support. We thank colleagues Igal Szleifer and Neelesh Patankar for many helpful discussions.



(1) George, C. B.; Ratner, M. A.; Szleifer, I. Langmuir monolayers with internal dipoles: understanding phase behavior using Monte Carlo simulations. J. Chem. Phys. 2010, 132, 014703. (2) Coninck, J. D.; Blake, T. D. Wetting and molecular dynamics simulations of simple liquids. Ann. Rev. Mater. Res. 2008, 38, 1−22. (3) de Gennes., P. G. Wetting: statics and dynamics. Rev. Mod. Phys. 1985, 57, 827−863. (4) Bertrand, E.; Dobbs, H.; Broseta, D.; Indekeu, J.; Bonn, D.; Meunier, J. First-order and critical wetting of alkanes on water. Phys. Rev. Lett. 2000, 85, 1282−1285. (5) Mann, S. Self-assembly and transformation of hybrid nanoobjects and nanostructures under equilibrium and non-equilibrium conditions. Nat. Mater. 2009, 8, 781−792. (6) Grzybowski, B. A.; Wilmer, C. E.; Kim, J.; Browne, K. P.; Bishop, K. J. M. Self-assembly: from crystals to cells. Soft Matter 2009, 5, 1110−1128. (7) Whitesides, G. M.; Laibinis, P. E. Wet chemical approaches to the characterization of organic surfaces: self-assembled monolayers, wetting, and the physical-organic chemistry of the solid-liquid interface. Langmuir 1990, 6, 87−96. (8) Ma, S. Modern Theory of Critical Phenomena; Westview Press: Reading, MA, 2000. (9) Young, T. An essay on the cohesion of fluids. Philos. Trans. R. Soc. London 1805, 95, 65−87. (10) Brochard-Wyart, F.; Meglio, J. M. D.; Quéré, D.; de Gennes, P. G. Spreading of nonvolatile liquids in a continuum picture. Langmuir 1991, 7, 335−338. (11) de Gennes, P.-G.; Brochard-Wyart, F.; Quéré, D. Capillarity and Wetting Phenomena: Drops, Bubbles, Pearls, Waves; Springer: New York, 2004. (12) Rowlinson, J.; Widom, B. Molecular Theory of Capillarity; Clarendon Press: Oxford, Oxfordshire, 1982. (13) Yeh, I.; Berkowitz, M. L. Ewald summation for systems with slab geometry. J. Chem. Phys. 1999, 111, 3155−3162. (14) Zhang, X.; Neogi, P.; Ybarra, R. M. Stable drop shapes under disjoining pressure: I. A hierarchical approach and application. J. Colloid Interface Sci. 2002, 249, 134−140. (15) Diaz, M. E.; Fuentes, J.; Cerro, R. L.; Savage, M. D. An analytical solution for a partially wetting puddle and the location of the static contact angle. J. Colloid Interface Sci. 2010, 348, 232−239. (16) Wu, Q.; Wong, H. A slope-dependent disjoining pressure for non-zero contact angles. J. Fluid Mech. 2004, 506, 157−185. (17) Yi, T.; Wong, H. Theory of slope-dependent disjoining pressure with application to Lennard−Jones liquid films. J. Colloid Interface Sci. 2007, 313, 579−591. (18) Allen, M. P.; Tildesley, D. J. Computer Simulation of Liquids; Oxford University Press: New York, 1989. (19) Lundgren, M.; Allan, N. L.; Cosgrove, T.; George, N. Wetting of water and water/ethanol droplets on a non-polar surface: a molecular dynamics study. Langmuir 2002, 18, 10462−10466. (20) Werder, T.; Walther, J. H.; Jaffe, R. L.; Halicioglu, T.; Koumoutsakos, P. On the water−carbon interaction for use in

ASSOCIATED CONTENT

S Supporting Information *

Parameters for all simulations and results for orientational entropy. This material is available free of charge via the Internet at http://pubs.acs.org.



REFERENCES

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]; s-lichter@ northwestern.edu. 15292

dx.doi.org/10.1021/la303252c | Langmuir 2012, 28, 15286−15293

Langmuir

Article

molecular dynamics simulations of graphite and carbon nanotubes. J. Phys. Chem. B 2003, 107, 1345−1352. (21) Hirvi, J. Y.; Pakkanen, T. A. Molecular dynamics simulations water droplets on polymer surface. J. Chem. Phys. 2006, 125, 144712. (22) Mognetti, B. M.; Kusumaatmaja, H.; Yeomans, J. M. Drop dynamics on hydrophobic and superhydrophobic surfaces. Faraday Discuss. 2010, 146, 153−165. (23) Xu, H.; Shirvanyants, D.; Beers, K. L.; Matyjaszewski, K.; Dobrynin, A. V.; Rubinstein, M.; Sheiko, S. S. Molecular Visualization of conformation-triggered flow instability. Phys. Rev. Lett. 2005, 94, 237801. (24) Gannon, G.; Greer, J. C.; Larsson, J. A.; Thompson, D. Molecular dynamics study of naturally occurring defects in selfassembly monolayer formation. ACS Nano 2010, 4, 921−932. (25) Ahn, Y.; Saha, J. K.; Schatz, G. C.; Jang, J. Molecular dynamics study of the formation of a self-assembled monolayer on gold. J. Phys. Chem. C 2011, 115, 10668−10674. (26) Wei, T.; Carignano, M. A.; Szleifer, I. Lysozyme adsorption on polyethylene surfaces: why are long simulations needed? Langmuir 2011, 27, 12074−12081. (27) Xue, B.; Wang, J.; Wang, W. Collapse of homopolymer chains with two fixed terminals. J. Chem. Phys. 2003, 119, 7534−7542. (28) Saletra, W. J.; Xu, H.; Vosch, T.; Feyter, S. D.; Amabilino, D. B. Chiral non-periodic surface-confined molecular nanopatterns revealed by scanning tunneling microscopy. CrystEngComm 2011, 13, 5578− 5583. (29) Pollock, E. L.; Alder, B. J. Static dielectric properties of Stockmayer fluids. Physica A 1980, 102, 1−21. (30) Adams, D. J.; Adams, E. M. Static dielectric properties of the Stockmayer fluid from computer simulation. Mol. Phys. 1981, 42, 902− 926. (31) Leeuw, S. W. D.; Perram, J. W.; Smith, E. R. Computer simulation of the static dielectric constant of systems with permanent electric dipoles. Annu. Rev. Phys. Chem. 1986, 37, 245−270. (32) Stevens, M. J.; Grest, G. S. Phase coexistence for the Stockmayer fluid in an applied field. Phys. Rev. E 1995, 51, 5976. (33) Payne, V. A.; Forsyth, M.; Kolafa, J.; Ratner, M. A.; de Leeuw, S. W. Dipole time correlation functions of the Stockmayer fluid in the microcanonical and canonical ensembles. J. Phys. Chem. 1993, 97, 10478−10485. (34) Kremer, K.; Grest, G. S. Dynamics of entangled linear polymer melts: a molecular-dynamics simulation. J. Chem. Phys. 1990, 92, 5057−5086. (35) Priezjev, N. V.; Troian, S. M. Influence of periodic wall roughness on the slip behaviour at liquid/solid interfaces: molecularscale simulations versus continuum predictions. J. Fluid Mech. 2006, 554, 25−46. (36) Priezjev, N. V. Relationship between induced fluid structure and boundary slip in nanoscale polymer films. Phys. Rev. E 2010, 82, 051603. (37) Thompson, P. A.; Troian, S. M. A general boundary condition for liquid flow at solid surfaces. Nature 1997, 389, 360−362. (38) Plimpton, S. Fast parallel algorithms for short-range molecular dynamics. J. Comput. Phys. 1995, 117, 1−19. (39) Humphrey, W.; Dalke, A.; Schulten, K. VMD - visual molecular dynamics. J. Mol. Graphics 1996, 14, 33−38. (40) Berendsen, H. J. C.; Postma, J. P. M.; van Gunsteren, W. F.; DiNola, A.; Haak, J. R. Molecular dynamics with coupling to an external bath. J. Chem. Phys. 1984, 81, 3684−3690. (41) Priezjev, N. V. Shear rate threshold for the boundary slip in dense polymer films. Phys. Rev. E 2009, 80, 031608. (42) Jones, R. A. L. Soft Condensed Matter; Oxford University Press: New York, 2003. (43) Yong, W.; Schmittmann, B.; Zia, R. K. P. Two-dimensional polymer networks near percolation. J. Phys. A: Math. Theor. 2008, 41, 025004. (44) Butt, H. J.; Graf, K.; Kappl, M. Physics and Chemistry of Interfaces; Wiley-VCH: Weinheim, Germany, 2006.

(45) Silberzan, P.; Leger, L. Evidence for a new spreading regime between partial and total wetting. Phys. Rev. Lett. 1991, 66, 185−188. (46) Xu, H.; Shirvanyants, D.; Beers, K.; Matyjaszewski, K.; Rubinstein, M.; Sheiko, S. S. Molecular motion in a spreading precursor film. Phys. Rev. Lett. 2004, 93, 206103. (47) Sedighi, N.; Murad, S.; Aggarwal, S. K. Molecular dynamics simulations of spontaneous spreading of a nanodroplet on solid surfaces. Fluid Dyn. Res. 2011, 43, 015507. (48) Wu, C.; Qian, T.; Sheng, P. Droplet spreading driven by Van der Waals force: a molecular dynamics study. J. Phys.: Condens. Matter 2010, 22, 325101. (49) Douezan, S.; Dumond, J.; Brochard-Wyart, F. Wetting transitions of cellular aggregates induced by substrate rigidity. Soft Matter 2012, 8, 4578−4583. (50) Popescu, M. N.; Oshanin, G.; Dietrich, S.; Cazabat, A.-M. Precursor films in wetting phenomena. J. Phys.: Condens. Matter 2012, 24, 243102. (51) Gelman-Constantin, J.; Carignano, M. A.; Szleifer, I.; Marceca, E. J.; Corti, H. R. Structural transitions and dipole moment of water clusters (H2O)(n=4−100). J. Chem. Phys. 2010, 133, 024506. (52) Gregory, J. K.; Clary, D. C.; Liu, K.; Brown, M. G.; Saykally, R. J. The water dipole moment in water clusters. Science 1997, 275, 814− 817. (53) Pal, S. K.; Peon, J.; Zewail, A. H. Biological water at the protein surface: dynamical solvation probed directly with femtosecond resolution. Proc. Natl. Acad. Sci. U.S.A. 2002, 99, 1763−1768. (54) Burnham, C. J.; Li, J.; Xantheas, S. S.; Leslie, M. The parametrization of a Thole-type all-atom polarizable water model from first principles and its application to the study of water clusters (n=2− 21) and the phonon spectrum of ice Ih. J. Chem. Phys. 1999, 110, 4566−4581. (55) Verma, A.; Sharma, A. Self-organized nano-lens arrays by intensified dewetting of electron beam modified polymer thin-films. Soft Matter 2011, 7, 11119−11124.

15293

dx.doi.org/10.1021/la303252c | Langmuir 2012, 28, 15286−15293