Molecular Dynamics of Spherical Nanoparticles in Dense Polymer Melts

Article pubs.acs.org/JPCB

Molecular Dynamics of Spherical Nanoparticles in Dense Polymer Melts Alessandro Patti*,† Institute of Advanced Chemistry of Catalonia (IQAC−CSIC) and CIBER de Bioingeniería, Biomateriales y Nanomedicina (CIBER-BBN), Jordi Girona 18-26, 08034 Barcelona,Spain ABSTRACT: By performing molecular dynamics simulations, we investigate the structural and dynamical properties of polymer melts containing probe spherical nanoparticles. Generally speaking, the behavior of these polymer nanocomposites is strongly affected by the interaction strength established between the nanoparticles and the chain monomers and by the nanoparticle sizes. We highlight that this dependence is not always evident and some intriguing properties, such as the heterogeneous dynamics of both polymer chains and nanoparticles and their nonGaussian behavior at short and long timescales, are not particularly influenced by the degree of attraction between nanoparticles and polymer for the range of interactions we study (up to 6 kBT). We find the existence of weakly ordered interdigitated structures with sequential arrangements of particles and polymer chains, which separate each other and hence inhibit the formation of nanoparticle clusters. This is especially evident with big nanoparticles, being less prone to aggregate than small ones, even when their interaction with the polymer chain is as low as 0.5 kBT. Moreover, by integrating the stress-tensor autocorrelation functions, we estimate the shear viscosity and determine its dependence on the strength of the polymer−nanoparticle interactions and on the nanoparticle size. By acting as plasticizers, small nanoparticles decrease the viscosity, especially at low-to-moderate interactions with the polymer. By contrast, big nanoparticles that establish strongly attractive interactions with the polymer chains behave as thickening agents and significantly increase the viscosity. This complex and perhaps still scantily understood balance between the geometry of nanoparticles and their interaction with the polymer is key to predict and fully control the macroscopic response of nanocomposite materials and hence suitably tailor their mechanical properties. storage, electronics, optics, lubricants, and photovoltaics.9−14 This widespread involvement is justified by a number of notable properties, such as tensile strength,15 conductivity,16,17 viscosity,18 refractive index,19 flame retardation,20 etc. that can be achieved or enhanced upon uniform distribution of NPs in the polymer matrix.21 Due to the tendency of host particles to form clusters, achieving a uniform dispersion is, in many practical cases, a rather challenging task. The interaction established with the polymer chains, the NP size with respect to the polymer radius of gyration, and the average interparticle distance, are some among the key ingredients affecting the NP distribution and the performance of a PNC, with temperature, pH and NP concentration other fundamental parameters.22−24 PNCs have been part of our everyday life since the very beginning of the last century (e.g., with the synthesis of Bakelite, a clay-reinforced resin that found massive applications as household material).25 A few decades before Bakelite, successful car tires produced by blending carbon black, zinc oxide, and/or magnesium sulfate NPs in vulcanized rubber, had been assembled.26 These early discoveries, however, did not motivate the scientific community enough to further investigate NP−polymer mixtures until a couple of decades ago, when

I. INTRODUCTION Polymer nanocomposite (PNC)-based engineering represents nowadays one of the most promising driving forces in nanotechnology by proposing stimulating debates aiming to enhance both fundamental research and practical applications. PNCs are hybrid materials incorporating, in polymer melts or solutions, organic or inorganic particles with at least one dimension on the submicrometer scale. A more restrictive definition of PNCs would imply the genesis of new distinctive characteristics not observed in their macroscopic counterparts or constituent pure components. If this is the case, the hybrid nanocomposite can be regarded as a fundamentally new material, in which the presence of properly synthesized nanoparticles (NPs) induces intrinsically new mechanical and thermal properties.1−5 Generally speaking, the behavior of a polymer and its macroscopic properties are significantly altered in the proximity of NPs and consequently depend strongly on the NP surface area and degree of exfoliation (dispersion).6 If a good dispersion would optimize the amount of NPs needed to obtain a true PNC with well-defined properties,7,8 a bad dispersion, despite the nanoscale of the particles, would simply result in a conventional composite. The degree of NP dispersion in the polymer media as well as the NP diffusion is a topic of remarkable scientific and technological interest in materials science and, especially, in electrochemical energy © 2014 American Chemical Society

Received: December 19, 2013 Revised: February 11, 2014 Published: March 12, 2014 3731

dx.doi.org/10.1021/jp412440g | J. Phys. Chem. B 2014, 118, 3731−3742

The Journal of Physical Chemistry B

Article

embedded in a dense melt of 1000 unentangled Lennard−Jones polymer chains. More specifically, by performing coarse-grain molecular dynamic (MD) simulations, we aim at shading light on how NPs, their architecture, and energy interactions influence the structural and dynamical properties of PNCs.

researchers at Toyota significantly improved physical, mechanical, and thermal properties of nylon by adding very moderate amounts of mica.27 This achievement, along with the discovery that it is possible to intercalate polymer melts in clays without organic solvents,28 sparked off a broad involvement among researchers to boost the synthesis of novel intriguing NPs of precise shape and size. Simultaneously, developing new experimental equipment, such as scanning force, laser scanning fluorescence, and electron microscopes, triggered further research efforts into understanding the influence of these NPs on the macroscopic properties of the polymer. Nowadays, mixtures containing polymers and NPs have a crucial relevance in several technological applications, such as paints, rheological modifiers, and coatings. Regardless of the application for which they are thought, any significant improvement toward rigorously designed PNCs can be accomplished only through a complete understanding of the physics behind their formulation. In the last ten years, several studies, mainly focused on the rheological behavior and dynamics of NPs at the filler−polymer interface, substantially contributed to the understanding and control of the basic laws to properly modify and enhance several fundamental material properties.18,29,30 Tsagaropoulos and Eisenberg showed that at particularly strong NP−monomer interactions, polymer chains tend to stick permanently on the NP surfaces and the resulting PNC exhibits two glass-transition temperatures: one strictly due to the chains surrounding the NPs and the other, lower, to those far from them.31 However, at much more moderate interactions, where a permanent attachment of the chains onto the NP surfaces is not observed, this hybrid system exhibits a single glass transition and a rather long relaxation time.32,33 The relaxation time is especially influenced by the relative NP−monomer interaction strength, as molecular dynamics simulations performed by Starr and coworkers showed.34 These authors demonstrated that the relaxation dynamics in PNCs containing weakly interacting NPs are very similar to that of PNCs with attractive NP− monomer interactions, with the relaxation times becoming indistinguishable at increasing temperature. Smith and coworkers investigated the effect of the NP−polymer interface on the viscoelastic properties of PNCs.35 Their simulations pointed out the deep dependence of the viscoelastic behavior upon the nature of the NP−monomer interactions. As a general trend, micrometer-sized fillers always increase the viscosity of a pure polymer melt, as predicted by the Einstein−Batchelor law: η = ηp(6.2x2f + 2.5xf + 1), with ηp as the viscosity of the pure polymer and xf as the volume fraction of the filler. By contrast, when nanosized particles are used, the rheological behavior of PNCs is less obvious, as at least NP size, shape, and polymerNP interaction contribute to modify its viscosity. Very recently, Kumar and co-workers showed that the shear viscosity of a polymer melt can be substantially reduced when filled with small NPs establishing weak interactions with it, regardless of the polymer chain length.36 Jain and co-workers showed that adding a relatively small amount of silica NPs reduces the viscosity of polypropylene melt significantly, hence improving the processability of the polymer.37 In this work, we aim to understand how the interaction strength between spherical NPs of different size and polymer chains affects their degree of dispersion in the polymer melt, their long-time relaxation dynamics, and hence the macroscopic response of the nanocomposite material. To this end, we investigate the motion and diffusion process of 100 model NPs

