4.1. Disc shape

We will first model the evolution of the disc shape. As shown in figures 2(a) and 2(b), at the onset of the transition, discs tend to wrap and form cylinder like particles. The exact shape of the disc is given by the equation of elasticity

(4.1)\begin{equation} \sigma \partial_{tt} r + B \nabla^2\nabla^2 r = \xi, \end{equation}

where $r$ is the position of a point on the disc, $\partial _{tt}$ the second temporal derivative, $\sigma =\rho e$ the surface density and $B=Ee^3/12(1-\nu )$ the bending modulus of the disc, with $\nu \simeq 0.5$ the Poisson coefficient of silicone, and $\xi$ the hydrodynamic stress. The disc being freely advected and as no external torque is applied at the boundary, the boundary condition is $\kappa (r=R)=0$.

To model the evolution of the shape, we consider a simple solution where a flat disc lies in the $x$–$z$ plane and can only deform in the $y$ direction according to

(4.2)\begin{equation} y_d(x)=\frac{K_0}{2}x^2\left(1-\frac{x^2}{6x_R^2}\right), \end{equation}

where $K_0$ is the curvature at the centre of the disc and $x_R$ the lateral extension of the disc. For an inextensible disc, $x_R$ is defined by the arc length which should be equal to the disc radius

(4.3)\begin{equation} R=\int_0^{x_R}\textrm{d}s=\int_0^{x_R}\sqrt{{\textrm{d} x}^2+{\textrm{d} y}^2}=\int_0^{x_R}\sqrt{1+K_0^2x^2\left(1-\frac{x^2}{3x_R^2}\right)^2}\,{\textrm{d} x}. \end{equation}

Representations of a flat and a bent disc are shown in figure 3. The $(x,y,z)$ axes being the principal axes of the object, the eigenvalues of the moment of inertia of the bent disc are given by

(4.4)$$\begin{gather} I_x=\int_{x_R}^{x_R} {\textrm{d} x} \int_0^{\infty} {\textrm{d} y}\int_{{-}z_m}^{z_m}\textrm{d}z ((y-y_G)^2+z^2) \delta\left(y-y_d(x)\right), \end{gather}$$

(4.5)$$\begin{gather}I_y=\int_{x_R}^{x_R} {\textrm{d} x} \int_0^{\infty} {\textrm{d} y}\int_{{-}z_m}^{z_m}\textrm{d}z (x^2+z^2) \delta\left(y-y_d(x)\right), \end{gather}$$

(4.6)$$\begin{gather}I_z=\int_{x_R}^{x_R} {\textrm{d} x} \int_0^{\infty} {\textrm{d} y}\int_{{-}z_m}^{z_m}\textrm{d}z ((y-y_G)^2+x^2) \delta\left(y-y_d(x)\right), \end{gather}$$

where $z_m(x)$ is the maximum height of the disc at the position $x$ and is given by

(4.7)\begin{equation} z_m^2(x)=R^2-\left(\int_0^x\,\textrm{d}s\right)^2, \end{equation}

and $y_G=3 K_0x_R^2 /20$ is the $y$-coordinate of the centre of gravity. The exact expressions of the moments of inertia are relatively tricky to compute analytically. We then compute the evolution of the shape factors numerically for different $K_0$, cf. figure 4(a). The theoretical evolution of the bent disc in the $(\lambda _1/\lambda _3, \lambda _2/\lambda _3)$ plane is represented by the plain line in figure 2(b). As one can see, the general trend is well captured by this model, but there is an offset between the measurements and the theoretical prediction. Two sources of discrepancy can be identified. First, the number of voxels forming the disc is much smaller than the number of voxels used to compute the theoretical moment of inertia. A second source of error is related to the convex hull volume method. As detailed previously, for voxels projecting onto the same pixel only the voxel closest to the camera is stored. The reconstructed disc is then, in general, not a disc but a partially filled cylinder, as illustrated by figure 1(d). To estimate the order of magnitude of this last source of error, we can compute the shape factor of an object made of touching two coaxial cylinders: one filled of radius $R$ and thickness $h$, and one empty cylinder of inner/outer radii $R_1$ and $R$ with the same height $h$. Here, $h$ is given by the voxel size and we take $R_1=R-h$. The evolution of $\lambda _2/\lambda _3$ for this object is shown in figure 4(b). The value of the shape factor for $h=0$ is equal to the one of a flat disc and increases when $h$ is increased. As mentioned previously, the exact value of the offset cannot be determined rigorously and has been considered as a constant free parameter to fit the data. The dashed line in figure 2(b) is obtained by shifting the theoretical prediction for a disc by 0.1 in the vertical direction ($\langle \lambda _2\rangle /\langle \lambda _3\rangle$) and 0.01 in the horizontal direction ($\langle \lambda _1\rangle /\langle \lambda _3\rangle$). These coefficients will remain fixed in the following sections.

