Tubes, Topology, and Polymer Entanglement - American Chemical

Aug 21, 2014 - purely topological approach can be used to compute the tube diameter by ... (A) Polymer chain confined by entanglements to a tube-like...
0 downloads 0 Views 2MB Size
Article pubs.acs.org/Macromolecules

Tubes, Topology, and Polymer Entanglement Jian Qin† and Scott T. Milner*,‡ †

Institute for Molecular Engineering, University of Chicago, Chicago, Illinois 60637, United States Department of Chemical Engineering, The Pennsylvania State University, University Park, Pennsylvania 16802, United States



ABSTRACT: Long polymer chains in concentrated solutions or dense melts interpenetrate extensively. This results in molecular entanglement and topological constraints that severely restrict the motion of polymers and dominate their flow properties. The effect of entanglements on chain motion has long been understood in terms of a confining tube. The tube is characterized by its diameter, which measures the strength of confinement. We show that a purely topological approach can be used to compute the tube diameter by establishing a direct link to the statistics of topological states in simulations of topologically equilibrated ring polymers. We determine how the tube diameter varies with solvent dilution and chain stiffness. We use these results to choose among several phenomenological proposals that relate the tube diameter to the packing length, which sets the scale for close encounters between nearby chain segments. he unique flow properties of melts and solutions of long polymer moleculeslong stress relaxation times, rubberlike response at shorter times, shear thinning, and normal stresses in shear flowultimately arise from the constraints on local motion imposed by the uncrossability of entangling chains surrounding a given chain. To account for entanglement effects on polymer motion, Edwards and de Gennes introduced the notion of a confining tube.1,2 They asserted that the effects of neighboring chain uncrossability on the dynamics of a test chain can be captured by positing a soft confining tube that encapsulates the test chain, inside which the test chain moves freely and beyond which the constraints are felt (Figure 1). Thus, to describe entangled polymer dynamics, one does not deal explicitly with interactions between chains, but only the deformation, motion, and relaxation of a given chain and its tube.3 By carefully examining and combining the various tube relaxation mechanisms, contemporary tube-based theories successfully account for many rheological properties of entangled polymers in both linear and nonlinear flow regimes.4,5 The tube is characterized by two geometrical properties: its diameter a and contour length L. The tube diameter a measures the transverse spatial range explored by monomers. It is a mesoscopic length scale, much larger than a single monomer but much smaller than the size of a polymer coil. For typical flexible polymers, the tube diameter ranges between 30 and 100 Å (33 Å for polyethylene and 79 Å for poly(dimethylsiloxane)6). The tube contour length L is the arc length of tube centerline or primitive path, which grows linearly with the mass of the polymer. The tube and its primitive path have been modeled in a simple way as a flexible random walk of linear steps of size a, with each step consisting of Ne monomers.3 Each step is an entanglement strand; Ne is the chain arc length between entanglements (“entanglement length”, for short). The total

T

© XXXX American Chemical Society

Figure 1. (A) Polymer chain confined by entanglements to a tube-like region of diameter a, spanned by a random walk of entanglement length Ne. A chain confined to a tube may also be regarded as a string of random-walk “blobs”of size a3.3 (B) Tube and Gaussian-shaped bead transverse spreading revealed by applying the isoconfigurational ensemble to a chain of 400 monomers.

contour length L is given by L = (N/Ne)a. The entire chain can be regarded either as a random walk of N monomers, with mean-square end-to-end distance R2 = Nb2 (in which b is the statistical segment length),7 or as a random walk of of N/Ne tube segments, with mean-square end-to-end distance (N/Ne) Received: April 12, 2014 Revised: July 11, 2014

A

dx.doi.org/10.1021/ma500755p | Macromolecules XXXX, XXX, XXX−XXX

Macromolecules

Article

a2. Equating these implies a2 = Neb2, which may be understood as saying that an entanglement strand of Ne monomers is a random walk with size equal to the tube diameter.

volume pervaded by one entanglement strand a3 is a constant (C3), with the characteristic volume of a close contact being p3. Its prediction for the tube diameter in semidilute solutions differs from that of ansatz I and is consistent with the concentration dependence of the rubber-like modulus for entangled theta solutions.13 Ansatz III was recently introduced14 in the study of tube diameter for oriented and stretched polymer melts, relevant to melts in strong aligning flows. It states that the number of close contacts between an entanglement strand and its neighboring chain segments is a constant (C2). Qualitatively, chains stretched out along some direction have a lesser tendency to fill the space surrounding their own monomers, hence more room for strands from other chains nearby, and more contacts with other chains. Thus, we expect Ne and the tube diameter to decrease for chains under tension. Careful scaling arguments reveal that ansatz I, II, and III all differ in their prediction of how strongly chain tension reduces the tube diameter. Ansatz III provides a better description of how the tube diameter decreases with chain tension in oriented melts, but the limited range of simulation data prevents us from completely ruling out the alternatives. The details of these three formulations have been discussed in ref 14, and the key predictions are summarized in Table 1.



WHAT CONTROLS Ne? The tube diameter is a material parameter which depends on local chain properties such as monomer bulkiness and chain flexibility, but not on large-scale chain architecture (length or presence of long branches). Together with the friction factor for local motion, the tube diameter is the key microscopic parameter for contemporary molecular theories of polymer rheology.4 Experimentally, the tube diameter is determined by comparing the measured rheological behavior of well-entangled polymers to predictions of theoretical models. The simplest such model predicts the rubber-like modulus G of wellentangled polymers on intermediate time scales short compared to the stress relaxation time for chains to explore wholly new conformations, by analogy to the shear modulus of a cross-linked rubber, as G = kBT/a3.3 No microscopic theory has been developed to predict the value of a based solely on local properties of chain monomers. Empirically, the most successful approach to correlate the entanglement length with other material properties is the Lin− Noolandi ansatz,8−11 which consistently describes experimental results for the tube diameter in a wide range of polymers.10 The Lin−Noolandi ansatz asserts that the number of chain segments of length Ne cohabiting the volume pervaded by one such strand is a universal constant for all flexible polymers,10,11 i.e. (Ne1/2b)3 Ne Ω

