- Email: [email protected]
Do bending and torsional potentials affect the unraveling dynamics of flexible polymer chains in extensional or shear flows?
Chemical Engineering Science 64 (2009) 4566 -- 4571
Contents lists available at ScienceDirect
Chemical Engineering Science journal homepage: w w w . e l s e v i e r . c o m / l o c a t e / c e s
Do bending and torsional potentials affect the unraveling dynamics of flexible polymer chains in extensional or shear flows? Semant Jain, Indranil Saha Dalal, Nicholas Orichella, Jeremy Shum, Ronald Gary Larson ∗ Department of Chemical Engineering, University of Michigan, Ann Arbor, MI 48109, USA
A R T I C L E
I N F O
Article history: Received 1 October 2008 Received in revised form 16 May 2009 Accepted 22 May 2009 Available online 2 June 2009 Keywords: Computational chemistry Polymers Simulation Rheology Rouse model Brownian dynamics
A B S T R A C T
We compare the flow-induced unraveling of a bead-spring model chain in which each stiff spring corresponds to a single Kuhn step to that of an atomistic representation of the backbone of the same polymer chain with realistic bending and torsional potentials. In our earlier work [Jain, S. and Larson, R.G. 2008. Effect of bending and torsional potentials on high-frequency viscoelasticity of dilute polymer solutions. Macromolecules 41(10), 3692–3700], we showed that in the linear viscoelastic regime, bending and torsional potentials suppress fast local diffusive modes. Now, in strong shear and extensional flows, we observe that bending and torsional potentials have only a slight effect on unraveling dynamics. However, in shear flow, we find that for relatively short coarse-grained chains having less than 20 or so Kuhn steps—representing polymers with fewer than around 150 backbone bonds—coarse-graining introduces periodic peaks in the probability distribution of steady-state stretch. We believe this occurs at high shear because of the dominance of a discrete set of stretch values each of which corresponds to an integral number of bonds that are nearly fully aligned between back-folds that are created during the chain's tumbling orbit. Even this difference between coarse-grained and realistic fine-grained chains disappears for more typical long chain lengths. © 2009 Elsevier Ltd. All rights reserved.
1. Introduction Since atomistic molecular dynamics are still too expensive to apply to long polymer chains (Larson et al., 1999), coarse-grained (CG) bead-spring or bead-rod chains, or constitutive equations derivable from such models are routinely used to determine the dynamics of polymers in flow fields (Larson, 2005). Seminal work by Morton Denn and co-workers (Fisher and Denn, 1977; Petrie and Denn, 1976), starting in the 1970s, has shown the value of this approach. Computer speeds prior to the 1990s limited much of the modeling of that period to dumbbell models (Kuhn, 1934). Nonetheless, many physical phenomena were included in those models including internal viscosity (Cerf, 1969; Manke and Williams, 1992), dissipative stress (Hinch, 1994; Orr and Sridhar, 1996; Rallison, 1997), self-entanglements (Armstrong et al., 1980; James and Sridhar, 1995; King and James, 1983), kinked or yo-yo conformations (Larson, 1990; Rallison and Hinch, 1988; Ryskin, 1987), conformation-dependent friction (de Gennes, 1974; Hinch, 1974), and straining inefficiencies (Dunlap and Leal, 1987; Hinch, 1994; Phan-Thien et al., 1984).
∗ Corresponding author. Tel.: +1 734 936 0772; fax: +1 734 763 0459. E-mail addresses: [email protected] (S. Jain), [email protected] (I.S. Dalal), [email protected] (R.G. Larson). 0009-2509/$ - see front matter © 2009 Elsevier Ltd. All rights reserved. doi:10.1016/j.ces.2009.05.035
In addition, localized effects were captured by using a finer-grained bead-spring chain (Rouse, 1953; Zimm, 1956) that still coarse-grains the polymer chain into a sequence of frictional beads connected by Hookean or non-Hookean springs. The finest-grain model typically used to describe the deformation of long polymers in dilute solution is the Kramer's freely jointed bead-rod chain (Doyle et al., 1997; Larson, 2005; Underhill and Doyle, 2004). This chain is resolved at the level of a Kuhn step wherein the flexibility that is in reality distributed over several backbone bonds of the real polymer is lumped into each free joint. As a consequence, bead-spring chains and even bead-rod chains are not expected to capture very fast dynamics of synthetic polymers like polystyrene (Larson, 2004a). These fast dynamics are controlled by local bond rotations that are influenced by bending and torsional potentials—potentials that are left out of conventional bead-spring and even bead-rod models. While these local chain details do have an effect on the overall coil size, for very long polymers, their effect on chain dynamics at low frequencies or low shear rates is insignificant. In this limit, a number of computational methods have been developed that make similar predictions as long as both hydrodynamic interactions and excluded volume effects are appropriately accounted for. One method is the renormalization group method, which when applied to a long chain in a shearing flow, yields shear thinning that becomes stronger as excluded volume effects increase (Baldwin and Helfand,
S. Jain et al. / Chemical Engineering Science 64 (2009) 4566 -- 4571
1990). A more commonly applied method is to solve for the moments of a Fokker–Planck equation for the chain configuration probability distribution using normal mode decompositions (Fixman, 1966; ¨ Herrchen and Ottinger, 1997; Magda et al., 1988; Ottinger, 1987; ¨ Ottinger, 1989; Rouse, 1953; Wedgewood, 1989). Finally, in more recent years, Brownian dynamics simulations have become conven¨ tional (Hsieh et al., 2003; Kroger et al., 2000; Kumar and Prakash, 2003; Liu, 1989; Lyulin et al., 1999; Somasi et al., 2002; Underhill and Doyle, 2004; Underhill and Doyle, 2007). To obtain universal results both the Fokker–Planck equation and especially the Brownian dynamics method require care to be sure that the asymptotic longchain limit is reached (Kumar and Prakash, 2003). At shorter distance and/or time scales, the local effects such as bending and torsional potentials become important once again. In principle, determination of these effects require using atomistically resolved molecular dynamics simulations that include local details such as bending and torsional potentials and perhaps even interactions between solvent molecules and the polymer chain. While chains of modest length—containing dozens to hundreds of backbone bonds—are not long enough for their properties to converge on the long-chain universal results, they remain too long to be simulated by atomistic molecular dynamic simulations except at great cost. Hence, polymer chains containing a modest number (10–150) of Kuhn step lengths have been simulated by the Kramer's freely jointed bead-rod chain (Larson, 2005) which is resolved at the level of a Kuhn step. What is not clear is if such an approach (Hinch, 1994) can capture very fast dynamics that might be controlled by local bond rotations that are influenced by bending and torsional potentials. Hence, we are interested in determining the accuracy of the freely-jointed chain model in describing the dynamics of modestlength chains at high frequency or high shear rate. Our earlier work (Jain and Larson, 2008) first explored the influence of these bending and torsional potentials on the linear viscoelastic relaxation of fine-grained (FG) chains. We showed that their presence slowed down the local relaxation modes of the chain so that their contributions to viscoelastic relaxation overlapped those of some of the slower collective diffusive modes. This study helped explain why these local modes are not evident in viscoelastic spectra of dilute polymer solutions (Lodge et al., 1982), why coarse-grained bead-spring models are so successful in describing the linear viscoelasticity of dilute polymer solutions over such a wide frequency range, and why there exists an optimal degree of fine-graining in such models. For polystyrene, the optimal fine-graining corresponds to choosing one coarse-grained spring to represent a molecular mass of 5000 Da (Amelar et al., 1991), so that a chain of molecular mass 100,000 Da is represented by 20 such springs. However, our earlier work was limited to dilute polymers in the linear viscoelastic regime. In the present work, we extend the previous methods to strong extensional and shear flows. At increasingly high deformation rates, the portions of the chain that can be assumed to remain well equilibrated become progressively smaller (Larson, 2005). Therefore, local details of the chain—especially the bending and torsional potentials—might interfere with the ability of the chain to deform and unravel in such strong flows. Previously, we also tested the ability of coarse-grained beadspring chains (which lack bending and torsional potentials) to capture the effects of fine-grained chains with torsional and bending potentials by mapping CG chains onto the FG chains and performing Brownian dynamics simulations of both. In order to carry out this mapping, the overall chain length, radius of gyration, and total drag of the two chains were matched. We used Fraenkel springs possessing a preferred spring length for both the CG and FG chains. For the FG chains, the preferred spring length represented the length of a backbone bond which is 1.53 Å for a carbon–carbon bond. The preferred spring length of the CG chain was that of a Kuhn step
4567
length. Since the stiff springs in our CG chains could neither stretch nor compress much, our CG chains were similar to Kramer's freely jointed chains that are commonly used as coarse-grained models of real chains. Thus, the coarse-grained chain had an order of magnitude fewer springs and beads but each spring was longer and each bead exerted more drag than their fine-grained counterparts. This allowed us to assess quantitatively the effect of the bending and torsional potentials on the dynamics of dilute polymers and to determine the degree to which the CG chain captures the dynamics of the FG chain. In what follows, we extend this approach to strong shear and extensional flows using the same FG and CG models used in the earlier study of viscoelastic relaxation. As in the earlier study, we neglect both hydrodynamic interactions and excluded volume interactions between the beads. 2. Simulation 2.1. Force calculations We model the polymer chain as a series of beads connected by stiff Fraenkel springs and use the same nomenclature and force calculations as before (Jain and Larson, 2008) except for the presence of a non-zero solvent velocity field as described below (Jain and Larson, 2008; Larson, 2005) F iD = (V − x˙ i )
(1)
where F iD is the drag force on the ith bead, is the bead frictional drag constant, x˙ i is the velocity of the ith bead, V = ri · ∇V is the imposed solvent flow velocity at the bead's location, and ri is the position vector of the ith bead. In planar extensional flow, this velocity is ⎛ ⎞ ˙ 0 0 ⎝ (2) ∇V ≡ 0 −˙ 0 ⎠ 0 0 0 where ˙ is the extension rate. Similarly, under simple shearing flow, the velocity gradient is ⎛ ⎞ 0 0 0 (3) ∇V ≡ ⎝ ˙ 0 0 ⎠ 0 0 0 where ˙ is the shear rate. Since a CG chain has fewer beads than its equivalent FG chain, to maintain the same total frictional drag coefficient for the whole chain, T , the frictional drag per bead, , of the CG chain must be made proportionally higher than that of the FG chain, i.e., T |CG = T |FG (NS + 1)|FG (4) ⇒ |CG = |FG (NS + 1)|CG where NS is the number of springs in the FG or CG chain. Neglecting inertial forces, to conserve momentum on the ith bead, the drag force, F iD , must be set equal to the sum of the stretching, F iS , bending, F iB , torsional, F iT , and random Brownian, F iR , forces yielding (Jain and Larson, 2008; Larson, 2005) F iD = F iS + F iB + F iT + F iR (5) i i i i ⇒ x˙ i = F S + F B + F T + F R + V In the above equation, the stretching force (stretching constant,
S /m = 2.5×1027 s−2 ) maintains a separation distance of around 1.53 Å between two beads through stiff Fraenkel springs, the bending force (bending constant, /m = 1.3×107 J/kg) maintains the 109◦ 47 tetrahedral bending angle, the torsional force (torsional constant, /m = 6.6344×105 J/kg) retains the desired distribution of trans and gauche conformations and barriers between them through
4568
S. Jain et al. / Chemical Engineering Science 64 (2009) 4566 -- 4571
the Rychaert–Bellman potential (1975), the random Brownian force mimics the random motion of the solvent molecules, and the drag force captures the effect of extensional or shear flow field. In these force constants, m is the mass of the bead. Since we are only concerned here with tracking the chain configurations and the stress, we only need to pick a time step small enough for these properties to converge. We showed earlier that this was accomplished when we used a time step of t = 5×10−16 s, which is one-tenth of the time-step value chosen by Helfand and coworkers (1980)—from whose work we had taken the force constants for our simulations. Similar to our previous work, at this value of the time step, we find that the springs in the present case remain within 5% of their equilibrium length even at the highest extension and shear rates used here. 2.2. Coarse-grained chain To construct a coarse-grained chain whose coil size and fully stretched length are equal to those for a fine-grained polymer chain with bending and torsional potentials, we equate the mean square end-to-end vector, R2 0 , and the fully extended length, L, of the two chains via the relationships R2 0 = C∞ nl2
(6)
L = 0.82nl = NK bK
(7)
where NK is the number of bonds (Kuhn steps), bK is the length of a bond in the equivalent freely jointed chain model, n is the number of backbone bonds, l is the length of each backbone bond, and C∞ is the characteristic ratio of the chain which incorporates the expansion in the chain due to bending and torsional forces. Thus, for a FG chain with bending and torsional potentials typical of real chains (with the torsional potential at 100% “base-case” as defined earlier (Jain and Larson, 2008)), and having 60 springs of equilibrium length, l, of 1.53 Å, we construct an equivalent CG chain having eight springs each with an equilibrium length of 8.9938 Å. The longest relaxation time, 1 , can be computed from the Rouse formula (Doi and Edwards, 1986; Larson, 1988; Larson, 2004b):
1 =
2
8kB T S sin2 [ /2(NS + 1)]
(8)
where the FG bead drag coefficient is = 1.4×1012 kg/mol s (Helfand 2 et al., 1980), and S is the mean square end-to-end distance of 2 B C T 2 = 2 × 2.41(1.53 × a spring computed as 3/2 S = R2 0,S = C∞ ∞ 0 10−10 )2 m2 = 1.1283 × 10−19 m2 (Jain and Larson, 2008) which gives us 1 |FG = 6.426 × 10−9 s. Using Eq. (4) to ensure the same total drag on FG and CG chains, we compute the bead drag coefficient of the CG chain which, along with Eq. (8), results in 1 |CG = 6.865 × 10−9 s. The unraveling of a polymer in a strong flow is dependent on the initial configuration of the chain. To begin a detailed comparison of the unraveling dynamics of CG and FG chains, we first map the initial bead positions of the FG chain onto the CG chain. We would like to choose a corresponding initial configuration of the CG chain that matches as closely as possible that of the FG chain. Since the CG chain has roughly seven times fewer springs than the FG chain, this initial mapping is accomplished by selecting the bead position of the first bead, every seventh subsequent bead, and the last bead of the FG chain as the initial guess bead positions for the CG chain. These bead positions usually yield CG spring lengths that differ considerably from their equilibrium length. So, to begin the simulation at the same roughly equivalent equilibrated state, we then carry out a no-flow simulation with just the CG chain using a very small time step of 10−20 s—a value which is 1/50,000th of the time step used in the flow simulations. Since spring forces at such an non-equilibrated state can
Fig. 1. Mean dimensionless extension, x , as a function of Hencky strain, , at Wi = 2, 5, and 10 for CG and FG chains. The x values of FG and CG chain data at Wi = 10 have been augmented by 10 units to separate these results from the results at Wi = 5.
be large, a time step any larger could lead to instabilities. So, to enable all CG springs to relax to their preferred lengths, this relaxation process is carried out for 2×107 time steps which is ∼3×10−3 times the longest relaxation time of the chain—a time long enough to relax the springs to near-equilibrium length but not long enough for the overall chain conformation to otherwise change appreciably from the original starting configuration of the FG chain. This completes the mapping of a given starting configuration of a FG chain onto a roughly equivalent starting configuration of a CG chain enabling us to begin a direct comparison between the behavior of individual FG and CG chains with similar starting configurations. 3. Results Here we compare the average dynamics using an ensembleaverage over 100 FG chains (each with 60 springs) with an ensembleaverage over 100 corresponding CG chains (each with eight springs) for the extensional flow study. The ensemble was averaged over 50 chains for the shear flow study. We carry out each simulation at various Weissenberg numbers, Wi = ˙ 1 or Wi = ˙ 1 . The dimen2 sionless mean fractional extension, x = X / , is plotted against Hencky strain units ( = ˙ t or = ˙ t). X 2 is the root-mean-square farthest distance between any two beads of a chain in the flow direction and is the equilibrium spring length of the FG chain. 3.1. Extensional flow In Fig. 1, we compare the average dynamics of fine-grained and coarse-grained chains at three Weissenberg numbers. Simulations are run for 9–12 Hencky strain units which translate into 7.5×107 , 2.5×107 , and 1.5×107 time steps when Wi = 2, 5, and 10, respectively. At low Weissenberg numbers, the FG chains stretch a bit less than the CG chains because of the internal rigidity produced by the bending and torsional potentials. This difference diminishes at higher Weissenberg numbers where the large drag forces seemingly overwhelm any dynamic effect of these potentials. The error bars in Fig. 1 show that the standard error decreases at higher Weissenberg numbers. The ensemble-averaged dynamics of unraveling of CG chains are in reasonably good agreement with those of FG chains with the agreement becoming better with increasing Weissenberg number (Wi ⱖ 10). Since we matched the initial configuration of each FG chain to an equivalent CG chain, we can also compare the unraveling behavior of the individual FG and CG chains. In doing so, we find
S. Jain et al. / Chemical Engineering Science 64 (2009) 4566 -- 4571
Fig. 2. Starting configurations of the FG and CG chains in three dimensions. The FG chain is shown as a gray line with closed circles and the CG chain is represented as a black line with closed diamonds.
4569
Fig. 4. Two-dimensional projections (onto the XY plane) of the configurations of the CG and FG chains whose stretch versus time dynamics are shown in Fig. 3. The FG chain is shown as a gray line with filled circles and the CG chain as a black line with filled triangles. The configurations are shown at strains, , of 0, 10, 15, 20, 30, and 40. For clarity, the y-coordinates of the FG and CG chains shown at = 10, 15, 20, 30, and 40 have been shifted up vertically by 3, 6, 9, 12, and 15 units, respectively.
Fig. 3. Mean dimensionless extension, x , as a function of Hencky strain, , for the FG and CG chains with similar initial configurations as shown in Fig. 2.
that while most of the time the CG chain shows similar unraveling behavior as the corresponding FG chain, there are large differences in some cases. One such case for Wi = 10 is shown in Figs. 2–4 where the FG chain unravels slower than the mapped CG chain. Fig. 2 shows the starting configuration of the FG chain and of the CG chain that was mapped onto it. The dimensionless mean extension as a function of strain is given in Fig. 3 and the projections of the FG and CG chain onto the X−Y plane as a function of strain are shown in Fig. 4. Figs. 3 and 4 depict a simulation where the FG chain gets temporarily trapped in a folded state that is avoided by the CG chain. In some other cases, which have not been shown here, the CG chain gets trapped instead and unravels slower than the corresponding FG chain. Thus, we find that one of the chains gets trapped temporarily in a folded state in all cases where the CG chain and FG chain dynamics differ greatly. To determine if the large differences between unraveling of FG and CG chains that occur in the few cases are because of non-linear amplifications of small differences in starting configurations, rather than due to the differences between FG and CG chains, we run pairs of simulations of eight-spring CG chains. The two chains in each pair have very similar but not identical starting configurations. These pairs are created by relaxing a CG chain for 4.00 and 4.01 relaxation times which are then subject to an extensional flow at Wi = 10 using identical Brownian force histories for both chains. Despite the fact
Fig. 5. Unraveling dynamics of two different chains with slightly different initial configurations and identical Brownian force histories plotted as mean dimensionless extension, x , as a function of strain, . The first chain had its starting configuration relaxed up to four relaxation times (4.00, solid black) and the second to 4.01 relaxation times (4.01, solid gray).
that the starting configurations are very similar to each other, and the random forces on the beads are identical, dissimilar behavior is obtained for some pairs of chains. One example of this behavior is shown in Fig. 5 which proves that minor differences in starting configurations can be amplified to produce widely different unraveling behavior in a few cases. Thus, the large differences in some cases between the unraveling dynamics of CG and FG chains can be attributed simply to small differences in starting configurations resulting from the mapping of multiple FG bonds onto a single CG bond and not to intrinsic differences between FG and CG chains themselves. Given the unavoidable small differences in starting conditions, especially at high Weissenberg numbers, the CG model captures the dynamics of FG chains as faithfully as possible. 3.2. Shear flow The mean stretch of the FG and CG chains under shear flow at Wi = 10, 50, and 100, is averaged over 50 runs with the data collected
4570
S. Jain et al. / Chemical Engineering Science 64 (2009) 4566 -- 4571
Table 1 Comparison of mean stretch of fine-grained and coarse-grained chains under shear flow conditions. Simulation no.
1 2 3
Wi
10 50 100
Mean stretch (with standard error) √ √ x¯ |CG ± / N − 1 x¯ |FG ± / N − 1
Error (%)
20.66 ± 0.15 25.51 ± 0.13 27.19 ± 0.13
−7.06 −2.66 1.44
22.12 ± 0.26 26.19 ± 0.23 26.79 ± 0.17
1 − x¯ |CG / x¯ |FG
The mean stretch is computed using 50 chains (N = 50).
Fig. 7. Probability density as a function of mean fractional extension, x /L, at Wi = 100 for a FG chain having 60 springs and CG chains having 8, 16, and 32 springs. The maximum extension, L, of the FG chain is 60×0.82×1.53 Å = 75.276 Å.
the peaks are only a problem for short chains. Chains represented by 32 or greater CG springs—which correspond to 250 or greater backbone bonds—show negligible artificial peaks in stretch distribution even at Wi = 100. 4. Conclusions Fig. 6. Probability density as a function of mean dimensionless extension, x , at Wi = 10, 50, and 100 for CG and FG chains.
in the first 25 strain units not included in the analysis. Since all shear flow simulations are carried out for 5×108 time steps, this translates into ∼400, ∼1000, and ∼1900 strain units at Wi = 10, 50, and 100, respectively. This analysis has been shown in Table 1. Next, in Fig. 6, we compare the probability distribution as a function of mean dimensionless extension (Hur et al., 2000) at various Weissenberg numbers for the CG and FG chains. With 60 fine-grained springs, the maximum mean fractional extension is 60×0.82 = 49.2 which is rounded up to 50 and then divided into 100 equally spaced bins. The factor of 0.82 arises on account of the zig-zag contour of the FG chain because of the bending angles. We notice that the overall shape of the distributions for FG and CG chains are similar but the CG chain distribution contains peaks or oscillations. These oscillations are due to the high degree of orientation often exhibited by the springs at high Weissenberg numbers and the small number (8) of stiff springs. Counting the small peak when x = 10, we observe seven equally spaced peaks which are one less than the eight-coarse grained springs in the chain. Evidently, any peak corresponding to the length of one spring is negligible. In a folded state of well-oriented springs, as the stiff Fraenkel springs each maintain a fixed length, only a discrete number of springs is aligned between folds and the preferred stretch values appear to correspond to these discrete states. To demonstrate that the oscillations are due to the small number of CG springs, we increase the number of springs from eight to 16 and then to 32. As in our earlier simulations, we show the probability distribution averaged over 50 molecules at Wi = 100 in Fig. 7. To retain the same Weissenberg number as the spring length is increased, in accordance with the Rouse theory, the strain rate is reduced as ˙ = Wi/(N|CG )2 where N|CG is the number of springs in the CG chain. While the chain with eight springs has seven peaks, these peaks become more numerous and diminish in magnitude as the number of springs in the chain is increased to 16 and 32. Thus,
We compared the non-linear response of a locally realistic finegrained bead-spring model with stretching, bending, torsional, and drag forces to that of an equivalent coarse-grained chain with only stretching and drag forces. At low extension rates, the CG chain deforms slightly more than the FG chain. This difference disappears at higher extension rates indicating that a CG model should be sufficient when the Weissenberg number is greater than ten. Our individual chain comparisons reveal that even minute differences in the starting configuration of chain can lead to drastically different unraveling behavior. This behavior is especially pronounced if any chain gets trapped during the unraveling process. In shear flow, the CG chain has a similar steady-state distribution of stretch as the FG chain does except that the former exhibits oscillations, or multiple values of stretch, resulting from the small number of springs used in the CG chain. We verify this by using longer chains that correspond to 250 or greater backbone bonds and find that the oscillatory behavior disappears almost entirely. Thus, for chains as long or longer than this, especially at high Weissenberg numbers, coarse-grained chains can accurately capture the fine-grained chain dynamics. Acknowledgments We gratefully acknowledge the financial support provided by Petroleum Research Fund through Grant number ACS PRF # 45882AC7, which partly funded this research. Furthermore, we would like to thank the Undergraduate Research Opportunity Program (UROP) at the University of Michigan which made the assistance of our undergraduate research assistants—Nicholas Orichella and Jeremy Shum—possible by providing them academic credit and financial support for working on this project. References Amelar, S., Eastman, C.E., Morris, R.L., Smeltzly, M.A., Lodge, T.P., Meerwall, E.D.V., 1991. Dynamic properties of low and moderate molecular weight polystyrenes at infinite dilution. Macromolecules 24, 3505–3516.
S. Jain et al. / Chemical Engineering Science 64 (2009) 4566 -- 4571
Armstrong, R.C., Gupta, S.K., Basaran, O., 1980. Conformational changes of macromolecules in transient elongational flow. Polymer Engineering and Science 20 (7), 466–472. Baldwin, P.R., Helfand, E., 1990. Dilute polymer solution in steady shear flow: nonNewtonian stress. Physical Review A 41 (12), 6772–6785. Cerf, R.J., 1969. Sur la theorie des proprietes hydrodynamiques de solutions de macromolecules en chaines. Journal de Chimie Physique et de Physico-Chimie Biologique 66, 479–488. Doi, M., Edwards, S.F., 1986. The Theory of Polymer Dynamics. University Press, Oxford. Doyle, P.S., Shaqfeh, E.S.G., Gast, A.P., 1997. Dynamic simulation of freely draining flexible polymers in steady linear flows. Journal of Fluid Mechanics 334, 251–291. Dunlap, P.N., Leal, L.G., 1987. Dilute polystyrene solutions in extensional flows: birefringence and flow modification. Journal of Non-Newtonian Fluid Mechanics 23, 5–48. Fisher, R.J., Denn, M.M., 1977. Mechanics of nonisothermal polymer melt spinning. A.I.Ch.E. Journal 23 (1), 23–28. Fixman, M., 1966. Polymer dynamics: non-Newtonian intrinsic viscosity. Journal of Chemical Physics 45 (3), 793–803. de Gennes, P.G., 1974. Coil-stretch transition of dilute flexible polymers under ultrahigh velocity gradients. Journal of Chemical Physics 60, 5030–5042. Helfand, E., Wasserman, Z.R., Weber, T.A., 1980. Brownian dynamics study of polymer conformational transitions. Macromolecules 13, 526–533. Herrchen, M., Ottinger, H.C., 1997. A detailed comparison of various FENE dumbbell models. Journal of Non-Newtonian Fluid Mechanics 68 (1), 17–42. Hinch, E.J., 1974. Mechanical models of dilute polymer solutions for strong flows with large polymer deformations. Colloques Internationaux du C.N.R.S., Polymeres et Lubrification 233, 241–247. Hinch, E.J., 1994. Uncoiling a polymer molecule in a strong extensional flow. Journal of Non-Newtonian Fluid Mechanics 54, 209–230. Hsieh, C.C., Li, L., Larson, R.G., 2003. Modeling hydrodynamic interaction in Brownian dynamics: simulations of extensional flows of dilute solutions of DNA and polystyrene. Journal of Non-Newtonian Fluid Mechanics 113, 147–191. Hur, J.S., Shaqfeh, E.S.G., Larson, R.G., 2000. Brownian dynamics simulations of single DNA molecules in shear flow. Journal of Rheology 44 (4), 713–742. Jain, S., Larson, R.G., 2008. Effect of bending and torsional potentials on highfrequency viscoelasticity of dilute polymer solutions. Macromolecules 41 (10), 3692–3700. James, D.F., Sridhar, T., 1995. Molecular conformation during steady-state measurements of extensional viscosity. Journal of Rheology 39, 713–724. King, D.H., James, D.F., 1983. Analysis of the Rouse model in extensional flow. II. Stresses generated in sink flow by flexible macromolecules and by finitely extended macromolecules. Journal of Chemical Physics 78 (7), 4749–4754. ¨ ¨ Kroger, M., Alba-Pérez, A., Laso, M., Ottinger, H.C., 2000. Variance reduced Brownian simulation of a bead-spring chain under steady shear flow considering hydrodynamic interaction effects. Journal of Chemical Physics 113 (11), 4767–4773. Kuhn, W., 1934. The shape of fibrous molecules in solution. Kolloidzschr 68, 2–15. Kumar, K.S., Prakash, J.R., 2003. Equilibrium swelling and universal ratios in dilute polymer solutions: Exact Brownian dynamics simulations for a delta function excluded volume potential. Macromolecules 36 (20), 7842–7856. Larson, R.G., 1988. Constitutive Equations for Polymer Melts and Solutions. Butterworth Publishers, Stoneham, MA. Larson, R.G., 1990. The unraveling of a polymer chain in a strong extensional flow. Rheological Acta 29 (5), 371–384.
4571
Larson, R.G., 2004a. An explanation for the high-frequency elastic response of dilute polymer solutions. Macromolecules 37 (13), 5110–5114. Larson, R.G., 2004b. Principles for coarse-graining polymer molecules in simulations of polymer fluid mechanics. Molecular Physics 102 (4), 341–351. Larson, R.G., 2005. The rheology of dilute solutions of flexible polymers: progress and problems. Journal of Rheology 49 (1), 1–70. Larson, R.G., Hu, H., Smith, D.E., Chu, S., 1999. Brownian dynamics simulations of a DNA molecule in an extensional flow field. Journal of Rheology 43 (2), 267–304. Liu, T.W., 1989. Flexible polymer chain dynamics and rheological properties in steady flows. Journal of Chemical Physics 90 (10), 5826–5842. Lodge, T.P., Miller, J.W., Schrag, J.L., 1982. Infinite-dilution oscillatory flow birefringence properties of polystyrene and poly(a-methylstyrene) solutions. Journal of Polymer Science: Polymer Physics Edition 20, 1409–1425. Lyulin, A.V., Adolf, D.B., Davies, G.R., 1999. Brownian dynamics simulations of linear polymers under shear flow. Journal of Chemical Physics 111 (2), 758–771. Magda, J.J., Larson, R.G., Mackay, M.E., 1988. Deformation-dependent hydrodynamic interaction in flows of dilute polymer solutions. Journal of Chemical Physics 89 (4), 2504–2513. Manke, C.W., Williams, M.C., 1992. Stress jumps predicted by the internal viscosity model with hydrodynamic interaction. Journal of Rheology 36 (7), 1261–1274. Orr, N.V., Sridhar, T., 1996. Stress relaxation in uniaxial extension. Journal of NonNewtonian Fluid Mechanics 67, 77–103. ¨ Ottinger, H.C., 1987. Generalized Zimm model for dilute polymer solutions under theta conditions. Journal of Chemical Physics 86 (6), 3731–3749. ¨ Ottinger, H.C., 1989. Gaussian approximation for Rouse chains with hydrodynamic interaction. Journal of Chemical Physics 90 (1), 463–473. Petrie, C.J.S., Denn, M.M., 1976. Instabilities in polymer processing. A.I.Ch.E. Journal 22 (2), 209–236. Phan-Thien, N., Manero, O., Leal, L.G., 1984. A study of conformation-dependent friction in a dumbbell model for dilute solutions. Rheological Acta 23 (2), 151–162. Rallison, J.M., 1997. Dissipative stresses in dilute polymer solutions. Journal of NonNewtonian Fluid Mechanics 68 (1), 61–83. Rallison, J.M., Hinch, E.J., 1988. Do we understand the physics in the constitutive equation? Journal of Non-Newtonian Fluid Mechanics 29, 37–55. Rouse, P.R., 1953. A theory of the linear viscoelastic properties of dilute solutions of coiling polymers. Journal of Chemical Physics 21 (7), 1272–1280. Ryckaert, J.P., Bellemans, A., 1975. Molecular dynamics of liquid n-Butane near its boiling point. Chemical Physics Letters 30 (1), 123–125. Ryskin, G., 1987. Calculation of the effect of polymer additive in a converging flow. Journal of Fluid Mechanics 178, 423–440. Somasi, M., Khomani, B., Woo, N.J., Hur, J.S., Shaqfeh, E.S.G., 2002. Brownian dynamics simulations of bead-rod and bead-spring chains: numerical algorithms and coarse-graining issues. Journal of Non-Newtonian Fluid Mechanics 108, 7–255. Underhill, P.T., Doyle, P.S., 2004. On the coarse-graining of polymers into bead-spring chains. Journal of Non-Newtonian Fluid Mechanics 122, 3–31. Underhill, P.T., Doyle, P.S., 2007. Accuracy of bead-spring chains in strong flows. Journal of Non-Newtonian Fluid Mechanics 145, 109–123. Wedgewood, L.E., 1989. A Gaussian closure of the second-moment equation for a Hookean dumbbell with hydrodynamic interaction. Journal of Non-Newtonian Fluid Mechanics 31, 127–142. Zimm, B.H., 1956. Dynamics of polymer molecules in dilute solution: viscoelasticity, flow birefrigence and dielectric loss. Journal of Chemical Physics 24 (2), 269–278.
Contents lists available at ScienceDirect
Chemical Engineering Science journal homepage: w w w . e l s e v i e r . c o m / l o c a t e / c e s
Do bending and torsional potentials affect the unraveling dynamics of flexible polymer chains in extensional or shear flows? Semant Jain, Indranil Saha Dalal, Nicholas Orichella, Jeremy Shum, Ronald Gary Larson ∗ Department of Chemical Engineering, University of Michigan, Ann Arbor, MI 48109, USA
A R T I C L E
I N F O
Article history: Received 1 October 2008 Received in revised form 16 May 2009 Accepted 22 May 2009 Available online 2 June 2009 Keywords: Computational chemistry Polymers Simulation Rheology Rouse model Brownian dynamics
A B S T R A C T
We compare the flow-induced unraveling of a bead-spring model chain in which each stiff spring corresponds to a single Kuhn step to that of an atomistic representation of the backbone of the same polymer chain with realistic bending and torsional potentials. In our earlier work [Jain, S. and Larson, R.G. 2008. Effect of bending and torsional potentials on high-frequency viscoelasticity of dilute polymer solutions. Macromolecules 41(10), 3692–3700], we showed that in the linear viscoelastic regime, bending and torsional potentials suppress fast local diffusive modes. Now, in strong shear and extensional flows, we observe that bending and torsional potentials have only a slight effect on unraveling dynamics. However, in shear flow, we find that for relatively short coarse-grained chains having less than 20 or so Kuhn steps—representing polymers with fewer than around 150 backbone bonds—coarse-graining introduces periodic peaks in the probability distribution of steady-state stretch. We believe this occurs at high shear because of the dominance of a discrete set of stretch values each of which corresponds to an integral number of bonds that are nearly fully aligned between back-folds that are created during the chain's tumbling orbit. Even this difference between coarse-grained and realistic fine-grained chains disappears for more typical long chain lengths. © 2009 Elsevier Ltd. All rights reserved.
1. Introduction Since atomistic molecular dynamics are still too expensive to apply to long polymer chains (Larson et al., 1999), coarse-grained (CG) bead-spring or bead-rod chains, or constitutive equations derivable from such models are routinely used to determine the dynamics of polymers in flow fields (Larson, 2005). Seminal work by Morton Denn and co-workers (Fisher and Denn, 1977; Petrie and Denn, 1976), starting in the 1970s, has shown the value of this approach. Computer speeds prior to the 1990s limited much of the modeling of that period to dumbbell models (Kuhn, 1934). Nonetheless, many physical phenomena were included in those models including internal viscosity (Cerf, 1969; Manke and Williams, 1992), dissipative stress (Hinch, 1994; Orr and Sridhar, 1996; Rallison, 1997), self-entanglements (Armstrong et al., 1980; James and Sridhar, 1995; King and James, 1983), kinked or yo-yo conformations (Larson, 1990; Rallison and Hinch, 1988; Ryskin, 1987), conformation-dependent friction (de Gennes, 1974; Hinch, 1974), and straining inefficiencies (Dunlap and Leal, 1987; Hinch, 1994; Phan-Thien et al., 1984).
∗ Corresponding author. Tel.: +1 734 936 0772; fax: +1 734 763 0459. E-mail addresses: [email protected] (S. Jain), [email protected] (I.S. Dalal), [email protected] (R.G. Larson). 0009-2509/$ - see front matter © 2009 Elsevier Ltd. All rights reserved. doi:10.1016/j.ces.2009.05.035
In addition, localized effects were captured by using a finer-grained bead-spring chain (Rouse, 1953; Zimm, 1956) that still coarse-grains the polymer chain into a sequence of frictional beads connected by Hookean or non-Hookean springs. The finest-grain model typically used to describe the deformation of long polymers in dilute solution is the Kramer's freely jointed bead-rod chain (Doyle et al., 1997; Larson, 2005; Underhill and Doyle, 2004). This chain is resolved at the level of a Kuhn step wherein the flexibility that is in reality distributed over several backbone bonds of the real polymer is lumped into each free joint. As a consequence, bead-spring chains and even bead-rod chains are not expected to capture very fast dynamics of synthetic polymers like polystyrene (Larson, 2004a). These fast dynamics are controlled by local bond rotations that are influenced by bending and torsional potentials—potentials that are left out of conventional bead-spring and even bead-rod models. While these local chain details do have an effect on the overall coil size, for very long polymers, their effect on chain dynamics at low frequencies or low shear rates is insignificant. In this limit, a number of computational methods have been developed that make similar predictions as long as both hydrodynamic interactions and excluded volume effects are appropriately accounted for. One method is the renormalization group method, which when applied to a long chain in a shearing flow, yields shear thinning that becomes stronger as excluded volume effects increase (Baldwin and Helfand,
S. Jain et al. / Chemical Engineering Science 64 (2009) 4566 -- 4571
1990). A more commonly applied method is to solve for the moments of a Fokker–Planck equation for the chain configuration probability distribution using normal mode decompositions (Fixman, 1966; ¨ Herrchen and Ottinger, 1997; Magda et al., 1988; Ottinger, 1987; ¨ Ottinger, 1989; Rouse, 1953; Wedgewood, 1989). Finally, in more recent years, Brownian dynamics simulations have become conven¨ tional (Hsieh et al., 2003; Kroger et al., 2000; Kumar and Prakash, 2003; Liu, 1989; Lyulin et al., 1999; Somasi et al., 2002; Underhill and Doyle, 2004; Underhill and Doyle, 2007). To obtain universal results both the Fokker–Planck equation and especially the Brownian dynamics method require care to be sure that the asymptotic longchain limit is reached (Kumar and Prakash, 2003). At shorter distance and/or time scales, the local effects such as bending and torsional potentials become important once again. In principle, determination of these effects require using atomistically resolved molecular dynamics simulations that include local details such as bending and torsional potentials and perhaps even interactions between solvent molecules and the polymer chain. While chains of modest length—containing dozens to hundreds of backbone bonds—are not long enough for their properties to converge on the long-chain universal results, they remain too long to be simulated by atomistic molecular dynamic simulations except at great cost. Hence, polymer chains containing a modest number (10–150) of Kuhn step lengths have been simulated by the Kramer's freely jointed bead-rod chain (Larson, 2005) which is resolved at the level of a Kuhn step. What is not clear is if such an approach (Hinch, 1994) can capture very fast dynamics that might be controlled by local bond rotations that are influenced by bending and torsional potentials. Hence, we are interested in determining the accuracy of the freely-jointed chain model in describing the dynamics of modestlength chains at high frequency or high shear rate. Our earlier work (Jain and Larson, 2008) first explored the influence of these bending and torsional potentials on the linear viscoelastic relaxation of fine-grained (FG) chains. We showed that their presence slowed down the local relaxation modes of the chain so that their contributions to viscoelastic relaxation overlapped those of some of the slower collective diffusive modes. This study helped explain why these local modes are not evident in viscoelastic spectra of dilute polymer solutions (Lodge et al., 1982), why coarse-grained bead-spring models are so successful in describing the linear viscoelasticity of dilute polymer solutions over such a wide frequency range, and why there exists an optimal degree of fine-graining in such models. For polystyrene, the optimal fine-graining corresponds to choosing one coarse-grained spring to represent a molecular mass of 5000 Da (Amelar et al., 1991), so that a chain of molecular mass 100,000 Da is represented by 20 such springs. However, our earlier work was limited to dilute polymers in the linear viscoelastic regime. In the present work, we extend the previous methods to strong extensional and shear flows. At increasingly high deformation rates, the portions of the chain that can be assumed to remain well equilibrated become progressively smaller (Larson, 2005). Therefore, local details of the chain—especially the bending and torsional potentials—might interfere with the ability of the chain to deform and unravel in such strong flows. Previously, we also tested the ability of coarse-grained beadspring chains (which lack bending and torsional potentials) to capture the effects of fine-grained chains with torsional and bending potentials by mapping CG chains onto the FG chains and performing Brownian dynamics simulations of both. In order to carry out this mapping, the overall chain length, radius of gyration, and total drag of the two chains were matched. We used Fraenkel springs possessing a preferred spring length for both the CG and FG chains. For the FG chains, the preferred spring length represented the length of a backbone bond which is 1.53 Å for a carbon–carbon bond. The preferred spring length of the CG chain was that of a Kuhn step
4567
length. Since the stiff springs in our CG chains could neither stretch nor compress much, our CG chains were similar to Kramer's freely jointed chains that are commonly used as coarse-grained models of real chains. Thus, the coarse-grained chain had an order of magnitude fewer springs and beads but each spring was longer and each bead exerted more drag than their fine-grained counterparts. This allowed us to assess quantitatively the effect of the bending and torsional potentials on the dynamics of dilute polymers and to determine the degree to which the CG chain captures the dynamics of the FG chain. In what follows, we extend this approach to strong shear and extensional flows using the same FG and CG models used in the earlier study of viscoelastic relaxation. As in the earlier study, we neglect both hydrodynamic interactions and excluded volume interactions between the beads. 2. Simulation 2.1. Force calculations We model the polymer chain as a series of beads connected by stiff Fraenkel springs and use the same nomenclature and force calculations as before (Jain and Larson, 2008) except for the presence of a non-zero solvent velocity field as described below (Jain and Larson, 2008; Larson, 2005) F iD = (V − x˙ i )
(1)
where F iD is the drag force on the ith bead, is the bead frictional drag constant, x˙ i is the velocity of the ith bead, V = ri · ∇V is the imposed solvent flow velocity at the bead's location, and ri is the position vector of the ith bead. In planar extensional flow, this velocity is ⎛ ⎞ ˙ 0 0 ⎝ (2) ∇V ≡ 0 −˙ 0 ⎠ 0 0 0 where ˙ is the extension rate. Similarly, under simple shearing flow, the velocity gradient is ⎛ ⎞ 0 0 0 (3) ∇V ≡ ⎝ ˙ 0 0 ⎠ 0 0 0 where ˙ is the shear rate. Since a CG chain has fewer beads than its equivalent FG chain, to maintain the same total frictional drag coefficient for the whole chain, T , the frictional drag per bead, , of the CG chain must be made proportionally higher than that of the FG chain, i.e., T |CG = T |FG (NS + 1)|FG (4) ⇒ |CG = |FG (NS + 1)|CG where NS is the number of springs in the FG or CG chain. Neglecting inertial forces, to conserve momentum on the ith bead, the drag force, F iD , must be set equal to the sum of the stretching, F iS , bending, F iB , torsional, F iT , and random Brownian, F iR , forces yielding (Jain and Larson, 2008; Larson, 2005) F iD = F iS + F iB + F iT + F iR (5) i i i i ⇒ x˙ i = F S + F B + F T + F R + V In the above equation, the stretching force (stretching constant,
S /m = 2.5×1027 s−2 ) maintains a separation distance of around 1.53 Å between two beads through stiff Fraenkel springs, the bending force (bending constant, /m = 1.3×107 J/kg) maintains the 109◦ 47 tetrahedral bending angle, the torsional force (torsional constant, /m = 6.6344×105 J/kg) retains the desired distribution of trans and gauche conformations and barriers between them through
4568
S. Jain et al. / Chemical Engineering Science 64 (2009) 4566 -- 4571
the Rychaert–Bellman potential (1975), the random Brownian force mimics the random motion of the solvent molecules, and the drag force captures the effect of extensional or shear flow field. In these force constants, m is the mass of the bead. Since we are only concerned here with tracking the chain configurations and the stress, we only need to pick a time step small enough for these properties to converge. We showed earlier that this was accomplished when we used a time step of t = 5×10−16 s, which is one-tenth of the time-step value chosen by Helfand and coworkers (1980)—from whose work we had taken the force constants for our simulations. Similar to our previous work, at this value of the time step, we find that the springs in the present case remain within 5% of their equilibrium length even at the highest extension and shear rates used here. 2.2. Coarse-grained chain To construct a coarse-grained chain whose coil size and fully stretched length are equal to those for a fine-grained polymer chain with bending and torsional potentials, we equate the mean square end-to-end vector, R2 0 , and the fully extended length, L, of the two chains via the relationships R2 0 = C∞ nl2
(6)
L = 0.82nl = NK bK
(7)
where NK is the number of bonds (Kuhn steps), bK is the length of a bond in the equivalent freely jointed chain model, n is the number of backbone bonds, l is the length of each backbone bond, and C∞ is the characteristic ratio of the chain which incorporates the expansion in the chain due to bending and torsional forces. Thus, for a FG chain with bending and torsional potentials typical of real chains (with the torsional potential at 100% “base-case” as defined earlier (Jain and Larson, 2008)), and having 60 springs of equilibrium length, l, of 1.53 Å, we construct an equivalent CG chain having eight springs each with an equilibrium length of 8.9938 Å. The longest relaxation time, 1 , can be computed from the Rouse formula (Doi and Edwards, 1986; Larson, 1988; Larson, 2004b):
1 =
2
8kB T S sin2 [ /2(NS + 1)]
(8)
where the FG bead drag coefficient is = 1.4×1012 kg/mol s (Helfand 2 et al., 1980), and S is the mean square end-to-end distance of 2 B C T 2 = 2 × 2.41(1.53 × a spring computed as 3/2 S = R2 0,S = C∞ ∞ 0 10−10 )2 m2 = 1.1283 × 10−19 m2 (Jain and Larson, 2008) which gives us 1 |FG = 6.426 × 10−9 s. Using Eq. (4) to ensure the same total drag on FG and CG chains, we compute the bead drag coefficient of the CG chain which, along with Eq. (8), results in 1 |CG = 6.865 × 10−9 s. The unraveling of a polymer in a strong flow is dependent on the initial configuration of the chain. To begin a detailed comparison of the unraveling dynamics of CG and FG chains, we first map the initial bead positions of the FG chain onto the CG chain. We would like to choose a corresponding initial configuration of the CG chain that matches as closely as possible that of the FG chain. Since the CG chain has roughly seven times fewer springs than the FG chain, this initial mapping is accomplished by selecting the bead position of the first bead, every seventh subsequent bead, and the last bead of the FG chain as the initial guess bead positions for the CG chain. These bead positions usually yield CG spring lengths that differ considerably from their equilibrium length. So, to begin the simulation at the same roughly equivalent equilibrated state, we then carry out a no-flow simulation with just the CG chain using a very small time step of 10−20 s—a value which is 1/50,000th of the time step used in the flow simulations. Since spring forces at such an non-equilibrated state can
Fig. 1. Mean dimensionless extension, x , as a function of Hencky strain, , at Wi = 2, 5, and 10 for CG and FG chains. The x values of FG and CG chain data at Wi = 10 have been augmented by 10 units to separate these results from the results at Wi = 5.
be large, a time step any larger could lead to instabilities. So, to enable all CG springs to relax to their preferred lengths, this relaxation process is carried out for 2×107 time steps which is ∼3×10−3 times the longest relaxation time of the chain—a time long enough to relax the springs to near-equilibrium length but not long enough for the overall chain conformation to otherwise change appreciably from the original starting configuration of the FG chain. This completes the mapping of a given starting configuration of a FG chain onto a roughly equivalent starting configuration of a CG chain enabling us to begin a direct comparison between the behavior of individual FG and CG chains with similar starting configurations. 3. Results Here we compare the average dynamics using an ensembleaverage over 100 FG chains (each with 60 springs) with an ensembleaverage over 100 corresponding CG chains (each with eight springs) for the extensional flow study. The ensemble was averaged over 50 chains for the shear flow study. We carry out each simulation at various Weissenberg numbers, Wi = ˙ 1 or Wi = ˙ 1 . The dimen2 sionless mean fractional extension, x = X / , is plotted against Hencky strain units ( = ˙ t or = ˙ t). X 2 is the root-mean-square farthest distance between any two beads of a chain in the flow direction and is the equilibrium spring length of the FG chain. 3.1. Extensional flow In Fig. 1, we compare the average dynamics of fine-grained and coarse-grained chains at three Weissenberg numbers. Simulations are run for 9–12 Hencky strain units which translate into 7.5×107 , 2.5×107 , and 1.5×107 time steps when Wi = 2, 5, and 10, respectively. At low Weissenberg numbers, the FG chains stretch a bit less than the CG chains because of the internal rigidity produced by the bending and torsional potentials. This difference diminishes at higher Weissenberg numbers where the large drag forces seemingly overwhelm any dynamic effect of these potentials. The error bars in Fig. 1 show that the standard error decreases at higher Weissenberg numbers. The ensemble-averaged dynamics of unraveling of CG chains are in reasonably good agreement with those of FG chains with the agreement becoming better with increasing Weissenberg number (Wi ⱖ 10). Since we matched the initial configuration of each FG chain to an equivalent CG chain, we can also compare the unraveling behavior of the individual FG and CG chains. In doing so, we find
S. Jain et al. / Chemical Engineering Science 64 (2009) 4566 -- 4571
Fig. 2. Starting configurations of the FG and CG chains in three dimensions. The FG chain is shown as a gray line with closed circles and the CG chain is represented as a black line with closed diamonds.
4569
Fig. 4. Two-dimensional projections (onto the XY plane) of the configurations of the CG and FG chains whose stretch versus time dynamics are shown in Fig. 3. The FG chain is shown as a gray line with filled circles and the CG chain as a black line with filled triangles. The configurations are shown at strains, , of 0, 10, 15, 20, 30, and 40. For clarity, the y-coordinates of the FG and CG chains shown at = 10, 15, 20, 30, and 40 have been shifted up vertically by 3, 6, 9, 12, and 15 units, respectively.
Fig. 3. Mean dimensionless extension, x , as a function of Hencky strain, , for the FG and CG chains with similar initial configurations as shown in Fig. 2.
that while most of the time the CG chain shows similar unraveling behavior as the corresponding FG chain, there are large differences in some cases. One such case for Wi = 10 is shown in Figs. 2–4 where the FG chain unravels slower than the mapped CG chain. Fig. 2 shows the starting configuration of the FG chain and of the CG chain that was mapped onto it. The dimensionless mean extension as a function of strain is given in Fig. 3 and the projections of the FG and CG chain onto the X−Y plane as a function of strain are shown in Fig. 4. Figs. 3 and 4 depict a simulation where the FG chain gets temporarily trapped in a folded state that is avoided by the CG chain. In some other cases, which have not been shown here, the CG chain gets trapped instead and unravels slower than the corresponding FG chain. Thus, we find that one of the chains gets trapped temporarily in a folded state in all cases where the CG chain and FG chain dynamics differ greatly. To determine if the large differences between unraveling of FG and CG chains that occur in the few cases are because of non-linear amplifications of small differences in starting configurations, rather than due to the differences between FG and CG chains, we run pairs of simulations of eight-spring CG chains. The two chains in each pair have very similar but not identical starting configurations. These pairs are created by relaxing a CG chain for 4.00 and 4.01 relaxation times which are then subject to an extensional flow at Wi = 10 using identical Brownian force histories for both chains. Despite the fact
Fig. 5. Unraveling dynamics of two different chains with slightly different initial configurations and identical Brownian force histories plotted as mean dimensionless extension, x , as a function of strain, . The first chain had its starting configuration relaxed up to four relaxation times (4.00, solid black) and the second to 4.01 relaxation times (4.01, solid gray).
that the starting configurations are very similar to each other, and the random forces on the beads are identical, dissimilar behavior is obtained for some pairs of chains. One example of this behavior is shown in Fig. 5 which proves that minor differences in starting configurations can be amplified to produce widely different unraveling behavior in a few cases. Thus, the large differences in some cases between the unraveling dynamics of CG and FG chains can be attributed simply to small differences in starting configurations resulting from the mapping of multiple FG bonds onto a single CG bond and not to intrinsic differences between FG and CG chains themselves. Given the unavoidable small differences in starting conditions, especially at high Weissenberg numbers, the CG model captures the dynamics of FG chains as faithfully as possible. 3.2. Shear flow The mean stretch of the FG and CG chains under shear flow at Wi = 10, 50, and 100, is averaged over 50 runs with the data collected
4570
S. Jain et al. / Chemical Engineering Science 64 (2009) 4566 -- 4571
Table 1 Comparison of mean stretch of fine-grained and coarse-grained chains under shear flow conditions. Simulation no.
1 2 3
Wi
10 50 100
Mean stretch (with standard error) √ √ x¯ |CG ± / N − 1 x¯ |FG ± / N − 1
Error (%)
20.66 ± 0.15 25.51 ± 0.13 27.19 ± 0.13
−7.06 −2.66 1.44
22.12 ± 0.26 26.19 ± 0.23 26.79 ± 0.17
1 − x¯ |CG / x¯ |FG
The mean stretch is computed using 50 chains (N = 50).
Fig. 7. Probability density as a function of mean fractional extension, x /L, at Wi = 100 for a FG chain having 60 springs and CG chains having 8, 16, and 32 springs. The maximum extension, L, of the FG chain is 60×0.82×1.53 Å = 75.276 Å.
the peaks are only a problem for short chains. Chains represented by 32 or greater CG springs—which correspond to 250 or greater backbone bonds—show negligible artificial peaks in stretch distribution even at Wi = 100. 4. Conclusions Fig. 6. Probability density as a function of mean dimensionless extension, x , at Wi = 10, 50, and 100 for CG and FG chains.
in the first 25 strain units not included in the analysis. Since all shear flow simulations are carried out for 5×108 time steps, this translates into ∼400, ∼1000, and ∼1900 strain units at Wi = 10, 50, and 100, respectively. This analysis has been shown in Table 1. Next, in Fig. 6, we compare the probability distribution as a function of mean dimensionless extension (Hur et al., 2000) at various Weissenberg numbers for the CG and FG chains. With 60 fine-grained springs, the maximum mean fractional extension is 60×0.82 = 49.2 which is rounded up to 50 and then divided into 100 equally spaced bins. The factor of 0.82 arises on account of the zig-zag contour of the FG chain because of the bending angles. We notice that the overall shape of the distributions for FG and CG chains are similar but the CG chain distribution contains peaks or oscillations. These oscillations are due to the high degree of orientation often exhibited by the springs at high Weissenberg numbers and the small number (8) of stiff springs. Counting the small peak when x = 10, we observe seven equally spaced peaks which are one less than the eight-coarse grained springs in the chain. Evidently, any peak corresponding to the length of one spring is negligible. In a folded state of well-oriented springs, as the stiff Fraenkel springs each maintain a fixed length, only a discrete number of springs is aligned between folds and the preferred stretch values appear to correspond to these discrete states. To demonstrate that the oscillations are due to the small number of CG springs, we increase the number of springs from eight to 16 and then to 32. As in our earlier simulations, we show the probability distribution averaged over 50 molecules at Wi = 100 in Fig. 7. To retain the same Weissenberg number as the spring length is increased, in accordance with the Rouse theory, the strain rate is reduced as ˙ = Wi/(N|CG )2 where N|CG is the number of springs in the CG chain. While the chain with eight springs has seven peaks, these peaks become more numerous and diminish in magnitude as the number of springs in the chain is increased to 16 and 32. Thus,
We compared the non-linear response of a locally realistic finegrained bead-spring model with stretching, bending, torsional, and drag forces to that of an equivalent coarse-grained chain with only stretching and drag forces. At low extension rates, the CG chain deforms slightly more than the FG chain. This difference disappears at higher extension rates indicating that a CG model should be sufficient when the Weissenberg number is greater than ten. Our individual chain comparisons reveal that even minute differences in the starting configuration of chain can lead to drastically different unraveling behavior. This behavior is especially pronounced if any chain gets trapped during the unraveling process. In shear flow, the CG chain has a similar steady-state distribution of stretch as the FG chain does except that the former exhibits oscillations, or multiple values of stretch, resulting from the small number of springs used in the CG chain. We verify this by using longer chains that correspond to 250 or greater backbone bonds and find that the oscillatory behavior disappears almost entirely. Thus, for chains as long or longer than this, especially at high Weissenberg numbers, coarse-grained chains can accurately capture the fine-grained chain dynamics. Acknowledgments We gratefully acknowledge the financial support provided by Petroleum Research Fund through Grant number ACS PRF # 45882AC7, which partly funded this research. Furthermore, we would like to thank the Undergraduate Research Opportunity Program (UROP) at the University of Michigan which made the assistance of our undergraduate research assistants—Nicholas Orichella and Jeremy Shum—possible by providing them academic credit and financial support for working on this project. References Amelar, S., Eastman, C.E., Morris, R.L., Smeltzly, M.A., Lodge, T.P., Meerwall, E.D.V., 1991. Dynamic properties of low and moderate molecular weight polystyrenes at infinite dilution. Macromolecules 24, 3505–3516.
S. Jain et al. / Chemical Engineering Science 64 (2009) 4566 -- 4571
Armstrong, R.C., Gupta, S.K., Basaran, O., 1980. Conformational changes of macromolecules in transient elongational flow. Polymer Engineering and Science 20 (7), 466–472. Baldwin, P.R., Helfand, E., 1990. Dilute polymer solution in steady shear flow: nonNewtonian stress. Physical Review A 41 (12), 6772–6785. Cerf, R.J., 1969. Sur la theorie des proprietes hydrodynamiques de solutions de macromolecules en chaines. Journal de Chimie Physique et de Physico-Chimie Biologique 66, 479–488. Doi, M., Edwards, S.F., 1986. The Theory of Polymer Dynamics. University Press, Oxford. Doyle, P.S., Shaqfeh, E.S.G., Gast, A.P., 1997. Dynamic simulation of freely draining flexible polymers in steady linear flows. Journal of Fluid Mechanics 334, 251–291. Dunlap, P.N., Leal, L.G., 1987. Dilute polystyrene solutions in extensional flows: birefringence and flow modification. Journal of Non-Newtonian Fluid Mechanics 23, 5–48. Fisher, R.J., Denn, M.M., 1977. Mechanics of nonisothermal polymer melt spinning. A.I.Ch.E. Journal 23 (1), 23–28. Fixman, M., 1966. Polymer dynamics: non-Newtonian intrinsic viscosity. Journal of Chemical Physics 45 (3), 793–803. de Gennes, P.G., 1974. Coil-stretch transition of dilute flexible polymers under ultrahigh velocity gradients. Journal of Chemical Physics 60, 5030–5042. Helfand, E., Wasserman, Z.R., Weber, T.A., 1980. Brownian dynamics study of polymer conformational transitions. Macromolecules 13, 526–533. Herrchen, M., Ottinger, H.C., 1997. A detailed comparison of various FENE dumbbell models. Journal of Non-Newtonian Fluid Mechanics 68 (1), 17–42. Hinch, E.J., 1974. Mechanical models of dilute polymer solutions for strong flows with large polymer deformations. Colloques Internationaux du C.N.R.S., Polymeres et Lubrification 233, 241–247. Hinch, E.J., 1994. Uncoiling a polymer molecule in a strong extensional flow. Journal of Non-Newtonian Fluid Mechanics 54, 209–230. Hsieh, C.C., Li, L., Larson, R.G., 2003. Modeling hydrodynamic interaction in Brownian dynamics: simulations of extensional flows of dilute solutions of DNA and polystyrene. Journal of Non-Newtonian Fluid Mechanics 113, 147–191. Hur, J.S., Shaqfeh, E.S.G., Larson, R.G., 2000. Brownian dynamics simulations of single DNA molecules in shear flow. Journal of Rheology 44 (4), 713–742. Jain, S., Larson, R.G., 2008. Effect of bending and torsional potentials on highfrequency viscoelasticity of dilute polymer solutions. Macromolecules 41 (10), 3692–3700. James, D.F., Sridhar, T., 1995. Molecular conformation during steady-state measurements of extensional viscosity. Journal of Rheology 39, 713–724. King, D.H., James, D.F., 1983. Analysis of the Rouse model in extensional flow. II. Stresses generated in sink flow by flexible macromolecules and by finitely extended macromolecules. Journal of Chemical Physics 78 (7), 4749–4754. ¨ ¨ Kroger, M., Alba-Pérez, A., Laso, M., Ottinger, H.C., 2000. Variance reduced Brownian simulation of a bead-spring chain under steady shear flow considering hydrodynamic interaction effects. Journal of Chemical Physics 113 (11), 4767–4773. Kuhn, W., 1934. The shape of fibrous molecules in solution. Kolloidzschr 68, 2–15. Kumar, K.S., Prakash, J.R., 2003. Equilibrium swelling and universal ratios in dilute polymer solutions: Exact Brownian dynamics simulations for a delta function excluded volume potential. Macromolecules 36 (20), 7842–7856. Larson, R.G., 1988. Constitutive Equations for Polymer Melts and Solutions. Butterworth Publishers, Stoneham, MA. Larson, R.G., 1990. The unraveling of a polymer chain in a strong extensional flow. Rheological Acta 29 (5), 371–384.
4571
Larson, R.G., 2004a. An explanation for the high-frequency elastic response of dilute polymer solutions. Macromolecules 37 (13), 5110–5114. Larson, R.G., 2004b. Principles for coarse-graining polymer molecules in simulations of polymer fluid mechanics. Molecular Physics 102 (4), 341–351. Larson, R.G., 2005. The rheology of dilute solutions of flexible polymers: progress and problems. Journal of Rheology 49 (1), 1–70. Larson, R.G., Hu, H., Smith, D.E., Chu, S., 1999. Brownian dynamics simulations of a DNA molecule in an extensional flow field. Journal of Rheology 43 (2), 267–304. Liu, T.W., 1989. Flexible polymer chain dynamics and rheological properties in steady flows. Journal of Chemical Physics 90 (10), 5826–5842. Lodge, T.P., Miller, J.W., Schrag, J.L., 1982. Infinite-dilution oscillatory flow birefringence properties of polystyrene and poly(a-methylstyrene) solutions. Journal of Polymer Science: Polymer Physics Edition 20, 1409–1425. Lyulin, A.V., Adolf, D.B., Davies, G.R., 1999. Brownian dynamics simulations of linear polymers under shear flow. Journal of Chemical Physics 111 (2), 758–771. Magda, J.J., Larson, R.G., Mackay, M.E., 1988. Deformation-dependent hydrodynamic interaction in flows of dilute polymer solutions. Journal of Chemical Physics 89 (4), 2504–2513. Manke, C.W., Williams, M.C., 1992. Stress jumps predicted by the internal viscosity model with hydrodynamic interaction. Journal of Rheology 36 (7), 1261–1274. Orr, N.V., Sridhar, T., 1996. Stress relaxation in uniaxial extension. Journal of NonNewtonian Fluid Mechanics 67, 77–103. ¨ Ottinger, H.C., 1987. Generalized Zimm model for dilute polymer solutions under theta conditions. Journal of Chemical Physics 86 (6), 3731–3749. ¨ Ottinger, H.C., 1989. Gaussian approximation for Rouse chains with hydrodynamic interaction. Journal of Chemical Physics 90 (1), 463–473. Petrie, C.J.S., Denn, M.M., 1976. Instabilities in polymer processing. A.I.Ch.E. Journal 22 (2), 209–236. Phan-Thien, N., Manero, O., Leal, L.G., 1984. A study of conformation-dependent friction in a dumbbell model for dilute solutions. Rheological Acta 23 (2), 151–162. Rallison, J.M., 1997. Dissipative stresses in dilute polymer solutions. Journal of NonNewtonian Fluid Mechanics 68 (1), 61–83. Rallison, J.M., Hinch, E.J., 1988. Do we understand the physics in the constitutive equation? Journal of Non-Newtonian Fluid Mechanics 29, 37–55. Rouse, P.R., 1953. A theory of the linear viscoelastic properties of dilute solutions of coiling polymers. Journal of Chemical Physics 21 (7), 1272–1280. Ryckaert, J.P., Bellemans, A., 1975. Molecular dynamics of liquid n-Butane near its boiling point. Chemical Physics Letters 30 (1), 123–125. Ryskin, G., 1987. Calculation of the effect of polymer additive in a converging flow. Journal of Fluid Mechanics 178, 423–440. Somasi, M., Khomani, B., Woo, N.J., Hur, J.S., Shaqfeh, E.S.G., 2002. Brownian dynamics simulations of bead-rod and bead-spring chains: numerical algorithms and coarse-graining issues. Journal of Non-Newtonian Fluid Mechanics 108, 7–255. Underhill, P.T., Doyle, P.S., 2004. On the coarse-graining of polymers into bead-spring chains. Journal of Non-Newtonian Fluid Mechanics 122, 3–31. Underhill, P.T., Doyle, P.S., 2007. Accuracy of bead-spring chains in strong flows. Journal of Non-Newtonian Fluid Mechanics 145, 109–123. Wedgewood, L.E., 1989. A Gaussian closure of the second-moment equation for a Hookean dumbbell with hydrodynamic interaction. Journal of Non-Newtonian Fluid Mechanics 31, 127–142. Zimm, B.H., 1956. Dynamics of polymer molecules in dilute solution: viscoelasticity, flow birefrigence and dielectric loss. Journal of Chemical Physics 24 (2), 269–278.