{'fraction_non_alphanumeric': 0.03097713097713098, 'fraction_numerical': 0.012751212751212751, 'mean_word_length': 4.372486969471333, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 1, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 1, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Linear transient growth of particulate pipe flows

We tackle the question of whether the presence of particles in a pipe flow can influence the linear transient growth of infinitesimal perturbations, in view of better understanding the behaviour of particulate pipe flows in regimes of transition to turbulence. The problem is tackled numerically by means of a simple model where particles are modelled as a second fluid, that interacts with the fluid phase through a two-way Stokes drag only. The transient growth is found to be enhanced by the presence of particles, especially so if the particles are localised at a specific radial location of the pipe. At the same time, the mechanisms of transient growth themselves remain those of non-particulate flows. The effect is maximised for particles of intermediate size (somewhere between the ballistic limit of very light particles and the point where they are too heavy to be influenced by the flow). Most remarkably, the Segre-Silberberg radius (around 2/3 of the pipe radius), where particles naturally cluster in laminar flows, turns out to be close to the optimal location to enhance the transient growth.


Introduction
This paper focuses on the topic of transition to turbulence for particulate flows. Transition to turbulence has been a widely studied topic since Reynolds first documented the phenomenon experimentally (1; 2). While much of the research on the topic has focused on single phase flows, there is a growing interest for particulate flows, due to their many applications. Examples range from the precise determination of the volume fraction of oil in the oil-water-sand-gas mixture that is extracted from offshore wells, to needs in the food processing industry (3), and flows of molten metal carrying impurities during recycling processes (4).
Transition to turbulence is, even for single phase flows, a complex problem. In the case of the pipe flow, there is no clearly defined critical Reynolds number. The problem is even more complex in the case of the particulate flow, due to the large number of parameters to consider for the solid phase. To the complexity of the dynamics involved is added the inherent difficulty of considering, theoretically or numerically, a large number of independent objects. Nonetheless, experimental (5; 6) and theoretical (7; 8; 9) knowledge has been amassed on the topic of transition to turbulence in particulate pipe flows. For the particles influence on the pipe flow stability in particular, the effect on the transition to turbulence depends non-trivially on the size and volume fraction of the particles.
Matas observed that transition occurs at lower flow rates after the addition of particles for small particles, while the effect is reversed for large particles. In a similar fashion, particles tend to have a destabilising effect on the pipe flow stability for low particle volume concentrations but a stabilising one for high volume concentrations (6). Numerical simulations based on accurate modelling of individual solid particles recovered this phenomenology for pipe flows (10).
The present knowledge on this topic is mostly empirical and there is a need for a better understanding of the underlying mechanisms underlying the transition to turbulence of particulate pipe flows. A previous study of the linear stability of particulate pipe flows uncovered a mechanism for instability (11). However, even when the flow is stable to infinitesimal disturbances, interactions between disturbances and the underlying flow can lead to large distortions of the base flow due to the non-normality of the linearised equations. Perturbations can experience large growth at finite time (12; 13), a phenomenon generally referred to as transient growth. This paper aims to further this understanding using linear transient growth analysis in order to study the flow behaviour below the critical Reynolds number.
The paper starts by an introduction of the two-fluids model used and the assumptions it relies on, the details of the variational method used to obtain the transient growth as well as the numerical methods used (section 2). We then consider the envelope of the optimal gain as a function of the time in section 3. Section 4 focuses on the effect of the Reynolds number and particle size on the optimal gain for homogeneously distributed particles. It is expanded for the case of a nonhomogeneous particle distribution in section 5. Finally, the topology of the optimal perturbations is studied in section 6.

Model and governing equations
The complexity of particulate flows means they are usually studied through modelling assumptions and approximations in order to simplify the problem while keeping as much of the underlying dynamics as possible. In general, a trade-off has to be struck between how accurately the model represents the particulate flow and its complexity. Models with an accurate particle description, such as fully Lagrangian models (14; 15; 16) and immersed boundary methods (17; 18; 19), suffer from a high computational cost. In order to avoid the computational cost incurred when accounting for particles as individual solids, here we describe the particulate flow using the 'two-fluid' model first derived in (7). The fluid and solid phases are treated as two inter-penetrating media, with the particles being described as a continuous field rather than as discrete entities. It is a two-way coupled model, taking into account the feedback of the solid phase on the fluid. On the other hand, particle-particle interactions such as collisions or clustering, as well as the deflection of the flow around the particles, are neglected. The two-fluid model is therefore valid for lower concentrations and in the limit where particles are sufficiently smaller than the characteristic scale of the flow. This model has been used in the context of channel flow (9; 20) and boundary layer flow (21).

Two-fluid model
We consider the flow of a fluid of density ρ f and dynamic viscosity µ through a straight pipe with a constant circular cross-section of radius r 0 and driven by a constant pressure gradient. The fluid carries spherical particles of radius a.
The particles are treated as a continuous field with a spatially varying density N . Their motion is coupled to the fluid via a Stokes drag force, S d = 6πaµ|u p − u|, where u and u p are the fluid and particle velocities respectively. When working with an averaging method, one has to ensure that the system of equations is closed. If only the Stokes drag is considered, no specific correction is required (22). The Stokes drag force is proportional to a. On the other hand, other forces commonly considered (such as virtual mass force, buoyancy, Magnus force, Saffman force and the Basset history force) are all quadratic or above in particle radius and so can in general be neglected. The Stokes force, by contrast, may become significant if the background shear is large rather than O(1) (23). Similarly, buoyancy is proportional to ρ f a 3 , regardless of which exact definition is chosen, and becomes vanishingly smaller than the Stokes drag S d = 6πaµ|u p − u| in the limit a → 0. More details on the relevance of the drag-only, two-fluid model used in this paper can be found in (11).
We use, the standard cylindrical set of coordinates (r, θ, z) aligned with the pipe, with the respective velocity components: u = (u r , u θ , u z ) and u p = (u pr , u pθ , u pz ) for the fluid and particulate phases. Where relevant, we distinguish quantities associated with the particles from those associated with the fluid by means of a subscript p. The fluid velocity is described using the standard Navier-Stokes set of equations to which a Stokes drag force is added to account for the interaction between fluid and solid phases. The solid phase is characterised by the conservation of the particles momentum and density. Nondimensionalising by the centreline velocity, U 0 , the pipe radius, r 0 , and the fluid density ρ f yields the following set of governing equations: where p is the flow pressure and N the local particle concentration. This system is governed by three non-dimensional parameters: the Reynolds numbers Re = U 0 r 0 ρ f /µ, the aforementioned Stokes number, which expresses a dimensionless relaxation time S = 2a 2 ρ p /9r 2 0 ρ f and the mass concentration f = m p /m f , corresponding to the ratio between the particles and fluid mass over the entire pipe. N is normalised such that N dV = 1. For a given position x, N (x) > 1 implies that the local concentration of particles is higher than the pipe average.
These equations satisfy an impermeable and no-slip boundary condition for the fluid and a no penetration boundary condition for the radial particle velocity The stability of the flow is studied through the addition of a small perturbation to the steady solution, U = U p = (1 − r 2 )ẑ : where N 0 is the average local particle concentration. Linearising equations (1) -(4) around this base state yields: The boundary conditions for the perturbations u ′ and u ′ p are the same as for the full flow, u and u p . From here on the primes are dropped for the sake of readability.
The gain corresponds to the ratio between the maximal energy a perturbation can have at a time T and the perturbation initial energy: The perturbation u(0) is the one causing the largest amount of growth, and is often referred to as the optimal disturbance. By optimising G(T, Re) over T , one can find the maximum possible gain, or optimal gain, at a given Reynolds number. This paper focuses on the optimal gain and the associated time of occurrence. A variational method approach, adapted from the single phase flow problem (24), is used to solve this optimisation problem. The problem described by equations (7)-(10) can be characterised with the following functional L: where λ, Υ, Υ p , Γ and Π are the Lagrange multipliers enforcing the constraints of the problem: λ enforces that the energy is fixed, Υ, Υ p and Γ enforce that Equations (7) and (8) hold true over t ∈ [0, T ], Π and Γ enforces the incompressibility of the flow and the conservation of the total number of particles. The brackets represent a normalised volume integral over the pipe, given any Finding the initial perturbation that will maximise energy growth is equivalent to maximising L, done here through finding the root of its variational derivative δL. By reordering δL, one can obtain the adjoint system of equations of our problem, with an additional set of conditions. The adjoint system of equation is: where Υ and Υ p are the adjoint fluid and particles velocities respectively; Γ is the adjoint particle local concentration while Π is the adjoint pressure. The adjoint equations must be true for δL to be equal to 0. Enforcing δL = 0 yields another set of conditions: In this paper we consider homogeneous and nonhomogeneous particle distributions. In the case of a homogeneous particle distribution, N 0 is held constant spatially. However, particles are not necessarily uniformly distributed in practice. In particular, they tend to aggregate in the radial direction, around r = 0.6 − 0.8 (5; 25). We parametrise this phenomenon by mean of a particle distribution of the form withÑ chosen such that 1 0 N 0 (r)rdr = 1. The distribution is then, in the radial direction, a Gaussian centred around radius r d with a standard deviation σ. N 0 is still homogeneous in the axial and azimuthal directions.
A point of note is that, as opposed to the single phase pipe flow which is well-known to be linearly stable, particulate flow can, within our theoretical framework, be linearly unstable in the case of nonhomogeneous particle distributions (11). However, only linearly stable cases are considered in this work.

Iterative variational method
We use an iterative procedure to minimise δL, akin to the one used in (26). Initially, a first guess of the initial perturbation is made for the fluid velocity u (0) (t = 0) = u 0 (0) and the particles velocity u p (0) (t = 0) = u p0 (0) . The initial perturbations of the fluid and solid phases for iteration (i + 1) are: for the fluid velocity, and: for the particle velocity. It entails that Υ(0) and Υ p (0) have to be computed for each iteration. To that effect, the iteration process is as follows: • At the i-th iteration, equations (7)-(10) are advanced from t = 0 until a target time t = T is reached in order to obtain u (i) (T ) and u p (i) (T ).
• The adjoint system of Equations (13) give the gradient in u 0 and u p0 and the initial conditions are updated as where ǫ is the step size.
The process is repeated until the norms of ∂L/∂u 0 and ∂L/∂u p0 are less than a threshold chosen for convergence.

Computational method
The code is derived from a standard DNS code (27). Temporal discretisation is done through a predictor-corrector scheme. The spatial discretisation is done using a fourth order finite difference method in the radial direction and Fourier spatial discretisation with 128 mesh points used in the azimuthal and streamwise directions. Any field g can then be written as: where α and m are the wavenumbers in the streamwise and azimuthal directions respectively. The numerical code has been modified in order to incorporate the solid phase, using a fully Eulerian method. First, we add a set of equations for the particles velocity, u p , for both the standard and adjoint problems (Equations (8) and (14) respectively). Initial and boundary conditions for the particle velocity (equation 6) are added as well. The initial fluid velocity is obtained from a previously saved state.  Figure 1: Optimal gain G f (green) and optimal time of gain T f (red) for the single phase flow as a function of the Reynolds number. The points correspond to values obtained using our code. The lines corresponds to the scaling given in (28):

Code validation and convergence
The code has been first verified against the literature on the single phase pipe flow, which is simulated by fixing the particle mass concentration f to 0. (29) found that the time of the peak in energy increases linearly with the Reynolds number while the optimal gain scales with Re 2 for all modes, with G f = 72.40Re 2 × 10 −6 and T f = 48.77Re × 10 −3 (28). These scalings have been recovered with our code, as illustrated in figure 1. Second, the growth rate of the leading eigenvalue obtained through a linear stability analysis of the system of equations (7)-(10) is proportional to the energy decay rate of the linear DNS at large times and therefore offers a convenient way to test the long term evolution of individuals modes in the DNS code. Table 1 and Table 2 show the leading eigenvalue found with linear stability analysis and LDNS simulation, for a single phase and particulate flow respectively. The normalised error is always below 10 −3 . The difference between the linear stability analysis and linear DNS results is not increased by the addition of particles. Figure 2 shows the difference between the values for optimal gain obtained for a given number of iterations, and a fully converged value, G (500) . The growth G is shown to converge as the process is iterated, reaching fully converged values after a sufficient number of iterations in the three cases considered in figure 2. The number of iterations needed to fully converge depends on the case, convergence is significantly faster in the case of a single phase flow, where 30 iterations are typically needed to reach machine precision; whereas in the case of particulate flows, this number varies between 80 and 100. The number of iterations to reach machine precision can also be decreased by choosing initial velocity profiles closer to the ones leading to optimal growth.
The optimal gain also converges as the time step decreases following a power law as illustrated in Figure 2. The time step chosen in this study is, unless otherwise specified, ∆t = 10 −3 , to obtain a good compromise between accuracy and computational cost. As we observe asymptotic behaviours for extreme values of S, these are less relevant, we therefore use in this work values of S ranging form 10 −4 to 10 −1 as it is the region where interesting behaviour is observed. We chose to keep f constant at f = 0.1, as f was not shown to significantly impacts the results found, similarly to what has been observed in the case of the linear stability analysis with the same model (11). Reynolds

Growth envelope
The value of the maximum transient growth depends on the target time chosen. While we are mostly interested in optimising for T , it is still interesting to see how G depends on T . Figure    the particle distribution is homogeneous or not. This result, plus the similar shapes of the envelopes, suggests that the mechanisms producing growth are the same for single phase and particulate flows. The mode (α, m) = (1, 1) is the most affected by the addition of particles, especially for nonhomogeneous particle distributions as illustrated in figure 3. However, (α, m) = (1, 0) is still the mode for which the gain is the strongest. From now on, G is optimised over T when studying the optimal gain and the value of the mass fraction f is also kept constant, to f = 0.1.

Homogeneous particle distribution
We first examine the effect of adding homogeneously distributed particles in the flow. In order to illustrate the effect of particles, we define the ratio between the growth for the particulate flow with a given set of parameters and the single phase flow with the same Reynolds number: where G p and G f , are the optimal gains for particulate and single phase flows respectively, both maximised over all values of T . A similar ratio is chosen between the times of optimal growth for particulate and single phase flow, where T p and T f are the target time associated to G p and G f .

Effect of the Stokes numbers on the gain
We first consider variations of these quantities as a with the Stokes number. Figure 4 shows G ′ and T ′ as a function of S, for different values of Re.
The addition of particles increases the optimal gain for all values of S. The curves are non monotonic, with G ′ increasing until it reaches a peak defined as G ′ peak , for an associated Stokes number S G . For large values of S, the ratio G ′ seems to decreases towards 1. In the limit of S → ∞, the particles are so heavy that they are effectively decoupled from the flow and have no effect on it. When S → 0, G ′ ≈ 1.21. The difference between single phase and particulate flows is due to the modification of the average density of the flow caused by the particles. With f = 0.1, the modified Reynolds number is Re ′ = (1 + f )Re = 1.1Re. Since for the single phase pipe flow, G f ∝ Re 2 , G p ∝ (1+f ) 2 Re 2 , as observed in Figure  4. It follows that there is an optimal Stokes number S for which the influence of homogeneously distributed particles on the optimal gain is greatest. A similar behaviour is observed for T ′ . For all values of Re and S considered, the growth is delayed for particulate flows compared to single phase flows. As S → 0, the time for which the growth is maximised increases by 10% compared to the single phase flow. This corresponds to the time of optimal growth for the modified Reynolds number Re ′ = (1 + f )Re since, as discussed in the previous section, the time for maximum growth increases linearly with Re. Similarly, a peak for the time ratio T ′ peak occurs at a Stokes number S T . The time ratio then decreases as the Stokes number continues to increase in a similar fashion as the ratio of growths. For all Reynolds numbers considered and f = 0.1, the peak Stokes number is around S T = 5 × 10 −2 .

Effect of the Reynolds number
The Reynolds number has little incidence on G ′ , as the envelopes have a very similar shape when Re is varied. For all Re considered, the curves of Figure 4 exhibit a peak at approximately the same Stokes number, S G ≃ 2.5 × 10 −2 . The value G ′ peak of the peak shows little change, varying by only 0.25% for Re ranging from 500 to 5000. Since G ′ peak is almost constant over the Reynolds number and the optimal gain for single phase flow G f (Re) scales with Re 2 (29), it follows that the optimal gain for particulate flows optimised over S also scales with Re 2 . Similarly, the ratio of delays T ′ peak varies little with the Reynolds number, with variations just under 1% for Re ranging from 500 to 5000. Moreover, T f scales linearly with the Reynolds number. Therefore T p optimised over S scales linearly with the Reynolds number as well.

Inhomogeneous particle distribution
Particles do tend to cluster in laminar pipe flows (6) such that considering homogeneously distributed particles is less realistic. Moreover, allowing for an inhomogeneous distribution dramatically increases the effect of the solid phase in the case of the linear stability (11). Figure 5 shows G ′ and T ′ for varying values of σ, centred around r d = 0.7 (figure 5a) and r d = 0.3 (figure 5b). The overall shape of the growth is the same as in the case of a homogeneous particle distribution, but the effect of the solid phase on the gain is significantly stronger. It also significantly varies with σ. The more the particles are concentrated, i.e. the smaller σ, the larger both G ′ and T ′ are, as illustrated in figure 5. Varying σ from 0.15 to 0.10, G ′ peak increases by 24% and T ′ peak by 18% for r d = 0.7 (table 3). The effect of σ is similar for r d = 0.3, as seen in figure 5b. However, the values of G ′ peak and T ′ peak are noticeably smaller for equivalent values of σ. This indicates that the position of the particles determine the amount of transient growth as well. However, both G ′ and T ′ still tends towards 1 as S → ∞ in all cases. Compared to homogeneous particles distribution, the values of S yielding the maximal gain shifts to a larger value for r d = 0.7. The time at which the optimal growth occurs is delayed for Gaussian particle distributions as well. On the other hand these effects are reversed for r d = 0.3. Moreover, while changing σ affected the growth ratio, it has little effect on the value of S G and S T for all cases observed.   Figure 6 shows the ratios of gains as a function of S for several average radii r d of the particle distribution, Re = 1000 and σ = 0.1. G ′ peak and T ′ peak show strong variations with r d . Indeed, the G ′ peak ranges from 1.95 to 3.50 while T ′ peak ranges from 1.40 to 2.20. The effect of particles on the flow is highest for r d in the range 0.5 − 0.6 both for the ratio of maximum growths and the ratio of optimal times as shown in Figure  6. This value is relatively close to the Segré-Silberberg radius where particles are known to naturally cluster (5), albeit a little closer to the pipe centre. Unlike the Reynolds number, r d has a very strong influence on the optimal Stokes number. Both S G and S T are larger than their counterpart in the case of a uniform particle distribution for all S, with the exception of r d = 0.3, which where this is only the case for S 0.07.

Topology of the optimal perturbations
In this section we study the topology of the optimal velocity fields. First, we consider the velocity fields at t = 0, subsequently called optimal perturbation and denoted u 0 for the fluid and u p0 for the particles' velocity. Second, the velocity fields at t = T , referred to as the velocity peak and denoted u T for the fluid and u pT for the particles' velocity. Two target times are shown here, T = 14 for which the mode (α, m) = (1, 1) is dominant, and T = 90 for which (α, m) = (0, 1) is dominant. The radial sections of the optimal perturbations are very different depending on whether the dominant mode is (α, m) = (1, 1) or (α, m) = (0, 1). Figures 7a-7f shows the contours of streamwise velocity and sectional fluid velocity vectors, for a single phase flow and a particulate flow with a nonhomogeneous distribution. The optimal perturbations contours are weakly affected by the addition of particles in this case. For T = 14 (figures 7a and 7c), the profile of the optimal perturbations shows two symmetric rolls in the spanwise direction. The streamwise velocity has a peak in the shape of an antisymmetric annulus between r = 0.5 and r = 0.7. Streamwise and spanwise velocities are of the same order of magnitude in both cases. For a larger target time T , the streamwise independent mode is dominating as illustrated in figures 7b and 7d. In the spanwise direction we observe two rolls that are distinctive of the usual single phase transient growth. Figures 7e and 7f shows the contours of streamwise velocity and sectional particles velocity vectors: both the streamwise and sectional particles velocities are strongest in the region where the particle concentration is the highest. Examples of peak velocity contours are given in figure 8. In all cases considered the streamwise velocity is larger than the spanwise velocity, for both u T and u pT . The effect is even more pronounced for T = 90. u T (= 14) spanwise component is similar to u 0 (T = 14)'s but the streamwise component is different with two opposed currents, one close to the centre and the other around it. In the case T = 90 the streamwise velocity take the form of two antisymmetric rolls. Figure 8 shows that fluid and particles profile are, in the case of the peak velocity, almost identical, due to the strong coupling between the fluid and solid phases. In summary, the addition of particles does not alter the transient growth mechanisms nor the topology of the optimal modes, even if the particles are distributed inhomohenously. Particles tend to be accelerated where the flow is, however their effect on the growth itself is significant when they are inhomogenously distributed.

Conclusion and discussion
This paper presented a study of the linear transient growth of particulate pipe flow through a simple two-fluid, two-way model for the solid and the liquid phases. The addition of particles has been found to increase the amount of transient growth regardless of Stokes number. However the modes that are responsible for the transient growth remain the same as those in flows without particles. The growth itself still varies with the Stokes number, with a sweet spot for which it is maximised. Interestingly, the corresponding Stokes number is independent of the Reynolds number. Moreover, the ratio of growths optimised over time for the particulate to non-particulate flow, G ′ peak , is also independent of the Reynolds number, implying that the growth for the particulate flows scales as Re 2 as it does for the single phase flow (13). We also showed that the solid phase has a delaying effect on the transient growth, again regardless of the Stokes number considered. We observe that there is a value of S for which the delay is maximised, this Stokes number is independent of the Reynolds number as well. Here too, the ratio of times of optimal growths for the particulate to the non-particulate flow T ′ peak , is independent of the Reynolds number. This implies that the time for which growth is optimised scales as Re as it does for the single phase flow. The most important result is that allowing for particles to be inhomogeneously distributed can drastically increase their impact on the transient growth, which is increased by more than 200% depending on their size and the shape of the spatial distribution. The way in which the particles are distributed is important too. We have considered particles in a Gaussian distribution of standard deviation σ located in annulus located at radius r d . The transient growth increases monotonically as σ decreases, i.e. when the particles are more localised. The effect of the solid phase on the transient growth was also found to be weaker when the particles are localised close to the wall (r d close to 1) or at the pipe centre (r d close to 0), and strongest in the intermediate region (r d = 0.6 − 0.7). This region seems to play a particularly important role both in the laminar state and in the process of transition to turbulence: not only do particle tend to naturally cluster there in the laminar state (5; 25), but particulate pipe flows have been found linearly unstable when particles of intermediate size are added in that region (11). This raises the question of whether the actual pathway to turbulence is indeed sensitive to particles being present in that region, a question that could be answered in further analysis including fully nonlinear effects. Secondly, the question of the robustness of the model is also crucial: further studies with more physically accurate model for the solid phase and their interaction with the fluid phase are needed to confirm the nature of the role played by the mechanisms identified in this work.