=C

Table 1. Scaling of ne ≡ Ne/Ne,0 with Volume Fraction ϕ in Solutions and for Melts Oriented with Tension F, under Ansatz I, II, and IIIa ansatz I II III

(1)

Here Ne1/2b is the size of one entanglement strand, Ω is the monomer volume, and C is the universal constant that has been found from data to be about 22.4.10 Equation 1 can be written explicitly for Ne, as Ne = C2Ω2/b6, which reveals how Ne varies with local chain dimensions. Bulkier monomers with larger Ω, or more flexible chains with smaller b, increase the entanglement length. In both cases, the reason is packing constraints. To see this, we write eq 1 in another way, replacing Ne with a = Ne1/2b to obtain a = CΩ/b2 ≡ Cp, in which p ≡ Ω/b2 is the packing length. The packing length is the length below which a strand adopting random walk statistics overfills space.11 The packing length can be thought of as setting the size scale of a local encounter between neighboring chain segments. If a chain has bulky monomers (large Ω), or is flexible (small b) so that its strands adopt compact random walk configurations, other chain segments cannot approach very closely. A “packing strand” consisting of Np monomers has a random-walk dimension of p, i.e., p2 = Npb2 or Np = p2/b2, with displaced volume NpΩ = p3. The Lin−Noolandi ansatz amounts to saying that the tube diameter is proportional to the packing length.10 Besides the Lin−Noolandi ansatz (ansatz I), two alternative formulations have been put forward, both of which yield identical predictions for the tube diameter a in polymer melts but make different predictions for how the tube diameter varies with solvent concentration or chain stretching in oriented melts. Ansatz II was proposed in the study of concentration dependence of the shear modulus in semidilute solutions by Rubinstein and Colby12 and Milner.11 It asserts that the number of close contacts between chain segments in the

solution −2

ϕ ϕ−4/3 ϕ−1

F < FT 1 1 1

F > FT −1

F F−1/2 F−2/3

F ≫ FT F F1/2 F−1/3

a

Ne is the entanglement length in solutions or oriented melts and Ne,0 its unperturbed value. The thermal tension FT ≡ 3kBT/a is the force required to deform an entanglement strand. In solutions and for small tension, the corresponding ratio a/a0 scales as ne1/2; for large tension, a/a0 scales as (3ne3 (x)/x)1/2, where 3 (x) ≡ coth(x) − 1/x is the Langevin function, x is the scaled tension FlK/kBT, and lK is the polymer Kuhn length.

The goal of this work is to develop a coherent description of the phenomenological tube and its underlying topological origin, which explains why the tube diameter varies as it does with monomer size, chain stiffness, solvent dilution, and chain stretching. To do this, we exploit a new method for studying the variation of tube diameter, by analyzing the topological states of simulated ring polymers.15 We examine systems with varying stiffness and show that the dependence of tube diameter on stiffness is consistent with predictions of the original Lin−Noolandi ansatz. We also study systems diluted with ideal solvents and measure the dependence of tube diameter on polymer volume fraction ϕ. It is generally expected that the tube diameter increases with added solvent, since dilution decreases the concentration of entanglement constraints. On dimensional grounds, we expect a power-law dependence of the tube diameter on volume fraction, of the form a(ϕ) = a(1)ϕα, with some dilation exponent α.16 However, the value of α has been controversial; α has been variously argued to be −116−20 or −4/3.12,21−23 It is difficult to determine α experimentally, since the tube diameter is determined only indirectly, by adjusting its value to give best agreement between B

dx.doi.org/10.1021/ma500755p | Macromolecules XXXX, XXX, XXX−XXX

Macromolecules

Article

and their tubes around. The amplitudes and relaxation times for these modes depend on dilution, and so relating ICE studies of systems at varying dilution depends on modeling the concentration dependence of these modes.

some tube-dependent rheological model and experimental rheology data, often linear dynamic rheology. Developing a reliable method to determine the tube diameter from simulations enables us to bypass the ambiguity of choosing the appropriate rheological model to compare to data. We can then conduct a stringent test of the different versions of the entanglement ansatz, since (see Table 1) the values of α are quite distinct for the different versions.