Figure 3.

Three-dimensional representation of an undistorted disc $K_0=0$, (a), and of a bent disc, $K_0 R=1$, (b).

Figure 4.

(a) Theoretical evolution of the shape factors of a bent disc represented in figure 3 as a function of the normalised curvature $K_0R$. (b) Evolution of the mean shape factor $\langle \lambda _2 \rangle /\langle \lambda _3 \rangle$ for a partially empty cylinder.

Here, the agreement between our model and our measurements suggests that near the transition the disc bend and adopt a U-shape which can be described by (4.2). All the physics of the deformation are then hidden in the parameter $K_0$, which depends on the disc deformability and the turbulent intensity. The following sections aim to model these dependencies.

4.2. A power budget argument

As shown in figure 4(a), if $K_0R\lesssim 10^{-2}$ the deformations are weak and the disc can be seen as rigid, whereas if $K_0R>10^{-2}$ deformations are important and the disc is indeed flexible. The curvature $K_0$ then plays a role similar to the persistence length $\ell _p$ for the deformation of worm-like chain polymers. The stiffness of polymers is characterised by the ratio of their length $L_p$ to the persistence length $\ell _p$, defined by $\langle t(s)\cdot t(s+\ell ) \rangle = \textrm {e}^{-\ell /\ell _p}$, where $t(s)$ is the tangent vector at the position $s$ (Yamakawa Reference Yamakawa1971). If $L_p\gg \ell _p$, deformations are important. On the contrary, if $L_p\ll \ell _p$, the polymer remains straight and is considered stiff. In polymer theory, the order of magnitude of the persistence length is given by the balance of the thermal energy $k_B T$ with the elastic energy $\mathcal {E}_{el}\sim EI/\ell _p$, where $k_B$ is the Boltzmann constant, $T$ the temperature and $I$ the area moment of inertia. Brouzet et al. (Reference Brouzet, Verhille and Le Gal2014) draw an analogy between this system and flexible fibres distorted by a turbulent flow. They showed that the energy balance was not able to capture the transition from rigid to flexible fibres. The reason is that the correlation time of turbulence is, in general, of the same order of magnitude as the dynamical time scale of the deformations. Therefore, contrary to polymers, the forcing cannot be modelled by delta correlated noise. The energy budget has then to be replaced by a balance of power. Indeed, the deformation of a fibre in a stationary turbulent flow is an out-of-equilibrium stationary process. The mean elastic energy stored by the fibre $\mathcal {E}_{el}$ is then constant and

(4.8)\begin{equation} \frac{\textrm{d} \mathcal{E}_{el}}{\textrm{d}t}=P_{inj}-P_{dis}=0, \end{equation}

where $P_{inj}$ is the power injected by the turbulence and $P_{dis}$ is the power dissipated by viscosity. For fibres whose lengths are in the inertial range, deformations are due to eddies of similar sizes. The turbulent power $\mathcal {P}_{turb}$ is then given by $\mathcal {P}_{turb}\sim \rho _f L_p^3 \varepsilon$. The dissipative term should be proportional to the elastic energy $\mathcal {E}_{el}$ stored in the fibre and should depend on the typical time scale of the deformation $\tau _{def}$. The ratio of these two terms defines an elastic power $\mathcal {P}_{el}$. For fibres, they showed that the time scale of the deformation is given by the relaxation time defined by the ratio of the viscous stress to the bending stress.

Following a similar argument, one can assume that the transition from a rigid to a flexible disc is given by the balance of the turbulent power

(4.9)\begin{equation} \mathcal{P}_{turb}\sim \rho_f R^3 \varepsilon, \end{equation}

and the elastic power

(4.10)\begin{equation} \mathcal{P}_{el,v}=B \kappa^2 R^2/\tau_{B}, \end{equation}

where $\tau _{B}$ is the relaxation time scale of the deformation given by the balance of the viscous stress $\xi _v$ with the bending stress $B\nabla ^2\nabla ^2 r$ in (4.1). From dimensional analysis, at small Reynolds number, $\xi _v$ scales as

(4.11)\begin{equation} \xi_v \sim \frac{\mu}{R} (u_f-\partial_t r), \end{equation}

where $u_f$ is the fluid velocity and $\mu$ the dynamical viscosity of the fluid. The relaxation time is then

