Understanding changes to microstructural dynamics under nonlinear deformations is critical for designing flow processes of entangled polymeric fluids, motivating the development of experimental methods to probe strain- and rate- dependent modifications to relaxation mechanisms. Although orthogonal superposition rheometry (OSR) holds promise as such a probe, the ability to interpret the superposition moduli accessible by OSR in the context of entangled polymer dynamics remains an open question. To fill this gap, we report model OSR predictions using detailed microstructural models for both monodisperse and polydisperse entangled polymers, i.e., the Rolie-Poly and the Rolie-Double-Poly models, respectively, which account for reptation, chain retraction, and convective constraint release. By combining numerical calculations with a perturbation analysis, we demonstrate that for polymers that can be described by a single-mode model, the OSR superposition moduli at different shear rates and frequencies can generally be collapsed onto a single master curve, with rate-dependent shift factors that depend on the nonlinear rate-dependent modification of polymer conformation and relaxation rates without changing the dominant relaxation mechanisms. We systematically study how the OSR moduli are sensitive to the shape and dispersity of the molecular weight distribution. We discuss the generality of our results for a broad class of constitutive models and suggest an analogy to Laun’s rule to relate OSR moduli to the first normal stress difference. Our results provide a foundation to guide the design and interpretation of future experiments and demonstrate that orthogonal superposition rheometry often probes features in nonlinear dynamics more directly than conventional rheometry techniques.

Flow processing of soft materials and complex fluids usually involves nonlinear deformations, which can nontrivially modify both the associated microstructural configuration and dynamics. This is particularly important for entangled polymeric fluids, where entanglements produce nonlinear relaxation mechanisms that are significantly modified when polymers stretch and orient in flow, leading to various flow instabilities and complications in polymer processing. Various recent flow protocols have been proposed to highlight the nonlinear contributions to the relaxation: one is large amplitude oscillatory shear (LAOS) [1] and another is superposition rheology, which was originally proposed by Tanner and Simmons [2,3]. Although in LAOS the frequency dependence offers the possibility to probe different time scales, the method relies on subjecting the sample to a very complex kinematic history [4]. Alternatively, superposition rheology still enables frequency dependent measurements while retaining simpler kinematics. As a result, the superposition of a small-strain oscillatory motion onto a steady or transient shear flow can provide a clearer insight into the effects of flow on the mechanisms underlying the nonlinear response of rheologically complex fluids [5]. The oscillatory motion can be imposed either parallel (γ˙, parallel superposition, PSR) or perpendicular (γ˙, orthogonal superposition, OSR) to the direction of motion for the steady shear flow. Of these, orthogonal superposition offers the advantage that the two flow fields are not coupled and, as we show in this work, the corresponding moduli can be more directly related to the microstructural configuration under flow. OSR provides measurement of the frequency- and shear-rate-dependent nonlinear viscoelastic superposition moduli (storage and loss moduli, G and G), which are analogous to their linear viscoelastic (LVE) counterparts that characterize equilibrium dynamics. So, in principle, OSR is a very promising experimental method for studying the coupling of microstructural dynamics and nonlinear rheological response in polymeric liquids.

In practice, however, it has been difficult to achieve OSR measurements experimentally due to instrumental challenges and measurement sensitivity. The earlier devices for OSR by Simmons [3], Mewis and Schoukens [6], and Zeegers et al. [7] have a complex and delicate mechanical design and only cover a rather narrow viscosity range. Vermant et al. [8] described an orthogonal superposition technique based on a simple modification of the force rebalance transducer in a commercial rheometer and introduced the design of a double wall Couette cell with open bottom for the inner cylinder to avoid pumping flow. These new instrument developments finally enabled the implementation of sensitive orthogonal superposition on commercial rheometers. A recent experimental and simulation study thoroughly investigated the calibration procedures and corrections needed for the nonidealized flow field in the OSR setup on commercial rheometers [9]. Orthogonal superposition rheometry (OSR) has been used to study structural changes during flow in complex fluids, including polymeric fluids [10,11], starlike and wormlike micelles [5,12], and both in Brownian dynamics simulations [13] and experiments for colloidal suspensions [4,14–17].

Because OSR has only recently become experimentally available on commercial rheometers, relatively little is known about how to interpret the nonlinear viscoelastic results in the context of entangled polymer dynamics. Specifically, there are very few theoretical and computational studies that would provide a fundamental basis for interpreting orthogonal superposition measurements. To date, most of the experimental studies have interpreted the data by analogy to the interpretation of linear viscoelastic measurements. However, though proposed analogies between linear viscoelasticity and orthogonal viscoelasticity have been used (e.g., using empirical rules to estimate chain length [18–20]), these analogies have not been theoretically verified. This provides strong motivation for a deeper theoretical and computational study to investigate what is fundamentally being measured in the OSR experiments.

To better understand superposition experiments, the resulting moduli should be described by suitable rheological constitutive models. Lacking a generally valid nonlinear rheological constitutive equation, previous studies have used various phenomenological models. Yamamoto proposed a strain rate-dependent relaxation spectrum, which demonstrates that the parallel moduli reflect the coupling between shear and superimposed flow, while orthogonal superposition moduli do not [21]. Wong and Isayev predicted the orthogonal moduli with the Leonov model [22]; Kwon and Leonov later corrected the predictions [11]. The corrected predictions agree qualitatively well with the early experimental work by Simmons on polymer solutions [10]. Kim et al. used the Giesekus model to predict parallel and orthogonal superposition moduli and compared them with experiments on wormlike micelles [5].

For studying entangled polymers in superposition flows, only two molecularly based models have been used. The first involved a prediction of parallel superposition moduli by Unidad and Ianniruberto [23] using the differential constitutive equation accounting for convective constraint release (CCR) proposed by Marrucci and Ianniruberto: the double-convection-reptation model with chain stretch (DCR-CS model) [24]. However, the interpretation of parallel superposition results is complicated by the coupling of the steady shear flow and oscillatory flow. The second study by Mead [25] uses the monodisperse and polydisperse MLD models to predict PSR and OSR moduli [26–28]. A potential weakness of the MLD model is that the stretch and orientation dynamics are treated using separate dynamical equations, which can lead to problems even in monodisperse rheology predictions [29]. As pointed out by Boudara et al. [30], the decoupling approximation is even more problematic for polydisperse blends because couplings between constraint release and chain retraction must be readmitted in an ad hoc fashion and do not arise naturally from the model. Furthermore, the studies by both Mead as well as Unidad and Ianniruberto neglected chain stretch. So far, no studies have systematically studied the effect of polydispersity on OSR moduli using models that treat stretch and orientation in the same equation and probed a sufficiently wide range of shear rates and frequencies so that significant shifts in the moduli can be observed with increasing shear rate.

The current investigation addresses these limitations of previous studies by using molecularly based models that treat orientation and stretch in the same equation for predicting OSR moduli. Specifically, we use the Rolie-Poly model [31] for monodisperse polymers and the Rolie-Double-Poly model [30] for polydisperse polymers. These models account for the detailed nonlinear relaxation processes of entangled polymers including reptation, chain retraction, and convective constraint release. Furthermore, by combining a perturbation analysis and numerical calculations, we can obtain OSR moduli across a wide range of frequencies and shear rates and systematically vary the degree of chain polydispersity. After the detailed modeling study, we will discuss our results in a broader context to help design and interpret future OSR experiments. We find that the relationships governing OSR moduli for the Rolie-Poly model also work for a broader class of constitutive models. We suggest an analogy to Laun’s rule [32] to relate OSR moduli to the first normal stress difference. Additionally, we derive expressions to extract model parameters from experiments without the need for detailed model fitting. By making an analogy to the linear relaxation spectrum, we show the utility of OSR experiments to obtain material structural information that would be difficult to obtain from conventional rheological measurements. In this way, studying orthogonal superposition computationally not only provides an important comparison to future experiments but also provides better sensitivity for testing models for nonlinear polymer processing.

Storage and loss moduli have been used ubiquitously to describe viscoelastic materials, such as entangled polymers. The storage modulus G is proportional to the energy storage and, therefore, reflects the elasticity in the system. The loss modulus G is proportional to the dissipation or loss of energy, reflecting the viscous nature of the system [33]. The most common approach to obtain information about the deformation-induced structural rearrangement within a complex fluid is to measure the two moduli using small amplitude oscillatory shear (SAOS). In this regime, the material response is independent of the strain amplitude, and this is known as the linear viscoelastic behavior [34].

The storage and loss moduli at equilibrium can be used to estimate the average contour length of entangled polymers and wormlike micelles, as demonstrated in previous works [18]. Specifically, based on the polymer network theory, the number of entanglements of a chain is given by the ratio of the plateau modulus GN0 and the minimum of the loss modulus Gmin. The entanglement length, le, is evaluated in the context of a particular chain model (e.g., a flexible chain or a semiflexible chain) [19,20]. Thus, one can estimate the average chain length, Lc, as follows:

Lc=leGN0Gmin.
(1)

