Topological Methods for Polymeric Materials: Characterizing the Relationship Between Polymer Entanglement and Viscoelasticity

Abstract

We develop topological methods for characterizing the relationship between polymer chain entanglement and bulk viscoelastic responses. We introduce generalized Linking Number and Writhe characteristics that are applicable to open linear chains. We investigate the rheology of polymeric chains entangled into weaves with varying topologies and levels of chain density. To investigate viscoelastic responses, we perform non-equilibrium molecular simulations over a range of frequencies using sheared Lees–Edwards boundary conditions. We show how our topological characteristics can be used to capture key features of the polymer entanglements related to the viscoelastic responses. We find there is a linear relation over a significant range of frequencies between the mean absolute Writhe $$Wr$$ and the Loss Tangent $$tan\left(\delta \right)$$. We also find an approximate inverse linear relationship between the mean absolute Periodic Linking Number $$L{K}_P$$ and the Loss Tangent $$tan\left(\delta \right)$$. Our results show some of the ways topological methods can be used to characterize chain entanglements to better understand the origins of mechanical responses in polymeric materials.

1. Introduction

A central aim in polymeric material science is to understand the relationships between chemistry, molecular-level interactions, and bulk material properties. We consider polymeric materials and develop approaches for investigating the relationship between the topology of polymeric chain interactions and the resulting bulk viscoelastic mechanical responses. The collective configurations of the chains can result in complex entangled structures that greatly restrict chain motions and augment transmission of mechanical stresses within the material [1,2,3]. Characterizing the intuitive notion of entanglement in a way that is quantitatively precise poses many interesting challenges. For instance, entanglement in polymer melts or gels often involves kinetics over a broad range of time-scales with collective responses that often depend on both the frequency and duration of the mechanical perturbations. In addition, the polymer chains often exhibit friction with respect to one another or coupling from immersion within a solvent fluid. This can result in significant variations in the viscoelastic responses that may depend sensitively on frequency [1,4,5].

When investigating polymer chain entanglement, a common approach is to consider two length and time scales. First, there is the length and time scale of the entire chain, where global entanglement occurs as the chains get knotted and linked to each other [1,6,7,8,9,10]. Second, one has local entanglement arising from the local constraints, obstacles, that a chain feels over a small length and time scale. This entanglement length often plays a central role in models of entangled polymer dynamics and is closely related to the tube diameter in tube model theories [3,11,12,13,14]. Understanding the relation between these scales poses a number of challenges in practice. For instance, experimentally there are different methods for the determination of the plateau modulus and different entanglement molecular weights give different results [15]. Since Edwards’ tube model, several improvements of this theory have emerged, but inconsistencies still remain [16,17,18,19]. Moreover, an even shorter length scale, that of the packing length, has been shown to have an influence and has been incorporated in some theories of polymer viscoelasticity [20,21,22].

The relation between entanglement and viscoelastic properties of materials has been studied indirectly by varying the density or molecular weight of the chains thereby influencing the number of contacts between neighboring chains [23,24,25,26,27]. A measure of polymer entanglement that has been very helpful in such studies is the number of “kinks” per chain, which is derived from the application of a contour reduction algorithm on a polymer melt [25,26,27,28]. However, the global entanglement complexity is more subtle and cannot be assessed by only measuring the number of contacts. This has led to the use of tools from knot theory to study entanglement in polymers [6,7,25,29,30,31,32]. The difficulty in using tools from knot theory is that they are defined on simple closed curves in space (rings) while the polymers often have other architectures, such as being open linear chains. Toward dealing with this issue, a new statistical definition of knotting was introduced in [33]. This method can determine the principal knot type of a fixed configuration of an open chain. However, the polymer chains move in time and this method may have problems when applied to study entanglement in non-equilibrium conditions. Moreover, the tools from knot theory have not yet been developed sufficiently for polymeric systems with three-dimensional Periodic Boundary Conditions (PBC) which are often used in practice [34].

In our work, we take a complimentary approach building on mathematical ideas from topology and geometry to quantify the complexity of the polymer entanglements. We introduce two quantities referred to as the collective “Linking Number” and “Writhe” [35,36,37]. The advantage of the Gauss linking integral is that it can be applied to both linear and ring polymers and it is a continuous function of the chain coordinates. Moreover, in [35] it has been shown that it can be extended to systems with 1,2 or 3 PBC to provide a new continuous measure of entanglement. These measures have been applied to study polymer entanglement in both equilibrium and non-equilibrium conditions. More precisely, it was shown that the writhe in combination with the Z1 algorithm can provide a new estimator of the entanglement length with several advantages over other estimators and our proposed measures have also been used to understand the dis-entanglement of polymer chains in a melt under an elongational force [36,37,38]. Recently, the Gauss linking integral has been also used to study protein folding kinetics [39]. These results indicate the promise of such topology-based estimators in polymer theories.

Here, we develop further these topological approaches for characterizing polymer entanglements within non-equilibrium systems subject to external materials deformations, such as oscillatory shearing. To investigate viscoelastic responses, we develop methods for non-equilibrium three-dimensional molecular simulations with shearing Lees–Edwards periodic boundary conditions (LE-PBC) [40]. We use our approach to study the frequency dependent viscoelastic responses and relationships to the underlying topology of the polymer chain entanglements.

