Calculation of Noncontact Forces between Silica Nanospheres

Article pubs.acs.org/Langmuir

Calculation of Noncontact Forces between Silica Nanospheres Weifu Sun,† Qinghua Zeng,‡ and Aibing Yu*,† †

Laboratory for Simulation and Modelling of Particulate Systems, School of Materials Science and Engineering, The University of New South Wales, Sydney, NSW 2052, Australia ‡ School of Computing, Engineering and Mathematics, University of Western Sydney, Penrith, NSW 2751, Australia ABSTRACT: Quantification of the interactions between nanoparticles is important in understanding their dynamic behaviors and many related phenomena. In this study, molecular dynamics simulation is used to calculate the interaction potentials (i.e., van der Waals attraction, Born repulsion, and electrostatic interaction) between two silica nanospheres of equal radius in the range of 0.975 to 5.137 nm. The results are compared with those obtained from the conventional Hamaker approach, leading to the development of modified formulas to calculate the van der Waals attraction and Born repulsion between nanospheres, respectively. Moreover, Coulomb’s law is found to be valid for calculating the electrostatic potential between nanospheres. The developed formulas should be useful in the study of the dynamic behaviors of nanoparticle systems under different conditions. force.20 Yet, many efforts have been limited to metal nanoparticles.21 Experimentally, atomic force microscopes,22 surface force apparatus,23 and colloidal particle scattering24 have been applied to measure the vdW forces between spheres or between plates at a separation larger than 1.5 nm.25,26 However, at this stage of development, such experimental techniques are often hampered by the surface roughness of the sample considered and issues such as surface contamination, ambient vibrations, and limited distance resolution.27 On the other hand, molecular simulation could overcome the experimental difficulty and play a valuable role in the quantification of the interactions between nanoparticles. Molecular dynamics (MD) or Monte Carlo simulations28−34 have been used to study the solvation force, vdW attraction, or Born repulsion between nanoparticles in a solution or polymer composite. Yet, very few attempts have been made to calculate the interaction forces between nanoparticles in vacuum by MD simulations.34,35 Although the pairwise additivity inherent in both Hamaker approach and MD simulations neglects the many-body (MB) effects, it was reported that MB effects can be ignored for shape-isotropic silica nanocluster,36 thus MD simulations could provide accurate results of internanosphere forces. Using coarse-grained MD simulations, Lee et al.34 studied the interaction potentials between silica nanoparticles, but they did not take into account electrostatic interactions. Our recent work35 has also demonstrated that MD simulation is a useful approach to determining the interactions between nanoparticles because it would allow the flexible movement of surface atoms, differentiate the unsaturated surface atoms from the core atoms, and account for nonuniform distributions of

1. INTRODUCTION Understanding and prediction of the behaviors of particles, e.g., aggregation,1 adsorption,2 self-assembly,3 dispersion,4 packing5 and flowing,6 require a detailed knowledge of mutual interactions between particles, such as the long-range van der Waals (vdW) attraction, short-range Born repulsion, and electrostatic interaction. Quantitative description of such interactions is of particular importance in nanoparticles and nanostructured materials since they exhibit unique sizedependent properties different from their macroscopic counterparts. The most significant size effect occurs in nanoparticles of 1−10 nm in diameter, which has prompted researchers to develop novel nanostructured materials to meet the demands for the state-of-the-art electronic, optical, magnetic, and mechanical nanoscale devices. A large variety of nanoparticles with special shape and properties have been synthesized in the past years.7,8 Self-assembly of such nanoparticles offers many opportunities for generating a spectrum of structures and functional materials.9 One important factor that governs the behavior of these nanoparticle systems is the interaction forces between nanoparticles, which is, however, not well understood.10 Theoretically, the vdW attractive force between colloidal particles can be estimated by the well-established Hamaker approach,11 Lifshitz’s continuum theory,12,13 or Derjaguin’s approximation.14 For example, Hamaker11 and Feke et al.15 derived expressions for the vdW attractive and Born repulsive forces between macrospheres, respectively. However, the Hamaker approach does not generally apply to nanoparticles because it neglects the nature of discrete atomic structure and surface effects.16−19 Some efforts have been made to investigate factors that cause the deviation from the Hamaker approach and to propose an alternative approach to calculate the vdW © 2013 American Chemical Society

Received: August 1, 2012 Published: January 22, 2013 2175

dx.doi.org/10.1021/la305156s | Langmuir 2013, 29, 2175−2184

Langmuir

Article

Figure 1. Typical snapshots of two silica nanospheres during the head-on collision, captured at the times of (a) 11.50 ps, (b) 13.00 ps, (c) 14.60 ps, (d) 16.50 ps, (e) 17.15 ps and (f) 18.10 ps. Here, Vr is the relative velocity of the two nanospheres, with Vr > 0 indicating they approach toward each other and Vr < 0 indicating they depart from each other; the initial relative velocity Vr,0 = 300 m/s; d is the shortest surface-to-surface separation, i.e., d = r − 2R̅ , where r is the center separation, R̅ = 3.063 nm is the radius of nanosphere, with d < 0 indicating the overlapping of the two nanospheres; and the red beads denote oxygen atoms while the yellow ones represent silicon atoms.

wave functions. DFT semicore pseudopotential core treatment was used to include some degree of relativistic effects and reduce the computational cost. All the atoms were allowed to relax to their ground state configuration with the convergence criteria of 1 × 10−5 Ha and 0.002 Ha/Å for the total energy and forces, respectively. Since the atomic partial charges in COMPASS are derived from the first-principle calculation of electrostatic potential, to be consistent with the COMPASS force field, the electrostatic potential partial charges of surface atoms in the nanospheres were assigned using the population analysis implemented in DMol3 module. Then, MD simulations were performed on a pair of approaching silica nanospheres using the Materials Studio software package. During the MD simulations, the COMPASS force field within the package was used, which consists of a Lennard-Jones n−m potential (LJn−m, here n = 9 and m = 6) function for the Born repulsion and vdW attraction terms, respectively, and a Coulombic potential function for the electrostatic interaction, given by38