II. MODEL In our coarse-grained model, each system contains Nnp = 100 spherical NPs in a melt of Nch = 1000 unentangled polymer bead-necklace chains consisting of l = 10 beads and whose radius of gyration is ⟨Rg⟩ ≃ 2.3σp. Therefore, the total number of chain beads (or monomers) is constant and set to Nb = 10000. All the polymer beads have the same diameter, σp, which is also our unit length, and are constrained to stay connected to each other by bond length and bond bending energy terms. By contrast, the NP diameter is a simulation parameter and can take three possible values: σnp = (1,2,3)σp. It should be noted that the mass of each NP is proportional to its volume, but its mass density remains constant. In other words, the mass of polymer beads and that of NPs of the same size are set to M, our unit mass. By contrast, the mass of bigger NPs with diameter σnp = 2σp or 3σp is equal to Mσ3np/σ3p = 8M or 27M, respectively. A typical polymer chain with NPs of different diameters is given in Figure 1.

Figure 1. Typical polymer chain containing l = 10 Lennard−Jones sites and bead diameter σp. NPs of different diameters are also shown.

All the nonbonded interactions are modeled via a truncated and shifted Lennard−Jones (LJ) potential that reads ⎧ ⎡⎛ ⎞12 ⎛ ⎞6 ⎤ σij ⎥ ⎪ ⎢⎜ σij ⎟ c c ⎪ 4εij⎢⎜ ⎟ − ⎜⎜ ⎟⎟ ⎥ − U0(rij) rij ≤ rij Uij = ⎨ ⎣⎝ rij ⎠ ⎝ rij ⎠ ⎦ ⎪ ⎪0 rij > rijc ⎩

(1)

where the subscripts i and j refer to chain beads or NPs, Uij is the LJ potential, σij = (σi + σj)/2 and εij are adjustable LJ parameters, and rcij = Crσij is the cutoff radius, with Cr = 2.5 for monomer−NP interactions and Cr = (2)1/6 for all the remaining pair interactions. Therefore, the interactions established in our systems are in general purely repulsive, but those between NPs and chain beads are attractive at large distances. Finally, the term U0(rcij) makes the LJ potential vanish at rij = rcij and beyond. Numerical details for all the systems studied are given in Table 1. The strength of interaction between NPs and chain monomers, here defined as εnp,m, is a simulation parameter assuming values between 0.5 and 6 kBT, where kB is the Boltzmann’s constant and T the absolute temperature. By 3732



Article

Initially, equilibration runs have been performed in the NPT ensemble, at constant Nb, T, and pressure Pσ3p/kBT = 1.0, to compress a completely random configuration up to the desired chain density ρb ≡ Nb/V = 0.7/σ3p, where V = L3 is the volume of the simulation box. At density ρb, which yields a box side L ≃ 24σp, we let the system relax in the NVT ensemble until the total energy achieved a steady value within the statistical fluctuations. In both ensembles, chains have been moved by reptation and configurational bias moves,38 whereas NPs have been moved by random displacements of maximum length L/2 in the three spatial directions. Configurational bias moves have been implemented to attempt the partial or complete regrowth of polymer chains. In accordance with the scheme proposed by Frenkel and Smit, first a trial configuration is generated and the corresponding Rosenbluth weight, Wtrial R , computed, and then the old configuration with its Rosenbluth weight, Wold R , is considered.38 The move is accepted with probability

Table 1. Details of the Systems Studied in this Paper, Each Consisting of Nch = 1000 Polymer Chains with Length l = 10 Beads and Nnp = 100a system

(εnp,m)/ (kBT)

(σnp)/ (σp)

(Dnpτ)/(σ2p) × 103

tR/τ

ηsτ−1σ3/ kBT

S1 S2 S3 S4 S5 S6 S7 S8 S9 S10 S11 S12 S13 S14 S15 S16 S17 S18 S19 S20 S21

0.5 0.5 0.5 1.0 1.0 1.0 2.0 2.0 2.0 3.0 3.0 3.0 4.0 4.0 4.0 5.0 5.0 5.0 6.0 6.0 6.0

1.0 2.0 3.0 1.0 2.0 3.0 1.0 2.0 3.0 1.0 2.0 3.0 1.0 2.0 3.0 1.0 2.0 3.0 1.0 2.0 3.0

161.0 28.6 10.1 149.0 27.5 10.0 132.0 25.7 9.3 116.0 24.8 9.0 103.0 23.9 8.9 91.7 23.6 8.6 82.8 22.5 8.5

6.34 24.30 40.57 8.02 24.41 52.08 9.70 30.20 69.33 10.35 33.83 65.59 12.57 43.75 96.52 13.39 43.02 101.67 15.80 43.32 92.69

0.93 1.25 4.64 0.96 1.11 4.59 0.87 1.48 5.65 0.88 1.30 4.77 0.99 1.58 6.24 0.90 1.17 5.90 0.85 1.56 7.42

Π = min[1, WRtrial /W Rold]

To compute a given number of trial orientations should be produced. In a lattice model, such a number cannot be larger than the lattice coordination number, but in an off-lattice model, it is infinite. Clearly, it is computationally opportune to calculate the Rosenbluth weight only for a subset of trial orientations. We selected k = 10, being a good compromise between the CPU time required by each trial move and the number of MC cycles required to equilibrate the energy of the system. In the NPT ensemble, volume changes have been attempted once every MC cycle, with one cycle equivalent to Nb + Nnp trial moves. Moves have been accepted according to the standard Metropolis algorithm. The maximum displacement of NPs and chain beads has been set to give an acceptance rate between 0.4 and 0.5, whereas the box volume was updated with acceptances between 0.3 and 0.4. For each system, an equilibrium configuration was then used as the starting configuration for the MD simulations with a time step δt = 3.96 · 10−4τ, where τ = σp(M/kBT)1/2 and M are the time and mass units, respectively. We introduced a Nosé− Hoover thermostat to set the velocity of the beads and hence the kinetic energy of the system, according to the imposed temperature. The thermostat was removed at t = 39.6τ, and with the distribution of velocities resulting from running MD simulations at constant temperature, we started MD simulations in the NVE ensemble, that is at constant energy, to collect the dynamical properties of the systems during t ≈ 200τ. In order to characterize the dynamics of the polymer chains and NPs, and the long-time structural relaxation of the systems, we estimated (i) the mean square displacement (MSD), (ii) the self-part of the van Hove autocorrelation function (VHF), (iii) the nonGaussian parameter (NGP), and (iv) the self-intermediate scattering function (s-ISF). As far as the polymer is concerned, each dynamical property was separately estimated from the motion of each single chain bead and also from the motion of the chains’ centers of mass, which gives information on the diffusion of the whole polymer chain. More specifically, the self-part of the VHF measures the probability distribution for the displacements of the beads at time t0 + t, given their positions at t0. It is defined as40

σnp are the NP diameters, whereas εnp,m is the well depth of the LJ potential. Additionaly, we include the NP diffusion coefficients, Dnpτ/ σ2p, the NP relaxation time, tR/τ, and the shear viscosity, η*s . a

tuning the value of εnp,m, we can investigate how (i) the dynamics of NPs and polymer chains, and (ii) the NP distributions in the polymer melt are influenced. By contrast, the remaining interaction terms, namely εm,m and εnp,np, are set to 1kBT for the whole set of systems studied. This choice provides NP−monomer interactions that are generally stronger than monomer−monomer interactions, as usually expected for many different NPs. In addition to the LJ interactions, a simple harmonic potential, UVIB, is used to model the bonded interactions describing the fluctuations of the bond length and angle: UVIB = K ϑ[cos(ϑ) − cos(ϑ0)]2 + K m(m − m0)2

(3)

Wtrial R ,

(2)

