1. Introduction
The dynamics of stretchable particles suspended in a fluid, such as droplets, macromolecules or elastic fibres, depends strongly on the nature of the flow that transports them. Indeed, several classification schemes for laminar flows are based on the flow’s ability to distort a test particle (Tanner Reference Tanner1976; Olbricht, Rallison & Leal Reference Olbricht, Rallison and Leal1982; Cunha et al. Reference Cunha, de Souza Mendes, Thompson, Rodrigues and Quintella2023). For turbulent flows, the study of particle stretching relies on a statistical approach. In particular, the deformation of a tiny stretchable particle is closely related to that of a material line element and, thus, may be analysed using techniques from the theory of Lagrangian chaos (Crisanti et al. Reference Crisanti, Falcioni, Vulpiani and Paladin1991; Bohr et al. Reference Bohr, Jensen, Paladin and Vulpiani1998; Cecconi, Cencini & Vulpiani Reference Cecconi, Cencini and Vulpiani2010). This approach has been used, in the context of incompressible turbulent flows, to predict the statistics of stretching and breakup of polymers (Balkovsky et al. Reference Balkovsky, Fouxon and Lebedev2000, Reference Balkovsky, Fouxon and Lebedev2001; Chertkov Reference Chertkov2000; Vincenzi et al. Reference Vincenzi, Watanabe, Ray and Picardo2021) and droplets (Biferale, Meneveau & Verzicco Reference Biferale, Meneveau and Verzicco2014; Ray & Vincenzi Reference Ray and Vincenzi2018). Our goal is to extend these results to compressible random flows, while considering a generic stretchable particle.
One scenario of interest is that of deformable particles confined to the free surface of an incompressible fluid in turbulent motion. The chaotic two-dimensional flow at the surface has compressible small-scale dynamics – the contraction of area elements at the surface produces strong inhomogeneities in the spatial distribution of floating particles (Sommerer & Ott Reference Sommerer and Ott1993; Bec, Gawȩdzki & Horvai Reference Bec, Gawȩdzki and Horvai2004; Cressman et al. Reference Cressman, Davoudi, Goldburg and Schumacher2004; Sanness Salmon et al. Reference Sanness Salmon, Baker, Kozarek and Coletti2023). In three dimensions, our study would be relevant to tiny dense particles whose motion departs weakly from that of fluid tracers; such particles can be assigned an effective velocity field that is compressible (Boffetta et al. Reference Boffetta, Celani, De Lillo and Musacchio2007; Haller & Sapsis Reference Haller and Sapsis2008), with a negative (positive) divergence in regions where the underlying fluid flow is strain-dominated (vortical) (Bec, Laenen & Musacchio Reference Bec, Laenen and Musacchio2017; Ravichandran, Deepu & Govindarajan Reference Ravichandran, Deepu and Govindarajan2017). A situation that is not of interest here is high-Mach-number turbulence. In this case, indeed, the small-scale dynamics is effectively incompressible (Yu et al. Reference Yu, Fouxon, Wang, Li, Yuan, Mao and Mond2023): while the instantaneous velocity gradient has a non-zero trace, the three Lyapunov exponents sum to zero, implying that volume elements maintain their volume on the average.
To appreciate how compressibility at small scales affects the deformation of particles, it is helpful to first recall the key concepts of the general theory of line-element stretching in chaotic flows. The deformation of a material line element
$\boldsymbol{\ell }(t)$
that is transported in a flow with velocity field
$\boldsymbol{u}$
is determined by the velocity gradient
$\boldsymbol{\nabla }{\boldsymbol{u}}$
sampled along a Lagrangian trajectory
$\boldsymbol{x}(t)$
:
${d}_t \ell _i = \ell _{\!j}\,\partial _{\!j} u_i(\boldsymbol{x}(t),t)$
. Let us denote the instantaneous rate of exponential stretching or contraction as
$\gamma (t)$
, so that
$\ell (t)/\ell (0) = \exp (\gamma t)$
, where
$\ell (t)=\vert \boldsymbol{\ell }(t)\vert$
. The long-time average of
$\gamma (t)$
is the principal Lyapunov exponent
$\lambda = \lim _{t\to \infty } \langle \gamma (t)\rangle$
, where
$\langle \boldsymbol{\cdot }\rangle$
denotes the average over the statistical realisations of the velocity field. Information on the fluctuations of
$\gamma (t)$
, and hence on the statistical deviations from the average dynamics, is provided by the generalised Lyapunov exponents
$L(q)$
, which determine the long-time evolution of moments of
$\ell (t)$
:
$\langle \ell ^q(t)\rangle \sim \ell ^q(0)\mathrm{e}^{L(q) t}$
as
$t\to \infty$
(Zel’dovich et al. Reference Zel’Dovich, Ruzmaikin, Molchanov and Sokoloff1984; Benzi et al. Reference Benzi, Paladin, Parisi and Vulpiani1985; see also Crisanti et al. Reference Crisanti, Falcioni, Vulpiani and Paladin1991; Bohr et al. Reference Bohr, Jensen, Paladin and Vulpiani1998; Cecconi et al. Reference Cecconi, Cencini and Vulpiani2010). The function
$L(q)$
is convex and satisfies
$L(0)=0$
and
$L'(0)=\lambda$
. Moreover, if
$\lambda \gt 0$
, then
$L(-d_c)=0$
with
$d_c$
being the correlation dimension of the attractor of fluid–particle trajectories (Baxendale Reference Baxendale1991). The generalised Lyapunov exponents are related, via a Legendre transform, to the rate function (also known as the Cramér function), which determines the form of the large deviation approximation to the probability density function of
$\gamma (t)$
(Bohr et al. Reference Bohr, Jensen, Paladin and Vulpiani1998). In general,
$L(q)$
is difficult to calculate analytically and is known explicitly only for a few random flows, including incompressible linear renewing flows (Young Reference Young1999), telegraph noise (Falkovich & Martins Afonso Reference Falkovich and Martins Afonso2007) and isotropic Gaussian flows with vanishing correlation time (Chertkov, Kolokolov & Vergassola Reference Chertkov, Kolokolov and Vergassola1998). A quartic polynomial approximation of
$L(q)$
valid for general random flows is derived in Fouxon, Ainsaar & Kalda (Reference Fouxon, Ainsaar and Kalda2019).
In incompressible flows,
$\lambda$
is positive and all positive-order moments of
$\ell (t)$
grow exponentially in time. Further, the correlation dimension
$d_c$
equals the spatial dimension of the flow, i.e. material points or tracers sample the fluid uniformly. The function
$L(q)$
or, equivalently, the rate function has been computed in incompressible turbulent flows, specifically in numerical simulations of isotropic turbulence (Bec et al. Reference Bec, Biferale, Boffetta, Cencini, Musacchio and Toschi2006; Johnson & Meneveau Reference Johnson and Meneveau2015, Reference Johnson and Meneveau2016) and turbulent channel flows (Bagheri et al. Reference Bagheri, Mitra, Perlekar and Brandt2012; Johnson et al. Reference Johnson, Hamilton, Burns and Meneveau2017).
Now, let us consider compressible flows. A measure of the Eulerian degree of compressibility of a velocity field
$\boldsymbol{u}$
is provided by the ratio
where
$\Vert \boldsymbol{\nabla }\boldsymbol{u}\Vert =\sqrt {\sum _{i,j=1}^d (\partial _{\!j} u_i)^2}$
is the Frobenius norm of the velocity gradient with
$d$
being the spatial dimension of the flow. The degree of compressibility satisfies
$0\leqslant \wp \leqslant 1$
, with
$\wp =0$
for an incompressible flow and
$\wp =1$
for a potential flow
$\boldsymbol{u}=\boldsymbol{\nabla }\phi$
, where
$\phi$
is a differentiable scalar field. For free-surface turbulent flows – a case of particular relevance to our study – numerical simulations (Boffetta et al. Reference Boffetta, Davoudi, Eckhardt and Schumacher2004; Cressman et al. Reference Cressman, Davoudi, Goldburg and Schumacher2004) and experiments (Cressman et al. Reference Cressman, Davoudi, Goldburg and Schumacher2004; Li et al. Reference Li, Sanness Salmon, Hassaini, Chang, Mucignat and Coletti2025; Sanness Salmon et al. Reference Sanness Salmon, Chang, Mucignat and Coletti2025) report values of
$\wp$
between 0.34 and 0.5.
Lagrangian dynamics in compressible flows does not exhibit the universality of its incompressible counterpart. Specifically, the transport of scalar fields and tracers in fluctuating compressible flows has been found to depend strongly on the temporal correlation of the flow. For example, the inhomogeneity of the spatial distribution of tracer particles, in a compressible flow, can be enhanced or depleted by temporal correlations, depending on the spatial structure of the flow (Vergassola & Avellaneda Reference Vergassola and Avellaneda1997; Boffetta et al. Reference Boffetta, Davoudi, Eckhardt and Schumacher2004; Gustavsson & Mehlig Reference Gustavsson and Mehlig2013; Dhanagare, Musacchio & Vincenzi Reference Dhanagare, Musacchio and Vincenzi2014).
Stretching statistics in compressible flows and, in particular, the Cramér rate function have been studied in numerical simulations of free-surface turbulence (Boffetta, Davoudi & De Lillo Reference Boffetta, Davoudi and De Lillo2006) and two-dimensional chaotic flows (Pérez-Muñuzuri Reference Pérez–Muñuzuri2014). In general,
$\lambda$
decreases with increasing
$\wp$
, i.e. the mean rate of exponential stretching is reduced. For some random flows, there is even a critical degree of compressibility,
$\wp _c$
, beyond which
$\lambda$
becomes negative (Chertkov et al. Reference Chertkov, Kolokolov and Vergassola1998; Gustavsson & Mehlig Reference Gustavsson and Mehlig2013; Dhanagare et al. Reference Dhanagare, Musacchio and Vincenzi2014). It is thus customary to distinguish between a weakly compressible regime (
$0\lt \wp \lt \wp _c$
) and a strongly compressible one (
$\wp _c \lt \wp \leqslant 1$
) (Falkovich, Gawȩdki & Vergassola Reference Falkovich, Gawȩdki and Vergassola2001). In the latter case, line elements contract in typical realisations of the flow. However, because of the convexity of
$L(q)$
, high-order moments of
$\ell (t)$
continue to grow in time, which indicates strong deviations from the mean contraction behaviour, i.e. rare events of intense stretching.
As mentioned earlier, the stretching of small particles is closely related to that of material line elements. Therefore,
$\lambda$
determines the mean dynamics of a particle (stretching or contraction). However, the full statistics of the deformation, including the large deviations that can lead to breakup, are determined by
$L(q)$
. Hence, flows with the same
$\lambda$
but different
$L(q)$
yield different deformation statistics. This fact was demonstrated, in the incompressible case, by Picardo, Plan & Vincenzi (Reference Picardo, Plan and Vincenzi2023) who compared the stretching statistics of polymers in isotropic turbulence with those in a Gaussian flow having the same
$\lambda$
and the same correlation times of strain and vorticity. Small but systematic differences in the probability of large extensions were observed: stiff polymers stretched more in the turbulent flow (owing to extreme-straining events), while highly elastic polymers stretched more in the Gaussian flow (owing to persistent mild straining).
Here, we study how changes in
$L(q)$
, caused by flow compressibility, affect the statistics of particle deformation and breakup. Since the influence of compressibility on Lagrangian dynamics is flow-dependent, a universal relationship between particle stretching and flow compressibility is unlikely. We therefore study particular compressible flows and determine how the flow-specific changes in
$L(q)$
impact particle stretching; we select random model flows for which
$L(q)$
either is known analytically or may be computed accurately.
To gain insight into the effects of flow compressibility for a wide class of deformable bodies, we adopt the vector model of Olbricht et al. (Reference Olbricht, Rallison and Leal1982), which generalises the elastic dumbbell model to include the effects of internal viscosity and non-affine deformation. This model was initially used to identify strong-flow conditions for spatially uniform steady flows (Olbricht et al. Reference Olbricht, Rallison and Leal1982). In a series of subsequent papers, Szeri, Wiggins & Leal (Reference Szeri, Wiggins and Leal1991) and Szeri & Leal (Reference Szeri and Leal1992, Reference Szeri and Leal1993, Reference Szeri and Leal1994) used the same model to reveal the rich dynamics of orientation and stretching of deformable particles in inhomogeneous and unsteady flows. More recently, Cunha et al. (Reference Cunha, de Souza Mendes, Thompson, Rodrigues and Quintella2023) have used this model to test a new flow classification scheme. Obviously, a vector model is, in principle, only applicable to axisymmetric particles and is ill-suited to capture the high distortion that precedes a breakup event. Nevertheless, this minimal particle model, in conjunction with model random flows, provides an ideal setting for understanding the qualitative effects of flow compressibility on particle stretching and breakup.
The rest of this article is organised as follows. We first describe the vector model for a stretchable particle in § 2. Next, in § 3, we extend the theory of Balkovsky et al. (Reference Balkovsky, Fouxon and Lebedev2000, Reference Balkovsky, Fouxon and Lebedev2001) to general stretchable particles in compressible random flows and, thereby, derive a relationship between the stationary statistics of the particle size and
$L(q)$
. In § 4, we present analytical results for the compressible Batchelor–Kraichnan model, where the flow is Gaussian and has zero correlation time. In particular, we obtain the exact expression for the slope of the probability of intermediate particle sizes and calculate the mean time for breakups as a function of the degree of compressibility of the flow. An intriguing conclusion of the analysis is that stiff particles can break up more rapidly in a compressible flow even as the mean straining is suppressed. We then address the case of a random flow with non-zero correlation time in § 5, by studying a compressible renewing flow, wherein the velocity gradient transforms randomly at fixed time instants. Importantly, the analytical predictions, obtained for the delta-correlated flow, anticipate all the qualitative features of the numerical results, computed for the time-correlated flow. We end in § 6 with a summary of our key findings and a discussion of their implications.
2. Vector model of a stretchable particle
The vector model of Olbricht et al. (Reference Olbricht, Rallison and Leal1982) describes the size and orientation dynamics of a stretchable axisymmetric particle. The particle is assumed to be neutrally buoyant and smaller than the smallest length scale of the flow. In turbulent flows, this means that the size of the particle lies in the viscous dissipation range. Consequently, the centre of mass of the particle moves as a material point, and the particle stretching and orientation dynamics is fully determined by the velocity gradient. The configuration of the particle is specified by a single state vector
$\boldsymbol{R}$
, which provides its length
$R=\vert \boldsymbol{R}\vert$
and orientation vector
$\boldsymbol{n}=\boldsymbol{R}/R$
. The time evolution of
$\boldsymbol{R}$
is governed by
where
$\unicode{x1D734}$
and
$\unicode{x1D64E}$
are the antisymmetric and symmetric parts, respectively, of the velocity gradient at the position of the centre of mass of the particle. Equation (2.1) involves three parameters. The parameter
$\tau _{\!p}$
is the relaxation time of the elastic restoring force that describes the resistance of the particle to stretching. The shape factor
$g$
, which may take values in the range
$0\leqslant g\leqslant 1$
, accounts for the non-affine deformation of particles that are not perfect rods: values of
$g$
smaller than unity model the inefficiency of strain in rotating and stretching a particle of finite aspect ratio. Finally,
$\epsilon$
is the internal viscosity parameter, which models internal dissipation mechanisms associated with the deformation of the particle. It can assume non-negative values, ranging from zero to infinity;
$\epsilon =0$
represents an elastically deformable particle, while the limit
$\epsilon \rightarrow \infty$
models a rigid particle. Thus, (2.1) involves the interplay of rotational and extensional effects, induced by the fluid, as well as internal dissipation and relaxation of the particle. As described in Olbricht et al. (Reference Olbricht, Rallison and Leal1982), this equation encompasses several microscopic bodies for different values of the parameters, including a linear elastic dumbbell (
$\epsilon =0$
,
$g=1$
,
$\tau _{\!p}\gt 0$
), a linear elastic dumbbell with internal viscosity (
$\epsilon \gt 0$
,
$g=1$
,
$\tau _{\!p}\gt 0$
), a solid spheroidal particle (
$\epsilon =\tau _{\!p}=\infty$
,
$0\leqslant g\leqslant 1$
) and a material line element (
$\epsilon =0$
,
$g=1$
,
$\tau _{\!p}=\infty$
).
Equation (2.1) can be rewritten as coupled evolution equations for the size and orientation (Olbricht et al. Reference Olbricht, Rallison and Leal1982):
The latter equation is simply Jeffery’s equation for the orientation vector of a rigid spheroidal particle. Furthermore, the orientation dynamics depends only on the parameter
$g$
, while the size dynamics depends on the two parameters
$g/(\epsilon +1)$
and
$1/2\tau _{\!p} (\epsilon +1)$
. The limitations of this vector model become apparent on considering the case of
$g=0$
. While the equation for the orientation (2.3) reduces to that of a sphere, the equation for the extension (2.2) exhibits a non-physical loss of flow-induced stretching (this breakdown of the model is not unexpected since an isotropic sphere cannot be described by a vector). Here, we are primarily concerned with elongated soft particles; in this context, the vector model provides a simple way to generalise beyond the elastic dumbbell and account for internal viscosity and the effects of a finite aspect ratio.
Another consequence of the simplicity of (2.2) is that it only describes moderate deviations from equilibrium. On the one hand, if the action of the strain rate is weak compared with the elastic force, the size of the particle collapses to zero. On the other hand, if the straining is sufficiently intense, the linearity of the elastic force allows unbounded growth of the particle size. To resolve these issues, the model of Olbricht et al. (Reference Olbricht, Rallison and Leal1982) can be augmented in various ways. A non-zero equilibrium length
$R_{{eq}}$
can be imposed by including Brownian fluctuations in (2.1), while a maximum length
$R_{\textit{max}}$
can be imposed by replacing the linear force with a nonlinear one, such as that in the finitely extensible nonlinear elastic (FENE) dumbbell model (Bird et al. Reference Bird, Curtiss, Armstrong and Hassager1987). For
$\epsilon =0$
and
$g=1$
, the resulting evolution equation for
$\boldsymbol{R}(t)$
is (the general equation for
$\epsilon \gt 0$
and
$g\lt 1$
is provided in Appendix A)
\begin{align} \dot {R}_i=R_{\!j}\partial _{\!j} u_i-\frac {f(R)R_i}{2\tau _{\!p}} + \sqrt {\dfrac {R_{{eq}}^2}{\tau _{\!p} d}}\,\xi _i(t), \end{align}
where
$d$
is the spatial dimension of the flow and
$\boldsymbol{\xi }(t)$
is a
$d$
-dimensional white noise. In the limit
$R_{\textit{max}}\gg R_{{eq}}$
, the quantity
$R_{{eq}}$
is the root-mean-square size of the particle at equilibrium. Equation (2.4) corresponds to the FENE dumbbell model, well studied in polymer physics (Bird et al. Reference Bird, Curtiss, Armstrong and Hassager1987).
As mentioned in § 1, our aim is to understand the qualitative effects of compressibility on particle stretching, and these are expected to manifest more prominently in the intermediate range
$R_{{eq}}\ll R\ll R_{\textit{max}}$
, where the statistics of the particle size result from the balance between the opposing actions of the strain rate and the elastic force. Instead, sizes near to
$R_{{eq}}$
correspond to events where the straining is weak and molecular noise dominates the dynamics, while for
$R$
near to
$R_{\textit{max}}$
the nonlinear elastic force completely overcomes the stretching action of the flow to keep the length of the particle finite. It is expected, therefore, that the physically interesting range is
$R_{{eq}}\ll R\ll R_{\textit{max}}$
, and that the effects of compressibility on the statistics of the particle size will be largely independent of the specific small-
$R$
and large-
$R$
regularisations of (2.1).
We therefore focus on the range of intermediate particle sizes
$R_{{eq}}\ll R\ll R_{\textit{max}}$
and use (2.1) for the derivation of our main predictions. Of course, the resulting predictions should be interpreted only within the same intermediate range
$R_{{eq}}\ll R\ll R_{\textit{max}}$
. A regularised version of the model with non-zero
$R_{{eq}}$
and finite
$R_{\textit{max}}$
will be used, for practical reasons, in numerical simulations (see § 4) and in the exact solution of the Batchelor–Kraichnan model (Appendix B). These calculations with the regularised model are found to be in excellent agreement with the predictions of the theory based on (2.1), for the statistics of intermediate extensions.
3. Stationary size statistics in chaotic flows
3.1. The case of
$\epsilon =0$
and
$g=1$
We first consider the case
$\epsilon =0$
,
$g=1$
,
$\tau _{\!p}\gt 0$
, i.e. a Hookean dumbbell. We will then generalise the results to arbitrary values of the parameters in § 3.2.
Balkovsky et al. (Reference Balkovsky, Fouxon and Lebedev2000, Reference Balkovsky, Fouxon and Lebedev2001) studied the dynamics of a Hookean dumbbell in an incompressible flow, in the context of single-polymer dynamics in turbulence. In the absence of breakup events, they showed that the stationary probability density function (PDF) of the size,
$P_{\textit{st}}(R)$
, displays a power-law behaviour:
The exponent
$\alpha$
is determined by the rate function of the straining rate
$\gamma (t)$
. Given the large deviation form of the PDF of
$\gamma (t)$
as
where
$S$
is the rate function, one obtains
$\alpha$
from
where
$S'$
denotes the derivative of
$S$
and
$\beta$
satisfies
(Here we follow the definition of Balkovsky et al. (Reference Balkovsky, Fouxon and Lebedev2000, Reference Balkovsky, Fouxon and Lebedev2001), where
$S$
is such that
$S(0)=0$
. This definition is different from that used, for instance, in Crisanti et al. (Reference Crisanti, Falcioni, Vulpiani and Paladin1991), Bohr et al. (Reference Bohr, Jensen, Paladin and Vulpiani1998) and Cecconi et al. (Reference Cecconi, Cencini and Vulpiani2010), where the rate function vanishes at
$\lambda$
. The two definitions only differ by a translation of the argument of the rate function.)
Alternatively,
$\alpha$
can be expressed in terms of the generalised Lyapunov exponents (Boffetta, Celani & Musacchio Reference Boffetta, Celani and Musacchio2003) as the solution of
where the
$q$
th-order generalised Lyapunov exponent
$L(q)$
is defined as
In what follows, we use the latter formulation, wherein the connection between
$\alpha$
and the statistics of the straining rate is more readily apparent.
The preceding predictions hold for general random flows, the only assumption being that the correlation time of the velocity gradient,
$\tau _{\!f}$
, is finite (in the limit where the Kubo number of the flow,
$\mathit{Ku}=\lambda \tau _{\!f}$
, tends to infinity,
$\alpha$
approaches zero independently of the value of
$\tau _{\!p}$
; see Musacchio & Vincenzi (Reference Musacchio and Vincenzi2011)). In particular, though Balkovsky et al. (Reference Balkovsky, Fouxon and Lebedev2000, Reference Balkovsky, Fouxon and Lebedev2001) considered incompressible turbulence, the theory developed therein holds for compressible flows as well. Indeed, the derivation of (3.1) to (3.5) relies on estimating the probability of strong stretching events, and, as mentioned in § 1, such events are also present in compressible chaotic flows.
The solution of (3.5) depends on the form of
$L(q)$
and therefore on the details of the random flow that transports the particle. Nevertheless, some aspects of the dependence of
$\alpha$
on the system parameters can be deduced by using the general properties of
$L(q)$
, namely the convexity, the behaviour at the origin (
$L(0)=0$
and
$L'(0)=\lambda$
) and the fact that the second zero of
$L(q)$
is located at
$q=-d_c\geqslant -d$
when
$\lambda \gt 0$
. For a given relaxation time
$\tau _{\!p}$
, (3.5) has two solutions, which correspond to the two points of intersection between the straight line
$\alpha /2\tau _{\!p}$
and the graph of the function
$L(\alpha )$
(see figure 1). One of the two points is the origin; the other gives the value of
$\alpha$
that corresponds to
$\tau _{\!p}$
. Since the slope of the straight line
$\alpha /2\tau _{\!p}$
decreases with
$\tau _{\!p}$
, an immediate consequence is that
$\alpha$
is a decreasing function of
$\tau _{\!p}$
, i.e. the probability of large sizes increases as the particle becomes more stretchable. Other features of
$\alpha$
depend on the degree of compressibility of the flow.
We begin by recalling the incompressible case (
$\wp =0$
,
$\lambda \gt 0$
and
$d_c=d$
) (Balkovsky et al. Reference Balkovsky, Fouxon and Lebedev2000, Reference Balkovsky, Fouxon and Lebedev2001). The solution of this case is illustrated in figure 1(a). For small values of
$\tau _{\!p}$
(stiff particles), the straight line
$\alpha /2\tau _{\!p}$
is steep, and its intersection with the graph of
$L(\alpha )$
is located at a positive value of
$\alpha$
. Thus,
$P_{\textit{st}}(R)$
decays rapidly as
$R$
deviates from the equilibrium size. As
$\tau _{\!p}$
increases, the straight line
$\alpha / {2 \tau _{\!p} }$
intersects the graph of
$L(\alpha )$
at smaller and smaller values of
$\alpha$
. Consequently,
$P_{\textit{st}}(R)$
becomes less steep and the probability of large deformations increases, although, as long as
$\alpha \gt 0$
, the configurations near the equilibrium one remain the most likely. When
$\tau _{\!p}$
equals
$ ( 2 \lambda )^{-1}$
the slope of
$\alpha /(2\tau _{\!p})$
equals the slope of
$L(\alpha )$
at the origin. Therefore, the two solutions of (3.5) coincide, and
$\alpha$
vanishes. Now
$P_{\textit{st}}(R)$
varies as
$R^{-1}$
and is no longer normalisable in the limit
$R_{max } \rightarrow \infty$
; in the absence of a nonlinear cutoff, the size of the majority of particles undergoes unbounded growth. Following Balkovsky et al. (Reference Balkovsky, Fouxon and Lebedev2000, Reference Balkovsky, Fouxon and Lebedev2001), the value
$\tau _{\!p}=(2\lambda )^{-1}$
can therefore be regarded as the threshold for the transition from a shrunk state, close to equilibrium, to a highly stretched state. In the context of polymer dynamics, such a transition is known as ‘the coil–stretch transition’ (Gerashchenko, Chevallard & Steinberg Reference Gerashchenko, Chevallard and Steinberg2005); we here refer to it as the shrink–stretch transition.
Graphical construction of the solution of (3.5) in (a) the incompressible regime (
$\wp =0$
), (b) the weakly compressible regime (
$0\lt \wp \lt \wp _c$
) and (c) the strongly compressible regime (
$\wp _c\lt \wp \leqslant 1$
). The convex black curve is the graph of
$L(\alpha )$
, while the straight lines are
$\alpha /(2\tau _{\!p})$
for small (green line) and large (orange line)
$\tau _{\!p}$
.

