Polymers in turbulence: stretching statistics and the role of extreme strain-rate fluctuations

We primarily model polymer molecules using the finitely extensible nonlinear elastic (FENE) dumbbell model (Bird et al., 1987; Öttinger, 1996; Graham, 2018). With $\tau$ as the relaxation time of the polymer, $R_{\mathrm{eq}}$ as the equilibrium root-mean-square (r.m.s.) end-to-end length, and $R_{m}$ as the maximum length, the dynamics of the end-to-end separation vector $\bm{R}$ in $d$ dimensions satisfies

where $\kappa_{ij}=\nabla_{j}u_{i}$ is the velocity gradient at the location of the centre of mass of the polymer, $f(R)=\big{(}1-R^{2}/R_{m}^{2}\big{)}^{-1}$ defines the FENE spring force, and $\bm{\xi}(t)$ is $d$-dimensional white noise that accounts for thermal fluctuations. This noise would also make an appearance in the equation of motion for the centre of mass. However, its effect on the transport of the dumbbell is very small compared to the advection by the turbulent carrier flow and thus may be neglected. The motion of the centre of mass of the polymer is therefore treated like that of a tracer.

Given a Lagrangian time series of $\bm{\kappa}(t)$, (1) is integrated using the Euler–Marujama method supplemented with the rejection algorithm proposed by Öttinger (1996), which rejects those time steps that yield extensions greater than $R_{m}\big{(}1-\sqrt{\mathrm{d}t/10}\big{)}^{1/2}$. Since the velocity gradient fluctuates, the numerical integration of (1) does not present the same difficulties as in the case of a laminar extension flow, and more sophisticated integration methods are not necessary. We have indeed checked that only a negligible fraction of time steps is rejected over a simulation.

Throughout this paper, we study the dumbbell model with $R_{\mathrm{eq}}=1$ and $R_{m}=110$. The elastic relaxation time $\tau$ defines the non-dimensional Weissenberg number, $\mathit{Wi}\equiv\lambda\tau$, where $\lambda$ is the maximum Lyapunov exponent of the carrier flow. The Weissenberg number provides a non-dimensional measure of the elasticity of the polymer, with high-$Wi$ polymers being easily extensible. We consider a wide range of $\mathit{Wi}$, from 0.3 to 40.

Note that while the dumbbell model (1) is used in most of this study, we do show in § 5 that our results also hold true for bead-spring chains.

To study polymer stretching in a turbulent flow, we use a database of Lagrangian trajectories from a DNS of homogeneous isotropic incompressible turbulence (at Taylor-microscale Reynolds number $\mathit{Re}_{\lambda}\approx 111$), generated at ICTS, Bangalore (James & Ray, 2017). The DNS solves the incompressible Navier–Stokes equations, discretized on a periodic cube, using a standard fully-dealiased pseudo-spectral method with $512^{3}$ grid points. Time integration is performed using a second-order slaved Adams–Bashforth scheme. The motion of $9\times 10^{5}$ tracers is calculated using a second-order Runge–Kutta method for time integration; the fluid velocity at the location of a tracer is obtained from its value on the grid using trilinear interpolation. The velocity gradient $\bnabla\bm{u}$ is calculated along these trajectories and stored at intervals of $0.11\tau_{\eta}$, where $\tau_{\eta}$ is the Kolmogorov time-scale (given below). This Lagrangian data provides the values of $\bm{\kappa}(t)$ for the integration of (1) along $9\times 10^{5}$ trajectories, allowing us to obtain good statistics of single-polymer stretching dynamics.