where m and m0 = σp are the bond length and its equilibrium value, respectively, whereas ϑ and ϑ0 = π are the bond angle and its equilibrium value, respectively. Kϑ and Km are force constants determining the most probable bond length and bending properties (stiffness) of the chains. In all the simulations, we set Kmσ2p/kBT = 110.7 and Kϑ/kBT = 10.0 rad−2.

III. SIMULATION METHODOLOGY To equilibrate the systems, Monte Carlo (MC) simulations have been performed at T = 350 K in cubic boxes of side L. Periodic boundary conditions have been applied, and the cell list (or linked list) method was implemented.38,39 The initial configurations consisted of chains and NPs sequentially placed as in a lattice network that were allowed to relax at a very high temperature by keeping the bond length constant to avoid chain breaking. This generated a completely random (isotropic) distribution of beads and provided the starting configuration for the MC simulations.

N

1 Gs(r, t ) = ⟨∑ δ[r − (ri(t0 + t ) − ri(t0))]⟩ N i=1 3733

(4)



Article

with N = Nb, Nch, or Nnp, depending on whether chain beads, chains’ centers of mass, or NPs, respectively, are considered; ⟨...⟩ is the ensemble average over N and initial time t0, and δ is the Dirac-delta. Note that Gs(r,t) would be a Gaussian distribution of r for freely diffusive particles. If the VHF is not Gaussian, then deviations from a normal Gaussian dynamical behavior should be expected. The nonGaussian parameter quantifies such deviations:41 5⟨Δr(t )4 ⟩ −1 3⟨Δr(t )2 ⟩2

α2,r(t ) =

(5)

where Δr(t) = r(t0 + t) − r(t0) denotes the displacement of the polymer chains (beads) or NPs in the time interval t starting at t0. Heterogeneous dynamics occurs on a timescale t, if the NGP is nonvanishing. Finally, the s-ISF contributes to quantify the structural relaxation and, in particular, the decay of the density fluctuations: N

1 ⟨∑ exp{iq·[rj(t + t0) − rj(t0)]}⟩ N j=1

Fs(q, t ) =

(6)

at the wave vector q = (qx, qy, qz), which corresponds to the first peak of the structure factor. The diffusion of the NPs and the shear viscosity of the system have been estimated by integrating the velocity and the stress tensor autocorrelation functions, respectively, through the following Green−Kubo formulas: D(t ) =

∫0

∞

⟨vi(t )vi(0)⟩dt

Figure 2. Radial distribution functions (RDFs) relative to systems with 1000 chains of 10 beads and NP diameter σnp = σp. The monomer− NP interaction parameter, εnp,m, varies between 0.5 and 6 kBT, as indicated in the panels. The RDF in the top left panel refers to a pure polymer melt. ○, □, and ◇ refer to monomer−monomer, NP− monomer, and NP−NP RDFs, respectively. Solid lines are guides for the eyes.

(7)

where vi(t) is the velocity of particle i at time t. The shear viscosity (ηs) has been estimated by integrating the time autocorrelation function of the off diagonal elements of the stress tensor Jxy:42,43 ∞ 1 ηs(t ) = ⟨J xy (0)J xy (t )⟩dt VkBT 0 (8)

∫

short distances and fast decay to 1, are practically undistinguishable (○ in Figure 2). By contrast, the gnp,m and gnp,np distributions (□ and ◇ in Figure 2, respectively) change substantially with εnp,m. In particular, at increasing interaction strengths, the NPs are on average much closer to the chain beads, resulting in a rather strong structural correlation at short distances as the height of the first peak of gnp,m indicates. This tendency is also observed from the analysis of the NP−NP RDFs, which display a gradual segregation between the NPs, whose mutual distance increases from r/σp ≃ 1 at εnp,m/kBT = 0.5 to r/σp ≃ 2.3 at εnp,m/kBT = 6, with the chain monomers separating them and creating a sort of interdigitated structure at short distances. This arrangement, which is promoted by the strong attraction between NPs and monomers, favors a more homogeneous distribution of NPs in the polymer melt and inhibits their tendency to aggregate in big clusters. In Figure 3, the RDFs of systems containing NPs with diameter σnp = 3σp are shown. Again, the monomer−monomer RDFs are identical in all systems, regardless of the interaction between chain beads and NPs. Interestingly enough, the first peak of gnp,np (green ◇ in Figure 3) is very small, indicating that big NPs are less prone to aggregate than small NPs, even when their interaction with the polymer is as low as 0.5 kBT. We note a clear tendency toward monomer−NP interdigitation with the NPs most likely located at one or two chain beads away from each other, as indicated by the double peak between r/σp = 4 and r/σp = 5. In other words, the polymer chains form a sort of single or double-layered shell surrounding the NPs,

To improve the statistics of the ensemble average in eq 8, we used all three independent off diagonal elements, Jxy, Jxz, and Jyz. For a two-component system, these tensors read44 2

J xy =

Nμ

∑ ∑ Mμ(vixviy)μ − μ=1 i=1

1 2

2

2

Nμ

Nν

∑ ∑ ∑ ∑ rμνx ij μ=1 ν=1 i=1 j>i

∂ϕ(rμνij) y ∂rμν ij

(9)

where i and j indicate the number of beads, Nμ and Nν are the number of beads in species μ and ν, respectively, (vxi )μ is the velocity along the x direction of a bead i belonging to species μ, and ∂ϕ(rμvij)/∂ryμvij is the force that particle i exerts on particle j along the y direction.

IV. RESULTS The dispersion of polymer chains and NPs in our systems can be investigated by computing the monomer−monomer, NP− monomer, and NP−NP radial distribution functions (RDFs), or more concisely, gm,m, gnp,m, and gnp,np, respectively. In Figure 2, we show several RDFs for systems containing NPs with diameter σnp = σp. The RDF of a pure polymer melt is also shown for comparison. From top to bottom, the interaction εnp,m between monomers and NPs increases from 0.5 to 6kBT. The strength of the NP−monomer attraction does not have any tangible effect on the distribution of the polymer chains, and consequently, all the gm,m functions, which show the typical behavior of a dense liquid with relatively strong correlations at 3734



Article

Figure 4. Snapshots of system S3 (top) and S21 (bottom). Each chain has a randomly selected color, whereas all NPs are in blue. In the right panels, only NPs are shown.

contrast, at high NP−monomer interactions, the polymer chains form an enthalpy-driven bridging structure linking the NPs in a sort of network. By decreasing the NP size, the solubility of NPs in the polymer melts increases and the miscibility region gets broader. The strength of the NP−monomer interaction influences not only the degree of exfoliation of NPs in the polymer melt but also their dynamical properties. In the main panel of Figure 5,

Figure 3. Radial distribution functions (RDFs) relative to systems with 1000 chains of 10 beads and NP diameter σnp = 3σp. The monomer− NP interaction parameter, εnp,m, varies between 0.5 and 6 kBT, as indicated in the panels. The RDF in the top left panel refers to a pure polymer melt. ○, □, and ◇ refer to monomer−monomer, NP− monomer, and NP−NP RDFs, respectively. Solid lines are guides for the eyes.

keeping them apart from each other. It should be noted that this interdigitated arrangement is basically independent from the interaction strength, as the close similarity between the gnp,np distributions from 3 to 6kBT demonstrate. To better appreciate the distribution of NPs in systems S3 (σnp/σp = 3, εnp,m/kBT = 0.5) and S21 (σnp/σp = 3, εnp,m/kBT = 6), two snapshots are given in Figure 4. The presence of interdigitated structures in polymer melts was also observed by Vacatello, with Monte Carlo simulations.45 This author studied NPs with diameter comparable with the radius of gyration of the polymer chains (similar scenario to that of the present work) and found the formation of shells made by polymer chains separating (dispersing) the particles. In general, systems of NPs and polymer melts show a complex variety of phase behaviors, very much depending on several properties, including at least temperature, NP sizes and concentrations, and interparticle interactions. Tuning these factors, Glotzer and co-workers showed that the same system can phase separate (a single big aggregate of NPs form) or display a very homogeneous dispersion in the polymer matrix.46 In their theoretical work, Hooper and Schweizer have shown that systems containing NPs with diameter between 10 and 32 times bigger than the polymer chain beads, display two distinct phase separations of clustered NPs, divided by a miscibility region with the NPs well dispersed in the polymer melt.47 At low NP−monomer interactions, they observe thin layers of polymer chains around the fillers as a result of an entropydriven phase separation mediated by depletion attraction. By

