Microscopic Structure, Conformation, and Dynamics of Ring and

Article pubs.acs.org/Macromolecules

Microscopic Structure, Conformation, and Dynamics of Ring and Linear Poly(ethylene oxide) Melts from Detailed Atomistic Molecular Dynamics Simulations: Dependence on Chain Length and Direct Comparison with Experimental Data Dimitrios G. Tsalikis,† Thanasis Koukoulas,† Vlasis G. Mavrantzas,*,†,‡ Rossana Pasquino,§,∥ Dimitris Vlassopoulos,§,⊥ Wim Pyckhout-Hintzen,# Andreas Wischnewski,# Michael Monkenbusch,# and Dieter Richter# †

Department of Chemical Engineering, University of Patras and FORTH-ICE/HT, GR 26504, Patras, Greece Particle Technology Laboratory, Department of Mechanical and Process Engineering, ETH-Z, CH-8092 Zürich, Switzerland § FORTH, Institute for Electronic Structure and Laser, Heraklion 71110, Greece ∥ Department of Chemical, Materials and Industrial Engineering, University of Napoli Federico II, P.le Tecchio 80, 80125 Napoli, Italy ⊥ Department of Materials Science & Technology, University of Crete, Heraklion 71003, Greece # Jülich Centre for Neutron Science (JCNS-1) and Institute for Complex Systems (ICS-1), Forschungszentrum Jülich GmbH, 52425 Jülich, Germany ‡

S Supporting Information *

ABSTRACT: We present results from very long (on the order of several microseconds) atomistic molecular dynamics (MD) simulations for the density, microscopic structure, conformation, and local and segmental dynamics of pure, strictly monodisperse ring and linear poly(ethylene oxide) (PEO) melts, ranging in molar mass from ∼5300 to ∼20 000 g/mol. The MD results are compared with recent experimental data for the chain center-of-mass self-diffusion coefficient and the normalized single-chain dynamic structure factor obtained from small-angle neutron scattering, neutron spin echo, and pulse-field gradient NMR, and remarkable qualitative and quantitative agreement is observed, despite certain subtle disagreements in important details regarding mainly internal ring motion (loop dynamics). A detailed normal-mode analysis allowed us to check the degree of consistency of ring PEO melt dynamics with the ring Rouse model and indicated a strong reduction of the normalized mode amplitudes for the smaller mode numbers (compared to the Rouse model scaling), combined with an undisturbed spectrum of Rouse relaxation rates. We have further measured the zero-shear rate viscosity η0 of the PEO-5k and PEO-10k rings at several temperatures and extracted their activation energies. These were compared with the activation energies extracted from the MD simulations via analysis of the temperature dependence of the corresponding Rouse relaxation times of the two rings in the same temperature range.

1. INTRODUCTION

and constraint release that require the presence of chain ends

Because of the absence of chain ends and their loop topology, polymer rings exhibit dynamic and viscoelastic properties that cannot be explained by the reptation theory.1,2 Because of their looped structure, in particular, processes such as reptation through an effective confining tube, contour length fluctuations

cannot directly apply to melts of polymer rings. The absence of

© XXXX American Chemical Society

Received: November 17, 2016 Revised: February 18, 2017

A

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

Article

Macromolecules

for the description of interatomic interactions, the use of graphics processing units (GPUs) to speed up code execution, and the availability of extremely powerful supercomputing systems nowadays. MD allows the study of well-defined systems, i.e., of a precise chain length, architecture, and composition.19−21 This is particularly important for rings for which purity is, by default, an issue. The simulation results are free of any ambiguities regarding contamination by linear analogues. A second advantage of the MD method is that it can provide direct information about the temporal evolution of the system, a feature which is particularly important for polymers that are characterized by a wide spectrum of relaxation times and lengths. For example, in an MD method, one can easily distinguish between the dynamics of atomistic segments and the dynamics of entire chains, can track the dynamics of normal modes, etc. Of course, an inherent limitation of the atomistic MD method is that the time span tracked in the course of the simulation can be short compared to the spectrum of relaxation times characterizing molecular motions in these systems. To overcome this serious limitation, coarse-grained models are typically invoked.22−25 If, however, the molecular length of the simulated polymer is relatively short (corresponding, e.g., to ten entanglement lengths), the method can track the evolution of the system, thus providing very accurate information about its dynamics at both short and long scales. Given, in particular, the current emphasis on the dynamics of even very short polymer ring melts, atomistic MD can prove the method of choice for the study of their physical properties. In addition, and in conjunction with state-of-the-art topological algorithms for analyzing the entanglement structure of the dense phases of polymers, the atomistic MD method can provide molecularlevel information about the type and role of topological constraints that govern dynamics in polymer rings. Through such a study, for example, threading events in pure PEO ring melts and in blends of ring−linear PEO melts have been recently analyzed providing a detailed picture for the effect of topological constraints on ring dynamics, in particular on the observed slow modes at long times.26−29 Furthermore, and similar to experiments, molecular simulations can serve the purpose of testing the predictive capabilities of analytical models30,31 for the conformational and dynamic properties of polymer rings. Despite the above-mentioned advances in analytical, experimental, and simulation methods, a universal description of both structure and dynamics of pure rings is still missing. We seek such a description here where, motivated by the experimental measurements reported in refs 11, 12, and 16, we have carried out detailed atomistic MD simulations of linear and ring PEO melts in order to obtain predictions for their microscopic structure, overall conformation and equilibrium dynamic properties, and their dependence on molar mass. The results have been obtained with a very accurate force field for PEO,32,33 and this will be reflected both in the quality of the results and in their excellent comparison with the experimental data, for practically all properties of interest. In fact, to enable a one-to-one comparison with the reported experimental data, we have simulated PEO melts of exactly the same molecular length as those used in the SANS, NSE, and PFG-NMR studies of refs 11 and 12. In all cases, long simulations were executed (especially for the longest PEO systems), on the order of several microseconds, to reliably track long-time (terminal) dynamics which governs the viscoelastic response.

end-segment relaxation renders them ideal systems for testing fundamental theories of polymer dynamics3 or for understanding the dynamics of important ring structures in nature such as mitochondrial and plasmic DNA.4,5 This explains the great interest in the (static, dynamic, and viscoelastic) properties of this unique class of polymers in the past decade. From an experimental point of view, a major difficulty to overcome in the measurements is the availability of highly purified and monodisperse polymer rings in sufficient quantities, since small contamination of the melt by linear chains can have a dramatic effect on the dynamics and molecular rheology of the melt.6−10 In recent studies,11−18 highly pure hydrogenated and deuterated PEO samples (purity of cyclic product larger than 99.5%) of three different molar masses Mw (equal to 5330, 10 044, and 20 086 g/mol) have been synthesized, characterized by a very small polydispersity index (∼1.01 for the hydrogenated samples and ∼1.03 for the deuterated ones). These have been used to perform small-angle neutron scattering (SANS), neutron spin echo (NSE), pulse field gradient NMR (PFG-NMR), and rheological measurements.11−18 From the experimental data, it has thus been possible to extract the dependence of the mean-squared radius of gyration ⟨Rg2⟩ of a ring PEO molecule on the number of monomers N, the dynamic structure factor S(q,t)/S(q,0) over a very large regime of time scales, and the long-time diffusive behavior. It was found that ⟨Rg2⟩ ∼ N2ν with v = 0.43 ± 0.015. Regarding the N-dependence of the ring center-of-mass selfdiffusion coefficient Dcm, the PFG-NMR measurements indicated that for the two higher Mw samples Dcm scales as Dcm ∼ N−2. The shorter ring (5330 g/mol) has further been reported to be characterized by a significant no-Gaussian behavior. The crossover time in the mean-square displacement (msd) of the rings centers-of-mass quantifying the passage from the subdiffusive regime to Fickian diffusion was also quantified in the three melts. It was argued that ring dynamics is characterized by a fast Rouse relaxation of subunits (loops) and a slower dynamics displaying a lattice animal-like loop displacement.11,12 The corresponding rheological measurements (which probed the linear rheology of several, as pure as currently possible, ring polymers of different chemistry and molar mass)16−18 indicated in all cases a power-law stress relaxation. As far as their zeroshear rate viscosity is concerned, this was represented in the experiments in the form η0,ring/η0,linear versus number of entanglements Z in the linear samples (where η0,ring and η0,linear denote the zero-shear rate viscosity of the ring and of its linear analogue, respectively).16 Although in the unentangled regime linear and ring polymers were observed to follow the same scaling with molar mass, in the entangled regime the experimental data confirmed not only the weaker molecular weight dependence of ring viscosity with respect to the linear analogue but also the universality of the scaling (independence on chemistry). Besides experiments, molecular simulations can tremendously help understand the molecular mechanisms underlying dynamics in melts of polymer rings because of their capability to provide direct information about the temporal evolution of the system at all scales.19−24 Atomistic-level molecular dynamics (MD) simulations, in particular, can provide a wealth of information about the structural, conformational, and dynamic properties of this extremely important class of materials thanks to the use of very accurate potential functions B


Article

Macromolecules

a few polymer rings immersed in a host matrix of longer linear polymers were studied.

The paper is organized as follows: Section 2 reviews the details of the simulated PEO systems and of the molecular model adopted to carry out the simulations and provides some technical details about the MD simulation algorithm itself. Results from the MD simulations for a variety of properties (density, pair correlation functions, mean-squared radius of gyration, normal-mode analysis, dynamic structure factor, segmental and chain center-of-mass msd’s, chain center-ofmass self-diffusion coefficient, and zero shear viscosity) are presented and directly compared with recent experimental data and previous simulation results in section 3. Finally, in section 4 we summarize the main conclusions of the present work and outline future plans.

3. RESULTS A. Density. MD predictions for the density ρ of the simulated ring and linear PEO melts at T = 413 K and P = 1 atm are shown in Table 1. For completeness, in the table we Table 1. MD Results for the Density of Several Short Linear PEO Systems at Various Temperatures along with Some Earlier Simulations26 with the Same Force Field and Experimental Data,40,41 Fully Validating the Atomistic Model Employed in the Present PEO Simulations ρ (g/cm3)

2. MOLECULAR MODEL AND SIMULATED SYSTEMS The ring PEO chains considered in our simulations are represented by the formula -CH2-O-(CH2-CH2-O)N-CH2while the linear analogues by the formula CH3-O-(CH2-CH2O)N-CH3, with N in both cases denoting the number of monomers per chain (the same for both species). Three different chain sizes were considered, equal to N = 120, 227, and 455, respectively (implying strictly monodisperse melts). The corresponding molar masses (they differ slightly between ring and linear molecules) were equal to 5330, 10 044, and 20 086 g/mol, and we will be referring to them in the following simply as the PEO-5k, PEO-10k, and PEO-20k melts, respectively. All MD simulations were carried out with the modified TrAPPE united-atom (UA) force-field proposed by Fischer et al.32,33 using GROMACS.34 The corresponding potential functions describing bond stretching, bond bending, torsional, Lennard-Jones, and electrostatic interactions can be found elsewhere.26,32 Initial configurations were prepared with MAPS.35 The model was thoroughly validated by comparing its predictions for several properties of a number of relatively short pure linear and pure ring PEO systems (oligomers) against reported experimental data. It reproduced experimental data for the density of several short linear PEO systems at various temperatures extremely accurate (within less than 1%). Additional simulation predictions for the temperature dependence of the zero-shear rate viscosity of the linear CH3-O-(CH2CH2-O)10-CH3 PEO system (at P = 1 atm) (obtained by inserting in the corresponding Rouse theory expression simulation predictions for the density and the longest relaxation time of the chain) were in excellent agreement26 with reported experimental data as well as with simulation results through a direct application of the Green−Kubo equation. All MD runs were conducted in the isothermal−isobaric (NPT) statistical ensemble by making use of the Nosé−Hoover thermostat36,37 coupled with the Parrinello−Rahman38 barostat to maintain temperature T and pressure P fixed at their prescribed values of T = 413 K and P = 1 atm. Large cubic simulation cells (subject to periodic boundary conditions) were employed in all cases (containing up to 125 chains in some cases) to ensure that the results are free of any finite system size effects and to reduce the statistical uncertainty of the predicted properties. For example, in the simulations with the largest PEO melt studied here (20k), the simulation cell contained more than 5 × 104 atoms and the MD simulation was carried out for about 2.2 μs, which is about 1.5 times the corresponding longest relaxation time of linear PEO-20k (∼1.4 μs) at T = 413 K. Some of the simulation results for the PEO-20k ring melt presented in this work have already appeared in a recent publication,39 where the microscopic dynamics and topology of