Before discussing the model for polydisperse polymers, we first take a moment to review the classic Rolie-Poly model of Likhtman and Graham [31], which has its origin in the full GLaMM model. The GLaMM model was originally proposed by Graham et al. [35] as a full chain molecular theory for entangled monodisperse linear polymer chains under fast deformation. The GLaMM model accounts for an accumulation of stress via affine deformation and relaxation of stress by reptation, chain retraction, and convective constraint release (CCR). Reptation was included in the original tube model introduced by Doi and Edwards [19] and describes the curvilinear motion of a polymer chain in a tube formed by entanglements with surrounding chains. When the deformation becomes much faster than the inverse of the Rouse time, the chain can become stretched, and as a result, the length of the chain and the occupied tube exceed their equilibrium configuration. When the strain stops, the chain retracts along the deformed tube until it regains its equilibrium contour. Convective constraint release (CCR) is the release of entanglement constraints due to chain retraction from affine deformation in nonlinear flows. The effect of CCR diminishes for shear rates larger than the inverse of the Rouse time.

Although the GLaMM model is successful in predicting the rheology of fast flows, it is computationally prohibitive for nonviscometric flow calculations and for OSR calculations that span several orders of magnitude of frequency and shear rate. Hence, as a simplified version of the GLaMM model, Likhtman and Graham derived a one-mode approximation of the GLaMM model to produce a differential constitutive model for entangled monodisperse polymer chains: the Rolie-Poly model (for Rouse linear entangled polymers) [31]. The time evolution of the contribution to the stress tensor due to the polymer, σ, is calculated via the conformation tensor for an entanglement segment, A, as

σ=GN0A
(2)

and

dAdt=κA+AκTconvection1τd(AI)reptation2(1λ1)τRAretractionβCCR2(1λ1)τRλ2ρ(AI)convectiveconstraintrelease.
(3)

Here, GN0 is the plateau modulus, κ is the velocity gradient tensor, and the stretch ratio of the polymer chains, i.e., the ratio of the current to equilibrium-with-no-flow chain contour lengths, is defined as

λ=trA3.
(4)

The reptation term (with τd) is a single time relaxation toward equilibrium. It is also possible to include the effect of contour length fluctuation (CLF) with the reptation term by modifying the reptation time [30,36]. However, as shown by Boudara et al. [30], the inclusion or absence of CLF only affects the quantitative agreement of the model with experiment without affecting the qualitative trends. To keep the physical picture simple, we will not include CLF in our calculations. The CCR term is also a relaxation toward equilibrium but with the rate dependent on the amount of stretch λ and the Rouse time τR. Retraction is relaxation to zero stress, also with the rate dependent on λ and τR. βCCR specifies the CCR contribution and is analogous to cν in the full GLaMM model. cν determines the number of retraction events necessary to result in one tube hop of a tube diameter with a value of cν1. [35] To fit the GLaMM model to experiments, the authors used cν=0.1, which approximately corresponds to βCCR=1. ρ is an additional fitting parameter. βCCR=1, ρ=0.5 are the optimal values to fit the Rolie-Poly model prediction to the full GLaMM theory, thus these are generally used in the Rolie-Poly model and also in this work.

The Rolie-Poly model incorporates changes in the conformation due to orientation and stretch in the same equation. Thus, the model avoids anomalous shear thickening behavior at shear rates corresponding to chain stretch predicted by models that decouple stretching and orientation [29]. Furthermore, since Rolie-Poly model is a single mode approximation of the full GLaMM model, it has sound microstructural basis while still making the calculation of OSR moduli computationally feasible.

To develop a constitutive model for polydisperse entangled polymer, Boudara and co-workers used ideas from both the Rolie-Poly model [31] and the double reptation approximation [37] to formulate the RDP model [30]. The “RDP model” is used to refer the “Rolie-Double-Poly” model, where “Double” signifies both “double reptation” and “double poly” as in “POLYdisperse POLYmers.” The RDP model can be formally derived as a greatly simplified approximation of the detailed molecular theory for bidisperse polymer melts of Read and co-workers [38].

In the RDP model, the total polymer contribution to the stress is the sum of the contributions coming from each species i, weighted by their volume fraction, ϕi,

σ=GN0i=1NϕiAi,
(5)

with

λi=trAi3.
(6)

The conformation tensor Ai accounts for the stresses that come from the interaction of the species i with itself and the other n1 species,

Ai=j=1NϕjAij,
(7)

where Aij is the stress conformation tensor on the i-chains coming from their entanglements with the j-chains,

dAijdt=κAij+AijκTconvection12(1τd,i(AijI)reptation+βthτd,j(AijI)constraintrelease)2(1λi1)τR,iAijretractionβCCR2(1λj1)τR,jλi2ρ(AijI)convectiveconstraintrelease.
(8)

In the original Rolie-Poly model, τd specifies the terminal stress relaxation time of the polymer. However, in the RDP model, the effects of reptation are explicitly separated out from thermal constraint release (CR) given by the parameter βth. Due to the interaction of different chains, reptation of a smaller chain causes the constraint release of a longer chain. Thus, CR needs to appear in a separate term. To ensure that the RDP model converges to the Rolie-Poly model in the monodisperse limit, the reptation term and the CR term are divided by two. This choice is physically sound because the reptation term in the Rolie-Poly model is the effective reptation term that actually includes both reptation and CR.

1. Monodisperse calculation

To evaluate the orthogonal superposition moduli, it is necessary to calculate the transient stress response under a combination of axial steady shear flow and orthogonal oscillation. The velocity gradient tensor is

κ=[0γ˙00000γ˙(t)0],
(9)

where γ˙ is the steady state shear rate and γ˙(t)=γ0,ωcosωt is the oscillatory superposition term. In this work, we kept γ0,=0.05 to remain in the linear viscoelastic region, as determined in the strain sweep (see supplementary material Fig. S1 [39]). The six chain conformation tensor components of Eq. (3) need to be calculated numerically by transiently solving the coupled differential equations. First, the equations are nondimensionlized using the Rouse time τR and the plateau modulus GN0. Specifically, the Deborah number De is defined as the dimensionless frequency (De=ωτR=ω~), the Weissenberg number is the dimensionless shear rate (Wi=γ˙τR), the dimensionless time is t=t/τR, and nondimensionalized stress tensor is σ~=σ/GN0. Furthermore, the ratio of the reptation time and the Rouse time is τd/τR=θ=3Z, where Z is the number of entanglement points per chain. In this work for the monodisperse calculations, we take Z=20 for a well entangled polymer system.

The six transient equations for the conformation tensor components of Eq. (3) were solved using the MATLAB built-in solvers “ode45” and “ode15s.” “ode15s” is applied for the low frequency range to take advantage of the significant acceleration of the computation using adaptive mesh refinement. However, for the high frequency range, the adaptive mesh refinement caused large fluctuations in the moduli for adjacent frequency values. Consequently, for the high frequency range the very robust solver “ode45” is applied, which allows a user-defined time discretization. For this work, a logarithmic time spacing enables a sufficient sampling of the oscillations, across the entire frequency range simulated (De=104 to 100).

The shear stress in the orthogonal direction () is used for the determination of G and G using

GGN0=σ~0,γ0,1+(tanδ)2,
(10)
GGN0=σ~0,tanδγ0,1+(tanδ)2+ηsω~,
(11)

where σ~0, and γ0, are the amplitudes of the dimensionless shear stress and the applied strain, whereas δ defines their phase shift. Furthermore, the term ηsω~ denotes the solvent contribution, which is purely viscous and, therefore, only appears in G, with the scaled solvent viscosity ηs=0.001 in this work.

To determine the values of δ and σ~0, from the simulation, we first need to identify the start of the steady state in the main shear flow. The shear stress σ~xy is plotted as a function of time [Fig. 1(a)]. Within the plateau, the stress continues oscillating with a very small amplitude due to the superposed oscillation. We specify the plateau criterion as the smaller value of 103 and 102ω~ to account for the dependence of the amplitude of the stress oscillation on the frequency. The start of the plateau is defined as the first point in the transient shear stress with a moving standard deviation that is smaller than the specified plateau criterion. Once the time for the start of the plateau for σ~xy is specified, the same time is used to identify the start of steady state for the stress in the orthogonal direction (σ~yz), which is plotted in Fig. 1(b) along with the applied strain rate (γ˙). The gray section of the figure illustrates the transient start-up before the harmonic oscillation. After neglecting the start-up phase, σ~yz and γ˙ are fitted by sine functions. Since the strain rate is a cosine function, it needs to be shifted by π/2 for the fitting of a sine function. Equations (12) and (13) define the two functions, including their fitting parameters c1, and c2,

σ~yz=σ~0,sin(ω~t+c1),
(12)
γ˙=γ0,ω~sin(ω~t+c2+π2).
(13)

In the equations, σ~0, and γ0, specify the amplitudes of the orthogonal stress and the applied strain in the orthogonal direction, respectively. Furthermore, the difference of the resulting coefficients defines the phase shift between the two wave functions: δ=c1c2. With the phase shift and the amplitude values, the OSR moduli G and G are calculated for each combination of shear flow strength Wi and superposed oscillation frequency De according to Eqs. (10) and (11).