We consider in our investigations polymeric weaves which have well-controlled topological properties that we can vary. We consider polymeric systems that have local topologies arising from short linear polymer chains and those with global topologies arising from long linear polymer chains. We also vary density to consider systems ranging from weakly entangled to strongly entangled. Our topological approaches allow us to investigate both the local and global entanglement effects. We find an approximate linear relation over a large range of frequencies between the mean absolute Writhe $Wr$ and the loss tangent $tan\left(\delta \right)$. We also find an approximate inverse linear relationship between the mean absolute Periodic Linking Number $L{K}_P$ and the loss tangent $tan\left(\delta \right)$. We expect our topological approaches could be useful in gaining additional quantitative information relevant to understanding the mechanics of polymeric materials.

We organize the paper as follows. In Section 2 we introduce approaches from knot theory to precisely characterize the topology of the polymeric chains. In Section 3 we describe a class of polymeric materials having a weave-like topology. In Section 4, we discuss the details of our computational methods and simulation approaches. In Section 5, we present our results and discuss relationships revealed by our methods between the topology of the polymer entanglements and bulk material responses.

2. Characterizing Polymer Entanglement

We measure the degree to which polymer chains interwind and attain complex configurations using the Gauss Linking Integral [41]. We define the Gauss Linking Number as

(1)$\begin{array}{c}\hfill L(l_1,l_2)=\frac{1}{4\pi }{\int }_{[0,1]}{\int }_{[0,1]}\frac{({\dot{\gamma }}_1\left(t\right),{\dot{\gamma }}_2\left(s\right),{\gamma }_1\left(t\right)-{\gamma }_2\left(s\right))}{\left|\right|{\gamma }_1\left(t\right)-{\gamma }_2\left(s\right){\left|\right|}^3}dtds.\end{array}$

We use this for two disjoint (closed or open) oriented curves $l_1$ and $l_2$ whose arc-length parametrizations are respectively ${\gamma }_1\left(t\right),{\gamma }_2\left(s\right)$. The double integral is over $l_1$ and $l_2$. In this notation $({\dot{\gamma }}_1\left(t\right),{\dot{\gamma }}_2\left(s\right),{\gamma }_1\left(t\right)-{\gamma }_2\left(s\right))$ denotes the scalar triple product of ${\dot{\gamma }}_1\left(t\right),{\dot{\gamma }}_2\left(s\right)$ and ${\gamma }_1\left(t\right)-{\gamma }_2\left(s\right)$. The Gauss Linking Number is a topological invariant for closed chains and a continuous function of the chain coordinates for open chains (see Figure 1 for an illustrative example).

We also define a one chain measure for the degree of intertwining of the chain around itself.

We define the Writhe of a chain as

(2)$\begin{array}{c}\hfill Wr\left(l\right)=\frac{1}{4\pi }{\int }_{[0,1]}{\int }_{[0,1]}\frac{(\dot{\gamma }\left(t\right),\dot{\gamma }\left(s\right),\gamma \left(t\right)-\gamma \left(s\right))}{\left|\right|\gamma \left(t\right)-{\gamma \left(s\right)\left|\right|}^3}dtds.\end{array}$

For a curve ℓ with arc-length parameterization $\gamma \left(t\right)$ is the double integral over l. The Writhe is a continuous function of the chain coordinates for both open and closed chains.

For systems employing Periodic Boundary Conditions (PBC), the linking that is imposed from one simulated chain on another chain propagates in three dimensional space by the images of the other chain. In other words, for a system with PBC each simulated chain gives rise to a free chain in the periodic system which consists of an infinite number of copies of the simulated chain. We call each copy of a chain an image of the free chain. It has been shown that a measure of entanglement that can capture the global linking in a periodic system is the periodic linking number $L{K}_P$ [35]:

We define the Periodic Linking Number as

(3)$\begin{array}{c}\hfill L{K}_P(I,J)=\sum _{v}L(I_u,J_v),\end{array}$

The I and J denote two (closed, open or infinite) free chains in a periodic system. Suppose that $I_u$ is an image of the free chain I in the periodic system. The sum is taken over all the images $J_v$ of the free chain J in the periodic system. We say that the Periodic Linking Number $L{K}_P$ is between two free chains I and J.

The Periodic Linking Number is a topological invariant for closed chains and a convergent series for open chains that changes continuously with the chain coordinates. For its computation, we use a cutoff, the local Periodic Linking Number [35,38].

3. Polymeric Materials with Weave-like Entanglements

We study the role of entanglement topology in the mechanical responses of polymeric materials. We consider both the case when the global topology is fixed and when the topology can change over time for a few different weave entanglements, see Figure 2 and Figure 3. We investigate how mechanical responses depend on the topology, chain density, or whether the polymers are to be considered open chains or closed (infinite) chains.

We refer to chains as closed when they are very large relative to the length-scale of the entanglements. In practice, we can think of these as effectively infinitely long chains. These chains have the important property that in computations no end-points occur within the simulation box. These chains are treated as extending periodically to create a toopology corresponding to an infinite weave.