Mw (g/mol)

T (K)

669 669 1046 1806 5330 10044 20086

303.5 343.5 343.35 413 413 413 413

MD ring melts

1.015 1.018 1.018 1.018

± ± ± ±

0.007 0.006 0.006 0.006

MD linear melts

exptl data

± ± ± ± ± ± ±

1.054 1.045 1.067

1.046 1.043 1.076 1.000 1.015 1.017 1.018

0.006 0.006 0.006 0.009 0.009 0.009 0.008

have included some additional results from previous simulations with shorter PEO systems with the same force field as well as several experimental data from various references. The error bars on the simulation data are very small due to the use of large simulation cells and the long duration of the MD runs. The results indicate that our MD simulations predict the density of linear PEO melts very satisfactorily; the maximum deviation from the experimental data is less than 1%. We also see that as Mw increases, ρ increases monotonically but for Mw > 5K a plateau is reached beyond which ρ does not change anymore. This tendency is seen both in the simulation and in the experimental data.40,41 B. Microscopic Structure. In Figures 1 and 2 as well as Figures S.I.1 and S.I.2 we show the simulation predictions for the pair correlation functions gαβ(r) (α, β ∈ {C, O}) in the three PEO melts (linear and ring) which provide direct information about local packing and structure over short and long length scales in the simulated systems. In general, a pair correlation function is the sum of an intramolecular part gintra αβ (from atom pairs that belong to the same chain) and of an intermolecular part ginter (from atom pairs that belong to αβ different chains). The intramolecular contribution is actually reported in terms of the intrachain pair density function intra intra wintra αβ (r) defined through gαβ = wαβ (r)/(n/V) where n denotes the total number of atomistic units in the system and V the volume of the system. Figure 1 shows the simulation results for inter the quantity wintra CC (r) and Figure 2 for the quantity gCC . The very sharp peaks observed in the intramolecular contribution, wintra CC (r), are due to chain connectivity which constrains pair distances within certain limits (depending on the stiffness of the bond stretching and bond bending potentials as well as on the temperature and pressure conditions considered). Of course, by definition, all (total) pair distribution functions should approach unity at long distances, a direct consequence of the absence of long-range structure in melt systems. The first Dirac peak in all plots of Figure 1 and Figure S.I.1 at small distances is due to pairs of atoms (e.g., C−C and C−O) separated by one bond. The second Dirac peak is due to pairs of atoms two bonds apart, while the other sharp peaks reflect conformational C


Article

Macromolecules

Figure 2. Same as with Figure 1 but for the intermolecular pair distribution function ginter αβ (r) of C−C pairs.

Figure 1. Predicted intrachain pair density function w(r) of C−C pairs for all ring (a) and linear (b) PEO systems studied in the present work.

Regarding the intermolecular functions, ginter CC (r), the graphs in Figure 2 show that for both ring and linear PEO melts, these rise more or less smoothly with increasing distance r exhibiting some subtle structural features only up to a certain distance, beyond which they approach the asymptotic value of one. Exactly the same behavior is observed for the C−O and O−O pairs shown in Figure S.I.2. In all intermolecular pair distribution functions, several maxima and minima appear at relatively short distances. The correlation hole effect (humps in the corresponding ginter αβ (r) plots) is evident in all intermolecular correlations and can be explained by considering the first coordination shell defined by the sum of the van der Waals radii of the atoms involved in the correlation. At least three pronounced humps (corresponding to the first, second, and third coordination shells) are discerned in all intermolecular correlations at distances approximately equal to 4.7, 9.4, and 13.7 Å for the C−C pairs (Figure 2), equal to 5.2, 10.2, and 14.7 Å for the O−O pairs (Figure S.I.2), and equal to 5.5, 9.9, and 14.2 Å for the C−O pairs (Figure S.I.2). We further see that for a given chain length the correlation hole effect is more pronounced in the ring than in the corresponding linear system. That is, at distances significantly smaller compared to the mean radius of gyration, the value of ginter αβ (r) for a given ring PEO melt is suppressed relative to its value for the corresponding linear analogue as the segment cloud of the cyclic chain excludes mers of other chains from coming close to the reference mer on the chain. Typically, the distance over which the correlation hole effect persists scales with the square root of the chain length42,43 and, in fact, is commensurate with the mean radius of gyration of the chain. The peaks in ginter αβ (r) at short distances are typically used to characterize local packing

preferences (trans and gauche) of the chains. Because the locations of these peaks reflect chain-length-independent intramolecular characteristics, they are practically the same for all linear and ring PEO systems (5k, 10k, and 20k) studied here. Where the systems differ appreciably is in the overall decay of wintra αβ (r) with distance r which reflects the finite size of the chain segment cloud: the longer the chain, the slower the decay of wintra αβ (r) at larger r values. Furthermore, for the same chain length N, the wintra αβ (r) curves for ring chains extend to shorter distances r than for linear ones which is a direct manifestation of the more compact structure of rings due to their looped conformation. This is more pronouncedly observed in the wintra αβ (r) vs r graphs between ring and linear chains in the case of the longest melt (PEO-20k). How wintra αβ (r) changes with r is analyzed in Figure S.I.3, where the data of Figure 1 for distances r greater than 10 Å and up to twice the mean radius of gyration Rg (Rg ≡ ⟨Rg2⟩1/2) for each melt are plotted in a log−log plot. For the three linear PEO melts, wintra αβ (r) drops with r following a power-law −b behavior wintra with values of b equal to 1.65, 1.39, and αβ (r) ∼ r 1.25 for linear PEO-5k, linear PEO-10k, and linear PEO-20k, respectively. Ring molecules, on the other hand, exhibit a different behavior: for distances r between 10 Å and half Rg, −b intra wintra with values αβ (r) follows a power-law behavior wαβ (r) ∼ r of b equal to 1.40, 1.20, and 0.99 for ring PEO-5k, ring PEO10k, and ring PEO-20k, respectively, but for longer distances (r > 0.5Rg), wintra αβ (r) drops much faster due to the more compact conformation of rings and their limited spatial extent compared to the corresponding linear molecules. D


Article

Macromolecules in the melt. We observe then that the positions of these peaks are the same in all three melts (5k, 10k, and 20k), which suggests that all of them are characterized by the same density (or, equivalently, that in this regime of chain lengths the density does not change with N). A direct comparison of the intrachain pair density functions wintra αβ (r) and of the intermolecular pair distribution functions ginter αβ (r) between linear and ring melts for all atom pairs (C−C, C−O, and O−O) is shown in Figure S.I.4. The graphs have been obtained from the simulations with the PEO-20k melt, but very similar curves were obtained for the rest of the PEO melts studied here. All intrachain pair density functions wintra αβ (r) for the rings are above those for the corresponding linear melts (except from the very long distances), which is due to the looped structure of ring chains that keeps atoms within the same chain closer to each other. At very long distances (beyond that corresponding to the mean radius of gyration), on the other hand, the wintra αβ (r) curves for ring chains fall below those for the corresponding linear ones, which reflects the smaller spatial extent of ring chains. For exactly the same reason (the more compact structure of rings that keeps different ring molecules relatively apart one from the other), the curves for all ginter αβ (r) functions for the rings are below those for the corresponding linear melts. An additional reason for this is that linear analogues tend to interpenetrate heavily. The characteristic distances, however, where the maxima and minima in the ginter αβ (r) vs r curves are observed are identical between ring and linear melts. That linear chains interpenetrate heavily irrespective of their size was confirmed by computing a quantity that can characterize structure at the level of entire chains; as such, we chose the chain center-of-mass pair distribution function, gcm(r) .44 The dependence of gcm(r) on r and its scaling with chain length N for the simulated PEO melts can be seen in Figure 3. The results are noisy because of the relatively small number of chains present in the simulation box. However, they clearly show the absence of any strong correlations between the centers-of-mass of the linear chains. The fact that gcm(r → 0) assumes a nonzero value is a direct proof that different linear chains interpenetrate profusely, which allows their centers-ofmass to come close enough without violating the nonoverlap requirement for atoms.44 The same but to a much lesser degree is observed for the ring melts: indeed at r = 0, gcm(r) for the PEO-5k and PEO-10k ring melts assumes a nonzero value (but well below one); at longer distances, this rises smoothly and monotonically with ring−ring center-of-mass distance r toward its asymptotic value of one. From this, we conclude that ring molecules are not as penetrable by other ring molecules as their linear analogues are by other linear chains or that they are penetrable to some degree but their centers-of-mass cannot easily overlap (this is more obvious in the case of the PEO-20k ring melt). Furthermore, and similar to the ginter αβ (r) plots discussed in Figure 2 and Figure S.I.2, the gcm(r) vs r curves for the ring melts are systematically below those for the linear analogues, which is another indication of the more compact structure of ring chains. From a topological point of view, a quantity that is of interest is the average number K1(r) of neighboring rings whose centerof-mass is within a distance r equal to the ring mean radius of gyration Rg (Rg ≡ ⟨Rg2⟩1/2) from the center-of-mass of a reference ring.22 The variation of K1(r) with ring chain length N (number of monomers per chain) is shown in Figure 4. We see that K1(r = Rg) attains values between 1.75 and 2.75 for the

Figure 3. Chain center-of-mass pair distribution function gcm(r) for the simulated PEO-5k (a), PEO-10k (b), and PEO-20k (c) ring and linear melts.

PEO ring melts addressed here, but these increase as N increases. The simulations reported in ref 22 with a generic bead−spring chain model yielded similar numerical results. For comparison, also shown in Figure 4 is the corresponding result for the linear PEO analogues. For the same value of chain length N, K1(r = Rg) comes out to be about 4 times larger in the linear melts than in the rings. This is also consistent with the results reported in ref 22 where the K1(r) data for the linear melts were shown for a distance r equal to the root-meansquare chain end-to-end distance (and not to the ring mean radius of gyration). C. Static Structure Factor. By taking the Fourier transform of the total pair distribution function gαβ(r) (α, β) ∈ {C, O, H}), we can compute the static structure factor F(q). To E


Article

Macromolecules