In the vicinity of the transition,
$\alpha$
can be calculated exactly. For small
$q$
, the function
$L(q)$
is well approximated by a quadratic polynomial,
$L(q) \approx \lambda q + \varDelta q^2 / 2$
, where
$\varDelta$
is the variance of
$\gamma (t)$
and is related to the rate function as
$S''(0) = \varDelta ^{-1}$
(Crisanti et al. Reference Crisanti, Falcioni, Vulpiani and Paladin1991; Bohr et al. Reference Bohr, Jensen, Paladin and Vulpiani1998; Cecconi et al. Reference Cecconi, Cencini and Vulpiani2010). The quadratic behaviour of
$L(q)$
near
$q=0$
is consistent with the central limit theorem: the low-order moments of
$\ell$
are indeed dominated by the typical fluctuations of the stretching rate. Plugging this approximation into (3.5) yields (Balkovsky et al. Reference Balkovsky, Fouxon and Lebedev2000, Reference Balkovsky, Fouxon and Lebedev2001)
This expression holds for
$\alpha$
close to zero irrespective of the statistics of the straining rate.
As
$\tau _{\!p}$
exceeds
$( 2 \lambda )^{-1}$
, the exponent
$\alpha$
becomes negative, which means that
$P_{\textit{st}}(R)$
is less steep than
$R^{-1}$
and can even increase with
$R$
(of course, nonlinear effects will eventually cause
$P_{\textit{st}}(R)$
to decrease rapidly as
$R$
nears
$R_{\textit{max}}$
). The exponent saturates to
$-d$
in the limit of
$\tau _{\!p}\to \infty$
: the straight line
$\alpha / 2\tau _{\!p}$
approaches the horizontal axis, in this limit, so that it intersects the graph of
$L(\alpha )$
at
$\alpha \to -d$
. The scenario described thus far is well established in the context of polymer dynamics in incompressible turbulent flows (see e.g. Watanabe & Gotoh Reference Watanabe and Gotoh2010; Bagheri et al. Reference Bagheri, Mitra, Perlekar and Brandt2012; Picardo et al. Reference Picardo, Plan and Vincenzi2023).
Let us now see how the scenario changes in the weakly compressible case, wherein
$\wp \gt 0$
and
$\lambda$
decreases relative to its incompressible value while still remaining positive (hence,
$d_c$
will also be positive but smaller than its incompressible value). The solution of (3.5) in the weakly compressible case is depicted in figure 1(b). Applying the same arguments as in the incompressible case, we find that
$\alpha$
is positive for
$\tau _{\!p} \lt ( 2 \lambda )^{-1}$
, vanishes for
$\tau _{\!p} = ( 2 \lambda )^{-1}$
, takes negative values
$-d_c\lt \alpha \lt 0$
for
$\tau _{\!p} \gt ( 2 \lambda )^{-1}$
and tends to
$-d_c$
as
$\tau _{\!p}\to \infty$
. It is not possible to make general statements regarding the dependence of
$\alpha$
on
$\wp$
for small values of
$\tau _{\!p}$
, because this relationship is sensitive to the specific form of
$L(q)$
. However, since
$\lambda$
and
$d_c$
generally decrease with
$\wp$
, we expect the shrink–stretch transition to shift to larger values of
$\tau _{\!p}$
. And, for large
$\tau _{\!p}$
, we expect the probability of large deformations to be subdued by compressibility (indeed
$\alpha \geqslant -d_c\gt -d$
).
In a strongly compressible flow (
$\wp \gt 0$
and
$\lambda \lt 0$
), the statistics of the particle size differ significantly from the two previous cases. Since
$L'(0)=\lambda \lt 0$
, the function
$L(q)$
is positive for all
$q\lt 0$
. Consequently,
$\alpha$
cannot become negative. Its value is actually bounded below by the positive zero of
$L(q)$
, i.e.
$\alpha \geqslant \alpha _\star$
, where
$\alpha _\star$
is such that
$L(\alpha _\star )=0$
and
$\alpha _\star \gt 0$
(see figure 1
c). The shrink–stretch transition is therefore not observed in a strongly compressible flow and, at large
$\tau _{\!p}$
, the probability of strong deformations is depleted by compressibility. It is interesting to note, however, that the PDF of the particle size continues to behave as a power law, and therefore significant deviations from the equilibrium size remain possible, despite the fact that in typical realisations of the flow particles contract (
$\lambda \lt 0$
).
3.2. The general case
Up to now, we have examined the special case where the shape factor
$g$
is unity and the internal viscosity coefficient
$\epsilon$
is zero, which corresponds to a Hookean dumbbell. The behaviour of the exponent
$\alpha$
for
$0 \lt g \lt 1$
and
$\epsilon \gt 0$
can be obtained from the
$g=1$
,
$\epsilon =0$
case by using a scaling argument, as in Vincenzi (Reference Vincenzi2021). From (2.2), we deduce that the size dynamics of a particle with
$0 \lt g \lt 1$
and
$\epsilon \gt 0$
is equivalent to that of a Hookean dumbbell having relaxation time
$(1+\epsilon )\tau _{\!p}$
in a flow with rate-of-strain tensor
$ g\unicode{x1D64E} / (1+\epsilon )$
. To understand the effect of rescaling the rate-of-strain tensor
$\unicode{x1D64E}$
, we return to the time evolution of an infinitesimal line element along a trajectory
$\boldsymbol{x}(t)$
:
The length of the line element evolves according to the equation
whose formal solution is
Consider now the rescaled rate-of-strain tensor
$\unicode{x1D64E}_\xi =\xi\unicode{x1D64E}$
and denote the
$q$
th generalised Lyapunov exponent associated with
$\unicode{x1D64E}_\xi$
as
$\mathcal{L}(\xi , q)$
which satisfies
$\mathcal{L}(1, q) = L(q)$
for
$\xi =1$
by definition. From (3.10) and the definition of the generalised Lyapunov exponents in (3.6), it is clear that
$\mathcal{L}(\xi , q) = \mathcal{L}(1, \xi q) = L(\xi q)$
, i.e. a rescaling of the rate-of-strain tensor results in a rescaling of the argument of
$L(q)$
by the same factor. Therefore, (3.5) yields
\begin{align} { \begin{aligned} \alpha ( \epsilon , g, \tau _{\!p} ) &= 2 (1+\epsilon ) \tau _{\!p} \, \mathcal{L}\left ( \frac {g}{1+\epsilon }, \alpha ( \epsilon , g, \tau _{\!p} ) \right ) \\ &= 2 (1+\epsilon ) \tau _{\!p} \, \mathcal{L}\left ( 1, \frac {g}{1+\epsilon }\,\alpha ( \epsilon , g, \tau _{\!p} ) \right ) \end{aligned} } \end{align}
or, equivalently, in terms of
$L(q)$
,
Now we observe that (3.12) is identical to (3.5) with relaxation time
$g\tau _{\!p}$
for the variable
${g}\,\alpha ( \epsilon , g, \tau _{\!p} ) / ({1+\epsilon })$
. In other words, the quantity
${g}\,\alpha ( \epsilon , g, \tau _{\!p} ) / ({1+\epsilon })$
is the solution of (3.5) where the internal viscosity, shape factor and relaxation time are
$0$
,
$1$
and
$g\tau _{\!p}$
, respectively. Thus we find
The dependence of
$\alpha$
on the shape factor and the internal viscosity is encapsulated in this relation. In particular, internal viscosity increases the absolute value of
$\alpha$
and therefore steepens the PDF of the particle size at both small
$\tau _{\!p}$
, where
$\alpha$
is positive, and large
$\tau _{\!p}$
, where
$\alpha \lt -1$
, thus making the shrink–stretch transition sharper (Vincenzi Reference Vincenzi2021). Reducing
$g$
has the effect of increasing
$\alpha$
, both because the prefactor increases and because the effective relaxation time decreases (recall that
$\alpha$
is a decreasing function of
$\tau _{\!p}$
). At small
$\tau _{\!p}$
, this results in a depletion of the tail of
$P_{\textit{st}}(R)$
, consistent with the interpretation of
$g$
as a measure of the efficiency of strain in stretching particles. At large
$\tau _{\!p}$
, the effect is more subtle and depends on the details of how
$\alpha$
varies with
$\tau _{\!p}$
. Indeed, as long as
$\tau _{\!p}$
is such that
$\alpha (0,1,\tau _{\!p})\lt 0$
, decreasing
$g$
amplifies the prefactor
$1/g$
, but multiplying the relaxation time by
$g\lt 1$
reduces the absolute value of
$\alpha$
.
Because of (3.13), we can understand the influences of internal viscosity and the shape factor using just the
$\epsilon =0$
,
$g=1$
case. So, we focus on this case in our subsequent investigation of the effects of compressibility. The general properties of
$L(q)$
have already yielded some qualitative insights into the dependence of
$\alpha$
on
$\wp$
, especially for large
$\tau _{\!p}$
. For a more detailed understanding, we study specific examples of compressible random flows: the next section discusses the Batchelor–Kraichnan flow, which admits an analytical solution, while § 5 presents numerical simulations for a renewing flow. We shall see that compressibility can have counterintuitive effects, particularly for the breaking of stiff particles.
4. Delta-correlated Batchelor–Kraichnan flow
4.1. The compressible Batchelor–Kraichnan model
In the Batchelor regime of the Kraichnan model, the velocity field is Gaussian, has zero correlation time and is statistically homogeneous, isotropic, parity invariant and linear in space (Falkovich et al. Reference Falkovich, Gawȩdki and Vergassola2001). The second-order correlation of the velocity components is
where
$i,j=1,\ldots ,d$
, the separation vector is defined as
$\boldsymbol{r}=\boldsymbol{x}-\boldsymbol{y}$
,
$D_0$
and
$D$
are positive constants (in particular,
$D$
has the dimension of inverse time) and the degree of compressibility
$\wp$
has been defined in (1.1). The corresponding Eulerian correlation of the velocity gradient is (Falkovich et al. Reference Falkovich, Gawȩdki and Vergassola2001)
where
The Lagrangian correlation of the velocity gradient, relevant for line-element stretching, is the same as the Eulerian one because of the
$\delta$
-correlation in time of the Batchelor–Kraichnan model.
Now, to study the evolution of a line element, it is convenient to replace the Lagrangian velocity gradient
$\boldsymbol{\nabla }\boldsymbol{u}(\boldsymbol{x}(t),t)$
in (3.8) with a random tensor having the same statistics, i.e. to consider the stochastic differential equation
where
$\unicode{x1D642}(t)$
is a Gaussian
$d\times d$
tensorial noise such that
$\langle\unicode{x1D642}(t)\rangle =0$
and
It is important to note, however, that while the Itô and Stratonovich interpretations of (3.8) coincide, (4.4) and (3.8) are statistically equivalent only if the former is interpreted in the Itô sense (see Gawȩdzki Reference Gawȩdzki2008, §§ 2.4.2 and 2.9.3, for a detailed discussion). This proviso is particularly relevant in the compressible case, wherein the Itô and Stratonovich interpretations of (4.4) differ, with the Stratonovich counterpart of (4.4) being
where
$\widetilde {{\unicode{x1D63E}}}_{il}={\unicode{x1D63E}}_{\textit{ijjl}}=\wp D(d-1)(d+2)\delta _{il}$
(clearly, the additional drift term vanishes if the flow is incompressible). The difference between the two interpretations becomes important when trying to approach the delta-correlated limit by decreasing the correlation time of a time-correlated random flow – such a limiting procedure yields the Stratonovich interpretation provided the correlation function of the non-white noise is time-reversible (Bo & Celani Reference Bo and Celani2013). In § 5, we compare the stretching statistics in the Batchelor–Kraichnan flow with those in a time-correlated renewing flow; in order for the statistics in the two flows to agree, in the small correlation-time limit, we will have to account for the drift term of (4.6) while prescribing the Lagrangian velocity gradient of the renewing flow.
We come now to the generalised Lyapunov exponents of the compressible Batchelor–Kraichnan flow. As a consequence of the Gaussian statistics and the temporal decorrelation of the velocity field,
$L(q)$
is quadratic for all
$q$
and not just near
$q=0$
(Chertkov et al. Reference Chertkov, Kolokolov and Vergassola1998; Falkovich et al. Reference Falkovich, Gawȩdki and Vergassola2001):
with
(In the above expressions,
$d\gt 1$
; the case
$d=1$
is recovered from the same formula by taking
$D\propto 1/(d-1)$
; see Falkovich et al. (Reference Falkovich, Gawȩdki and Vergassola2001).) It follows that: (i) the critical degree of compressibility, above which
$\lambda \lt 0$
, is
$\wp _c=d/4$
; since
$\wp \leqslant 1$
, the strongly compressible regime only exists for
$d\leqslant 4$
; (ii)
$L(q)$
vanishes at
$q=0$
and
$q=-{2\lambda } /{\varDelta }=- ( d-4\wp ) / ( 1+2\wp )$
; therefore, for
$\wp \lt d/4$
, the correlation dimension is
$d_c= ( d-4\wp ) / ( 1+2\wp )$
; (iii)
$L(2)$
is independent of
$\wp$
for all
$d$
, i.e. the growth rate of the mean square of the length of a line element is independent of the degree of compressibility.
The function
$L(q)$
is plotted in figure 2(a) for
$d=2$
and various values of
$\wp$
(the qualitative behaviour of
$L(q)$
is similar for all
$d \lt 4$
). Two properties of
$L(q)$
are worth highlighting. First, even in the strongly compressible regime, where
$\lambda \lt 0$
, the moments
$\langle \ell ^q(t)\rangle$
grow in time as soon as
$q\gt ( 4 \wp - d ) / ( 1+2 \wp )$
, where the maximum value of the lower bound is
$(4-d)/3$
(for
$\wp =1$
). This growth in time of relatively low-order moments of
$\ell$
implies that line elements experience episodes of strong stretching, which deviate significantly from the mean shrinking behaviour. Second,
$L(q)$
increases with
$\wp$
for
$q\geqslant 2$
, i.e. strong-stretching events, caused by extreme fluctuations of the the strain rate (
$\hat {\boldsymbol{\ell }} \boldsymbol{\cdot }\unicode{x1D64E} \boldsymbol{\cdot } \hat {\boldsymbol{\ell }}$
in (3.9)), not only exist in the compressible regime, but become more prominent with increasing compressibility. Next, we examine the consequences of these features of
$L(q)$
on the statistics of particle stretching, in particular, on the exponent
$\alpha$
that determines the power-law form of the stationary PDF of the particle size.
Plots of (a)
$L(q)/D$
as a function of
$q$
, (b)
$g \alpha /(1+\epsilon )$
as a function of
$gD\tau _{\!p}$
and (c)
$g \alpha /(1+\epsilon )$
as a function of
$g \lambda \tau _{\!p}$
, for different values of
$\wp$
in the Batchelor–Kraichnan flow. In all panels,
$d=2$
. The values of
$\wp$
considered in (a,b) cover three regimes of compressibility: incompressible (
$\wp =0$
), weakly compressible (
$0\lt \wp \lt d/4$
) and strongly compressible
$(d/4 \lt \wp \leqslant 1)$
; only the incompressible and weakly compressible regimes are considered in (c).

