Using Theory and Simulations To Calculate Effective Interactions in

Dec 15, 2016 - ... between sites i and j, and the angle brackets represent ensemble averaging over uncorrelated snapshots in the simulation trajectory...
2 downloads 0 Views 3MB Size
Article pubs.acs.org/Macromolecules

Using Theory and Simulations To Calculate Effective Interactions in Polymer Nanocomposites with Polymer-Grafted Nanoparticles Tyler B. Martin† and Arthi Jayaraman*,†,‡ †

Department of Chemical and Biomolecular Engineering, Colburn Laboratory, and ‡Department of Materials Science and Engineering, University of Delaware, Newark, Delaware 19716, United States S Supporting Information *

ABSTRACT: Using theory and large-scale simulations, we demonstrate how one can program structure and thermodynamics into polymer-grafted particles filled polymer nanocomposites (PNCs). We simulate varying graft (G) and matrix (M) polymer compositions for varying model graft−matrix bead pairwise interactions, χGM, and calculate structural features and the effective graft−matrix interaction parameter, χeff GM, in the PNC. Varying the graft (G) and matrix (M) polymer compositions provides tunability of morphology (particle dispersion/aggregation) and graft−matrix interpenetration at each χGM. Thermodynamically, for all composites the χeff GM exhibits negative values (effective attraction) at low values of χGM, with a sharp transition to positive values (effective repulsion) at large values of χGM. The sharp transition in χeff GM coincides with the structurally characterized particle dispersion−aggregation transition marked by the onset of upturn in the matrix−matrix structure factor at zero wavenumber. Strikingly, regardless of the composition of the graft and matrix chains or the dispersion−aggregation transition point, universally, the effective interactions in the PNC at the dispersion−aggregation transition is identical to the analogous athermal PNC.

I. INTRODUCTION Polymer nanocomposites are ubiquitous materials that are composed of a matrix polymer and a typically minority filler material that modifies and improves the macroscopic properties of the matrix polymer. While the chemistry of the matrix polymer and filler material has a large impact on the properties of the nanocomposite, equally important is the composite morphology, i.e., the order and arrangement of the matrix chains and the nanoscale fillers.1−4 Driven by the need to control and tune composite morphology, much work has gone into designing and understanding composites where the nanoscale fillers are grafted with polymers of the same chemistry as the matrix polymer.5−8 These studies have established how varying grafting density, graft to matrix chain length ratio, filler fraction, filler particle chemistry, filler particle shape and size, graft and matrix polymer dispersity, and flexibility, etc., lead to a range of composite morphologies dispersed or aggregated states with isotropic and anisotropic assemblies of the filler particles in the matrix.9−19 In most of these studies, the particle dispersion to aggregation transition is linked to the mixing/demixing of the graft and matrix chains, i.e., wetting/dewetting of the grafted layer by matrix chains. The established understanding is that the graft−matrix mixing (i.e., grafted layer wetting) is synonymous with particle dispersion while graft−matrix demixing (i.e., grafted-layer dewetting) is synonymous with particle aggregation. Recently, using both simulations and experiments, we challenged this established understanding for the case of polymer nano© XXXX American Chemical Society

composites with chemically dissimilar graft and matrix chains9,10,20−27 by demonstrating that the wetting−dewetting and dispersion−aggregation transitions are in fact distinct.28 We showed that in polymer nanocomposites with chemically dissimilar graft and matrix chains the dispersion−aggregation transition is a sharp, first-order transition, while the wetting− dewetting transition is a gradual transition occurring over a broad temperature range. Furthermore, using simulations, we showed that at the dispersion−aggregation transition the extent of wetting in the chemically dissimilar composites is identical to the extent of wetting of an equivalent chemically identical (i.e., athermal) composite. In a follow-up simulation study, we demonstrated that by tuning the attractive composition of the graft and matrix chains, we have the capability to tailor the extent of wetting in the composites and, as a result, program into the composite the desired dispersion−aggregation transition.29 Universally, we saw that irrespective of the graft−matrix chain composition, the dispersion−aggregation transition occurred when the extent of wetting equaled that of an analogous athermal system. Additionally, at the dispersion− aggregation transition, the graft and matrix chain conformations of the chemically dissimilar composites were also identical to the equivalent chemically identical (athermal) composite. Received: August 31, 2016 Revised: November 22, 2016

A

DOI: 10.1021/acs.macromol.6b01920 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules

interaction parameter, χeff GM, that we calculate from the simulations and theory (described in section II.C), this model enthalpic interaction parameter χGM is not a function of the graft or matrix design, grafting density, particle size, or composite composition. χGM simply captures how strong the attractive graft−matrix interactions are compared to the attractive graft−graft or matrix−matrix interactions, while χeff GM characterizes the effective monomer−monomer interactions due to the overall enthalpic and entropic driving forces in the composite. All simulations in this paper are of polymer grafted particles with diameter D = 5d, graft length Ngraft = 10, and grafting density 0.76 chains/d2 in a chemically dissimilar matrix of length Nmatrix = 50. The filler fraction of the nanocomposites in this study is maintained at ϕG = 0.13. The filler fraction ϕG is defined as

While all of the above experimental and simulation characterization of the wetting−dewetting and dispersion− aggregation transitions were done based on the structure of the polymer nanocomposites (e.g., wet monomer fraction, structure factors at low wave vector), to the best of our knowledge, there has not been a study establishing the universal thermodynamic signatures of the wetting−dewetting and dispersion−aggregation transitions in chemically dissimilar composites. Such an understanding is valuable not only from a fundamental perspective but also in enabling programmable structure and thermodynamics into the polymer nanocomposites via molecular design and engineering of their components. In this paper, we go beyond these previous studies to show that the thermodynamics quantified by the effective graft− matrix interaction parameter, χeff GM, calculated using trajectories from coarse-grained molecular dynamics simulations along with polymer reference interaction site model (PRISM) theory and random phase approximation (RPA) theory, captures a dispersion−aggregation transition that is consistent with the dispersion−aggregation transition found purely from structure eff factors. The χGM obtained from both PRISM and RPA formalisms correctly predicts an effective graft−matrix attraction for the dispersed morphologies and a sharp transition to effective graft−matrix repulsion for the aggregated morphologies. The effective particle−particle interactions captured in the second virial coefficient and the particle− particle potential of mean force present a picture consistent with effective graft−matrix interactions and also show a signature of the wetting−dewetting transition. All of these effective interaction parameters show that universally at the, dispersion−aggregation transition, the chemically dissimilar graft−matrix composites display effective interactions that are identical to those of the equivalent athermal composite, irrespective of the graft and matrix chain composition or the dispersion−aggregation transition point. The paper is organized as follows. In section II we present the model and the details of the simulations and the theoretical formalisms for calculating effective interactions. In section III we present the results for homopolymer chemically dissimilar composites and random copolymer chemically dissimilar composites. In the last section of the paper we present the major conclusions of this work and a brief discussion of the implications of the results presented in this paper.