atoms and the small anisotropy of nanospheres at small separation. The interaction potentials of two nanoparticles can be calculated by manually assigning a set of separation between two static and rigid nanoparticles. Yet, when two nanoparticles interact with each other, their relative position and orientation are subjected to changing dynamically. The present work is to extend that work to use fully atomistic MD simulations to develop formulas for the calculation of the noncontact potentials between two dynamically interacting silica nanospheres. The aim is to generate equations to calculate the interaction forces between nanospheres, which can be directly used in the granular dynamics simulation3,5,6 of nanoparticle systems.

2. SIMULATION METHOD AND CONDITIONS In this work, nanospheres were generated from α-quartz, which has a stable low-energy structure of hexagonal symmetry (a = 4.91 Å, c = 5.40 Å) and a density of 2.65 g/cm3. Each unit cell has three silicon and six oxygen atoms, two types of Si−O bond lengths, and four types of O−O distances.37 Each nanosphere was constructed by cutting out atoms within a specified radius from such silica crystal as done in our previous work.35 Then, the isolated silicon or oxygen atoms on the surface were removed, which produced silica nanospheres with a radius in the range of 0.975 to 5.137 nm. As we will discuss in subsection 3.1, the atoms in a nanosphere will be classified into core atoms and surface atoms. The atomic partial charges of core atoms with full coordination were assigned based on the COMPASS force field, while those of surface atoms that do not bear full coordination were achieved by density functional theory (DFT) calculation. Specifically, DFT calculation was carried out for a neutral silica nanosphere with a radius of 0.60 nm using DMol3 module as implemented in the Materials Studio software package. Prior to geometry optimization, the total charge of the nanosphere was set to 0. During this optimization, the generalized gradient approximation scheme and Perdew−Burke−Ernzerhof exchange correlation functional were employed. Spin unrestricted geometry optimization calculations were performed using double numerical plus polarization basis set to expand the

Eij =

∑ i,j

qiqj rij

n m⎤ ⎡ m ⎛⎜ σij ⎞⎟ n ⎛⎜ σij ⎞⎟ ⎥ ⎢ + ∑ εij − ⎢⎣ n − m ⎜⎝ rij ⎟⎠ n − m ⎜⎝ rij ⎟⎠ ⎥⎦ i,j

(1)

where qi and qj are the charges of atoms i and j, separated by a distance rij, respectively; ε is the potential well depth for the interactions between atoms; and σij is the collision diameter of atoms i and j at which the potential becomes zero. In the COMPASS force field, the short-range repulsion models the short-range stabilizing effect of the Born repulsion, which originates from the strong repulsive forces between atoms as their electron shells interpenetrate each other. Such Born repulsion is explained according to Pauli’s exclusion principle, that is, no two electrons can have the same four quantum numbers. For a pair of dissimilar atoms, σij is given by eq 2 in terms of the σi and σj values of the individual atoms based on the Waldman−Hagler combination rule.39 2176


Langmuir

Article

⎡ (σ )6 + (σ )6 ⎤1/6 i j ⎥ σij = ⎢ ⎢⎣ ⎥⎦ 2

atomic structure on the surface of nanoparticles is expected to make a significant impact on their behaviors. From this point of view, it is necessary to define the surface and core atoms of a nanoparticle and hence the size of a nanoparticle. In this work, atoms on the surface of a nanosphere with unsaturated coordination are considered surface atoms while those in the core with full coordination are considered core atoms.42 Accordingly, the following structural parameters are defined: nanosphere core radius (Rcore), surface roughness (rms), and surface thickness (δ) (Figure 2). These parameters can be determined by analyzing the equilibrium structure obtained from the MD simulations (Table 1).

(2)

The εij value is given by ⎛ σ 3σ 3 ⎞ i j ⎟ εij = 2 εiεj ⎜⎜ 6 6⎟ σ ⎝ i + σj ⎠

(3)

The values of σi and σj for silicon and oxygen are 4.45 and 3.30 Å, respectively,38,40 thus the σ value for Si and O atom interactions σSi−O = 4.07 Å; the values of εi and εj for silicon and oxygen are 0.198 and 0.080 kcal/mol, respectively, giving rise to an εSi−O value of 0.088 kcal/mol. Additionally, the Hamaker constant (A = π2C/v2, where C = 3εSi−Oσ6Si−O is the vdW attraction interaction parameter, and v = 4/3π(σSi−O/2)3 is the atom volume of about 35.3 Å3) is estimated to be 6.6 × 10−20 J, which is consistent with the experimental value of 6.5− 8.86 × 10−20 J, given by Bergstrom.41 MD simulations were performed using the NVE ensemble (i.e., constant number of atoms, constant volume, and constant energy) at an initial temperature of 298.0 K, which allows both the temperature and stress of the system to change under the external load. Note that for the purpose of comparison, some cases using the NVT ensemble (i.e., constant number of atoms, constant volume, and constant temperature) were also carried out in this work. Prior to performing MD simulations, each nanosphere was fully relaxed using the NVT ensemble at 298.0 K, then two fully relaxed silica nanospheres were allowed to move toward each other at an equal but opposite initial velocity. Numerical integration was performed using the velocity Verlet integration algorithm and a time step of 1 fs (1.0 × 10−15 s). Energies and other particle information were recorded every 100 steps as a function of the shortest surfaceto-surface separation distance (d) between the two nanospheres in an output trajectory file. The initial d value between the surfaces of the two silica nanospheres was set to 10.0 nm. Figure 1 shows a typical sequence of approaching, colliding, and departing of a pair of nanospheres in a MD simulation. Initially, the two nanospheres move close to each other (Figure 1a) at an initial relative velocity of Vr,0 =300 m/s. Then, they continue moving toward each other and reach a position of close proximity (Figure 1b). After that, the overlapping of two nanospheres occurs (Figure 1c). Figures 1d, e, f demonstrate the following dynamic departing process. It is shown that the two nanoparicles can rebound out of the reach of each other after the overlap if a very large initial velocity is applied. Changing the initial relative velocity Vr,0 will produce different collision dynamics, as discussed in subsection 3.2. The interaction energy, EInter, between two nanospheres takes the form35 E