4.2. Stationary PDF of the particle size
Since
$L(q)$
is quadratic in
$q$
, the exponent
$\alpha$
is obtained by inserting
$\lambda$
and
$\varDelta$
from (4.8) into (3.7) and (3.13):
(In Appendix B,
$\alpha$
is calculated in an alternative way by solving the Fokker–Planck equation for the PDF of
$R(t)$
.) The rescaled exponent is plotted as a function of
$gD\tau _{\!p}$
, in figure 2(b), for various compressibilities and for
$d=2$
(the behaviour is analogous for any
$d\lt 4$
). In the incompressible and weakly compressible cases (
$0\leqslant \wp \lt \wp _c$
),
$\alpha$
decreases from positive to negative values as
$gD\tau _{\!p}$
increases (vanishing at
$g\tau _{\!p}=(2\lambda )^{-1}$
, which marks the shrink–stretch transition). In contrast,
$\alpha$
is always positive in the strongly compressible case (
$\wp _c\lt \wp \leqslant 1$
). The curves
$\alpha$
versus
$gD\tau _{\!p}$
cross at
$gD\tau _{\!p}=[2(d-1)(d+2)]^{-1}$
. Therefore, two different situations can be identified. At large
$gD\tau _{\!p}\gt [2(d-1)(d+2)]^{-1}$
, the effect of compressibility is to increase
$\alpha$
; in particular, the asymptotic value of
$\alpha$
is
These features of the variation of
$\alpha$
with
$\wp$
and
$\tau _{\!p}$
are entirely consistent with the predictions, of § 3, that follow from the general properties of
$L(q)$
. The behaviour of
$\alpha$
at small
$\tau _{\!p}$
, however, depends on the flow-specific behaviour of
$L(q)$
at large
$q$
. Figure 2(b) shows that, for small
$gD\tau _{\!p} \lt [2(d-1)(d+2)]^{-1}$
, compressibility decreases the value of
$\alpha$
. This can be confirmed by comparing
$\alpha$
with its value for
$\wp =0$
(denoted as
$\alpha _0$
) in the limit of vanishing relaxation times:
As noted earlier, compressibility simultaneously decreases
$\lambda$
and increases the probability of large fluctuations of the line-element stretching rate (
$L(q)$
increases with
$\wp$
for
$q\gt 2$
). To isolate the consequences of the latter effect on
$\alpha$
, figure 2(c) presents the rescaled exponent as a function of
$g \lambda \tau _{\!p}$
(in the incompressible and weakly compressible regimes). This plot accounts for the displacement of the zero of
$\alpha$
due to the variation of
$\lambda$
with
$\wp$
(all curves now pass through zero at the same point); the large decrease of
$\alpha$
with compressibility, seen in figure 2(c) for small
$\tau _{\!p}$
, is entirely caused by strong fluctuations of the strain rate.
Plots of the stationary distribution
$P_{\textit{st}}(R)$
as a function of
$R$
for (a) small
$\tau _{\!p} \ ( \tau _{\!p} = 0.05)$
and (b) large
$\tau _{\!p} \ (\tau _{\!p}=2)$
for different values of the compressibility. For these distributions, we used the analytical expression (B5) for particles with
$\epsilon = 0$
and
$g = 1$
in the Batchelor–Kraichnan flow (we set
$R_{{eq}} =1.0$
and
$R_{\textit{max}}=100$
). For small
$\tau _{\!p}$
(stiff particles), all distributions decay as a power of
$R$
, corresponding to positive values of
$\alpha$
. Moreover, the particles clearly stretch more as the compressibility
$\wp$
increases. For large
$\tau _{\!p}$
(highly elastic particles), in contrast, particles stretch less as
$\wp$
increases.