thus, in Figure 5 only the curve for the PEO-20k system is depicted. Despite a number of complicating factors and the different weightings at different q regimes in the experiments, simulated F(q) curves can be directly compared against experimental F(q) data obtained either from X-ray diffraction or neutron scattering (NS) measurements. To this, in Figure 5 we also plot experimental points from NS measurements at 1 atm and 363 K reported by Johnson et al.,45 which provide a means for validating the predicted intermolecular pair inter distribution functions gαβ (r). The agreement between simulation and experimental F(q) vs q curves is very satisfactory. It is particularly encouraging that the first peak in F(q) is captured quite well by the simulation. This first peak represents intermolecular correlations while the second and third ones at large q values reflect intramolecular packing. Furthermore, the MD predictions for the location and magnitude of the major intermolecular peak near 1.5 Å−1 and of the intramolecular peaks near 3.2 and 6.3 Å−1 agree with the F(q) pattern reported in ref 43 with the OPLS force field for a linear PEO melt at 500 K and 1 atm. D. Conformational Properties. The conformational properties of the simulated ring and linear PEO melts have been analyzed in terms of the mean-squared chain radius of gyration ⟨Rg2⟩, and the results are shown in Figure 6 as a

Figure 4. Average number of neighboring linear or ring molecules K1(r = Rg) versus number N of monomers per linear or ring PEO molecule.

account for the different atomic species (carbons, oxygens, and hydrogens) present in the simulated systems, the static structure factor F(q) is calculated as the Fourier transform of the function H(r) defined as 3

H (r ) =

3

xaxβfa fβ

∑∑

3

a=1 β=1

(∑γ = 1 xγfγ )2

(gαβ (r ) − 1) (1)

where xα and xβ denote the number fraction of α-type and βtype atoms in the system and fα and fβ the respective scattering factors. Values of xα and fα for the three different elements (C, O, H) in a PEO chain are reported in Table S.I.1. By taking the Fourier transform of the function H(r) F(q) = 1 +

n V

∫0

∞

4πr 2

sin(qr ) H (r ) d r qr

(2)

we computed the X-ray diffraction pattern shown in Figure 5. In eq 2, n/V denotes the total number density of atomistic units in the simulation system. F(q) was seen to be very similar for all three systems studied here (PEO-5k, PEO-10k, and PEO-20k); Figure 6. Dependence of the mean-squared chain radius of gyration ⟨Rg2⟩ on chain length N as obtained from the present MD simulations for ring and linear PEO melts, together with experimental data for the rings from SANS measurements. The scaling of ⟨Rg2⟩ with N for the two types of melts (ring versus linear) is also indicated.

function of chain length N. The values obtained for ⟨Rg2⟩ for all simulated ring and linear PEO melts are reported in Table S.I.2. A direct comparison shows that for the same chain length the average dimension of a ring polymer is always smaller than that of its linear counterpart, clearly pointing to a more compact overall structure of rings. This is best reflected in the values of the characteristic exponents quantifying the power-law dependence of ⟨Rg2⟩ on N for ring and linear PEO melt chains. We find that ⟨R g2⟩ Figure 5. MD simulation predictions for the X-ray diffraction pattern of ring PEO (red line) and comparison with the experimental pattern (black squares) taken from the NS measurements reported by Johnson et al.45 The MD results have been obtained from the simulations with the PEO-20k system.

0.85 ± 0.03 ⎧ , for ring molecules ⎪N ∼⎨ ⎪ 1.0 ± 0.04 , for linear chains ⎩N

⟨Rg2⟩

(3)

That ∼N for the simulated linear PEO melts is a direct manifestation of the Gaussian structure of linear PEO chains. On the other hand, the scaling ⟨Rg2⟩ ∼ N0.85±0.03 for the F

1.0±0.04


Article

Macromolecules rings is consistent with previous simulation,19−22,25 experimental,12,46 and theoretical works.31 We further computed the values of the mean-squared chain end-to-end distance ⟨R2⟩ for the linear PEO melts, and we found that for all systems the results obey the random-coil relation, namely that ⟨R2⟩/⟨Rg2⟩ = 6 to an excellent approximation. In Figure 6, we have also plotted the ⟨Rg2⟩ values obtained from the recent experimental study12 based on SANS data. Since some data in the figure overlap between actual measurements and simulations, the corresponding best fit curve to the experimental points is shown with the correct slope but shifted a little downward for better visualization. Admittedly, the agreement between simulation and measured data is extremely good especially as the chain length increases, which can also be considered as an indirect verification of the purity of the ring PEO samples used in the SANS and NSE measurements. To check the influence of ring closure on the microscopic structure of PEO chains, in Figure S.I.5 we report simulation results for the distribution functions of the two torsional angles (C−O−C−C and O−C−C−O) along the backbone of a PEO molecule from the MD runs with the ring and linear PEO-20k melts (similar plots were obtained for the rest of the systems). For both dihedrals, the resulting probability distribution functions between ring and linear PEO melts are identical (indistinguishable). E. Kratky Plots. The simulation data for the 5k, 10k, and 20k rings are further compared with SANS experimental data in two different representations: (a) in the so-called Kratky representation, i.e., as I(q)q2 versus q (Figure 7), and (b) as I(q)q3 versus qRg (Figure 8).12,14 For the latter figure, we remind the reader that the chain length dependence of Rg was just discussed in Figure 6. As noted by Gooßen at al.,12 in the low qRg regime the data fall on one master curve if properly normalized for volume fraction and number of monomers. At intermediate values of qRg a plateau develops while in the high qRg regime the curves rise rapidly. The simulation data follow the SANS measurements quite accurately, especially for the 10k and 20k melts, except for qRg values greater than approximately 5 where they somehow underestimate the rapid increase of I(q) q3 with qRg. These are due to some minor quantitative differences between simulation and measured data for I(q) which are magnified when we multiply by q2 or q3. F. Local and Terminal Relaxation. Local Dynamics. Local dynamics in the simulated PEO systems was quantified in terms of the torsion autocorrelation function, defined as P(φ(t )) =

⟨cos(φ(t )) cos(φ(0))⟩ − ⟨cos(φ(0))⟩2 ⟨cos(φ(0)) cos(φ(0))⟩ − ⟨cos(φ(0))⟩2

(4)

Results for the decay of the function P(φ(t)) for the two different torsional angles along a PEO chain (C−O−C−C, O−C−C−O) for the three different systems (PEO-5k, PEO10k, and PEO-20k) are reported in Figure S.I.6. The simulation curves can be accurately fitted with a stretched exponential or β KWW function of the form P(φ(t)) = e−(t/tKWW) , where tKWW denotes the characteristic time and β the characteristic stretching exponent. The corresponding correlation time can be calculated from the integral below the P(φ(t)) vs t curve and is equal to

where Γ denotes the gamma function. The values obtained for tKWW, β, and tc for all simulated ring and linear PEO melts are reported in Table S.I.3. The following conclusions can be drawn: (a) for all systems, the values of the stretching exponent β are rather scattered and different for the two angles, i.e., of ϕ1 describing the torsion of bond C−C connecting successive monomers along a PEO chain and of ϕ2 describing the torsion of bond O−C; (b) the correlation times tc characterizing the relaxation of the two dihedrals are practically of the same order of magnitude (a few picoseconds) in all systems; (c) for each

()

Γ tc = t KWW

Figure 7. Comparison of the Kratky plots between MD simulations and SANS measurements for the three PEO rings (5k, 10k, and 20k).

1 β

β

(5) G


Article

Macromolecules

equation with suitable boundary conditions reflecting the absence of chain ends. As a result, odd modes disappear and only even modes contribute to the dynamics:20,48 ∞

∑

R n(t ) = X 0 +

p;even

⎡ ⎛ pπ n ⎞ ⎛ pπ n ⎞⎤ ⎟ + 2Y sin⎜ ⎟ ⎢⎣2X p cos⎜⎝ p ⎠ ⎝ N ⎠⎥⎦ N (6a)

X p(t ) =

1 N

∫0

N

⎛ pπ n ⎞ ⎟ dn R n(t ) cos⎜ ⎝ N ⎠

for p = 0, 2, 4, 6, ... (6b)

Yp(t ) =

1 N

∫0

N

⎛ pπ n ⎞ ⎟ dn R n(t ) sin⎜ ⎝ N ⎠

for p = 2, 4, 6, ...

(6c)

where p;even in the summation represents only the even modes (i.e., p = 0, 2, 4, 6, ...), Rn(t) denotes the position of the nth bead along the chain, Xp(t) is the pth cosine mode for p = 0, 2, 4, 6, ..., and Yp(t) is the pth sine mode for p = 2, 4, 6, .... The normal coordinates satisfy the following autocorrelation functions (for details, see refs 20 and 48): ⟨X pα(t )Xqβ(0)⟩ = δpqδαβ

kBT exp( − t /τp) kp

for p , q = 2, 4, 6, ... (7a)

kT ⟨Ypα(t )Yqβ(0)⟩ = δpqδαβ B exp( − t /τp) kp

for p , q = 2, 4, 6, ... (7b)

⟨Χ pαYqβ⟩ = 0

for all p and q

(7c)

where kp =

6π 2kΒT Nb2

p2

for p = 2, 4, 6, ...

(8)

denotes the magnitude (strength) of the normal mode and τp =

4τ 2Νζ ζ Ν2b2 1 = = 22 2 2 kp 3π kBT p p

for p = 2, 4, 6, ... (9)

the corresponding relaxation time. The last equation shows that the characteristic relaxation spectrum of a ring molecule bears exactly the same form as its linear analogue, the only difference being that while all modes are permitted in the case of linear chains (p = 1, 2, 3, 4, ...) only even modes are allowed in the case of rings (i.e., p = 2, 4, 6, 8, ...). The longest relaxation (Rouse) time for a ring molecule is therefore

Figure 8. Comparison of MD simulation predictions and SANS measurements data for the three PEO rings (5k, 10k, and 20k) plotted as I(q)q3 versus qRg.

τR,ring = τ2 =

dihedral, the differences in the correlation times between ring and linear melts are small. Normal Mode Relaxation. In the Rouse model, a polymer molecule is considered as a sequence of coarse-grained beads, characterized by an effective hydrodynamic friction coefficient ζ, that are connected by Gaussian and massless springs (connectors) with a modulus (or spring constant) K equal to 3kBT/b2, where kB is Boltzmann’s constant, T the temperature, and b2 the average mean-squared distance between adjacent beads at equilibrium. Structural and dynamical correlations between beads in the Rouse model are thus essentially determined by ζ and K only.2,47 For ring melts, the dynamics of such an N-bead Rouse chain under equilibrium conditions (no flow applied) is described by the corresponding Langevin

ζ Ν2b2 12π 2kBT

(10)

which is equal to one-fourth of the corresponding Rouse time for a linear chain comprising the same number of beads N (τR,linear = τ1 = ζN2β2/3π2kBT); therefore, τR,ring = τR,linear/4, which is consistent with the experimental results of refs 17 and 18. The Rouse model for ring polymers presented in refs 20 and 48 was formulated in continuum space using the position vectors of the beads. In the literature, the effect of ring closure on the rheological behavior has also been studied in the framework of the kinetic theory of dilute solutions of bead− spring rings by Wiest et al.49 and Liu and Ö ttinger.50 In particular, Liu and Ö ttinger50 obtained exact expressions for the time constants of a Hookean-spring ring with an arbitrary number of beads N in the presence of hydrodynamic H


Article

Macromolecules