FIG. 1.

Schematic showing the steps for calculating superposition moduli numerically (for Wi=1 and ω~=1). (a) Identification of steady state from shear stress response as a function of time. (b) Applied strain rate and yzcomponent of stress response as a function of time, the inset shows the phase shift.

FIG. 1.

Schematic showing the steps for calculating superposition moduli numerically (for Wi=1 and ω~=1). (a) Identification of steady state from shear stress response as a function of time. (b) Applied strain rate and yzcomponent of stress response as a function of time, the inset shows the phase shift.

Close modal

Finally, before using the numerical method above to make predictions of G and G, we verified that γ0,=0.05 is small enough to be in the linear viscoelastic region by conducting a strain sweep at a fixed oscillation frequency [see supplementary material Fig. S1(a) [39]]. Furthermore, we validated the numerical method by comparing steady state predictions with predictions reported with the original Rolie-Poly model formulation [31] and comparing orthogonal superposition predictions using a Giesekus model with analytical results from the work by Kim and co-workers [5] (see supplementary material Fig. S2 [39]).

2. Polydisperse calculation

The G and G calculations for polydisperse polymers follow the same workflow as the monodisperse calculation. The only difference is that several different conformation tensors Aij need to be calculated before the summation, given by Eqs. (5) and (7), is computed. For the nondimensionalization, the Rouse time (τR) of the weight-average chain length (Z¯=20) is used to define the dimensionless shear rate Wi, the dimensionless frequency ω~, and the dimensionless time t. Furthermore, the continuous length distribution must be discretized into different species, each with its own chain length and volume fraction. However, the computational expense significantly increases with the number of species by a factor of N2. For this study, a maximum N=32 is used. To systematically study the effect of polydispersity on the OSR moduli, we assumed a lognormal distribution for the volume fraction distribution with varying degrees of polydispersity while keeping the weight-average molecular weight (Mw) the same. The lognormal distribution is frequently seen in polymers used for rheological studies [40]. The continuous distribution is defined as

ϕ(Z)=1Zs2πexp[(ln(Z)μ)22s2],
(14)

with s2=ln(Ip) and μ=ln(Z¯/Ip). ϕ(Z)δZ gives the volume fraction of chains with entanglement number falling in the range Z±δZ/2 for asymptotically small δZ. We also define Z¯=ϕZdZ as the mean number of entanglements per chain and Ip=Mw/Mn=Z¯ϕ(Z)/ZdZ as the polydispersity index. We then partition the continuous distribution into N sections, where the i=1,2,3,N section is bracketed by chains with entanglement numbers Z^i1 and Z^i. In the approach used in this work and first developed by Peterson [41,42], the partitions are drawn at lines of equal volume fraction. In our view, this is the simplest discretization scheme that guarantees convergence to the correct cumulative distribution when N. Written in terms of the cumulative distribution function Φ0(Z), the distribution scheme is formally given by

Φ0(Z)=0Zϕ(Z)dZ=12[1+erf[ln(Z)μ)s2]],
(15)
Φ0(Z^i)=i/N.
(16)

For the lognormal distribution, Φ0(Z) can be inverted to obtain the set of Z^i,

Z^i=exp[μ+s2erf1[2iN1]].
(17)

Now defining Zi as the weight average molecular weight of chains in the range between Z^i1 and Z^i, we approximate the continuous distribution by a discrete distribution of N components with entanglement numbers Z1,Z2,Z3,,ZN and equal volume fractions ϕ1=ϕ2=ϕ3==1/N. This strategy ensures that the mean entanglement number Z¯=iϕiZi is not changed by discretization. The end result for the discretization scheme is

Zi=12Nexp[μ+12s2][erf(erf1(2iN1)s2)erf(erf1(2(i1)N1)s2)].
(18)

For orthogonal superposition, the oscillatory flow is assumed to be in the linear viscoelastic regime. Thus, the strain amplitude γ0,1 and γ0,ω~Wi. Wi is the shear rate of the mean flow nondimensionalized by the Rouse time τR. In this case, we can use a perturbation analysis with a small parameter ϵ being the ratio of the magnitudes of the strain amplitudes of the orthogonal and main flow: ϵ=γ0,|γ˙|. We shall assume a regular asymptotic solution structure and retain only the first order perturbation to the steady state solution. The superscript “(0)” denotes the constant, order-one steady state solution and the superscript “(1)” denotes the order-ϵ perturbation to the steady state value. With these notations, the stress tensor components and the chain stretch can be expressed as an expansion about the steady state levels: σij(t)=σij(0)+ϵσij(1)(t) and λ(t)=λ(0)+ϵλ(1)(t). Furthermore, as in the numerical calculation, we nondimensionalize the stress by the plateau modulus GN0: σ~ij=σij/GN0, which also equals the conformation tensor components Aij according to Eq. (2). Thus, we have Aij(t)=Aij(0)+ϵAij(1)(t). We will first write the mathematical framework for the perturbation analysis using the monodisperse Rolie-Poly model. A similar procedure is used for the RDP model. For all six equations of the stress components in Eq. (3), we collect the O(1) terms (0) and the O(ϵ) terms (1). For the purpose of calculating the OSR moduli, we only need to solve the O(ϵ) equation for Ayz(1),

dAyz(1)dt=Ayy(0)γ0ω~cos(ω~t)(1θ+2(11λ(0))+2(11λ(0))λ(0))Ayz(1).
(19)

We can write this equation completely in terms of λ(0) by obtaining the steady state solution for Ayy(0),

Ayy(0)=λ(0)+2(11λ(0))θλ(0)+2(11λ(0))θ+2(λ(0)1)θ.
(20)

The equation for Ayz(1) can then be rewritten as

dAyz(1)dt+BAyz(1)=Ayy(0)γ0ω~cos(ω~t),
(21)

with

B=1θ+2(11λ(0))+2(11λ(0))λ(0).
(22)

Equation (21) can be solved analytically for the OSR storage and loss moduli,

GGN0=Ayy(0)ω~2B2+ω~2,
(23)
GGN0=Ayy(0)Bω~B2+ω~2.
(24)

Here, only the polymeric contribution is included. To include the solvent contribution, we can add ηsω~ as in Eq. (11). In Eqs. (23) and (24), Ayy(0)(γ˙) and B(γ˙) are evaluated numerically from steady state Rolie-Poly equations without superposition. The values are then used in the equations to calculate G and G analytically. The perturbation analysis gives a speedup of more than 100 times while offering the same accuracy compared to the numerical calculation of OSR moduli (see supplementary material Fig. S3 [39] for a comparison of the two results).

A similar procedure is used for the perturbation analysis using the RDP model for polydisperse polymers. In this case, Aijyy(0)(γ˙) and Bij(γ˙) of each entanglement pair need to be determined from the steady state calculation. The OSR moduli are then sums of contributions from all entanglement pairs,

GGN0=i=1nj=1nϕiϕjAijyy(0)ω~2Bij2+ω~2,
(25)
GGN0=i=1nj=1nϕiϕjAijyy(0)Bijω~Bij2+ω~2.
(26)

At this point, we reflect on the way in which we derived these relationships and their potential generality. Results for the orthogonal moduli in terms of the yy-component of the polymer conformation tensor come about from the convection of the normal stress through the upper-convected time derivative to the yz-component of the stress, which is the direction of the oscillation. Regardless of the constitutive model that is used, this is the procedure and the mathematical propagation that gives us a relationship similar to Eqs. (25) and (26). In fact, the linear proportionality between the plateau value of G and Aijyy(0)(γ˙) not only holds for Rolie-Poly and RDP models but also for any model with the upper-convected time derivative of the conformation tensor A and stress relaxation terms that are linear in the conformation tensor and any linear superposition of these models. The results are valid, regardless of how the chain conformation is quantified, i.e., whether the end-to-end vector or the orientation tensor is used. A similar analysis on the DCR-CS model by Marrucci and Ianniruberto [24] (see supplementary material Sec. V [39]) and the MLD model by Mead [25] produces the same relationship.

We first present numerical calculations to obtain the orthogonal superposition moduli using the Rolie-Poly model. Figure 2 shows the computational results of a frequency De sweep for the oscillatory shear at different shear rates of the steady flow, i.e., Weissenberg numbers, Wi=γ˙τR, ranging from 1 to 100. Initial studies of OSR on entangled polymers were limited to relatively low Wi, but recent advances in instrumentation have extended the accessible range to cover the range of Wi we study. These results are for a monodisperse polymer with an entanglement number of Z=20, which is the same as what accompanied the original Rolie-Poly model formulation [31]. Consequently, the ratio of the two relaxation time scales θ=τd/τR=60. We find two distinctive trends in the moduli behavior with increasing Wi. First of all, a horizontal shift is observed, which can be seen best in the rightwards shifting of the GG crossover point. Second, G decreases monotonically with increasing main flow shear rates Wi. This effect can be seen best by the plateau values of the G curves. Comparing G and G at Wi=0.001 to the linear viscoelastic moduli G and G in the absence of flow [Fig. S1(b) in supplementary materials [39]] verifies that results at low flow rates approximate the result at equilibrium.