TOPOLOGICAL DEFINITION OF Ne To avoid the uncertainties that accompany the various heuristic approaches for visualizing and measuring the tube in simulations, in this work we relate the tube diameter directly to the statistics of topological states for simulated ring polymers. The basic idea is simple: the more entangled a polymer melt is, the more likely it is to be knotted, and the more different knots it can tie. For every additional entanglement strand we add to an entangled melt or solution, that strand can make one more random decision, to wrap this way or that way around some neighboring strand. This multiplies the total number of different random knots tied by roughly a factor of 2, as the new random decision leads to a new pair of knots for every original knot we had before. We again study ring polymers because (unlike a melt of linear chains) their topological states are well-defined and cannot be altered without breaking and rejoining chains or otherwise allowing chains to cross each other. For our simulations, we used a bead−spring model15 similar to that introduced by Kremer and Grest.31 The number density of beads is fixed at 0.7/σ3 (σ is the Lennard-Jones diameter) to represent a melt. The simulation box is cubic, with various boundary conditions imposed (see below). As for our earlier ICE studies, we generate an ensemble of topologically equilibrated ring polymers by implementing a family of ringrebridging Monte Carlo moves, which effectively allows rings to cross.28 From this ensemble, we count how often different topologically distinct states appear. States are topologically distinct if they cannot be continuously deformed to each other without cutting and rejoining chain segments, i.e., without chains crossing. To identify the topological states of ring polymers, we exploit the elegant and powerful mathematical theory of knot invariant polynomials. Each ring configuration can be regarded as a knot, and a characteristic polynomial can be computed for each knot, which is invariant with respect to all deformations that do not involve cutting and reconnecting. These deformations can be built from elementary “Reidemeister moves”, which correspond to untwisting a simple overhand twist, sliding a simple loop out from under a crossing strand, or sliding a strand over a crossing of two other strands. The invariant polynomial serves as a (nearly unique) label for the knot and can be used to identify the ring configurations (see Appendix for sample knots and their polynomials). Historically, several knot invariant polynomials have been developed, including the Alexander, Jones, and HOMFLY polynomials.34,35 We use the Jones polynomial since it is known to be an excellent knot identifier and, more importantly, since it appears to be the only knot invariant polynomial that has been generalized to periodic knots.36 Periodic knots result from simulating ring polymers in periodic boundary conditions (aperiodic boundary conditions, with impenetrable walls, gives ordinary knots). Just as for many other simulation studies, using periodic boundary conditions for knotted rings helps to minimize finite size effects when collecting statistics of topological states. Confined rings are less entangled within a tube diameter or so from an impenetrable wall, because of the absence of chain segments on the other side of the wall with



VISUALIZING ENTANGLEMENTS Previous simulation work on the tube diameter relies on a family of chain-shrinking algorithms.24−27 The algorithms convert equilibrated polymer configurations into networks of piecewise linear segments between binary contact points, by either energy minimization24,27 or geometrical length reduction,25,26 subject to the constraints that the chain ends are anchored and chain segments cannot cross. The resulting piecewise linear strands are interpreted as the tube primitive path, and the mean distance between successive contact points along a chain is taken to be the entanglement length. These methods capture the essence of chain uncrossability but are somewhat unphysical in that the system configurations are distorted in an ad hoc manner, and the resulting primitive paths contain sharp corners. We have introduced a noninvasive isoconfigurational ensemble (ICE) averaging method for investigating the properties of the tube.28 In this approach, tube confinement is revealed by examining an ensemble of molecular dynamics trajectories derived from the same equilibrated configuration with different random initial velocities. The collection of monomer positions at different times for different trajectories gives rise to a monomer “cloud” manifesting the shape and extent of the tube. The cloud center indicates the tube contour, and the transverse spreading gives the tube diameter. Using this approach, we investigated the effects of free surface on tube diameter,28 the effects of stretch and compression on the tube diameter in entangled networks,29 the effects of chain orientation on tube diameter,14 and the effects of tube semiflexibility on chain relaxation dynamics.30 In all these investigations, we studied melts of one or more long ring polymers, topologically equilibrated by allowing rings to occasionally cross each other, as a proxy for a melt of long entangled linear chains. For simulation efficiency, we employ a simple well-studied bead−spring model for polymers, with harmonic springs of nonzero rest length between bonded monomer “beads”, and repulsive Lennard-Jones interactions between beads.31,15 Long rings offer advantages for studying the tube because in simulations of rings (with the chain-crossing moves turned off), entanglements and hence the tube are permanent. Similar methods have been described by Read32 and Likhtman,33 who simulated linear chains stiffened to decrease Ne, used time averaging of chain conformations without isoconfigurational averaging to determine primitive paths, and investigated the flexibility of the tube32 and the survival probability of entanglement points.33 Unfortunately, the ICE method suffers from a shortcoming that renders it less than ideal for precise quantitative measurements of how the tube diameter varies with concentration. Namely, the tube diameter depends to some extent on the duration of the multiple simulations used to generate the monomer cloud, in part because the entangled melt has thermally excited, longwavelength, slowly relaxing “gel modes” that move the chains C

dx.doi.org/10.1021/ma500755p | Macromolecules XXXX, XXX, XXX−XXX

Macromolecules

Article

The observed power-law distribution of knot frequencies for our largest rings is qualitatively consistent with the expected asymptotic result for the topological entropy St(N). Suppose we assume a power law for the rank-ordered knot probability p(k) of the form p(k) = (β − 1)k−β, with β some function of ring length N. Now we estimate St from St = −kB∫ ∞ 1 dk p(k) log p(k), to obtain St/kB = β/(β − 1) − log(β − 1). This is consistent with the asymptotic limit St = (3/2)NkB if the exponent β for large N approaches 1 + 2/(3N). For our largest rings, the observed exponents β are indeed approaching unity from above. (We emphasize however that in computing St from our histograms, we evaluated the sum −∑kp(k) log p(k) explicitly, without approximating it as a continuous power law.) The behavior of the individual knot probabilities {pi} with ring length N is also interesting (see Figure 3). For any specific

which the ring would have otherwise been entwined. With periodic boundary conditions, a simulated melt of one or more rings can entangle across the periodic boundary with the periodic image of the system. The details for generalizing the Jones polynomial, and an efficient algorithm to compute it for both aperiodic and periodic knots, are described in our previous work.15 We use the Jones polynomials to label the topological states of a single topologically equilibrated ring, simulated in aperiodic, 1D periodic, or 2D periodic boundary conditions. In the nonperiodic directions, a repulsive wall potential is used.15 We have not studied the topological states of rings in 3D periodic boundary conditions because our method for computing the extended Jones polynomial requires a projection of the ring configuration along some axis to obtain a planar “crossing diagram”; it is not evident how to extend this approach to 3D periodic boundary conditions (see Appendix and ref 15). We interpret the set of distinct Jones polynomial values and their frequency of occurrence as the histogram for topologically distinct states, from which we calculate the probability of the system visiting the ith topological state pi. Figure 2 displays our