The Lyapunov exponent of the flow, required for defining $\mathit{Wi}$, is computed from these trajectories via the continuous $QR$ method, implemented using an Adams–Bashforth projected integrator along with the composite trapezoidal rule (for further details see Dieci et al., 1997). We find $\lambda=0.136/\tau_{\eta}$ which is compatible with previous simulations of isotropic turbulence (Bec et al., 2006). We also estimate the Lagrangian correlation time-scales of strain-rate and vorticity in the turbulent flow, $\tau_{S}$ and $\tau_{\varOmega}$, which serve as inputs to the Gaussian random model described in the next subsection. The rate-of-strain and rotation tensors are defined as $\mathsfbi{S}=(\bnabla\bm{u}+\bnabla\bm{u}^{\top})/2$ and $\mathsfbi{\varOmega}=(\bnabla\bm{u}-\bnabla\bm{u}^{\top})/2$, respectively. The autocorrelation functions of $S_{11}$ and $\varOmega_{12}$ are calculated and found to display an approximately exponential decay. Integrating these functions yields the integral time scales $\tau_{S}=2.20\,\tau_{\eta}$ and $\tau_{\varOmega}=8.89\,\tau_{\eta}$, in agreement with previous numerical simulations at comparable $R_{\lambda}$ (Yeung, 2001). The Kolmogorov time-scale $\tau_{\eta}$ is determined from $S_{11}$, using isotropy, as $\tau_{\eta}=\left({15\langle S_{11}^{2}\rangle}\right)^{-1/2}=3.72\times 10^{-2}$.

One of the goals of this study is to compare the stretching of polymers in a turbulent flow to that in a flow with Gaussian statistics, in order to determine the effect of extreme-valued fluctuations of the turbulent velocity gradient. For this, we use a Gaussian random velocity-gradient model to generate a time series of $\bm{\kappa}(t)$ for each polymer trajectory. Following Brunk, Koch & Lion (1997), we take $\bm{\kappa}(t)=\mathsfbi{S}(t)+\mathsfbi{\varOmega}(t)$ with

where $A$ determines the magnitude of the velocity gradient and $\zeta_{i}(t)$ ($i=1,\dots,5$) and $\varpi_{i}(t)$ ($i=1,2,3$) are independent zero-mean unit-variance Gaussian random variables with exponentially decaying autocorrelation functions and integral times $\tau_{S}$ and $\tau_{\varOmega}$, respectively. Therefore, $S_{ij}$ and $\varOmega_{ij}$ are Gaussian variables such that $\langle S_{ij}\rangle=\langle\varOmega_{ij}\rangle=0$,

and

As a consequence, $\langle\varOmega_{ij}\varOmega_{ij}\rangle=\langle S_{ij}S_{ij}\rangle$, which reproduces the relation $\nu\langle\omega^{2}\rangle=\langle\epsilon\rangle$, where $\omega$ is the vorticity and $\epsilon$ is the energy dissipation rate (Frisch, 1995). We take $\tau_{S}$ and $\tau_{\varOmega}$ to be the same as in the turbulent flow (§ 2.2). Furthermore, we set the coefficient $A=2.538$ so as to obtain approximately the same Lyapunov exponent $\lambda$ as in the turbulent flow. As a consequence, the Kubo number $\mathit{Ku}=\lambda\tau_{S}$ is also nearly the same in the turbulent and Gaussian flows.

The theory of Balkovsky et al. (2000) and Chertkov (2000) is recalled here in terms of the generalized Lyapunov exponents, rather than the Cramér function of the strain rate; the two properties are equivalent and related via a Legendre transform (see Boffetta, Celani & Musacchio, 2003). If $\bm{\ell}(t)$ is a line element in a random flow, the Lyapunov exponent is defined as

where $\langle\cdot\rangle$ denotes the average over the statistics of the velocity field. The $q$-th generalized Lyapunov exponent,

gives the asymptotic exponential growth rate of the $q$-th moment of $\ell(t)$. The function $\mathscr{L}(q)$ is positive and convex and satisfies $\mathscr{L}^{\prime}(0)=\lambda$ (see, e.g., Cecconi, Cencini & Vulpiani, 2010, for more details). If the flow is incompressible, then we also have $\mathscr{L}(-d)=\mathscr{L}(0)=0$, where $d$ is the space dimension (Zel’dovich et al., 1984).

For the elastic dumbbell (1), the p.d.f of the extension has a power-law form, $p(R)\sim R^{-1-\alpha}$ for $R_{\mathrm{eq}}\ll R\ll R_{m}$, with an exponent that satisfies