ϕG =

Model. We use a generic coarse-grained model to represent polymer grafted spherical nanoparticles in a polymer matrix. The model used in this study is identical to the models we described in refs 28 and 29. Briefly, the nanoparticle (denoted as P) is modeled as a spherical rigid-body of several 1d diameter beads (where d ≈ 1 nm), and both the graft (G) and matrix (M) polymers are modeled as fully flexible bead−spring chains30,31 with beads of diameter d. All attractive graft or matrix monomer−monomer interactions are modeled using the Lennard-Jones (LJ) potential with a cutoff of 2.5d while all athermal monomer−monomer interactions are modeled via a purely repulsive Weeks−Chandler−Andersen (WCA) potential with a cutoff of 21/6d. In all composites, the particle−particle and particle−(graft or matrix) monomer interactions are modeled using the WCA potential. We describe the model graft−matrix interactions via a purely enthalpic interaction parameter, χGM, defined as 1 (εGG + εMM) 2

(2)

where Vgrafts is the total volume occupied by graft monomers across all graft chains and grafted particles and Vmatrix is the total volume of all matrix monomers across all matrix chains. The graft and matrix chains are either homopolymers or random copolymers of attractive and athermal monomers; each random copolymer graft or matrix chain has an independently randomized sequence of attractive and athermal monomers. The compositions of the graft and matrix polymers, described via the fraction of attractive beads per chain, are denoted as f G and f M, respectively, both ranging from 0.0 to 1.0. For example, f G = f M = 1.0 corresponds to the fully attractive homopolymer graft and matrix case (i.e., chemically dissimilar homopolymer graft−matrix case) and f G = f M = 0.0 corresponds to the fully athermal homopolymer case (i.e., chemically identical graft−matrix case). Figure 1 shows schematically the model and the types of graft−matrix composites studied in this paper. Simulation. We conduct NVT Brownian dynamics simulations in the canonical ensemble using the HOOMD-blue and LAMMPS packages.32−35 Briefly, we initialize each system at an extremely low density and high temperature by randomly placing polymer grafted particles and matrix chains in a cubic simulation box. After isotropic compression of the simulation box to the desired total volume fraction, we geometrically anneal the composite for 50 × 106 timesteps from T* = 5.0 to T* = 1.0. All data presented are calculated from sampling configurations every 0.5 × 106 timesteps over 50 × 106 timesteps after the annealing stage, and we estimate the error in these calculations via the standard deviation across five block averages. The complete details of our initialization, equilibration, and production methodologies are described in our recent work.28,29 We also summarize the details of our equilibration procedure in the Supporting Information. To explore the parameter space while ensuring the absence of any finite size effects in our calculations, we run sets of simulations at both a “small” and a “large” system size. Using the HOOMD-blue package for the small system size, we simulate 15 grafted particles in 1200 matrix chains (70 875 total beads) which results in a simulation box length of L ≈ 47.3d at a total volume fraction of η = 0.35. The large systems are simulated using the LAMMPS package and have 100 grafted particles in 8000 matrix chains (472 500 total beads) with L ≈ 88.9d at η = 0.35. We highlight that conducting all the simulations presented in this paper with the large system sizes would be extremely prohibitive due to computational intensity of both carrying out those large simulations and analyzing their large trajectory files. Therefore, after we confirmed that the larger systems and smaller systems produced quantitatively similar results, we limited the next set of calculations to simulations on the smaller system size. Lastly, we used two different simulation packages for the small and large systems for the following reasons. We began this study using HOOMD-blue for the small system size, but we were unable to scale our simulations past 1 GPU because, at that time, HOOMD-blue did not support rigid bodies for multi-GPU or multi-CPU simulations. Since our particle model requires rigid bodies to model fixed graft points, we require a package that could provide rigid body features with the possibility of scaling up to multi-CPU or multi-GPU configurations. To this end, we

II. METHODS

χGM = εGM −

Vgrafts Vgrafts + Vmatrix

(1)

where εGM, εGG, and εMM are the well depths of the graft−matrix, graft−graft, and graft−matrix LJ potentials. In contrast to the effective B

DOI: 10.1021/acs.macromol.6b01920 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules

correlation function. We note that the caret above each variable denotes it is a Fourier space variable, so the (k) may not be explicitly written in equations in later sections for brevity. Each of the variables in eq 3 represents a matrix of values corresponding to each pair of “sites” or bead types in a model, with the “·” representing matrix multiplication and “+” the element-wise addition. In this study, using the G, M, and P labels for graft, matrix, and particle, respectively, the 3site Ĉ (k), Ω̂(k), and Ĥ (k) matrices are written as

⎡ c ̂ (k) c ̂ (k) c ̂ (k) ⎤ GM GP ⎥ ⎢ GG Ĉ(k) = ⎢ cGM ̂ (k) cMM ̂ (k) cMP ̂ (k)⎥ ⎥ ⎢ ⎢⎣ cGP ̂ (k) cMP ̂ (k) c PP ̂ (k) ⎥⎦

(4a)

⎡ ω̂ (k) ω̂ (k) ω̂ (k) ⎤ GM GP ⎢ GG ⎥ ̂ ⎢ ̂ (k) ω̂ MM(k) ω̂ MP(k)⎥ Ω(k) = ωGM ⎢ ⎥ ⎢⎣ ωGP ̂ (k) ω̂ MP(k) ω̂ PP(k) ⎥⎦

(4b)