Figure 3. A log−linear plot of knot probability versus length of rings in aperiodic boundary conditions, for the simplest and most common knots in the “periodic table” of knots (31 (overhand), 41 (figure eight), 51,52,61, and the composite knot (31)2, from top to bottom).

knot k, the probability of finding that knot at first increases as the ring becomes long enough to tie a knot as complicated as k and then decreases as the ring becomes long enough that it almost always ties more complicated than k. The exception to this rule is the “unknot” probability p0, i.e., the probability of not being knotted, which decreases monotonically to zero as N increases, and the ring becomes long enough to be almost always knotted. We observe that the decrease beyond the peak of pk(N) for some given k is roughly exponential in N. This is consistent with the intuitive argument that for any given knot there is a constant probability per unit length of “mistying” the knot by making a mistake. Unlike previous work on polymer knots, which focused for example on the coil dimensions of a knot of a given type,37 we emphasize that here we average over all chain conformations consistent with melt density and the system boundary conditions, as our ring is topologically equilibrated by allowing chain crossings, which can change the knot. In our simulations, the ring conformation (and in particular its fluctuating gyration radius) will influence the probability of tying various knots. Our simulation results for the unknot probability and the topological entropy as a function of ring length N are plotted in Figure 4A,B. The three curves in the figures correspond to aperiodic, 1D periodic, and 2D periodic boundary conditions. Figure 4A shows that the unknot probability indeed decreases with N, with short rings almost always unknotted, and long rings almost always knotted. For intermediate ring length N there is a well-defined transition region, which shifts to the smaller N as the number of periodic boundary directions increases. This result is consistent with the idea that nonperiodic boundaries reduce the entanglement of nearby chain segments, which was our motivation for using periodic

Figure 2. A log−log plot of rank-ordered knot probability versus rank order, for rings of different length (N = 100, 200, 400, 600, and 800 from left to right) in aperiodic boundary conditions. “Stair steps” at large rank order correspond to knots observed a small number (1, 2, 3, ...) of times.

results for the probability of occurrence P of different knots, plotted log−log versus the rank order of the observed knot (i.e., the knot probabilities are sorted, and plotted in decreasing order), for rings of increasing size (N = 100, 200, 400, 600, and 800 from left to right). We note that the distributions tend to a power law form, with an exponent that varies with ring length. From the set of knot probabilities {pi}, we can compute the topological entropy St, defined by St = −kB∑ipi log pi, in which the sum is with respect to all topologically distinct states. The entropy St measures the number of topologically distinct states available to the system; St is extensive, so it should grow linearly with the system size N for sufficiently large N. In ref 15, it was shown that for large N the topological entropy approaches (3kB/2)(N/Ne)the number of topological states grows by 3kB/2 for each additional entanglement strand (consistent with our qualitative argument above). Physically the factor 3kB/2 measures the configuration entropy loss associated with confining one entanglement strand into a tube segment; however, this is a matter of bookkeeping, since a topologically equilibrated melt consists of ideal random walks, unhindered by entanglement constraints or tubes. This entropy loss on confinement of the chain to a tube is exactly compensated by the entropy of topologically distinct tube configurations. D

dx.doi.org/10.1021/ma500755p | Macromolecules XXXX, XXX, XXX−XXX

Macromolecules

Article

Figure 4. Length dependence of (A) unknot probability and (B) topological entropy for ring polymer simulations using a bead−spring model described in ref 15 for aperiodic (blue), 1D (red), and 2D (yellow) boundary conditions. Solid curves in (A) are phenomenological fits (see text). Dashed lines in (B) are linear fits to terminal slopes of simulation results.

Figure 5. (A) Length dependence of the unknot probability p0 from 2D periodic simulations of stiffened and diluted systems. Solid curves are tanh fits; inset is master curve. (B) Shift factors to obtain master curve from unknot probability at different solvent concentrations. Slopes −4/3 and −1 are predictions of ansatz II and III, respectively.

In general, we expect the most important corrections to the asymptotic linear scaling of St(N) for long chains to arise from the effects of impenetrable walls. In a transversely infinite slab of finite thickness made from a melt of long entangled topologically equilibrated rings, the chains will of course be confined by the hard top and bottom walls. The resulting conformations will be ideal random walks that reflect at the walls; the effects of the walls on local chain conformations do not extend beyond a distance of order the persistence length from the walls. This confinement has no effect on the degree of entanglement in the bulk of the slab because it does not matter which segments are locally entangled. Intuitively, we expect the “entangledness” of the slabNe and the rubberlike modulus to be independent of chain length for chains much longer than Ne. More precisely, we expect St/N to become independent of N. However, we have shown previously that a slab is “less entangled” near a repulsive wall because there are no chains on the other side of the wall with which to entangle. This effect extends of order a tube diameter from the wall28 and gives rise to the finite-size corrections which we minimize by adopting 2d periodic boundary conditions, as previously discussed. For our simulations of a ring of N = 800 monomer beads, the box dimension is L = 10.5σ; the tube diameter as determined by ICE analysis is about 4.9σ. Thus, our simulations are marginally large enough that surface effects should begin to be small, in the sense that the system is about two tube diameters across. It is worth noting that limited as they are, our simulation results look at a very large number of distinct knots. The entropy is as large as 12, which if all knots appeared with equal probability would mean there were e12 = 1.63 × 105 different knots present. In fact, the more commonly observed knots are vastly more common than the rarely observed knots, which makes getting good statistics on the rarer knots very