Figure 3. Long description
Panel A: A line graph shows the stationary distribution P_st(R) as a function of R for small values of R. The x-axis is labeled R and spans from 10^-1 to 10^2. The y-axis is labeled P_st(R) and spans from 10^-12 to 10^0. Different lines represent various values of compressibility, denoted by different colors: blue for 0, orange for 0.25, green for 0.50, red for 0.75, and purple for 1.00. All distributions decay as a power of R, corresponding to positive values of R. The particles stretch more as the compressibility increases. Panel B: Another line graph shows the stationary distribution P_st(R) as a function of R for large values of R. The x-axis is labeled R and spans from 10^-1 to 10^2. The y-axis is labeled P_st(R) and spans from 10^-5 to 10^0. Different lines represent various values of compressibility, denoted by different colors: blue for 0, orange for 0.25, green for 0.50, red for 0.75, and purple for 1.00. For large values of R, particles stretch less as the compressibility increases.
The above results imply that flow compressibility displaces the shrink–stretch transition to higher values of
$D\tau _{\!p}$
and eventually suppresses the transition once the strongly compressible regime is attained. Moreover, compressibility depletes the probability of large extensions for highly elastic particles (large
$D\tau _{\!p}$
) and increases it for stiff particles (small
$D\tau _{\!p}$
). Figure 3 illustrates these opposing effects of compressibility on the stretching of stiff (figure 3
a) and highly elastic (figure 3
b) particles, for the case of
$\epsilon = 0$
and
$g = 1$
. Here, the PDF of
$R$
, calculated in Appendix B, is plotted for different values of
$\wp$
. Clearly, increasing compressibility leads to more stretching of stiff particles and less stretching of highly elastic particles. Consequently, the shrink–stretch transition becomes evermore gradual as
$\wp$
increases, until it disappears completely in the strongly compressible regime.
The contrasting effects of compressibility on elastic and stiff particles arise from the differing response of these particles to fluctuations of the strain rate (
$\boldsymbol{n}\boldsymbol{\cdot }\unicode{x1D64E}\boldsymbol{\cdot }\boldsymbol{n}$
in (2.2), see also (3.9)). Elastic particles integrate the straining action of the rate-of-strain tensor
$\unicode{x1D64E}$
over long times; therefore, the large but rare fluctuations of the strain rate are filtered out, and the particle dynamics is dominated by mild and frequent fluctuations, which become weaker as compressibility increases. Stiff particles have a very short response time to the fluctuations of the strain rate; their dynamics is therefore sensitive to large fluctuations even if these are short-lived, and the frequency of large fluctuations of the strain rate increases with
$\wp$
. Therefore, compressibility increases the probability of highly stretched stiff particles.
4.3. Mean breakup time
We have thus far analysed the stationary PDF of particle size while ignoring the possibility of breakup under the action of strong velocity gradients. We now assume that a particle breaks when the elastic force in (2.1) exceeds a threshold, which is equivalent to assuming that breakup occurs when the particle size exceeds a threshold extension
$R_{{br}}$
(with
$R_{{eq}}\ll R_{{br}}\ll R_{\textit{max}}$
). We aim to study the mean breakup time
$T_{{br}}$
and its dependence on flow compressibility.
As shown in Appendix B, the PDF of
$R$
at time
$t$
, denoted as
$P(R,t)$
, satisfies the following Fokker–Planck equation:
where the expressions for the drift and diffusion coefficients are given in Appendix B for
$\epsilon =0$
and
$g=1$
. Here, we obtain asymptotic results for the general case, using just the power-law behaviour of the stationary PDF for
$R_{{eq}}\ll R\ll R_{\textit{max}}$
and the fact that
$\mathfrak{D}_2\sim R^2$
in the same range of extensions (the diffusion coefficient must scale as
$R^2$
in view of the linearity of the velocity field in the neighbourhood of the particle).
In the absence of breakups, (4.12) would be complemented by reflecting boundary conditions at
$R=0$
and
$R=R_{\textit{max}}$
(i.e. the probability current
$\partial _{R} (\mathfrak{D}_2 P) - \mathfrak{D}_1 P$
is set to zero at both ends of the interval
$0\leqslant R\leqslant R_{\textit{max}}$
). The corresponding steady solution is the stationary PDF
$P_{\textit{st}}(R)\propto e^{-\varPhi _a(R)}$
, where
and the value of the positive constant
$a$
is fixed by the normalisation of
$P_{\textit{st}}(R)$
(Gardiner Reference Gardiner1985). Appendix B shows that
$P_{\textit{st}}(R)\sim R^{-1-\alpha }$
for
$R_{{eq}} \ll R \ll R_{\textit{max}}$
(with
$\alpha$
given by (4.9)), in agreement with the general theory of § 3.
Now, to account for particle breakups, we impose an absorbing boundary condition on the right-hand side, at the threshold value
$R_{{br}}$
, i.e. we set
$P(R_{{br}},t)=0$
for all
$t$
. Given a particle with initial size
$R(0)=\rho$
, such that
$R_{{eq}}\leqslant \rho \lt R_{{br}}$
, we obtain the average time it takes to break,
$T_{{br}}$
, by calculating the mean first-passage time through the boundary at
$R=R_{{br}}$
. We thus have (see Gardiner Reference Gardiner1985, § 5.2.7)
(We have rewritten the original expression in Gardiner (Reference Gardiner1985) so as to avoid the divergence of
$\varPhi _a(R)$
for
$a\to 0$
.)
We can estimate the asymptotic behaviour of
$T_{{br}}$
as a function of
$R_{{br}}/\rho$
, in the limit of
$R_{{br}}\gg R_{{eq}}$
, by ignoring the contribution of the range
$0\leqslant z\leqslant R_{{eq}}$
to the second integral in (4.14); this amounts to ignoring the small amount of time spent by the particle near
$R_{{eq}}$
compared with the time it takes to reach
$R_{{br}}$
. Since
${\rm e}^{-\varPhi _a(R)}\sim R^{-1-\alpha }$
and
$\mathfrak{D}_2(R)\sim R^2$
for
$R_{{eq}}\ll R\ll R_{\textit{max}}$
, we have
On inserting the latter expression into (4.14), we find
\begin{align} \textit{DT}_{{br}} \sim \operatorname {sgn}(\alpha )\int _\rho ^{R_{{br}}} \mathrm{d}y\, y^{\alpha -1}\big(R_{{eq}}^{-\alpha }-y^{-\alpha }\big) \propto \begin{cases} (R_{{br}}/\rho )^\alpha , & \alpha \gt 0, \\ \ln (R_{{br}}/\rho ), & \alpha \lt 0. \end{cases} \end{align}
This result implies that, in the incompressible and weakly compressible regimes (
$0\leqslant \wp \lt \wp _c$
), the mean breakup time increases as
$(R_{{br}}/\rho )^\alpha$
for
$\tau _{\!p}\lt (2\lambda )^{-1}$
and as the logarithm of
$R_{{br}}/\rho$
for
$\tau _{\!p}\gt (2\lambda )^{-1}$
. In the strongly compressible regime (
$\wp _c\lt \wp \leqslant 1)$
, the breakup time exhibits a power-law variation with
$R_{{br}}$
for all
$\tau _{\!p}$
.
The dependence of
$T_{{br}}$
on
$\epsilon$
and
$g$
cannot be deduced entirely from (4.16). However, we can use the scaling relation in (3.13) to predict that internal viscosity increases the slope of the power law when
$\alpha \gt 0$
, consistent with the fact that internal viscosity opposes sudden size changes that could result in particle breakups. The shape factor can have a more significant effect: unlike
$\epsilon$
, a variation of the shape factor can cause a change in the sign of
$\alpha$
, leading to a transition between a logarithmic and a power-law growth of
$T_{{br}}$
.
(a) Contour plot of the rescaled mean breakup time
$\textit{DT}_{{br}}$
as function of
$D\tau _{\!p}$
and
$\wp$
. The contours are spaced logarithmically and the horizontal axis is in logarithmic scale. (b,c) Plots of
$\textit{DT}_{{br}}$
versus
$D\tau _{\!p}$
and
$\wp$
for small and large values of
$D\tau _{\!p}$
, respectively; the horizontal axis is now in linear scale (also note the different range of values of
$\textit{DT}_{{br}}$
in the two panels). In all panels,
$d=2$
,
$\epsilon =0$
,
$g=1$
,
$R_{{br}}/\rho =10$
and the elastic force is linear.