Inter

12

1

=E −E −E

2

Figure 2. A schematic structure of a spherical nanoparticle with the defined parameters. Note that not all atoms in the core are shown.

Table 1. Structural Parameters of Silica Nanospheres R0 (Å)

10

20

30

40

50

Rcore (Å) R̅ (Å) rms (Å) rms/R δ (Å) δMax (Å)

8.78 9.75 0.587 0.067 1.557 3.03

18.90 20.06 0.677 0.034 1.837 3.42

29.60 30.63 0.643 0.021 1.679 3.30

40.20 41.20 0.632 0.016 1.632 3.80

50.30 51.37 0.674 0.013 1.774 3.72

Note: R0 is the radius of an arbitrarily cut nanosphere.

In this work, the core radius (Rcore) of a nanosphere is defined as Rcore ≡ max{R jcore}

(5)

The surface roughness (rms) of a nanosphere is defined as the root mean squared:43 n

rms =

(4)

where E1, E2 and E12 are the potential energies of particle 1, particle 2, and the two particles, respectively. The individual potential contributions (e.g., the vdW attraction, Born repulsion, and electrostatic interaction) can also be obtained in a similar way.

∑i = 1 (R isurf − R̅ )2 n

(6)

The maximum surface thickness (δMax) is defined as δMax = max{R isurf } − max{R jcore}

(7)

The effective surface thickness (δ) is written as δ ≡ R̅ + rms − Rcore

3. RESULTS AND DISCUSSION 3.1. Characterization of Atomic Structure of Nanospheres. Many properties of materials depend on the characteristics of their atomic structure. In particular, the

(8)

The average radial distance (R̅ ) of surface atoms is defined as the radius of a nanosphere (R) and used in the Hamaker approach. 2177


Langmuir

Article n

R̅ =

∑i = 1 R isurf

Rsurf i

produced, which indicates that interparticle potentials can also be considered to be independent of temperature. 3.2. Effect of Initial Relative Velocity on Interparticle Potentials. In this subsection, the effects of initial relative velocity Vr,0, as well as other factors, will be studied, and the relationships among parameters such as Vr,0, surface separation d, and the minimum surface separation dmin will be examined. The aim is to identify variables that should be taken into account in the calculation of interparticle potentials. Figure 4 shows the interparticle potentials as a function of d for silica nanospheres of R = 3.063 nm for various initial relative

(9)

n

Rcore j

where and are the radial distances of surface atom i and core atom j from the nanosphere center, respectively. n is the number of surface atoms. The effects of nanosphere size on the above structural parameters are shown in Table 1. rms, δ and δMax vary within a certain limit, which can be regarded as independent of nanosphere radius. rms is about 0.64 ± 0.04 Å, which is significantly smaller than that of atomic force microscopic measurements for microparticles or nanoparticles.22 This can be attributed to the fact that for the nanosphere attained in the present study, the surface roughness at the atomic scale arises only from its discrete atomic structure,19 rather than from the nonspherical surface profile. δ is about 1.70 ± 0.11 Å, and δMax is about 3.45 ± 0.31 Å. However, the relative surface roughness decreases (rms/R) significantly with the increase in nanosphere size (Figure 3).

Figure 4. Interparticle potentials as a function of the surface separation of two neutral silica nanospheres (R = 3.063 nm) approaching each other at different initial relative velocities, Vr,0.

velocities. The five curves are almost identical except for the minor discrepancy at small separations, indicating that interparticle potentials are independent of the initial relative velocity prior to the contact of the two nanospheres. The relative velocity Vr changes when both nanospheres are approaching to each other (Figure 5). In particular, for a given Figure 3. Relative surface roughness as a function of the radius of nanospheres R.

Moreover, the influence of temperature on particle’s structural parameters was also explored using the NVT ensemble for comparison. The results (Table 2) indicate that Table 2. Structural Parameters of Silica Nanospheres of R = 3.063 nm at Different Temperature T T (K)

1

150

298

450

600

Rcore (Å) R̅ (Å) rms (Å) rms/R δ (Å) δMax (Å)

29.78 30.71 0.554 0.018 1.483 2.57

29.76 30.70 0.572 0.019 1.508 2.84

29.60 30.63 0.643 0.021 1.679 3.30

29.56 30.60 0.654 0.021 1.697 3.27

29.47 30.56 0.686 0.022 1.772 3.63

Figure 5. Relative velocity as a function of surface separation of two neutral nanospheres (R = 3.063 nm) during the approaching process at different initial relative velocities.

the change of temperature from 1 to 600 K only leads to a slight change of particle size (R) and surface roughness (rms). Therefore, the influence of temperature on the structure of a silica nanosphere can be ignored, and the defined structural parameters are considered to be independent of temperature, at least in the temperature range considered. As mentioned in section 2, the present MD simulations were performed using a NVE ensemble. However, if changed to a NVT ensemble, as will be discussed in subsection 3.2, almost identical results were

initial relative velocity Vr,0, as d decreases, Vr increases gradually to a peak because of the interparticle vdW attractive force, followed by a sharp decrease because of the interparticle Born repulsion force. The peak becomes less obvious when the initial relative velocity Vr,0 is higher than a certain value, about 200 m/ s. Moreover, the minimum gap, dmin, where the relative velocity becomes almost zero and starts to change from positive to 2178


Langmuir

Article

Figure 6. Dependence of dmin on (a), Vr,o and (b) R.