boundary conditions: for periodic boundaries, the rings can become entangled with their periodic images, whereas near repulsive walls, chain segments lose half of their neighbors to be knotted with. The observed unknot probability is well described by a tanh function p0 (log N ) = (1/2)( −tanh(log(N /N0)/w))

(2)

(solid curves in Figures 4A and 5). For single rings in free space, the unknot probability is well-known to decrease exponentially.38 In contrast, our results for the unknot probability of a ring in melt conditions entangled with itself (and its periodic images when periodic boundary conditions are applied) are not well described by an exponential. Indeed, single rings in free space have been heavily studied; for example, swelling due to self-avoidance causes single rings in free space to be more likely unknotted, as one might expect.39 Whereas, because our rings are in melt conditions, selfavoidance is screened away; our ring conformations are ideal random walks, confined of course by the boundary conditions to reflect from the hard walls. Figure 4B shows that the topological entropy St increases monotonically with N for all three boundary conditions and exhibits a substantial finite size effect, similar to that for the unknot probability. We have argued that St should increase linearly with N for large N, approaching St = (3kB/2)(N/Ne) independent of boundary conditions, so that the terminal slopes dSt/dN for different boundary conditions should all approach 3kB/(2Ne). The aperiodic simulation data clearly have not reached this asymptotic regime. If we fit the terminal data in the 1D and 2D boundary conditions to straight lines, we obtain Ne = 73 and 65, which values are consistent with Ne values obtained from the chain shrinking algorithm. E

dx.doi.org/10.1021/ma500755p | Macromolecules XXXX, XXX, XXX−XXX

Macromolecules

Article

challenging indeed. The histogram for the largest N (N = 1600, aperiodic boundary conditions) contained over 12 million entries, with over 4 million distinct knots observed. Given the considerable difficulty of reaching the linear regime for the entropy data, it would be prohibitive to obtain the variation in Ne with dilution or other parameters from direct observation of St. Finally, we should point out a limitation of our method for obtaining St: even the extended Jones polynomial is not a perfect knot identifier. It has been proven that the Jones polynomial never misidentifies two identical knots as distinct but sometimes misidentifies two distinct knots as identical. This becomes more prevalent for more complicated knots, which would contribute increasingly for large N. The quantitative effect of this limitationstatistically, how often the polynomial makes a mistakeis not known; the first mistake is for two knots of 10 crossings. We speculate that this problem is made less severe by using periodic boundary conditions because the generalized Jones polynomial has substantially more degrees of freedom. In any case, our results for St should be understood as lower bounds to the true topological entropy. The results for St(N) for aperiodic systems converge more slowly to the asymptotic scaling that for systems with periodic boundary conditions because the boundary effects are larger (six unentangled sides rather than two for 2d PBC). It may also be that our estimate of the entropy fails at large N for the aperiodic system because too many knots are misidentified as being the same. Intuitively, this should be less of a problem for 2d PBC, for which there are many more different knots and knot polynomials of a given complexity (minimal number of crossings needed to draw the knot), because of the additional toroidal winding numbers.15 Thus, our knot identifier “knows more different knots” and can perhaps provide a better lower bound estimate of St for larger N.

consistent with the value obtained from the chain shrinking algorithm (≃8015). Without a theoretical prediction for the shape of in terms of Ne, it is not obvious what is the best approach to obtain a definite numerical value for Ne. However, as implied by eq 3, we may use Ne to scale the shape of p0(N), requiring that the data for different values of parameters such as chain stiffness or solvent concentration overlap to form a master curve. This does not tell us directly the value of Ne but can be used to determine how Ne varies with other physical parameters. This approach is sufficient to determine the dilution scaling exponent α, as we now show. To examine the effects of solvent dilution and chain stiffness on entanglement length, we performed two sets of simulations. For the dilution effect, we introduced explicit ideal solvent beads, identical in their nonbonded Lennard-Jones interactions to the bonded beads on the polymers. At sufficient dilution, this would result in a good solution of ring polymers, in which the polymers would exhibit swelling from self-avoidance. However, we studied three polymer volume fractions ϕ = 0.71, 0.50, and 0.33, which are much higher than the typical crossover ϕ‡ = 0.1 for the Edwards or concentrated regime (Table 7.2 and eq 7.26 in ref 40). Thus, the excluded volume interactions are screened, and the rings adopt ideal random-walk configurations. To study the effect of chain stiffness, we add to the simulation Hamiltonian an angle potential Ha, of the form Ha(θ) = κ(1− cos θ), where θ is the angle spanned by the two neighboring bonds and κ is the stiffness parameter. Increasing κ stiffens the polymer. In this work, we used κ = 0.75 and κ = 1.5 (the unperturbed system has κ = 0). The statistical segment lengths of the resulting models are determined by observing the length dependence of the radius of gyration for both ring and linear polymers in 3D periodic boundary conditions. The unknot probability data under 2D periodic boundary conditions from both sets of simulations are shown in Figure 5A. For the diluted systems, Figure 5A shows that the unknot probability curves are shifted progressively to larger N as ϕ decreases, implying that the entanglement length increases as we dilute. To quantify the fractional change in Ne, we use eq 3. We fit the unknot probability p0(N) for each ϕ value by a hyperbolic tangent (solid curves in Figure 5A). Then we shift the unknot probability results for diluted systems horizontally (dividing the N values by a shift factor N*(ϕ) for each ϕ) to overlap all the data onto a master curve (Figure 5A inset). The shift factor N*(ϕ) equals Ne(ϕ)/Ne(1.0), displayed in Figure 5B, in which the lines of slope −4/3 and −1 are predictions from ansatz II and III, respectively (ansatz I predicts a much stronger slope of −2). The tube dilation exponent in the Edwards regime is consistent with −1, in support of ansatz III, consistent with our previous results on the effect of chain orientation14 and with estimates obtained by applying the chain-shrinking algorithm to diluted systems.41 In future work, we may explore a greater range of concentration dependence by using an implicit solvent model for our simulations, in which we dispense with the explicit solvent beads, and adjust the range of the attractive LennardJones tail to control solvent quality. In this way the theta condition can be reached,42 and scaling predictions about the concentration dependence of Ne in theta solutions11,12 may be tested. With implicit solvent, the computational demands on our simulations and knot counting would likely be manageable. We would dilute the system to a roughly constant number of