these calculations for several p modes (the case p = 2 corresponds to half of the molecule) for the three PEO ring melts are shown in Figure 10. The figure shows two sets of curves: the first is represented by the symbols and corresponds to the normalized probability distribution functions for the distances Rp computed directly from the simulation trajectories; the second is shown by the continuous curves and corresponds to best fits of the simulation curves with the Gaussian function:

interactions. These showed that the only effect of ring closure was to alter the spectrum of time constants and not the form of the constitutive equation (compared to that of a dilute solution of Rouse chains). In the special case that hydrodynamic interactions are neglected, their results were found to reduce to those of Wiest et al.49 for Hookean-spring rings obtained by a different approach. The long trajectories accumulated in the course of the MD simulations conducted here serve as an excellent starting point for testing the range of chain lengths for which the above version of the Rouse model for ring polymers provides a reliable description of their dynamics. There are several tests that one can think of. The first is that the mean-squared amplitude ⟨Xp(0)2⟩ of the pth Rouse mode scales with N and p as N/p2 (see eqs 7a, 7b, and 8). Thus, if we plot the MD data in the form (p2/N)⟨Xp⟩(0)2 vs p, then, according to the Rouse model, these should collapse on the same horizontal line for all chain lengths N. We show this test in Figure 9. The prevailing

Figure 9. Normalized squared amplitudes of the Rouse normal modes ⟨Xp(0)2⟩p2/N plotted separately for each of the three PEO ring melts simulated here as a function of mode number p.

observation from Figure 9 is that these normalized amplitudes are significantly reduced for small p modes and this effect becomes more pronounced with increasing size of the ring molecule, which implies that longer-range motions are more and more hindered. For the PEO-20k, for example, the amplitude of the lowest p modes is about half that of its highest p modes. Our MD results suggest that we can distinguish between two types of dynamics: unrestricted Rouse motion at short length scales (high-p modes) where the amplitudes of the normal modes follow the Rouse model and restricted dynamics at larger scales (low-p modes) which most likely should be associated with strong ring−ring threading events.29 The effect of restricted dynamics at long length scales observed in our simulations is even more pronounced in the experiments where it relates to the suppression of Rouse modes. We also observe that the regime of Rouse modes p for which the Rouse scaling holds changes as the chain length N increases. We remind the reader that similar conclusions to those drawn here have been reported in the literature for several linear polymers (e.g., PE51 and 1,4-PB52) and other rings.20 We have further computed the probability distribution function of the characteristic molecular distances Rp in a ring molecule that correspond to segments spanning N/p monomers along the molecule. Representative results from

Figure 10. Normalized probability distribution function for the distance Rp between two backbone atoms along a PEO ring molecule corresponding to a segment of size N/p, for several values of p. The continuous lines correspond to best fits of the simulation data with eq 11. Results are shown from the simulations with the three ring PEO chain lengths studied here: (a) PEO-5k, (b) PEO-10k, and (c) PEO20k. I


Article

Macromolecules f (R p) =

4πR p 2 (π ⟨R p 2⟩0 )3/2

⎛ R 2 ⎞ p ⎟ exp⎜⎜ − 2 ⎟ ⎝ ⟨R p ⟩0 ⎠

(11)

The departure of ring conformations from the Gaussian shape is significant for all three rings and increases with p value. In fact, at large p numbers, Figure 10a indicates that the Gaussian distribution cannot anymore fit the simulation data. Similar observations are drawn from the graphs shown in Figure S.I.7 depicting how the values of the squared distances ⟨Rp2⟩0 that provide the best fits of eq 11 to the plots of Figure 10 compare with the values computed directly from the MD trajectories for the mean-squared end-to-end distance of a part of ring molecule of length N/p segments long along a PEO chain. To facilitate the comparison between the various sets of data, in Figure S.I.7 we show the ratio ⟨Rp2⟩/⟨Rp2⟩0 of the two values ⟨Rp2⟩ and ⟨Rp2⟩0 as a function of p for the three PEO ring melts. We see that the best-fit values ⟨Rp2⟩0 come closer to the true ⟨Rp2⟩ values computed by directly averaging the simulation results for ⟨Rp2⟩ over all segments of molecular length equal to N/p as the mode number p increases, but we have to keep in mind that as p increases considerably, the Gaussian function cannot describe at all the simulation data (see Figures 10a and 10b). In particular, the large deviation of the ratio ⟨Rp2⟩/⟨Rp2⟩0 from one for the smallest normal mode (p = 2) observed in Figure S.I.7 (by ∼45−60%) for all three PEO rings simulated here (5k, 10k, and 20k) is a direct manifestation of the nonGaussian behavior of rings at the level of entire molecules. We also observe that for a given p mode the deviation from the Gaussian result (i.e., the deviation of the ratio ⟨Rp2⟩/⟨Rp2⟩0 from one) is larger for the longer PEO melts. A second test of the Rouse model for melts of ring chains refers to the characteristic times describing the decay of the normalized autocorrelation function ⟨Xp(t)·Xp(0)⟩/⟨Xp(0)2⟩ for p = 2, 4, 6, 8, .... According to eq 7a, a log−linear plot of ⟨Xp(t)·Xp(0)⟩/⟨Xp(0)2⟩ versus t for each of the simulated melts should yield a straight line. Numerical results for the relaxation of the first six modes (p = 2, 4, 6, 8, 10, and 12) from our simulations are shown in Figure S.I.8. We see that to a good degree, and except from the very short times, the decay of normal modes seems to be exponential-like, i.e., to exhibit Rouse-like characteristics. As mentioned above, given that the pth mode probes the dynamics of a subchain of N/p segments long, this fair agreement with the Rouse scaling indicates only a mild departure of the dynamics of these subchains from the Gaussian behavior, resembling to a very good degree what is already known in the literature for linear unentangled melts.51−53 It is also very interesting that the deviations from the Rouse scaling are similar for all chain lengths analyzed. In Figure 11 we show the scaling of the characteristic times τp describing the relaxation of the pth mode with node number p and ring length N after fitting the curves of Figure S.I.8 with a single-exponential function (which does not seem to be exactly the case at short times). For all three PEO rings addressed in the present work, our MD data show that to an excellent approximation (with the exception of the smallest modes where the data show some bending toward the low p values)

⎛ N ⎞2.0 τp ∝ ⎜ ⎟ ⎝ p⎠

Figure 11. Scaling of the characteristic times τp describing the relaxation of the Rouse normal modes p with chain length N. According to the Rouse model, τp ∼ (N/p)2.

time τ2 (also known as the Rouse or ring relaxation time, τR,ring) varies with PEO chain length N as τR,ring ∼ N2 exactly as predicted by the ring Rouse theory. This is another indication that the Rouse theory can describe many of the salient features of ring polymers. Overall, our MD simulation data indicate a strong reduction of the normalized mode amplitudes for the smaller mode numbers compared to the Rouse model scaling combined with an undisturbed Rouse relaxation rate. How well such a model can describe experimentally measured data for the dynamic structure factor S(q,t) is discussed in great detail in section G. Normal Mode Cross-Correlation. To test whether or not the Rouse modes are reasonable eigenmodes of the ring PEO chains,54 we computed the amplitudes of the normalized mode correlators Φpq =

X p(0)Xq(0) |X p(0)||Xq(0)|

(13)

for the three PEO melts studied here. To give the full picture, we have plotted these correlators versus the variable x defined as x = (N − 1)(p − 1) + q. The computed Φpq vs x plots are shown in Figure S.I.9 where we have ignored the cases where p = q, since then the correlator is equal to one by construction. For all three PEO systems, all the Φpq values are below 0.15, thereby confirming the orthogonality of the Rouse modes (at least within small deviations). Terminal Relaxation. Terminal relaxation is quantified by monitoring the decay of the orientational autocorrelation function ⟨u(t)·u(0)⟩ of a unit vector u that is a measure of the longest chain dimension. For linear chains, u is taken to be equal to the unit vector along the chain end-to-end vector R. For rings, u is taken to be the unit vector along the ring diameter Rd, namely the vector that connects atoms 1 and N/2 along the ring. Simulation results for all simulated systems are shown in Figure S.I.10. The rate with which ⟨u(t)·u(0)⟩ approaches the zero value is a measure of how fast the chain “forgets” its initial configuration, i.e., of the rate of the overall orientational relaxation of the chain. According to the ring Rouse model, the time autocorrelation function of the ring diameter unit vector is described by the equation

(12)

which is identical to the Rouse scaling (see eq 9). In particular, our MD data indicate that the characteristic longest relaxation J


Article

Macromolecules ⟨u(t ) ·u(0)⟩ =

32 π2

∞

∑ p = 2,6,10,14,...

1 exp( −t /τp) p2

(14)

with τp given by eq 9. Figure S.I.10a presents the simulation data for the three PEO ring systems studied here along with the corresponding Rouse model predictions in a linear−log plot; the latter have been computed by keeping only the first four terms in the summation in eq 14 and by using the values of the characteristic relaxation times τp obtained directly from the normal-mode analysis discussed in the previous paragraphs. We find that for all three PEO melts ⟨u(t)·u(0)⟩ relaxes faster than what is predicted by the Rouse model. Then, to describe terminal relaxation, and similar to the torsional orientational relaxation functions, we resorted to a description of the ⟨u(t)· u(0)⟩ curves through KWW functions of the form described above, namely β0 ⎞ ⎛ ⎛ t ⎞ ⎟ ⎜ ⟨u(t ) ·u(0)⟩ = exp⎜ −⎜ 0 ⎟ ⎟ ⎝ ⎝ t KWW ⎠ ⎠

Figure 12. N dependence of the time t0 quantifying terminal orientational relaxation in the simulated (ring and linear) PEO melts.

for the regime of chain lengths addressed in our atomistic simulations. Given that the entanglement molecular weight of PEO is Me = 2020 g/mol, the scaling t0 ∼ Nb with b = 3.1 for the linear melts is consistent with the predictions of the reptation theory for entangled linear melts.1,2 For rings, on the other hand, the corresponding scaling t0 ∼ N1.9 is consistent with a recent experimental study17 according to which ring polystyrene (PS) melts with Mw smaller than 5Me show the Mw dependence of η0 predicted by the Rouse model for rings whereas ring PS melts with Mw larger than 5Me show a drastic deviation from the Rouse scaling. From a computational point of view, the fast increase of t0 with chain length reveals how difficult it is for MD to relax the conformational characteristics of polymer chains (ring, linear, branched, etc.) as their chain length increases. Given, however, that rings are characterized by relaxation times significantly smaller than their linear counterparts, we understand that atomistic MD simulations (either equilibrium or nonequilibrium) with truly long ring melts are not impossible, and this (atomistic simulation of even longer PEO ring melts) will be the subject of a future study. Mean-Square Displacement of Atomistic Segments. We discuss next the simulation results obtained for the msd of chain segments, g1(t) = ⟨(Rn(t) − Rn(0))2⟩, averaged over all atomistic units n, with time t. Figure S.I.11 shows logarithmic plots of the function g1(t) for the three simulated ring melts (PEO-5k, PEO-10k, and PEO-20k). For PEO-5k and PEO-10k, the g1(t) vs t curves display two distinct breaks that are marked by the arrows in Figure S.I.11, corresponding to three characteristic diffusive regimes. [For the longest system, PEO20k, significantly longer simulation times are needed in order to unambiguously observe the various breaks.] The values of the characteristic exponents quantifying the scaling of g1(t) with t in these three regimes are equal to 0.57, 0.49, and 1.00 for the PEO-5k melt and equal to 0.51, 0.42, and 1.00 for the PEO-10k melt. We recall here that for a linear polymer undergoing reptation the corresponding g1(t) vs t curves exhibit three (and not two) characteristic breaks, with corresponding exponents equal to 0.50, 0.25, 0.50, and 1.00, respectively. Our results for the segmental msd, g1(t), are in excellent agreement with the atomistic MD simulation findings of Hur et al.21 with ring PE melts, who reported that for small ring melts (up to N = 100 carbon atoms per chain) g1(t) shows a crossover from a t0.5 power-law Rouse regime to the free diffusion regime t1. With increasing chain length, however, an intermediate time regime develops, with a slope significantly