It is known that due to the different bonding state and chemical environment, with a decrease in particle size, the surface atoms represent an increased contribution and tend to affect the behavior of particles. The effect of the surface of particles under the static circumstance has been tested in the previous work,35 showing that this could lead to the fluctuation of interaction potentials, particularly at small separations. During the present dynamic head-on impact process, the surface atoms may be subject to change dynamically. Therefore, as part of the present study, the effects of surface atoms due to different relative orientations were also studied, as done in the previous study.35 Figure 8 shows the interparticle potentials as a function of surface separation, where the two nanospheres collided with

negative, depends on the initial relative velocity and demonstrates a linear relationship (Figure 6a): at a small Vr,0, dmin > 0 indicating that the two particles cannot be in contact because of the strong Born repulsion force; at large Vr,0, dmin < 0 becomes possible because the inertial force overcomes the repulsion force leading to the collision between the two particles. Not surprisingly, therefore, dmin is also dependent on particle size. As shown in Figure 6b, for a given Vr,0, dmin decreases linearly with the increase of R in the range of 0.975 to 5.137 nm. Furthermore, the effects of simulation ensembles, temperature, and relative orientations of two nanospheres on the interaction potentials were examined. In Figure 7, the first four

Figure 7. Interparticle potentials as a function of surface separation of two neutral silica nanospheres (R = 3.063 nm) approaching each other at an initial relative velocity of 300 m/s under different ensembles or temperatures.

Figure 8. Interparticle potentials as a function of the surface separation d of two neutral silica nanospheres (R = 3.063 nm) approaching each other at different configurations.

curves were obtained under the NVT ensemble and the last one from the NVE ensemble at an initial relative velocity of 300 m/ s. The results indicate that both NVE and NVT give almost the same results, and the effect of temperature on the interparticle potentials can be ignored. Note that temperature is a thermodynamic quantity and is largely meaningful only at equilibrium. However, the interaction between two nanospheres in the present MD simulation is a nonequilibrium dynamic process in which the NVE ensemble is usually recommended.44 Furthermore, the effect of temperature ranging from 1 to 600 K on the structural parameters of silica nanospheres can be ignored as shown in Table 2. Therefore, the NVE ensemble is chosen in the following MD simulations.

each other with different (relative) orientations. In doing so, prior to head-on dynamic impact, one particle remained unchanged while the other rotated clockwise by an angle of 0°, 90°, 180° or 270° along its central axis. Similar to the effect of relative velocity, the results show that the four curves are almost identical except for a small discrepancy at small separations, indicating that interparticle potentials are independent of interacting configurations before they are close to contact. In fact, changing relative orientation is equivalent to altering interacting surfaces and hence surface roughness. As observed from Figure 8, the influence of surface roughness could be largely ignored. It seems that there is a conflict with the results reported by Qin et al.29 This is because in their work, there is 2179


Langmuir

Article

on noncontact forces. The forces after contact will be considered in a separate study. Although there are several continuum models available to describe interparticle forces, such as the Lifshitz theory12 and Derjaguin approximation,14 the Hamaker approach11 is widely used due to its convinence. Therefore, the interparticle forces obtained from the MD simulations are compared with those from the conventional Hamaker approach. According to the Hamaker approach, for two solid spherical particles of radii R1 and R2, their vdW attraction interaction potential at a separation distance d is given by15

an additional force, i.e., the solvation force due to the ordered solvents on the particle surface,9 and this force is very sensitive to the surface roughness. In fact, because the nanospheres produced are very symmetrical, different orientations should not change the surface atoms much. 3.3. Interparticle vdW Attraction and Born Repulsion Forces. The results in subsection 3.2 suggest that the interparticle potentials between silica nanospheres are mainly affected by two variables: surface separation d and particle size R. To be quantitative, MD simulations were performed to obtain the interparticle potentials between neutral silica nanospheres with a radius in the range of 0.975 to 5.137 nm as a function of the surface separation d (Figure 9). It can be

A⎡ 2λ 2λ E HvdW = − ⎢ 2 + 2 6 ⎣ r′ − (1 + λ)2 r′ − (1 − λ)2 + ln

r′2 − (1 + λ)2 ⎤ ⎥ r′2 − (1 − λ)2 ⎦

(10)

Similarly, integrating the repulsive (or ‘n to the power’) part of LJn−m potential of eq 1 over two spherical colloidal particles with n = 9 and m = 6, the Hamaker approach yields the Born repulsion interaction potential as follows:15 ⎛ σ ⎞n − 6 (n − 8)! 1 E HBorn = 4A⎜ ⎟ ⎝ R1 ⎠ (n − 2)! r′ ⎡ − r′2 − (n − 5)(λ − 1)r′ − (n − 6)[λ 2 − (n − 5)λ + 1] ×⎢ ⎣ (r′ − 1 + λ)n − 5

Figure 9. Interparticle potentials as a function of separation d between two identical silica nanospheres of radius ranging from 0.975 to 5.137 nm. Note that the radii 2.475 and 2.875 nm are the geometric mean radii of a silica nanosphere of 2.006 nm and its interacting nanosphere of 3.063 and 4.120 nm, respectively.

+

− r′2 + (n − 5)(λ − 1)r′ − (n − 6)[λ 2 − (n − 5)λ + 1] (r′ + 1 − λ)n − 5

+

r′2 + (n − 5)(λ + 1)r′ + (n − 6)[λ 2 + (n − 5)λ + 1] (r′ + 1 + λ)n − 5

+

r′2 − (n − 5)(λ + 1)r′ + (n − 6)[λ 2 + (n − 5)λ + 1] ⎤ ⎥ ⎦ (r′ − 1 − λ)n − 5

(11)