The dependence of
$T_{{br}}$
on the degree of compressibility is more difficult to estimate analytically. We therefore insert the expressions of
$\mathfrak{D}_1(R)$
and
$\mathfrak{D}_2(R)$
for a Hookean dumbbell, given in (B4), into (4.14) and calculate
$T_{{br}}$
numerically (since
$R_{{br}}\ll R_{\textit{max}}$
, we have used a linear elastic force, i.e.
$f(R)=1$
in (B4)). The variation of
$T_{{br}}$
with
$D\tau _{\!p}$
and
$\wp$
is depicted in figure 4. Unsurprisingly, the breakup time decreases as the particle becomes more elastic, i.e.
$T_{{br}}$
is a decreasing function of
$D\tau _{\!p}$
(figure 4
a). With regard to compressibility, its contrasting effects on the stretching of stiff and elastic particles, uncovered in § 4.2, naturally produce contrasting effects on the breakup time:
$T_{{br}}$
decreases with
$\wp$
at small
$D\tau _{\!p}$
(see figure 4
b) and increases with
$\wp$
at large
$D\tau _{\!p}$
(figure 4
c), with the transition between the two behaviours occurring at
$D\tau _{\!p}=[2(d-1)(d+2)]^{-1}=1/8$
. So, stiff particles not only stretch more but also break up faster in flows that are more compressible.
5. Time-correlated renewing flow
We now check whether the effects of flow compressibility, identified in the previous section, remain the same in a flow with a finite time correlation. Specifically, we consider a two-dimensional renewing (or renovating) flow, where the velocity gradient along a Lagrangian trajectory is a sequence of independent identically distributed matrices, each of which remains constant for a time
$\tau _{\!f}$
(see e.g. Zel’dovich et al. Reference Zel’Dovich, Ruzmaikin, Molchanov and Sokoloff1984; Childress & Gilbert Reference Childress and Gilbert1995; Young Reference Young1999). The parameter
$\tau _{\!f}$
, thus, plays the role of the Lagrangian correlation time. We also account for the drift term of (4.6), so that the results in the renewing flow approach those in the Batchelor–Kraichnan flow as
$\tau _{\!f} \to 0$
. To be precise, the velocity gradient along a Lagrangian trajectory
$\boldsymbol{x}(t)$
in the renewing flow is
where
$D$
is a positive constant and
Here,
$\unicode{x1D642}_n=\unicode{x1D64E}_n+\unicode{x1D734}_n$
with
\begin{align}\unicode{x1D64E}_n= \sqrt {\frac {2D}{\tau _{\!f}}}\, \begin{pmatrix} \sigma _{1,n} +\sqrt {2\wp }\,\sigma _{2,n} & \sigma _{3,n} \\[5pt] \sigma _{3,n} & {-\sigma _{1,n}} +\sqrt {2\wp }\,\sigma _{2,n} \end{pmatrix}, \quad \unicode{x1D734}_n=2\sqrt {\frac {D}{\tau _{\!f}}(1-\wp )}\, \begin{pmatrix} 0 & \omega _n \\[3pt] -\omega _n & 0 \end{pmatrix}\!. \end{align}
In the above expressions,
$\sigma _{i,n}$
(
$i=1,2,3$
) and
$\omega _n$
are independent, zero-mean, Gaussian random variables such that
Hence,
$\unicode{x1D642}(t)$
is Gaussian and zero-mean. Its two-time correlation is
with
\begin{align} F(t,t')= \begin{cases} \tau _{\!f}^{-1} & \text{if $\lfloor t/\tau _{\!f}\rfloor = \lfloor t'/\tau _{\!f}\rfloor $}, \\[5pt] 0 & \text{otherwise}, \end{cases} \end{align}
where
$\lfloor \boldsymbol{\cdot }\rfloor$
denotes the floor function and
${\unicode{x1D63E}}_{\textit{ijkl}}$
has been defined in (4.3). Thus, in the limit
$D\tau _{\!f}\to 0$
, the matrix
$\unicode{x1D642}(t)$
has the same statistics as in the two-dimensional Batchelor–Kraichnan model (note, indeed, that
$\int _{-\infty }^{\infty } F(t,t')\,{\rm d}t'=1$
). The form of
$\unicode{x1D642}_n$
is obviously not unique, and other choices that reduce to the Batchelor–Kraichnan flow are also possible; however, the choice in (5.3) has a minimal number of noises and is therefore advantageous in numerical simulations. Finally, the additional term
$-4\wp D\unicode{x1D644}$
in (5.1) ensures that (3.8) becomes equivalent to the Itô equation (4.4) as
$D\tau _{\!f}\to 0$
. Hence, in the limit
$D\tau _{\!f}\to 0$
, the dynamics of a line element in the renewing flow tends to that in the Batchelor–Kraichnan flow with the same coefficient
$D$
and the same degree of compressibility
$\wp$
. This convergence is illustrated in figure 5 via a comparison of
$L(q)$
(the calculation of
$L(q)$
in the renewing flow is discussed below).
Plots of
$L(q)/\vert \lambda \vert$
versus
$q$
for the two-dimensional Batchelor–Kraichnan model and the renewing flow with either
$D\tau _{\!f}=0.01$
or
$D\tau _{\!f}=0.2$
. The degree of compressibility is (a)
$\wp =0$
, (b)
$\wp =0.25$
and (c)
$\wp =1$
. (d) The variation of
$\lambda /D$
with
$\wp$
in the two flows.

Our first goal is to determine
$\alpha$
, the power-law exponent of
$P_{\textit{st}}(R)$
, and to verify the dynamical-systems theory of § 3. To this end, we compare the values of
$\alpha$
obtained from simulations of the stochastic differential equation for
$\boldsymbol{R}(t)$
with those calculated from (3.5). This comparison requires, first, the computation of
$L(q)$
, which is eased by the temporal properties of the renewing flow. Since the velocity gradient in each interval
$\mathcal{I}_n$
remains fixed, the solution of (3.8) is
where
with
\begin{align}\unicode{x1D648}_n=\begin{pmatrix} c_n+\textstyle \sqrt {\frac {2D\tau _{\!f}}{\varSigma _n}}\,\sigma _{1,n} s_n & \sqrt {\frac {2D\tau _{\!f}}{\varSigma _n}}(\sigma _{3,n}+\sqrt {2}\omega _n\sqrt {1-\wp }) s_n \\[15pt] \sqrt {\frac {2D\tau _{\!f}}{\varSigma _n}}(\sigma _{3,n}-\sqrt {2}\omega _n\sqrt {1-\wp }) s_n & c_n-\textstyle \sqrt {\frac {2D\tau _{\!f}}{\varSigma _n}}\,\sigma _{1,n} s_n \end{pmatrix}\!. \end{align}
Here,
$\varSigma _n = 2D\tau _{\!f}[\sigma _{1,n}^2+\sigma _{3,n}^2-2(1-\wp )\omega ^2_n]$
, and
$c_n=\cosh (\sqrt {\varSigma _n})$
,
$s_n=\sinh (\sqrt {\varSigma _n})$
if
$\varSigma _n\geq 0$
or
$c_n=\cos (\sqrt {\vert \varSigma _n\vert })$
,
$s_n=\sin (\sqrt {\vert \varSigma _n\vert })$
if
$\varSigma _n\lt 0$
. The statistical isotropy of
$\unicode{x1D642}(t)$
leads to the following formula for the Lyapunov and generalised Lyapunov exponents:
where
$\boldsymbol{e}$
is any constant unit vector and the average is taken over the statistics of
$\unicode{x1D642}_n$
(see e.g. Childress & Gilbert (Reference Childress and Gilbert1995, Ch. 11) and Young (Reference Young1999)). The averages in (5.10) are calculated using the resampled Monte Carlo method of Vanneste (Reference Vanneste2010); this method accelerates statistical convergence which would otherwise be extremely slow in the compressible case. (For large values of
$\vert q\vert$
and values of
$\wp$
close to unity, statistical convergence required up to
$N=10^3$
time iterations and
$K=16\times 10^6$
realisations of the renewing flow.)
The possibility of calculating
$L(q)$
via (5.10) is a significant simplification that is specific to renewing flows. For more complex chaotic flows,
$L(q)$
must generally be calculated from direct numerical simulations of the Lagrangian dynamics; this can be an extremely challenging computation, especially for large
$\vert q\vert$
, because it requires accurate statistics of rare stretching events. Accurate knowledge of
$L(q)$
for large positive
$q$
is important while applying (3.5) to stiff particles; we have found it necessary to go up to
$q=4$
, at least, in order to solve (3.5) in the small-
$D\tau _{\!p}$
regime.
Plots of (a)
$L(q)/D$
as a function of
$q$
and (b)
$g \alpha /(1+\epsilon )$
versus
$gD\tau _{\!p}$
for different values of
$\wp$
in the renewing flow with
$\tau _{\!f}=0.2$
. In (b), the lines correspond to the numerical roots of (3.5), or equivalently (3.12), whereas the symbols are estimates of the slope of
$P_{\textit{st}}(R)$
for
$R_{{eq}}\ll R\ll R_{\textit{max}}$
obtained from numerical solutions of the stochastic differential equation for
$\boldsymbol{R}(t)$
. To obtain a non-zero equilibrium length and a finite maximum length, we modified (2.1) by adding Brownian noise with an amplitude such that
$R_{{eq}}=1$
and by replacing the linear elastic force with a nonlinear force that enforces
$R_{\textit{max}}=10^3$
. For simplicity, we also took
$\epsilon =0$
and
$g=1$
. The resulting stochastic differential equation for
$\boldsymbol{R}(t)$
is given in (2.4) and (2.5). We integrated (2.4) numerically by using the Euler–Maruyama method with time step
${\rm d}t=5\times 10^{-4}$
in combination with Öttinger’s rejection algorithm, to preserve the constraint
$R\lt R_{\textit{max}}$
(see Öttinger Reference Öttinger1996, § 4.3.2). These simulations confirm the validity of (3.5) in the compressible case.

Figure 6. Long description
The image contains two line graphs labeled as Panel A and Panel B. Panel A: This line graph shows the relationship between L divided by q and q for different values of phi. The x-axis is labeled q, and the y-axis is labeled L divided by D. The graph includes four lines representing different values of phi: 0, 0.25, 0.50, and 1.00. Each line shows a distinct curve, with the lines for higher phi values being more steeply curved. Panel B: This line graph illustrates the relationship between g alpha divided by 1 plus epsilon and g D tau p for different values of phi. The x-axis is labeled g D tau p and uses a logarithmic scale, while the y-axis is labeled g alpha divided by 1 plus epsilon. The graph includes four lines with symbols representing different values of phi: 0, 0.25, 0.50, and 1.00. The lines show a decreasing trend, with the lines for higher phi values being less steep.
In incompressible renewing flows,
$L(q)$
has been shown to vary continuously as a function of
$\tau _{\!f}$
(e.g. Childress & Gilbert Reference Childress and Gilbert1995; Musacchio & Vincenzi Reference Musacchio and Vincenzi2011). We have observed the same behaviour for the compressible renewing flow considered here. Therefore,
$\alpha$
and
$T_{{br}}$
are continuous functions of
$\tau _{\!f}$
. In subsequent calculations, the correlation time is set to
$D\tau _{\!f}=0.2$
(unless stated otherwise), which for
$\wp =0$
yields a Kubo number of
$\mathit{Ku} = \lambda \tau _{\!f}\approx 0.3$
. For comparison,
$\mathit{Ku}\approx 0.6$
in three-dimensional incompressible homogeneous isotropic turbulence (Watanabe & Gotoh Reference Watanabe and Gotoh2010).
Let us first examine the generalised Lyapunov exponents. Figure 5 shows that the stretching-rate statistics in the renewing flow, for
$D\tau _{\!f}=0.2$
, deviate significantly from those in the Batchelor–Kraichnan flow, when both flows are compressible. (The function
$L(q)$
remains quadratic, however, at least for the range of
$q$
computed here.) We find that the critical degree of compressibility of the renewing flow, for which
$\lambda =0$
, is
$\wp _c\approx 0.38$
. The variation of
$L(q)$
with
$\wp$
is depicted in figure 6(a); interestingly,
$L(q)$
increases with
$\wp$
for
$q\gtrsim 2$
, just as it does in the Batchelor–Kraichnan flow (see figure 2
a). So, even in the renewing flow, high-order moments of
$\ell$
grow faster in time as the flow becomes more compressible, i.e. large-stretching events become more prominent with increasing compressibility. Consequently, the dependence of the exponent
$\alpha$
on
$D\tau _{\!p}$
and
$\wp$
will be analogous to that in the Batchelor–Kraichnan model.
Figure 6(b) presents the predictions of
$\alpha$
from the theory of § 3 (the lines correspond to the numerical roots of (3.5), or equivalently (3.12)). As expected, compressibility has contrasting effects at small and large
$gD\tau _{\!p}$
: the exponent
$\alpha$
decreases with
$\wp$
at small
$gD\tau _{\!p}$
, while it increases with
$\wp$
at large
$gD\tau _{\!p}$
; the transition between the two behaviours occurs at
$gD\tau _{\!p}\approx 0.18$
. The consequences of such variations of
$\alpha$
for the stationary statistics of the particle size have already been discussed in § 4. Figure 6(b) also presents measurements of
$\alpha$
from simulations of the regularised stochastic differential equation for
$\boldsymbol{R}(t)$
in the renewing flow (symbols), considering particles with
$\epsilon = 0$
and
$g = 1$
(see the caption of figure 6 for details). We find excellent agreement between the stochastic simulations and the theory of § 3 (symbols and lines in figure 6
b, respectively). As anticipated, the small- and large-
$R$
regularisations do not affect the power-law behaviour of
$P_{\textit{st}}(R)$
over intermediate extensions.
Next, we measure the dependence of the mean breakup time on the critical size and the degree of compressibility, by simulating the stochastic differential equation for
$\boldsymbol{R}(t)$
(see the caption of figure 6 for details). As in the Batchelor–Kraichnan model, opposing variations of
$\textit{DT}_{{br}}$
with
$\wp$
appear on either side of
$gD\tau _{\!p}\approx 1.8$
. At small
$gD\tau _{\!p}$
(stiff particles), the mean breakup time reduces as the compressibility of the flow increases (figure 7
a). As explained in § 4, this surprising outcome can be understood by noting that (i) the dynamics of stiff particles is more sensitive to rare but intense straining events rather than to persistent but mild straining and (ii) large fluctuations of the strain rate become more prominent as
$\wp$
increases. At large
$gD\tau _{\!p}$
(highly elastic particles),
$\textit{DT}_{{br}}$
is an increasing function of
$\wp$
, since breakup events are essentially determined by the mean dynamics (see figure 7
b,c).
Regarding the dependence of
$\textit{DT}_{{br}}$
on the critical size for breakup, we find a transition from a logarithmic to a power-law dependence as
$\alpha$
changes from negative to positive values (figure 7
b,c). This result agrees with the prediction of the Batchelor–Kraichnan model in (4.16). In particular, the exponents of the power laws in figure 7(c) are close to the values of
$\alpha$
that determine the slope of
$P_{\textit{st}}(R)$
in the absence of breakups (figure 6
b).
(a) Rescaled mean breakup time as a function of the degree of compressibility for
$R_{{br}}/\rho = 10$
and two different values of
$D\tau _{\!p}$
. (b,c) Rescaled mean breakup time versus
$R_{{br}}/\rho$
for
$D\tau _{\!p}=1$
and different values of
$\wp$
. In (c), the constant
$c$
is determined by the power-law fit of the data (
$c =-3.4$
for
$\wp =0.5$
and
$c =-7.2$
for
$\wp =1.0$
). In all panels,
$\epsilon =0$
and
$g=1$
.

6. Concluding remarks
We have examined the extensional dynamics of tiny deformable particles in compressible fluctuating flows, using a generic particle model and simple, synthetic, random flows. In these flows, the mean stretching of line elements is suppressed as compressibility increases, as evidenced by the decreasing value of the principal Lyapunov exponent
$\lambda$
. This mean behaviour, however, does not determine entirely the effect of compressibility on particle stretching and breaking. Indeed, a key insight of this work is that the fluctuations of the line-element stretching rate play a critical role, giving rise to several important effects: (i) a broad, power-law distribution of extension, which persists even when the flow is strongly compressible (with
$\lambda \lt 0$
); this implies that, regardless of the degree of compressibility, large deviations of the stretching rate can produce highly extended particles; (ii) a counterintuitive increase with compressibility of the stretching and breakup rate of stiff particles, provided the large positive fluctuations of the stretching rate – quantified by high positive-order generalised Lyapunov exponents – increase with compressibility; and (iii) a reduction in the sharpness of the shrink–stretch transition; this outcome follows from the aforementioned increase in the stretching of stiff particles, as well as a contrasting decrease in the stretching of highly elastic particles.
The opposing effects of compressibility on stiff and highly elastic particles are due to (i) the increase in the relative frequency of intense-straining events as the compressibility of our model flows is increased and (ii) the differing impact of such events on stiff and highly elastic particles. The long elastic relaxation time of highly elastic particles averages out the strain rate and renders the extension insensitive to short bursts of straining. In contrast, stiff particles with a short relaxation time are more responsive to sudden strong straining; in fact, such events are necessary to overcome the strong elastic restoring force and stretch out the particle; so an increase in the frequency of strong-straining events, even at the cost of suppressing the mean straining, will produce more extended particles.
These effects of compressibility are reminiscent of the effects of intermittency, revealed by Picardo et al. (Reference Picardo, Plan and Vincenzi2023) in the context of elastic dumbbells and chains in incompressible turbulent flows. Increasing either property increases the proportion of strong straining relative to mild straining, which as discussed above produces contrasting changes in the extension of stiff and highly elastic particles; these opposing effects result in a shallower shrink–stretch transition. Now, if the model flows considered in this work were rendered intermittent, then we would expect the effects of intermittency to augment those of compressibility, so that the shrink–stretch transition will be even more shallow than what is seen in our Gaussian compressible flows. The simultaneous action of intermittency and compressibility will be relevant in compressible turbulent flows, especially at high Reynolds numbers, and both these effects should be borne in mind when considering particle deformation in such flows.
A summary of our main technical results is provided below for the convenience of the reader:
-
(i) In § 3, we have extended the large deviations theory of Balkovsky, Fouxon & Lebedev (Reference Balkovsky, Fouxon and Lebedev2000) to compressible random flows and to stretchable axisymmetric particles with internal viscosity and undergoing non-affine deformation. For intermediate sizes, the stationary PDF of the particle size,
$P_{\textit{st}}(R)$
, behaves as a power-law
$R^{-1-\alpha }$
. Some of the effects of compressibility on the exponent
$\alpha$
are deduced from the general properties of the generalised Lyapunov exponents (see figure 1). The dependence of
$\alpha$
on the shape factor and on internal viscosity is expressed by a simple scaling relation given in (3.13). -
(ii) For the Bachelor–Kraichnan flow, the exact expression of
$\alpha$
is given in (4.9); its dependence on the degree of compressibility
$\wp$
and the rescaled relaxation time of the particle,
$D\tau _{\!p}$
, is presented in figure 2. The PDF
$P_{\textit{st}}(R)$
is calculated explicitly in Appendix B, see (B5) and (B6), and is plotted in figure 3. -
(iii) For the Batchelor–Kraichnan flow, we have also estimated the mean breakup time of particles,
$T_{{br}}$
, as a function of the threshold size
$R_{{br}}$
, see (4.16). A transition between a power-law and a logarithmic dependence is observed as
$\wp$
and
$D\tau _{\!p}$
are varied. Unsurprisingly, stiff particles break slower than highly elastic particles, in a flow with a fixed degree of compressibility. Interestingly, on varying the compressibility of the flow, we find that the breakup time of stiff particles decreases as a function of compressibility, while that of highly elastic particles increases (figure 4). -
(iv) Our theoretical predictions are confirmed by numerical calculations for a time-correlated renewing flow. The results for
$\alpha$
and
$T_{{br}}$
in the renewing flow are shown in figures 6(b) and 7, respectively.
The two model flows considered here, namely the time-decorrelated Batchelor–Kraichnan flow and the time-correlated renewing flow, both exhibit an increase with flow compressibility of the magnitude of the high positive-order generalised Lyapunov exponents. There is no reason why this should be true in general, and it is quite possible that other flows show the opposite trend; if so, then the frequency of strong-straining events will decrease with compressibility and both stiff and highly elastic particles will stretch less in flows of increasing compressibility. This qualitative dependence of particle stretching on the full statistics of straining is the crucial insight of the present work.
Predicting the effect of compressibility on deformable particles in a given flow will require knowledge of the variation of the generalised Lyapunov exponents with compressibility. Since computing these exponents is, in general, very challenging, it would be expedient to directly simulate the dynamics of deformable particles in the flow of interest. Even so, the physical connection between the generalised Lyapunov exponents and particle stretching provides a rationale for understanding the behaviour of particles in different kinds of compressible flows.
We have limited our study to flows with just two values of the temporal correlation, corresponding to
$\mathit{Ku} = 0$
(delta-correlated) and
$\mathit{Ku} = 0.3$
. Our goal here was just to check whether the analytical predictions made using the delta-correlated flow continued to hold for a flow with a finite correlation time. The quantitative details of the particle dynamics will vary with
$\mathit{Ku}$
and an investigation of this dependence could be taken up in future work. (Our limited explorations suggest that the results will change significantly only for
$\mathit{Ku} \gg 1$
; this case, however, corresponds to a near steady flow and so is not of interest in the context of fluctuating flows.) A task that demands more immediate attention is the study of particle deformation in direct numerical simulations, say of a surface flow, to see how our qualitative predictions manifest in a turbulent suspension. Such simulations could also consider the feedback force exerted by the particle on the fluid, which will be relevant if the suspension is not ultra-dilute.
Another important question concerns the applicability of our study to complex particles that possess multiple deformation modes. We expect our results for
$\epsilon =0$
and
$g=1$
(elastic dumbbells) to extend to bead–spring chains with negligible interbead hydrodynamic interactions (Bird et al. Reference Bird, Curtiss, Armstrong and Hassager1987). Indeed, by applying an appropriate rescaling of the parameters (Jin & Collins Reference Jin and Collins2007), it has been shown that, in incompressible turbulence, the stationary PDF of the end-to-end extension of such chains is independent of the number of elastic links and can therefore be described by a vector model like the one used here (Watanabe & Gotoh Reference Watanabe and Gotoh2010; Vincenzi et al. Reference Vincenzi, Watanabe, Ray and Picardo2021). We expect this equivalence between elastic chains with different numbers of links to carry over to compressible random flows, since the rescaling of the chain parameters on which the equivalence is based is independent of flow compressibility. However, if interbead hydrodynamic interactions are significant, then the dynamics of multibead chains becomes qualitatively different from that of dumbbells, even in incompressible turbulence (Ganesh et al. Reference Ganesh, Vincenzi, Prabhakar and Picardo2026). In this case, the standard vector model will fail and bead–spring chain models will have to be considered; though it may be possible to capture the stretching behaviour of chains by including an extension-dependent relaxation time in the vector model (Hinch Reference Hinch1977; Celani, Puliafito & Vincenzi Reference Celani, Puliafito and Vincenzi2006; Prabhakar, Sevick & Williams Reference Prabhakar, Sevick and Williams2007). If detailed bead–spring or bead–rod models are used, then one can also refine the treatment of polymer scission by using a stochastic scission model (Sim, Khomami & Sureshkumar Reference Sim, Khomami and Sureshkumar2007).
Acknowledgements
J.R.P. acknowledges his Associateship with the International Centre for Theoretical Sciences (ICTS), Tata Institute of Fundamental Research, Bengaluru, India. D.R. and D.V. are grateful to the Université Côte d’Azur’s Center for High-Performance Computing (OPAL infrastructure) for providing resources and support.
Funding
This work was supported by the Indo-French Centre for the Promotion of Advanced Research (IFCPAR/CEFIPRA, project no. 6704-1) and by the French government through the France 2030 investment plan managed by the National Research Agency (ANR), as part of the Initiative of Excellence Université Côte d’Azur under reference number ANR-15-IDEX-01. D.V. acknowledges the support of Agence Nationale de la Recherche through project no. ANR-21-CE30-0040-01. D.R. acknowledges support from the Indian Institute of Science Education and Research Thiruvananthapuram (IISER TVM).
Declaration of interests
The authors report no conflict of interest.
Appendix A. A regularised version of the model of Olbricht et al. (Reference Olbricht, Rallison and Leal1982)
To study the statistics of
$R$
over the entire interval
$[0,R_{\textit{max}}]$
, it is necessary to modify (2.1) in such a way as to regularise the behaviour at small and large extensions. This can be achieved by introducing thermal noise, which sets the equilibrium size
$R_{{eq}}$
, and by assuming a nonlinear elastic force that enforces a maximum extension
$R_{\textit{max}}$
. Here we consider the FENE force, which is commonly used in bead–spring polymer models (Bird et al. Reference Bird, Curtiss, Armstrong and Hassager1987). The resulting evolution equation for
$\boldsymbol{R}(t)$
is
\begin{eqnarray} \nonumber \dot {\boldsymbol{R}} &= & \unicode{x1D734}\boldsymbol{\cdot }\boldsymbol{R} + g \left (\unicode{x1D64E}\boldsymbol{\cdot }\boldsymbol{R}- \dfrac {\epsilon }{\epsilon +1}\,\dfrac {\boldsymbol{R}\boldsymbol{\cdot }\unicode{x1D64E}\boldsymbol{\cdot }\boldsymbol{R}}{R^2}\,\boldsymbol{R}\right ) -\dfrac {f(R)\boldsymbol{R}}{2\tau _{\!p}(\epsilon +1)} \\ &&-\frac {\epsilon R_{{eq}}^2}{d(\epsilon +1)\tau _{\!p}}\,\frac {\boldsymbol{R}}{R^2} +\sqrt {\frac {R_{{eq}}^2}{\tau _{\!p} d}} \left [\unicode{x1D644} - \left (1-\sqrt {\frac {1}{1+\epsilon }}\right ) \dfrac {\boldsymbol{R}\otimes \boldsymbol{R}}{R^2} \right ]\boldsymbol{\cdot }\boldsymbol{\xi }(t), \end{eqnarray}
where
$\unicode{x1D644}$
is the identity matrix,
$d$
is the spatial dimension of the flow and
$\boldsymbol{\xi }(t)$
is a
$d$
-dimensional white noise. Equation (A1) is a stochastic differential equation, and the noise term is interpreted in the Itô sense; the coefficient of the noise is such that it satisfies the fluctuation–dissipation theorem (Schieber Reference Schieber1992; Hua & Schieber Reference Hua and Schieber1995). Here
$R_{{eq}}$
is the root mean square of
$\boldsymbol{R}(t)$
at equilibrium of an equivalent linearly elastic particle (
$R_{\textit{max}}\gg R_{{eq}}$
). For
$\epsilon =0$
and
$g=1$
, (A1) reduces to (2.4).
Appendix B. The PDF of the particle size in the compressible Batchelor–Kraichnan model
Here, we calculate the exact PDF of
$R$
, over the interval
$[0,R_{\textit{max}}]$
, and see how the power-law behaviour (discussed in § 4) appears over intermediate extensions
$R_{{eq}}\ll R\ll R_{\textit{max}}$
. We restrict ourselves to the case of
$\epsilon =0$
and
$g=1$
, for simplicity, and consider (2.4) and (2.5). As discussed in § 4, (2.4) is equivalent to the Itô stochastic differential equation:
\begin{align} \dot {\boldsymbol{R}}=\unicode{x1D642}(t)\boldsymbol{R}-\frac {f(R)\boldsymbol{R}}{2\tau _{\!p}} + \sqrt {\dfrac {R_{{eq}}^2}{\tau _{\!p} d}}\,\boldsymbol{\xi }(t), \end{align}
where
$\unicode{x1D642}(t)$
has been defined in (4.5). The PDF of
$\boldsymbol{R}$
at time
$t$
, denoted as
$\mathscr{P}(\boldsymbol{R},t)$
, satisfies the Fokker–Planck equation (Plan et al. Reference Plan, Ali and Vincenzi2016):
\begin{align} &\hspace{3pc} =D(d+1-2\wp )\dfrac {\partial ^2}{\partial R_i\partial R_i}(R^2\mathscr{P})+2D(\wp d-1) \dfrac {\partial ^2}{\partial R_i\partial R_{\!j}}(R_i R_{\!j}\mathscr{P})\\ \nonumber & \hspace{4pc} + \dfrac {1}{2\tau _{\!p}}\,\dfrac {\partial }{\partial R_i}[f(R)R_i\mathscr{P}] + \dfrac {R_{{eq}}^2}{2\tau _{\!p} d} \,{\nabla} ^2_{\boldsymbol{R}} \mathscr{P}, \end{align}
where summation over repeated indices is implied. The correlation tensor
$\unicode{x1D63E}$
is defined in (4.3) and depends on the space dimension of the flow,
$d$
, and the degree of compressibility
$\wp$
. Assuming that the initial condition of (B2) is independent of the orientation of
$\boldsymbol{R}$
, we expect the solution of (B2) to retain this property at all times, since the velocity gradient is statistically isotropic. To take advantage of this fact, we transform from Cartesian variables
$\boldsymbol{R}$
to azimuthal variables
$(R,\varPsi )$
, where
$\varPsi$
denotes a suitable set of angles in
$d$
dimensions. The PDF
$\widetilde {\mathscr{P}}(R,\varPsi ,t)$
is connected to
$\mathscr{P}(\boldsymbol{R},t)$
via the relation
$\widetilde {\mathscr{P}}(R,\varPsi ,t)=J \mathscr{P}(\boldsymbol{R},t)$
, where
$J$
is the Jacobian of the transformation between the two sets of variables and is of the form
$J=R^{d-1}h(\varPsi )$
with
$h(\varPsi )$
being a function of the angles. Owing to statistical isotropy, the dependence of
$\widetilde {\mathscr{P}}(R,\varPsi ,t)$
on the angles is contained entirely in
$J$
; therefore,
$\widetilde {\mathscr{P}}(R,\varPsi ,t)=P(R,t)h(\varPsi )$
, where
$P(R,t)$
is the marginal PDF of
$R$
. This implies that
$\mathscr{P}(\boldsymbol{R},t)=R^{1-d}P(R,t)$
. The Fokker–Planck equation for
$P(R,t$
) is obtained from (B2) by changing the variables in the partial derivatives, discarding the derivatives with respect to the angular variables and replacing
$\mathscr{P}$
with
$R^{1-d}P$
:
where
where
$\lambda$
and
$\varDelta$
are defined in (4.8) and depend on both
$d$
and
$\wp$
. With reflecting boundary conditions at both ends of the interval
$0\leqslant R \leqslant R_{\textit{max}}$
(which set the probability current
$\partial _{R} (\mathfrak{D}_2 P) - \mathfrak{D}_1 P$
to zero at
$R=0$
and
$R=R_{\textit{max}}$
), (B3) has the stationary solution
$P_{\textit{st}}(R) = {\rm e}^{-\varPhi _a(R)}$
, where
$\varPhi _a(R)$
has been defined in (4.13). Using the expressions of the drift and diffusion coefficients in (B4) yields
\begin{align} P_{\textit{st}}(R) \propto R^{d-1} \left (1+ {\rm d}\varDelta \tau _{\!p}\dfrac {R^2}{R_{{eq}}^2}\right )^{- ({\alpha +d})/{2}}\left (1-\frac {R^2}{R_{\textit{max}}^2}\right )^{ ({\alpha }/{2})+ ({\lambda }/{\varDelta })} \end{align}
with
The stationary PDF of
$R$
consists of three factors, each of which dominates in a different range of sizes. For
$R\ll R_{{eq}}$
, the dynamics is dominated by molecular noise and
$P_{\textit{st}}(R)\sim R^{d-1}$
. For sizes near to
$R_{\textit{max}}$
, the third factor on the right-hand side of (B5) makes
$P_{\textit{st}}(R)$
rapidly decay to zero. The nonlinear elastic force indeed diverges as
$R$
approaches
$R_{\textit{max}}$
, so as to overcome the action of the strain rate and keep the size of the particle bounded. In the intermediate range
$R_{{eq}}\ll R\ll R_{\textit{max}}$
, the stationary PDF behaves as
$R^{-1-\alpha }$
with
$\alpha$
given in (B6). It is interesting to note that, in the limit
$b\to \infty$
, the exponent
$\alpha$
becomes
as predicted by the large deviations theory with a linear elastic force; see (3.7) (indeed, (B7) becomes (4.9) upon substituting
$\lambda$
and
$\varDelta$
from (4.8) and setting
$\epsilon =0$
and
$g=1$
). This fact shows that, in the range
$R_{{eq}}\ll R\ll R_{\textit{max}}$
, the statistics of the particle size are determined by the balance between the actions of the strain rate and the elastic restoring force, which can be assumed to be linear in this range. Molecular noise and the nonlinear elasticity of the particle do not affect the statistics of moderate extensions, and ultimately they only play the role of small- and large-
$R$
regularisations of the dynamics.
Here, we have restricted the derivation to the case of
$\epsilon =0$
and
$g=1$
, for simplicity. The general case could be solved in a similar manner by considering the Fokker–Planck equation associated with (A1), with
$\unicode{x1D64E}$
and
$\unicode{x1D734}$
having the statistics of the symmetric and anti-symmetric parts of the Batchelor–Kraichnan velocity gradient.
Appendix C. Zero
$\boldsymbol{\tau} _{\!\boldsymbol{f}}$
limit of the renewing flow
Consider the equation
with
$\unicode{x1D642}(t)$
as in (5.2). In each time interval
$\mathcal{I}_n$
, its solution is
Hence,
As
$D\tau _{\!f}\to 0$
,
$\unicode{x1D642}(t)$
tends to the
$d\times d$
white noise defined in (4.5), while
$\boldsymbol{\ell }(t)$
tends to the time-ordered exponential
which is the solution of the Stratonovich equation
To recover the Stratonovich equation equivalent to (4.4) in the limit
$D\tau _{\!f}\to 0$
, it is therefore necessary to subtract the drift term
$4\wp D \boldsymbol{\ell }$
from the right-hand side of the evolution equation for
$\boldsymbol{\ell }(t)$
(see (4.6)).


℘=0
0<℘<℘c
℘c<℘⩽1
L(α)
α/(2τp)
τp
L(q)/D
q
gα/(1+ϵ)
gDτp
gα/(1+ϵ)
gλτp
℘
d=2
℘
℘=0
0<℘
(d/4<℘⩽1)
Pst(R)
R
τp (τp=0.05)
τp (τp=2)
ϵ=0
g=1
Req=1.0
Rmax=100
τp
R
α
℘
τp
℘
DTbr
Dτp
℘
DTbr
Dτp
℘
Dτp
DTbr
d=2
ϵ=0
g=1
Rbr/ρ=10
L(q)/|λ|
q
Dτf=0.01
Dτf=0.2
℘=0
℘=0.25
℘=1
λ/D
℘
L(q)/D
q
gα/(1+ϵ)
gDτp
℘
τf=0.2
Pst(R)
Req≪R≪Rmax
R(t)
Req=1
Rmax=103
ϵ=0
g=1
R(t)
dt=5×10−4
R
Rbr/ρ=10
Dτp
Rbr/ρ
Dτp=1
℘
c
c=−3.4
℘=0.5
c=−7.2
℘=1.0
ϵ=0
g=1