FIG. 2.

Orthogonal superposition moduli for monodisperse polymers at six different shear rates (Wi).

FIG. 2.

Orthogonal superposition moduli for monodisperse polymers at six different shear rates (Wi).

Close modal

Furthermore, evaluating the contribution of each stress relaxation mechanism (reptation, retraction and CCR) in the yz-component of the stress, σ~yz, provides a better understanding of the causes of the shifts in the OSR moduli. Figure 3 shows the fractional contribution of relaxed orthogonal stress terms separated by relaxation mechanism, at a frequency of De=0.001. The shape of the plots is independent of the choice of the frequency. As Wi increases, contributions of the stretch-based phenomena (retraction and CCR) increase, with CCR developing a maximum and then decreasing. In contrast, reptation starts at 100% for very low shear rates and decreases with Wi to a value close to 0%. This is because the stress relaxed by reptation approaches a constant when the inverse of the disentanglement time τd is reached. Consequently, the reptation contribution in the Rolie-Poly model remains constant, while retraction and CCR start dominating relaxation for increasing Wi. In the first part of their increase, the same slope is observed. However, as Wi continues to increase the chains can retract only partially from affine stretch, and convective constraint release is less effective and the relative magnitude of the CCR contribution decreases. At the same frequencies, the retraction term develops a shoulder. Finally, at sufficiently large Wi, retraction dominates the relaxed stress with a value close to 100%.

FIG. 3.

Fractional stress contribution of different relaxation mechanisms in yz-component of polymer stress.

FIG. 3.

Fractional stress contribution of different relaxation mechanisms in yz-component of polymer stress.

Close modal

By examining the numerical results alone, it is difficult to pinpoint what physical phenomena directly cause the horizontal and vertical shifts. If these were purely linear viscoelastic measurements under different temperatures and concentrations, researchers would typically perform a time-temperature superposition or time-concentration superposition and force all of these curves to collapse. Here, we show that we can achieve such a collapse with shift factors. Furthermore, these shift factors can be predicted based on the perturbation analysis results in Eqs. (23) and (24). By assuming the superposed oscillation is a small perturbation to the main shear flow (γ0,1 and γ0,DeWi), we obtained analytical expressions for G and G in Eqs. (23) and (24). Using the two parameters Ayy(0)(γ˙) and B(γ˙) identified in the analytical expressions, we scaled the numerical results from Fig. 2 after subtracting out the solvent contribution (ηsω~). This procedure completely collapses all the Wi-dependent moduli onto a single set of curves. Specifically, we divided the OSR moduli by Ayy(0)(γ˙) and divided the dimensionless frequency by B(γ˙) [Fig. 4(a)]. The collapse of the curves suggests that the steady shear flow only speeds up the relaxation and decreases the effective oscillatory driving force, without changing the relaxation mechanism.

FIG. 4.

(a) Rescaling superposition moduli by Ayy(0) and frequency by B collapses superposition moduli at different shear rates onto one set of curves. (b) Changes in plateau modulus (Ayy(0)) as a function of Wi. (c) Changes in crossover frequency as a function of Wi.

FIG. 4.

(a) Rescaling superposition moduli by Ayy(0) and frequency by B collapses superposition moduli at different shear rates onto one set of curves. (b) Changes in plateau modulus (Ayy(0)) as a function of Wi. (c) Changes in crossover frequency as a function of Wi.

Close modal

We can further quantify the magnitude of the vertical shift by plotting the plateau value of dimensionless G as a function of Wi [Fig. 4(b)]. In fact, in the plateau, we have limω~1GGN0=Ayy(0)(γ˙). With increasing Wi, G monotonically decreases. Around Wi=1, a shoulder shape appears in the plateau [also in Ayy(0)(γ˙)], similar to the shoulder shape in the retraction contribution in Fig. 3. This shoulder again arises due to the competition of retraction and convective constraint release. Since Wi is defined in terms of the Rouse time τR, Wi=1 represents the start of significant chain stretch, which results in increased contribution of the stress relaxation by chain retraction. We can also quantify the horizontal shift by plotting the frequency at the crossover of G and G as a function of Wi [Fig. 4(c)]. With increasing Wi, the crossover frequency monotonically increases. In linear viscoelastic measurements, the inverse of the crossover frequency gives a characteristic relaxation time of the system: τ=1/ωc. So, Fig. 4(c) indicates the effective relaxation time of the polymers is decreasing with increasing Wi. By examining Eqs. (23) and (24), we find the crossover frequency equals B. To get a better understanding of what Ayy(0)(γ˙) and B(γ˙) represent physically, we can return the G and G results in Eqs. (23) and (24) to dimensional form,

G=σyy(0)(τeffω)21+(τeffω)2andG=σyy(0)τeffω1+(τeffω)2,
(27)

with

τeff=τd1+τdτR(22λ2).
(28)

The effective relaxation time τeff represents the speedup of the stress relaxation in the flow. Its asymptotic value changes from τeff=τd at low Wi before the onset of chain stretch to τeff=τR/2 for very large chain stretch. So, Ayy(0)(γ˙) represents the nonlinear shear-dependent chain conformation in the yy-direction, which is orthogonal to the main flow direction, and B(γ˙) represents the nonlinear shear-dependent relaxation rate.

As discussed in Sec. II A, for conventional linear viscoelastic measurements, the plateau modulus GN0 is used to estimate the entanglement length while the ratio of the plateau modulus and Gmin is used to estimate the number of entanglements. Thus, LVE G and G moduli are used to estimate average chain length, and previous investigators have assumed that the corresponding OSR moduli can be interpreted in the same way, thus providing a way to explore the effect of a steady shear flow on chain scission for wormlike micelles and other “living” polymers [5,43]. However, we see that the decrease in G with increasing Wi in Fig. 4(b) and the decrease in the ratio of the plateau and Gmin in Fig. 2 occur without any modification to the properties of the chain; for example, those that might be extracted from a scaling analysis of the linear viscoelasticity (i.e., the number of entanglements is fixed at 20 per chain). The plateau value and the crossover frequency change due to chain stretch and orientation, as illustrated in the contribution of stress relaxation in Fig. 3. Therefore, changes in the plateau modulus and relaxation time cannot be interpreted simply as due to changes in the scaling variables (e.g., entanglement length, entanglement density, etc.). This is consistent with a recent Brownian dynamics simulation for star polymers [44], in which Metri and Briels showed that the shear relaxation modulus from orthogonal superposition is not the same as the stress relaxation modulus from step strain (as is the case for linear viscoelastic experiment at equilibrium), even for weak shear flow in the OSR experiment.

Although the data interpretation from the OSR measurement is complicated even for a monodisperse polymer, the perturbation analysis nevertheless offers physical insights to pinpoint the material properties that cause the change in the OSR moduli. Before moving to the polydisperse calculations, we first demonstrate the generality of the perturbation analysis by applying it to a different model, the DCR-CS (double-convection-reptation model with chain stretch) by Marrucci and Ianniruberto [24]. As shown in supplementary material Sec. V [39], the resulting expressions for the OSR moduli are the same as those from the Rolie-Poly model [Eqs. (23) and (24)]. The differences are in the expressions for Ayy(0)(γ˙) and B(γ˙). Even though these expressions are different, the similarity in the functional dependence of G and G on Ayy(0)(γ˙) and B(γ˙) indicates that the dependence of the plateau shift on the yy-component of the chain conformation and the horizontal shift on the effective relaxation time are model-independent and represent the physical phenomena that determine the OSR moduli for linear polymers of fixed chain length.

To systematically study the effect of polydispersity on the OSR moduli, we used a discretized lognormal molecular weight distribution, as discussed in Sec. III A 2. Specifically, we used 32 bins to represent the continuous MW distribution with equal volume fraction in each bin. The polydispersity index was varied from 1.01, which is nearly monodisperse, to 5, which has a long tail. The weight-average molecular weight (weight-average number of entanglements per chain, Z¯) is fixed at 20. Figure 5(a) plots all of the molecular weight distributions used in this work. First, we compare the numerical calculations at two representative cases of the degree of polydispersity: PDI=1.01 and PDI=2 in Figs. 5(b) and 5(c), respectively. For PDI=1.01, the results are nearly identical to the monodisperse results in Fig. 2. For PDI=2, the effects of polydispersity are more evident. Specifically, the shoulder in the storage modulus becomes broader at all Wi and the plateau value decreases less. The gradual evolution of the polydisperse moduli reflects the broad relaxation spectrum in the highly polydisperse system (i.e., many relaxation times). As we will show below, the weakened softening of the plateau for the polydisperse case compared to the monodisperse case arises from the nonlinear dependence of Aijyy(0) on the chain length Zi. The presence of the short chains causes the plateau to decrease less than the monodisperse case.

FIG. 5.