Here A is the Hamaker coefficient, λ = R2/R1, and r′ = (R1 + R2 + d)/R1 is the center-to-center separation, made dimensionless by R1. As shown in Figure 10, the ratios of the interparticle forces obtained from the MD simulations and the Hamaker approach indicate that with the increase in surface separation, the ratios for both vdW attraction and Born repulsion forces first increase sharply and then decrease drastically. The maximum occurs at d ≈ 0.4 nm. The results clearly demonstrate that the ratio of the interparticle forces obtained from the MD simulation and the Hamaker approach is a function of particle size and surface separation. In particular, in most cases, the vdW attraction and Born repulsion forces from the MD simulations are much larger than those from the Hamaker approach. This can be attributed to the surface effect and the nature of discrete atomic structures that have been ignored in the Hamaker approach.16−19 On one hand, nanoparticle should not be always treated as rigid (incompressible) as in the Hamaker approach, because the surface atoms, which is relatively high for nanoparticles, are unsaturated having a certain degree of flexibility. This will result in a significant contribution to the interaction between nanoparticles, especially at small separations.45 On the other hand, in the Hamaker approach, a particle is assumed to behave as a continuum medium with uniform (density) distribution of atoms. Yet, such an assumption results in a large error for nanoparticles because the amount of surface atoms is relatively high.16,18 In the more rigorous Lifshitz theory, the expressions of interparticle forces (i.e., eqs 10 and 11) remain unchanged,

observed that both the interparticle vdW attraction and Born repulsion potentials vary with particle size in a systematic manner. That is, for a given separation, the magnitude of each of the interparticle potentials increases with particle size, particularly at small interparticle separations. At the atomic scale, it is often assumed that two atoms are in contact or bonding if their separation is less than a certain distance.45−47 As a result of the competing forces of attraction and repulsion between individual atoms of two interacting bodies, two solid surfaces may have an equilibrium separation.46,48 In fact, a separation ranging from 1.6547 to 5.0 Å49 is recognized as an important parameter to calculate friction forces between particles in the Maugis−Dugdale model.50 This treatment can be supported by the present MD results to some degree. As shown in Figure 5, particle collision occurs only when the initial relative velocity is very high. Figure 9 further suggests that for the silica nanospheres considered, when surface separation d is smaller than about 2 nm, the attraction becomes effective. When d = 0.2 nm, the interparticle potential reaches a minimum where the vdW attraction force equals the Born repulsion force, giving an equilibrium separation for the two nanospheres when there is no external disturbance. After that point, repulsion gradually becomes effective. However, real repulsion may occur only when two nanospheres collide and overlap. At this stage, contact forces such as adhesion51 and elastic restoring force due to deformation52 should be considered. As mentioned earlier, the present work is focused 2180


Langmuir

Article

vdW FMD