Given the properties of $\mathscr{L}(q)$, equation (7) implies that $\alpha$ is a decreasing function of $\mathit{Wi}$ that crosses zero at the critical value $\mathit{Wi}_{\rm cr}=1/2$ and saturates to $-d$ for very large $\mathit{Wi}$. In the limit $R_{m}\to\infty$, the p.d.f. of $R$ ceases to be normalizable for $\mathit{Wi}\geq\mathit{Wi}_{\rm cr}$, i.e. highly stretched configurations predominate, and so $\mathit{Wi}_{\rm cr}$ is taken to mark the coil–stretch transition.

In principle (7) may be used to determine $\alpha$ as a function of $Wi$. However, in general, calculating the function $\mathscr{L}(q)$ is very challenging. A useful approximation can be obtained in the vicinity of the coil–stretch transition by expanding about $q=0$: $\mathscr{L}(q)=\lambda q+\varDelta q^{2}/2+O(q^{3})$ with $\varDelta=\int\left(\left\langle\zeta(t)\zeta(t^{\prime})\right\rangle-\lambda^{2}\right)dt^{\prime}$ and $\zeta(t)=\bm{R}\cdot\bm{\kappa}(t)\cdot\bm{R}/R^{2}$. Substituting this quadratic expansion into (7) yields

Interestingly, in the limiting case of a time-decorrelated flow, (8) is accurate for all $\mathit{Wi}$, since $\mathscr{L}(q)$ is quadratic for all $q$ (Falkovich, Gawȩdki & Vergassola, 2001). Moreover, in this case $\lambda/\Delta=d/2$, which further simplifies (8) to yield:

To obtain $\alpha$ for polymers in a general chaotic flow and for arbitrary $\mathit{Wi}$, we must measure $\mathscr{L}(q)$; due to statistical errors this is especially difficult for values of $q$ that are negative or large and positive (Vanneste, 2010). Thus, past studies have been restricted to values of $\mathit{Wi}$ sufficiently near $\mathit{Wi}_{\rm cr}$ so that $\alpha$ does not deviate far from zero (Gerashchenko et al., 2005; Bagheri et al., 2012). In order to illustrate the validity of (7) over a wider range of $\mathit{Wi}$, we now consider a renewing (or renovating) random flow, for which $\mathscr{L}(q)$ can be calculated easily (Childress & Gilbert, 1995). To generate this flow, the time axis is divided into intervals $\mathscr{I}_{n}=[t_{n},t_{n+1})$ with $t_{n}=n\tau_{c}$ and $n=1,2,\dots$. The velocity field changes randomly at the beginning of each interval and then remains frozen for the rest of the time interval. Thus, the parameter $\tau_{c}$ sets the velocity correlation time. The velocity field is chosen to be a Couette flow, i.e. a two-dimensional linear shear flow, with a direction that is randomly rotated by an angle $\theta_{n}$ at the beginning of each time interval $\mathscr{I}_{n}$ (Young, 1999, 2009). For $t\in\mathscr{I}_{n}$ the velocity gradient takes the form

where $\sigma$ is the magnitude of the shear and $\theta_{n}$ is distributed uniformly over $[0,2\pi]$. The Lyapunov exponent as well as the generalized Lyapunov exponents for this flow can be calculated exactly (Young, 1999, 2009):

where $P_{q/2}$ is the Legendre function of order $q/2$.

The solution of (7) for the renewing Couette flow is presented in figure 1(a) for several values of $\sigma\tau_{c}$ (this non-dimensional group is the only free parameter that remains after substituting (11) in (7)). As $\tau_{c}$ is decreased towards zero the results are seen to approach the prediction for a delta-correlated flow (9), shown by the dashed (black) line. We also expect the results for all cases of $\sigma\tau_{c}$ to be well approximated by (9) in the vicinity of the coil–stretch transition, $\mathit{Wi}=\mathit{Wi}_{\rm cr}=1/2$. This is indeed the case. In addition, for this renewing flow,  (9) is seen to provide an excellent approximation for all $\mathit{Wi}$ greater than $\mathit{Wi}_{\rm cr}$. In general, one may expect the deviation of $\alpha$ from the prediction of (9) to be higher for small $\mathit{Wi}$ than for large $\mathit{Wi}$. This is because $\alpha$ is negative for large $\mathit{Wi}$ and the form of $\mathscr{L}(q)$ for negative arguments is strongly constrained by its convexity and the general properties $\mathscr{L}(0)=\mathscr{L}(-d)=0$ and $\mathscr{L}^{\prime}(0)=\lambda$.