Polydisperse calculations: (a) Lognormal molecular weight distributions used in this work. (b) OSR moduli for PDI=1.01 (c) OSR moduli for PDI=2.

FIG. 5.

Polydisperse calculations: (a) Lognormal molecular weight distributions used in this work. (b) OSR moduli for PDI=1.01 (c) OSR moduli for PDI=2.

Close modal

To better quantify the horizontal and vertical shifts in the OSR moduli, we can plot the plateau value and the location of the apparent crossover frequency as a function of Wi for different PDI. However, having to account for N2 number of chain pairs for all combinations of De and Wi significantly increases the computational time for the polydisperse case compared to the monodisperse case. Thus, once again, we turn to a perturbation analysis to obtain analytical expressions for the OSR moduli, which makes prediction of OSR moduli at finely resolved Wi and De possible. The procedure for the perturbation analysis and the main results are described in Sec. III B. We verified that the perturbation analysis produced identical results as the full numerical calculations. So, subsequently, we only report results from the perturbation analysis.

Figure 6 quantifies changes in the OSR moduli and crossover frequency as a function of Wi for varying degrees of polydispersity: PDI=1.01,1.05,1.2,2,5 and for the monodisperse result (shown as the solid black curve). The crossover frequency G and G increases monotonically with increasing Wi for all cases except for PDI=5 [Fig. 6(a)]. For PDI=5 and Wi>10, the decrease of the crossover frequency at large Wi is a result of significant stretching of the longest chains in the system. Such extreme stretching is unphysical for finitely extensible chains and can be removed with the inclusion of a Warner spring [51] (see supplementary material Fig. S7(a) [39]). Figure S7 compares the linear spring calculation with the Warner spring calculation. Subfigures (a)–(d) are organized similarly to those in Fig. 6 here. Predictions using the Warner spring only start to differ significantly from predictions for the linear spring for Wi>10 and PDI2. On the one hand, across most of the parameter space (spanning different shear rates and PDI) that is accessible experimentally, the predictions for the linear and Warner spring are very similar. On the other hand, the Warner spring adds extra complications associated with an additional parameter λmax, the maximum chain stretch. Thus, we will only discuss results for the linear spring in the main paper. At low dispersity, the moduli value at the crossover decreases monotonically with Wi [Fig. 6(b)]. However, for PDI=2 and 5, the moduli at the crossover first increase then decrease with Wi. The G plateau decreases monotonically with increasing Wi for all cases [Fig. 6(c)]. With increasing PDI, the curves deviate more and more from the monodisperse curve. Additionally, the shoulder shape in the monodisperse result that arises from the competition between retraction and CCR becomes smoothed out with increasing polydispersity. Interestingly, the monodisperse result for the moduli at the crossover has the same shape as the monodisperse plateau value, with the former being half of the magnitude. This effect is better illustrated in Fig. 6(d), in which the ratio of the crossover moduli and the plateau moduli is plotted as a function of Wi for the monodisperse case and all polydisperse cases. The monodisperse case is a flat line at 0.5 for all Wi which indicates that for a monodisperse polymer, the moduli value at the crossover is exactly half of the plateau value for all Wi. Examining the perturbation result for the monodisperse case [Eqs. (23) and (24)] illustrates this point more clearly. At the crossover, G=G and B(γ˙)=ω~. Thus, the moduli value at the crossover equals Ayy(0)(γ˙)/2. However, the Wi-independent, constant moduli ratio is not preserved in the polydisperse cases. With increasing polydispersity, the ratio decreases and becomes more nonmonotonic [Fig. 6(d)].

FIG. 6.

Quantify changes in OSR moduli as a function of Wi for varying degrees of polydispersity: (a) crossover frequency, (b) moduli at crossover, (c) plateau modulus, (d) ratio of moduli at crossover and plateau.

FIG. 6.

Quantify changes in OSR moduli as a function of Wi for varying degrees of polydispersity: (a) crossover frequency, (b) moduli at crossover, (c) plateau modulus, (d) ratio of moduli at crossover and plateau.

Close modal

The reason for all the changes in the polydisperse results compared to the monodisperse result, as illustrated in Figs. 6(a)6(d), arises from the fact that the OSR moduli for the polydisperse case are summations of the contributions from all entanglement pairs and the OSR moduli are affected by the interactions between chains of different length [Eqs. (25) and (26)]. In particular, instead of being solely determined by the Ayy(0)(γ˙) component of a single chain length as in the monodisperse case, the polydisperse plateau modulus is determined by i=1nj=1nϕiϕjAijyy(0)(γ˙)=Ayy(0)(γ˙), which represents the dimensionless weight-average chain conformation in the yy-direction of the entire population of chains. However, a weight-average of Bij(γ˙) for all entanglement pairs does not yield the overall effective relaxation rate (i.e., the crossover frequency) for the polydisperse case. This is because Bij(γ˙) appears in the denominator of the individual moduli terms that are summed together. Thus, for polydisperse systems, although there is a qualitative correlation of the data with the degree of polydispersity, as shown in Fig. 6 especially in (d) for the ratio of the moduli, no simple scaling exists to explain the shifts in the moduli, as in the monodisperse case.

One may wonder to what extent the changes in the OSR moduli in the polydisperse case can be understood by simple superposition of results for the discrete chains as opposed to needing the full double reptation physics included in the RDP model. Here, we make a comparison of the RDP results for PDI=2 shown in Fig. 6 and the results of the simple superposition of the monodisperse results for each discrete species in the polydisperse calculation. Thus, for the linear superposition case, we have GGN0=i=1nϕiAiyy(0)ω~2Bi2+ω~2 and GGN0=i=1nϕiAiyy(0)Biω~Bi2+ω~2, where Aiyy(0)(γ˙) and Bi(γ˙) are calculated from the perturbation analysis for the monodisperse case for chains of length Zi and the values for Zi are obtained from the discretization of the lognormal distribution. The subfigures (a)–(d) in Fig. 7 are organized similarly to the subfigures in Fig. 6. Although the two sets of predictions share similar qualitative trends, the simple superposition cannot quantitatively capture the changes in the crossover frequency and crossover moduli at low Wi. Thus, the full double reptation physics is needed to quantitatively describe the OSR moduli. We also note that the amount of discrepancy between the simple superposition results and the RDP results increases with increasing PDI. Nevertheless, since the plateau value at large Wi (Wi=100) is nearly identical for the simple superposition and the RDP model, the simple superposition results can be used to explain the weakened softening of the plateau value for the polydisperse case compared to the monodisperse case, as observed in Figs. 5(c) and 6(c). As shown in supplementary material Fig. S7 [39], the yy-component of the conformation tensor Aiyy(0) depends nonlinearly on the chain length Zi. Thus, for a discretized distribution of the lognormal distribution with a weight-average the same as the chain length in the monodisperse case, the weight-average of Aiyy(0) is not the same as Ayy(0) for the monodisperse case.

FIG. 7.

Comparison of results for RDP model with a simple linear superposition (LSP) of monodisperse results for PDI=2; (a) crossover frequency, (b) modulus at crossover, (c) plateau modulus, (d) ratio of moduli at crossover and plateau.

FIG. 7.

Comparison of results for RDP model with a simple linear superposition (LSP) of monodisperse results for PDI=2; (a) crossover frequency, (b) modulus at crossover, (c) plateau modulus, (d) ratio of moduli at crossover and plateau.

Close modal

In conclusion, the presence of chains of different length results in the smoother transition to the plateau. The presence of short chains and the nonlinear dependence of Aijyy(0) on the chain length results in the weakened softening of the plateau compared to the monodisperse case. Although a simple superposition of the monodisperse results can capture the qualitative trends observed in the polydisperse OSR moduli, to obtain quantitative predictions, the double reptation physics is needed.

Given the results just presented, we will now provide perspective on how this work could potentially help aid the design and interpretation of orthogonal superposition experiments.

By using a combination of numerical calculations and perturbation analysis, we showed for both the monodisperse case and the polydisperse case, the plateau of the orthogonal modulus is completely determined by the yy-component of the polymer conformation in the steady imposed shear, namely,

limω~1G(γ˙)GN0=Ayy(0)(γ˙),
(29)

where Ayy(0)(γ˙)=i=1nj=1nϕiϕjAijyy(0)(γ˙) is the weight-averaged result for the polydisperse case. As we discussed in the theory Sec. III B, Eq. (29) is generally valid for any model with the upper-convected derivative of the conformation tensor A and relaxation terms that are linear in the conformation tensor or any linear superposition of these models. Given this result and its generality, we can use orthogonal superposition measurements to obtain estimates of Ayy(0)(γ˙). Experimentally, this is accomplished by shifting G and G at different shear rates to achieve a complete overlap, in the same manner performed routinely for linear viscoelastic measurements for the purpose of time-temperature or time-concentration superposition, one can obtain shift factors α(γ˙)=G(ω,γ˙)G(ω,γ˙=0)=G(ω,γ˙)G(ω,γ˙=0). We can use the vertical shift factors to obtain Ayy(0)(γ˙):