(4.12)\begin{equation} \tau_B \sim \mu R^3/B, \end{equation}

and, assuming that at the transition $\kappa$ can be replaced by $1/R$ in (4.10), the elastic power scales as

(4.13)\begin{equation} \mathcal{P}_{el,v}\sim \frac{B^2}{\mu \rho_f R^3}. \end{equation}

The balance of power $\mathcal {P}_{turb}=\mathcal {P}_{el,v}$ defines a critical turbulent dissipation rate $\varepsilon _v$

(4.14)\begin{equation} \varepsilon_v = \frac{B^2}{\rho_f \mu R^6}. \end{equation}

If this model were able to capture the physics of the deformation, the evolution of the mean shape factors should only depend on the ratio $\varepsilon /\varepsilon _v$. To compare the theoretical evolution of the shape factor with this prediction one needs to relate the curvature $K_0$ to $\varepsilon$. Using (4.9), (4.10) and (4.12), which define the turbulent and elastic power, one can show that, near the threshold, the typical curvature $\kappa$, and therefore $K_0$, should be proportional to $\varepsilon ^{1/2}$. The comparison of the theoretical and measured evolution of the mean shape factors as a function of $\varepsilon /\varepsilon _v$ is shown in figures 5(a) and 5(b). The dashed line is the best fit of the experimental measurements considering $K_0$ as a free parameter. Here, the global trend is well captured by the theoretical prediction. However, the estimation of the threshold with this model is not valid as evidenced by the scattering of the experimental points around the theoretical prediction. In particular, the threshold of the discs having the highest densities (red stars and blue circles) seems to be at least one order of magnitude lower than the model prediction.

Figure 5.

Evolution of the mean shape factors $\langle \lambda _1\rangle /\langle \lambda _3\rangle$ (a) and $\langle \lambda _2\rangle /\langle \lambda _3\rangle$ (b) as a function of the normalised turbulent dissipation rate $\varepsilon /\varepsilon _v$. The dashed line represents the best fit of an ideal bent disc; see text for more details.

In this model, inertia is always neglected. It is then not surprising that the dependency of the transition on the disc density is not captured for the higher densities. Another model, derived by Rosti et al. (Reference Rosti, Banaei, Brandt and Mazzino2018), considers the density of fibres. The transposition of their model to discs is presented in the following section.

4.3. A temporal argument

The second model used to quantify fibre deformation is based on a temporal argument. If the forcing time scale is smaller than the deformation time scale, deformations are weak as fibres do not have the time to adapt their shapes to the flow. On the contrary, if the forcing time scale is larger than the deformation time scale, the fibre is distorted by the coherent structures of the flow.

As for the model presented in § 4.2, the transposition of this idea to the current study requires an estimation of the deformation time scale of the discs. Two time scales are relevant depending on the influence of particle inertia. The first one is the relaxation time scale, as defined by (4.12). It is used to define the Weissenberg number, which generally quantifies polymer deformations in turbulent flows (Vincenzi et al. Reference Vincenzi, Watanabe, Sankar Ray and Picardo2021). The second one is the resonant frequency, which has been used by Rosti et al. (Reference Rosti, Banaei, Brandt and Mazzino2018) to quantify fibre deformations in the underdamped regime, i.e. when viscous stress is negligible as compared with particle inertia. As shown previously, the relaxation time is independent of the particle density. Therefore, a model solely based on this time scale should not be able to describe our measurements.

The resonant frequency of the disc $\omega _d$ is given by the equation of elasticity (4.1) by balancing the inertial term $\sigma \partial _{tt}r$ and the bending term $B\nabla ^2\nabla ^2r$

(4.15)\begin{equation} \omega_d^2 = c_\omega \frac{B}{\sigma R^4}, \end{equation}

where $c_\omega$ is a constant which depends on the excited mode. In the following we will consider the case where $c_\omega =1$. When the disc radius $R$ is in the inertial range of turbulence, the relevant time scale for the forcing $\tau _{turb}$ is the typical time of eddies of similar size

(4.16)\begin{equation} \tau_{turb}\sim R/u_R\sim R^{2/3}\varepsilon^{{-}1/3}, \end{equation}

where $u_R\sim (\varepsilon R)^{1/3}$ is the typical velocity at scale $R$. From these two time scales, one can define an inertial Weissenberg number $Wi=\omega _d \tau _{turb}$ which should be of order unity at the transition. This defines a critical turbulent dissipation rate $\varepsilon _{r}$ above which a disc can be bent