We refer to chains as open when they are finite in length. These chains are open in the sense they always have end-points within the simulation box. These chains can still cross the periodic boundary where they interact with the periodic image points generated by the unit cell. In this case, the topology of the material simulated can change over time in response to the deformations and stresses of the material.

We consider in our studies the specific weave topologies referred to as (w0) for aligned, (wI, wII) for orthogonal non-interlaced at different densities and (wIII) alternating interlacing. We show the base-line chain density and lengths in Table 1. We show example configurations of each of these weaves in Figure 2 and Figure 3.

We start with the weave w0 which has a relatively simple global topology. The polymers are simply arranged parallel to one another without any entanglement, providing a good reference case for topology and responses. We then consider the weaves wI and wII that arrange polymers in a regular criss-cross patterm. Both wI and wII have the same global topology, but we take these to have different densities. We take weave wIII to be an inter-woven topology alternating in-out entanglements. We take wIII to have the same density as wII. The weave wIII has a non-trivial global topology. As a result of the polymer sterics preventing the crossing of chains, the global topology of the infinite systems cannot change without breakage of the bonds. It should be noted that the local configurations of the chain interactions can lead to local entanglements that change over time under the global constraints of the topology.

We create open chain systems for w0, wI, wII, wIII by deleting the same bond from each chain in the simulation box. Given the periodicity this creates a standard procedure for obtaining an initial open-chain configuration. As a result of the polymer sterics preventing the crossing of chains, the global topology of the infinite systems cannot change without breakage of the bonds. In contrast to the infinite chain systems, the global topology of the open chain system can change over time by slippage of the chains past the entanglements. This creates the possibility in response to mechanical deformations for topology rearrangements over time.

We investigate the mechanical responses of the material by subjecting the polymeric chains to an oscillatory shear through deformation of the simulation box in the style of Lees–Edwards [40]. This provides us with shear stresses for the material which we can correlate with simultaneous measurements of the chain density and topology of the materials. We remark that similat to our polymeric weave configurations, there have also been studies using weaves for investigating metal organic frameworks and crystals [42,43].

4. Simulation of the Polymers

We investigate entangled polymeric chains using three-dimensional molecular simulations. The polymers are treated as elastic macromolecules modeled with harmonic bond potential of energy $E=K_b{(r-r_0)}^2$, $K_b=250$, r denotes the length of extension of the polymer bonds and $r_0=1$ denotes the rest length of the bond. The polymer bending stiffness is controlled with a harmonic angle potential with energy $E=K_{\theta }(1-cos(\theta -{\theta }_0))$, with $K_{\theta }=8$, where $\theta$ is the angle between two consecutive bonds. The rest angle is ${\theta }_0=\pi$. The length of the simulation box is approximately $20\sigma$. Each polymer chain has approximately 20 monomers inside the simulation box for the densities and parameters given in Table 1 and Table 2. With this choice of $K_{\theta }$, the chains have persistence length of approximately $1/5$ of the length of the simulation box. With these potentials, there is no maximum permitted length or bond angle constraints, but there is a high energy penalty for large deviations from the rest length. This does not exclude the possibility of chains crossing through each other, especially for large deformations. Our results however, show that chain crossings are rare enough so as to not influence the qualitative effects of entanglement observed here (see Section 5). The beads of our polymers interact through the Lennard–Jones (LJ) potential with energy

(4)$E_{LJ}=4ϵ\left({\left(\frac{\sigma }{r}\right)}^{12}-{\left(\frac{\sigma }{r}\right)}^6\right)$

We use a cutoff of $2.5\sigma$. We simulate the finite temperature and kinetics of the polymer chain dynamics using the Langevin Thermostat

(5)$m\frac{dV}{dt}=-\mathsf{{\rm Y}}V-\nabla \Phi \left(X\right)+\sqrt{2k_{B}T\gamma }\frac{d{B}_t}{dt}.$

The X denotes the position of the atoms, $V=dX/dt$ is the velocity, $-\nabla \Phi \left(X\right)$ denotes the interaction forces, $\mathsf{{\rm Y}}$ denotes the friction coefficient, and $\sqrt{2k_{B}T\gamma }d{B}_t/dt$ denotes the random force accounting for thermal fluctuations [44]. We perform all simulations using the LAMMPS molecular dynamics package and our custom extension packages [45,46].

To study the bulk mechanics of the polymeric system, we perform rheological studies using oscillatory shearing motions based on Lees–Edwards boundary conditions [40]. In Lees–Edwards conditions, periodic boundary conditions are used with shifted image interactions (see Figure 4). We use a sinusoidal oscillation of the displacement $L\left(t\right)=L_0+Asin(2\pi t/T_p)$ with amplitude A and time periodicity $T_p$. This corresponds to a cosine oscillation of the strain with rate $\dot{\gamma }={\dot{\gamma }}_{0}cos\left(\omega t\right)$ where $\omega =2\pi /T_p$ and ${\dot{\gamma }}_0=A\omega$.