Ayy(0)(γ˙)=α(γ˙)Ayy(0)(γ˙=0)=α(γ˙).
(30)

It is worth noting here that in time-temperature superposition, the conventional nomenclature for the vertical shift factor is b. Since we already used B to represent the dimensionless crossover frequency, which is related to the horizontal shift, we are using α to represent the vertical shift factor.

Since a rheological constitutive equation for an incompressible fluid only specifies the stress to within an arbitrary isotropic component, we cannot use Ayy(0)(γ˙) to directly determine the normal stress component in the yy-direction, σyy. However, we can attempt to make predictions of the first normal stress difference, N1. In particular, for the Rolie-Poly model, we can solve for the steady state value of Axx(0)(γ˙) using the O(1) equation in the perturbation analysis (see Sec. VIII in supplementary material for the detailed derivation [39]). Then, using results for Ayy(0)(γ˙) and B(γ˙), we obtain

N1(γ˙)=GN0(Axx(0)Ayy(0))=2GN0Ayy(0)[WiB]2.
(31)

If we recall that Wi=γ˙τR and B(γ˙) represents the dimensionless crossover frequency, B(γ˙)=ω~c,γ˙=ωc,γ˙τR, we get

N1(γ˙)=2GN0Ayy(0)[γ˙ωc,γ˙]2.
(32)

Experimentally, this means we can use the plateau value of G(γ˙), i.e., GN0Ayy(0)(γ˙), the applied shear rate, γ˙, and the shear-rate dependent crossover frequency, ωc,γ˙, to completely determine the first normal stress difference N1(γ˙). In the limit where a single-mode model accurately describes the experimental data, one can use this relation to determine N1(γ˙). We can further rearrange Eq. (32) to get

limω[G(ω,γ˙)ωc,γ˙2]=N1(γ˙)2γ˙2.
(33)

Previously, several relations have been proposed to connect the first normal stress difference or the yy-component of the stress to more easily measurable quantities in shear flow. Specifically, Holroyd et al. [45] derived analytical solutions of the Rolie-Poly model and found an expression to relate the recoverable strain to the yy-component of the stress, γrec=σxy0/σyy0. This theoretical result has been verified in molecular dynamics simulations of entangled polymers.[46] In LAOS measurements, Lee et al. found N1(t)=2G(t)γrec2(t), where G(t) is the recoverable strain-dependent elastic modulus and γrec is the recoverable strain [47]. Tanner and Williams [48] used the phenomenological K-BKZ model [49,50] to derive a relationship between the OSR moduli and the first normal stress difference limω0[G(ω,γ˙)ω2]=N1(γ˙)2γ˙2. This relation only works in the low-frequency limit and is specific to the K-BKZ model. The well-known Laun’s rule [32] relates the first normal stress difference to the linear viscoelastic moduli and is predicted to work for any general viscoelastic model,

N1(γ˙)=2G(ω)[1+(G(ω)G(ω))2]0.7|ω=γ˙.

Our result is an alternative set of relations that may prove more useful depending on where the sensitivity of the measurement lies. Previous rules provide relationships that are generally valid at asymptotically small frequencies and shear rates. We can now make a similar comparison at asymptotically large frequencies, which is usually the region with better instrumental sensitivity. Notably, the normal stress component is obtained as a by-product of measuring the plateau of the orthogonal modulus without having to build extra instrumentation to measure the normal stress.

From the perturbation analysis, we found the horizontal shift of the OSR moduli is determined by the effective relaxation rate B(γ˙). In the context of the Rolie-Poly model (i.e., for monodisperse polymers), since both Ayy(0)(γ˙) and B(γ˙) only depend on the ratio of relaxation time θ, the CCR parameter βCCR, and chain stretch due to the steady shear flow λ(0), OSR experiments can be used to estimate chain stretch and the CCR parameter βCCR. Specifically, one can first use linear viscoelastic measurements to estimate the Rouse time τR and the reptation time τd. Then, one can subtract the solvent contribution from the orthogonal superposition measurement (in the values for G). Next, G and G at different shear rates are shifted numerically to achieve a complete overlap with associated shift factors α(γ˙) and B(γ˙). Finally, we can use the crossover frequencies to obtain the relaxation rate B(γ˙). For the simplest case of βCCR=1 as we assumed in the previous calculations, we can obtain chain stretch λ(γ˙) from the crossover frequency

λ(γ˙)=22(ωc,γ˙τRωc,eqmτR).
(34)

In order to obtain both λ(γ˙) and βCCR, we need to use expressions for Ayy(0)(γ˙) and B(γ˙) that include βCCR,

Ayy(0)=λ(0)+2βCCR(11λ(0))θλ(0)+2βCCR(11λ(0))θ+2(λ(0)1)θ,
(35)
B=1θ+2(11λ(0))+2βCCR(11λ(0))λ(0).
(36)

By fitting the shear rate dependent experimental values of Ayy(0)(γ˙) and B(γ˙) to Eqs. (35) and (36), we can obtain λ(γ˙) and βCCR. Furthermore, using Eq. (30) and converting B(γ˙) to dimensional form [B(γ˙)=ωc,γ˙τR], we obtain expressions for λ(γ˙) and βCCR in terms of the vertical shifting factor α(γ˙) and crossover frequency measured in experiments,

λ(γ˙)=22ωc,γ˙τR(1α(γ˙))
(37)

and

βCCR=α(γ˙)ωc,γ˙τRτRτdωc,γ˙τR(1α(γ˙))ωc,γ˙2τR2(1α(γ˙))22.
(38)

Typically, the CCR parameter is assumed to be a constant and independent of shear rate. In the case that model predictions are used to generate the OSR moduli, as in the present analysis, Eq. (38) can be used to test the self-consistency of βCCR. More interestingly, in the case of experimental data where the apparent value(s) of the nonlinear constitutive model parameter(s) is unknown and in the limit where a single mode model accurately describes the experimental data, one can use this shifting procedure in order to extract and fit model parameters of the desired constitutive equation. In order to successfully estimate model parameters from OSR experiments, we, therefore, give the following recommendations: (1) The shear rates need to span a large enough range to cover both the case of no chain stretch and the case of significant chain stretch (γ˙τR>1). (2) The frequency range needs to span from values below the crossover frequency at equilibrium to frequency values corresponding to the plateau of G in flow.

We now discuss the implications of our analysis on the more general issue of the physical interpretation of OSR moduli measured in experiments in the absence of a particular constitutive model. In light of Fig. 6, the plateau modulus for both the monodisperse and polydisperse cases decreases with increasing Wi and the moduli shifts horizontally without any modification to the properties of the chain; e.g., those that might be extracted from a scaling analysis of the linear viscoelasticity. For example, in previous orthogonal superposition experiments on wormlike micelles by Vermant and co-workers [5] and Rothstein and co-workers [43], the authors assumed by analogy with linear viscoelasticity ξm(γ˙)=(kBT/G(γ˙)|ω)1/3 to estimate the apparent, shear-rate dependent mesh size, ξm, from the plateau modulus of orthogonal superposition moduli, G(γ˙)|ω. Vermant and co-workers linked the systematic decrease of the plateau modulus to shear-induced changes in the mesh size of the micellar network. Rothstein and co-workers went further to suggest the change in the mesh size is related to the breakage of micelles in flow. However, our analysis of OSR for the Rolie-Poly model clearly shows, even in the context of regular polymers, the plateau modulus decreases with increasing shear flow due to stretch-based relaxation mechanisms (CCR and retraction), in the absence of any chain breakage. A comparison between the experiments on WLMs by Rothstein and co-workers and our RDP model calculation with a lognormal distribution of PDI=2 shows a similar change in the plateau modulus (20%) over the comparable range of Wi, suggesting that changes in entanglement dynamics are important and must be considered in this range. Additionally, the effect of polydispersity adds another complication for any analysis that tries to extract chain length in flow from the plateau of OSR moduli. As such, one cannot simply use empirical relationships for estimating length that are developed for equilibrium linear viscoelastic measurements, such as those described in Sec. II A, to interpret OSR results. This knowledge from simulations is very useful for interpreting OSR experiments.

In addition to the applicability of linear viscoelastic scaling analysis, one may wonder about the potential for using OSR moduli to infer or understand polydispersity effects given the present analysis. We identified new parameterization of the OSR moduli to qualitatively relate to the degree of polydispersity. Specifically, the ratio of the moduli at the crossover to the plateau modulus as a function of shear rate becomes increasingly nonmonotonic and deviates from the constant value of 0.5 with increasing degree of polydispersity [Fig. 6(d)]. Although Fig. 6 focuses on the result for the lognormal distribution, the qualitative trend also holds for other continuous, broad distributions, except for a bi-modal distribution with sharp peaks. In particular, for PDI=1.01, the shape of the molecular weight distribution is nearly identical for lognormal and Gaussian; similarly, for PDI=2, the shape of the molecular weight distribution is very similar for the lognormal and exponential. Although in our analysis we primarily focused on the lognormal distribution, by changing the degree of polydispersity, we are also changing the shape of the molecular weight distribution.

