Turbulence statistics in a negatively buoyant multiphase plume

TL;DR

This study examines turbulence in a negatively buoyant multiphase particle plume, revealing similarities and key differences with bubble plumes, including turbulence structure, energy budget, and spectral characteristics, advancing understanding of multiphase flow dynamics.

Contribution

It provides detailed turbulence statistics of a heavy particle plume using advanced imaging, highlighting structural differences from bubble plumes and extending classical plume theories to multiphase flows.

Findings

Turbulence exhibits self-similarity in the plume centerline.

The particle plume shows a non-reversed skewness compared to bubble plumes.

TKE budget includes significant shear production, unlike bubble plumes.

Abstract

We investigate the turbulence statistics in a {multiphase plume made of heavy particles (particle Reynolds number at terminal velocity is 450)}. Using refractive-index-matched stereoscopic particle image velocimetry, we measure the locations of particles {whose buoyancy drives the formation of a multiphase plume,} {together with the local velocity of the induced flow in the ambient salt-water}. {Measurements in the plume centerplane exhibit self-similarity in mean flow characteristics consistent with classic integral plume theories.} The turbulence characteristics resemble those measured in a bubble plume, {including strong anisotropy in the normal Reynolds stresses. However, we observe structural differences between the two multiphase plumes. First, the skewness of the probability density function (PDF) of the axial velocity fluctuations is not that which would be predicted by simply reversing the direction of a bubble plume. Second, in contrast to a bubble plume, the particle plume has a non-negligible fluid-shear production term in the turbulent kinetic energy (TKE) budget. Third, the radial decay of all measured terms in the TKE budget is slower than those in a bubble plume.} Despite these dissimilarities, a bigger picture emerges that applies to both flows. The TKE production by particles (or bubbles) roughly balances the viscous dissipation, except near the plume centerline. The one-dimensional power-spectra of the velocity fluctuations show a -3 power-law that puts both the particle and bubble plume in a category different from single-phase shear-flow turbulence.

Figures

Tables

Table 1. Table 1 : Physical properties of particles and surrounding salt-water solution.

Parameter	Particle	Fluid
	(Teflon PFA)	(0.4% salt water solution)
Length (mm)	2.0	–
Bulk volume (mL)	110	2.28 $\times {10}^6$
Density ($g/{\mathrm{cm}}^3$)	2.1	1.04
Dynamic viscosity (kgm-1s-1)	–	1.059$\times {10}^{-3}$
Refractive index	1.34	1.338

Table 2. Table 2 : The second moment of the six available components of the velocity gradient tensor averaged axially and radially. The quantities within the bracket indicate the 95% uncertainty bounds obtained with the bootstrap method

$i$
1	1555 [1536, 1577]	1509 [1487, 1530]
2	624 [619, 629]	565 [561, 568.5]
3	399 [396, 402]	368.5 [366, 371]

Table 3. Table 3 : Simplified similarity expressions and values of various integrals used in equation 6 . The shape functions F 4 ​ ( η ) subscript 𝐹 4 𝜂 F_{4}(\eta) and F 5 ​ ( η ) subscript 𝐹 5 𝜂 F_{5}(\eta) are summarized in Appendix A .

Term in eq. 5	Full form	Similarity form	Value
Flux of $K$, $I_1$		$2\pi {\int }_0^{\infty }\eta e^{-3{\eta }^2}𝑑\eta$	1.047
Flux of $k$, $I_2$		$2\pi {\int }_0^{\infty }\eta e^{-{\eta }^2}{F}_4\left(\eta \right)𝑑\eta$	0.095 – 0.105
Dissipation, $I_3$		$2\pi {\int }_0^{\infty }\eta F_5\left(\eta \right)𝑑\eta$	0.021 – 0.027
$K\approx U_1^2/2$	$k=\left(u_1^2+u_2^2+u_3^2\right)/2$

Table 4. Table 4 : Qualitative comparison summary of three TKE budget terms in our particle plume and those in a bubble plume (Lai & Socolofsky, 2019 ) , a single-phase jet (Darisse et al. , 2015 ; Lai & Socolofsky, 2018 ) , and a variable density single-phase (VDSP) jet (Charonko & Prestridge, 2017 ) .

TKE budget term	Multiphase		Single phase
	Particle	Bubble	Jet	VD jet
Advection, $A$	\textcolorgreen✓	\textcolorgreen✓	\textcolorgreen✓	\textcolorgreen✓
Shear production, $P_s$	\textcolorgreen✓	\textcolorred✗	\textcolorgreen✓	\textcolorred✗
Turbulent transport, $T$	\textcolorgreen✓	\textcolorgreen✓	\textcolorred✗	\textcolorgreen✓
\textcolorgreen✓: similar profile shapes	\textcolorred✗: dissimilar profile shapes

Table 5. Table 5 : Coefficients of double-Gaussian function in equation 9 used to fit various parameters measured in the experiment

Function ($F_i\left(\eta \right)$)	Parameter	$f_1$	$f_2$	$f_3$
$F_1$		0.069	0.74	0.63
$F_2$		0.011	0.57	1.00
$F_3$		0.009	0.62	0.84
$F_4$	$k/U_c^2$	0.045	0.71	0.68
$F_5$		0.004	0.76	0.72

Equations

$${\theta }_1=atan\left(\frac{{⟨U_2^{\prime }⟩}_0}{{⟨U_1^{\prime }⟩}_0}\right);{\theta }_2=atan\left(\frac{{⟨U_3^{\prime }⟩}_0}{{⟨U_1^{\prime }⟩}_0}\right).$$

$$\frac{⟨U_2⟩}{U_c}=-\frac{\alpha }{\eta }\left((1-e^{-{\eta }^2})-\frac{\beta }{\alpha }{\eta }^2{e}^{-{\eta }^2}\right).$$

$$\text{centering}M(x_1)=(3\sqrt{2\pi }\beta F_0/4{)}^{2/3}(x_1-x_0{)}^{4/3}.\text{@add@centering}$$

$$⟨ϵ⟩=\nu$$

$$+⟨{\left(\frac{\partial u_2}{\partial x_1}\right)}^2⟩+⟨{\left(\frac{\partial u_2}{\partial x_2}\right)}^2⟩+⟨{\left(\frac{\partial u_2}{\partial x_3}\right)}^2⟩$$

$$\displaystyle+\left<\left(\frac{\partial u_{3}}{\partial x_{1}}\right)^{2}\right>+\left<\left(\frac{\partial u_{3}}{\partial x_{2}}\right)^{2}\right>+\left<\left(\frac{\partial u_{3}}{\partial x_{3}}\right)^{2}\right>\Biggr{]}.$$

$$\text{centering}{\iiint }_V\frac{D}{Dt}\left(⟨K⟩+⟨k⟩\right)dV=-f_c{\iiint }_V⟨ϵ⟩dV+{\iiint }_V\frac{\Delta \rho }{\rho }⟨U_1⟩gdV,\text{@add@centering}$$

$$\text{centering}{\left[\left(0.5I_1+I_2\right)\left(b_g^2{U}_c^3\right)\right]}_{x_1}^{x_2}=f_c{\int }_{z_1}^{z_2}{I}_3{b}_g{U}_c^{3}dx+1.44F_0\Delta x_1\text{@add@centering}$$