In figure 1(b) the prediction of (7) is compared with the results of Brownian dynamics (BD) simulations of the dumbbell model (§ 2.1), with $\bm{\kappa}$ given by (10), for the case of $\tau_{c}=0.1$ and $\sigma=10$. The decorrelated-flow approximation is also shown for comparison (black dashed line). To obtain $\alpha$ from the simulations, we fit the stationary p.d.f of the extension $p(R)$ with a power law over a range $R_{\mathrm{eq}}\ll R\ll R_{m}$. Figure 1(b) shows excellent agreement between the large deviations theory and the simulations of the dumbbell model.

It is interesting to note that (7) may also be regarded as a tool for measuring the generalized Lyapunov exponents of a turbulent flow from the statistics of polymer extensions. Since $\alpha$ is a monotonic function of $\mathit{Wi}$, one could invert the $\alpha$ vs $\mathit{Wi}$ relation obtained from simulations and express $\mathit{Wi}$ in terms of $\alpha$ on the left-hand-side of (7), so as to obtain an explicit formula for $\mathscr{L}(\alpha)$. However, even with this strategy, measuring the generalized Lyapunov exponents for large negative and positive orders will remain challenging. On the one hand, because $\alpha\geq-d$ this approach cannot yield the generalized Lyapunov exponents of order less than $-d$. On the other hand, it is computationally intensive to construct the tail of the p.d.f. of $R$ for small $\mathit{Wi}$, which corresponds to large positive values of $\alpha$ and hence to generalized Lyapunov exponents of large positive order.

Let us now examine the p.d.f. of the extension for FENE dumbbells in homogeneous isotropic turbulence (§ 2.2). These results are presented in figure 2(a) (solid lines). A power-law range is apparent for $R$ between $R_{\mathrm{eq}}$ and $R_{m}$ (vertical dashed lines) and the corresponding exponents, $-1-\alpha$, are seen to increase past $-1$ as $\mathit{Wi}$ increases beyond $\mathit{Wi}_{\mathrm{cr}}$. (The $R^{2}$ behaviour for $R<R_{\mathrm{eq}}$ is a consequence of thermal fluctuations.) This panel also shows results for a Gaussian random flow (dotted lines), constructed so as to match the turbulent flow in terms of its integral correlation times of vorticity and strain as well as its Lyapunov exponent $\lambda$ (see § 2.3). Of course, the higher-order generalized Lyapunov exponents of these two flows will not be the same (see Biferale, Meneveau & Verzicco, 2014), and this is reflected in the differing values of $\alpha$ shown in figure 2(b) (markers). The thin solid line in this panel corresponds to the result (9) for a three-dimensional ($d=3$) time-decorrelated random flow. Though small, the differences between these results show a systematic dependence on $\mathit{Wi}$ which warrants further scrutiny.

Consider first the effect of the non-Gaussian statistics of turbulence. Figure 2(b) shows that low-$\mathit{Wi}$ stiff polymers stretch more in a turbulent flow than in a Gaussian flow: the values of $\alpha$ are less negative in the turbulent flow (compare the filled and open markers) which implies a greater power-law exponent for $p(R)$. Therefore, encountering a higher frequency of extreme-valued velocity gradients aids in stretching stiff polymers. Surprisingly, this is no longer true when $\mathit{Wi}$ is increased beyond $\mathit{Wi}_{\rm cr}$. Rather, these moderately high-$\mathit{Wi}$ polymers which are relatively easy to stretch are seen to be more extended in a Gaussian flow than in a turbulent flow. On further increasing $\mathit{Wi}$, all three flows in figure 2(b) eventually exhibit nearly identical values of $\alpha$; this is to be expected since $\alpha$ must attain the limiting value of $-3$ regardless of flow statistics (see § 3.1).