Figure 5. MSD of NPs with diameter σnp = σp at increasing strength of NP−monomer interaction, as indicated in the legend. In the inset, the MSDs of systems with the same NP−monomer interaction strength, εnp,m = 0.5kBT, but different NP sizes are shown. The dashed lines indicate the dependence of ⟨Δr2⟩ on time at short and long timescales.

we exemplarily show the MSDs of NPs with unit diameter and interactions in the range 0.5 ≤ εnp,m/kBT ≤ 6.0. NPs with diameter 2σp or 3σp produce MSDs with similar tendencies, as can be noted from the inset of Figure 5, where the MSD of NPs with diameter σp (black ○), 2σp (red □), and 3σp (green ◇) are shown for comparison at the same NP−monomer interaction strength. On very short timescales, below ∼0.5τ, 3735



Article

diffusive regime is very fast, with the NPs starting to interact with their surrounding neighbors, mostly chain monomers, which try to trap them in a sort of cage. Such a temporary confinement is rather weak, and its effect on the diffusion of NPs, which succeed to escape from it in less than a time decade, is actually rather feeble. A similar trend is also observed in the dynamics of the polymer chains and monomers. In Figure 7, we

the dynamics of NPs are Newtonian, practically independent from the strength of their interactions with the chain monomers, and ri(t0 + t) − ri(t0) = vi(t0 + t) + 6 (t2), where vi is the velocity of particle i at time t. In this time interval, the diffusion is ballistic (the MSD changes quadratically with time), and the NPs are free to move, being their mutual collisions or those with the chain monomers virtually absent. Soon after an almost unsignificant transition time, the NPs undergo a different dynamical behavior and enter the so-called diffusive regime, with the MSD changing linearly with time. At such longer timescales, the strength of NP−monomer interactions unveils its major influence on the NP dynamics by slowing down their diffusion as soon as the ballistic regime has been exceeded. This long time behavior can be better quantified by computing the NP diffusion coefficients, Dnp, displayed in Figure 6 as a function of the interaction strength for different

Figure 7. MSD of chain monomers (empty symbols) and polymer chains (full symbols) at different NP−monomer interactions and NP sizes. In the inset, the MSDs are plotted in a double linear scale. The dashed lines in the main panel indicate the dependence of the MSD on time at short and long timescales.

show the MSDs of the chain monomers (empty symbols) and the MSDs of the chains’ centers of mass (solid symbols), which are associated with the diffusion of the chain as a whole. Also in this case, the transition from the ballistic to the diffusive regime is rather fast, being the cage effect not particularly significant. The polymer chains display slower dynamics compared to that of their monomers, which, in principle, can move even though the center of mass of their chain does not. At long timescales, this difference still holds, as it is more easily detected from the double linear plot in the inset of Figure 7. Unlike the welldefined cause-effect in the NP dynamics, the dynamics of the polymer chains and monomers lack an unambiguous dependence on the NP size and NP−monomer interaction strength, at least in the range of interactions and diameters considered in this work. Although the MSDs seem to decrease with increasing NP size and NP−monomer interaction strength, this tendency is rather ambiguous and surely less appreciable than that observed in Figure 5, where the long-time diffusion slows down at higher interactions. Beside the remarkable fact that the polymer dynamics do not undergo a deep change with NP sizes and interaction, we believe that this effect should be much more evident at higher NP concentrations, which however would imply to simulate rather larger systems. Very recently, Meth et al. proposed an interesting theory to investigate the polymer center of mass diffusion into nanocomposites at low volume loading of spherical NPs.49 Their results are in good agreement with experiments,50,51 but in some cases overestimate the slowing down of the diffusion, possibly due, as the same authors suggest, to the nonGaussian dynamics of the polymer, not included in their model and that we highlight here.

Figure 6. Diffusion coefficient of NPs with diameter σp (black ●), 2σp (red ■), and 3σp (green ▲), as a function of the NP−monomer interaction strength. In the inset, dependence of the diffusion coefficient on the relative NP diameter in systems with different interaction strengths. In both panels, the solid lines are a guide for the eyes and error bars are smaller than the symbol size.

NP sizes. Regardless of the NP diameter, the diffusion coefficients decrease with increasing attraction between NPs and chain monomers. In the range of interactions studied here, the diffusion coefficients of small NPs (black ● in Figure 6) are clearly larger than those of NPs with a bigger size. However, the long time diffusion of big NPs (red ■ and green ▲ in Figure 6) are less influenced by a change in the NP−monomer interaction strength (see also Table 1). It is therefore evident that by properly tuning the NP size and NP−monomer interaction strength, we can modify the dynamics of NPs, determine their diffusion, and, hence, control their degree of distribution in the polymer melt. As recently pointed out by Liu and co-workers, the dependence of the diffusion coefficient on the interaction between NPs and chain monomers, not considered in the Stokes−Einstein equation, which underestimates the particle diffusion coefficient, can be included by introducing different reasonably defined effective hydrodynamic radii of the nanoparticle.48 As a general tendency, and regardless of the value assumed by εnp,m, the dynamics of NPs is characterized by the presence of two different regimes. The transition from the ballistic to the 3736



Article

Despite its relatively short time extension, the abovementioned cage effect has an interesting influence on the dynamical behavior of the whole system by determining small, but not negligible, deviations from the normal Gaussian behavior detected at short and long times, as can be appreciated in Figures 8 and 9, where the NGPs are displayed for polymer

Figure 9. NonGaussian parameters relative to the diffusion of NPs with (a) σnp/σp = 3, (b) σnp/σp = 2, and (c) σnp/σp = 1. The two vertical dashed lines at t = 0.1τ and t = 2τ emphasize the delayed nonGaussian deviations in the diffusion of NPs as a function of NP size. The legend of the symbols is given in Figure 7.

Figure 8. NonGaussian parameters relative to the diffusion of chain beads (top) and to the chains’ centers of mass (bottom) in all the systems studied. The two vertical dashed lines at t = 0.1τ and t = 10τ emphasize the delayed nonGaussian deviations in the diffusion of the chains’ centers of mass as compared to the diffusion of the chain beads. The legend of the symbols is given in Figure 7.

the chains as a whole (bottom panel of Figure 8), is characterized by higher and delayed nonGaussian deviations than those observed in the dynamics of their monomers. In other words, the polymer chains escape their cage, composed by other neighboring chains, and hence, enter the long-time diffusive regime rather slow compared to the time taken by their own monomers and display more pronounced deviations from a normal diffusion as the height of the peaks in the bottom panel of Figure 8 shows. Such nonGaussian deviations indicate the existence of heterogeneous dynamics determined by the simultaneous presence of fast and slow chains and monomers displacing distances of different length in the same lapse of time. Fast and slow objects play a crucial role in the structural relaxation of the systems. However, as generally observed in supercooled liquids,63 their influence on the dynamics is limited to a relatively short time interval, fading out at long timescales when the deviations from a normal behavior go exponentially to zero. The NGPs of Figure 9, relative to the dynamics of NPs, present a more defined dependence on NP size but are not particularly influenced by the NP−monomer interaction strength. In particular, the NGPs of NPs with diameters 3σp, 2σp, and σp, overlap eachother in panels (a), (b), and (c), respectively, by making it rather difficult to determine an evident tendency unambiguously correlating εnp,m to the extension and magnitude of these deviations. More specifically, for a given NP diameter, the interaction energies established with the polymer melt do not affect (i) the existence of fast and