$$0$$

$$\stackrel{\text{T}}{\overbrace{-\left(\frac{1}{r}\frac{\partial ⟨u_{r}k⟩}{\partial r}+\frac{\partial ⟨u_{z}k⟩}{\partial z}\right)}}$$

$$\stackrel{\text{Tp}}{\overbrace{-\frac{1}{\rho }\left(\frac{1}{r}\frac{\partial ⟨u_{r}p⟩}{\partial r}+\frac{\partial ⟨u_{z}p⟩}{\partial z}\right)}}$$

$$\stackrel{\text{Ps}}{\overbrace{-\left(⟨u_z^2⟩\frac{\partial ⟨U_z⟩}{\partial z}+⟨u_r^2⟩\frac{\partial ⟨U_r⟩}{\partial r}+\frac{⟨U_r⟩⟨u_{\theta }^2⟩}{r}+⟨u_z{u}_r⟩\left(\frac{\partial ⟨U_z⟩}{\partial r}+\frac{\partial ⟨U_r⟩}{\partial z}\right)\right)}}$$

$$-⟨ϵ⟩$$

$$+P_b$$

$$\text{centering}⟨{\left(\frac{\partial u_1}{\partial x_1}\right)}^2⟩\approx \frac{2⟨u_1^2⟩}{{\lambda }_f^2}.\text{@add@centering}$$

$$\text{centering}{F}_i(\eta )=f_1\left[\exp \left(-{\left(\frac{\eta -f_3}{f_2}\right)}^2\right)+\exp \left(-{\left(\frac{\eta +f_3}{f_2}\right)}^2\right)\right]\text{@add@centering}$$

$$\text{centering}S=⟨u_z{u}_r⟩/U_c^2=\left(0.0124\eta +0.0286{\eta }^3-0.0113{\eta }^5\right)\exp \left(-1.47{\eta }^2\right)\text{@add@centering}$$

$$\text{centering}{T}_1=⟨u_{z}k⟩/U_c^3=\left(-0.0016-0.0252{\eta }^2-0.0181{\eta }^4+0.0746{\eta }^6\right)\exp \left(-2.62{\eta }^2\right)\text{@add@centering}$$

$$\text{centering}{T}_2=⟨u_{r}k⟩/U_c^3=\left(-0.0008\eta 0.0048{\eta }^3-0.0139{\eta }^5-0.0049{\eta }^7\right)\exp \left(-1.76{\eta }^2\right)\text{@add@centering}$$

$$\text{centering}A=-[⟨U_z⟩\frac{\partial k}{\partial z}+⟨U_r⟩\frac{\partial k}{\partial r}]\text{@add@centering}$$

$$\text{centering}A\times \frac{b_g}{U_c^3}=\left[\frac{2}{3}{F}_4(\eta )+\eta \frac{d{F}_4(\eta )}{d\eta }\right]\beta e^{-{\eta }^2}-\left[G(\eta )\frac{d{F}_4(\eta )}{d\eta }\right]\text{@add@centering}$$

$$\text{centering}P=-\left[\underset{P_1}{\underbrace{⟨u_z^2⟩\frac{\partial ⟨U_z⟩}{\partial z}}}+\underset{P_2}{\underbrace{⟨u_r^2⟩\frac{\partial ⟨U_r⟩}{\partial r}}}+\underset{P_3}{\underbrace{\frac{⟨U_r⟩⟨u_{\theta }^2⟩}{r}}}+\underset{P_4}{\underbrace{⟨u_z{u}_r⟩\left(\frac{\partial ⟨U_z⟩}{\partial r}+\frac{\partial ⟨U_r⟩}{\partial z}\right)}}\right]\text{@add@centering}$$

$$P_1\times \frac{b_g}{U_c^3}$$

$$P_2\times \frac{b_g}{U_c^3}$$

$$P_3\times \frac{b_g}{U_c^3}$$

$$P_4\times \frac{b_g}{U_c^3}$$

$$\text{centering}T=-\left[\frac{\partial ⟨u_{z}k⟩}{\partial z}+\frac{1}{r}\frac{\partial ⟨u_{r}k⟩}{\partial r}\right]\text{@add@centering}$$

$$\text{centering}T\times \frac{b_g}{U_c^3}=\beta \left(T_1(\eta )+\eta \frac{dT1}{d\eta }\right)-\left(\frac{T_2(\eta )}{\eta }+\frac{\partial T_2(\eta )}{\partial \eta }\right)\text{@add@centering}$$

Full text

1 Introduction

Plumes containing bubbles, particles and droplets are present in both environmental and industrial applications. A few examples of environmental interest are settling sediment, volcanic eruption columns, CO2 ocean sequestration plumes, and rising oil droplets and gas bubbles from oil well blowouts (Freeth, 1987; Baines & Sparks, 2005; Baines, 2008; Socolofsky & Bhaumik, 2008; Woods, 2010; Socolofsky et al., 2011; Huppert & Neufeld, 2014; Wang et al., 2016). Suspension-flow plumes differ from traditional single-fluid plumes in that the energy due to buoyant forcing is transmitted indirectly from the suspended phase to the continuous phase. The relative motion between the two phases introduces additional length and time scales; these must be included in the model formulation employed to predict plume behavior, for example, when extending classic single-phase integral plume models (Morton et al., 1956) to multiphase plumes (Milgram, 1983; Sun & Faeth, 1986). Identifying such scales is nontrivial due to the complexity of particle-particle and particle-fluid interactions.

For the special case of air bubbles in water, empirical data collection has allowed accurate closure of predictive schemes (Lance, 1991; Risso, 2018; Risso & Ellingsen, 2002; Mercado et al., 2010; Almeras et al., 2017). However, the bubble-in-water plume can be quite different from other suspension-flow plumes of interest, because bubbles are far less dense than the ambient fluid (specific gravity of order 10*-3*) and have negligible inertia. Other plumes of interest are droplets in water (specific gravity of order 1), solid particles in water (specific gravity of order 2 to 10), or liquid droplets in air (specific gravity of order 103). Each of these different suspensions can behave quite differently than the others in terms of the interphase interactions. Empirical data are relatively scarce for these other suspension-flow plumes.

Recent progress in suspension flows (especially turbulent flows) offers the hope that predictive techniques will eventually describe overall plume behavior from a direct description of the internal dynamics of particle-fluid coupling. To help support this strategy and address some of the open questions that remain about turbulent suspension flows (Balachandar & Eaton, 2010; Guazzelli & Morris, 2012), we take an observational approach herein, measuring the particle and fluid behaviors in a particle-laden plume. Because particles are dynamically quite different from air bubbles in water, we are curious to see how the two-phase kinematics and the interstitial fluid turbulence behave.