As a measure of material response, we consider the dynamic complex modulus $G\left(\omega \right)=G_1\left(\omega \right)+i{G}_2\left(\omega \right)$. The components are defined from measurements of the stress as the least-squares fit of the periodic stress component ${\sigma }_{xy}$ by the function $g\left(t\right)=G_1\left(\omega \right){\gamma }_{0}sin\left(\omega t\right)+G_2\left(\omega \right){\gamma }_{0}cos\left(\omega t\right)$. This offers a characterization of the response of the material to oscillating applied shear stresses and strains as the frequency $\omega$ is varied. The $G_1$ is referred as the Elastic Storage Modulus and $G_2$ is described as the Viscous Loss Modulus. These dynamic moduli are motivated by considering the linear response of the stress components ${\sigma }_{xy}$ to applied stresses and strains. At low frequency the polymer stresses appear to have sufficient time to equilibrate to the applied shear stresses. At high frequencies, the polymer stresses do not appear to have sufficient time to equilibrate to the applied shear stresses. This is manifested in ${\sigma }_{xy}\left(t\right)$ which is seen to track the applied stress very closely. A phase lag 0 is representative of solids and $\pi /2$ is representative of liquids. This delay is caused by propagating the stress through the domain via the chain topology. The increase of $G_2$ indicates that the mechanics arises effectively from chains’ resistance to more rapid motions, such as sliding, while the increase of $G_1$ indicates in the mechanics a resistance to direct deformation represented by increases in the elastic bond lengths or from the bending stiffness of chains.

To estimate the dynamic complex modulus in practice, the least-squares fit is performed for ${\sigma }_{xy}\left(t\right)$ over the entire stochastic trajectory of the simulation after some transient period of approximately $10T$ (see Table 4), which is of the order of the diffusion time of the open chains under study. In our simulations, the maximum strain over each period was chosen to be half the periodic unit cell in the x-direction, corresponding to a strain amplitude ${\gamma }_0=\frac{1}{2}$. A description of the parameters and specific values used in the simulations can be found in Table 2 and Table 3. We notice that the applied strain is large and would imply a nonlinear regime for polymer melts. However, the systems considered in this study are polymer solutions of very low molecular weight and our results indicate a linear regime so we can neglect higher harmonic contributions to stress [47].

The effective stress tensor associated with the polymers at a given time is estimated using the Irving-Kirkwood method [1,48]

(6)${\sigma }_{l,k}=\frac{1}{V}\sum _n\sum _{j=1}^{n-1}⟨f_j^{\left(l\right)}·(x_{q_n}^{\left(k\right)}-x_{q_j}^{\left(k\right)})⟩$

where V is the volume of the periodic box, $x_{q_v}^{\left(k\right)}$ is the k-th coordinate of the $q_v$-th atom (the minimum image convention is used for the difference) and $f_j^{\left(l\right)}$ is the l-th coordinate of the pairwise interaction between the two atoms.

Table 4 shows how the simulation time and oscillation period range compare to characteristic times in our systems (computed using the parameters used in our simulations, shown in the previous tables). The advection time is computed as ${\tau }_A=m/\mathsf{{\rm Y}}$. The Rouse time, ${\tau }_R$, is computed for an ideal chain of length 20. The critical time ${\tau }_0$, refers to the characteristic time where cross-overs are observed in our simulations, shown in Section 5. Notice that we do not report an entanglement time because our chains are short (with number of discrete local topological constraints $Z<2$ in many cases) and the notion of entanglement length does not apply to them.

5. Bulk Mechanical Responses

5.1. Complex Modulus

We show the log-log plot of the Elastic Storage Modulus $G_1$ and Viscous Loss Modulus $G_2$ for all the infinite and open weave polymeric materials as the shear response frequency is varied in Figure 5. The frequency of oscillation is normalized by ${\omega }_0=2\pi /{\tau }_0$ where ${\tau }_0=943\tau =598ps=1.98{\tau }_D$ is a time-scale on the order of the diffusion time ${\tau }_D$ (see Table 4).

Comparing $G_1,G_2$ for the infinite weaves we see that, in the range of frequencies studied, we have $G_1>G_2$ for all the simple weaves. The crossover of $G_1$ and $G_2$ is absent for those systems within this range of frequencies, which indicates no behavioral transition in the samples which exhibit solid properties. When $G_1$ is larger than $G_2$ the elastic response is dominant indicating there is relatively few polymer rearrangements (reptation) within the network structure. This indicates that energy is mainly stored elastically in the stretching and bending of bridging polymeric chains. This can be verified by our Writhe quantity for the chains as it reaches a minimum at the extrema of the oscillatory strain period within this regime (see Section 5.2). The systems with large $G_1$ behave like stiff materials having strong entanglements similar to imperfect networks having transient covalent crosslinking [51,52,53,54,55,56,57]. This indicates that polymer solutions of long linear semiflexible chains can behave like crosslinked networks, even in the absence of explicit crosslinks. A similar observation was reported in [58].

Initially, $G_1$ and $G_2$ are independent of the frequency of oscillation and we see a crossover at frequency ${\omega }_0$ that corresponds to period ${\tau }_0$. At frequencies higher than ${\omega }_0$ (period times shorter than ${\tau }_0$) there is a significant dependence of moduli on the frequency which increases with increasing topological complexity. This is in agreement with predictions for polymeric networks [49]. The line segments shown in the figure indicate a scaling between ∼${\omega }^{1/4}$ and ∼${\omega }^{1/2}$, respectively, to be compared with that of Rouse chains.