chains and NPs, respectively. Deviations from the Gaussian behavior have been extensively investigated in liquid,52,53 glassy,54,55 and liquid crystalline56−62 systems. In close analogy to these systems, the NGPs of Figures 8 and 9 exhibit a timedependent behavior due to the transient cage effect of the neighbors surrounding the chain monomers and NPs. As a general trend, α2,r is essentially zero at very short timescales, where the beads are still very close to their original position and the presence of their nearest neighbors has not yet been perceived. In Figure 8, we display the evolution of α2,r associated to the diffusion of the chain beads (top panel) and of the chains’ centers of mass (bottom panel) in all systems studied. In both panels of Figure 8, we only detect a rather feeble dependence of the NGP on the NP−monomer interaction strength and on the NP size, as all the systems show a very similar nonGaussian behavior with maximum deviations at t ≃ 10τ (chain monomers, top panel) and t ≃ 20τ (chains’ centers of mass, bottom panel). The highest deviations are actually determined in the dynamics of systems containing big NPs strongly interacting with the chain monomers, but, as previously underlined for the MSDs in Figure 7, it is not a very marked effect, but only a rather fuzzy tendency. Furthermore, the dynamics of the chains’ centers of mass or, equivalently, of 3737



Article

Figure 10. Self-part of the van Hove autocorrelation function of chain monomers (left) and chains’ centers of mass (right) in systems containing NPs with diameter σnp = σp and εnp,m/kBT = 0.5, at t/τ ≃ 0.04 (black ○), 0.12 (red □), 0.28 (green ◇), and 0.36 (blue △). In the insets, the VHFs are computed at τ ≃ 7.88 (pink ○) and 118.76 (green ◇). In all panels, the thick dashed lines are Gaussian fits.

Figure 11. Self-part of the van Hove autocorrelation function of NPs with diameter σp (left), 2σp (middle), and 3σp (right), with εnp,m/kBT = 0.5 at t/ τ ≃ 0.04 (black ○), 0.12 (red □), 0.28 (green ◇), and 0.36 (blue △). In the insets, the VHFs are computed at t/τ ≃ 7.88 (pink ○) and 118.76 (green □). In all panels, the thick dashed lines are Gaussian fits.

slow particles (even at low εnp,m, the NGPs are not zero), (ii) the number distribution of their population (the peak of α2,r is always located at roughly the same height), (iii) the average lifetime of the cage effect (the extension of the deviations is practically the same). Nevertheless, at σnp/σp = 2 (panel b), the magnitude of the deviations decreases, as compared to those given in panel (c), suggesting the appearance of a less heterogeneous scenario where the existence of extremely fast (or slow) NPs is still possible but not probable. In other words, the dynamics of big NPs are on average more uniform than that of small NPs, and, as confirmed by the inset of Figure 5, it is

also slower. This scenario is further confirmed by the NGPs of NPs with diameter 3σp, given in panel (a) of Figure 9, where it is quite arduous to distinguish physically meaningful nonGaussian deviations from a simple statistical noise. In this last case, the dynamics of the systems are rather homogeneous, making a distinction among fast and slow NPs almost completely worthless. Deviations from Gaussian dynamics in our systems can also be investigated by computing the VHFs, defined in eq 4 and given in Figure 10 for monomers and chains and in Figure 11 for NPs of different sizes. In both figures, the NP−monomer 3738



Article

interaction strength is set to εnp,m/kBT = 0.5. At higher interactions, these autocorrelation functions display a very similar tendency, analogous to what was already noticed in Figures 7 and 8 and are not shown here. The VHFs shed light on the distribution of displacements of monomers, polymer chains, and NPs, after different intervals of time. In particular, at short timescales, chains and NPs are still very close to their original positions, having displaced only a small fraction of their unit length, σp. As time passes, the distribution of displacements acquires a broader shape, indicating a wide spectrum of dynamical behaviors, where the fastest objects are determined by the functions’ tail. At very long timescales, the peak of each distribution moves toward longer displacements and the probability of finding polymer chains, monomer, or NPs, still in their original locations goes gradually to zero. For comparison, we also plot the Gaussian approximation which, in general, underestimates the presence of fast-moving objects as well as their contribution to the structural relaxation of the whole system, especially at those time intervals in which the impact of the cage effect is particularly pronounced (i.e., the transition from the ballistic to the diffusive regime). Therefore, the analysis of the VHFs confirms the above conclusions on the dynamics in polymer melts: (i) chain monomers displace faster than their chains, especially in the ballistic regime, being this difference less evident at long times and (ii) the size of NPs do influence their diffusion, as small NPs displace significantly faster than big NPs. The lack of a quantitative agreement with the Gaussian fits is also observed for the pure polymer melt (not shown here) and can be taken as an estimation of the deviations from the Rouse model. The Rouse model is generally applied to predict the long-time diffusion of unentangled polymer chains whose beads interact with each other only by bonded interactions. In particular, the chains are modeled as Gaussian random coils of monomers following a Brownian dynamics, with the diffusion coefficient (D) given by the Einstein relation D=

kBT Nζ

Figure 12. Time dependence of the s-ISF of NPs (main panels) with diameters (a) σp, (b) 2σp, and (c) 3σp and interaction strengths set to εnp,m/kBT = 0.5 (○), 1 (□), 2 (◇), 3 (△), 4 (◁), 5 (×), and 6 (+), at σpq* ≃ 7. In the insets, the s-ISFs relative to the relaxation dynamics of chain beads are given. Lines are exponential (main panel) and stretched exponential (inset) fits.

of Fs(t) is well-fitted by an exponential function of the form exp[−(t/tR)β], where tR is the time at which Fs(t) decays to e−1, usually referred to as relaxation time. The decay of the s-ISF relative to the NP dynamics is practically exponential, with 0.9 < β < 1 for NPs with diameter σp or 2σp, in panels (a and b), respectively, and β ≃ 0.8 for the biggest NPs in panel (c). The NP relaxation times vary between tR ≃ 6τ and 100τ, depending on their size and interaction with the chain monomers. In Figure 13, we can appreciate the unequivocal interplay

(1)

where N and ζ are the number of beads per chain and the bead friction coefficient, respectively.64 The chains modeled in the present work are made by 10 beads, that is short enough to be considered unentangled, but nonbonded interactions are taken into account and nonGaussian deviations, not predicted by the Rouse model, observed. In the diffusive regime, these deviations decay to zero and the MSD of the polymer chains in Figure 7 changes linearly with time. Therefore, at long timescales, despite the presence of nonbonded interactions, the Gaussian assumption of the Rouse model can be accepted. However, at short to moderate timescales, when the cage effect delays the motion of our polymer chains, the Rouse picture would deviate from the simulation results and hence show a simple diffusive behavior. Deviations from the Rouse scenario have been reported in several works.65−69 In addition, we study the structural relaxation of NPs and chain monomers by computing the self-intermediate scattering functions, Fs(q,t). Since most of the relevant structural information is contained at the first maximum of the static structure factor, S(q*), with q*σp ≃ 7, we evaluate the s-ISF at the wave vector q*, so that Fs(t) = Fs(q*,t). In Figure 12, we show the s-ISF relative to the relaxation dynamics of chain monomers (inset) and NPs (main panel). The long-time decay

Figure 13. Relaxation times, tR/τ, of NPs as a function of the NP− monomer interaction strength. ○, □, and △ refer to NP diameters σp, 2σp, and 3σp, respectively. The solid lines are power-law fits of the type tR/τ ∝ εωnp,m, with ω ≃ 0.3.