Why are moderately high-$\mathit{Wi}$ polymers stretched more in a Gaussian flow? To answer this, it is helpful to examine the p.d.f. of the rate of strain $s=\sqrt{S_{ij}S_{ij}}$ sampled by polymers. Figure 3(a) compares the distributions of $s$ for the turbulent and Gaussian flows. We see that the Gaussian flow has a comparative abundance of mild strain-rate events to compensate for its lack of extreme-valued events. This suggests that high-$\mathit{Wi}$ polymers that have long relaxation times are stretched more effectively by mild persistent straining rather than by strong but short-lived straining. In contrast, stretching small-$\mathit{Wi}$ polymers which have short relaxation times requires strong strain-rate events. These observations are consistent with a previous result by Terrapon et al. (2004) for a turbulent channel flow. It was shown that stretching events at low Wi are typically preceded by a burst of the strain rate; such bursts were not seen at high Wi.

If it is true that high-$\mathit{Wi}$ polymers are stretched primarily by mild persistent straining, then they should not only stretch more in a Gaussian flow but also remain in an extended configuration for a longer duration of time. To detect this behaviour, we carry out a persistence time analysis and quantitatively compare how long polymers stay stretched in the turbulent and Gaussian random flows. Interestingly, the concept of persistence, which arose out of problems in non-equilibrium statistical physics (Majumdar, 1999; Bray et al., 2013), has recently been used to study the turbulent transport of particles and filaments (Perlekar et al. 2011; Kadoch et al. 2011; Bhatnagar et al. 2016; Singh, Picardo & Ray 2022).

We begin by defining a polymer to be in a ‘stretched’ state if $R>\ell$, where the threshold $\ell=R_{m}/2$ is set well-within the power-law range. The non-stretched state is then defined as $R\leq\ell$. We have verified that varying the value of $\ell$ within the power-law range does not change our conclusions. With these states defined, we examine the Lagrangian history of each polymer and detect the time intervals $t_{\mathrm{str}}$ over which the polymer remains in a stretched state (see figure 3b). The distribution of this persistence time, $p(t_{\mathrm{str}})$, is presented in figure 3(c) for various $\mathit{Wi}$ and for both the turbulent and Gaussian flows. At large $t_{\mathrm{str}}$ the distribution displays an exponential tail, $p(t_{\mathrm{str}})\sim\mathrm{exp}(-t_{\mathrm{str}}/T_{\mathrm{str}})$, from which we extract the persistence time scale $T_{\mathrm{str}}$.

Figure 3(d) presents $T_{\mathrm{str}}$ for both flows and for various values of $\mathit{Wi}$. Clearly, polymers with $\mathit{Wi}>\mathit{Wi}_{\rm cr}$ typically remain stretched for a significantly longer time in the Gaussian flow as compared to the turbulent flow. This is true even for very large $Wi$ for which the exponent of the power-law of $p(R)$ is nearly the same in both flows (see the results for $\mathit{Wi}\geq 2$ in figures 2b and 3d). So though the probability of large extensions is the same, the nature of stretching is different: Polymers experience many short-lived episodes of large extension in the turbulent flow, whereas such episodes are fewer but last longer in the Gaussian flow.

Thus far we have been studying the stationary p.d.f. of $R$, attained after the polymers have spent enough time in the flow for their statistics to equilibrate. We now examine how the p.d.f. of $R$ evolves with time. Past work has shown that the p.d.f. relaxes exponentially to its stationary form, with a $\mathit{Wi}$-dependent time scale that exhibits a pronounced maximum at the coil–stretch transition (Celani et al., 2006; Watanabe & Gotoh, 2010). But what is the shape of the p.d.f. as it evolves? And how, if at all, does this evolution depend on $\mathit{Wi}$?

To answer these questions, we use our Lagrangian simulations to construct the p.d.f. of $R$ as a function of time, $p(R,t)$. We consider an initial state in which the polymers are in equilibrium with a static fluid. And so we first evolve the polymers with $\bm{\kappa}=0$ (in (1)) until the p.d.f. of $R$ attains the equilibrium distribution (Bird et al., 1987). We then ‘turn on’ the turbulent flow.