⎡ ρ ρ h ̂ (k) ρ ρ h ̂ (k) ρ ρ h ̂ (k) ⎤ G M GM G P GP ⎢ G G GG ⎥ ⎢ ̂ (k) ρ ρ h ̂ (k) ρ ρ h ̂ (k)⎥ Ĥ (k) = ⎢ ρG ρM hGM ⎥ M M MM M P MP ⎢ ⎥ ̂ (k) ρ ρ h ̂ (k) ρ ρ h ̂ (k) ⎥⎦ ⎢⎣ ρG ρP hGP M P MP P P PP

(4c)

where ρX is the number density of site X, ĉXY(k) is the X−Y direct correlation function, ω̂ XY(k) is the X−Y intramolecular correlation function, and ĥXY(k) is the X−Y total correlation function. The ω̂ GG, ω̂ GP, and ω̂ MM curves are calculated from the trajectories of the Brownian dynamics simulations, using the Debye scattering relation41,42

ω̂XY (k) =

NM NX NY

∑∑∑ m

i

j

sin(krij) krij

(5)

where X and Y represent site types, NM is the total number of molecules containing type X and type Y sites, NX is the total number of sites of type X in molecule m, NY is the total number of sites of type Y in molecule m, Ntot = (NX + NY) if X ≠ Y otherwise Ntot = NX, rij is the distance between sites i and j, and the angle brackets represent ensemble averaging over uncorrelated snapshots in the simulation trajectory. We calculate ω̂ GG, ω̂ GP, and ω̂ MM from three types of simulation trajectories: (Ω̂_1) chemically distinct graft and matrix composite (attractive, f G = f M = 1.0) at the same conditions as the “small” system (described in the previous section); (Ω̂_2) chemically identical graft and matrix composite with athermal WCA interactions between all pairs of sites and f G = f M = 0.0, at the same conditions as the “small” system; (Ω̂_3) a single grafted particle in vacuum or a single matrix chain in vacuum and with purely WCA interactions between the nonbonded pairs of sites. The purpose of comparing these different simulation sources for Ω̂ is primarily to identify how the intramolecular structure of the grafted particles and matrix chains contributes to the overall effective interactions in the composite. For the particle−particle intramolecular correlation function, there is only one particle per “grafted particle sin(x) molecule”; therefore, ω̂ PP (k) = 1 due to lim x = 1. All pairs of site

Figure 1. Schematic of the three primary classes of composites considered in this study. The graft monomers are isotropically grafted to the particle surface; in these images, some of the graft chains in the front are hidden for clarity.

x→0

use the LAMMPS package, which supports rigid bodies for multi-CPU simulations. Analyses. The primary analysis used in this study is the graft− matrix effective interaction parameter, χeff GM, which we calculate via both the polymer reference interaction site model (PRISM)36−39 and random phase approximation (RPA)37,38,40 formalisms. The fundamental PRISM equation is written in Fourier space as

Ĥ (k) = Ω̂(k)· Ĉ(k)· [Ω̂(k) + Ĥ (k)]

1 NMNtot

types that do not exist within the same molecule have ω̂ XY(k) = 0; therefore, ω̂ MP(k) = 0 and ω̂ GM(k) = 0. The total correlation function ĥXY(k) is also calculated from the molecular dynamics trajectories, as it is equivalent to the translated intermolecular pair correlation function ginter XY (r) in real space: inter hXY (r ) = gXY (r ) − 1

(3)

(6)

̂ (k) using a Fourier Once hXY(r) is calculated, it is inverted to hXY transform. The details of this Fourier space inversion are presented in the Supporting Information.

with Ĉ (k) representing the direct correlation function, Ω̂(k) the intramolecular correlation function, and Ĥ (k) the total intermolecular C

DOI: 10.1021/acs.macromol.6b01920 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules Using the Ω̂(k) and Ĥ (k) calculated from the simulation trajectories as described above, eq 3 can be used to calculate the direct correlation function matrix Ĉ (k). Using Ĉ (k), we calculate the PRISM graft− matrix effective interaction parameter: ρ [cGG ̂ (k) + cMM ̂ (k) − 2cGM ̂ (k)] 2

eff χGM (k) =

(7)

where ρ is the total system density and ĉXY(k) represents component XY of the Ĉ (k) matrix in Fourier space.37−39 In this study, we report the low-wavenumber value of the effective interaction parameter χeff GM(k → 0) as it represents the closest analogue to the traditional Flory− Huggins χGM parameter.37−39 Furthermore, rather than extrapolating our low-k data that already contain significant uncertainty, we report the lowest wavenumber in a data set as the k → 0 value; that lowest k value depends on the box size of the simulation. We also calculate the RPA χeff GM using the intramolecular correlation functions: eff 2χGM (k) =

1 1 1 + − ̂ (k) ̂ (k) ϕGωGG ϕM ω̂ MM(k) SMM

(8)

where ϕG and ϕM are the volume fractions of graft and matrix chains and ϕG + ϕM = 1.0. To calculate the matrix−matrix structure factor, ŜMM(k), we take the Fourier transform of the full (i.e., both inter- and intramolecular) pair correlation function:

̂ (k) = 1 + SMM

4πρM k