Both bubbles and particles have been shown to modulate the turbulent properties of multiphase flows in a manner sensitive to volume fraction ($\phi$), particle Reynolds number ($Re_{d}$), and particle Stokes number ($St$). One starting point for understanding such turbulence modulation is the agitation of a quiescent fluid by bubbles or particles distributed homogeneously in space. This is often referred to as pseudo-turbulence (e.g. Lance (1991); Cartellier & Rivière (2001)) as its characteristic scales are set by the bubble or particle wakes. Pseudo-turbulence is fundamentally anisotropic, with stronger velocity fluctuations aligned in the direction of suspended-phase motion ($x_{1}$, typically vertical). Numerical simulations of pseudo-turbulence by Riboux et al. (2013) show that the velocity fluctuation energy spectrum has a slope of $k^{-3}$ for wavelengths smaller than the integral length scale; they also find that the integral length scale is the ratio of the bubble or particle diameter and the drag coefficient, Cd.

A recent topical review article by Risso (2018) provides a comprehensive summary of pseudo-turbulence due to a homogeneous swarm of bubbles rising in a quiescent medium. First, the turbulence intensity is anisotropic and is dominated by the vertical component. The probability density function (PDF) of the vertical component is positively skewed due to the wake immediately behind each rising bubble, whereas the other two components are symmetrically distributed with exponential tails. Because the turbulent kinetic energy (TKE) production in suspension flows occurs at two length scales (first at the scale of bubble wakes and second at the scale of overall bubble population), the combination of the two results in a $k^{-3}$ spectral subrange.

To our knowledge, there are no existing studies examining turbulence in negatively buoyant particle plumes released into an unstratified, initially quiescent fluid. However, there are many studies of turbulence in bubble plumes with notable contributions coming from Soga & Rehmann (2004); Wain & Rehmann (2005); García & García (2006); Seol et al. (2009); Simiano et al. (2009); Lai & Socolofsky (2019). These investigations reveal different velocity characteristics at different axial distance ($x_{1}$) from the origin. Simiano et al. (2009) examined the near-field characteristics ($x_{1}/D<2$, where $D$ is the dynamic length scale to be defined in section 3.1) and showed that the meandering nature of bubbles influences the effective spreading, volume fraction, and mean velocity profiles. Seol et al. (2009) and Lai & Socolofsky (2019) focused on the far field ($2<x_{1}/D<11$) turbulence characteristics, such as Reynolds-stress, turbulent transport and TKE budget; their results will be used to compare with our particle plume data.

Characterizing the turbulence statistics in a multiphase plume is challenging primarily due to the difficulty in simultaneously measuring velocity in both phases. Traditional intrusive techniques such as hotwire anemometry suffer from the possibility of being damaged by the solid particles. Optical measurements from techniques such as particle image velocimetry (PIV) and laser doppler anemometry are usually distorted due to the difference in optical properties of the two media. We overcome the distortion herein by carefully choosing two media with matched refractive indices.

This paper is organized as follows: section 2 provides a detailed description of the plume-generating facility and method for characterizing the two-phase flow. We report and discuss our experimental observations in section 3. Finally, section 4 provides the key conclusions from this work.

2 Experimental setup and methodology

2.1 Plume facility

The plume experiments were conducted in a rectangular tank (80 cm deep $\times$ 80 cm wide $\times$ 365 cm long) as described in Bordoloi & Variano (2017). The tank was filled with tap water filtered to 5 microns and maintained using an UV purifier. 91.1 kg of commercial sodium chloride (Cargill Top-Flo) was then mixed to yield a salt concentration of 0.04 g/mL. The resulting solution has density $\rho_{s}$ = 1.04 $\mathrm{g/cm^{3}}$ and kinematic viscosity $\mu_{s}$ = 1.059$\times 10^{-3}$ Pa.s (see Table 1).

A schematic of the experimental setup is shown in figures 1(a). The negatively buoyant plume was created by releasing 110 mL of cylindrical PFA pellets (Neoflon AP-202) from a height of 56.5 cm via a screw-feed particle release mechanism (see figure 1(b)). The pellets are right circular cylinders with length=diameter=2 mm (see figure 1(b)-A). The physical properties of the solid particles and the surrounding salt solution are summarized in Table 1. Because of the hydrophobic nature of PFA, the particles tend to trap and hold air bubbles when added to water. To prevent these air bubbles from entering the experiment, particles were presoaked in water in a separate container and rapidly stirred to dislodge all bubbles. Once the particles were free of air bubbles, they were placed in a funnel for eventual release into the quiescent salt-water mixture. The funnel was kept partially submerged 23 cm below the free surface through a nozzle with internal diameter $d_{0}$ = 11.25 mm so that the particles did not contact air (see figure 1(b)). Particle release was governed by a motor-driven helical screw. Prior to release, the particles were held inside the funnel by the blades of the screw. Upon release, the motor rotated the screw at a constant rate of 0.5 RPM operated via Lego Mindstorms software. The particle flux, $Q_{0}$ = 7.5 cm3/s was measured by video-recording the release of a known quantity of particles and measuring the time difference between the exit of the first and last particle from the funnel nozzle. Before each experiment, the tank fluid was seeded with optical tracers, specifically 13 $\mu$m silver-coated hollow glass spheres (SH400S20, Conduct-O-Fil, Potters Industries). Twenty plume releases provided a sufficient number of samples for the analysis described in Section 3.

2.2 Refractive index matching

To measure the fluid velocity inside and around the plume, we matched the refractive index (n) of the particle and fluid phases. A target refractive index of PFA (n $\approx$ 1.34) was achieved from an aqueous solution of commercial sodium chloride with a concentration of 0.04 g/mL. The refractive index of the solution was measured to be 1.338 using an ATAGO refractometer and found to match with the salinity vs. n prediction given in Tan & Huang (2015). Although the refractive indices of the two phases were matched, because of the scattering properties of PFA the particles appeared bright under the laser illumination. To limit the intensity of particles below the saturation threshold of the camera’s sensor we used a circular polarizer on each camera lens. Figure 2 shows example images of Teflon particles in different-salinity water samples illuminated by a laser sheet (wavelength=532 nm). The high-intensity sections in the images are the particles intersected by the laser sheet, and the low-intensity sections are the out-of-plane particles. Fluid tracers are much more visible when particles and fluid have matched refractive indices (figure 2(b) compared to figure 2(a)), especially behind or in front of out-of-plane particles as shown in the highlighted regions.

2.3 Multiphase velocimetry

We performed stereoscopic particle image velocimetry (SPIV) to find the two-dimensional three-component (2D3C) velocity field of the fluid phase of the plume. The origin of our coordinate system is located at the tank center. The plume is axisymmetric about the vertical axis, but we use Cartesian coordinates ($x_{1}$, $x_{2}$, $x_{3}$) for our measurements as described in figure 1(b). Velocity measurement of the particle phase is currently underway and beyond the scope of this paper.

The $x_{3}$ = 0 plane was illuminated with a laser sheet (1 mm thick at beam waist; Quantel/Big Sky Lasers, 532-nm double-pulse Nd-YAG). Two charge-coupled device (CCD) cameras (Imager PRO-X, 1600 pixels $\times$ 1200 pixels) were synchronized with the laser pulses. They were placed $\pm 55^{\circ}$ to the laser sheet (c.f. 90∘ for standard 2D2C PIV). To minimize distortion due to this angle, water-filled acrylic prisms were placed between the camera lenses and the tank walls. The cameras were mounted with Nikkor 105 mm lenses, circular polarizers, and Scheimpflug adapters. The interframe delay ($\Delta$t) was optimized to 0.5 ms. The PIV images were collected at a frequency of 14.0 Hz. Both fluid and particle motions were decorrelated at timescales $>$ 1/14 s. This setup gave approximately 50 independent samples during the steady-state phase of each experimental run.