The evolution of $p(R,t)$ is illustrated in figure 4. Interestingly, we find two qualitatively distinct regimes. At small or moderate $\mathit{Wi}$ (figure 4a), the p.d.f. is seen to quickly attain a power-law form with an exponent $\beta(t)$ that increases in time until it reaches its stationary value $\beta_{\infty}=-1-\alpha$. By fitting the distributions with a power law in the range $R_{\mathrm{eq}}\ll R\ll R_{m}$, we find that $\beta(t)$ relaxes exponentially, i.e. $\beta_{\infty}-\beta\sim\mathrm{exp}(-t/T_{\beta})$, as demonstrated in figure 4(b). The time-scale $T_{\beta}$ is analyzed later in § 4.3.

The evolution at high $\mathit{Wi}$ is quite different and is shown in figure 4(c)-(d). Here, the polymers stretch rapidly and quickly produce a local peak close to the maximum extension $R_{m}$. Thus, although a transient power law appears at very early times, it is quickly lost as the local peak at $R_{m}$ begins to dominate the distribution (figure 4c). The long-time equilibration of the p.d.f. occurs by the peak near $R_{m}$ first approaching its stationary value; then the stationary power law $R^{\beta_{\infty}}$ gradually emerges, starting near $R_{m}$ and then extending its range down-scale towards $R_{\mathrm{eq}}$ (figure 4d).

These two regimes of equilibration are termed the evolving power-law regime and the rapid-stretching regime. The former occurs for $\mathit{Wi}\lesssim 3/4$, while the latter occurs for larger $\mathit{Wi}$. Note that the cross-over point, $\mathit{Wi}=3/4$, is marked by a stationary p.d.f. with $\beta_{\infty}=0$ that has neither a decaying tail at large $R$ nor a local peak near $R_{m}$; the evolution near $\mathit{Wi}=3/4$ is a blend of the two regimes.

We now lend credence to the above characterization of the evolving power-law regime by using a stochastic model to show that, for $0\leq\mathit{Wi}\lesssim 3/4$, an evolving power-law is a natural consequence of scale separation between $R_{\mathrm{eq}}$ and $R_{m}$. The associated analysis also reveals the $\mathit{Wi}$ dependence of the relaxation time-scale $T_{\beta}$.

We consider the Batchelor–Kraichnan flow wherein the velocity gradient $\bm{\kappa}(t)$ is a statistically isotropic, time-decorrelated, $3\times 3$ Gaussian tensor (Falkovich et al., 2001). In terms of the scaled variables $r=R/R_{\mathrm{eq}}$ and $\tilde{t}=t/2\tau$, the p.d.f. of the extension $p(r,\tilde{t})$ is governed by a Fokker–Planck equation (Martins Afonso & Vincenzi, 2005; Plan et al., 2016):

where $f(r)=1/(1-r^{2}/r_{m}^{2})$, $r_{m}=R_{m}/R_{\mathrm{eq}}$, and

Note that while (12) is exact for the Batchelor–Kraichnan flow, an equation of the same form may be obtained for general turbulent flows by considering the Fokker–Planck equation associated with (1), assuming statistical isotropy, and modelling the stretching term à la Richardson, i.e. as a diffusive process. The associated extension-dependent eddy diffusivity should scale as $R^{2}$ since the extension remains below the viscous scale.

We assume a wide separation between the scales of the equilibrium and maximum extensions, i.e. $R_{m}\gg R_{\mathrm{eq}}$, or $r_{m}\gg 1$. Further, we focus on the long-time evolution of $p(r,\tilde{t})$ which is characterized by distinct behaviours for small, intermediate, and very large extensions: (i) In the range of extensions where thermal fluctuations dominate ($r<1$), the p.d.f. of $r$ assumes a quadratic profile (figure 2); (ii) For intermediate extensions $1\ll r\ll r_{m}$, a power-law solution $r^{\beta(\tilde{t})}$ forms with an exponent $\beta(\tilde{t})$ that converges to its stationary value $\beta_{\infty}$; (iii) For extreme extensions, $p(r,\tilde{t})$ decreases rapidly as $r$ approaches $r_{m}$. We therefore introduce the ansatz:

Here, is a threshold value of $r/r_{m}$, less than unity, beyond which the nonlinear part of the FENE spring coefficient $f(r)$ becomes non-negligible and terminates the power-law behaviour. The function $g(r,\tilde{t})$ is such that, for $r=\ratio r_{m}$, we have $g(\ratio r_{m},\tilde{t})=c(\tilde{t})(\ratio r_{m})^{\beta(\tilde{t})}$; while for $r>\ratio r_{m}$, $g(r,\tilde{t})$ decreases faster than $r^{\beta(\tilde{t})}$ as $r$ approaches $r_{m}$. The latter assumption only applies if $\mathit{Wi}<3/4$, for which $\beta_{\infty}<0$. For larger $\mathit{Wi}$, our Lagrangian simulations show that the p.d.f. of $r$ develops a pronounced time-dependent peak between $\gamma r_{m}$ and $r_{m}$ (see figures 4c and 4d). Thus, the subsequent analysis applies only for $\mathit{Wi}<3/4$, which in fact corresponds to the regime in which the Lagrangian simulations exhibit an evolving power law.

The coefficient $c(\tilde{t})$ is determined in terms of the exponent $\beta(\tilde{t})$ through the normalization condition

Clearly, $I_{1}=c/3$, while $I_{2}=c[(\gamma r_{m})^{1+\beta}-1]/[1+\beta]$ for $\beta\neq 1$ and $I_{2}=c\ln(\ratio r_{m})$ for $\beta=-1$. For $r_{m}\gg 1$, $I_{3}$ is subdominant with respect to $I_{1}$ and $I_{2}$ and will be ignored. Thus, solving (15) for $c$ yields

To calculate $\beta(\tilde{t})$, we observe that, for $1\ll r\ll r_{m}$, (12) can be simplified as follows:

At steady-state we know that $p(r)\propto r^{\beta_{\infty}}$, which when substituted into (17) yields $\beta_{\infty}=2-3/2\mathit{Wi}$, in accordance with the large deviations theory for a time-decorrelated flow (see (9) and recall that $\beta_{\infty}=-1-\alpha$).

Now, substituting $p$ from (14b) into (17) results in

with

The left hand side of (18) contains $\mathrm{d}c/{\mathrm{d}\beta}$ which is evaluated using (16):

On taking the limit $r_{m}\gg 1$, the right-hand-side of (20a) simplifies to $-3c/(\beta^{2}-\beta-2)$ for $\beta<-1$ and to $-c\ln{(\ratio r_{m})}$ for $\beta>-1$. By substituting these expressions for $\mathrm{d}c/{\mathrm{d}\beta}$ in (18) and considering that $r_{m}\gg r\gg 1$, we obtain the following leading-order equation for $\beta$:

Because $\ln{(x)}$ is a weak function of $x$ for $x\gg 1$, we approximate $\ln r$ in (21a), as well as $\ln(\ratio r_{m})$ in (21b), as $\ln r_{m}$ to obtain

with $a=1$ for $\beta\leq-1$ and $a=-1$ for $\beta>-1$.

Being independent of $r$, (22) shows that $p(R,t)$ may be approximated (up to logarithmic corrections) by a time-dependent power-law in the range $R_{\mathrm{eq}}\ll R\ll R_{m}$.

Next, to reveal the long-time relaxation of $\beta$ to $\beta_{\infty}$, we substitute $\beta=\beta_{\infty}+\beta^{\prime}$ in (22) and linearize for small $\beta^{\prime}$. Noting that $F(\beta_{\infty},Wi)=0$ and that $\beta_{\infty}<-1$ for $\mathit{Wi}<Wi_{\mathrm{cr}}$, we obtain

In terms of dimensional time $t$, this implies an exponential relaxation of the power-law exponent,

with a time scale

that diverges as $\mathit{Wi}$ approaches $\mathit{Wi}_{\mathrm{cr}}$.

Note that (22) actually has two fixed points: $\beta_{\infty}$ and $-1$. The latter, though, can be shown to be dynamically unstable, leaving $\beta_{\infty}$ as the only stable fixed point.

