J. Phys. Chem. B 2001, 105, 9503-9508
9503
Molecular Dynamics Study of Partial Monolayer Ordering of Chain Molecules Raymond D. Mountain* Physical and Chemical Properties DiVision, 100 Bureau DriVe Stop 8380, National Institute of Standards and Technology, Gaithersburg, Maryland 20899-8380
Joseph B. Hubbard, Curtis W. Meuse, and Vernon Simmons Biotechnology DiVision, 100 Bureau DriVe Stop 8313, National Institute of Standards and Technology, Gaithersburg, Maryland 20899-8313 ReceiVed: June 11, 2001; In Final Form: July 20, 2001
We present the results of a molecular dynamics study of a partial monolayer of self-assembled octadecanethiol molecules. The correlations between various statistical measures of surface induced chain ordering are examined. These include the density profile, chain morphology, tilt angle distributions, and gauche defect distributions. Particular attention is focused on the significance of the strength of the alkane chain surface interaction, as well as the role of temperature, on the type and degree of disorder which we observe. Our simulations are in accord with experimental evidence which indicates that the quality, in the sense of reproducibility in the laboratory, of dense films increases as the coupling of the chain molecules with the surface decreases. We believe that our findings have general implications for the establishment of experimental protocols for selfassembled surface films of organic molecules with varying degrees of prescribed disorder.
1. Introduction Molecular dynamics simulations have an important role in the investigation of condensed systems. One significant feature of simulation studies is the ability to examine features that are not readily accessible to experimental methods. In this paper we use that ability of molecular dynamics to examine structural features of self-assembled monolayers of long chain molecules on a surface. These molecules could be alkanethiols or alkanesilanes. In particular, we investigate the structure of partially filled layers as the attraction of the surface for the methylene and methyl groups is varied. There is experimental evidence to indicate that the quality, in the sense of reproducibility in the laboratory, of dense films increases as the coupling of the chain molecules with the surface is decreased.1 Our simulations support that experimental observation and show that there is a range in the strength of the interaction between the methyl and methylene groups and the surface where the character of the structure of the partially filled films changes from islands of nearly upright chains coexisting with empty regions of the surface (the weak coupling case) to less structured regions of partially ordered chains coexisting with highly disordered chains that completely cover the surface (the strong coupling case). The simulations also show that as the system temperature is increased above ambient, the structure in the intermediate coupling case (defined below) becomes less ordered. This also is consistent with experimental results.2-4 The model employed in the simulations is described in section 2. Section 3 contains a description of the simulations and also describes the measures of order used to characterize the structure of the monolayers. Three types of order are considered. The first is the collective arrangement of the molecules as indicated by the density profile of the chains relative to the tethering * Corresponding author.
10.1021/jp012193y
surface. The second feature that measures order is the distribution of tilt angles of the chains. The final measure of order is the fraction of gauche defects found along the chains. We shall see that the defect fraction, a quantity that is measurable in an average sense,5-7 correlates closely with the density profile and tilt angle distribution. These quantitative measures of order are complemented by images of the chains. The images provide a qualitative, but informative picture of how the chains are arrayed. Section 4 is devoted to the results obtained at ambient conditions. Section 5 reports the results of the temperature variation study. Section 6 contains a summary of what has been learned here. 2. The Model A united atom representation of a 19-site chain molecule (octadecanethiol) is used in these simulations. Other long chain molecules such as octadecyltriclorosilane could be considered with appropriate changes in the headgroup interactions. Figure 1 provides an illustration of the geometrical arrangement of a chain and some notation that will be used below. The intermolecular and intramolecular interaction potentials are those of model I of Hautman and Klein8 with the addition of harmonic stretch interactions between adjacent intramolecular sites. The details of the model are included here. The coefficients for the potential parameters are listed in Tables 1-4. Unless otherwise noted, the values of the parameters are those of Model I for Hautman and Klein. This model has been shown to result in dense-packed films (number density of 4.65 molecules per nm2) that are consistent with the structure observed experimentally.8 The coefficients in these potentials of course depend on the sites involved. The extra indices needed to indicate the site dependence have been supressed to simplify the notation. The
This article not subject to U.S. Copyright. Published 2001 by the American Chemical Society Published on Web 09/01/2001
9504 J. Phys. Chem. B, Vol. 105, No. 39, 2001
Mountain et al. TABLE 3: Lennard-Jones Parameters: Unlike Site Interaction Parameters Determined Using the Geometric Mean of the Like Site Parametersa sites
σ, Å
, K
SI-SI CH3-CH3 CH2-CH2 CH3-CH2 CH3-S CH2-S
4.25 3.905 3.905 3.905 3.723 3.723
200.0 88.1 59.4 72.3 105.4 86.5
a For the purpose of determining the S-CH3 and the S-CH2 parameters, the S is 126 K and the σS is 3.55 Å rather than the values listed in this table. For example, σ for the CH3-S interaction is [3.55 × 3.905]1/2 ) 3.723 Å
TABLE 4: Parameters for the Surface Interactions8
Figure 1. Here is a single chain molecule containing 19 sites in an all trans configuration. The molecule makes an angle of ψ ) 30° with the z-axis which is normal to the tethering surface. The thiol headgroup site is indicated by a large filled circle and is referred to as site number 1. Sites 2 through 18 are CH2 united atom sites indicated by small shaded circles. Site 19 is the terminal methyl united atom site and is indicated as the small filled circle. The bonds between adjacent sites where torsion motions are possible are indicated by bold lines.
TABLE 1: Stretch Parameters kS and d0a parameter 107
kS, d0, Å a
K/nm2
S-C
C-C
4.529 1.82
2.265 1.54
The CH3 and CH2 sites are equivalent for the stretch interaction.13
103
kθ, θ0, deg
K/rad2
C-C-C
S-C-C
62.5 109.5
62.5 114.4
intramolecular interaction consists of stretch, bend, torsion, and Lennard-Jones terms. Explicitly,
Vintramolecular ) V2 + V3 + V4 + VL-J The bond stretch term is
V2 )
kS N-1 2
(|ri+1 - ri| - d0)2 ∑ i)1
where d0 is the equilibrium bond length. The sum runs over the N - 1 bonds and ri is the position of the ith site in the chain. The stretch parameters kS and d0 are listed in Table 1. The bend term is
V3 )
kθ N-2 2
C12, 107 KÅ12
C 3, K Å 3
z 0, Å
CH3 CH2 S
3.41 2.80 4.089
20800 17100 180600
0.860 0.860 0.269
N-3 5
V4 )
(cosθi+1 - cosθ0)2 ∑ i)1
where θi+1 is the angle with its vertex at site i + 1 formed by bonds connecting sites i, i + 1 and sites i + 1, i + 2, and q0 is the equilibrium angle at site i + 1. The bend parameters are listed in Table 2. The torsion term involves adjacent quadruples of sites with the form9
al cosl(φ) ∑ ∑ i)1 i)0
where φ is the dihedral angle between the two planes formed by the four adjacent sites. The vertexes of the dihedral angles are indicated as bold lines in Figure 1. There is no site dependence for these coefficients with values (in units of K) of a0 ) 1116, a1 ) 1462, a2 ) -1578, a3 ) -368, a4 ) 3156, and a5 ) -37788. The final intramolecular interaction term is a Lennard-Jones interaction between sites separated by three or more sites.
VL-J )
TABLE 2: Bend Parameters kθ and θ0 parameter
site
4[(σ/r)12 - (σ/r)6] ∑ sites
The intermolecular interaction is the same Lennard-Jones interaction used in VL-J and acts between sites on different molecules. The Lennard-Jones potential parameters are listed in Table 3. There is also a surface interaction of the form
Vsurface )
C12 (z - z0)
12
RC3 (z - z0)3
The interaction of the chain sites with the surface depends only on the distance of the site above the surface. The factor R multiplying the coefficient C3 is used to vary the attraction of the surface for the methyl and methylene groups. For the thiol group, R is always unity. The parameters for the surface potential are listed in Table 4. The values of C3 in Table 4 for the CH3 and CH2 sites are for the strong coupling case mentioned in the Introduction, that is for R ) 1. The weak coupling case has R values for the methyl and methylene sites of 1/10. We also consider scaling factors of 1/ (intermediate coupling) and of 3/ and 7/ to examine the 2 4 8 change from weak and intermediate coupling to strong coupling to the surface. 3. The Simulations Our simulations involve 225 chains tethered to a smooth surface by the surface potential. The chains are located in an L × L square planar region subject to periodic boundary conditions in the planar directions. The z-axis is normal to the smooth surface that is located in the z ) 0 plane. The size of
Monolayer Ordering of Chain Molecules the square L is 9.48 nm so that the average density is 2.50 molecules/nm2, about 1/2 the close packed density. The simulations were started by placing the chains in a vertical orientation (ψ ) 0) at random positions with no overlaps. The thiol sites were placed at the stable position of 0.24 nm above the z ) 0 surface. Then the system was allowed to evolve, subject to a Nose´-Hoover thermostat10 that maintains the temperature of the system at 297 K (ambient conditions). The equations of motion were integrated using the Beeman algorithm with a time step of 1 fs.11,12 The system was allowed to evolve until “stabilized”. Here stability means that the various order parameters described below were unchanged over sequential, 100 ps duration runs. We do not imply that thermal equilibrium was obtained in each case discussed. In fact, different initial configurations lead to significant differences in the measures of order when R ) 1. In addition to the simulations at 297 K, a second set of simulations were performed to investigate the temperature variation of the amount of order in the system for intermediate coupling, R ) 1/2. The details of this set of simulations are described in the section 5. Three measures of order are considered. The first is the density profile. During the simulation, the number of sites located at a distance between z and z + dz above the surface is monitored. The resolution in z, dz, is 0.00316 nm. The second measure of order is the distibution of tilt angles of the chains. As indicated in Figure 1, the tilt angle ψ is taken to be the angle between the normal to the surface and the vector connecting the headgroup to the methyl group. The cosines of the tilt angles are binned with a resolution of 0.01. Finally, the mean number of gauche defects fG for each of the bonds is determined. A gauche defect is said to be present at a bond between two sites in a chain if the torsion angle about the bond is greater in magnitude than 66°, the position of the maximum in the torsion potential.9 4. Results for Ambient Conditions In this section we describe the results for the three measures of order, the density profiles, the distribution of tilt angles ψ, and the mean number of gauche defects per site fG. The variation of the density profile n(z) of the chains with the change in the strength of the surface attraction is indicated in Figure 2. The distance of the sites above the surface is z. Figures 2a-d are for values of the coupling parameter R of 1/ , 1/ , 7/ , and 1, respectively. 10 2 8 There are several features of these profiles to be noted. When the surface coupling is weak, R ) 1/10, the profile is quite structured. There are 19 distinct peaks in Figure 2a suggesting that the chains are in nearly identical configurations. This structure becomes less pronounced as the coupling increases. The peaks for larger values of z, the distance above the surface, become blurred into less pronounced features. Also, the maximum extension of the chains decreases. For R ) 1/2 (and also for 3/4 which is not shown as the profile is essentially that for 1/2), the main effect is the blurring. The profile remains above n(z) ) 1/2 for z < 2 nm. As R increases from 3/4 to 7/8, the profile structure decreases significantly and n(z) > 1/2 occurs only for z < 1.6 nm. The further increase of R to 1, the strong coupling case, produces an almost linear decrease in density for z > 1.2 nm. The maximum height of a chain remains at about 2.2 nm as the other structural features become fainter. Next we examine the distribution of tilt angles ψ. The distributions of cos(ψ) are shown in Figure 3 for R’s of 1/10 (solid line), 1/2 (short dashed line), 7/8 (long dashed line), and 1
J. Phys. Chem. B, Vol. 105, No. 39, 2001 9505
Figure 2. The density profiles are shown for (a) R ) 1/10, (b) R ) 1/2, (c) R ) 7/8, and (d) R ) 1. The profile for R ) 3/4 is basically the same as for R ) 1/2 and is not displayed.
Figure 3. The tilt angle distributions, P[cosψ], are shown for R ) 1/10 (solid line), R ) 1/2 (dashed line), R ) 7/8 (long dashed line), and R ) 1 (long short dashed line). The distribution for R ) 3/4 is not very different from the distribution for R ) 1/2 and is not displayed.
(long-short dashed line). The distributions for 1/2 and 3/4 are not very different, so the distribution for R ) 3/4 is not shown. The chains in the weak coupling case are tilted about 20° to the normal with a narrow distribution of angles. As the coupling increases to R ) 1/2, the average tilt increases to about 35° with a broader distribution about the maximum in P[cos(ψ)]. Further increases in the coupling to the surface lead to a significant fraction of the chains with tilt angles greater than 50°. As is the case for the density profiles, the major change in the tilt angle distributions occurs when R > 3/4. The average fractions of gauche defects fG for each site on the chains are shown in Figure 4. The term “site” refers to the bond between two united atoms. For example, the defects on site 2.5 are associated with the bond between united atoms 2 and 3 as shown in Figure 1. The weak coupling case has few defects although the fraction increases for the last two bonds. As the coupling increases in strength, the fraction of defects increases, particularly close to the surface and close to the free end of the chain. As in the other two order measures, the difference between R ) 1/2 and R ) 3/4 is small. When R increases from 3/4 to 7/8, the major increase in defects is in the bottom half of the chain. The further increase of R to 1 results in an overall increase in defects over the entire length of the chain.
9506 J. Phys. Chem. B, Vol. 105, No. 39, 2001
Mountain et al.
Figure 4. The mean number of gauche defects for the chains are indicated here. The filled circles are for the weak coupling case, R ) 1/ , the filled squares are for R ) 1/ , the open diamonds are for 3/ , 10 2 4 the open trangles are for 7/8, and the filled triangles are for the strong coupling case, R ) 1. Figure 6. A snapshot view of the chains with R ) 1/2 shows that the disordered chains are located mostly along the interface between the island and the cavity. The tilt of the chains is more evident than in Figure 5.
Figure 5. A snapshot view of the chains with R ) 1/10 shows that most of the chains are in a nearly upright position with local spatial order as well as indicated by the rows of chain molecules. The periodic boundary conditions mean that the island of chains and the cavity have dimensions on the order of, or greater than, the simulation box size of 9.48 nm.
In addition to the quantitative measures of order that have been discussed, a more qualitative measure of order is available in the form of snapshots of the chain configurations. Snapshots of the final configurations for the R ) 1/10, 1/2, 7/8, and 1 cases are displayed in Figures 5-8, respectively. A view of the R ) 3/ case is not included as it is not particularly different from 4 the view of the 1/2 case. The view is looking down on the film along the z-axis. In these figures, the methyl groups are the large, yellow spheres, the methylene groups are smaller, blue spheres, and the thiol groups are small red spheres. The thiol groups are mostly not visible as they are obscurred by the rest of the chain. The surface to which the chains are tethered is white. These views of the chain configurations provide an indication of the morphology of the films. An indication of spatial order in these figures is the presence of well-defined rows of chain molecules. These are prominent in Figures 5 and 6, and less so in Figures 7 and 8. The tilt angle ψ is small when much of the methylene part of the chain is obscurred by neighboring chains. This is most pronounced
Figure 7. A snapshot view of the chains with R ) 7/8 shows that the chains have nearly filled in the cavity.
in Figure 5. As the tilt angle increases, more of an individual chain becomes visible. This is certainly the case in Figures 7 and 8 where according to Figure 3, there is a broad distribution of tilt angles. A well-defined cavity in the snapshots occurs when ordered rows of chains are present in the stable system. 5. Temperature Variation of Order A separate set of simulations were made to examine the temperature variation of the order in the system. For this set, the size of the simulation cell was reduced slightly to L ) 8.85 nm as part of a more extensive survey of film structures. To be certain that the ordered state evolves from the disordered one rather than being a consequence of initial conditions, we started with R ) 1 at 450 K so that a disordered, filled in configuration of chains was realized. Next, the temperature was reduced to 297 K and the system was run for 320 ps. There was no observable change in the disorder of the system in terms of the
Monolayer Ordering of Chain Molecules
Figure 8. A snapshot view of the chains with R ) 1 shows the ordered regions are smaller than the size of the simulation cell.
Figure 9. A snapshot view of the chains with R ) 1, just before the coupling with the surface was changed.
density profile, the tilt angle distribution, and the gauche defects, and as indicated in Figure 9. The format of the snapshots displayed in Figures 9-12 is described in the previous section. Snapshots are not a substitute for the quantative order measures. We monitored the order measures described in section 4 during the sequence of runs described below. The snapshots and our description of the state of the system are consistent with the order measures. The snapshots are used only to illustrate the trends in the order of the system. The surface coupling was reduced by setting R ) 1/2, and the system then evolved from the configuration shown in Figure 9 for 180 ps. At 100 ps the system was disordered, but at 180 ps, the system had developed a cavity as a portion of the surface was uncovered. This cavity evolved into a well-defined object surrounded by nearly upright chains over the next 170 ps. The configuration at the end of this part of the simulation is shown in Figure 10. It is evident from Figure 10 that the formation of the cavity precedes the development of spatial order in the form of ordered rows of molecules. The development of ordered rows
J. Phys. Chem. B, Vol. 105, No. 39, 2001 9507
Figure 10. A snapshot view of the chains with R ) 1/2, after 350 ps. A well-defined cavity has formed, but the spatial ordering of the chains is not well-developed.
Figure 11. A snapshot view of the chains for T ) 329 K, after 300 ps. Note the ordering of the chains into well-defined rows.
required and additional 200 ps. This evolution process gave us confidence that it was not necessary to achieve a fully stabilized system in order to infer the trends in the order of the system as the temperature was increased. As noted previously, a cavity is associated with ordered rows of chains in the stabilized system. The presence or absence of a cavity is an indicator of how the system will evolve. The temperature of the system shown in Figure 10 was then increased to 329 K and run for 300 ps. The cavity remained well-defined and the spatial ordering of the chains increased with the formation of well-defined rows of chains. This is shown in Figure 11. For this slightly higher density system (2.87 molcules/nm2), the cavity appears to be embedded in the chains rather than coexisting with an island. Next, the temperature was increased to 360 K and run for 220 ps. The cavity decreased slightly in size as the “edge” of the cavity became more irregular. Then the temperature was increased to 391 K and run for 320 ps. The size of the cavity decreased during this interval, as is shown in Figure 12. Also, the spatial extent of well-defined rows has decreased.
9508 J. Phys. Chem. B, Vol. 105, No. 39, 2001
Figure 12. A snapshot view of the chains for T ) 360 K after 220 ps. The size of the cavity has decreased, the edge of the cavity is less sharp, and the spatial extent of the ordered rows is reduced from Figure 11.
Finally the temperature was increased to 422 K. During the next 100 ps, the cavity became smaller and considerable disorder was evident. After another 85 ps, the cavity disappeared as the surface was completely covered by the chains, quite similar in appearance to Figure 9. While there are some moderately long simulations in this set, no attempt was made to demonstrate stability of the states as was done in the previous set of simulations. 6. Discussion The types of order that develop in tethered chain molecules when the coverage is on the order of 1/2 of the close packed coverage have been examined in terms of the strength of the attraction between the tethering surface and the methyl and methylene groups of the chains. In the weak coupling case, R ) 1/10, the chains form ordered islands of nearly upright chains coexisting with regions of the surface that contain no chains (a cavity). The distribution of tilt angles is narrowly peaked about 20° and the density profile is highly structured. The snapshot of this case suggests that the island consists of ordered rows of nearly upright chains arranged in an approximately hexagonal arrangement. Consistent with this view is relatively low fraction of gauche defects over most of the chain sites. The fraction increases to about 1/10 near the methyl end group. On the other hand, for the strong coupling case, R ) 1, the chains form small islands or ordered regions coexisting with disordered chains that completely cover the tethering surface. The distribution of tilt angles is broad and there are a significant fraction of gauche defects over all sites in the chains. The density profile is correspondingly less structured and decreases in magnitude as the height above the surface exceeds 1.2 nm.
Mountain et al. The case of intermediate coupling, 1/2 e R e 3/4, moderates the order of the weak coupling case by increasing the mean tilt angle to about 35°. The island of chains remains in coexistence with the cavity. In the narrow range of coupling, 3/4 < R e 7/8, the arrangement of the chains changes. The rows of ordered chains are not as extensive and the cavity begins to fill in with chains that are strongly tilted. It is not possible to say if this change in ordering is a phase transition on the basis of the results obtained here. The temperature variation study does indicate that for intermediate coupling, the ordered coexistence of upright chains and a cavity is the stable arrangement, for temperatures from ambient to 360 K. At higher temperatures, the system reverts to a disordered arrangement with no cavity and no clear ordering of rows of chains. Our simulations of submonolayer alkanethiol films are consistent with experimental findings. In particular, we find that the degree of order in the films correlates with the strength of the surface interaction and that increasing temperature promotes chain disorder. In addition, we find that several detailed statistical measures of surface film order that are not directly experimentally observable correlate nicely with more readily accessible averaged properties. These include density profiles, domain morphology, tilt angle distributions, and gauche defect distributions. Our results help understand the experimental finding that weak surface coupling promotes reproducibility of laboratory films.1 With weak coupling, the tendency of the chains to form ordered islands leaves room for the addition of chains with minimal rearrangement of already present molecules. With strong coupling, it becomes necessary to move one or more already attached molecules so that the thiol group of an added molecule can reach the equilibrium point 0.24 nm above the surface. Clearly, the growth of well-ordered films will be easier when the surface coupling is weak. Insofar as our simulations employ generic interaction potentials, we believe that our findings have general implications for the establishment of experimental protocols for self-assembled surface films of organic molecules with varying degrees of prescribed disorder. References and Notes (1) Allara, D. L.; Parikh, A. N.; Rondelez, F. Langmuir 1995, 11, 2357. (2) Parikh, A. N.; Allara, D. L.; Azouz, I. B.; Rondelez, F. J. Phys. Chem. 1994, 98, 7577. (3) Brzoska, J. B.; Azouz, I. B.; Rondelez, F. Langmuir 1994, 10, 4367. (4) Sung, M. M.; Carraro, C.; Yauw, O. W.; Kim, Y.; Maboudian, R. J. Phys. Chem. B 2000, 104, 1556. (5) Parikh, A. N.; Allara, D. L. J. Chem. Phys. 1992, 96, 927. (6) Terrill, R. H.; Tanzer, T. A.; Bohn, P. W. Langmuir 1998, 14, 845. (7) Meuse, C. W.; Connor, J. T.; Richter, L. J.; Plant, A. L. Quantitative Analysis of the Molecular Structure of Thin Films Using Infrared Reflection Spectroscopy. In preparation. (8) Hautman, J.; Klein, M. L. J. Chem. Phys. 1989, 91, 4994. (9) Clarke, J. H. R.; Brown, D. Mol. Phys. 1986, 58, 815. (10) Martyna, G. J.; Klein, M. L.; Tuckerman, M. J. Chem. Phys. 1992, 97, 2635. (11) Schofield, P. Comput. Phys. Commun. 1973, 5, 17. (12) Beeman, D. J. Comput. Phys. 1976, 20, 130. (13) Brooks, B. R.; Bruccoliri, R. E.; Olafson, B. D.; States, D. J.; Swaminathan, S.; Karplus, M. J. Comput. Chem. 1983, 4, 187.