EXPLOITING THE UNKNOT PROBABILITY In contrast to the topological entropy, we can with modest effort obtain very accurate and reliable results for the unknot probability p0, since no known knot shares the same Jones polynomial as an unknot (whether or not the Jones polynomial is truly a perfect unknot identifier is an open problem). We expect that any statistical property regarding the topological states, in particular pi, should be a function of N/Ne alone, because formation of entanglements is how polymers change topology, and the number of entanglement strands in the system is N/Ne. That is, we postulate that pi only depends on parameters such as the monomer size, chain stiffness, and solvent concentration through their effect on Ne: pi = pi ̃ (N /Ne)

(3)

Here p̃i are some universal functions of the ratio N/Ne. In particular, we assert that the unknot probability p0 depends only on N/Ne, with no separate dependence on solvent concentration or chain stiffness. The location of the inflection point in p0 as a function of log N provides an approximate value of the entanglement length. Extrapolating our results for rings in aperiodic, 1D, and 2D periodic boundary conditions, we estimate the inflection point in p0 for 3D periodic boundary conditions at Ne = 90, which is again F

dx.doi.org/10.1021/ma500755p | Macromolecules XXXX, XXX, XXX−XXX

Macromolecules

Article

entanglement strands; assuming Ne(ϕ) scaled as strongly as ϕ−4/3, to reach ϕ = 0.01 would require about 450 times as many polymer beads. But the knot complexities would be about the same (which is the point of diluting at fixed Ne), so the computational cost of counting the knots that dominates our efforts would still be manageable. Likewise, the expected good-solvent scaling behavior of Ne(ϕ), as ϕ1/(1−3ν) or approximately as ϕ−5/4 (with Flory value ν = 3/5 for the swelling exponent), could be investigated with implicit “athermal solvent” (purely repulsive Lennard-Jones), for concentrations sufficiently below the “swelling concentration” ϕ‡, at which chain segments of increasing length exhibit self-avoidance before they overlap.11 For flexible bead−spring chains, ϕ‡ can be estimated to be about 0.25 or so.43 Hence, to observe 1 decade or so of concentration scaling in the good solvent regime, we must dilute to a volume fraction of 0.025, which again would be practical with present methods. Figure 5A shows that the entanglement length decreases with stiffness κ. By shifting the data for different κ onto a master curve, we found that the entanglement lengths for κ = 0.75 and κ = 1.50 decrease by factors of 1.5 and 2.4, respectively, as compared to the values 1.6 and 3.35 predicted by the Lin− Noolandi ansatz (Table 2; Ne ∝ b−6 from all three formulations).

Figure 6. Variation of Ne with chain stiffness, plotted as the ratio /Ne versus a/lK (stiffer chains have smaller Ne and smaller a/lK). N(max) e Squares, present results; circles, results from chain-shrinking methods;24,44 solid curves, Lin−Noolandi scaling prediction; dashed curves, semiflexible scaling prediction.

flexible chains, with a larger ratio a/lK ≈ 8 than for typical “flexible” polymers [typical values of a/lK range from 2.4 (for polyethylene) to 7.6 (for polyisoprene)].14 It seems likely that the stiffest system we considered is so rigid that the assumption of flexible chains implicit in the Lin− Noolandi ansatz begins to fail. Entanglement of “semiflexible” (stiff but bendable) chains, with a much less than lK, is described by a different scaling regime,45,46 for which Ne (obtained from the scaling result for the plateau modulus) should scale as (a/lK)3/5. In Figure 6, the dashed lines indicate this scaling behavior (normalized to match the stiffest chains studied). The simulation results, and indeed the values of a/lK for typical flexible polymers, appear to reside in the crossover region between flexible and semiflexible entangled chains.

Table 2. Shift Factors for Stiffened Systemsa κ

b

(b/b0)6

N*

a/lK

0 0.75 1.50

1.41 1.53 1.73

1 1.60 3.35

1 1.5 2.4

8.0 5.8 3.6



CONCLUSIONS In summary, we have introduced a new, noninvasive, purely topological method for investigating the tube and entanglement length in polymer melts and solutions, exploiting knot invariant polynomials to obtain a purely topological way to determine Ne. This method exploits simulations of long topologically equilibrated rings as a proxy for entangled linear chains. With rings, entanglements are permanent and topological states are well-defined. Our recent work on isoconfigurational ensemble (ICE) averaging also provides a nondestructive way to visualize the tube and measure its properties. ICE averaging in principle can be used to estimate the tube diameter in a wide range of systems but in practice is limited by the fact the results show a weak dependence on the length of averaging time.28 With considerable effort, the topological entropy St can be estimated and found to grow proportional to N/Ne, but in practice is at present too time-consuming to apply broadly. Much simpler is to focus on how the unknot probability vanishes as N grows larger than Ne, when rings first become long enough to typically tie a knot. Our new approach sheds light on longstanding questions about the nature of the tube and entanglements. These questions are evident in the distinctions between the various forms of a phenomenological relation between the tube diameter a and packing length p, which all give the same predictions for unperturbed polymer melts, but different predictions for oriented melts under tension14 and concentrated solutions. Is an entanglement formed when a sufficient number of strands cohabit the same volume (ansatz I)? Or when a volume a3 contains a sufficient number of close contacts between chains