(15)

where t0KWW and β0 denote the characteristic KWW time and stretching exponent, respectively. The corresponding correlation time t0 quantifying terminal relaxation is computed from the integral below the ⟨u(t)·u(0)⟩ vs t curve and is given by

() 1 β0

Γ t0 =

0 t KWW

β0

(16)

The values of β0, and t0 for the simulated (linear and ring) PEO melts are reported in Table 2. The dependence of t0 on chain length is further shown graphically in Figure 12. t0KWW,

Table 2. MD Predictions for the Values of the Characteristic Stretching Exponent β0, KWW Relaxation Time t0KWW, and Total Correlation Time t0 Describing the Time Decay of the Function ⟨u(t)·u(t)⟩ Quantifying Terminal Relaxation in the Three PEO Melts (5k, 10k, and 20k) and Comparison between Ring and Linear Melts (T = 413 K, P = 1 atm) terminal relaxation system

β0

t0KWW (ns)

t0 (ns)

ring PEO-5k ring PEO-10k ring PEO-20k linear PEO-5k linear LEO-10k linear PEO-20k

0.73 0.63 0.64 0.53 0.40 0.36

12 38 126 46 251 1320

17 40 200 85 653 6290

Our simulation data (see Table 2) show that (a) the values of the stretching exponent β0 characterizing terminal dynamics are systematically larger for the ring melts, indicative of a smaller degree of cooperativity; (b) the relaxation time t0 is significantly shorter in the ring melts than in the corresponding linear analogues; and (c) for both types of melts (rings and linear), t0 increases rapidly with chain length. More precisely, the data of Figure 12 are in support of the following power law dependence of t0 on N 1.9 ⎧ , for ring molecules ⎪N t0 ∼ ⎨ ⎪ 3.1 ⎩ N , for linear chains

(17) K


Article

Macromolecules smaller than 1/2. For the ring system with N = 1500, Hur et al.21 found the dependence t0.5 reflecting the presence of entanglements. On larger time scales, g1(t) crosses over to a linear dependence corresponding to free diffusion. The g1(t) vs t curves are very important for another reason: they can be used to compare the monomer friction coefficient between ring and linear melts. To this, in Figure 13, we compare the segmental msd in the simulated ring PEO melts with the msd of the most central atoms in the corresponding linear analogues. We have chosen to work only with the most central atoms in linear melts, since end groups are characterized

by a higher mobility due to the excess free volume at the two ends of a linear chain. To make the comparison between linear and ring systems more clear in Figure 13, we show the two msd’s only up to the corresponding Rouse time for each melt. Clearly, for all three different molecular weight PEO melts studied (5k, 10k, and 20k), segments in rings move faster than central atoms in the corresponding linear melts, which demonstrates that the monomer friction factor is smaller in a ring than in the corresponding (equivalent) linear melt. This is consistent with a previous analysis by Tsolou et al.,20 who extracted the segmental friction coefficient ζ for PE indirectly from MD simulation data for the chain self-diffusion coefficient and compared its values between ring and linear PE melts of the same chain length N. It was found that for N > 50 ring PE melts are characterized by smaller ζ values than linear analogues, a direct consequence of the faster dynamics of ring chains compared to linear ones as the chain length increases. Tsolou et al.20 have further reported an additional significant difference between ring and linear PE melts regarding ζ: for linear PE melts and for chain lengths above approximately N = 156, ζ was found to increase considerably which was taken as a direct consequence of the onset of entanglements and of the passage to the reptation-like type of dynamics. In contrast, for ring PE melts, ζ was observed to continue having the same (approximately) value even for chain lengths as long as N = 400 (the longest PE ring melt simulated in ref 20) which was taken as evidence that for rings the Rouse dynamics extends over a significantly broader regime of chain lengths than for linear chains. Loop Dynamics. In a recent publication, Gooßen et al.12 have reconstructed the time evolution of the msd for the 20K ring PEO system from the NSE measurements for the dynamic structure factor. They constructed three curves: the first referred to the time-dependent loop relaxation at short times, the second displayed the dynamics of the loops as derived from the motion of the slowly decaying low p modes, and the third was the sum of the first two which yielded a slope of 0.32 at intermediate times. We check this in Figure S.I.12a showing the time evolution of the segmental msd after having subtracted out the motion of the ring center-of-mass, g2(t) = ⟨(Rn(t) − Rcm(t) − Rn(0) + Rcm(0))2⟩. Our results indicate indeed a slope equal to 0.32 at intermediate times, exactly as reconstructed experimentally from the NSE data. How simulated and reconstructed curves compare with each other is shown in Figure S.I.12b. We observe that while the msd’s reveal a power law regime with a similar exponent, the MD simulations result in segmental msd’s that are significantly larger than those derived from the experimental measurements (more than a factor of 3 in the power-law regime). This may correspond to the much weaker suppression of the lower-indexed modes in the MD simulations discussed above (see Figures 9 and 11). G. Dynamic Structure Factor. A quantity that can be easily computed by the MD simulations and is experimentally accessible by techniques such as dynamic light scattering (DLS) and neutron spin echo (NSE) is the single-chain dynamic structure factor S(q,t). For an isotropic medium (such as a polymer melt), S(q,t) is calculated in an MD simulation through55 S(q , t ) =

Figure 13. Comparison between the segmental msd in the PEO ring melts studied here (PEO-5k, PEO-10k, and PEO-20k) and the msd of the most central atoms in the corresponding linear analogues.

1

Nch

N

N

∑∑ ∑

N N ∑i =ch1 ∑n = 1 fi , n 2 i = 1 n = 1 m = 1

fi , n fi , m

sin(qR i , nm(t )) qR i , nm(t ) (18)

L


Article

Macromolecules where Nch denotes the total number of chains in the system, q the magnitude of the scattering vector q, f i,n and f i,m the scattering factors of atoms n and m along the same chain i, and Ri,nm(t) the magnitude of the displacement vector Ri,nm(t) = Ri,n(t) − Ri,m(t) between atoms n and m on chain i, and the brackets denote a configurational average. What is actually measured in NSE experiments on melts containing a few labeled (usually protonated) chains in a system of chemically equal but deuterated chains is the normalized single-chain dynamic structure factor S(q,t)/S(q,0). Computing S(q,t)/ S(q,0) therefore by recording the time correlation function Ri,n(t) − Ri,m(0) provides a means for comparing MD simulation results to experimental measurements to gain additional insight into the polymer melt dynamics. Parts a−c of Figure 14 present results for S(q,t)/S(q,0) for various values of the scattering vector q as obtained directly from the MD simulations (open symbols) by using eq 18 for all systems studied in this work. Also reported in Figure 14 (full symbols) are the experimentally measured NSE spectra.12,13 For all systems studied and all wavelengths analyzed, the agreement between MD simulations and experimental data is remarkable. The comparison with the Rouse model confirmed the validity of the normal-mode analysis, with the low-q dynamics representing basically translational diffusion. The MD data describe the dynamics at low q values reasonably well despite some slight quantitative deviations in the actual values of S(q,t)/S(q,0). The largest deviations between simulation data and experimentally measured S(q,t)/S(q,0) spectra are observed at intermediate q values. In this regime, the shape of the predicted S(q,t)/S(q,0) curves differs from that of the NSE measurements: at short times the experimental data fall below the simulation curves while at longer times the experimental lines seem to cross the simulation curves. The general impression from the curves shown in Figure 14, however, is that the agreement between predicted and actually measured S(q,t)/S(q,0) spectra is remarkable. The agreement between simulated and experimentally measured S(q,t)/S(q,0) spectra deserves some further investigation. In the NSE measurements, the fast decay of the S(q,t)/ S(q,0) curves at very short times is considered as a signature of loop dynamics; this ceases at longer times, since it is replaced by a slower process, that of loop migration. At high q values, the decay is very rapid and the long-time loop migration is not easily distinguishable. Our MD analysis suggests a different interpretation. From Figures 9 and 11, we recall that the low-p Rouse modes are reduced in their amplitude while the relaxation times are unchanged with respect to the ring Rouse model. In order to compare in detail the MD results with the NSE data for S(q,t)/S(q,0), we made fits of both with the same type of mode amplitude modified ring Rouse model, combined with sublinear center-of-mass diffusion. The long-time Fickian diffusion with ⟨r2(t)⟩F = 6Dcmt is obtained by PFG-NMR and plays only a minor role at the long time end of the spectra; the

Figure 14. Comparison of the normalized single-chain dynamic structure factor S(q,t)/S(q,0) as obtained from the MD simulations (open symbols) for various q values and as measured experimentally by NSE spectroscopy (full symbols) for the three ring PEO melts studied here: (a) PEO-5k, (b) PEO-10k, and (c) PEO-20k.

f u n c t i o n i s g i v e n ( G a u s s i a n a p p r o x i m a t i o n) b y 1 S(q , t )diff,cm = exp⎡⎣ − 6 q2⟨r 2(t )⟩⎤⎦. In the context of the ring Rouse model, the contribution of the Rouse modes to the dynamic structure factor is given by20

v

( ) . In

sublinear diffusion is described by ⟨r 2(t )⟩sub = r0 2

t tc

order to use it in a combined fit of the Rouse model with these types of center-of-mass diffusion regimes, a (smooth) transition is imposed by letting ⟨r2(t)⟩ = [(⟨r2(t)⟩F)a + (⟨r2(t)⟩sub)a]1/a where a controls the sharpness of the transition between the two regimes; for the following model fits, the value of a was fixed to 8 (a = 8). The diffusion part of the scattering

Sring(q , t )

1 = Sring(q , 0) N

N

N

∑ ∑ exp(−Φn,m(t )) n=1 m=1

(19)

with M


Article

Macromolecules Φn , m(t ) = Φn , m(0) + × {1 − e

4Nb2q2 6π 2

−t / τp

∞

∑ p;even

}

⎡ p ⎤ 1 cos⎢π (n − m)⎥ ⎣ Ν ⎦ p2 (20)

(note that only even modes are to be summed up). The Rouse model expression for the dynamic structure factor is obtained by substituting eq 20 in eq 19 and multiplying the result with 1 S(q,t)diff,cm = exp⎡⎣ − 6 q2⟨r 2(t )⟩⎤⎦ . A more general expression for the scattering function Sring(q,t) is provided by using for Φn,m(t) the following expression:12

(

)

2 2 1 q | n − m | ( N − | n − m | )b Ν 6 ∞ 2 2 a 4Nb q ⎡ p ⎤ ∑ p /22 cos⎢⎣π (n − m)⎥⎦{1 − e−t /τp} + 2 Ν 6π p

Φn , m(t ) =

p;even

(21)