Conventionally, viscoelastic measurements of polydisperse systems have been used to extract information about the distribution of relaxation times, i.e., the discrete linear response relaxation spectrum. In Figs. 7 and S9 [39], we show that the nonlinear decomposition of the orthogonal modulus into its relative contributions Aijyy(γ˙) is sensitive to the spectrum of the nonlinear relaxation times. In this way, these measurements could be used as a nonlinear analog to linear viscoelastic measurements in terms of obtaining the spectra of nonlinear relaxation rates. The idea of extracting a shear rate-dependent relaxation spectrum from OSR measurement was proposed by Yamamoto [21]. However, this approach has a drawback because one has to assume a functional form of the relaxation spectrum beforehand. For the Rolie-Poly and RDP models, we can analytically derive the results for the superposition moduli that naturally contain a sum of the individual chain pairs in the system [Eqs. (25) and (26)]. In this case, the relaxation rate of individual modes is proportional to Aijyy(γ˙). Thus, if the experimental system can be accurately described by these relationships, then one can fit experimental data to these two equations and just leave Aijyy(γ˙) and Bij(γ˙) as adjustable parameters for some number of modes. Then, what one measures is the discrete nonlinear relaxation spectrum and the resulting Aijyy(γ˙) and Bij(γ˙) in the nonlinear case are equivalent to G and τ for the multimode Maxwell model in the linear viscoelastic limit. The advantage of this approach is that we do not need to presuppose anything about what Aijyy(γ˙) and Bij(γ˙) are, as long as a general multimode model with the upper-convected time derivative of the conformation tensor and relaxation terms linear in the conformation tensor describes the material of interest. For the specific case of the RDP model, we actually know what Aijyy(γ˙) and Bij(γ˙) are, whereas if we presuppose a form for the relaxation spectrum, we would not necessarily know what Aijyy(γ˙) and Bij(γ˙) are. In other words, our results reaffirm the approach by Yamamoto that one can define some effective nonlinear relaxation spectrum, but beyond this, it gives proof to the idea that the relaxation spectrum should be related to the underlying constitutive behavior of the fluid. Going a step further, for the case of the Rolie-Poly and the RDP models, we can obtain analytical expressions for the discrete nonlinear relaxation spectrum from the expressions for the OSR moduli [Eqs. (25) and (26)].

To summarize the previous discussions, our study emphasizes the critical importance of comparing OSR experiments with this type of detailed OSR predictions using rheological models to understand the physical meaning of changes in orthogonal superposition moduli as is now routinely done for other types of nonlinear rheological measurement protocols.

In this work, we combine numerical calculations and a perturbation analysis using detailed microstructural models (Rolie-Poly and Rolie-Double-Poly models) to study orthogonal superposition for monodisperse and polydisperse entangled linear polymers. A perturbation analysis allows identification of key parameters (Ayy(0)(γ˙) and B(γ˙)) that determine the shifts of the OSR moduli. By quantifying the shifts in an OSR experiment and using the expressions for Ayy(0)(γ˙) and B(γ˙), we can estimate the average chain stretch λ(γ˙), the CCR parameter βCCR, and the first normal stress difference N1.

By systematically varying the degrees of polydispersity while keeping the weight-average molecular weight constant in a lognormal distribution of the molecular weight, we find qualitative signatures in the OSR moduli that correlate with the degree of polydispersity. Furthermore, the presence of short chains and the nonlinear dependence of Aijyy(0) on the chain length results in the weakened softening of the plateau compared to the monodisperse case. Although a simple superposition of the monodisperse results can capture the qualitative trends observed in the polydisperse OSR moduli, to obtain quantitative predictions, the double reptation physics is needed. By making an analogy with the discrete relaxation spectrum obtained from LVE moduli, we showed that OSR moduli can be used to extract a nonlinear relaxation spectrum, which reveals a potential significant advantage of OSR measurement over conventional rheological measurements.

We also found that the contributions of different relaxation mechanisms and polydispersity to the shifts of OSR moduli complicate the analysis of estimating changes in chain length in flow. Thus, we cannot simply use conventional scaling analyses developed for LVE measurements to interpret OSR results.

The perturbation analysis and simulation scheme developed in this work are quite general and can be applied for other rheological models and material systems. In fact, we found the result for yy-component of the chain conformation determining the plateau value of the orthogonal moduli to be a general result for any model with the upper-convected time derivative of the conformation tensor and relaxation terms that are linear in the conformation tensor. Our result can be used as an analogy to Laun’s rule for relating OSR moduli to the first normal stress difference.

In summary, we found that orthogonal superposition gives very useful information about nonlinear material moduli under flow, which can provide better sensitivity for testing constitutive models for nonlinear polymer processing. Results in our work have important implications for the design and interpretation of future orthogonal superposition experiments.

This research was primarily supported by U.S. Department of Energy, Office of Science, Office of Basic Energy Sciences under No. DE-SC0020988. The authors made use of computational facilities purchased with funds from the National Science Foundation (No. CNS-1725797) and administered by the Center for Scientific Computing (CSC). The CSC is supported by the California NanoSystems Institute and the Materials Research Science and Engineering Center (MRSEC; NSF No. DMR 1720256) at UC Santa Barbara. The authors wish to thank Dr. Joseph Peterson, Dr. Patrick Corona, and Michael Burroughs for helpful discussions.