The alternate interlacing weave, wIII, is the only infinite weave for which $G_1$, $G_2$ intersect and for which $G_1$ and $G_2$ both seem to scale as ${\omega }^{1/2}$ in the intermediate frequencies. Moreover, for wIII, $G_1\approx G_2$, with $G_1<G_2$ for low frequencies. We find that the original configuration of wIII is not favorable to the stiffness of the chains and the chains need to stretch resulting in a larger $G_1$. This causes extra collisions with other chains which results in larger values of $G_2$ as well. At high frequencies, we notice a shift from filament bending to stretching which results in higher values of $G_1$. Such transitions have also been observed in networks of actin filaments [59,60].

Comparing $G_1,G_2$ for the open systems we find that both $G_1$ and $G_2$ are initially constant up to $\omega \approx {\omega }_0$ and then increase and intersect at $\omega \approx 10{\omega }_0$. We have $G_2<G_1$ for $\omega <10{\omega }_0$ and $G_1<G_2$ for $\omega >10{\omega }_0$. This suggests two critical times in the polymer chain dynamics. The first is ${\tau }_e={\tau }_0/10$ and the second is ${\tau }_0$ at which we find have a trend of slightly increasing with decreasing density of the systems as predicted in [61]. We find that with increasing frequency the response tends to become dissipative.

At low frequencies $G_1$∼${\omega }^{1/2}$ we find the trends follow the Rouse model. For larger frequencies we find that $G_1$∼${\omega }^{3/2}$, $G_2$∼${\omega }^{3/4}$ and then $G_1,G_2$ tend to a plateau value. Similar scalings were reported in [62] and in [47] for polymer solutions of linear FENE chains of similar molecular weight, which further shows that the use of harmonic bonds does not significantly influence our findings. We find a decreasing slope of $G_1$ for increasing entanglement which suggests a slower relaxation mechanism.

We show the loss tangent as a function of the frequency of oscillation in Figure 6. We remark that $tan\delta$ can be interpreted as reflecting the strength of what is sometimes called “colloidal forces”. In other words, if $tan\delta <1$ then the particles are highly associated and sedimentation could occur. If $tan\delta >1$, the particles are highly unassociated. The loss tangent is almost constant, close to 0, for all the simple infinite weaves (w0,wI,wII). The values of the open weaves are greater than one and then decrease to the values of the corresponding infinite weaves. The asymptotic ordering of the phase lag of the systems is $w0<wI<wII<wIII$.

We find that, our data at larger frequencies that all the materials behave like elastic solids, as is often seen in large frequency responses. The inset graph shows the corresponding log-log plot only for the open systems. It reveals a cross-over at approximately ${\omega }_0$, which corresponds to times on the order ${\tau }_0$. This time-scale could be related to the entanglement time as in [49,63]. This characteristic time-scale, seems to decrease with the topological complexity of the weave. The large frequency tail of $tan\delta$ decreases more slowly with increasing topological complexity and density indicating a substantial dissipation effect related to entanglement.

5.2. Conformational Analysis

We show configurations of the polymer weaves at different times during deformation in Figure 7. For the small oscillation frequencies, the infinite chains follow the deformation of the defining box, attaining an s-shape conformation. The open chains, significantly rearrange in time and tend to avoid the boundary by aligning with the orientation of the deformation. This process happens more slowly for wI and even more slowly for wII and wIII systems due to topological obstacles. We note that the chains tend to form bundles of chains, giving an inhomogeneous material, suggesting that the inhomogeneity decreases with increasing density and entanglement complexity. Similar phase separation of polymer solutions in oscillatory shear has been observed experimentally in [64].

We find the transition from bundle-dominated structures to entanglement dominated structures is related to the entanglement of the chains as has been also reported in [65]. A possibility for the bundle formation is finite-size effects. To examine that, we performed similar simulations in equilibrium in the NVE and in the NPT ensemble for the w0 infinite system. In both cases, bundle formation occured rapidly in the simulation. This indicates that the bundle formation is not a finite size effect.

We propose that the chains bundle together in order to decrease their deformation which occurs due to thermal fluctuations and due to the deformation of the cell. As the chains bundle, they form tubes of chains. A larger diameter tube resists the deformation stronger than the individual chains. This larger tube structure is apparent in the w0 infinite system. The same happens in two and three directions respectively for wI, wII and wIII. This also happens transiently at initial times for the open systems. At first the open chains form these bundles then they keep rearranging and entangling further until finally forming globules.

Our topological methods can be used to more precisely characterize the conformational features of the chains. We do this by measuring the Writhe and the Periodic Linking Number topological quantities that we introduced in Section 2. We remark that the chains in our system are loosely entangled relative to highly knotted systems yielding, as a consequence, Periodic Linking Numbers and Writhe that are less than one. While the quantities appear small relative to those of highly knotted systems, they still provide a useful characterization of the collective configurations of the polymeric chains of the material and their rearrangements. As our results indicate these are useful in understanding the connection between topology and mechanical responses even in systems of low molecular weight (below that of the entanglement length).