allowing for possibly modified mode amplitudes ap/2 (the ringRouse default is ap/2 = 1). We may specify the first few normalized amplitudes ap/2 explicitly, or we can impose an additional, more systematic normalized amplitude variation by using a modulation (Fermi) function: ap/2 = [1 + exp(−{p/2 − pmim}/pwidth)]−1. The mode relaxation time is given by τp = 2W[1 − cos(pπ/N)] ≅ Wp2π2/N2 where W follows from the well-established Rouse parameter Wb4 (= 17180 Å4/ns at T = 413 K) for PEO. Typically, if we use this, the only parameters for the modeling of Sring(q,t) that are left are pmin and pwidth. Finally, we have S(q,t) = Sdiff,cm(q,t)Sring(q,t). Using this model for S(q,t), we performed fits of the NSE and MD data for the 20K ring melt. The long time center-ofmass diffusion constant was fixed to the PGF-NMR value D0 = 0.36 Å2/ns for both (the NSE and the simulation) data fits. The covered time range hardly extends to the transition (at 200− 400 ns) between Fickian and sublinear center-of-mass diffusion. The sublinear diffusion has been fitted by adjusting the exponent ν and the prefactor DS = (r02/tcν). The Rouse rate and the corresponding segment length b were fixed to literature values. Using then S(q,t) = Sdiff,cm(q,t)Sring(q,t) as the model function, simultaneous fits for all available q values with a unique parameter set were performed. To arrive at a good data description by the model, it was absolutely necessary to reduce the contribution of the lowest modes, as described in eq 21. The comparisons illustrated in Figure 15 have been obtained using a systematic mode amplitude cutoff with a “Fermi distribution” type of function. The plain Rouse model for rings without mode suppression but using the established molecular friction value (as expressed in terms of Wb4) grossly fails to describe the data as well as the simulation predictions (see Figure S.I.13). Also, fits (see Figure S.I.14) allowing for a variation of the molecular friction coefficient yielded very unsatisfactory results and in addition would infer a factor of 5 or more deviating (larger) friction coefficient. Since this is a local segmental property, it would be highly improbable that such a deviation is caused by ring topology. Figure 16 displays the resulting mode amplitudes for both fits (to the MD and NSE data). In both cases, the lowest modes are heavily suppressed. However, for the simulation results, the full amplitude is in effect after four modes, while for the NSE data the approach to the full mode spectrum is more gradual and affects about 10 modes. The exact transition shape is difficult to infer from the fits, but it is highly significant that a more

Figure 15. NSE data compared to the Rouse model fit with low mode suppression (a) and fits of the same model to the MD simulation results (b).

Figure 16. Mode amplitudes obtained from a fit with systematic low mode suppression. Blue and red dots correspond to the MD simulation and NSE data fits, respectively.

extensive mode suppression is seen in the experiments than in the simulations. Further, the sublinear diffusion observed in the NSE experiments is somewhat different from that obtained by simulation. In both cases, the scattering function at the lowest q values is dominated by diffusion and thus determines the obtained values. We further clarify that quantitatively the lowN


Article

Macromolecules mode suppression amplitude scheme that comes out from the fitting procedure is somewhat different from the suppression scheme detected directly from the MD simulations (reported in Figure 9). This was necessary in order for the above theoretical model, namely S(q,t) = Sdiff,cm(q,t)Sring(q,t), to provide the best fits to the experimental (and not to the MD-based) NSE data. However, the results are internally consistent since in all cases it is found that the amplitude of the smaller modes should be strongly reduced (compared to the value suggested by the ring Rouse theory) in order for the Rouse model to describe the measured NSE data, with the mode amplitude suppression being stronger in the case of the experimental data. The simulation data yield a subdiffusion exponent ν = 0.756 close to theoretical expectation, while the NSE data for the 20K PEO rings yield ν = 0.567. The corresponding msd’s are shown in Figure 17; the transition between Fickian diffusion with the

Figure 17. Center-of-mass msd as inferred from the fitted diffusion model parameters. Solid line: NSE experiment; dashed line: MD simulation. In both cases, the long t-limit was fixed to the PFG-NMR value.

PFG-NMR value and short time subdiffusion occurs at 360 ns for the NSE data and at 230 ns for the simulation parameters. In order to have a defined common fitting procedure for NSE and simulation results, the present approach is slightly different from the data analysis procedure described in ref 12; however, the salient features of the NSE results corroborate the previous results and have further allowed for a direct comparison with the simulation data. The simulation exhibits similar mode suppression as the experimental data but to a lesser extent. This may indicate that there are still some chain interaction effects that are underestimated in the simulation. It would be interesting to further analyze the origin of this in a future study. H. Self-Diffusion Coefficient. The diffusive behavior of the PEO rings is examined in Figure 18 showing the MD results for the temporal evolution of the msd of the rings centers-ofmass, φcm(t) = ⟨(Rcm(t) − Rcm(0))2⟩, for the simulated PEO systems in a log−log plot. All systems display a similar behavior, characterized by the appearance of two characteristic regions. The first shows up at early times where φcm(t) is observed to scale with t as φcm(t) ∼ tb with b in the interval [0.74, 0.80] for all systems. The second is associated with normal diffusion where φcm(t) scales linearly with t. The chain center-of-mass self-diffusion coefficient Dcm is evaluated by the slope of the msd curve in this long time behavior of φcm(t). The time it takes for the passage from the subdiffusive regime to free diffusion is the so-called crossing time tcross; according to Figure 18 and Table 3, this depends strongly on chain size. For

Figure 18. MD simulation predictions for the msd of the chain centerof-mass φcm(t) as a function of time t in a log−log plot, for all ring PEO chain lengths studied in the present work: (a) PEO-5k, (b) PEO10k, and (c) PEO-20k. Also shown in the figure are the scalings of φcm(t) with t in the subdiffusive and diffusive regimes. The selfdiffusion coefficient is extracted from the slope of the msd curve in the Fickian (long-time) regime (slope equal to 1).

comparison, in Table 3 we also report the characteristic crossing times obtained via the NSE data.12 The values reported in Table 3 indicate that the distances traveled by rings to reach free diffusion are larger than the ring dimensions. For example, in PEO-20k, the transition from the subdiffusive regime to Fickian diffusion occurs at ∼1150 ns, and the corresponding msd of the ring center-of-mass is ∼1350 Å2, implying that the ring has traveled a distance larger than its mean radius of gyration. O


Article

Macromolecules Table 3. MD Predictions for the Crossover Time tcross of the Ring Center-of-Mass from the Subdiffusive Regime to Normal Diffusion for the Simulated PEO Ring Melts and Comparison with the Corresponding Estimates Based on the NSE Data12 tcross (ns) system

MD data

exptl data (ref 12)

PEO-5k PEO-10k PEO-20k

∼60 ∼240 ∼1150

∼100 ∼180 ∼1725

An important point to check here is if the scaling of Dcm with N, together with the scaling of Rg with N discussed in section 3D, is consistent with the scaling of the longest ring relaxation time with N discussed in section 3F. A simple scaling argument suggests that this time should scale as N2.0+0.85, i.e., as N2.85 which differs from the scaling t0 ∼ N1.9 discussed in section 3F. However, if we inspect more carefully Figure 12, we will see that as we go from PEO-10k to PEO-20k, t0 grows faster than N1.9 (roughly N2.30), which is not far from the above scaling. In Figure 19 we show representative atomistic snapshots of three ring PEO chains, one from each of the three simulated systems (5k, 10k, and 20k), and their evolution in time. In all cases, the chains appear to be rather flexible with their shape exhibiting large fluctuations in the course of the simulation. They are seen to assume structures that are quite extended in space and which consist of many small and large loops. This picture is in excellent agreement with the findings of ref 12 where it was argued that dynamics in PEO rings of molecular weight equal to 20k is dominated by loop relaxation at short times, by loop migration and global loop motion at intermediate times, and by the translational diffusivity of entire rings at long times (beyond the crossing time tcross). The dependence of the extracted self-diffusion coefficients Dcm on chain length N is shown in Figure 20, together with experimental data.12 Our conclusion from the data sets displayed in Figure 20 is the excellent agreement between

Figure 20. Predicted and experimental self-diffusion coefficients Dcm versus chain length N in a log−log plot for the simulated ring PEO melts. The scaling of Dcm with N is also shown.

simulation predictions and experimental measurements for Dcm; for all PEO melts, the MD results practically coincide with the experimentally measured values. Overall, ring dynamics in PEO melts with Mw > 5K is described quite accurately by a power law of the form Dcm ∼ N−b with b = 2. In contrast, linear PEO melts show a different scaling: for melts with Mw < 5K, Dcm ∼ N−b with b ≈ 1 (i.e., chain motion is approximately Rouse-like), whereas for melts with Mw > 5K, Dcm ∼ N−b with b ≈ 2.1 (i.e., chain motion becomes reptation-like). Very similar scaling exponents have been reported by Hur et al.21 through atomistic MD simulations of ring and linear melts of a different polymer, polyethylene (PE), with carbon atom numbers per chain up to 500 for the linear system and up to 1500 for the rings. I. Zero-Shear Rate Viscosity η0 and Flow Activation Energies Ea. The zero-shear rate viscosity η0 of long-chain systems is a very difficult property to be predicted directly from equilibrium (through the Green−Kubo relation using the time integral of the shear-stress autocorrelation function) or nonequilibrium (through the response to an imposed steady

Figure 19. Atomistic configurations of single ring chains at different time intervals from the simulations with all PEO melts studied here (PEO-5k, PEO-10k, and PEO-20k). P


Article

Macromolecules shear field) molecular dynamics simulations. In this section, we explore a new path, namely the indirect prediction of the viscosity η0 and its dependence on temperature T from the longest ring relaxation time τR,ring (see eq 10) describing the relaxation of the first even mode (p = 2) which is obtained rather straightforwardly and reliably from the MD simulations. Then, to map τR,ring values to η0 values, one can make use of the Rouse theory for rings22,48 according to which η0 =

π 2 NchkBTτR,ring 6 V

cone-and-plate geometry (25 mm diameter−1° angle) and parallel plates with 8 mm diameter and small gap thickness (between 0.4 and 0.6 mm) to test reproducibility. Flow curves in the shear rate range 1−1000 s−1 were performed in order to measure the viscosity at each temperature. At each temperature, a Newtonian behavior for both melts was recovered (Figure S.I.16) in the entire shear rate range. Figure 21 shows the mean viscosity values obtained at each temperature for the two measured PEO systems in the form of

(22)

where Nch is the number of chains in a system of volume V. According to eq 22, in the context of the Rouse model, η0 is determined by the Rouse relaxation time τR,ring of the ring and the density of the melt. Given that the density does not vary much with temperature, eq 22 was used in the present study to carry out an indirect comparison between simulation and experimental data for the viscosity at various temperatures as follows: instead of using eq 22 to compute η0 from the simulation data for τR,ring at various temperatures, we analyzed the temperature variation of τR,ring itself as obtained directly from our MD simulations with the PEO-5k and PEO-10k ring melts at various temperatures to get the activation energy Ea; and this was subsequently compared with the corresponding activation energy values obtained from direct rheological measurements of the viscosity of the two melts at similar temperatures. Simulating the two melts at the very low temperatures accessed in the rheological measurements (T = 363, 373, and 393 K) necessitated very long MD runs, much longer than 1 μs, for the complete relaxation of their conformational properties and the subsequent reliable extraction of their Rouse relaxation times. This can be seen in Figure S.I.15 showing the time autocorrelation function (ACF) of the first even normal mode (p = 2) in the two systems at T = 363, 373, 393, and 413 K. The longest relaxation time at each temperature was then obtained by the integral below each ACF curve, and the results obtained are reported in Table 4.