The prediction (24) of an exponentially relaxing power-law exponent is certainly in agreement with the results of our Lagrangian simulations (see figure 4b). The $\mathit{Wi}$-dependence of the relaxation time-scale will be examined in light of (25) in the next section. But first, it is instructive to compare (24)-(25) against numerical simulations of the Fokker-Planck equation (12). Beginning from a distribution of coiled dumbbells, we solve (12) using second-order central differences and the LSODA algorithm (Petzold, 1983) with adaptive time-stepping (Wolfram Research, 2019). We take $R_{m}=10^{3}R_{\mathrm{eq}}$, thereby realizing a much larger scale separation than is practical in our Lagrangian simulations. The results are presented in figure 5.

The stochastic model exhibits the same two regimes of evolution as the Lagrangian simulations: the evolving power-law regimes for $\mathit{Wi}<3/4$ (cf. figures 4a and 5a) and the rapid-stretching regime for larger $\mathit{Wi}$ (cf. figures 4c-4d and 5b). Within the evolving power-law regime, we extract the exponent $\beta(t)$ and find that it does in fact evolve exponentially, in accordance with (24). The time-scale $T_{\beta}$, obtained from a least-squares fit, is presented as a function of $\mathit{Wi}$ in figure 5(c). The predicted divergence of $T_{\beta}$ at $\mathit{Wi}_{\rm cr}$ (see (25)), manifests in the simulations (which have a finite scale separation) as a prominent maximum of $T_{\beta}$ at $\mathit{Wi}_{\rm cr}$. Such a slowing down of the stretching dynamics at the coil–stretch transition was demonstrated by Martins Afonso & Vincenzi (2005) and Celani et al. (2006), but in regard to the equilibration of the entire p.d.f. of $R$. The behaviour of $T_{\beta}$ vs $\mathit{Wi}$ gives an alternative characterization of the same phenomenon.

We now quantitatively compare the time scales of equilibration of $p(R,t)$ in the turbulent and Gaussian flows (Lagrangian simulations) as well as in the Batchelor–Kraichnan decorrelated flow (numerical solution of the stochastic model).

One time scale of interest is that associated with the evolving power-law tail, $T_{\beta}$. However, this time scale is only relevant for $\mathit{Wi}<3/4$. So, to characterize the equilibration time for all $\mathit{Wi}$, we use the time scale associated with the entire p.d.f. of $R$. Martins Afonso & Vincenzi (2005) and Celani et al. (2006) showed, for random flows, that $p(R,t)$ approaches its asymptotic stationary form exponentially, with a time-scale $T_{P}$ that displays a maximum near $\mathit{Wi}_{\rm cr}$. Watanabe & Gotoh (2010) confirmed this prediction in a DNS of isotropic turbulence. We obtain $T_{P}$ from our simulations by fitting the evolution of the first moment of $p(R,t)$ (approximately the same value is obtained by fitting the higher moments).

Figure 6 compares the equilibration time scales $T_{\beta}$ and $T_{P}$ in panels (a) and (b), respectively, for the turbulent flow (filled markers) and the Gaussian flow (open markers). The results for the time-decorrelated Batchelor–Kraichnan flow are also presented (line). We see that the qualitative variation of these time scales with $\mathit{Wi}$ is the same in all flows and exhibits a maximum near $\mathit{Wi}_{\rm cr}$. In the case of $T_{P}$, a $\mathit{Wi}^{-1}$ asymptotic behaviour is visible (see the inset of figure 6b), in agreement with Martins Afonso & Vincenzi (2005) and Celani et al. (2006).

On comparing the results for the Gaussian and decorrelated flows, we see that the time correlation of the flow slows down the equilibration of the p.d.f. of $R$: both $T_{\beta}$ and $T_{P}$ are higher for the Gaussian flow (figures 6a and  6b). The equilibration is even slower for the turbulent flow (DNS) with its non-Gaussian velocity gradient statistics.

Qualitatively, chaotic random flows are seen to provide a good approximation for studying polymer stretching in homogeneous isotropic turbulence. Our results also show that the accuracy of the predictions, especially for finite-time statistics, are improved by incorporating temporal correlation in the random flow.