full (r ) − 1) sin(kr ) dr ∫ r 2(gMM

(9)

Although eq 9 is written in the continuous space format for clarity, we use the discrete version presented in the Supporting Information, when calculating the transform. To characterize the particle−particle interactions in these composites, we also calculate the second virial coefficient, B2, and the particle−particle potential of mean force, WPP(r). The second virial coefficient is calculated from the particle− particle pair correlation function B2 = − 2π

∫ r 2(gPP(r) − 1) dr

(10)

while potential of mean force is calculated using

where nM is the total number of matrix monomers in the system and nM,wet is the number of matrix monomers that wet the grafted layers. This calculation is described in detail in our recent work.28,29

Figure 2. (a) Wet monomer fraction, ϕW, and matrix−matrix structure eff function ŜMM(k → 0). (b) PRISM χeff GM and (c) RPA χGM versus the model enthalpic interaction parameter χGM. In all plots the vertical dashed lines mark the dispersion−aggregation transition point for each composite, while the horizontal dashed lines mark the value of y-axis parameter for an equivalent athermal composite. In subplots b and c, the symbol indicates the type of simulation trajectory used to calculate intramolecular structure factor Ω̂(k): Ω̂_1 (diamonds), Ω̂_2 (squares), and Ω̂_3 (crosses), as described in the Methods section. The data for all the plots in this figure are from the small system size as described in the Methods section.

III. RESULTS Homopolymer-Grafted Nanoparticles in a Homopolymer Matrix. In Figure 2a, we present the low wavenumber (low-k) value of the matrix−matrix structure factor ŜMM(k → 0) and wet monomer fraction ϕW as a function of the enthalpic interaction parameter χGM. These two structure-based parameters characterize the dispersion−aggregation transition and extent of grafted layer wetting by matrix chains, respectively. As we previously reported,28,29 ŜMM(k → 0) is small at the low χGM values, goes through a sharp increase at an intermediate χGM value, and reaches a plateau at higher χGM. Supporting ̂ (k) which does not Information Figure S1 presents the SMM diverge at k → 0 in simulations due to finite system sizes. Because of the lack of divergence, we call the onset of upturn in ŜMM(k → 0) the dispersion to aggregation transition, even though the experimental results in our previous work confirm

that it is indeed a first-order (spinodal) phase transition.28 We also visualize the dispersion−aggregation transition χGM through simulation snapshots as shown in Figure S2. The dispersion−aggregation transition χGM based on the location of the sharp increase in ŜMM(k → 0) is marked with a vertical dashed line in Figure 2a. We note that while the composite system discussed in this section has a dispersion−aggregation transitions which occurs around χGM = 0, for other composites (discussed in following sections), this transition occurs at nonzero values of χGM. The wet monomer fraction, ϕW, in Figure 2a is a measure of the overall extent of wetting in the composite. The ϕW smoothly decreases across the χGM range below and above the dispersion−aggregation transition χGM. Comparing this gradual wetting−dewetting transition to the sharp dispersion aggregation transition lead us to report, for the first time, that these are two distinct transitions in attractive

WPP(r ) = − kBT log(gPP(r ))

(11)

In addition to the above thermodynamic analyses, we also calculate the wet matrix fraction. It quantifies the extent of wetting of the grafted layer by matrix monomers: nM,wet wet monomer fraction = nM (12)

D

DOI: 10.1021/acs.macromol.6b01920 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules

Figure 3. Graft−graft (a−c) and matrix−matrix (d−f) intramolecular pair correlation function ω̂ vs wavenumber k (units of d−1) for Ω̂_1 simulations (diamonds), Ω̂_2 simulations (squares), and Ω̂ _3 simulations (crosses). For the Ω̂_1 case, the color of the line indicates the enthalpic interaction parameter as described by the color bar. See the Methods section in the main text for definitions of Ω̂_1, Ω̂_2, and Ω̂_3.

composites happens to occur at χ GM = 0 for these homopolymer composites, we show later that some random copolymer composites display this “athermal critical point” at χGM < 0. The same behavior was noted in Figure 2a in SMM(k → 0) and ϕW, where the crossover between the attractive and athermal data marked the dispersion−aggregation transition point. We have previously reported that the conformations of the graft and matrix chains in the attractive systems at the transition point are also nearly identical to those in the equivalent athermal composites.28,29 The coincidence of the attractive and athermal χeff GM at the transition further highlight that the competing thermodynamic driving forces bring about an effectively athermal state at the dispersion−aggregation transition, leading to the composite morphology that is driven purely by entropy. It is important to note that the athermal eff composites have a nonzero and positive χGM , which is consistent with the fact that these systems have Nmatrix = 5Ngraft, a regime where entropy can drive the polymer grafted particles to aggregate in a chemically identical polymer matrix.5−8 We now discuss why the choice of Ω̂ has little to no effect on the effective interaction parameter χeff GM. One would expect that variations in Ω̂, the intramolecular structure factor which captures the structure and conformations of the graft and matrix chains, would affect the effective graft−matrix interactions. In recent work we have shown that with increasing χGM, the distribution of graft chains’ radius of gyration shifts to lower values as the graft chain conformations are impacted by increasing dewetting of the grafted layer. Since Ω̂ is included in both the PRISM and RPA calculations of χeff GM, one would assume that variations in the graft and matrix chain conformations should directly impact the effective graft−matrix interactions. Despite this reasoning, we do not see a large impact of the choice of Ω̂ on χeff GM. To understand this, we first present in Figure 3 the graft−graft and matrix−matrix intramolecular correlation functions (ω̂ GG(k) and ω̂ MM, respectively) used in the calculation of χeff GM. Each row shows three “views” of the same data (at low-, mid-, and long-range k) in order to clearly display the full behavior of each ω̂ (k) curve. For the data associated with the attractive composite (Ω̂_1)

graft−matrix polymer nanocomposites.28 In Figure 2a, the horizontal red and black dashed lines are the ŜMM(k → 0) and ϕW value for the equivalent athermal composite (i.e., f G = 0.0, f M = 0.0), where all monomer−monomer interactions are modeled as purely repulsive using the Weeks−Chandler− Andersen (WCA) potential. We see that the SMM(k → 0) and ϕW for the attractive composite equals that of the athermal composite at the dispersion−aggregation transition point. In Figures 2b and 2c we show the graft−matrix effective interaction parameter χeff GM calculated using the PRISM and RPA formalisms, respectively, as a function of the model enthalpic interaction parameter, χGM. The three lines in each plot correspond to the three types of simulations used to calculate the intramolecular correlation function, Ω̂, as described in the Methods section. Overall, we see only small effects of the choice of intramolecular correlation function Ω̂ on the value of χeff GM and discuss the reasons underlying this trend later in this section where we present Figure 3. In Figures 2b and 2c, the PRISM and RPA χeff GM exhibit qualitatively the same behavior. The χeff GM increases at low values of χGM, undergoes a sharp increase at intermediate χGM, and plateaus at higher values of χGM. For both methods, the graft−matrix interactions are effectively attractive before the dispersion−aggregation transition (vertical dashed lines) and then become effectively repulsive at higher χGM. The quantitative differences between the methods arise from the basic fact that RPA is a mean-field formalism while PRISM includes concentration fluctuations at all length scales. In order to verify the absence of any system size effects, we repeated these χeff GM calculations with larger simulation sizes using a different simulation package, with ∼5 times the number of simulation beads (∼475 000 beads), and observe quantitatively nearly identical χeff GM (see Figure S3). Furthermore, in Figure S4, we found that varying the density (volume fraction) of our simulation had only minor effects on the calculated χeff GM for systems that are in the “melt-like” regime (η ≈ 0.35−0.45). Strikingly, at the dispersion−aggregation transition point, the eff χeff GM values of the attractive systems are equal to the χGM of the equivalent athermal composites (horizontal dashed lines). While the crossover point between the attractive and athermal E

DOI: 10.1021/acs.macromol.6b01920 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules simulations (diamond symbols), the color of the line indicates the value of the enthalpic interaction parameter χGM. For ω̂ GG (Figure 3a−c), the midrange peak in ω̂ GG shrinks and the entire curve shifts to lower k as χGM increases. This matches our observation of the graft chains stretching to accommodate the increasing number of matrix chains wetting the grafted layer with increasing χGM.29 ω̂ GG of the attractive systems matches most closely with the athermal melt composite (Ω̂_2) at high χGM (triangles) and the athermal single grafted particle/matrix in vacuo (Ω̂_3) system (crosses) at low χGM. These data show that as χGM increases, the graft chain conformations transition from having collapsed melt-like statistics to extended solutionlike statistics. This is in-line with the graft Rg2 distributions presented previously for this system.29 In that same study, we also observed that the chain conformations of the matrix chains were relatively constant at all χGM, and this observation is reflected in the ω̂ MM data (Figure 3d−f). All melt cases with either attractive or athermal interactions (Ω̂_1 and Ω̂_2) show nearly identical behavior at all k. Only the Ω̂_3 system shows any difference for the ω̂ MM data, displaying expanded chain statistics at low-k and mid-k. Despite these differences between the various ω̂ GG and ω̂ MM, the χeff GM data in Figures 2b and 2c show little variation with changing Ω̂. This is explained as follows. First, the largest ̂ differences in χeff GM when comparing the different Ω correspond to the Ω̂_3 case (Figure 2b,c crosses). The Ω̂_3 ω̂ GG (Figure 3a−c, crosses) are similar to the Ω̂_1 at low χGM, while the Ω̂_3 ω̂ MM (Figure 3d−f, crosses) deviate from the other Ω̂ at low-k. It follows that if chain structure has any effect on the effective graft−matrix interactions, variations in the matrix chain structure should have the largest effect on the effective graft−matrix interactions due to the matrix being the majority component in the system. In addition, as discussed above, the eff presented χeff GM values are the k → 0 limit of the full χGM(k) curves and are therefore more likely connected to the longwavelength dispersion−aggregation trends (SMM(k → 0)) rather the local wetting−dewetting (ϕW) and mixing−demixing trends. Based on Figure 3, it is clear that the largest differences in the ω̂ GG are in the mid-k range, which would not contribute to the χeff GM value at the k → 0 limit. In contrast, the differences between the various Ω̂ at lowest-k are negligibly small for most of the systems. Going back to Figure 2, the sharp increase in χeff GM mimics the sharp increase in SMM(k → 0) (Figure 2a), indicating that χeff GM describes the dispersion−aggregation transition rather than the wetting−dewetting transition. This result may not be intuitive as one would expect the graft−matrix interaction parameter to mirror the data measuring the graft−matrix mixing−demixing (ϕW). Since we report the traditional, k → 0 limit of the PRISM and RPA χeff GM, this interaction parameter is linked to the longwavelength, macrophase separation (aggregation) behavior of the composite rather the local mixing−demixing (wetting− eff , the particle−particle dewetting). However, unlike χGM potential of mean force, WPP, in Figure 4a contains the signature of both the wetting−dewetting and dispersion− aggregation behavior of the composites. In Figure 4a at the low χGM = −0.4 (diamonds), WPP has a single attractive well at interparticle distance r ∼ 16d surrounded by a large repulsive barrier at lower r and a small barrier at higher r. The large, low-r repulsion is due to the loss of both the attractive graft−matrix interactions and graft−matrix mixing entropy as the particles are brought closer together. At the dispersion−aggregation transition, WPP (downward triangles) maintains the same shape

Figure 4. (a) Particle−particle potential of mean force WPP (units of kBT) as a function of center-to-center distance r (units of d) and (b) particle−particle second virial coefficient B2 as a function of the enthalpic interaction parameter χGM. In part a, the colors and symbols indicate χGM: 0.4 (red upward triangles), 0.0 (green downward triangles), −0.4 (blue diamonds), athermal (black circles). In part b, the vertical line marks the dispersion−aggregation transition, and the horizontal line marks the second virial coefficient for the equivalent athermal composite. The data for all the plots in this figure is from the large system size as described in the Methods section.

as WPP at lower χGM, but the entire curve is shifted such that the well is at r ∼ 12d. The shift in the attractive well of WPP to smaller interparticle distance with increasing χGM is due to the reduced wetting of the grafted layers by matrix chains and the collapsing grafted layer. Yet again at the dispersion−aggregation transition, the attractive simulations have the same WPP as the analogous athermal system (circles). At the highest χGM values (upward triangles), WPP has multiple wells which are stronger/ deeper than the lower χGM data. The shift in WPP to lower r indicates the complete dewetting and collapse of the grafted layer while the presence of multiple wells is indicative of strong particle aggregation and ordering. Lastly, the second virial coefficient, B2 (Figure 4b), also predicts the same qualitative picture seen with χeff GM from both the PRISM and RPA formalisms and the particle−particle potential of mean force WPP. The B2 is positive (i.e., effective particle−particle repulsion) at the lowest χGM in agreement with SMM(k → 0) which predicts effective attraction between the graft and matrix chains promoting grafted layer wetting by F

DOI: 10.1021/acs.macromol.6b01920 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules matrix chains and, in turn, particle dispersion. With increasing χGM, B2 becomes negative, i.e., effective particle−particle attraction, again agreeing with the aggregated morphologies seen visually and as predicted by effective graft−matrix and effective particle−particle interactions. At the dispersion− aggregation transition point (dashed vertical line), B2 is slightly negative and is close to the value for the equivalent athermal composite (dashed horizontal line). As was discussed above, the interactions for the athermal system are expected to be indicative of the onset of particle aggregation, and these athermal B2 values support that assertion. Even though the overall trend in B2 fits well with the trends in χeff GM in describing the dispersion−aggregation behavior of the composite, there is a discrepancy between the athermal (WCA interactions) composite and attractive (LJ interactions) composite’s B2 values at the transition point. The fluctuations in WPP, which are indicative of the fluctuations in gPP(r) used to calculate B2, could be large enough to explain the differences in Figure 4b between the two. So far, for a chemically distinct attractive homopolymer graft−matrix composite, χeff GM and B2 describe the dispersion aggregation transition only, and WPP clearly exhibits signatures of both wetting−dewetting and dispersion−aggregation. The comparison of these analyses display the importance of considering multiple microscopic and macroscopic (shortand long-wavelength) analyses when evaluating the effective interactions in phase separating soft matter systems. Treating the chemically distinct attractive homopolymer composites and the chemically similar athermal homopolymer composites as the two extremes cases of varying enthalpic driving forces with constant entropic driving forces, next we investigate composites with varying fractions of attractive and athermal monomers in the random copolymer graft and matrix chains. This allows us to tune the enthalpic driving forces that compete with the entropic driving forces governing wetting− dewetting and dispersion−aggregation. Random Copolymer Grafted Nanoparticles in a Random Copolymer Matrix. Structurally, we have showed that in polymer nanocomposites with random copolymers of attractive and athermal monomers in the graft and matrix chains, the dispersion−aggregation and wetting−dewetting transitions can be tuned via the graft and matrix chains’ attractive monomer compositions (f G and f M, respectively).29 In Figures 5a and 5b we show ŜMM(k → 0) and ϕW for four different chain compositions, as described in the legend. If the f G/f M ratio is kept constant (compare f G = 1.0, f M = 1.0 to f G = 0.7, f M = 0.7), the dispersion−aggregation transition stays constant at χGM = 0 (Figure 5a), but the extent of wetting of the grafted layer in the dispersed state is reduced (Figure 5b). If the f G/f M ratio is varied asymmetrically (compare f G = 0.7, f M = 0.7 to f G = 0.7, f M = 0.3 and f G = 0.7, f M = 0.1), not only is the wetting of the grafted layer reduced, but the dispersion− aggregation transition is also shifted to χGM < 0. The composites with asymmetric f G/f M ratio require a more negative χGM to stabilize the dispersed morphology because the ratio of attractive graft and matrix monomers are significantly out of balance. For such asymmetric f G/f M ratio, at χGM = 0, even though the graft and matrix monomers are equally attractive to themselves (G−G and M−M) and each other (G−M), there are far more attractive graft monomers than attractive matrix monomers when f G/f M ratio is greater than 1. At χGM = 0, in order to maximize the overall favorable interactions (and entropy) of the system, the composite takes

Figure 5. (a) ŜMM(k → 0), (b) wet monomer fraction ϕW, (c) PRISM eff χeff GM, and (d) RPA χGM versus the enthalpic interaction parameter χGM. In all plots the vertical dashed lines mark the dispersion−aggregation transition point for each composite, while the horizontal dashed lines mark the value of that parameter for an equivalent athermal composite. In all subplots the symbol indicates the graft−matrix composition as shown by the legend. The data for all the plots in this figure is from the small system size. The χeff GM in parts c and d were calculated using using Ω̂_1, as described in the Methods section.

on an aggregated morphology as it can make many more attractive G−G contacts than G−M contacts. It is only at significantly reduced (negative) χGM that the enthalpic driving forces for graft−matrix mixing are strong enough to stabilize dispersion. Overall, the attractive compositions of the graft and matrix chains enable tunability in both the extent of wetting of the grafted layer and the dispersion−aggregation transition. For all f G and f M compositions, the χeff GM calculated using PRISM formalism begin at approximately the same value (≈ −0.0075) and remain constant with increasing χGM until the dispersion−aggregation transition of that composite. At the transition χGM, the PRISM χeff GM sharply increase and all compositions plateau again at similar values (≈0.006). The RPA χeff GM data in Figure 5d display qualitatively similar behavior to those in Figure 5c but do not show a plateau at low χGM, instead displaying rapidly decreasing magnitude with increasingly large error bars. Based on χeff GM for the f G = 1.0 and f M = 1.0 case with the larger simulations (Figure S3), these plateau values are not dependent on the system size but are characteristic of the effective interactions for a dispersed or aggregated composite at these Ngraft, Nmatrix, particle sizes, and grafting densities. For all chain compositions, the PRISM and RPA methods display the correct switch from negative (effective graft−matrix attraction) to positive (effective graft− matrix repulsion) values with increasing χGM. Most importantly, all of the chain compositions and both methods for calculating χeff GM display that at the dispersion−aggregation transition, irrespective of the value of the graft and matrix chain compositions or the value of transition χGM, the χeff GM of the attractive composite systems match those of the equivalent athermal composites (indicated by horizontal dashed line). The fact that this behavior occurs across all graft and matrix G

DOI: 10.1021/acs.macromol.6b01920 Macromolecules XXXX, XXX, XXX−XXX

Macromolecules



compositions studied, even those where the dispersion− aggregation transition occurs at nonzero χGM, highlights that the effectively athermal nature of these attractive composites at the dispersion−aggregation transition is a universal feature of these systems.

Article

ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.macromol.6b01920. Some simulation details and additional result figures (PDF)

IV. DISCUSSION AND CONCLUSIONS In summary, we have shown in this study that at the dispersion−aggregation transition point, the extent of wetting of the grafted layer, effective graft−matrix interaction, and effective particle−particle interaction are all identical to those of an athermal composite, regardless of the graft or matrix chain composition or the χGM of the dispersion−aggregation transition. At this athermal “critical” point, the competing enthalpic interactions are in perfect balance, bringing about a state that is driven purely by entropy. The athermal “critical” point serves as a tipping point where making the effective graft−matrix interactions any more (less) repulsive or the particle−particle interactions any less (more) repulsive will drive the composite toward aggregation (dispersion). These results combined with our previous studies identify chemically dissimilar graft−matrix composites as highly configurable alternatives to more traditional chemically identical composites. We have shown that by simply varying the composition of the graft and matrix chains, the thermodynamics and, in turn, the equilibrium phase behavior and wetting/dewetting of these composites can be tuned to fit specifications of a particular application. While there have been previous theoretical and experimental studies of chemically dissimilar graft−matrix composites, to the best of our knowledge, no prior studies have taken advantage of the thermodynamic and structural programmability that our results have identified. Bockstaller and co-workers reported greatly improved fracture toughness in chemically dissimilar graft− matrix composites when compared their chemically identical counterparts.22 Archer and co-workers studied the dynamics of a chemically dissimilar composite system and found that even for long, entangled polymer matrices, the addition of chemically dissimilar polymer grafted particles significantly altered the short and long time scale dynamics of the matrix chains.27 In the nanocomposite community, it is well-known that the interface between the filler material and matrix in a composite plays a role in determining its mechanical and rheological properties, so tuning the wetting−dewetting in a composite is essentially equivalent to tuning its macroscopic properties.1,2,4,7 So while these studies demonstrate that increased wetting can yield modified or improved material properties, our results demonstrate that chemically dissimilar composites go beyond this idea and further offer the ability to directly program in the desired material properties via direct control of the wetting/ dewetting in the composite. In addition, composites with chemically dissimilar graft and matrix chains present the possibility of engineering thermo- or solvoresponsive composites with precisely programmed transitions. While chemically identical composites generally require physical modifications (i.e., resynthesis) in order to switch from dispersed to aggregated states, even with fixed graft and matrix compositions, the chemically dissimilar composite morphology can be tuned via varying χGM (e.g., via temperature, solvent quality, pH, etc.). Overall, our results in this paper present strong evidence supporting the advantages in the use of structurally and thermodynamically programmable chemically dissimilar graft and matrix polymer nanocomposites.



AUTHOR INFORMATION

Corresponding Author

*E-mail [email protected] (A.J.). ORCID

Arthi Jayaraman: 0000-0002-5295-4581 Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS T.B.M. acknowledges partial financial support from NSF graduate fellowship (Award DGE 1144083). This work used the Farber supercomputer at University of Delaware.



REFERENCES

(1) Kumar, S. K.; Krishnamoorti, R.; Prausnitz, J. M.; Doherty, M. F.; Segalman, M. A. Annu. Rev. Chem. Biomol. Eng. 2010, 1, 37−58. (2) Jancar, J.; Douglas, J. F.; Starr, F. W.; Kumar, S. K.; Cassagnau, P.; Lesser, A. J.; Sternstein, S. S.; Buehler, M. J. Current Issues in Research on Structure-Property Relationships in Polymer Nanocomposites. Polymer 2010, 51, 3321−3343. (3) Krishnamoorti, R.; Vaia, R. A. Polymer Nanocomposites. J. Polym. Sci., Part B: Polym. Phys. 2007, 45, 3252−3256. (4) Crosby, A. J.; Lee, J. Y. Polymer Nanocomposites: The “Nano” Effect on Mechanical Properties. Polym. Rev. 2007, 47, 217−229. (5) Karatrantos, A.; Clarke, N.; Kröger, M. Modeling of Polymer Structure and Conformations in Polymer Nanocomposites from Atomistic to Mesoscale: A Review. Polym. Rev. 2016, 56, 385−428. (6) Ganesan, V.; Jayaraman, A. Theory and Simulation Studies of Effective Interactions, Phase Behavior and Morphology in Polymer Nanocomposites. Soft Matter 2014, 10, 13−38. (7) Kumar, S. K.; Jouault, N.; Benicewicz, B.; Neely, T. Nanocomposites with Polymer Grafted Nanoparticles. Macromolecules 2013, 46, 3199−3214. (8) Green, P. F. The Structure of Chain End-Grafted Nanoparticle/ Homopolymer Nanocomposites. Soft Matter 2011, 7, 7914−7926. (9) Ferrier, R. C.; Koski, J.; Riggleman, R. A.; Composto, R. J. Engineering the Assembly of Gold Nanorods in Polymer Matrices. Macromolecules 2016, 49, 1002−1015. (10) Ferrier, R. C.; Huang, Y.; Ohno, K.; Composto, R. J. Dispersion of Pmma-Grafted, Mesoscopic Iron-Oxide Rods in Polymer Films. Soft Matter 2016, 12, 2550−2556. (11) Martin, T. B.; Jayaraman, A. Effect of Matrix Bidispersity on the Morphology of Polymer-Grafted Nanoparticle-Filled Polymer Nanocomposites. J. Polym. Sci., Part B: Polym. Phys. 2014, 52, 1661−1668. (12) Lin, B.; Martin, T. B.; Jayaraman, A. Decreasing Polymer Flexibility Improves Wetting and Dispersion of Polymer-Grafted Particles in a Chemically Identical Polymer Matrix. ACS Macro Lett. 2014, 3, 628−632. (13) Li, Y.; Krentz, T. M.; Wang, L.; Benicewicz, B. C.; Schadler, L. S. Ligand Engineering of Polymer Nanocomposites: From the Simple to the Complex. ACS Appl. Mater. Interfaces 2014, 6, 6005. (14) Tao, P.; Li, Y.; Siegel, R. W.; Schadler, L. S. Transparent Luminescent Silicone Nanocomposites Filled with Bimodal PdmsBrush-Grafted Cdse Quantum Dots. J. Mater. Chem. C 2013, 1, 86−94. (15) Martin, T. B.; Dodd, P. M.; Jayaraman, A. Polydispersity for Tuning the Potential of Mean Force between Polymer Grafted Nanoparticles in a Polymer Matrix. Phys. Rev. Lett. 2013, 110, 018301. H

DOI: 10.1021/acs.macromol.6b01920 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules

(36) Schweizer, K. S.; Curro, J. G. Prism Theory of the Structure, Thermodynamics, and Phase-Transitions of Polymer Liquids and Alloys. Atomistic Modeling of Physical Properties 1994, 116, 319−377. (37) Curro, J. G.; Schweizer, K. S. An Integral-Equation Theory of Polymer Blends - Athermal Mixtures. Macromolecules 1990, 23, 1402− 1411. (38) Schweizer, K. S.; Curro, J. G. Micropscopic Theory of the Structure, Thermodynamics, and Apparent X-Parameter of Polymer Blends. Phys. Rev. Lett. 1988, 60, 809−812. (39) Curro, J. G.; Schweizer, K. S. Theory for the Chi Parameter of Polymer Blends - Effect of Attractive Interactions. J. Chem. Phys. 1988, 88, 7242−7243. (40) de Gennes, P.-G. Scaling Concepts in Polymer Physics; Cornell University Press: 1979. (41) McCoy, J. D.; Honnell, K. G.; Curro, J. G.; Schweizer, K. S.; Honeycutt, J. D. Single-Chain Structure in Model Polyethylene Melts. Macromolecules 1992, 25, 4905−4910. (42) Flory, P. J. Statistical Mechanics of Chain Molecules; John Wiley & Sons, Inc.: 1969.

(16) Sunday, D.; Ilavsky, J.; Green, D. L. A Phase Diagram for Polymer-Grafted Nanoparticles in Homopolymer Matrices. Macromolecules 2012, 45, 4007−4011. (17) Meng, D.; Kumar, S. K.; Lane, J. M. D.; Grest, G. S. Effective Interactions between Grafted Nanoparticles in a Polymer Matrix. Soft Matter 2012, 8, 5002−5010. (18) Hore, M. J. A.; Frischknecht, A. L.; Composto, R. J. Nanorod Assemblies in Polymer Films and Their Dispersion-Dependent Optical Properties. ACS Macro Lett. 2012, 1, 115−121. (19) Akcora, P.; Liu, H.; Kumar, S. K.; Moll, J.; Li, Y.; Benicewicz, B. C.; Schadler, L. S.; Acehan, D.; Panagiotopoulos, A. Z.; Pryamitsyn, V.; Ganesan, V.; Ilavsky, J.; Thiyagarajan, P.; Colby, R. H.; Douglas, J. F. Anisotropic Self-Assembly of Spherical Polymer-Grafted Nanoparticles. Nat. Mater. 2009, 8, 354−U121. (20) Zhang, R.; Lee, B.; Bockstaller, M. R.; Kumar, S. K.; Stafford, C. M.; Douglas, J. F.; Raghavan, D.; Karim, A. Pattern-Directed Phase Separation of Polymer-Grafted Nanoparticles in a Homopolymer Matrix. Macromolecules 2016, 49, 3965−3974. (21) Dang, A.; Ojha, S.; Hui, C. M.; Mahoney, C.; Matyjaszewski, K.; Bockstaller, M. R. High-Transparency Polymer Nanocomposites Enabled by Polymer-Graft Modification of Particle Fillers. Langmuir 2014, 30, 14434−14442. (22) Ojha, S.; Dang, A.; Hui, C. M.; Mahoney, C.; Matyjaszewski, K.; Bockstaller, M. R. Strategies for the Synthesis of Thermoplastic Polymer Nanocomposite Materials with High Inorganic Filling Fraction. Langmuir 2013, 29, 8989−8996. (23) Heo, K.; Miesch, C.; Emrick, T.; Hayward, R. C. Thermally Reversible Aggregation of Gold Nanoparticles in Polymer Nanocomposites through Hydrogen Bonding. Nano Lett. 2013, 13, 5297− 5302. (24) Hore, M. J. A.; Composto, R. J. Using Miscible Polymer Blends to Control Depletion−Attraction Forces between Au Nanorods in Nanocomposite Films. Macromolecules 2012, 45, 6078−6086. (25) Borukhov, I.; Leibler, L. Enthalpic Stabilization of Brush-Coated Particles in a Polymer Melt. Macromolecules 2002, 35, 5171−5182. (26) Borukhov, I.; Leibler, L. Stabilizing Grafted Colloids in a Polymer Melt: Favorable Enthalpic Interactions. Phys. Rev. E: Stat. Phys., Plasmas, Fluids, Relat. Interdiscip. Top. 2000, 62, R41−R44. (27) Mangal, R.; Srivastava, S.; Archer, L. A. Phase Stability and Dynamics of Entangled Polymer-Nanoparticle Composites. Nat. Commun. 2015, 6, 7198. (28) Martin, T. B.; Mongcopa, K. I. S.; Ashkar, R.; Butler, P.; Krishnamoorti, R.; Jayaraman, A. Wetting−Dewetting and Dispersion−Aggregation Transitions Are Distinct for Polymer Grafted Nanoparticles in Chemically Dissimilar Polymer Matrix. J. Am. Chem. Soc. 2015, 137, 10624−10631. (29) Martin, T. B.; Jayaraman, A. Tuning the Wetting−Dewetting and Dispersion−Aggregation Transitions in Polymer Nanocomposites Using Composition of Graft and Matrix Polymers. Mater. Res. Express 2016, 3, 034001. (30) Grest, G. S.; Kremer, K. Molecular-Dynamics Simulation for Polymers in the Presence of a Heat Bath. Phys. Rev. A: At., Mol., Opt. Phys. 1986, 33, 3628−3631. (31) Ceperley, D.; Kalos, M. H.; Lebowitz, J. L. Computer Simulation of the Dynamics of a Single Polymer Chain. Phys. Rev. Lett. 1978, 41, 313−316. (32) Glaser, J.; Nguyen, T. D.; Anderson, J. A.; Lui, P.; Spiga, F.; Millan, J. A.; Morse, D. C.; Glotzer, S. C. Strong Scaling of GeneralPurpose Molecular Dynamics Simulations on Gpus. Comput. Phys. Commun. 2015, 192, 97−107. (33) Nguyen, T. D.; Phillips, C. L.; Anderson, J. A.; Glotzer, S. C. Rigid Body Constraints Realized in Massively-Parallel Molecular Dynamics on Graphics Processing Units. Comput. Phys. Commun. 2011, 182, 2307−2313. (34) Anderson, J. A.; Lorenz, C. D.; Travesset, A. General Purpose Molecular Dynamics Simulations Fully Implemented on Graphics Processing Units. J. Comput. Phys. 2008, 227, 5342−5359. (35) Plimpton, S. Fast Parallel Algorithms for Short-Range Molecular Dynamics. J. Comput. Phys. 1995, 117, 1−19. I

DOI: 10.1021/acs.macromol.6b01920 Macromolecules XXXX, XXX, XXX−XXX