1.
Hyun
,
K.
,
M.
Wilhelm
,
C. O.
Klein
,
K. S.
Cho
,
J. G.
Nam
,
K. H.
Ahn
,
S. J.
Lee
,
R. H.
Ewoldt
, and
G. H.
McKinley
, “
A review of nonlinear oscillatory shear tests: Analysis and application of large amplitude oscillatory shear (LAOS)
,”
Prog. Polym. Sci.
36
,
1697
1753
(
2011
).
2.
Tanner
,
R.
, and
J.
Simmons
, “
Combined simple and sinusoidal shearing in elastic liquids
,”
Chem. Eng. Sci.
22
,
1803
1815
(
1967
).
3.
Simmons
,
J. M.
, “
A servo-controlled rheometer for measurement of the dynamic modulus of viscoelastic liquids
,”
J. Sci. Instrum.
43
,
887
892
(
1966
).
4.
Colombo
,
G.
,
S.
Kim
,
T.
Schweizer
,
B.
Schroyen
,
C.
Clasen
,
J.
Mewis
, and
J.
Vermant
, “
Superposition rheology and anisotropy in rheological properties of sheared colloidal gels
,”
J. Rheol.
61
,
1035
1048
(
2017
).
5.
Kim
,
S.
,
J.
Mewis
,
C.
Clasen
, and
J.
Vermant
, “
Superposition rheometry of a wormlike micellar fluid
,”
Rheol. Acta
52
,
727
740
(
2013
).
6.
Mewis
,
J.
, and
G.
Schoukens
, “
Mechanical spectroscopy of colloidal dispersions
,”
Faraday Discuss. Chem. Soc.
65
,
58
64
(
1978
).
7.
Zeegers
,
J.
,
D.
van den Ende
,
C.
Blom
,
E. G.
Altena
,
G. J.
Beukema
, and
J.
Mellema
, “
A sensitive dynamic viscometer for measuring the complex shear modulus in a steady shear flow using the method of orthogonal superposition
,”
Rheol. Acta
34
,
606
621
(
1995
).
8.
Vermant
,
J.
,
P.
Moldenaers
,
J.
Mewis
,
M.
Ellis
, and
R.
Garritano
, “
Orthogonal superposition measurements using a rheometer equipped with a force rebalanced transducer
,”
Rev. Sci. Instrum.
68
,
4090
4096
(
1997
).
9.
Tao
,
R.
, and
A. M.
Forster
, “
End effect correction for orthogonal small strain oscillatory shear in a rotational shear rheometer
,”
Rheol. Acta
59
,
95
108
(
2020
).
10.
Simmons
,
J. M.
, “
Dynamic modulus of polyisobutylene solutions in superposed steady shear flow
,”
Rheol. Acta
7
,
184
188
(
1968
).
11.
Kwon
,
Y.
, and
A. I.
Leonov
, “
Remarks on orthogonal superposition of small amplitude oscillations on steady shear flow
,”
Rheol. Acta
32
,
108
112
(
1993
).
12.
Jacob
,
A. R.
,
A. S.
Poulos
,
A. N.
Semenov
,
J.
Vermant
, and
G.
Petekidis
, “
Flow dynamics of concentrated starlike micelles: A superposition rheometry investigation into relaxation mechanisms
,”
J. Rheol.
63
,
641
653
(
2019
).
13.
Lee
,
Y. J.
,
H.
Jin
,
S.
Kim
,
J. S.
Myung
, and
K. H.
Ahn
, “
Brownian dynamics simulation on orthogonal superposition rheology: Time–shear rate superposition of colloidal gel
,”
J. Rheol.
65
,
337
354
(
2021
).
14.
Mobuchon
,
C.
,
P. J.
Carreau
,
M.-C.
Heuzey
,
N. K.
Reddy
, and
J.
Vermant
, “
Anisotropy of nonaqueous layered silicate suspensions subjected to shear flow
,”
J. Rheol.
53
,
517
538
(
2009
).
15.
Jacob
,
A. R.
,
A. S.
Poulos
,
S.
Kim
,
J.
Vermant
, and
G.
Petekidis
, “
Convective cage release in model colloidal glasses
,”
Phys. Rev. Lett.
115
,
218301
(
2015
).
16.
Sung
,
S. H.
,
S.
Kim
,
J.
Hendricks
,
C.
Clasen
, and
K. H.
Ahn
, “
Orthogonal superposition rheometry of colloidal gels: Time-shear rate superposition
,”
Soft Matter
14
,
8651
8659
(
2018
).
17.
Moghimi
,
E.
,
J.
Vermant
, and
G.
Petekidis
, “
Orthogonal superposition rheometry of model colloidal glasses with short-ranged attractions
,”
J. Rheol.
63
,
533
546
(
2019
).
18.
Granek
,
R.
, and
M. E.
Cates
, “
Stress relaxation in living polymers: Results from a poisson renewal model
,”
J. Chem. Phys.
96
,
4758
4767
(
1992
).
19.
Doi
,
M.
, and
S. F.
Edwards
,
The Theory of Polymer Dynamics
(
Oxford University
,
New York
,
1988
), Vol.
73
.
20.
De Gennes
,
P.-G.
,
Scaling Concepts in Polymer Physics
(
Cornell University
,
Ithaca
,
1979
).
21.
Yamamoto
,
M.
, “
Rate-dependent relaxation spectra and their determination
,”
Trans. Soc. Rheol.
15
,
331
344
(
1971
).
22.
Wong
,
C.
, and
A.
Isayev
, “
Orthogonal superposition of small and large amplitude oscillations upon steady shear flow of polymer fluids
,”
Rheol. Acta
28
,
176
189
(
1989
).
23.
Unidad
,
H. J.
, and
G.
Ianniruberto
, “
The role of convective constraint release in parallel superposition flows of nearly monodisperse entangled polymer solutions
,”
Rheol. Acta
53
,
191
198
(
2014
).
24.
Marrucci
,
G.
, and
G.
Ianniruberto
, “
Flow-induced orientation and stretching of entangled polymers
,”
Philos. Trans. R. Soc. London, Ser. A
361
,
677
688
(
2003
).
25.
Mead
,
D. W.
, “
Small amplitude oscillatory shear flow superposed on parallel or perpendicular (orthogonal) steady shear of polydisperse linear polymers: The MLD model
,”
J. Non-Newtonian Fluid Mech.
195
,
99
113
(
2013
).
26.
Mead
,
D. W.
, “
Development of the ‘binary interaction’ theory for entangled polydisperse linear polymers
,”
Rheol. Acta
46
,
369
395
(
2007
).
27.
Mishler
,
S.
, and
D.
Mead
, “
Application of the MLD ‘toy’ model to extensional flows of broadly polydisperse linear polymers: Part II. Comparison with experimental data
,”
J. Non-Newtonian Fluid Mech.
197
,
80
90
(
2013
).
28.
Mead
,
D. W.
,
R. G.
Larson
, and
M.
Doi
, “
A molecular theory for fast flows of entangled polymers
,”
Macromolecules
31
,
7895
7914
(
1998
).
29.
Wapperom
,
P.
, and
R.
Keunings
, “
Impact of decoupling approximation between stretch and orientation in rheometrical and complex flow of entangled linear polymers
,”
J. Non-Newtonian Fluid Mech.
122
,
33
43
(
2004
).
30.
Boudara
,
V. A. H.
,
J. D.
Peterson
,
L. G.
Leal
, and
D. J.
Read
, “
Nonlinear rheology of polydisperse blends of entangled linear polymers: Rolie-double-poly models
,”
J. Rheol.
63
,
71
91
(
2019
).
31.
Likhtman
,
A. E.
, and
R. S.
Graham
, “
Simple constitutive equation for linear polymer melts derived from molecular theory: Rolie–poly equation
,”
J. Non-Newtonian Fluid Mech.
114
,
1
12
(
2003
).
32.
Laun
,
H.
, “
Prediction of elastic strains of polymer melts in shear and elongation
,”
J. Rheol.
30
,
459
501
(
1986
).
33.
Eckstein
,
A.
,
J.
Suhm
,
C.
Friedrich
,
R.-D.
Maier
,
J.
Sassmannshausen
,
M.
Bochmann
, and
R.
Mülhaupt
, “
Determination of plateau moduli and entanglement molecular weights of isotactic, syndiotactic, and atactic polypropylenes synthesized with metallocene catalysts
,”
Macromolecules
31
,
1335
1340
(
1998
).
34.
Macosko
,
C. W.
, Rheology: Principles, Measurements, and Applications, Advances in Interfacial Engineering Series (VCH, New York, 1994).
35.
Graham
,
R. S.
,
A. E.
Likhtman
,
T. C. B.
McLeish
, and
S. T.
Milner
, “
Microscopic theory of linear, entangled polymer chains under rapid deformation including chain stretch and convective constraint release
,”
J. Rheol.
47
,
1171
1200
(
2003
).
36.
Likhtman
,
A. E.
, and
T. C. B.
McLeish
, “
Quantitative theory for linear dynamics of linear entangled polymers
,”
Macromolecules
35
,
6332
6343
(
2002
).
37.
des Cloizeaux
,
J.
, “
Double reptation vs simple reptation in polymer melts
,”
Europhys. Lett.
5
,
437
442
(
1988
).
38.
Read
,
D. J.
,
K.
Jagannathan
,
S. K.
Sukumaran
, and
D.
Auhl
, “
A full-chain constitutive model for bidisperse blends of linear polymers
,”
J. Rheol.
56
,
823
873
(
2012
).
39.
See supplementary material at https://www.scitation.org/doi/suppl/10.1122/8.0000272 for (i) identification of linear viscoelastic region, (ii) validation of numerical method, (iii) comparison of perturbation analysis and numerical calculation, (iv) comparison of SRDP and RDP models for bidisperse polymers, (v) generalization of perturbation analysis, (vi) comparison of linear spring and nonlinear spring, (vii) dependence of Aijyy(0) on chain length Zi, and (viii) derivation of the relationship between first normal stress difference and OSR moduli.
40.
Nichetti
,
D.
, and
I.
Manas-Zloczower
, “
Viscosity model for polydisperse polymer melts
,”
J. Rheol.
42
,
951
969
(
1998
).
41.
Peterson
,
J. D.
, “Shear induced demixing in polymer melts and solutions,” Ph.D. thesis, University of California, Santa Barbara, 2018.
42.
Peterson
,
J. D.
,
G. H.
Fredrickson
, and
L.
Gary Leal
, “
Shear induced demixing in bidisperse and polydisperse polymer blends: Predictions from a multifluid model
,”
J. Rheol.
64
,
1391
1408
(
2020
).
43.
Khandavalli
,
S.
,
J.
Hendricks
,
C.
Clasen
, and
J. P.
Rothstein
, “
A comparison of linear and branched wormlike micelles using large amplitude oscillatory shear and orthogonal superposition rheology
,”
J. Rheol.
60
,
1331
1346
(
2016
).
44.
Metri
,
V.
, and
W.
Briels
, “
Brownian dynamics investigation of the Boltzmann superposition principle for orthogonal superposition rheology
,”
J. Chem. Phys.
150
,
014903
(
2019
).
45.
Holroyd
,
G. A.
,
S. J.
Martin
, and
R. S.
Graham
, “
Analytic solutions of the Rolie Poly model in time-dependent shear
,”
J. Rheol.
61
,
859
870
(
2017
).
46.
Anwar
,
M.
, and
R. S.
Graham
, “
Nonlinear shear of entangled polymers from nonequilibrium molecular dynamics
,”
J. Polym. Sci., Part B: Polym. Phys.
57
,
1692
1704
(
2019
).
47.
Lee
,
J. C.-W.
,
K. M.
Weigandt
,
E. G.
Kelley
, and
S. A.
Rogers
, “
Structure-property relationships via recovery rheology in viscoelastic materials
,”
Phys. Rev. Lett.
122
,
248003
(
2019
).
48.
Tanner
,
R.
, and
G.
Williams
, “
On the orthogonal superposition of simple shearing and small-strain oscillatory motions
,”
Rheol. Acta
10
,
528
538
(
1971
).
49.
Kaye
,
A.
, “Non-Newtonian flow in incompressible fluids,” College of Aeronautics Note 134 & 149s, 1962.
50.
Bernstein
,
B.
,
E.
Kearsley
, and
L.
Zapas
, “
A study of stress relaxation with finite strain
,”
Trans. Soc. Rheol.
7
,
391
410
(
1963
).
51.
Warner
,
H. R.
, “
Kinetic theory and rheology of dilute suspensions of finitely extendible dumbbells
,”
Ind. Eng. Chem. Fundam.
11
,
379
387
(
1972
).

Supplementary Material