2.3.1 Fluid-phase processing

We computed the fluid-phase velocity field using DaVis 8.2 software (LaVision GmbH). Before stereo-PIV processing, we performed stereo self-calibration to reduce disparity in the alignment between the laser sheet and the measurement plane to below 0.1 pixel. Tracers and particles in an image were separated by intensity thresholding; erosion and dilation were applied to isolate individual particles. Fluid velocity fields were obtained by masking particles and correlating tracer locations through multipass stereoscopic cross-correlation. The cross-correlation was applied with an initial interrogation window of 64 $\times$ 64 pixels and final window of 48 $\times$ 48 pixels with 50% overlap, yielding a spacing between vectors of 0.7 mm. Interrogation windows were weighted according to the symmetric two-dimensional Gaussian function. Vectors were discarded based on universal outlier detection (Westerweel & Scarano, 2005) and left as data gaps. Figure 3 shows a representative field of fluid velocity fluctuations superposed with the corresponding particle image (transformed into the lab coordinate system). The grid cells are color-coded with the out-of-plane velocity component while the vectors show the in-plane velocity field.

2.3.2 Particle-phase processing

To estimate the particle number density across the plume, we conducted a series of additional image-processing steps (see figure 4). First, the raw images were transformed from the camera coordinate system to the lab coordinate system using the stereoscopic mapping function. Small-size background noise features (including tracers) were removed using a simple median filter.

In stereo-mapping, after transforming into the lab coordinate system, the section of the particle intersected by the laser sheet should overlap in both cameras. Figure 4(a) shows a sample instant with the two camera images overlayed. The high intensity regions in the image represent the intersection between the two camera images, whereas the low intensity background is from the out-of-plane particles captured in only one of the two cameras. We remove the background noise by setting an intensity threshold and convolving the two images. The result is shown in figure 4(b).

After inspecting our dataset, we find many instances where particles are nearly touching as exemplified in the lower left pair in the image. The nonuniform intensity gradients in the bordering regions among different particle pairs yielded different intensity gradients. A traditional intensity-gradient-based image segmentation therefore could not differentiate between two particles, and so adopted a segmentation technique based on the watershed transform. The watershed transform treats an image as a topographic map with brighter intensity pixels as heights, and finds the separating line that runs along the ridges (Gonzalez & Woods, 2007). The result, with individual particles identified with separate colors, is shown in figure 4(c). We applied an area threshold to remove partially illuminated particles and clusters that our method could not isolate. The centroids of particles identified from the final binary image are shown in figure 4(d).

2.4 Coordinate transformation

Typically in jets and plumes the cross-stream velocity components ($U_{2}$ and $U_{3}$) are smaller than the axial velocity component ($U_{1}$) near the centerline. Thus, any small misalignment between the plume axis and the measurement axis could lead to a systematic bias in the measured $U_{2}$ and $U_{3}$. The scenario is schematically shown in inset B in figure 1(b) in a simplified $x_{1}-x_{2}$ plane in which the plume axis ($x_{1}$) makes an angle $\theta_{1}$ with the measurement axis $x_{1}^{\prime}$. For our data, if $\theta_{1}$ is defined the same way and if $\theta_{2}$ is the angle between the plume axis and the measurement axis on the $x_{1}-x_{3}$ plane, using axisymmetry and assuming that $\left<U_{2}\right>,\left<U_{3}\right>=0$ at the plume centerline, the measured velocities can be transformed into the plume co-ordinate system via correction angles,

[TABLE]

The measurement bias in $U_{2}$ and $U_{3}$ at the plume centerline are captured in figure 5(a) that shows the respective mean radial profiles, $\left<U_{2}^{\prime}\right>$ and $\left<U_{3}^{\prime}\right>$. We perform this transformation for the 20 replicate datasets independently with correction angles ($\theta_{1}$, $\theta_{2}$) reported in figure 5(c). The bias in the non-zero centerline velocities are reflected in $\theta_{1}$ and $\theta_{2}$ which are small and less than 1.5∘. The two velocity components transformed into the plume coordinate systems are shown in figure 5(b). The analysis that follows will be in the plume coordinate system.

3 Results and discussion

3.1 Mean flow characteristics

We first characterize the mean interstitial fluid velocity based on approximately 1000 independent PIV snapshots. The mean axial velocity field $\left<U_{1}(x_{1},x_{2})\right>$ in dimensional form is shown in figure 6(a). The radial ($x_{2}$) variations of $\left<U_{1}\right>$ normalized by the local mean centerline velocity $\left<U_{c}(x_{1},0)\right>$ (written as $U_{c}$ for simplicity from here onward) are shown in figure 6(b) at representative axial locations ($12.9d_{0}<x_{1}<15.1d_{0}$) across our measurement region. The Gaussian plume radius, $b_{g}$, defined as the $x_{2}$ location where $\left<U_{1}\right>=e^{-1}U_{c}$, is obtained via a nonlinear least-squares fit of each measured profile to a normalized Gaussian function. The data closely follow the Gaussian curve (solid line in figure 6(b)), which is also seen in single- and multi-phase jets and plumes (Darisse et al., 2012; Milgram, 1983; Lai & Socolofsky, 2019).

Next, we examine the axial decay of centerline velocity ($U_{c}$) and the axial growth of Gaussian plume radius ($b_{g}$) for a particle plume by comparing them to existing data from bubble plumes (see figure 7). For the purpose of comparison, we use two nondimensional parameters: a velocity scale ($u_{s}$) and a length scale ($D=gQ_{0}/4\pi\alpha_{0}^{2}u_{s}^{3}$). These are based on an integral model of bubble plumes (Bombardelli et al., 2007), in which a constant $\alpha_{0}=0.083$ is used as the entrainment coefficient. The computed entrainment coefficient of our particle plume to be presented later is different from this a priori constant. The velocity scale ($u_{s}$ = 0.23 m/s) is the terminal velocity of a single particle in quiescent fluid and is computed using a simple drag model (Clift et al., 1978). Based on $Q_{0}$ = 7.5 $\mathrm{cm^{3}/s}$, the length scale $D$ is approximately 6.8 cm. This situates our axial measurement region between 2-2.5$D$.

The overall trends in both centerline velocity and plume radius show similarity between particle and bubble plumes. The axial variation of $U_{c}/u_{s}$ in a particle plume is not measurable from our data, but we can say that it falls above the curve for bubble plumes (Lai & Socolofsky, 2019) (see figure 7(a)). The dashed and solid lines in figure 7(a) are the $\mathrm{-1/3^{rd}}$ power-law fits, with $A$=2.05 (current data, see inset) and $A$=1.6 (Lai & Socolofsky (2019)). The variation of $b_{g}$ for all the data (including ours) is captured by the solid line ($b_{g}/D=0.114x_{1}/D$) in figure 7(b). Within the limited axial extent ($\approx 0.5D$) of our measurement, the local spreading rate ($\beta=db_{g}/dx_{1}$) is difficult to estimate. However, justified based on the collapse of our data on the solid line in figure 7(b), we will use $\beta=0.114$ in the remaining analyses in this paper.