Figure 21. Temperature dependence of the longest relaxation time τR,ring (from the MD simulations) and of the zero shear rate viscosity η0 (from the rheological measurements). Dashed lines through the data are best fits using an Arrhenius-type of equation (see main text and Supporting Information for details).

an Arrhenius plot. Error bars are not shown along because they are smaller than the size of the symbols. Also, shown in Figure 21 (again in the form of an Arrhenius plot) are the MD simulation predictions for the temperature dependence of τR,ring (= τ2) for the two PEO ring melts (PEO-5k and PEO-10k). The dashed lines in Figure 21 denote the best fits to the data using the Arrhenius-type expression for η0 or τR,ring, η0 or τR,ring = A exp[Ea/RT], with T the absolute temperature, R the universal gas constant, and A and Ea fitting parameters; the parameter Ea is the activation energy. Practically perfect straight lines are obtained in both cases, with the relative slopes providing the values of the corresponding activation energies Ea/R. The fitted values for both samples are shown in Table 5.

Table 4. MD Predictions for the Rouse (Orientational) Relaxation Time τR,ring as a Function of Temperature T for the Simulated PEO-5k and PEO-10k Systems τR,ring (ns) T (K)

PEO-5k

PEO-10k

363 373 393 413

73 60 36 23

280 197 108 81

Table 5. Parameters Obtained from the Linear Regression (Dashed Lines in Figure 20) Quantifying the Temperature Dependence of η0 (Experimental Measurements) and τR,ring (MD Simulations)

The rheological measurements were carried out with a stresscontrolled Physica MCR-702 rheometer (Anton Paar, Austria). PEO samples were loaded into the rheometer plates at the highest temperature in a way to cancel the thermal history. The temperature was then decreased, and flow curves were performed at each temperature in the range [70−110 °C]. Thermal equilibrium was reached by waiting roughly 30 min at each temperature. The lower temperature limit was set by the glass transition temperature and the upper one by the degradation temperature of the sample. Temperature control of ±0.1 °C was achieved with a Peltier system under a nitrogen atmosphere to avoid thermal degradation at high temperature. Because of the limited quantities of samples available, we used a

system

A × 103 (ns)

PEO-5k (MD) PEO-10k (MD) PEO-5k (exptl data) PEO-10k (exptl data)

4.8 16.6

A × 105 (Pa·s)

2.1 12.9

Ea/R ×10−3 (K) 3.5 3.5 3.15 3.02

± ± ± ±

0.1 0.2 0.06 0.01

For PEO-5k, the predictions from the MD results and the experimental data are equal to Ea/R = (3.5 ± 0.1) × 103 K and Ea/R = (3.15 ± 0.06) × 103 K, respectively. For PEO-10k, the activation energy from the MD simulations is Ea/R = (3.5 ± 0.2) × 103 K whereas that extracted experimentally from the Tdependence of η0 is Ea/R = (3.02 ± 0.01) × 103 K. The Q


Article

Macromolecules

experimental data, the long-time translational ring motion is well depicted in the MD simulations although the lower modes are not suppressed enough (at least not as much as observed in the experiments). We also showed that the predicted internal segmental msd is at least a factor of 3 larger than in the experiments. Given that the MD simulations provide an excellent description of the measured NSE data for S(q,t)/ S(q,0), this implies that the subdiffusivity is slower in the MD simulations than what is recorded in the experimental measurements (in some sense, the faster internal motion compensates for the slower overall diffusive motion). The slowing-down of motion at larger scales is found to express itself in the correlator amplitudes which are suppressed for the lower modes compared to higher ones, while the characteristic times are not affected much. A direct comparison of the activation energies extracted from the simulation studies for the temperature dependence of the Rouse time τR,ring of the simulated ring PEO-5k and PEO-10k systems and from direct experimental measurements for the zero shear rate viscosity η0 showed excellent qualitative and acceptable quantitative agreement. The non-negligible differences found between computed (from the MD simulations) and measured (from the rheological measurements) activation energies should be taken as a signal that the dependence of η0 on τR,ring is not exactly linear as suggested by the Rouse model. We will further pursue this issue in the near future by carrying out detailed NEMD simulations in order to compute the entire viscosity curve (the viscosity as a function of applied shear rate) as well as the first and second normal stress differences. The agreement of the MD predictions with the experimental data (from SANS, NSE, and rheological measurements) proves that detailed atomistic simulations with a validated force field can provide a reliable description of ring melt dynamics. This opens up the way to further employ MD simulations in order to study the dynamics of ring PEO chains in a matrix of linear PEO chains of the same or longer length39 or to carry out detailed NEMD simulations to study the flow behavior of PEO ring melts.

predictions for the activation energies from rheometry for both PEO chain lengths are in a very good quantitative agreement within the statistical uncertainty of the calculations, and the same is true for the MD measurements. Moreover, the activation energies measured in this work are in agreement with previous measurements of a PEO ring of a lower molecular weight (∼2000 g/mol) in ref 15 which gave an activation energy Ea of ∼3.3 × 103 K. It is clear that the activation energy does not depend on the polymer molar mass, at least in the range of molecular weights addressed in this work. The activation energy has been demonstrated not to depend on architecture either;15 ring polymers and their linear counterparts have similar activation energies. Moreover, this behavior is independent of chemistry (see Figure S.I.17 showing the temperature-dependent viscosities of polystyrene architectures with different molecular weights). The comparison between MD and experimental rheometry suggests a slightly milder dependence of η0 (compared to τR,ring) on T, implying a different (nonlinear) dependence of η0 on τR,ring and, consequently, on chain length N (compared to the simple linear relationship, eq 22, suggested by the ring Rouse theory). We are currently trying to directly predict the zero shear viscosity η0 of the PEO ring polymers studied here using nonequilibrium molecular dynamics (NEMD) simulations at several shear rates and then extrapolating in the limit of zero shear rates.

4. CONCLUSIONS We have presented results from microseconds long, equilibrium atomistic MD simulations with model PEO melt systems (ring and linear) of molecular weights ranging from ∼5K up to ∼20K. Simulation predictions for the single-chain intermediate scattering factor S(q,t)/S(q,0), the scaling of ⟨Rg2⟩ and of the ring center-of-mass self-diffusion coefficient Dcm with molecular length N, and the crossover time tcross from subdiffusive to diffusive dynamics have been compared directly to recent experimental data and have been found to be in excellent agreement. Our exhaustive computational experiments provide unique insight for the mechanisms that govern the structural and dynamic properties of pure ring melts. For example, by monitoring changes in the conformational features of individual PEO rings in the course of the MD simulation, we have been able to observe the formation of loop structures and their evolution. This behavior is very similar to that described by Gooßen et al.12 for the time-dependent msd of 20k long PEO rings based on dynamic structure factor measurements through NSE. We have also examined how our simulation results for the conformation and dynamics of normal modes compare with the predictions of the Rouse model. Remarkably, the Rouse modes decay exponentially and to a great degree respect the orthogonality condition. However, deviations from the Rouse behavior are observed for the scaling of the normalized squared amplitudes of the normal modes at small p modes and the shape of the distribution function for the molecular distances Rp between ring atoms along a ring molecule spanning N/p segments. Overall, our MD simulations have demonstrated the validity of the Rouse model for the dynamics of the simulated PEO ring melts. The general agreement with the experimental measurements is also very satisfactory, despite certain deviations in important details referring mainly to internal ring motion. Judging from the comparison with the

■

ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.macromol.6b02495. Tables S.I.1−S.I.3 and Figures S.I.1−S.I.17 provide additional results for the mean-squared radius of gyration ⟨Rg2⟩, the parameters of the KWW function describing the time decay of the torsional autocorrelation functions, the intrachain C−O and O−O pair density functions w(r), the intermolecular C−O and O−O pair distribution functions ginter αβ (r), the dihedral angle distribution functions, the time autocorrelation function of dihedral angles, the mean-squared distances ⟨R2p⟩ for several p values, the normalized time autocorrelation functions of several Rouse modes, the mode cross-correlator Φpq, the decay of the time autocorrelation function ⟨u(t)·u(t)⟩ quantifying terminal orientational relaxation, the meansquare displacement g1(t) of atomistic segments, and the comparison of the plain ring Rouse model to NSE and simulation results from the MD simulations with the studied ring and linear PEO melts; they also provide additional experimental data for the rheological response R


Article

Macromolecules

■

(12) Gooßen, S.; Brás, A. R.; Krutyeva, M.; Sharp, M.; Falus, P.; Feoktystov, A.; Gasser, U.; Pyckhout-Hintzen, W.; Wischnewski, A.; Richter, D. Molecular Scale Dynamics of Large Ring Polymers. Phys. Rev. Lett. 2014, 113, 168302. (13) Richter, D.; Gooßen, S.; Wischnewski, A. Topology matters: structure and dynamics of ring polymers. Soft Matter 2015, 11, 8535. (14) Gooßen, S.; Krutyeva, M.; Sharp, M.; Feoktystov, A.; Allgaier, J.; Pyckhout-Hintzen, W.; Wischnewski, A.; Richter, D. Sensing Polymer Chain Dynamics through Ring Topology: A Neutron Spin Echo Study. Phys. Rev. Lett. 2015, 115, 148302. (15) Brás, A. R.; Pasquino, R.; Koukoulas, T.; Tsolou, G.; Holderer, O.; Radulescu, A.; Allgaier, J.; Mavrantzas, V. G.; Pyckhout-Hintzen, W.; Wischnewski, A.; Vlassopoulos, D.; Richter, D. Structure and Dynamics of Polymer Rings by Neutron Scattering: Breakdown of the Rouse Model. Soft Matter 2011, 7, 11169. (16) Pasquino, R.; Vasilakopoulos, T. C.; Jeong, Y. C.; Lee, H.; Rogers, S.; Sakellariou, G.; Allgaier, J.; Takano, A.; Brás, A. R.; Chang, T.; Gooßen, S.; Pyckhout-Hintzen, W.; Wischnewski, A.; Hadjichristidis, N.; Richter, D.; Rubinstein, M.; Vlassopoulos, D. Viscosity of Ring Polymer Melts. ACS Macro Lett. 2013, 2, 874. (17) Doi, Y.; Matsubara, K.; Ohta, Y.; Nakano, T.; Kawaguchi, D.; Takahashi, Y.; Takano, A.; Matsushita, Y. Melt Rheology of Ring Polystyrenes with Ultrahigh Purity. Macromolecules 2015, 48, 3140. (18) Yan, Z.-C.; Costanzo, S.; Jeong, Y.; Chang, T.; Vlassopoulos, D. Linear and Nonlinear Shear Rheology of a Marginally Entangled Ring Polymer. Macromolecules 2016, 49, 1444. (19) Hur, K.; Winkler, R. G.; Yoon, D. Y. Comparison of Ring and Linear Polyethylene from Molecular Dynamics Simulations. Macromolecules 2006, 39, 3975. (20) Tsolou, G.; Stratikis, N.; Baig, C.; Stephanou, P. S.; Mavrantzas, V. G. Melt Structure and Dynamics of Unentangled Polyethylene Rings: Rouse Theory, Atomistic Molecular Dynamics Simulation, and Comparison with the Linear Analogues. Macromolecules 2010, 43, 10692. (21) Hur, K.; Jeong, C.; Winkler, R. G.; Lacevic, N.; Gee, R. H.; Yoon, D. Y. Chain Dynamics of Ring and Linear Polyethylene Melts from Molecular Dynamics Simulations. Macromolecules 2011, 44, 2311. (22) Halverson, J. D.; Lee, W. B.; Grest, G. S.; Grosberg, A. Y.; Kremer, K. Molecular dynamics simulation study of nonconcatenated ring polymers in a melt. I. Statics. J. Chem. Phys. 2011, 134, 204904. (23) Halverson, J. D.; Lee, W. B.; Grest, G. S.; Grosberg, A. Y.; Kremer, K. Molecular dynamics simulation study of nonconcatenated ring polymers in a melt. II. Dynamics. J. Chem. Phys. 2011, 134, 204905. (24) Halverson, J. D.; Grest, G. S.; Grosberg, A. Y.; Kremer, K. Rheology of Ring Polymer Melts: From Linear Contaminants to RingLinear Blends. Phys. Rev. Lett. 2012, 108, 038301. (25) Halverson, J. D.; Kremer, K.; Grosberg, A. Y. Comparing the results of lattice and off-lattice simulations for the melt of nonconcatenated rings. J. Phys. A: Math. Theor. 2013, 46, 065002. (26) Tsalikis, D. G.; Koukoulas, T.; Mavrantzas, V. G. Dynamic, conformational and topological properties of ring−linear poly(ethylene oxide) blends from molecular dynamics simulations. React. Funct. Polym. 2014, 80, 61. (27) Tsalikis, D. G.; Mavrantzas, V. G. Threading of Ring Poly(ethylene oxide) Molecules by Linear Chains in the Melt. ACS Macro Lett. 2014, 3, 763. (28) Lee, E.; Kim, S.; Jung, Y. Slowing Down of Ring Polymer Diffusion Caused by Inter-Ring Threading. Macromol. Rapid Commun. 2015, 36, 1115. (29) Tsalikis, D. G.; Mavrantzas, V. G.; Vlassopoulos, D. Analysis of Slow Modes in Ring Polymers: Threading of Rings Controls LongTime Relaxation. ACS Macro Lett. 2016, 5, 755. (30) Obukhov, S. P.; Rubinstein, M.; Duke, T. Dynamics of a Ring Polymer in a Gel. Phys. Rev. Lett. 1994, 73, 1263. (31) Ge, T.; Panyukov, S.; Rubinstein, M. Self-Similar Conformations and Dynamics in Entangled Melts and Solutions of Nonconcatenated Ring Polymers. Macromolecules 2016, 49, 708.