(4.17)\begin{equation} \varepsilon_r = \frac{B^{3/2}}{\sigma^{3/2} R^4} . \end{equation}

Contrary to the previous model, this model cannot predict the evolution of the mean curvature $\kappa$ as a function of $\varepsilon$. There was relatively good agreement between the shape of the theoretical prediction and the global trend of the measurement in the previous scaling argument, so we will assume that $\kappa \sim \varepsilon ^{1/2}$. Once again, $K_0$ will be considered as a free parameter in order to fit the experimental data.

The measurements and the theoretical predictions are compared in figures 6(a) and 6(b). As for the previous case, there is a scattering of the experimental points around the theoretical prediction for both shape factors. The failure of this approach is not surprising as it was derived to describe the transition between two modes of deformation of fibres: the excitation of the first bending mode if $\omega _d\tau _{turb}\ll 1$ and the deformation dependent upon the coherent structures of the flow if $\omega _d\tau _{turb}\gg 1$ (Rosti et al. Reference Rosti, Banaei, Brandt and Mazzino2018). On the contrary, we focus here on the transition from rigid to flexible discs which should occur at smaller $\varepsilon$.

Figure 6.

Evolution of the mean shape factors $\langle \lambda _1\rangle /\langle \lambda _3\rangle$ (a) and $\langle \lambda _2\rangle /\langle \lambda _3\rangle$ (b) as a function of the normalised turbulent dissipation rate $\varepsilon /\varepsilon _r$. The dashed line represents the best fit of an ideal bent disc; see text for more details.

4.4. An inertial power balance argument

In the first model presented in § 4.2, the transition from rigid to flexible discs was modelled using a balance between the turbulent power $\mathcal {P}_{turb}$ and the elastic power $\mathcal {P}_{el}=B\kappa ^2R^2/\tau _{el}$ assuming that the time scale of the deformation $\tau _{el}$ was the relaxation time scale, cf. (4.12). In the second model (§ 4.3), we introduce a second time scale for the deformation: the frequency of the first bending mode, cf. (4.15). The relevance of each time scale depends on the importance of the disc inertia in the equation of elasticity (4.1). This is quantified by the ratio of the inertial term $\sigma \partial _{tt} r\sim \sigma \omega _f^2\zeta$, where $\omega _f$ is the forcing frequency and $\zeta$ the typical displacement, to the viscous term $\mu \partial _{t}r/R\sim \mu \omega _f\zeta /R$. This ratio defines an elastic Stokes number $St_d$

(4.18)\begin{equation} St_d \sim \frac{\rho}{\rho_f} \left(\frac{e}{\eta_K}\right)^{4/3} \left(\frac{R}{e}\right)^{1/3}. \end{equation}

Here, we have assumed that the forcing time scale is given by the eddies at the scale of the disc $\omega _f\sim u_R/R$. In this study, $St_d$ varies between 14 and 320. Inertia is then larger than the viscous dissipation and the time scales of the deformation should be given by $\omega _d$. This regime corresponds to the underdamped regime investigated by Rosti et al. (Reference Rosti, Banaei, Brandt and Mazzino2018) for fibres. Replacing the relaxation time scale by $\omega _d^{-1}$ in the definition of the elastic power leads to

(4.19)\begin{equation} \mathcal{P}_{el,w}\sim \frac{B^{3/2}\kappa^2}{\sigma^{1/2}}. \end{equation}

As for the first model, balancing this elastic power with the turbulent power $\mathcal {P}_{turb}=\rho R^3\varepsilon$, we find that $\kappa \sim \varepsilon ^{1/2}$ near the threshold. Moreover, the critical turbulent dissipation rate $\varepsilon _w$ can be estimated by defining the threshold of the transition by $\kappa \sim 1/R$

(4.20)\begin{equation} \varepsilon_w = \frac{B^{3/2}}{\rho_f \sigma^{1/2} R^5}. \end{equation}

The evolution of the shape factors as a function of the normalised turbulent dissipation rate $\varepsilon /\varepsilon _w$ is shown in figures 7(a) and 7(b). There is now a good agreement between the measurements and the theoretical predictions for all the discs except for the most flexible ones. In that case, deformations probably involved several modes of bending. This cannot be captured by our model and will be discussed in the following section.

Figure 7.

Evolution of the mean shape factors $\langle \lambda _1\rangle /\langle \lambda _3\rangle$ (a) and $\langle \lambda _2\rangle /\langle \lambda _3\rangle$ (b) as a function of the normalised turbulent dissipation rate $\varepsilon /\varepsilon _r$. The dashed line represents the best fit of an ideal bent disc; see text for more details.