Profiles of mean radial ($\left<U_{2}\right>$) and out-of-plane ($\left<U_{3}\right>$) velocities at various axial locations are shown in figures 8(a) and 8(b), respectively. The mean radial velocity exhibits typical jet/plume behavior: within the plume radius ($|x_{2}/b_{g}|<1$) the interstitial fluid flows away from the centerline, whereas outside the plume radius ($|x_{2}/b_{g}|>1$) the surrounding fluid is entrained into the plume. The mean out-of-plane velocity should be zero by design (swirling motions are not introduced by the particle feeder) but it shows some variations and asymmetry far from the plume axis. The solid line in figure 8(a) is a nonlinear least-squares fit to a shape function,

[TABLE]

The shape function in equation 2 is obtained from the radial integration of the continuity equation written in a cylindrical coordinate system after adopting the entrainment hypothesis for jets/plumes and the Gaussian-profile assumption for $\left<U_{1}\right>$ (Papanicolaou & List, 1988; Lee & Chu, 2003). Herein, the entrainment coefficient, ${\alpha}$=0.07 is computed as a fitting parameter in equation 2. The spreading vs. entrainment ratio, $\beta/\alpha$, has been measured as 2 for pure jets (Lee & Chu, 2003), and 1.2 for bubble plumes (Lai & Socolofsky, 2019). Our results show $\beta/\alpha$ = 0.114/0.07 = 1.63, situating the particle-laden plume between pure jets and bubble plumes.

We compute particle number density by counting the particles across each sample (detection method described in section 2.3.2). The normalized distribution of particle number density $\phi$ across the radius of the plume is shown in figure 9. Here, $\left<\phi\right>$ signifies the probability of finding a particle in the specified radial location, and $\phi_{c}$ is the centerline value, which measures as 0.11 at all axial locations. The distribution of $\phi$ allows us to fit a Gaussian profile and measure the half-width ($b_{\phi}$), which we use to designate the ‘particle-core’ from here onward. The axially averaged $b_{\phi}$ obtained from the particle concentration profiles, when normalized by the axially averaged $b_{g}$, shows that $b_{\phi}/b_{g}$=0.56. This ratio for the bubble plume in Lai & Socolofsky (2019) was not measured, but in earlier studies of bubble plumes, e.g. Milgram (1983), reported values of $b_{\phi}/b_{g}$ are in the range of 0.8-0.9. This value is somewhat higher than our value, implying that solid particles spread less rapidly than bubbles. This may be due to the fact that rising bubbles exhibit swirling motions and experience significant lateral lift force when compared to particles (Lai & Socolofsky, 2019). In addition, large-eddy simulations (LES) for a range of monodispersed and polydispersed bubble plumes in Fraga & Stoesser (2016) show that the ratio $b_{\phi}/b_{g}$ is very sensitive to the size distribution of bubbles across the plume. They showed that the size distribution across the plume resembles a reverse-Gaussian profile with larger bubbles populating away from the centerline. This reveals the complexity related to clearly defining a bubble-core which depends on the distribution of the number-density as well as bubble-size across the plume. A future investigation on bubble plumes examining the relationship between polydispersity and turbulence characteristics will help to better explain the above differences between the particle and bubble plumes.

3.2 Fluctuating flow characteristics

In this section, we examine the nature of the fluctuating components of velocity ($u_{i}=U_{i}-\left<U_{i}\right>$; $u_{i,rms}=\sqrt{\left<u_{i}^{2}\right>}$). The normalized probability density functions (PDF) of $u_{i}$ at the plume centerline is shown in figure 10. The respective PDF for a bubble plume (Lai & Socolofsky, 2019) is also included in each plot for comparison. Interestingly, the distribution of the axial velocity fluctuations is negatively skewed for a particle plume, such that axial fluctuations opposite to the direction of particle motion are more common than those along the direction of particle motion. This behavior is strikingly opposite to what is seen in bubble plumes and homogeneous bubble swarms (Riboux et al., 2010; Lai & Socolofsky, 2019; Prakash et al., 2016), for which the axial velocity fluctuations moving in the direction of bubble motion are more common than those moving in the opposite direction. Another way to see this effect is that the mode of the distribution is slightly positive for the particle plume, while it is slightly negative for a bubble plume (see figure 10(a)). One implication of this contrasting result is that the fluid flow in a multiphase plume is sensitive to the direction of the plume with respect to gravity; particle plumes are not a simple reversal of bubble plumes.

The PDFs of radial and out-of-plane velocity fluctuations are symmetric about their means (figure 10(b)), and closely follow a Gaussian curve. We do not observe the prominent signatures of intermittency in the cross-stream components typically observed in bubble plumes and bubble swarms (Lai & Socolofsky, 2019; Prakash et al., 2016).

In figure 11, we further examine the velocity fluctuations across the plume at different radial locations ($x_{2}/b_{g}$). For this purpose, we sort the data into eight equal bins (bin width = 0.16$b_{g}$) for $x_{2}\geq 0$ after checking that the choice of bin width does not change the results. Across the plume, the radial ($u_{2}$) and out-of-plane ($u_{3}$) velocity fluctuations are symmetric with nearly zero skewness (see figure 11(b)). The axial velocity fluctuations ($u_{1}$) switch from negatively skewed at the centerline (as discussed earlier) to positively skewed outside the half-radius ($x_{2}>0.5b_{g}$) of the plume. These characteristics are better captured in the statistical moments (rms, skewness and kurtosis) as shown in figure 12. Figure 12(b) shows that the change of sign in $S(u_{1})$ occurs at three-quarters of the plume width $b_{g}$, which also coincides with the maximum $u_{1,rms}$ (see figure 12(a)). Both $u_{2,rms}$ and $u_{3,rms}$ show their maximum at the centerline, and are consistently smaller than $u_{1,rms}$. The radial variation of kurtosis $K(u_{i})$ in figure 12(c) captures the increasing flatness of each distribution outside the plume half-radius.

3.3 Reynolds stresses

The Reynolds stresses across the plume are reported at all measured axial locations in figure 13(a) (normal stresses) and 13(b) (shear stress in the measurement plane). The turbulent kinetic energy ($k$) based on 13(a) is shown separately in figure 13(c). A nonlinear least-squares fit of the data to a shape-preserving function (see Appendix A) captures well each profile (see fitted lines in figure 13). Results show that the turbulent kinetic energy is primarily dominated by the axial velocity fluctuations, and that it increases from the plume centerline to its maximum value near the edge of the plume ($x_{2}\approx 0.75b_{g}$). The maximum TKE located at $x_{2}/b_{g}=0.75$ is about 44% of $U_{c}^{2}$ and is about 1.5 times the TKE at the plume center. The in-plane shear stress, $\left<u_{1}u_{2}\right>/U_{c}^{2}$, follows a trend similar to that of a turbulent jet and increases from zero at the centerline to a maximum located near $x_{2}/b_{g}=0.75$ (see figure 13(b)). The location of maximum shear is consistent with that of $\left<u_{1}^{2}\right>/U_{c}^{2}$, suggesting that the shear stress is dominated by the axial fluctuations $u_{1}$.