connecting the NP−monomer interaction strength and NP sizes with the relaxation time of the polymer nanocomposite. It can be observed that at given NP size, the relaxation time increases with εnp,m following a power law of the type tR/τ ∝ εωnp,m, with ω ≃ 0.3. By contrast, the s-ISFs of the chain beads (insets in Figure 12) do not vary significantly from one system to another and are better described by stretched exponential fits with β ≃ 0.7 and 20 < tR/τ < 30. 3739



Article

more pronounced in PNCs containing bigger particles (panel c), where the shear viscosity can be as big as 8η*s,0 at the highest interaction strength. In other words, while small NPs act as plasticizers, big NPs act as thickening agents, by lowering and increasing, respectively, the viscosity of the polymer melt. These results are in agreement with several experiments,70−72 recently systematically investigated by Kumar and co-workers,36 where the viscosity can either increase or decrease according to the NP size and the degree of entaglement of the polymer chains. The simulation results are fitted by an exponential function resembling the usual response to a stepwise change from a steady-state condition to another:

To address the effect of NP size and NP−polymer interaction on the macroscopic behavior of our systems, we have estimated the shear viscosity of the polymer melt, ηs, by integrating the stress-tensor time autocorrelation function (not shown here). We found similar tendencies previously pointed out by nonequilibrium MD simulations in ref 36 with the NP size and NP−monomer interaction playing a key role in determining the viscosity of the system. Due to the nonnegligible fluctuations in the decay of the stress-tensor autocorrelation functions, computing the shear viscosity through equilibrium MD simulations produces relatively noisy statistics. In order to reduce these random fluctuations, the stress-tensor autocorrelation functions have been smoothed out by an exponential fit of the type c1exp(−tc2/c3), with ci constants. This yields a clearer picture of their long time decay by leaving out noise while retaining all the significant information. The application of eq 8 to these exponential fits provides the results given in Figure 14, where the adimensional

ηs(t ) − η0 η∞ − η0

= 1 − exp( −t /ts*)

(10)

∞

where η and η are the values of the shear viscosity at t = 0 and t = ∞, respectively, and ts* is a time constant in the order of 10−1 τ. Such a time dependence of the shear viscosity is simply due to the time decay of the stress tensor autocorrelation function, which is also exponential. 0

V. CONCLUSIONS In summary, we have investigated the effect of NP size and NP−monomer interaction strength on the structural and dynamical properties of dense polymer melts as well as on their rheological behavior. By increasing the NP−monomer interaction, we observe an interesting redistribution of NPs, which get closer and closer to the neighboring chain beads and form with them a sort of local interdigitated arrangement preventing their aggregation. The tendency to form interdigitated structures is more tangible in systems containing bigger NPs, where the polymer chains are arranged in a layered shell-like fashion surrounding the NPs and hampering their clustering. The NP diffusion through the dense polymer matrix is slightly nonGaussian, with fast and slow particles escaping their cage on different timescales and hence generating a heterogeneous scenario, which is particularly evident when NPs as small as a single chain bead are present. By contrast, big NPs display a rather homogeneous and significantly slower dynamics, with diffusion coefficients up to 2 orders of magnitude smaller than those of small NPs. The existence of fast and slow NPs is further verified by the analysis of the probability distribution of their displacements (self-part of VHFs). By assuming a broader and broader shape as time passes, these distribution functions highlight a wide spectrum of dynamical behaviors, with their tail confirming the presence of fast particles. The dynamics of the polymer chains is also nonGaussian, with the monomers moving and relaxing faster than the polymer chains they belong to. The short-time dynamics is Newtonian and basically independent of the interaction established with the NPs. The transition from the Newtonian to the diffusive regime is rather fast, and it does not seem to be particularly influenced by the NP size or by the NP−monomer interaction strength, even if the long-time MSD appears to decrease with increasing both parameters. By integrating the stress-tensor autocorrelation function, we have also investigated the dependence of the shear viscosity on the NP size and NP−monomer interaction strength. Although both parameters modify the viscosity of the polymer melt and determine its macroscopic response, the size is clearly the most significant factor. More specifically, small particles can reduce the viscosity of the polymer melt and hence act as plasticizers,

Figure 14. Shear viscosity in systems containing NPs with diameter σnp = (a) σp, (b) 2σp, and (c) 3σp and interaction strengths set to εnp,m/kBT = 0.5 (□), 1 (◇), 2 (Δ), 3 (+), 4 (×), 5 (▷), and 6 (◁). The shear viscosity of the pure polymer melt is given in the three panels for reference (●). The solid lines are exponential fits obtained from eq 10.

shear viscosity, ηs* = ηsτ−1σ3/kBT, of PNCs containing NPs with diameters (panel a) σp, (panel b) 2σp, and (panel c) 3σp are shown. In the limits imposed by this computation technique, interesting quantitative and qualitative tendencies can be unambiguously detected and are worth mentioning here. The shear viscosity of the pure polymer melt (● in the three panels of Figure 14) is roughly η*s,0 = 1.0. The addition of spherical NPs alter this value as a result of the combined action of their size and interaction with the polymer chains. Our results indicate that the particle size has a more tangible effect on the macroscopic response of the systems, whereas the effect of the NP−monomer interaction strength is less evident. More specifically, NPs with unit diameters do not significantly alter the viscosity of the polymer melt, although a moderate reduction of up to 15% is clearly observed in panel (a). By contrast, adding NPs with σnp/σp = 2 increases the viscosity of the polymer up to 1.6ηs,0 * at εnp,m = 6kBT. This tendency is even 3740



Article