of ring PEO-10k and the viscosity as a function of the inverse of the absolute temperature for several ring and linear polystyrene (PS) melts (PDF)

AUTHOR INFORMATION

Corresponding Author

*Phone +30-2610-997398; e-mail [email protected] (V.G.M.). ORCID

Vlasis G. Mavrantzas: 0000-0003-3599-0676 Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS Financial support for this study has been provided by the project “ARISTEIA I” (Project: RINGS, Grant number: 2486) that is implemented under the NSRF−“OPERATIONAL PROGRAMME EDUCATION AND LIFELONG LEARNING” and is cofunded by the European Union (European Social Fund) and the Greek State (Ministry of Education and Religious Affairs−Greek General Secretariat for Research and Technology). We also acknowledge support by the Limmat Foundation, Zurich, Switzerland, through the project “Multiscale Simulations of Complex Polymer Systems” (MuSiComPS). The work was supported by computational time granted from the Greek Research & Technology Network (GRNET) in the National HPC facilityARISunder project pr002038. To this, we further feel indebted to Dr. Ioannis Liabotis and Dr. Dimitrios Ntelis from GR-NET, Greece, for their support on several technical aspects of the work.

■

REFERENCES

(1) de Gennes, P. G. Reptation of a Polymer Chain in the Presence of Fixed Obstacles. J. Chem. Phys. 1971, 55, 572. (2) Doi, M.; Edwards, S. F. The Theory of Polymer Dynamics; Oxford University Press: Oxford, 1986. (3) Qin, J.; Milner, S. T. Tube Dynamics Works for Randomly Entangled Rings. Phys. Rev. Lett. 2016, 116, 068307. (4) Cremer, T.; Cremer, C. Chromosome Territories, Nuclear Architecture and Gene Regulation in Mammalian Cells. Nat. Rev. Genet. 2001, 2, 292. (5) Rosa, A.; Everaers, R. Structure and Dynamics of Interphase Chromosomes. PLoS Comput. Biol. 2008, 4, e1000153. (6) Roovers, J. Viscoelastic Properties of Polybutadiene Rings. Macromolecules 1988, 21, 1517. (7) McKenna, G. B.; Hostetter, B. J.; Hadjichristidis, N.; Fetters, L. J.; Plazek, D. J. A Study of the Linear Viscoelastic Properties of Cyclic Polystyrenes Using Creep and Recovery Measurements. Macromolecules 1989, 22, 1834. (8) Lee, H. C.; Lee, H.; Lee, W.; Chang, T.; Roovers, J. Fractionation of Cyclic Polystyrene from Linear Precursor by HPLC at the Chromatographic Critical Condition. Macromolecules 2000, 33, 8119. (9) Kawaguchi, D.; Masuoka, K.; Takano, A.; Tanaka, K.; Nagamura, T.; Torikai, N.; Dalgliesh, R. M.; Langridge, S.; Matsushita, Y. Comparison of Interdiffusion Behavior between Cyclic and Linear Polystyrenes with High Molecular Weights. Macromolecules 2006, 39, 5180. (10) Kapnistos, M.; Lang, M.; Vlassopoulos, D.; Pyckhout-Hintzen, W.; Richter, D.; Cho, D.; Chang, T.; Rubinstein, M. Unexpected power-law stress relaxation of entangled ring polymers. Nat. Mater. 2008, 7, 997. (11) Brás, A. R.; Gooßen, S.; Krutyeva, M.; Radulescu, A.; Farago, B.; Allgaier, J.; Pyckhout-Hintzen, W.; Wischnewski, A.; Richter, D. Compact structure and non-Gaussian dynamics of ring polymer melts. Soft Matter 2014, 10, 3649. S


Article

Macromolecules (32) Fischer, J.; Paschek, D.; Geiger, A.; Sadowski, G. Modeling of Aqueous Poly(oxyethylene) Solutions: 1. Atomistic Simulations. J. Phys. Chem. B 2008, 112, 2388. (33) Fischer, J.; Paschek, D.; Geiger, A.; Sadowski, G. Additions and Corrections in Modeling of Aqueous Poly(oxyethylene) Solutions: 1. Atomistic Simulations. J. Phys. Chem. B 2008, 112, 8849. (34) van Der Spoel, D.; Lindahl, E.; Hess, B.; Groenhof, G.; Mark, A. E.; Berendsen, H. J. C. GROMACS: Fast, Flexible, and Free. J. Comput. Chem. 2005, 26, 1701. (35) Scienomics, MAPS Platform, Version 3.4.2, France. 2015. Available online: http://www.scienomics.com/ (accessed on 30 September 2016). (36) Nosé, S. Constant Temperature Molecular Dynamics Methods. Prog. Theor. Phys. Suppl. 1991, 103, 1. (37) Hoover, G. W. Canonical dynamics: Equilibrium phase-space distributions. Phys. Rev. A: At., Mol., Opt. Phys. 1985, 31, 1695. (38) Parrinello, M.; Rahman, A. Polymorphic transitions in single crystals: A new molecular dynamics method. J. Appl. Phys. 1981, 52, 7182. (39) Papadopoulos, G. D.; Tsalikis, D. G.; Mavrantzas, V. G. Microscopic Dynamics and Topology of Polymer Rings Immersed in a Host Matrix of Longer Linear Polymers: Results from a Detailed Molecular Dynamics Simulation Study and Comparison with Experimental Data. Polymers 2016, 8, 283. (40) Borodin, O.; Douglas, R.; Smith, G. D.; Trouw, F.; Petrucci, S. MD Simulations and Experimental Study of Structure, Dynamics, and Thermodynamics of Poly(ethylene oxide) and Its Oligomers. J. Phys. Chem. B 2003, 107, 6813. (41) Wick, C. D.; Theodorou, D. N. Connectivity-Altering Monte Carlo Simulations of the End Group Effects on Volumetric Properties for Poly(ethylene oxide). Macromolecules 2004, 37, 7026. (42) Brown, S.; Szamel, G. Structure and dynamics of ring polymers. J. Chem. Phys. 1998, 108, 4705. (43) Brown, S.; Szamel, G. Computer simulation study of the structure and dynamics of ring polymers. J. Chem. Phys. 1998, 109, 6184. (44) Mavrantzas, V. G.; Theodorou, D. N. Atomistic Simulation of Polymer Melt Elasticity: Calculation of the Free Energy of an Oriented Polymer Melt. Macromolecules 1998, 31, 6310. (45) Johnson, J. A.; Saboungi, M.-L.; Price, D. L.; Ansell, S.; Russell, T. P.; Halley, J. W.; Nielsen, B. Atomic structure of solid and liquid polyethylene oxide. J. Chem. Phys. 1998, 109, 7005. (46) Annis, B. K.; Borodin, O.; Smith, G. D.; Benmore, C. J.; Soper, A. K.; Londono, J. D. The structure of a poly(ethylene oxide) melt from neutron scattering and molecular dynamics simulations. J. Chem. Phys. 2001, 115, 10998. (47) Bird, R. B.; Curtiss, C. F.; Armstrong, R. C.; Hassager, O. Dynamics of Polymeric Liquids, 2nd ed.; John Wiley & Sons: New York, 1987; Vol. 2. (48) Watanabe, H.; Inoue, T.; Matsumiya, Y. Transient Conformational Change of Bead-Spring Ring Chain during Creep Process. Macromolecules 2006, 39, 5419. (49) Wiest, J. M.; Burdette, S. R.; Liu, T. W.; Bird, R. B. Effect of ring closure on rheological behavior. J. Non-Newtonian Fluid Mech. 1987, 24, 279−295. (50) Liu, T. W.; Ö ttinger, H. C. Bead-spring rings with hydrodynamic interaction. J. Chem. Phys. 1987, 87, 3131−3136. (51) Harmandaris, V. A.; Mavrantzas, V. G.; Theodorou, D. N. Atomistic Molecular Dynamics Simulation of Polydisperse Linear Polyethylene Melts. Macromolecules 1998, 31, 7934. (52) Tsolou, G.; Mavrantzas, V. G.; Theodorou, D. N. Detailed Atomistic Molecular Dynamics Simulation of cis-1,4-Poly(butadiene). Macromolecules 2005, 38, 1478. (53) 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. (54) Brodeck, M.; Alvarez, F.; Moreno, A. J.; Colmenero, J.; Richter, D. Chain Motion in Nonentangled Dynamically Asymmetric Polymer

Blends: Comparison between Atomistic Simulations of PEO/PMMA and a Generic Bead-Spring Model. Macromolecules 2010, 43, 3036. (55) Warren, B. E. X-ray Diffraction; Addison-Wesley: Reading, 1969.

T


Microscopic Structure, Conformation, and Dynamics of Ring and

Recommend Documents