Statistical segment length of unperturbed system with κ = 0 is b0 = 1.41. N* values are obtained by overlap of unknot probability to a master curve.

a

To quantify the chain flexibility relative to entanglement, we compute the ratio of tube diameter to the Kuhn length a/lK. We estimate the tube diameter a from a = Cp (with C = 22.4 the Lin−Noolandi constant) and obtain a value a/lK = 8.0. After the polymers are stiffened, the bond length l0 = 1.0 remains virtually unchanged; the Kuhn length is obtained from lK = b2/l0. The ratio a/lK can be written (CΩl0)b−4. For κ = 0.75 and 1.50, this gives us a/lK to 5.8 and 3.6, respectively. Replacing the variation of b in terms of a/lK, we find that under the Lin−Noolandi ansatz we have Ne scaling as (a/lK)3/2. Chain-shrinking methods have also been employed to determine Ne for stiffened bead−spring chains.24,44 Figure 6 compares our results from topological considerations to the results from primitive path analysis. Plotted are the ratio N(max) e versus the ratio a/lK, where N(max) denotes the largest value of e Ne from the chains with no added stiffness (squares are our results; circles are from chain shrinking). The two sets of results are clearly consistent with each other and not quantitatively consistent with the Ne ∼ (a/lK)3/2 scaling prediction of the Lin−Noolandi ansatz (solid curves). The weaker dependence of Ne on stiffness relative to the Lin−Noolandi scaling exhibited by both the chain shrinking and topological results appears to indicate the onset of a breakdown of the Lin−Noolandi ansatz, which assumes the polymers are flexible enough that entanglement strands are Gaussian random walks. When chains become so stiff that a/lK becomes of order unity, this assumption breaks down. Note that the unstiffened bead−spring model corresponds to quite G

dx.doi.org/10.1021/ma500755p | Macromolecules XXXX, XXX, XXX−XXX

Macromolecules

Article

torus axes; (3) the winding number depends on the choice of the unit cell, which redundancy is eliminated by the prescription that the unit cell be chosen to minimizing the winding number. The generalized Jones polynomial for the chain-mail pattern is V = (−a−4 + 3 − a4)(0,0) + (a2 − a6)(2,0) + (−a−6 − a−2)(0,2), where pairs in bold are winding numbers. The computational time grows exponentially with the number of crossings in the planar diagram. The details of our algorithm are given in ref 15. Previous work on knot theory in polymers is discussed in ref 52.

(ansatz II)? Or when a strand suffers a sufficient number of close contacts with other strands (ansatz III)? These hypotheses differ in a fundamental way, as to whether an entanglement emerges from the collective interaction of many chains (ansatz I and II) or from binary interactions between chains (ansatz III) Our results here and in previous work14 support the third assertion, consistent with a “binary view” of entanglements, as suggested by the binary contacts resulting from chain-shrinking algorithms.47,48 This binary view is also supported by recent simulation results on the local primitive path rearrangements that result when two chain segments in the melt are allowed to cross each other.49 The topological entropy as well can be interpreted in this way: each time a given chain closely approaches another chain, it can wrap around the other, this way or that the given chain then goes on to encounter other chains, each time randomly dodging, resulting in the randomly braided structure we call an entangled melt.



AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected] (S.T.M.). Notes

The authors declare no competing financial interest.

■ ■



ACKNOWLEDGMENTS The authors acknowledge financial support from NSF DMR0907370.

APPENDIX. GIVING KNOTS A GOOD NAME To give knots a “good name” (a nearly unique label), we make use of the Jones polynomial. To compute the Jones polynomial, first project the knot to create a planar “crossing diagram”, which is then oriented by choosing a direction for the contour.34 The unknot, a trefoil knot, and their oriented planar diagrams are shown in Figure 7A,B. Each oriented planar

REFERENCES