Retardant Polymer Nanocomposites; Morgan, A. B., Wilkie, C. A., Eds.; John Wiley & Sons, Inc.: New York, 2007; pp 31−66. (7) Lu, C.; Mai, Y.-W. Influence of Aspect Ratio on Barrier Properties of Polymer-Clay Nanocomposites. Phys. Rev. Lett. 2005, 95, 088303. (8) Mitchell, C. A.; Bahr, J. L.; Arepalli, S.; Tour, J. M.; Krishnamoorti, R. Dispersion of Functionalized Carbon Nanotubes in Polystyrene. Macromolecules 2002, 35, 8825−8830. (9) Croce, F.; Curini, R.; Martinelli, A.; Persi, L.; Ronci, F.; Scrosati, B.; Caminiti, R. Physical and Chemical Properties of Nanocomposite Polymer Electrolytes. J. Phys. Chem. B 1999, 103, 10632−10638. (10) Bockstaller, M. R.; Thomas, E. L. Optical Properties of PolymerBased Photonic Nanocomposite Materials. J. Phys. Chem. B 2003, 107, 10017−10024. (11) Balazs, A. C.; Emrick, T.; Russell, T. P. Nanoparticle Polymer Composites: Where Two Small Worlds Meet. Science 2006, 314, 1107−1110. (12) Akcora, P.; Liu, H.; Kumar, S. K.; Moll, J.; Li, Y.; Benicewicz, B. C.; Schadler, L. S.; Acehan, D.; Panagiotopoulos, A. Z.; Pryamitsyn, V.; et al. Anisotropic Self-Assembly of Spherical Polymer-grafted Nanoparticles. Nat. Mater. 2009, 8, 354−359. (13) Srivastava, S.; Agarwal, P.; Archer, L. A. Tethered NanoparticlePolymer Composites: Phase Stability and Curvature. Langmuir 2012, 28, 6276−6281. (14) Mackay, M. E.; Tuteja, A.; Duxbury, P. M.; Hawker, C. J.; Van Horn, B.; Guan, Z.; Chen, G.; Krishnan, R. S. General Strategies for Nanoparticle Dispersion. Science 2006, 311, 1740. (15) Coleman, J. N.; Khan, U.; Gun’ko, Y. K. Mechanical Reinforcement of Polymers Using Carbon Nanotubes. Adv. Mater. 2006, 18, 689−706. (16) White, S. I.; Vora, P. M.; Kikkawa, J. M.; Winey, K. I. Resistive Switching in Bulk Silver NanowirePolystyrene Composites. Adv. Funct. Mater. 2011, 21, 233−240. (17) Feng, Y.; Zou, H.; Tian, M.; Zhang, L.; Mi, J. Relationship between Dispersion and Conductivity of Polymer Nanocomposites: A Molecular Dynamics Study. J. Phys. Chem. B 2012, 116, 13081−13088. (18) Oh, H.; Green, P. F. Polymer Chain Dynamics and Glass Transition in Athermal Polymer/Nanoparticle Mixtures. Nat. Mater. 2009, 8, 139−143. (19) Nussbaumer, R. J.; Caseri, W. R.; Smith, P.; Tervoort, T. Polymer-TiO2 Nanocomposites: A Route Towards Visually Transparent Broadband UV Filters and High Refractive Index Materials. Macromol. Mater. Eng. 2003, 288, 44−49. (20) Kashiwagi, T.; Du, F.; Douglas, J. F.; Winey, K. I.; Harris, R. H., Jr.; Shields, J. R. Nanoparticle Networks Reduce the Flammability of Polymer Nanocomposites. Nat. Mater. 2005, 4, 928−933. (21) Paul, D. R.; Robeson, L. M. Polymer Nanotechnology: Nanocomposites. Polymer 2008, 49, 3187−3204. (22) Tung, W.-S.; Clarke, N.; Composto, R. J.; Winey, K. I. Temperature Dependence of Polymer Diffusion in MWCNT/PS Nanocomposites. Macromolecules 2013, 46, 2317−2322. (23) Wang, H.; Yang, H.; Liu, H.; Yu, Y.; Xin, H. A Mesoporous Silica Nanocomposite Shuttle: pH-Triggered Phase Transfer between Oil and Water. Langmuir 2013, 29, 6687−6696. (24) Basu, T. S.; Ghosh, S.; Gierlotka, S.; Ray, M. Collective Charge Transport in Semiconductor-metal Hybrid Nanocomposite. Appl. Phys. Lett. 2013, 102, 053107. (25) Baekeland, L. H. Bakelite, a New Composition of Matter: Its Syntehsis, Constitution, and Uses. Sci. Am. 1909, 68 (Suppl), 322. (26) Goodyear, C. Verfahren zur Fabrication von Gegenständen welche mehr oder weniger hart, biegsam oder elastisch sind, aus Kautschuk und Gutta-percha in Verbindung mit anderen Stoffen. Dingler’s Polytechnisches Journal 1856, 139, 376. (27) Usuki, A.; Kojima, Y.; Kawasumi, M.; Okada, A.; Fukushima, Y.; Kurauchi, T. Synthesis of Nylon 6-Clay Hybrid. J. Mater. Res. 1993, 8, 1179. (28) Vaia, R. A.; Ishii, H.; Giannelis, E. P. Synthesis and Properties of Two-Dimensional Nanostructures by Direct Intercalation of Polymer Melts in Layered Silicates. Chem. Mater. 1993, 5, 1694−1696.

whereas big particles, by acting as thickening agents, significantly increase it. With dependence on the strength of NP−monomer interactions, both the plasticizer and the thickening effects can be tuned to achieve the desired value of viscosity for applications spanning from paint and printing industries (thickeners) to the huge PVC market (plasticizers). Although a comparison with a specific experimental system would definitely be important, the model used here is rather simple to reproduce a particular PNC, which generally includes significantly longer polymer chains. A more sophisticated model would better reproduce a real system, but most likely could not map its dynamics along 5 or 6 time decades. Beside the fact that this would lead to more demanding computational effort for simulating properly scaled systems, an increase of the chain length should maintain the overall tendencies observed in the present work. As far as NP−monomer interactions are concerned, we have explored a rather broad range of attractions, which are responsible for a change in the NP distribution and/or NP dynamics. Unless very strong and localized NP−polymer attractions are investigated, which could deeply change the aggregation behavior and lead to phase separation, these results should be of generic validity.

■

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. Present Address

† School of Chemical Engineering and Analytical Science, The University of Manchester, Sackville Street, Manchester M13 9PL, U.K.

Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS The author acknowledges Grants JCI-2010-06943 and CTQ2011-29336-C03 offered by the Spanish Ministry of Science and Innovation. Allan Mackie and Josep Bonet from University Rovira i Virgili are gratefully acknowledged for having provided computational resources. CIBER-BBN is an initiative funded by the VI National R&D&i Plan 2008-2011, Iniciativa Ingenio 2010, Consolider Program, CIBER Actions and financed by the Instituto de Salud Carlos III with assistance from the European Regional Development Fund.

■

REFERENCES

(1) Tyan, H.-L.; Liu, Y.-C.; Wei, K.-H. Thermally and Mechanically Enhanced Clay/Polyimide Nanocomposite via Reactive Organoclay. Chem. Mater. 1999, 11, 1942−1947. (2) Ash, B. J.; Siegel, R. W.; Schadler, L. S. Mechanical Behavior of Alumina/Poly(methyl methacrylate) Nanocomposites. Macromolecules 2004, 37, 1358−1369. (3) Ramanathan, T.; Liu, H.; Brinson, L. C. Functionalized SWNT/ Polymer Nanocomposites for Dramatic Property Improvement. J. Polym. Sci., Part B: Polym. Phys. 2005, 43, 2269−2279. (4) Ramanathan, T.; Abdala, A. A.; Stankovich, S.; Dikin, D. A.; Herrera-Alonso, M.; Piner, R. D.; Adamson, D. H.; Schniepp, H. C.; Chen, X.; Ruoff, R. S.; et al. Functionalized Graphene Sheets for Polymer Nanocomposites. Nat. Nanotechnol. 2008, 3, 327−331. (5) Sánchez-Jiménez, P. E.; Pérez-Maqueda, L. A.; Perejón, A.; Criado, J. M. Nanoclay Nucleation Effect in the Thermal Stabilization of a Polymer Nanocomposite: A Kinetic Mechanism Change. J. Phys. Chem. C 2012, 116, 11797−11807. (6) Manias, E.; Polizos, G.; Nakajima, H.; Heidecker, M. J. Fundamentals of Polymer Nanocomposites Technology. In Flame 3741



Article