We find a very different behavior of the Writhe for the weave of type w0 between the open and infinite systems (the situation is similar for wI and wII). In the case of the infinite periodic systems, the mean absolute Writhe of the chains shows a sinusoidal behavior. This is seen most clearly when the frequency of oscillation is small. During the shearing cycle of the unit cell in our simulations, the Writhe is maximum when the shear deformation is the least and the Writhe reaches a minimum when the shear deformation is at its greatest. This behavior is indicative of the chains stretching at the maximum deformation and relaxing to a more entangled state when the shear deformation is relaxed.

In the case of open systems, the mean absolute Writhe of the chains also follows a sinusoidal behavior, but it changes significantly in time. This happens because the chains are free to attain any possible configuration and tend to disentangle and relax to configurations similar to those of random coils. Indeed, the final values are similar to those of a semiflexible random coil of comperable length as reported in [36]. This behavior becomes less pronounced as we increase the density and complexity of the weave because the disentanglement time increases and the chains do not have sufficient time to rearrange.

5.3. Topology and Mechanical Responses

We investigate both closed and open chain systems. The closed chain systems were found to be less interesting since the values of the Writhe and the Periodic Linking Number did not change significantly as a function of the frequency of oscillation. This is as expected, since without breakage of bonds the topology must remain close to that of the original configuration. Similarly, the loss tangent of the infinite systems does not change significantly throughout the experiment.

Here, we focus on the open chains systems where there is the potential for significant rearrangements of polymer chain configurations and topology. We show the mean absolute Writhe of open chains as a function of the loss tangent for small frequencies of oscillation in Figure 8.

We find a decay of the mean absolute Writhe with the loss tangent and a clustering of the data for each system. It is notable that the proposed Writhe measure groups together the systems of similar material response. The clustering observed for these systems indicates that the global entanglement of the chains imposed by the original conformation affects the response of the material significantly. The responses at small frequencies are clustered and ordered with their Writhe decreasing as $w0>wI>wII>wIII$.

We remark that as the frequency increases we found the Writhe of the open systems tended to meet the values of their corresponding closed (infinite) chain versions. This occurs since the chains cannot escape their original configuration as readily at large oscillation frequencies. For the $wIII$ case, approximately the same values for open and closed chains were found throughout.

Previous studies have shown that there is an almost linear relationship between the number of kinks per chain and the mean absolute Writhe of a chain in the case of a melt of linear FENE chains [36]. The relationship emerges as a function of the molecular weight in equilibrium conditions.

In non-equilibrium conditions, the relation becomes more complex [37]. The viscosity can be obtained as the limit $\eta ={lim}_{\omega \to 0}{G}_2\left(\omega \right)/\omega$. Because the frequencies studied here are relatively large, however, we can examine how the ratio $G_2/\omega$ depends on topology as a reference. Our results for the smallest frequencies indicate that the Writhe of the open chains decreases with $G_2/\omega$, suggesting a decrease of the Writhe with viscosity, see Figure 9. The number of kinks has been shown to increase with viscosity [9,10]. A similar inverse behavior between the Writhe of the chains and the number of kinks was observed in the initial time of an elongation of the chains [37].

We show the mean absolute Periodic Linking Number of the open systems at small frequencies as a function of the loss tangent in Figure 10. We find that the mean absolute Periodic Linking Number decreases with the loss tangent for $tan\delta >1$. We find significant clustering of the data corresponding to the different systems. The clustering of the responses indicates that the imposed global topology of the initial configuration significantly affects the response of the material. We remark that with increasing frequency we found the open systems tended to the values of the corresponding infinite systems. This occurs because, at large frequencies, the open chains cannot escape their original configurations.

The increase of the linking number implies the presence of persistent entanglements. Persistent entanglements are tight contacts between chains that significantly restrict their motion [66]. Such contacts are likely to cause significant bond stretching under deformation that is followed by a decrease of the Loss Tangent. These results indicate that interactions underlying mechanical responses can be effectively captured by the Periodic Linking Number.

From our results we obtain information about the equilibrium modulus. The equilibrium modulus is defined as the limit $G_{eq}={lim}_{\omega \to 0}{G}_1\left(\omega \right)$. We show the storage modulus for the smallest frequencies in Figure 9. Our results suggest that the equilibrium modulus increases with the linking of the chains. Interestingly, a similar linear relation between the entanglement density as obtained from the Gauss linking integral and the shear modulus for ring polymers was reported in [32].

Our results on the Writhe and Periodic Linking Number show a competing relation between these two with respect to the Loss Tangent. We provide a brief explanation for this effect. If the Writhe is large and the Periodic Linking Number is small, this can be interpreted as meaning the chains attain random conformations with no significant topological constraints. This would result in a behavior that is primarily dissipative. This suggests that inter-chain contributions to stress dominate through collisions of molecules induced by the Brownian motion.

In contrast when the Writhe is small and the Periodic Linking Number is large, we expect that the chains get stretched by the presence of persistent entanglements. In this case, intra-chain contributions to stress would dominate. As the ratio $⟨|L{K}_P|⟩/⟨\left|Wr\right|⟩$ increases, the persistence of entanglements increases. This implies that the bond stretching increases, which decreases the Loss Tangent. Therefore, we expect the loss tangent to increase with the decreasing ratio of Periodic Linking Number versus Writhe. In fact, this is confirmed in our results as seen in Figure 11.