(1) Edwards, S. F. Proc. Phys. Soc. 1967, 92, 9−16. (2) de Gennes, P. G. J. Chem. Phys. 1971, 55, 572−579. (3) Doi, M.; Edwards, S. F. The Theory of Polymer Dynamics; Clarendon Press: Oxford, 1986. (4) McLeish, T. C. B. Adv. Phys. 2002, 51, 1379−1527. (5) Graham, R. S.; Likhtman, A. E.; McLeish, T. C. B.; Milner, S. T. J. Rheol. 2003, 47, 1171−1200. (6) Rubinstein, M.; Colby, R. Polymer Physics; Oxford University Press: New York, 2003. (7) Graessley, W. W. Polymeric Liquids & Networks: Structure and Properties; Garland Science: New York, 2003. (8) Lin, Y.-H. Macromolecules 1987, 20, 3080−3083. (9) Kavassalis, T. A.; Noolandi, J. Phys. Rev. Lett. 1987, 59, 2674− 2677. (10) Fetters, L. J.; Lohse, D. J.; Richter, D.; Witten, T. A.; Zirkel, A. Macromolecules 1994, 27, 996−998. (11) Milner, S. T. Macromolecules 2005, 38, 4929−4939. (12) Colby, R. H.; Bubinstein, M. Macromolecules 1990, 23, 2753− 2757. (13) Adam, M.; Delsanti, M. J. Phys. (Paris) 1984, 48, 1513. (14) Qin, J.; Milner, S. T. Macromolecules 2013, 46, 1659−1672. (15) Qin, J.; Milner, S. T. Soft Matter 2011, 7, 10676−10693. (16) Marrucci, G. J. Polym. Sci., Polym. Phys. Ed. 1985, 23, 159−177. (17) Ball, R. C.; McLeish, T. C. B. Macromolecules 1989, 22, 1911− 1913. (18) Larson, R. G. Macromolecules 2001, 34, 4556−4571. (19) Park, S. J.; Larson, R. G. Macromolecules 2004, 37, 597−604. (20) Das, C.; Inkson, N. J.; Read, D. J.; Kelmanson, M. A. J. Rheol. 2006, 50, 207−234. (21) Milner, S. T.; McLeish, T. C. B. Macromolecules 1997, 30, 2159−2166. (22) Lee, J. H.; Archer, L. A. Macromolecules 2002, 35, 6687−6696. (23) Park, S. J.; Larson, R. G. J. Rheol. 2003, 47, 199−211. (24) Everaers, R.; Sukumaran, S. K.; Grest, G. S.; Svaneborg, C.; Sivasubramanian, A.; Kremer, K. Science 2004, 303, 823−826. (25) Kröger, M. Comput. Phys. Commun. 2005, 168, 209−232. (26) Tzoumanekas, C.; Theodorou, D. N. Macromolecules 2006, 39, 4592−4604. (27) Shanbhag, S.; Larson, R. G. Phys. Rev. Lett. 2005, 94, 076001. (28) Bisbee, W.; Qin, J.; Milner, S. T. Macromolecules 2011, 44, 8792−8980. (29) Qin, J.; Milner, S. T. Macromolecules 2012, 45, 9816−9822. (30) Qin, J.; Milner, S. T.; Stephanou, P. S.; Mavrantzas, V. G. J. Rheol. 2012, 56, 702−723.

Figure 7. Simplest knots: unknot and trefoil (A) and their planar diagrams (B). (C) Two-dimensional periodic knots, such as this “chain mail” pattern, can be represented as knots on a torus. (D) The Jones polynomial for a given knot can be related to that of simpler knots by using the skein relation.

diagram is associated with a Jones polynomial V(a) (here a is the polynomial argument, not the tube diameter). The unknot has V(a) = 1. The Jones polynomials for three diagrams differing by just one crossing are related by the skein relation (Figure 7D).34 By iteratively using the skein relation together with operations called Reidemeister moves to simplify diagrams,34 the Jones polynomial for an arbitrary knot can be computed. For the trefoil knot in Figure 7B, this gives V = −a−16 + a−12 + a−4. Figure 7C shows a chain-mail pattern, an example of a 2D periodic knot. To compute its Jones polynomial, the periodic unit cell is projected onto a torus, with outgoing and incoming paths connected across the “seams” of the torus.36,50 Then reductions using the skein relation are performed on the torus. There are three additional complexities for periodic knots: (1) the skein relation alone is not sufficient to simplify the planar diagramthe bracket summation51 is also needed; (2) the Jones polynomial is generalized to account for the “winding number” of the knot, i.e., the number of turns made around the H

dx.doi.org/10.1021/ma500755p | Macromolecules XXXX, XXX, XXX−XXX

Macromolecules

Article

(31) Kremer, K.; Grest, G. S. J. Chem. Phys. 1990, 92, 5057−5086. (32) Read, D. J.; Jagannathan, K.; Likhtman, A. E. Macromolecules 2008, 41, 6843−6853. (33) Likhtman, A. E.; Ponmurugan, M. Macromolecules 2014, 47, 1470−1481. (34) Cromwell, P. Knots and Links; Cambridge University Press: Cambridge, 2004. (35) Sossinsky, A. Knots: Mathematics with a Twist; Harvard University Press: Cambridge, MA, 2002. (36) Grishanov, S.; Meshkov, V.; Omelchenko, A. J. Knot Theory Ramifications 2009, 16, 779−788. (37) Grosberg, A. Y. Polym. Sci., Ser. A 2009, 51, 70−79. (38) Moore, N. T.; Grosberg, A. Y. J. Phys. A 2006, 39, 9081−9092. (39) Koniaris, K.; Muthukumar, M. Phys. Rev. Lett. 1991, 66, 2211− 2214. (40) Graessley, W. W. Polymeric Liquids & Networks: Dynamics and Rheology; Garland Science: New York, 2008. (41) Shanbhag, S.; Larson, R. G. Macromolecules 2006, 39, 2413− 2417. (42) Graessley, M.; Hayward, R.; Grest, G. Macromolecules 1999, 32, 3510−3517. (43) Milner, S. T.; Lacasse, M.-D.; Graessley, W. W. Macromolecules 2009, 42, 876−886. (44) Sukumaran, S. K.; Grest, G. S.; Kremer, K.; Everaers, R. J. Polym. Sci., Part B: Polym. Phys. 2005, 43, 917−933. (45) Semenov, A. N. J. Chem. Soc., Faraday Trans. 2 1986, 82, 317− 329. (46) Morse, D. C. Phys. Rev. E 2001, 63, 031502. (47) Everaers, R.; Kremer, K. Phys. Rev. E 1996, 53, R37−40. (48) Everaers, R. Phys. Rev. E 2012, 86, 022801. (49) Cao, J.; Qin, J.; Milner, S. T. Macromolecules 2014, 47, 2479− 2486. (50) Grishanov, S.; Meshkov, V.; Omelchenko, A. Text. Res. J. 2009, 79, 702−713. (51) Kauffman, L. H. Knots and Physics; World Scientific Publishing: Singapore, 2002. (52) Orlandini, L.; Whittington, S. G. Rev. Mod. Phys. 2007, 79, 611− 642.

I

dx.doi.org/10.1021/ma500755p | Macromolecules XXXX, XXX, XXX−XXX