(29) Modeling of Nanocomposites; Picu, C.; Keblinski, P., Eds.; WileyVCH: Weinheim, Germany, 2003. (30) Allegra, G.; Raos, G.; Vacatello, M. Theories and Simulations of Polymer-based Nanocomposites: From Chain Statistics to Reinforcement. Prog. Polym. Sci. 2008, 33, 683−731. (31) Tsagaropoulos, G.; Eisenburg, A. Direct Observation of Two Glass Transitions in Silica-filled Polymers. Implications to the Morphology of Random Ionomers. Macromolecules 1995, 28, 396− 398. (32) Kropka, J. M.; Putz, K. W.; Pryamitsyn, V.; Ganesan, V.; Green, P. F. Origin of Dynamical Properties in PMMA-C60 Nanocomposites. Macromolecules 2007, 40, 5424−5432. (33) Kropka, J. M.; Garcia Sakai, V.; Green, P. F. Local Polymer Dynamics in Polymer-C60 Mixtures. Nano Lett. 2008, 8, 1061−1065. (34) Starr, F. W.; Schrøder, T. B.; Glotzer, S. C. Molecular Dynamics Simulation of a Polymer Melt with a Nanoscopic Particle. Macromolecules 2002, 35, 4481−4492. (35) Smith, G. D.; Bedrov, D.; Li, L.; Byutner, O. A Molecular Dynamics Simulation Study of the Viscoelastic Properties of Polymer Nanocomposites. J. Chem. Phys. 2002, 117, 9478−9490. (36) Kalathi, J. T.; Grest, G. S.; Kumar, S. K. Universal Viscosity Behavior of Polymer Nanocomposites. Phys. Rev. Lett. 2012, 109, 198301. (37) Jain, S.; Goossens, J. G. P.; Peters, G. W. M.; van Duin, M.; Lemstra, P. J. Strong Decrease in Viscosity of Nanoparticle-filled Polymer Melts through Selective Adsorption. Soft Matter 2008, 4, 1848. (38) Frenkel, D.; Smit, B. Understanding Molecular Simulation: From Algorithms to Applications; Academic Press: London, 2nd ed., 2002. (39) Hockney, R. W.; Eastwood, J. W. Computer Simulations Using Particles; McGraw-Hill: New York, 1981. (40) Hansen, J. P.; McDonald, I. R. Theory of Simple Liquids; Academic: London, 1986. (41) Rahman, A. Correlations in the Motion of Atoms in Liquid Argon. Phys. Rev. 1964, 136, A405. (42) Gubbins, K. E. In Statistical Mechanics; Singer, K., Ed.; The Chemical Society Burlington House: London, 1972; Vol 1, pp 194− 253. (43) Steele, W. A. In Transport Phenomena in Fluids; Hanley, H. J. M., Ed.; Marcel Dekker: New York, 1969; pp 209−312. (44) Fernández, G. A.; Vrabec, J.; Hasse, H. Fluid Phase Equilib. 2004, 221, 157−163. (45) Vacatello, M. Monte Carlo Simulations of Polymer Melts Filled with Solid Nanoparticles. Macromolecules 2001, 34, 1946. (46) Starr, F. W.; Douglas, J. F.; Glotzer, S. C. Origin of Particle Clustering in a Simulated Polymer Nanocomposite and Its Impact on Rheology. J. Chem. Phys. 2003, 119, 1777. (47) Hooper, J. B.; Schweizer, K. S. Theory of Phase Separation in Polymer Nanocomposites. Macromolecules 2006, 39, 5133. (48) Liu, J.; Cao, D.; Zhang, L. Molecular Dynamics Study on Nanoparticle Diffusion in Polymer Melts: A Test of the StokesEinstein Law. J. Phys. Chem. C 2008, 112, 6653−6661. (49) Meth, J. S.; Gam, S.; Choi, J.; Lin, C.-C.; Composto, R. J.; Winey, K. I. Excluded Volume Model for the Reduction of Polymer Diffusion into Nanocomposites. J. Phys. Chem. B 2013, 117, 15675− 15683. (50) Gam, S.; Meth, J. S.; Zane, S. G.; Chi, C.; Wood, B. A.; Seitz, M. E.; Winey, K. I.; Clarke, N.; Composto, R. J. Macromolecular Diffusion in a Crowded Polymer Nanocomposite. Macromolecules 2011, 44, 3494−3501. (51) Lin, C.-C.; Gam, S.; Meth, J. S.; Clarke, N.; Winey, K. I.; Composto, R. J. Do Attractive PolymerNanoparticle Interactions Retard Polymer Diffusion in Nanocomposites? Macromolecules 2013, 46, 4502−4509. (52) Hurley, M. M.; Harrowell, P. NonGaussian Behavior and the Dynamical Complexity of Particle Motion in a Dense Twodimensional Liquid. J. Chem. Phys. 1996, 105, 10521.

(53) Kob, W.; Donati, C.; Plimpton, S. J.; Poole, P. H.; Glotzer, S. C. Dynamical Heterogeneities in a Supercooled Lennard−Jones Liquid. Phys. Rev. Lett. 1997, 79, 2827. (54) Chaudhuri, P.; Berthier, L.; Kob, W. Universal Nature of Particle Displacements close to Glass and Jamming Transitions. Phys. Rev. Lett. 2007, 99, 060604. (55) Vorselaars, B.; Lyulin, A. V.; Karatasos, K.; Michels, M. A. J. NonGaussian Nature of Glassy Dynamics by Cage to Cage Motion. Phys. Rev. E 2007, 75, 011504. (56) Bier, M.; van Roij, R.; Dijkstra, M.; van der Schoot, P. SelfDiffusion of Particles in Complex Fluids: Temporary Cages and Permanent Barriers. Phys. Rev. Lett. 2008, 101, 215901. (57) Grelet, E.; Lettinga, M. P.; Bier, M.; van Roij, R.; van der Schoot, P. Dynamical and Structural Insights into the Smectic Phase of Rod-like Particles. J. Phys.: Condens. Matter 2008, 20, 494213. (58) Patti, A.; El Masri, D.; van Roij, R.; Dijkstra, M. Stringlike Clusters and Cooperative Interlayer Permeation in Smectic Liquid Crystals Formed by Colloidal Rods. Phys. Rev. Lett. 2009, 103, 248304. (59) Matena, R.; Dijkstra, M.; Patti, A. NonGaussian Dynamics in Smectic Liquid Crystals of Parallel Hard Rods. Phys. Rev. E 2010, 81, 021704. (60) Belli, S.; Patti, A.; van Roij, R.; Dijkstra, M. Heterogeneous Dynamics in Columnar Liquid Crystals of Parallel Hard Rods. J. Chem. Phys. 2010, 133, 154514. (61) Patti, A.; El Masri, D.; van Roij, R.; Dijkstra, M. Collective Diffusion of Colloidal Hard Rods in Smectic Liquid Crystals: Effect of Particle Anisotropy. J. Chem. Phys. 2010, 132, 224907. (62) Patti, A.; Belli, S.; van Roij, R.; Dijkstra, M. Relaxation Dynamics in the Columnar Liquid Crystal Phase of Hard Platelets. Soft Matter 2011, 7, 3533−3545. (63) Darst, R. K.; Reichman, D. R.; Biroli, G. Dynamical Heterogeneity in Lattice Glass Models. J. Chem. Phys. 2010, 132, 044510. (64) Rubinstein, M.; Colby, R. H. Polymer Physics; Oxford University Press: Oxford, 2003. (65) Kremer, K.; Grest, G. S. Dynamics of Entangled Linear Polymer Melts: A Molecular Dynamics Simulation. J. Chem. Phys. 1990, 92, 5057. (66) Smith, G. D.; Paul, W.; Monkenbusch, M.; Richter, D. A Comparison of Neutron Scattering Studies and Computer Simulations of Polymer Melts. Chem. Phys. 2000, 261, 61−74. (67) Paul, W. Anomalous Diffusion in Polymer Melts. Chem. Phys. 2002, 284, 59. (68) Paul, W.; Smith, G. D. Structure and Dynamics of Amorphous Polymers: Computer Simulations Compared to Experiment and Theory. Rep. Prog. Phys. 2004, 67, 1117−1185. (69) Farago, J.; Meyer, H.; Semenov, A. N. Anomalous Diffusion of a Polymer Chain in an Unentangled Melt. Phys. Rev. Lett. 2011, 107, 178301. (70) Grabowski, C. A.; Adhikary, B.; Mukhopadhyay, A. Dynamics of Gold Nanoparticles in a Polymer Melt. Appl. Phys. Lett. 2009, 94, 021903. (71) Schmidt, R. G.; Gordon, G. V.; Dreiss, C. A.; Cosgrove, T.; Krukonis, V. J.; Williams, K.; Wetmore, P. M. A Critical Size Ratio for Viscosity Reduction in Poly(dimethylsiloxane) Polysilicate Nanocomposites. Macromolecules 2010, 43, 10143−10151. (72) Ndong, R. S.; Russel, W. B. Rheology of Surface-Modified Titania Nanoparticles Dispersed in PDMS Melts: The Significance of the Power Law. J. Rheol. 2012, 56, 27−43.

3742


Molecular Dynamics of Spherical Nanoparticles in Dense Polymer Melts

Recommend Documents