In figure 14, we compare the in-plane turbulent intensities normalized by the local mean axial velocity, ($\sqrt{\left<u_{1}^{2}\right>}/\left<U_{1}\right>$), with two earlier studies on bubble plumes (Duncan et al., 2009; Lai & Socolofsky, 2019). For this comparison, we use the shape functions (solid and dashed lines in figure 13(a). The present data and the data from Duncan et al. (2009) show reasonably similar trends with unbounded growth away from the centerline as $\left<U_{1}\right>$ approaches zero asymptotically outside the plume. The results from Lai & Socolofsky (2019) show somewhat higher turbulent intensity inside the plume core. Also, the growth in their normalized turbulent intensities is not unbounded. Lai & Socolofsky (2019) attribute this deviation to finite $\left<U_{1}\right>$ outside of the plume caused by flow recirculation set up by the tank walls. We do not observe such wall effects in our experiment. The centerline turbulent intensities for the present data are 18%, 11% and 9% for $u_{1}$, $u_{2}$, $u_{3}$, respectively. The axial turbulent intensity ($\sqrt{\left<u_{1}^{2}\right>}$) at the centerline for a bubble plume (Lai & Socolofsky, 2019) is also approximately twice the other two normal intensities ($\sqrt{\left<u_{2,3}^{2}\right>}$), suggesting stronger anisotropy in multiphase plumes when compared to a single-phase jets/plumes in which the ratio $\sqrt{\left<u_{1}^{2}\right>/\left<u_{2,3}^{2}\right>}$ is about 1.4 (Wang & Law, 2002).

3.4 Conservation of kinematic momentum flux of the plume

Using the method in Lai & Socolofsky (2019), we assess the conservation of the total kinematic momentum flux, $M=\left<{M}\right>+M_{t}$, of the induced liquid flow in our particle plume in figure 15. Here, the momentum flux contributed by the mean flow is $\left<{M}\right>=2\pi\int_{0}^{\infty}x_{2}U_{1}^{2}dx_{2}$. The momentum flux carried by the turbulence fluctuations is computed as $M_{t}=2\pi\int_{0}^{\infty}x_{2}\left(\left<u_{1}^{2}\right>-0.5\left(\left<u_{2}^{2}\right>+\left<u_{3}^{2}\right>\right)\right)dx_{2}$. The results are compared to the analytical expression for $M(x_{1})$ for pure vertical plumes (Lee & Chu, 2003),

[TABLE]

Here, $x_{0}=-5.6d_{0}$ is the virtual origin of the plume, and $F_{0}=Q_{0}g$ is the buoyancy flux of the particles. Our results show reasonably good agreement with this model, suggesting that the particle plume obeys the scaling law $M\sim x_{1}^{4/3}$, which is consistent with buoyancy-driven plumes. The total momentum flux is primarily contributed by the momentum flux due to the mean flow $\left<M\right>$. The remaining contribution comes from $M_{t}$, and it is customary to quantify the contribution by the local momentum amplification factor $\gamma=M/\left<{M}\right>$. For our particle plume, $\gamma$ is 1.2, averaged across all axial measurements. This result shows that $\gamma$ for a particle plume is larger than that for a single-phase jet/plume ($\gamma=1.07-1.09$ Wang & Law (2002)) and smaller than a bubble plume ($\gamma=1.4-1.6$, Lai & Socolofsky (2019)).

3.5 Velocity triple-correlation and turbulent transport

Because the PDF of axial velocity fluctuations in a particle and a bubble plume are oppositely skewed (see figure 10(a)), some interesting differences between the two flows can be identified in terms of the velocity triple-correlation terms. These terms contribute to the transport of turbulent kinetic energy and thus are important in the TKE budget (section 3.7). The radial profiles of the triple-correlation of in-plane velocity fluctuations normalized by $U_{c}^{3}$ are shown in figure 16. To compare the trends we also show the respective profiles for a single-phase jet (Darisse et al., 2012) and a bubble plume (Lai & Socolofsky, 2019). The profiles extracted from the two references are multiplied by a factor for visual comparison (see legends in figure 16).

Other than $\left<u_{1}^{3}\right>$, all triple-correlation profiles for the particle plume show trends similar to a turbulent single-phase jet. The triple-correlation profiles of the first two terms ($\left<u_{1}^{3}\right>/U_{c}^{3}$ and $\left<u_{1}u_{2}^{2}\right>/U_{c}^{3}$) for the bubble plume in Lai & Socolofsky (2019) are significantly different from the single-phase jet and our particle plume, and they show 10–20 times larger magnitude compared to our particle plume results (see figure 16(a) and 16(b)). The bubble plume shows a positive $\left<u_{1}^{3}\right>$ near the centerline as opposed to the negative $\left<u_{1}^{3}\right>$ observed near the core of the particle plume, due to the opposite skewness of their $u_{1}$ distributions discussed earlier (see figure 10(a)).

Based on the triple-correlation terms computed above, we estimate the axial and radial transport of $k$ in figure 17(a) and 17(b), respectively. Since the TKE in the particle plume is primarily contributed by the axial stress ($\left<u_{1}^{2}\right>$), these two profiles are nearly identical to those in figures 16(a) and 16(c). These profiles are fitted to two shape functions (see appendix A) and are shown as solid lines in figures 17(a) and 17(b), respectively. We use them for estimating the TKE budget in section 3.7.

3.6 Mean square gradients of velocity fluctuations

To test the conditions of local isotropy and local axisymmetry (George & Hussein, 1991), we compute the in-plane derivatives of all three components of velocity fluctuations employing the $\mathrm{2^{nd}}$-order-accurate central-difference scheme. Table 2 and figure 18 show the second moments of axial ($\partial/\partial x_{1}$) and radial ($\partial/\partial x_{2}$) derivatives of all three components of velocity. These moments clearly do not satisfy the conditions of local isotropy (e.g. $\left<\left(\frac{\partial u_{1}}{\partial x_{1}}\right)^{2}\right>\neq\left<\left(\frac{\partial u_{2}}{\partial x_{2}}\right)^{2}\right>$ and $\left<\left(\frac{\partial u_{1}}{\partial x_{2}}\right)^{2}\right>\neq 2\left<\left(\frac{\partial u_{1}}{\partial x_{1}}\right)^{2}\right>$) nor of local axisymmetry (e.g. $\left<\left(\frac{\partial u_{3}}{\partial x_{1}}\right)^{2}\right>\neq\left<\left(\frac{\partial u_{2}}{\partial x_{1}}\right)^{2}\right>$).

The violation of local axisymmetry was also observed for bubble plumes (Lai & Socolofsky, 2019), suggesting that this is a common characteristic of multiphase plume turbulence contrary to single-phase jets.

Computing the true viscous dissipation rate (stress power) requires all nine components of the velocity gradient tensor and each component has to be adequately resolved, which is not easy to achieve experimentally. Instead, we compute the pseudo-dissipation rate $\left<\epsilon\right>$:

[TABLE]

Here, $\nu=1.02\times 10^{-6}$ $\mathrm{m^{2}/s}$ is the kinematic viscosity of saline water at 23∘ C. In our calculations all but the three out-of-plane derivative terms in equation 4 are computed directly from our stereo-PIV data, and the missing out-of-plane mean-square derivatives ($\left<\left(\frac{\partial u_{i}}{\partial x_{3}}\right)^{2}\right>$) are replaced by $\left<\left(\frac{\partial u_{i}}{\partial x_{2}}\right)^{2}\right>$. The mean pseudo-dissipation rate $\left<\epsilon\right>$ across the plume is shown for various axial locations in figure 19.

To assess the accuracy in the dissipation calculation, we adopt a method based on conservation of total kinetic energy for multiphase plume in a Lagrangian framework as outlined in Lai & Socolofsky (2019). In essence, this method estimates the balance of the total kinetic energy as production minus dissipation:

[TABLE]

and calculates the correction factor $f_{c}$ to balance the two sides. The mean kinetic energy, $\left<K\right>\approx\left<U_{1}\right>^{2}/2$ because $U_{1}\gg U_{2},U_{3}$. This equation can be simplified in the cylindrical coordinate system ($r,\theta,z$) by invoking steady state as follows. The left hand side of the equation is computed using the shape functions of the mean radial profiles for $\left<U_{1}\right>$ and $\left<k\right>$ (see section 3.1 and 3.2). The first term on the right hand side is obtained by using $\left<\epsilon\right>$ at multiple axial locations. The last term in equation 5 refers to the kinetic energy production due to the particle buoyancy flux and equals $CF_{0}\Delta x_{1}$, where the prefactor $C=\frac{1+\lambda^{2}}{\pi\lambda^{2}}\sqrt{4/\gamma}\int_{0}^{\infty}2\pi\eta e^{-\left(\frac{2\lambda^{2}+1}{\lambda^{2}}\right)\eta^{2}}d\eta$ depends on $\lambda=\left<b_{\phi}\right>/\left<b_{g}\right>$ and the momentum amplification factor $\gamma$. Based on $\lambda=0.6$ and $\gamma=1.2$, we obtain $C=1.44$. The simplified equation for a volume bounded by two axial locations $x_{1}$ and $x_{2}$ is written as:

[TABLE]

Here, $I_{1}$, $I_{2}$ and $I_{3}$ are the axisymmetric surface integrals ($2\pi\int_{0}^{\infty}\left(.\right)rdr$) of the respective terms inside the volume integral in equation 5 (summarized in Table 3).

We test the balance of equation 6 at multiple axial locations separated by $1.5b_{g}$ and have found $f_{c}$ to vary between 2.8 and 5.1 with a mean value of 4.0. This method estimates the amount by which $\left<\epsilon\right>$ is underestimated by our evaluation of equation 7, but it is unable to correct $\left<\epsilon\right>$ locally along the radius of the plume. Because of this limitation in our measurement of $\left<\epsilon\right>$, we use the corrected dissipation profile only for a qualitative assessment of the TKE budget discussed in the following section.

3.7 Turbulent kinetic energy budget

In the cylindrical coordinate system, the time-averaged transport equation for $k$ of an axisymmetric, steady, turbulent plume without swirl is (Kataoka & Serizawa, 1989):

[TABLE]

In addition to the terms already present in single-phase jets/plumes ($A$: advection, $P_{s}$: shear-production, $T$: turbulent transport, and $T_{p}$: pressure transport), a multiphase plume has another source term ($P_{b}$) that represents the interfacial energy transfer at the boundaries between the dispersed phase and the surrounding fluid. Obtaining the complete TKE budget for a particle plume is challenging for multiple reasons. First, it is difficult to accurately measure experimentally the mean dissipation rate $\left<\epsilon\right>$ because of the under-resolved velocity gradients (as seen in section 3.6). Second, there is no existing model for pressure-velocity correlation in multiphase flow turbulence, and therefore the pressure transport term ($T_{p}$) requires an assumption. Finally, the interfacial energy transfer ($P_{b}$) is the most complex and inaccessible term, and it requires the local distribution of particles and their relative velocity with the surrounding fluid (Santarelli et al., 2016). In an ideal scenario with all other terms measured correctly, the term $P_{b}$ can be sought as the closing term of equation 7, such that $P_{b}=\left<\epsilon\right>-(A+P_{s}+T+T_{p})$. We will use this method to estimate $P_{b}$ from our data.

We directly compute $A$, $P_{s}$ and $T$ from our experimental data. The expressions for these in terms of similarity variables are summarized in Appendix B. Each plot contains an error-band corresponding to the 95% confidence interval in the constituent fits. Each uncertainty band is computed via the rule of uncertainty propagation described in Charonko & Prestridge (2017). The pressure transport term ($T_{p}$), not directly available from experiment, is replaced by Lumley’s model $T_{p}=-2/5T$, widely used for single-phase turbulence (Lumley, 1978). These four terms and the mean dissipation rate ($\left<\epsilon\right>$) multiplied by the correction factor $f_{c}=4$ are shown in figure 20(a). The remaining terms of the balance equation are the interfacial energy transfer ($P_{b}$) from particles to fluid and out-of-balance experimental error (OOBE). The budget terms for a bubble plume from Lai & Socolofsky (2019) are reproduced in figure 20(b) for comparison. Each budget term is normalized by $U_{c}^{3}/b_{g}$; the order-of-magnitude difference in the two budgets arises because of the difference in the downstream measurement locations with respect to the source, leading to the differences in $U_{c}$ and $b_{g}$ (see figures 7(a) and 7(b)). Below we discuss the comparison of the TKE budget for the two multiphase plumes in light of previous experiments on a single-phase jet (Lai & Socolofsky, 2018) and a variable-density single-phase (VDSP) jet (Charonko & Prestridge, 2017). The qualitative summary of this comparison is provided in table 4.

The shape of the advection ($A$) profiles are very similar across all four plumes with a positive central lobe inside the plume core. The shear production term ($P_{s}$) in the particle plume increases from zero at the centerline to a dominant peak at about $r/b_{g}=0.65$ near the edge of the plume. This behavior strongly resembles a single-phase jet both in shape and relative magnitude. The $P_{s}$ profile in the VDSP jet also matches this trend, except that a negative $P_{s}$ is observed at the jet centerline. By contrast, the bubble plume shows almost negligible $P_{s}$ across the radius of the plume. The turbulent transport term ($T$) in the particle plume exhibits a positive peak near the centerline and a negative lobe near the edge of the plume. The bubble plume and the VDSP jet are similar, but with a larger peak in the former. The single-phase jet does not match the others, i.e. near the centerline the turbulent transport is consistently observed to be negative or marginally positive (Darisse et al., 2015).

Because the energy production at the two-phase boundaries ($P_{b}$) is computed indirectly, the peak location in $P_{b}$ is sensitive to the experimental error in $\left<\epsilon\right>$. Despite this limitation, it is worth noting that in both particle and bubble plumes, the sum ($A+T+P_{s}+T_{p}$) is significantly smaller than the mean dissipation rate outside the plume core ($r>0.5b_{g}$). This results in the TKE budget being an approximate balance between $P_{b}$ and $\left<\epsilon\right>$.

3.8 Two-point correlation and energy spectra

We obtain the two-point statistics from the spatial autocorrelation function $R_{11}$ computed along the axial measurements of axial velocity fluctuations, $u_{1}$. Figure 21 shows $R_{11}$ inside the plume core ($x_{2}/b_{g}<0.5$) at various radial locations (shaded circles) and their mean $\bar{R}_{11}$ (red circles). The horizontal axis is normalized by the particle diameter $d_{p}$. We estimate the Taylor microscale ($\lambda_{f}$) of this flow by fitting an osculating parabola, $1-\frac{x_{1}^{2}}{\lambda_{f}^{2}}$ (Pope, 2000) to the first four points of $\bar{R}_{11}$ (see solid line in figure 21), yielding $\lambda_{f}=1.3d_{p}$, which is of the order of the particle size. The Taylor microscale ($\lambda_{f}$) is related to the mean square velocity gradient in homogeneous isotropic turbulence (Pope, 2000) and we use this relation to approximate

[TABLE]

Based on this method, $\left<\left(\frac{\partial u_{1}}{\partial x_{1}}\right)^{2}\right>\left(b_{g}/U_{c}\right)^{2}$ at $x_{2}=0$ is 3.2$\pm$0.6 and at $x_{2}=0.5b_{g}$ is 5.8 $\pm$0.6. Within the uncertainty bounds, these values match those obtained from direct computation of gradients in section 3.6 (see figure 18(a)).

Risso (2018) suggests that the inter-scale energy transfer in bubble-induced agitation (BIA) is quite different from the classical picture of shear-induced turbulence (SIT). As a first step toward understanding the energy cascade in particle-plume turbulence, we generate the one-dimensional energy spectrum ($E_{11}$) from the autocorrelation function $\bar{R}_{11}$ (Pope, 2000). The normalized spectra for our particle plume is compared with experimental results (Riboux et al., 2010) and DNS simulations (Lai et al., 2018) in figure 22(a). The result shows striking similarity in the spectra between homogeneous bubble swarms (both experiment and DNS simulation) and our particle plume, all following the earlier prediction of a $\kappa^{-3}$ power law for wave numbers smaller than $1/\lambda_{f}$. Beyond $\kappa_{1}>1/\lambda_{f}$, both experimental results (ours and Riboux et al. (2010)) show different slopes than the DNS results, which can be attributed to experimental uncertainty in resolving velocity fluctuations below $\lambda_{f}$. Nonetheless, this result is valuable and it suggests power-law universality of multiphase flow turbulence which is quite different from that for single-phase turbulent shear flows.

The one-dimensional dissipation spectra $D_{11}(\kappa_{1})=2\nu\kappa^{2}E_{11}(\kappa_{1})$ is shown in figure 22(b). Quite strikingly, the peaks of dissipation for our data and bubble plume/swarms reside at $\kappa_{1}d_{p}\approx 0.2-0.4$. This indicates that in both bubble– and particle–laden turbulence, production and dissipation occur near the scale of particle size. This lack of scale separation was first postulated for homogeneous bubbly flows in Lance (1991) and is observed in direct numerical simulations reported in Lai et al. (2018).

4 Conclusions

We report an experimental study characterizing the turbulence inside a heavy particle plume descending under gravity within a salt-water solution. We measure the three components (2D3C) of the interstitial fluid velocity and the spatial distribution of particles in the central plane of the plume using refractive-index-matched stereoscopic particle image velocimetry. Below we summarize our key findings primarily in the light of the results for a bubble plume (Lai & Socolofsky, 2019).

The induced liquid flow inside the plume evolves with the mean flow characteristics of a bubble plume (Lai & Socolofsky, 2019), that show Gaussian mean axial velocity ($\left<U_{1}\right>$) profile. The radial profiles of mean particle number density are also Gaussian with a half-width $b_{\phi}$ equals to 0.56$b_{g}$, where $b_{g}$ is the Gaussian half-width of the mean streamwise velocity profile in the fluid.

The turbulence inside the plume is highly anisotropic with maximum streamwise turbulence intensity measuring up to 2.7 times that of the other two components, a result consistent with that observed in a bubble plume. The PDF of the axial velocity fluctuation ($u_{1}$) at the centerline of the plume is skewed so that the strongest events are in the direction opposite to mean flow. This behavior is strikingly opposite to that observed previously for a bubble plume and a homogeneous swarm of bubbles. The turbulent kinetic energy (TKE) and the in-plane shear-stress peak near $\approx 0.75b_{g}$, locating the shear layer slightly outside the edge of the particle-core ($b_{\phi}$). Although these quantities for a bubble plume in Lai & Socolofsky (2019) peaked at a similar location ($\approx 0.7b_{g}$), the location of the shear layer with respect to their bubble-core ($b_{\phi}$) remains unknown.

Another distinction between the particle and the bubble plume is observed in the fluid shear production term ($P_{s}$) in the TKE budget. While the shear production in the bubble plume is negligible across the plume compared to the other terms, we observe a distinct peak in $P_{s}$ at the shear layer region of a particle plume, resembling a single-phase jet.

The other TKE terms including advection ($A$) and the turbulent transport ($T$) in the two flows exhibit similar profiles. Further, we show that the turbulent transport ($T$) term follows similar qualitative profiles for both plumes and also an earlier reported variable-density single-phase jet (Charonko & Prestridge, 2017), all of which differ from a typical single-phase jet. The difference in the relative magnitude in $T$ at the centerline between the two plumes is attributed to the difference in their third moment (skewness) of the axial velocity fluctuations.

Despite the above differences inside the core, both particle and bubble plumes show qualitatively that the TKE production by the particles ($P_{b}$) approximately balances the mean dissipation rate ($\left<\epsilon\right>$) away from the centerline. Further, the one-dimensional spectrum in the particle plume exhibits the $\mathrm{-3^{rd}}$ power-law consistent with bubble plume and homogeneous swarms of bubbles. These two results support the notion that there is a lack of separation between the scales of production and dissipation in multiphase turbulent flows like ours.

Acknowledgement

This research was supported by the National Science Foundation (Grant No. OCE- 1334788) and by The Gulf of Mexico Research Initiative. The authors also thank Michael Heath for his help in the design of the experimental setup and data collection.

Appendix A

The three components of normal stresses, the turbulent kinetic energy ($k$), and the mean dissipation rate $\left<\epsilon\right>$ are fitted into a shape preserving double-Gaussian curve, given in equation 9. All terms are expressed in similarity variable $\eta=r/b_{g}$.

[TABLE]

The in-plane shear-stress is fitted into a polynomial-Gaussian function in equation 10.

[TABLE]

The axial and radial transport of the turbulent kinetic energy ($k$) are fitted into the polynomial-Gaussian functions in equations 11 and 12, respectively.

[TABLE]

[TABLE]

Appendix B

Expressions for the advection ($A$), shear-production ($P_{s}$), and turbulent tranport ($T$) terms in cylindrical co-ordinate system in the turbulent kinetic energy budget are summarized below. Each term is expressed in terms of the similarity variable $\eta=r/b_{g}$.

• •

Advection term

[TABLE]

[TABLE]

• •

Production term

[TABLE]

[TABLE]

• •

Turbulent Transport

[TABLE]

[TABLE]