Our results show, for the open systems at small frequencies, such a trend of the ratio of the Periodic Linking Number over the Writhe as a function of the Loss Tangent. This indicates that one can control the viscoelastic properties of a material by controlling the ratio of the Writhe and Periodic Linking Numbers of the constituent chains. This also suggests that $⟨|L{K}_P|⟩/⟨\left|Wr\right|⟩$ is a measure of the inter-chain contribution versus the intra-chain contribution to the stress. This finding may contribute to our understanding of the interplay between these two contributions [22].

6. Conclusions

We have introduced topological measures for polymeric materials allowing for quantifying the relationship between polymer chain entanglement and viscoelastic responses. We demonstrated how our topological approaches can be utilized in practice by performing non-equilibrium molecular dynamics simulations of a few polymeric systems. We found that our topological approaches provide a measure of the molecular-level polymer chain entanglements that contributes to aggregate mechanical responses in the storage modulus and loss modulus. Our results indicate some of the ways topological tools can be used for characterization of the interplay between entanglement and mechanics.

In particular, in this study, we have introduced new rigorous mathematical measures of topological complexity of open or closed curves in space that: (a) can be computed for linear and ring polymer chains in any setting (3-space or systems employing PBC), (b) are continuous functions of the chain coordinates and (c) correlate to the viscoelastic properties of polymeric materials. Our measures provide an effectively new “definition” of entanglement complexity of polymers. We showed our measures can be relevant to the rheological characteristics of polymers without relying too strongly on a particular polymer theory. For example, our methods did not require the definition of a tube or a primitive path for its application.

Basing our measures on well-defined mathematical principles allows for computing our measures for systems which are unentangled in the sense of tube model theories and can be used to understand the transition from entangled to unentangled states. Our measures also have the advantage of being computed directly by the chain coordinates and do not require any chain deformation (such as contour reducing algorithms, PPA or end-to-end or other type of closure of the open chains). We want to point out that our measures do not contradict the successful tube model theories, but we hope serve to compliment them by providing a mathematically rigorous basis from which to extend and improve further these theories.

We expect our methods also to be particularly helpful in the case of polymer melts or polydisperse systems with solutions of mixtures of different topological structures, such as rings and linear chains, and even for polymers of different architectures, like bottlebrush polymers [67,68,69]. Our methods could be particularly useful in handling in a common way the heterogeneity of such systems. Overall, we expect our methods to be broadly applicable in the further characterization and study of the responses of polymeric materials.

Author Contributions

Conceptualization, P.J.A., K.C.M. and E.P.; methodology, P.J.A., K.C.M. and E.P.; software, P.J.A., K.C.M. and E.P.; validation, P.J.A., K.C.M. and E.P.; formal analysis, P.J.A., K.C.M. and E.P.; investigation, P.J.A., K.C.M. and E.P.; resources, P.J.A., K.C.M. and E.P.; data curation, P.J.A., K.C.M. and E.P.; writing–original draft preparation, P.J.A., K.C.M. and E.P.; writing–review and editing, P.J.A., K.C.M. and E.P.; visualization, P.J.A., K.C.M. and E.P.; supervision, P.J.A., K.C.M. and E.P.; project administration, P.J.A., K.C.M. and E.P.; funding acquisition, P.J.A., K.C.M. and E.P.

Funding

This research was funded by research grant NSF CAREER DMS-0956210, NSF DMS - 1616353, and DOE ASCR CM4 DE-SC0009254. We also acknowledge UCSB Center for Scientific Computing NSF MRSEC (DMR-1121053) and UCSB MRL NSF CNS-0960316.

Acknowledgments

The author P.J.A. acknowledges support from research grant NSF CAREER DMS-0956210, NSF DMS-1616353, and DOE ASCR CM4 DE-SC0009254. We also acknowledge UCSB Center for Scientific Computing NSF MRSEC (DMR-1121053) and UCSB MRL NSF CNS-0960316.

Conflicts of Interest

The authors declare no conflict of interest.

Figures and Tables

Figure 1. The Gauss linking integral between two curves measures the degree to which they entangle, it is a real number that varies continuously with the chain coordinates for open chains, and it is an integer topological invariant for closed chains. Similarly, the Writhe (the Gauss linking integral over one curve) measures the degree to which a curve entangles with itself and varies continuously with chain coordinates.

Figure 2. We consider polymeric chains entangled with weave-like topologies. The weave0 (w0) denotes case with aligned chains, weaveI (wI) the case with smaller density of orthogonal and non-interlaced chains, weaveII (wII) the case with larger density of orthogonal and non-interlaced chains, and weaveIII (wIII) the case with alternating interlaced chains.

Figure 3. We show the entaglements of the polymer chains of wIII. The weave wIII has a topology with alternating interlaced chains. We see the chains in the x direction which alternatingly go over and under chains in the perpendicular y direction. We show one chain in the x-direction (orange curve) that can be seen locally to meet with three chains in the perpendicular y direction.

Figure 4. In Lees–Edwards boundary conditions, a monomer that exits the cell in a direction perpendicular to the shear direction will re-enter from the opposite face at position shifted by ±L(t)modL0.