⎧ −4d 0.6 )k vdWFHvdW , d ≥ 0.4 nm ⎪(1 + 26e =⎨ ⎪ e3.2(0.4 − d)0.65F vdW MD, d = 0.4nm , 0 < d < 0.4 nm ⎩

(13)

It was found that similar treatments can also be applied to the Born repulsion force. Thus, by fitting the MD results, the following equation is obtained to estimate the Born repulsion forces: Born FMD

⎧ −6d 0.5 )kBornFHBorn , d ≥ 0.4 nm ⎪(1 + 701.6e =⎨ 0.75 ⎪ e 4.65(0.4 − d) F Born MD, d = 0.4nm , 0 < d < 0.4 nm ⎩ (14)

where kBorn, the asymptotic ratio of the Born repulsion force calculated by the MD simulation to that calculated by the Hamaker approach, is described by the following equation: ⎛ rms ⎞0.10 ⎟ ln kBorn = 2.09⎜ ⎝ R ⎠

as shown in the inset of Figure 10b, kBorn decreases with the increase of particle size. Both kvdW and kBorn become unity when R becomes infinitely large. This reflects the fact that if particle size is large enough, the Hamaker approach is valid. As noted earlier, the ratio of the interparticle forces obtained from the MD simulation FMD and the Hamaker approach FH reaches its maximum at d = 0.4 nm. The point is used as a reference in formulating the equation for d < 0.4 nm. The combination of eqs 10 and 13, or eqs 11 and 14, allows the calculation of the corresponding vdW interaction forces, or the Born repulsion forces between two silica nanospheres. As shown in Figure 11, the so calculated results can match the MD simulated results reasonably well. The forces of small nanospheres (R < 1 nm) somehow fluctuate when separation d is small, which may need to be investigated further in the future. On the other hand, it should be noted that an effort has been made to demonstrate that the above equations may be applied to nanospheres with a radius much larger than 10 nm. This is achieved by use of eq 12 or 15, respectively, for the vdW attraction and Born repulsion forces, so that when R is large enough, eqs 13 and 14 produce results the same as those given by the Hamaker approach, i.e., eqs 10 and 11. Note that the Hamaker approach can be used for large particles.17,47 Figure 12 shows the results, generated by the newly proposed equations. It is clear that decreasing separation d increases both the attraction and repulsion forces, and the larger the particle size, the higher the magnitudes of the forces. This is consistent with our general understanding of the two forces. While taking the advantage of the Hamaker approach, the proposed equations may have also taken its disadvantage. For example, it may not be accurate when the so-called retardation effect becomes significant at separations larger than the “characteristic wavelength”, which is usually assumed to be 100 nm.53 According to the conventional Hamaker approach, when d is approaching zero, the repulsion force becomes dominant, giving an infinitely large value, hence two particles cannot be in contact at all. This is obviously contrary to reality because collision between particles is known to occur. To overcome this problem, a cutoff separation has to be introduced when using eq 10 or 11 to calculate the vdW attraction or Born repulsion force.5,6 The treatment is however a bit arbitrary, although

Figure 10. The ratios of interparticle forces between two silica nanospheres obtained from the MD simulations to those from the Hamaker approach: (a) the vdW attraction; and (b) the Born repulsion. The solid symbols represent the results for two identical nanospheres, while the open ones represent the results for two nanospheres with different size. The insets show the variation of the asymptotic ratio kvdW and kBorn as a function of particle radius R.

but the Hamaker coefficient is calculated in terms of bulk properties such as dielectric constants and refractive indices; strictly speaking, this coefficient is not constant but varies with the surface separation.13,47 Figure 10 also shows that the ratios decrease sharply with the increase of surface separation, and finally become constant asymptotically. The asymptotic ratio of the forces, denoted here as kvdW, decreases with the increase of particle size as shown in the inset of Figure 10a. The results are related to relative surface roughness, given by ⎛ rms ⎞0.18 ⎟ ln k vdW = 3.16⎜ ⎝ R ⎠

(15)

(12)

where the averaged surface roughness of rms = 0.064 nm was assumed in formulating eq 12. On the basis of the MD results, an equation has been formulated to describe the vdW attraction force between two nanospheres, given as 2181


Langmuir

Article

Figure 11. Interparticle forces between silica nanospheres obtained from MD simulations: (a) the vdW attraction force and (b) the Born repulsion force. The solid symbols represent results for two identical nanospheres, while the open ones represent the results for two nanospheres of different size. The solid lines represent fitting curves that correspond to the solid symbol colors, and the dotted lines represent fitting curves that correspond to the open symbol colors.

Figure 12. Interparticle vdW attraction (a) and Born repulsion (b) forces as a function of d between silica spheres of different radii R. Note that “negative” represents attraction force here, and lines 1−6 correspond to particle radii of 1−5000 nm, respectively.

first determined. The interparticle ESPs are then calculated by the MD simulation and the Coulomb’s law, respectively. As described in section 2, DFT calculation was performed on a silica nanosphere of 1.20 nm in diameter, where partial charge was determined for each atom (Figure 13). However, for large particles containing more than 100 atoms, DFT calculation becomes expensive. Therefore, different atoms in such a large particle were defined according to their status of atomic

shown to be effective. This arbitrary treatment can now be removed by use of the proposed modified Hamaker expressions. The MD results also confirm that an infinitely large vdW attraction or Born repulsion force would not occur at the molecular scale. Moreover, the above expressions may be extended to two silica particles with different radii R1 and R2 by replacing R with the following geometrical mean R = (R1R2)1/2 and/or two silica particles interacting across a (third) medium:47 A132 ≈ ( A11 −

A33 )( A 22 −

A33 )

(16)

where A11, A22 and A33 are the Hamaker constants of particles 1, 2, and medium 3, respectively (in this case, A11 = A22). The validity of the geometrical mean approximation for silica nanospheres can be seen from the results in Figures 10 and 11 for the cases with mean radius of 2.475 and 2.875 nm. The treatmemts have been widely used for particles;47 here we just confirm their applicability to nanoparticles. Their validity probably requires further studies in the future. 3.4. Interparticle Electrostatic Potentials. Due to the significant amount of unsaturated atoms on the surface, nanoparticles usually possess a certain amount of electrostatic charge. As part of this work, the possibility of applying the Coulomb’s law to estimate the electrostatic potential (ESP) between two nanospheres is examined. Specifically, the partial charges of atoms and the net charge of each nanosphere are

Figure 13. DFT calculated partial charge of atoms in a silica nanosphere with a radius of 0.6 nm (side view). 2182


Langmuir

Article

Table 3. Average Partial Charges (q) of Atoms Calculated from DFT Si1 0.29 ± 0

type of atom charge q (e)

Si2 0.63 ± 0.40

Si3 0.27 ± 0.06

2.006 −18.94

3.063 −33.53

4.120 6.38

5.137 10.44

Once the partial charges of atoms are known, MD simulation can be performed to calculate the ESPs between two nanospheres. On the other hand, the ESPs can also be calculated based on the net charge of nanospheres. In this study, we assumed that the net charges are uniformly distributed within a nanosphere, and then ESPs between two nanospheres is given by Coulomb’s law, ES ECoulomb =

1 4πεε0

■

Q 1Q 2 r

O2 0.31 ± 0.15

4. CONCLUSIONS MD simulations have been performed to determine the interaction forces of two silica nanospheres including the vdW attraction, Born repulsion, and electrostatic interaction. The results show that the interaction potentials of two nanospheres are independent of their initial relative velocity, temperature, and relative orientation, although the overlap between nanospheres varies with initial relative velocity and particle size. The vdW attraction and Born repulsion forces obtained from the MD simulations are different from those determined by the conventional Hamaker approach. On the basis of the MD results, modified Hamaker expressions have been developed to calculate the vdW attraction and Born repulsion potentials between silica nanoparticles. Moreover, the Coulomb’s law is shown to be applicable to calculate the ESP of two nanoparticles. The modified Hamaker expressions, together with those for other contact forces, can be applied to discrete element simulations for studying the dynamic behaviors of silica nanoparticles, as done for coarse particles.5,6 Their application in the study of the packing and flow of silica nanospheres will be reported in the future.

Table 4. Total Charges (Q) of Silica Nanospheres with Different Radii 0.975 −10.98

O1 0.32 ± 0.12

long as the net charge of a particle can be determined experimentally or numerically, eq 17 can be used to calculate the ESP or force between two nanospheres.

coordination (i.e., bonding conditions). For instance, silicon atoms bonded with 1, 2, 3, and 4 oxygen atoms are denoted by Si1, Si2, Si3, and Si4, respectively. Similarly, oxygen atoms bonded with 1 and 2 silicon atoms are denoted by O1 and O2, respectively. For each type of atoms, the ESP partial charges were averaged and assigned to surface atoms according to their bonding conditions (Table 3). The partial charges of Si = 0.89 e, and O = −0.445 e, which were adopted from the COMPASS force field and originally derived based on ESP energies, were assigned to the core atoms. The net charge of each nanosphere could then be obtained based on the partial charges of all atoms as listed in Table 4.

particle radius R (nm) total charge Q (e)

Si4 0.85 ± 0.21

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. Phone: +61-2-9385-4429. Fax: +61-2-9385-5956.

(17)

where r = d + 2R is the center-to-center nanosphere separation, Q1 and Q2 are the net charges of both nanospheres, respectively, and ε0 and ε are the permittivities of the vacuum and intervening medium, respectively. Note that here ε0 ≈ 8.85418782 × 10−12 F/m and ε = 1. Figure 14 shows that the results obtained from the MD simulations agree well with those from the Coulomb’s law. Therefore, Coulomb’s law can still be applied to calculate the ESP of two nanospheres. That is, as

Notes

The authors declare no competing financial interest.

■ ■

ACKNOWLEDGMENTS The authors would like to thank the Australian Research Council (ARC) for the financial support for this work. REFERENCES

(1) Bharti, B.; Meissner, J.; Findenegg, G. H. Aggregation of Silica Nanoparticles Directed by Adsorption of Lysozyme. Langmuir 2011, 27, 9823−9833. (2) Wi, H. S.; Cingarapu, S.; Klabunde, K. J.; Law, B. M. Nanoparticle Adsorption at Liquid-Vapor Surfaces: Influence of Nanoparticle Thermodynamics, Wettability, and Line Tension. Langmuir 2011, 27, 9979−9984. (3) Yu, A. B.; An, X. Z.; Zou, R. P.; Yang, R. Y.; Kendall, K. SelfAssembly of Particles for Densest Packing by Mechanical Vibration. Phys. Rev. Lett. 2006, 97, 265501−1−4. (4) Liu, J.; Gao, Y. G.; Cao, D. P.; Zhang, L. Q.; Guo, Z. H. Nanoparticle Dispersion and Aggregation in Polymer Nanocomposites: Insights from Molecular Dynamics Simulation. Langmuir 2011, 27, 7926−7933. (5) Dong, K. J.; Yang, R. Y.; Zou, R. P.; Yu, A. B. Role of Interparticle Forces in the Formation of Random Loose Packing. Phys. Rev. Lett. 2006, 96, 145505−1−4. (6) Zhu, H. P.; Zhou, Z. Y.; Yang, R. Y.; Yu, A. B. Discrete Particle Simulation of Particulate Systems: A Review of Major Applications and Findings. Chem. Eng. Sci. 2008, 63, 5728−5770. (7) Zhang, Z. L.; Glotzer, S. C. Self-Assembly of Patchy Particles. Nano Lett. 2004, 4, 1407−1413.

Figure 14. Interparticle ESPs between two identical silica nanospheres with different radii in the range of 0.975 to 5.137 nm: solid symbols, data obtained from the MD simulations; solid lines, calculations according to Coulomb’s law. 2183


Langmuir

Article

(8) Xia, Y. N.; Halas, N. J. Shape-Controlled Synthesis and Surface Plasmonic Properties of Metallic Nanostructures. MRS Bull. 2005, 30, 338−344. (9) Min, Y. J.; Akbulut, M.; Kristiansen, K.; Golan, Y.; Israelachvili, J. The Role of Interparticle and External Forces in Nanoparticle Assembly. Nat. Mater. 2008, 7, 527−538. (10) Jiang, X. C.; Zeng, Q. H.; Chen, C. Y.; Yu, A. B. Self-Assembly of Particles: Some Thoughts and Comments. J. Mater. Chem. 2011, 21, 16797−16805. (11) Hamaker, H. C. The London−van der Waals Attraction between Spherical Particles. Physica 1937, 4, 1058−1072. (12) Dzyaloshinskii, I. E.; Lifshitz, E. M.; Pitaevskii, L. P. The General Theory of van der Waals Forces. Adv. Phys. 1961, 10, 165− 209. (13) Pashley, R. M. DLVO and Hydration Forces between Mica Surfaces in Li+, Na+, K+, and Cs+ Electrolyte Solutions - A Correlation of Double-Layer and Hydration Forces with Surface Cation-Exchange Properties. J. Colloid Interface Sci. 1981, 83, 531−546. (14) Schiller, P.; Kruger, S.; Wahab, M.; Mogel, H. J. Interactions between Spheroidal Colloidal Particles. Langmuir 2011, 27, 10429− 10437. (15) Feke, D. L.; Prabhu, N. D.; Mann, J. A. A Formulation of the Short-Range Repulsion between Spherical Colloidal Particles. J. Phys. Chem. 1984, 88, 5735−5739. (16) Girifalco, L. A.; Hodak, M.; Lee, R. S. Carbon Nanotubes, Buckyballs, Ropes, and a Universal Graphitic Potential. Phys. Rev. B 2000, 62, 13104−13110. (17) Parsegian, V. A. van der Waals Forces: A Handbook for Biologist, Chemists, Engineers, and Physicist; Cambridge University Press: Cambridge, U.K., 2005. (18) Bishop, K. J. M.; Wilmer, C. E.; Soh, S.; Grzybowski, B. A. Nanoscale Forces and Their Uses in Self-Assembly. Small 2009, 5, 1600−1630. (19) Luan, B. Q.; Robbins, M. O. The Breakdown of Continuum Models for Mechanical Contacts. Nature 2005, 435, 929−932. (20) Kudryavtsev, Y. V.; Gelinck, E.; Fischer, H. R. Theoretical Investigation of van der Waals Forces between Solid Surfaces at Nanoscales. Surf. Sci. 2009, 603, 2580−2587. (21) Klimov, V. V.; Lambrecht, A. van der Waals Forces between Plasmonic Nanoparticles. Plasmonics 2009, 4, 31−36. (22) Jallo, L. J.; Chen, Y. H.; Bowen, J.; Etzler, F.; Dave, R. Prediction of Inter-particle Adhesion Force from Surface Energy and Surface Roughness. J. Adhes. Sci. Technol. 2011, 25, 367−384. (23) Anderson, T. H.; Donaldson, S. H.; Zeng, H. B.; Israelachvili, J. N. Direct Measurement of Double-Layer, van der Waals, and Polymer Depletion Attraction Forces between Supported Cationic Bilayers. Langmuir 2010, 26, 14458−14465. (24) Wu, X.; vandeVen, T. G. M. Retarded van der Waals Forces between Triblock-Coated Latex Spheres. Langmuir 1996, 12, 6291− 6294. (25) Mohideen, U.; Roy, A. Precision Measurement of the Casimir Force from 0.1 to 0.9 mum. Phys. Rev. Lett. 1998, 81, 4549−4552. (26) Israelachvili, J. N.; Tabor, D. The Measurement of van der Waals Dispersion Forces in Range 1.5 to 130 nm. Proc. R. Soc. London A 1972, 331, 19−38. (27) Yaminsky, V. V.; Stewart, A. M. Interaction of Glass Surfaces in Air: Dispersion Forces in the Retarded Regime. Langmuir 2003, 19, 4037−4039. (28) Bedrov, D.; Smith, G. D.; Smith, J. S. Matrix-Induced Nanoparticle Interactions in a Polymer Melt: A Molecular Dynamics Simulation Study. J. Chem. Phys. 2003, 119, 10438−10447. (29) Qin, Y.; Fichthorn, K. A. Molecular-Dynamics Simulation of Forces between Nanoparticles in a Lennard-Jones Liquid. J. Chem. Phys. 2003, 119, 9745−9754. (30) Patel, N.; Egorov, S. A. Interactions between Nanocolloidal Particles in Polymer Solutions: Effect of Attractive Interactions. J. Chem. Phys. 2005, 123, 144916−1−10.

(31) Marla, K. T.; Meredith, J. C. Simulation of Interaction Forces between Nanoparticles: End-Grafted Polymer Modifiers. J. Chem. Theory Comput. 2006, 2, 1624−1631. (32) Qin, Y.; Fichthorn, K. A. Solvation Forces between Colloidal Nanoparticles: Directed Alignment. Phys. Rev. E 2006, 73, 020401−1− 4. (33) Lane, J. M. D.; Ismail, A. E.; Chandross, M.; Lorenz, C. D.; Grest, G. S. Forces between Functionalized Silica Nanoparticles in Solution. Phys. Rev. E 2009, 79, 1−4. (34) Lee, C. K.; Hua, C. C. Nanoparticle Interaction Potentials Constructed by Multiscale Computation. J. Chem. Phys. 2010, 132, 224904−1−9. (35) Zeng, Q. H.; Yu, A. B.; Lu, G. Q. Evaluation of Interaction Forces between Nanoparticles by Molecular Dynamics Simulation. Ind. Eng. Chem. Res. 2010, 49, 12793−12797. (36) Kim, H. Y.; Sofo, J. O.; Velegol, D.; Cole, M. W.; Lucas, A. A. van der Waals Forces between Nanoclusters: Importance of ManyBody Effects. J. Chem. Phys. 2006, 124, 74504−1−4. (37) Keskar, N. R.; Chelikowsky, J. R. Anomalous Elastic Behaviour in Crystalline Silica. Phys. Rev. B 1993, 48, 16227−16233. (38) Sun, H. COMPASS: An Ab Initio Force-Field Optimized for Condensed-Phase Applications - Overview with Details on Alkane and Benzene Compounds. J. Phys. Chem. B 1998, 102, 7338−7364. (39) Waldman, M.; Hagler, A. T. New Combination Rules for RareGas van der Waals Parameters. J. Comput. Chem. 1993, 14, 1077− 1084. (40) Sun, H.; Rigby, D. Polysiloxanes: Ab Initio Force Field and Structural, Conformational and Thermophysical Properties. Spectrochim. Acta, Part A: Mol. Biomol. Spectrosc. 1997, 53, 1301−1323. (41) Bergstrom, L. Hamaker Constants of Inorganic Materials. Adv. Colloid Interface Sci. 1997, 70, 125−169. (42) Van Hoang, V. Molecular Dynamics Simulation of Amorphous SiO2 Nanoparticles. J. Phys. Chem. B 2007, 111, 12649−12656. (43) Rabinovich, Y. I.; Adler, J. J.; Ata, A.; Singh, R. K.; Moudgil, B. M. Adhesion between Nanoscale Rough Surfaces - I. Role of Asperity Geometry. J. Colloid Interface Sci. 2000, 232, 10−16. (44) Zhang, Q.; Qi, Y.; Hector, L. G.; Cagin, T.; Goddard, W. A. Origin of Static Friction and Its Relationship to Adhesion at the Atomic Scale. Phys. Rev. B 2007, 75. (45) Cheng, S. F.; Robbins, M. O. Defining Contact at the Atomic Scale. Tribol. Lett. 2010, 39, 329−348. (46) Johnson, K. L. Contact Mechanics; Cambridge University Press: Cambridge, U.K., 1985; Chapter 5. (47) Israelachvili, J. N. Intermolecular and Surface Forces, 3rd ed.; Academic Press: New York, 2011; Chapter 13. (48) Carpick, R. W.; Ogletree, D. F.; Salmeron, M. A General Equation for Fitting Contact Area and Friction vs Load Measurements. J. Colloid Interface Sci. 1999, 211, 395−400. (49) Chandross, M.; Lorenz, C. D.; Stevens, M. J.; Grest, G. S. Simulations of Nanotribology with Realistic Probe Tip Models. Langmuir 2008, 24, 1240−1246. (50) Maugis, D. Adhesion of Spheres: The JKR-DMT Transition Using a Dugdale Model. J. Colloid Interface Sci. 1992, 150, 243−269. (51) Prokopovich, P.; Perni, S. Multiasperity Contact Adhesion Model for Universal Asperity Height and Radius of Curvature Distributions. Langmuir 2010, 26, 17028−17036. (52) Wang, X. D.; Shen, Z. X.; Zhang, J. L.; Jiao, H. F.; Cheng, X. B.; Ye, X. W.; Chen, L. Y.; Wang, Z. S. Contact between Submicrometer Silica Spheres. Langmuir 2010, 26, 5583−5586. (53) Gregory, J. Approximate Expressions for Retarded van der Waals Interaction. J. Colloid Interface Sci. 1981, 83, 138−145.

2184


Calculation of Noncontact Forces between Silica Nanospheres

Recommend Documents