Figure 5. Polymer Weave Frequency Response: Dynamic Moduli. The Elastic Storage and Viscous Loss Moduli of the infinite chain weaves are shown left and those of the open chain weaves right. The infinite weaves behave like crosslinked polymers with a primarily elastic behavior throughout the range of frequencies simulated. The exponents 1/2 and 1/4 are similar to those in the Rouse model [49]. The open weaves transition from an elastic to a viscous behavior as frequency is varied. The exponents 3/4 and 3/2 indicate the predicted scaling for semi-dilute solutions of semi-flexible chains and for the BEL model respectively [50]. The slope increases with decreasing topological complexity.

Figure 6. Polymer Weave Frequency Response: Loss Tangent. The infinite systems behave like crosslinked polymers with a loss tangent less than 1 at all frequencies. The open chains transition from a liquid-like behavior to that of a solid-like behavior as the frequency increases. The inset plot shows the log-log plot for open chains. These results show that the crossover frequency increases with decreasing topological complexity. Similarly, the slope of decrease increases with decreasing topological complexity.

Figure 7. Polymer configurations subject to oscillatory shear. Configurations at the end of simulations for the closed chain and open chain systems. In both cases, the chains tend to form bundles. By forming tubes of larger radius the chains decrease their individual deformation. For weave0, the open chains can significantly rearrange their conformations to those of random coils and the system becomes inhomogeneous, disconnected accross the periodic boundary. As the entanglement complexity and density increases, the open chains cannot fully escape their original conformations, forming tubular connected domains (weaveI), lamellar structures (weaveI) and even retain entanglement percolation in three dimensions (weaveIII).

Figure 8. Writhe and Loss Tangent of open chains for small frequencies of oscillation (frequencies corresponding to periods T>6τD). We find a linear behavior of the Writhe as a function of the Loss Tangent. The data points corresponding to the different weaves form clusters. As the topological complexity of the weave decreases both the Writhe and Loss Tangent increase.

Figure 9. (Left) Writhe and G2/ω for open weaves at small frequencies (log-log plot). The viscosity can be obtained as the limit η=limω→0G2(ω)/ω. We find an indication that the Writhe decreases with viscosity, while Z shows the opposite relation to viscosity [10]. (Right) Periodic Linking Number as a function of G1 for open chains for small frequencies (log-log plot). The equilibrium modulus can be obtained as the limit Geq=limω→0G1(ω). We find an indication that Geq increases with increasing linking, similar to what was reported for rings in [32].

Figure 10. Periodic Linking Number and Loss Tangent of open systems for small frequencies of oscillation (corresponding to periods T>6τD). We find that the responses corresponding to the different weaves form clusters. The Periodic Linking Number increases with weave complexity and decreases with loss tangent.

Figure 11. The ratio of the Periodic Linking Number over the Writhe as a function of the loss tangent for open systems at low frequencies. We see that the ratio ⟨|LKP|⟩/⟨|Wr|⟩ decreases with the loss tangent.

Table 1

Densities associated with the polymeric weaves shown in Figure 2.

Weave	Topology	Density	MW (open)
W0	parallel, non-interlaced	0.0625 (15 amu/nm${}^3$)	20 $m_0$
WI	orthogonal (non interl.)	0.1875 (45 amu/nm${}^3$)	20 $m_0$
WII	orthogonal (non interl.)	0.33 (80 amu/nm${}^3$)	15 $m_0$
WIII	alternating interlaced	0.35 (84 amu/nm${}^3$)	21–17 $m_0$

Table 2

Parameterization for the polymer weave models. These parametrizations (and Table 3), are chosen for relatively stiff polymers exhibiting worm-like chain responses for polymeric solutions and melts.

Parameter	Description	Value
$\sigma$	monomer radius	1.0 nm
$ϵ$	energy scale	2.5 amu·nm${}^2$/ps${}^2$
$m_0$	reference mass	1 amu
$w_c$	energy potential width	$2.5\sigma$
m	monomer mass	240 $m_0$
$\tau$	LJ-time-scale	$\sigma \sqrt{m_0/ϵ}$ = 0.6 ps
$k_{B}T$	thermal energy	1.0 $ϵ$
$\rho$	solvent mass density	39 $m_0/{\sigma }^3$
$\mu$	solvent viscosity	25 $m_0/\tau \sigma$
$\mathsf{{\rm Y}}$	drag coefficient	476 $m_0/\tau$

Table 3

Parameterization for the stiffness and connectivity of the polymer chains.

Parameter	Description	Value
$E_b$	harmonic bonds potential constant	619.5 amu/ps${}^2$
b	harmonic bonds rest length	1.0 nm
$E_{\theta }$	harmonic angle potential constant	19.8 amu·nm${}^2$/ps${}^2$
${\theta }_0$	harmonic angles’ rest length	180${}^{\circ }$

Table 4

Characteristic time scales.

Parameter	Description	value
${\tau }_A$ advection time	propagation in fluid	0.0013 ps
${\tau }_D$ diffusion time	monomer moves a dist. $\sigma$	302 ps
${\tau }_R$ Rouse time	ideal chain $N=20$	6937 ps
${\tau }_0$ critical time	cross-over reference time	598 ps
T	period of oscillation	6 ps < T < 3600 ps
t simulation time	longest simulation time	150 ps < t < 90,000 ps