1. Introduction
Fluid flows in both the natural world and in engineering applications are rarely homogeneous; rather, they often carry suspended particulate material. The dynamics of such multiphase flows can be extremely complex (Balachandar & Eaton Reference Balachandar and Eaton2010), even when the particles themselves are simple and passively advected by the flow. In the more general case where the particles are ‘active’ and can generate their own forces or modulate their interaction with the fluid in addition to being carried by the flow, we have little general predictive capability. A particular example of such active particles that is our focus here are those whose effective buoyancy can change in response to the local flow around them.
Actively buoyant particles are not uncommon, and include, for example, living ‘particles’ such as algae or plankton that modulate their buoyancy for a variety of reasons (Reynolds, Oliver & Walsby Reference Reynolds, Oliver and Walsby1987; Gemmell et al. Reference Gemmell, Oh, Buskey and Villareal2016; Arrieta et al. Reference Arrieta, Jeanneret, Roig and Tuval2020). Our primary motivation for studying actively buoyant particles, however, is the problem of so-called spot ignition by transported firebrands in wildfires. Firebrands are small pieces of burning or smouldering material that are ejected from actively burning regions. If firebrands are sufficiently lofted and are entrained into the atmospheric boundary layer, they can be transported far downstream of the burning front, potentially leading to the ignition of ‘spot fires’ in distant locations (Koo et al. Reference Koo, Pagni, Weise and Woycheese2010; Fernandez-Pello Reference Fernandez-Pello2017; Manzello et al. Reference Manzello, Suzuki, Gollner and Fernandez-Pello2020). Predicting the landing sites of firebrands is critically important for wildfire containment and mitigation, but remains an unsolved challenge (Koo et al. Reference Koo, Pagni, Weise and Woycheese2010; Thurston et al. Reference Thurston, Kepert, Tory and Fawcett2017; Or et al. Reference Or, Furtak-Cole, Berli, Shillito, Ebrahimian, Vahdat-Aboueshagh and McKenna2023).
Fire spotting is often broken into three logically distinct elements: generation of firebrands in the burning part of the fire, their transport downstream and the ignition of unburned fuel at the landing site (Koo et al. Reference Koo, Pagni, Weise and Woycheese2010). Firebrand generation and fuel ignition are primarily combustion problems, and are not our focus here. The transport of firebrands, however, is a complex multiphase flow problem. The current state of the art in modelling firebrand transport makes a number of assumptions that may not be valid (Tohidi & Kaye Reference Tohidi and Kaye2017). For example, it is common to assume that firebrands always travel at their terminal velocity relative to the wind speed (Thomas, Sharples & Evans Reference Thomas, Sharples and Evans2020), and their drag is estimated using empirical drag coefficients determined at isothermal conditions (Albini, Alexander & Cruz Reference Albini, Alexander and Cruz2012; Dal-Ri Dos Santos Reference Dal-Ri dos Santos and Yaghoobian2023). These assumptions have been shown to overestimate transport distances when compared with measurements (Thomas et al. Reference Thomas, Sharples and Evans2020). The firebrands relevant to spotting continue to undergo combustion reactions as they travel; however, the effect of these reactions is often only incorporated by a modelled mass loss and temperature change (Fernandez-Pello Reference Fernandez-Pello2017), without consideration for how these characteristics affect the particle motion aside from modifying the particle size. Importantly, heat transfer from the particle to the air and its consequences is largely neglected, even though fundamental studies of heated particles moving in flows have demonstrated that there can be strong two-way coupling effects on the local particle and fluid motion (Frankel et al. Reference Frankel, Pouransari, Coletti and Mani2016; Zamansky et al. Reference Zamansky, Coletti, Massot and Mani2016).
A potentially key factor that is missing from current models is the dynamic modification of the firebrand’s effective buoyancy due to heat transfer. The elevated temperature of the firebrand relative to the surrounding air far from the burning fire makes them more buoyant than would be expected since they can heat a ‘jacket’ of fluid around them through thermal conduction that is then less dense than the rest of the air. Additionally, this heated air will tend to rise, exerting an upward drag force on the particle while also changing the viscous drag exerted on the particle from the surrounding cold fluid. Convection, however, can act to suppress this enhanced effective buoyancy and modulated drag by advecting away the scalar that is its source and thinning the low-density jacket, but in a way that depends on the local flow. These density effects are illustrated schematically in figure 1.
Some theoretical results exist for the drag on a heated sphere in the limit of low particle Reynolds number (Mograbi & Bar-Ziv Reference Mograbi and Bar-Ziv2005a
,
Reference Mograbi and Bar-Zivb
; Ganguli & Lele Reference Ganguli and Lele2019, Reference Ganguli and Lele2020) or in the Boussinesq limit (Kotouč et al. Reference Kotouč, Bouchet and Sušek2009), and demonstrate that, even for these simplified cases, the problem is complex. Firebrands, however, are not in either of these limits. Typical firebrands in wildfires have linear sizes of a few centimetres (Manzello, Maranghides & Mell Reference Manzello, Maranghides and Mell2007) and terminal velocities of a few to tens of metres per second (Sánchez Tarifa et al. Reference Sánchez Tarifa, Pérez del Notario and García Morena1965; Oliveira et al. Reference Oliveira, Lopes, Baliga, Almeida and Viegas2014; Thurston et al. Reference Thurston, Kepert, Tory and Fawcett2017; Thomas et al. Reference Thomas, Sharples and Evans2020), suggesting particle Reynolds numbers of order
$10^4$
, far from the creeping-flow limit. By the same token, their typical temperatures of
$500{-}700$
$^\circ$
C (Anthenien, Tse & Fernandez-Pello Reference Anthenien, Tse and Fernandez-Pello2006) suggest that the heat transfer from a firebrand to the surrounding air is far from Boussinesq. Thus, appropriately modelling the behaviour of a transported firebrand requires working in this difficult regime where common simplifying assumptions are not valid.

Figure 1. Schematic of the density field near an actively buoyant particle. The particle density
$\rho _p$
is larger than the fluid density
$\rho _{\!f}$
. Due to the buoyant scalar attached to the particle that mixes into the fluid directly around the particle, there is a ‘jacket’ of fluid of density
$\rho _s \lt \rho _{\!f}$
around the particle that can change the effective density of the particle. Because of relative motion between the particle and the fluid and the tendency of the lightweight jacket fluid to rise, however, convection will strip this lighter fluid away into a wake behind the particle where the density
$\rho$
is spatially dynamic and lies between
$\rho _s$
and
$\rho _{\!f}$
, further changing the effective density and the drag on the particle.
Since a first-principles theoretical analysis of the dynamics of an actively buoyant particle in the relevant parameter range is therefore intractable, we conducted a series of experiments to measure the settling of actively buoyant particles in a quiescent fluid. Before presenting our results, we first identify the relevant non-dimensional parameter space that characterises this system. We then describe our experimental methods, including our choice of buoyant scalar and the fabrication of appropriate particles, before moving on to our results. We find that the settling velocity of our actively buoyant particles can be either suppressed or enhanced relative to what would be expected for an inert particle, but that there is no simple dependence of the results on a single control parameter. We also analyse the trajectories of the settling particles. We find that these trajectories are well described as single-mode sinusoids with an approximately constant Strouhal number even though the Reynolds number range would suggest that the trajectories ought to be more complex. Taken together, our results indicate that the effects of active buoyancy on the particle dynamics are complex and not small, and thus that understanding the behaviour of such particles is ripe for further study.
2. Methods
2.1. Relevant parameters
Even the simplest situation of an actively buoyant particle settling in a quiescent fluid is governed by a large number of parameters. The background fluid is characterised by its mass density
$\rho _{\!f}$
and dynamic viscosity
$\mu$
. Likewise, the particle has a density
$\rho _p$
and a diameter
$d$
. Because in general
$\rho _{\!f} \ne \rho _p$
, buoyancy effects will be relevant so that the gravitational acceleration
$g$
must enter the problem. Given our motivation, we consider the active buoyancy of the particle to arise due to an internal reservoir of a scalar associated with the particle that can mix into the surrounding fluid. This scalar is characterised by a density
$\rho _s$
, which we define as the density of the fluid–scalar mixture when the scalar takes its peak value (e.g. the density of air heated to the temperature of a firebrand in the wildfire context), and a diffusion constant
$\kappa$
. Note that we assume both that the scalar reservoir is effectively infinite and that the particle density is constant in time; although these approximations are never strictly valid, in many practical situations these properties change on much slower time scales than the particle velocity. In the settling problem, we seek to determine the terminal settling velocity
$w$
as a function of these input parameters. Application of the Buckingham Pi theorem indicates that the non-dimensional settling velocity can be written as a function of four non-dimensional input parameters.
Previous work on the settling of finite-size particles has indicated that mass effects and size effects are important, and are not interchangeable (Cabrera-Booman, Plihon & Bourgoin Reference Cabrera-Booman, Plihon and Bourgoin2024); we therefore choose two of these input parameters to be the specific gravity
$\varGamma = \rho _p/\rho _{\!f}$
and the Galileo number
\begin{equation} Ga = \sqrt {\frac {gd^3 \rho _{\!f} |(\rho _p - \rho _{\!f})|}{\mu ^2}} = \sqrt {\frac {gd^3 \rho _{\!f}^2 | \varGamma - 1|}{\mu ^2}}, \end{equation}
which compares the effects of the (passive) particle buoyancy and viscous effects. To capture effects of the buoyancy due to the scalar (that is, the active contribution to the particle buoyancy), we use the particle Grashof number
where we have assumed that the scalar density varies linearly with concentration. We refer to
${\textit{Gr}}_p$
as a particle Grashof number since we use the particle length scale
$d$
rather than a (larger) length scale defined by the buoyant plume attached to the particle. Finally, we also consider the Prandtl number
$Pr = \mu / (\rho _{\!f} \kappa )$
(or, equivalently, the Schmidt number
${\textit{Sc}}$
if the scalar is not heat).
We non-dimensionalise the settling velocity using what is sometimes referred to as the settling number
${\textit{Sv}} = w / w_{\textit{ex}}$
, where
$w_{\textit{ex}}$
is the expected settling velocity. If the particle were both small and close to neutrally buoyant, it could be appropriate to set
$w_{\textit{ex}}$
to the Stokes settling velocity. In the more general case where the particle inertia is non-negligible, it is more appropriate to define
\begin{equation} w_{\textit{ex}} = \sqrt {\frac {4gd(\rho _p-\rho _{\!f})}{3 \rho _{\!f} C_d}}, \end{equation}
where
$C_d$
is a drag coefficient. In general, for passive particles,
$C_d$
depends on both
${\textit{Ga}}$
and
$\varGamma$
. Here, we use the recent parameterisation reported by Cabrera-Booman et al. (Reference Cabrera-Booman, Plihon and Bourgoin2024) for
$C_d$
and therefore
$w_{\textit{ex}}$
, as it conveniently expresses
$C_d$
purely as a function of
${\textit{Ga}}$
. This parameterisation is given by
\begin{equation} C_d(Ga) = \frac {4}{3} \left ( \frac {0.0258 {\textit{Ga}}^{2.6973} + 2.81 {\textit{Ga}}^{2.0306} + 18 {\textit{Ga}}^{1.364} + 405}{Ga (22.5 + {\textit{Ga}}^{1.364})} \right )^2\!. \end{equation}
Other choices for the relevant non-dimensional parameters are, of course, possible. It is common, for example, to define a particle Reynolds number
${\textit{Re}}_p = \rho _{\!f} w\textit{d} / \mu$
. In our case,
$w$
is an output parameter and emerges from a steady-state balance between the forces on the particle rather than being set as the result of an externally imposed flow, and so
${\textit{Re}}_p$
should also be considered to be an output. Here, we prefer to use
${\textit{Sv}}$
as the non-dimensional output variable, as it directly indicates whether the settling velocity is expected or anomalous. Previous studies of heated particles in external flows have shown the importance of the parameter
${\textit{Gr}}_p / {\textit{Re}}_p^2$
in controlling the heat transfer between the particle and the fluid (Acrivos Reference Acrivos1958; Mograbi & Bar-Ziv Reference Mograbi and Bar-Ziv2005a
,
Reference Mograbi and Bar-Zivb
). This ratio, which can be thought of as a Richardson number, distinguishes the limits of heat transfer by natural convection (
${\textit{Gr}}_p / {\textit{Re}}_p^2 \gg 1$
) or by forced convection (
${\textit{Gr}}_p / {\textit{Re}}_p^2 \ll 1$
). When
${\textit{Gr}}_p / {\textit{Re}}_p^2 \sim O(1)$
, the case of mixed convection, the dynamics is expected to be especially complex (Mograbi & Bar-Ziv Reference Mograbi and Bar-Ziv2005a
,
Reference Mograbi and Bar-Zivb
; Ganguli & Lele Reference Ganguli and Lele2020). Although we return to this parameter and the convective regime below, we again note that it is not a control parameter as such since it mixes input (
${\textit{Gr}}_p$
) and output (
${\textit{Re}}_p$
) variables.
2.2. Buoyancy source and particles
By defining the non-dimensional parameter space for this problem, we can take advantage of dynamical similitude to choose a combination of particle, buoyancy source and carrier fluid that allows the application of convenient experimental techniques while still making it possible to examine the fundamental physics governing the behaviour of actively buoyant particles. We use water as the working fluid, both because it is compatible with a wide range of flow visualisation techniques and because its density and viscosity relative to those of air slow down settling time scales, allowing for more precise measurements. Water does make it difficult, however, to use temperature as the buoyancy-modulating scalar, for several reasons. The relatively high thermal conductivity of water means that heat transfer from a hot particle is fairly efficient even by conduction, so that it is difficult to fabricate a small particle with enough internal heat reserves to remain at an elevated temperature for any appreciable length of time. And second, water remains a liquid only up to a temperature of 100
$^\circ$
C, and the density of 100
$^\circ$
C water is only approximately 4 % less than that of room-temperature water. Thus, the excess buoyancy that can be generated by heating water is fairly small.
Instead of temperature, we use ethanol as our source of buoyancy. Like temperature, ethanol concentration is a diffusing scalar field, since ethanol is fully soluble in water. However, the diffusion of ethanol is slower than the diffusion of heat in water, so that elevated levels of ethanol concentration around a particle will persist longer than elevated temperatures. Pure ethanol also has a specific gravity of 0.79, so that the excess buoyancy that can be generated using ethanol is much larger than what can be generated using temperature in water. This choice means that we are essentially considering compositional convection rather than thermal convection near the particle (and so, for that reason, we report the Schmidt number
${\textit{Sc}}$
rather than the Prandtl number
$ \textit{Pr} $
below); however, the underlying physics is sufficiently similar that our results will still have bearing on the wildfire problem. And even though the density contrast between ethanol and water is somewhat smaller than what would be expected for a several-hundred-degree temperature difference in air, it is still sufficiently high to fall in the strongly non-Boussinesq regime. Using ethanol as the buoyancy source imposes requirements on the particles that can be used: we need particles that can hold a reservoir of ethanol that can diffuse into the surrounding water, ideally slowly enough that the reservoir does not deplete rapidly compared with the time scale of the experiments. In addition, as we wish to study settling, the particles must be dense enough that they remain negatively buoyant in water even when their ethanol reservoir is full.

Figure 2. (a) Photograph of the settling chamber used in these experiments. (b) An example 3D-printed mould used for casting the foam particles. The core ball bearing is held in place by alignment screws arranged in a tetrahedral configuration. The two halves of the mould are precisely aligned with pins. (c) An example completed particle. (d) An example image of a falling particle charged with ethanol. For qualitative visualisation of the scalar wake in this example, the ethanol has been dyed blue.
To meet these requirements, we developed a method to fabricate porous particles from a spongy open-cell polyurethane foam that is both water and ethanol compatible. By soaking the particle in a water/ethanol solution of fixed concentration, we could ‘charge’ the particle with a reservoir of a passive scalar with a specified
$\rho _s$
, which then diffused out of the particle when placed in water. As long as the particle was sufficiently large, this internal ethanol reservoir did not appreciably deplete over the course of each experiment. However, because the foam has only a small mass, an ethanol-charged foam particle is positively buoyant in water. We therefore used a core–shell design, casting the foam around a small steel core to add mass. By varying the size of the steel core, we could vary
$\rho _p$
(which we defined as the net mass of the particle divided by its total volume) independently of
$d$
.
Specifically, we used 440C stainless steel ball bearings of various sizes as the particle cores. The outer shell was fabricated from a castable polyurethane foam (Smooth-On FlexFoam-iT! 6). To ensure the ball bearing was centred within the foam shell, we designed a custom 3D-printed mould (figure 2 b). The mould used four screws arranged in a tetrahedral configuration, with each screw bottoming out to suspend the ball bearing precisely at the centre. Two alignment pins were included to ensure alignment between the top and bottom hemispheres of the mould, which were fit together while the initially liquid foam precursor expanded and cured. An example of the final particle is shown in figure 2(c), and an example of a particle settling while loaded with ethanol is shown in figure 2(d). As can be seen in figure 2(c), removing the screws used to hold the ball-bearing core during curing does create four small holes in the particle, and the two-part mould leaves an unavoidable ring around the equator of the particle. These features are small compared with the particle dimensions, however, and so we do not anticipate that they contribute in any significant way to the particle dynamics. Additionally, we did not observe significant differences in the settling behaviour of the particles when dropped in different orientations. Still, to control for any slight systematic differences as a function of these surface imperfections, we kept the initial orientation of the particles the same for all the results reported below.
2.3. Experimental set-up
To study actively buoyant particle settling, we conducted experiments using four particle types. The outer diameter of each particle mould was the same in all cases (2.19 cm), but we used four different ball-bearing cores (with diameters of 0.87, 0.95, 1.11 and 1.27 cm). We loaded the particles with five different water–ethanol solutions: 0 % (pure water), 25 %, 50 %, 75 % and 100 % ethanol by volume. The mass densities of these mixtures (corresponding to
$\rho _s$
in our definitions above) are 1.00, 0.97, 0.93, 0.87 and 0.79 g cm−
$^3$
, respectively. We note that the diffusion constant for a water–ethanol mixture is a function of the ethanol concentration; thus, the Schmidt number for each of the cases differed. In our experiments, however, we could not change
${\textit{Sc}}$
independently from
${\textit{Gr}}_p$
, since we set each by varying the concentration of the ethanol solution. When loaded with fluid, we found that the particles swelled somewhat, so that their diameter was larger than that of the mould. We also found that the swelling was weakly dependent on the ethanol concentration. We therefore measured the volume of each loaded particle using Archimedes’ principle; this method was more accurate than a direct measurement of particle diameter with callipers, given that the particles were soft. We also measured the mass of each loaded particle using a precision scale, allowing us to compute
$\rho _p$
for each experimental case. In table 1, we report the ranges of the control parameters spanned in our experiments.
Table 1. Control parameters for the different experimental cases considered. Here,
$d_{\textit{core}}$
is the diameter of the inner ball-bearing core in each particle, and
$d_{\textit{mold}} = 2.19$
cm is the diameter of the particle mould (and so the outer diameter of the dry particle).

To contextualise these values, and in particular the seemingly very high Grashof numbers, we note again that typical firebrands have sizes of a few centimetres (Manzello et al. Reference Manzello, Maranghides and Mell2007), temperatures of
$500{-}700$
$^\circ$
C (Anthenien et al. Reference Anthenien, Tse and Fernandez-Pello2006) and are made of wood or similar organic material, giving
${\textit{Gr}}_p \sim O(10^6)$
and
${\textit{Ga}} \sim O(10^4)$
. The values of
$\varGamma$
and
${\textit{Sc}}$
in our experiments are, however, quite different from the firebrand case.
In each experimental run, the particle was released from rest into a tank of quiescent water, as shown in figure 2(a). The tank had a square cross-section measuring 25.4 cm on a side, and was 127 cm tall. These dimensions were chosen to minimise potential wall effects on the settling (Miyamura, Iwasaki & Ishii Reference Miyamura, Iwasaki and Ishii1981; Chhabra Reference Chhabra1995) and to be large enough to allow the particles to reach their terminal velocity well before hitting the bottom (Khatmullina & Isachenko Reference Khatmullina and Isachenko2017; Wang et al. Reference Wang, Dou, Ren, Sun, Jia and Zhou2021). Particles were released through a funnel to ensure a consistent initial location, and the initial orientation of the particle was controlled in order to minimise any variability due to small surface deformities resulting from the screw locations in the casting mould. We observed that the particles did not tumble or rotate as they fell, suggesting both that any such deformities did not have a large effect on the particle dynamics and that the centre of mass of the particles was co-located with their geometric centre.
A LED light sheet illuminated the particle from behind, and a grey scale CMOS camera (FLIR Flea3) imaged the particles settling at 100 frames per second. The camera was positioned far enough below the funnel that the particles had reached their terminal velocity in the imaging region. Two-dimensional particle trajectories were reconstructed using a predictive tracking algorithm (Kelley & Ouellette Reference Kelley and Ouellette2011), and the terminal velocity was measured by linear fits to the time series of the vertical particle position (see figure 8 below for examples of these time series). We chose to fit the time series rather than differentiating the data to avoid the amplification of any experimental error due to numerical derivatives. We performed 10 trials for each of the 20 experimental conditions. The fluid in the tank was allowed to come to rest over a period of several minutes between runs; measurements of the fluid velocity after this time period indicated that any residual fluid motion was far smaller the run-to-run variability in our experiments.
3. Results and discussion
3.1. Settling velocity
Our primary goal in this work is to characterise the effects of active buoyancy, as indicated by a non-zero particle Grashof number
${\textit{Gr}}_p$
, on the quiescent settling velocity of a particle. In figure 3, we therefore plot the settling number
${\textit{Sv}} = w / w_{\textit{ex}}$
as a function of
${\textit{Gr}}_p$
. As discussed above, we measured
$w$
by fitting the time series of the vertical position. We computed
$w_{\textit{ex}}$
from (2.3).

Figure 3. The settling number
${\textit{Sv}}$
plotted as a function of the particle Grashof number
${\textit{Gr}}_p$
. The colour of each data point shows the specific gravity
$\varGamma$
, as indicated by the colour bar. The error bars correspond to the standard error over 10 trials for each set of control parameters.

Figure 4. The settling number
${\textit{Sv}}$
plotted as a function of the specific gravity
$\varGamma$
. The colour of each data point shows the particle Grashof number
${\textit{Gr}}_p$
, as indicated by the colour bar. The error bars correspond to the standard error over 10 trials for each set of control parameters.
The data shown in figure 3 have several intriguing features. First, for
${\textit{Gr}}_p=0$
, corresponding to the case with no active buoyancy,
${\textit{Sv}} \approx 0.8$
, indicating that even without active buoyancy our particles settle approximately 20 % more slowly than expected. The cause of this discrepancy is not immediately clear. We note that our particle Reynolds numbers
${\textit{Re}}_p$
are significantly higher than those considered by Cabrera-Booman et al. (Reference Cabrera-Booman, Plihon and Bourgoin2024), which may explain part of the deviation; they reported results for
${\textit{Re}}_p$
up to approximately 500, while our particles have
${\textit{Re}}_p \sim O(10^4)$
. Cabrera-Booman et al. (Reference Cabrera-Booman, Plihon and Bourgoin2024) also noted that
$C_d$
may have a systematic dependence on the specific gravity
$\varGamma$
that was not captured in their parameterisation, with lighter particles experiencing more drag than expected and heavier particles experiencing less drag. To investigate whether such an effect is present in our data as well, we both colour the data points in figure 3 by
$\varGamma$
and plot
${\textit{Sv}}$
as a function of
$\varGamma$
(coloured by
${\textit{Gr}}_p$
) in figure 4. Although we do observe spread in our data as a function of
$\varGamma$
, the overall trends are not so clear – and if anything appear to suggest the opposite of what Cabrera-Booman et al. (Reference Cabrera-Booman, Plihon and Bourgoin2024) observed, as
${\textit{Sv}}$
is larger for smaller values of
$\varGamma$
than for larger values.
An additional possibility is that the porosity or surface roughness of our particles or mass transfer effects may have an impact on their settling. Emadzadeh & Chiew (Reference Emadzadeh and Chiew2020) reported enhanced drag coefficients for porous particles for
${\textit{Re}}_p \gt 100$
, which they attributed to a combination of increased frictional drag, reduced pressure drag, a modified wake and boundary layer separation. The magnitude of increased drag observed in their experiments is consistent with the reduction in
${\textit{Sv}}$
that we observe for
${\textit{Gr}}_p = 0$
. However, they did not report a general parameterisation of the drag coefficient for porous particles, so we cannot directly extrapolate their results for our situation. Mass transfer also occurs on the surface of our particles, since as the ethanol mixture loaded into the particles flows out it is replaced the background pure water. This mass transfer could lead to a Stefan flow, which may affect the drag as well (Jayawickrama et al. Reference Jayawickrama, Haugen, Babler, Chishty and Umeki2019), although to an extent that is difficult to predict at our high particle Reynolds numbers. A potential option to at least partially account for these complications would be to define
${\textit{Sv}}$
using the measured settling velocity for the
${\textit{Gr}}_p = 0$
case as the velocity scale, so that
${\textit{Sv}}(Gr_p=0)$
would be unity. However, as we change
${\textit{Gr}}_p$
by varying the ethanol concentration loaded into the particles,
${\textit{Ga}}$
also changes – and we expect the drag to be a function of
${\textit{Ga}}$
. Thus, in estimating the expected settling velocity to define
${\textit{Sv}}$
, we choose to use the
$C_d$
model from Cabrera-Booman et al. (Reference Cabrera-Booman, Plihon and Bourgoin2024) and simply compare with the
${\textit{Gr}}_p = 0$
case as a reference point for assessing the relative influence of active buoyancy.

Figure 5. The settling number
${\textit{Sv}}$
plotted as a function of the particle Grashof number
${\textit{Gr}}_p$
. The colour of each data point shows the Schmidt number
${\textit{Sc}}$
, as indicated by the colour bar. The error bars indicate the standard error over 10 trials for each set of control parameters. Note that, although the data for each
${\textit{Gr}}_p$
have the same
${\textit{Sc}}$
,
${\textit{Sc}}$
in our experiments does not vary monotonically with
${\textit{Gr}}_p$
.
Independent of the lower-than-expected
${\textit{Sv}}$
at
${\textit{Gr}}_p=0$
, figure 3 shows that
${\textit{Sv}}$
does indeed vary with
${\textit{Gr}}_p$
, that the settling can be either slower or faster than expected, and that this variation is non-monotonic. For
$0.5 \times 10^7 \lt Gr_p \lt 1.5 \times 10^7$
,
${\textit{Sv}}$
increases slightly relative to the
${\textit{Gr}} = 0$
case. A slight decrease in
${\textit{Sv}}$
occurs around
${\textit{Gr}}_p = 2.5 \times 10^7$
, followed by a further increase in
${\textit{Sv}}$
dependent on
$\varGamma$
near
${\textit{Gr}}_p = 4 \times 10^7$
, which is more pronounced for particles with
$\varGamma$
closer to unity. A potential explanation for the non-monotonic behaviour we observe may be Schmidt number dependence. As described above, the primary way we vary
${\textit{Gr}}_p$
is by changing the concentration of the ethanol solution loaded into the particles, which changes
$\rho _s$
. Here,
${\textit{Gr}}_p$
increases monotonically as this concentration increases. However, the diffusion constant for ethanol–water mixtures does not vary monotonically with concentration (Pratt & Wakeham Reference Pratt and Wakeham1974), meaning that our Schmidt numbers do not vary monotonically with
${\textit{Gr}}_p$
, as shown in figure 5. Given the lack of theory for this inertial, non-Boussinesq regime and the co-variation of
${\textit{Sc}}$
and
${\textit{Gr}}_p$
in our experiments, we cannot at present clarify the distinct dependence of
${\textit{Sv}}$
on
${\textit{Sc}}$
.

Figure 6. The settling number
${\textit{Sv}}$
plotted as a function of the ratio of the particle Grashof number
${\textit{Gr}}_p$
to the square of the Galileo number
${\textit{Ga}}$
. This composite parameter is equivalent to the ratio of the density anomaly due to the scalar and the density anomaly due to the particle (see text). The colour of each data point shows the Schmidt number
${\textit{Sc}}$
, as indicated by the colour bar. The error bars correspond to the standard error over 10 trials for each set of control parameters.
This lack of theory also generally prevents us from exploring any dependence of
${\textit{Sv}}$
on combinations of our four control parameters. Nevertheless, we do consider the variation of
${\textit{Sv}}$
with the composite parameter
${\textit{Gr}}_p / {\textit{Ga}}^2$
, as shown in figure 6, because its limits are physically interpretable. As can readily be seen from their definitions
and can thus be interpreted as the ratio of the density anomaly due to the scalar attached to the particle (which contributes to the active, dynamical part of the particle’s effective buoyancy) and the density anomaly due to the intrinsic particle density (which is purely passive and static). As seen in figure 6, we observe an approximately linear relationship between
${\textit{Sv}}$
and
${\textit{Gr}}_p / {\textit{Ga}}^2$
in our parameter range, with a cross-over from
${\textit{Sv}}\lt 1$
to
${\textit{Sv}}\gt 1$
at
${\textit{Gr}}_p / {\textit{Ga}}^2 \approx 2$
. This result suggests that the enhancement of the settling velocity begins to occur when the active buoyancy dynamics begins to outweigh the passive contribution. However, this result also has other possible interpretations. We note that
${\textit{Gr}}_p / {\textit{Ga}}^2$
can also be written as
$\varLambda / (\varGamma - 1)$
, where
$\varLambda = (\rho _{\!f} - \rho _s) / \rho _{\!f}$
is a parameter that can be used to assess how non-Boussinesq the density perturbations due to the scalar are (Ganguli & Lele Reference Ganguli and Lele2019). As
$\varLambda$
, and therefore
${\textit{Gr}}_p/{\textit{Ga}}^2$
, grows, we might expect the scalar field shed from the particle to begin to affect its momentum directly via a momentum flux rather than only indirectly through a renormalisation of its effective density.

Figure 7. The settling number
${\textit{Sv}}$
plotted as a function of the ratio of the particle Grashof number
${\textit{Gr}}_p$
to the square of the particle Reynolds number
${\textit{Re}}_p$
. For
${\textit{Gr}}_p / {\textit{Re}}_p \ll 1$
, the scalar wake behind the particle is expected to result from free convection, while for
${\textit{Gr}}_p / {\textit{Re}}_p \gg 1$
it is expected to result from forced convection. The case
${\textit{Gr}}_p / {\textit{Re}}_p \approx 1$
corresponds to the complex, mixed-convection case (see text). The colour of each data point shows the Schmidt number
${\textit{Sc}}$
, as indicated by the colour bar. The error bars correspond to the standard error over 10 trials for each set of control parameters.
As a final way to characterise our settling results, in figure 7 we plot
${\textit{Sv}}$
as a function of the mixed-convection parameter
${\textit{Gr}}_p / {\textit{Re}}_p^2$
. As discussed above, this parameter indicates the relative strength of free and forced convection, and thus brings in some of the dynamics. Although the results should be interpreted with some caution, given that the measured settling velocity
$w$
appears in both
${\textit{Sv}}$
and
${\textit{Re}}_p$
, figure 7 suggests that the enhanced settling of our actively buoyant particles is most pronounced when
${\textit{Gr}}_p / {\textit{Re}}_p^2 \sim 1$
, indicating strong and potentially nonlinear coupling between free and forced convection (Ganguli & Lele Reference Ganguli and Lele2020).
3.2. Particle trajectories
In addition to characterising the particle settling velocity, we also investigated the particle paths themselves. As can be seen in figure 2(d), our particles do not follow straight paths; rather, we observe oscillatory behaviour that is nearly confined in a single two-dimensional plane (sometimes with some small out-of-plane motion as well) in all cases. These zig-zag paths are not unexpected, given the particle Reynolds numbers (
${\textit{Re}}_p \sim O(10^{4})$
) with which our particles fall. A significant amount of previous work has characterised the paths of falling objects (Ern et al. Reference Ern, Risso, Fabre and Magnaudet2012) as the particle Reynolds number increases, which can be described as arising from instabilities of the inertial wake shed from the falling particle (Uhlmann & Dušek Reference Uhlmann and Dušek2014; Dušek Reference Dušek2017). The resulting trajectories can be straight, angled, oscillatory or even chaotic. Our actively buoyant particles, however, are more complex: in addition to the typical inertial wake, they also leave behind a density wake. The interaction of this density anomaly with the fluid and particle inertia may potentially affect the wake instabilities that drive the deviation of a settling particle from a straight path (Ern et al. Reference Ern, Risso, Fabre and Magnaudet2012; Uhlmann & Dušek Reference Uhlmann and Dušek2014; Dušek Reference Dušek2017).

Figure 8. (a) Reconstructed trajectories of falling particles for ten experimental trials with
${\textit{Gr}}_p = 4.0 \times 10^7$
,
${\textit{Ga}} = 7.03\times 10^3$
,
$\varGamma = 1.26$
and
${\textit{Sc}} = 950$
. These particles had a diameter of 2.69 cm. Here,
$z$
is aligned with the gravity direction (with the acceleration due to gravity pointing in the
$-\hat {\mathbf{z}}$
direction) and
$x$
is an orthogonal direction. Each trajectory is shifted so that
$x=0$
corresponds to its mean horizontal position and
$z=0$
corresponds to its initial detected location; these origins have no other physical meaning. Additionally, for ease of visualisation, some trajectories have been flipped across the
$x$
axis so that all of the trajectories appear to move initially in the positive
$x$
direction. (b) Time series of
$z$
for each of the ten trials. (c) Time series of
$x$
for each of the ten trials. Each trajectory is plotted in the same colour in (a), (b) and (c). The black dashed lines in (c) show least-squares fits of single-mode sinusoids to the time series of the horizontal position (see text for more details).
In figure 8, we show the trajectories for the set of ten experiments with
${\textit{Gr}}_p = 4.0 \times 10^7$
,
${\textit{Ga}} = 7.03\times 10^3$
,
$\varGamma = 1.26$
and
${\textit{Sc}} = 950$
. Although we only show trajectories for this set of control parameters, they are qualitatively similar for all the other cases we considered. It is immediately clear from figure 8(a) that the particles do not fall along a straight path, but rather oscillate. This oscillation, however, occurs only the direction perpendicular to gravity; as shown in figure 8(b), the vertical position of the particle is linear in time (justifying our choice to extract particle velocities via linear fits rather than numerical differentiation). The oscillations we observe are very nearly single-mode sinusoids. This result is somewhat surprising, as our
${\textit{Ga}}$
and
$\varGamma$
values are in the regime where we would expect chaotic trajectories (Uhlmann & Dušek Reference Uhlmann and Dušek2014; Dušek Reference Dušek2017; Cabrera-Booman et al. Reference Cabrera-Booman, Plihon and Bourgoin2024). It is possible that the shed scalar wake behind the particle, which, like the momentum wake, moves in the direction anti-parallel to gravity, may have a stabilising effect on any path instabilities. We do, however, also observe similar behaviour for the case of
${\textit{Gr}}_p=0$
, suggesting (similar to the smaller-than-expected
${\textit{Sv}}$
for
${\textit{Gr}}_p=0$
discussed above) that other effects such as porosity may also be a factor in controlling the behaviour of the particles.

Figure 9. The Strouhal number
${\textit{St}}$
plotted as a function of the Galileo number
${\textit{Ga}}$
. The colour of each data point shows the particle Grashof number
${\textit{Gr}}_p$
, as indicated by the colour bar. The error bars correspond to the standard error over 10 trials for each set of control parameters.
To extract a frequency
$f$
for the oscillations, we fit the time series of the horizontal particle position
$x$
to functions of the form
$\cos (2\pi \textit{ft})$
; figure 8(c) shows that these fits capture the observed behaviour very accurately. We then non-dimensionalised the frequencies to define a Strouhal number
${\textit{St}} = \textit{fd} / w$
for each set of control parameters. In figure 9, we plot
${\textit{St}}$
as a function of
${\textit{Ga}}$
. Although there is scatter in the data and potentially different behaviour for
${\textit{Ga}}\lt 3000$
, our results are consistent with a constant Strouhal number of approximately
$0.05$
. There is also no apparent trend in the results as a function of
${\textit{Gr}}_p$
. Similar
${\textit{St}}$
values have been reported previously for inert settling spheres (Ern et al. Reference Ern, Risso, Fabre and Magnaudet2012), although again for lower
${\textit{Ga}}$
ranges. These
${\textit{St}}$
values are thus again consistent with the hypothesis that the active buoyancy in our experiments may act to suppress wake instabilities that would occur due to inertia alone.
4. Conclusions
To briefly summarise our results, we conducted a series of experiments to measure the settling behaviour of particles that released a buoyant scalar into the surrounding fluid. Because the density anomaly due to the introduction of this scalar was sourced from the particle and strongly non-Boussinesq, we refer to these particles as being actively buoyant. We found that the settling velocity could be both suppressed or enhanced relative to what would be expected for inert particles of the same size, shape and density. These deviations were not small; in particular, for our experimental case with the highest particle Grashof number and specific gravity closest to unity, we found that the settling velocity was enhanced by nearly 60 %. Given that several of our control parameters were similar to what would be expected for real firebrands, our results suggest that firebrand transport models that do not account for active buoyancy effects risk severe mis-prediction of their likely landing locations.
There is no theory currently available for how the settling velocity would be expected to depend on the governing non-dimensional parameters in the inertial, non-Boussinesq regime we studied. Indeed, our results indicate that all four control parameters (
${\textit{Gr}}_p$
,
${\textit{Ga}}$
,
$\varGamma$
and
${\textit{Sc}}$
) likely play a role in determining the settling velocity. Plotting our results as a function of the composite parameters
${\textit{Gr}}_p / {\textit{Ga}}^2$
and
${\textit{Gr}}_p / {\textit{Re}}_p^2$
, which both show monotonic increase of
${\textit{Sv}}$
, suggests that the buoyancy of the scalar and the nature of the resulting convective flow play key roles in determining the particle settling behaviour. When combined also with the strong effects for
$\varGamma$
close to unity and the highly non-Boussinesq nature of the density fluctuations, we hypothesise that the direct momentum fluxes imparted on the particle by the shed scalar dominate over any simple renormalisation of the effective density of the particle. Future studies of this complex problem using, for example, numerical simulation or experimental techniques such as combined particle image velocimetry and laser-induced fluorescence could probe this hypothesis in more detail.
In contrast to the settling velocity, where the addition of active buoyancy led to significantly different values relative to the case of inert particles, we did not observe significant modifications to the trajectories of our settling particles. As would be expected for the range of
${\textit{Ga}}$
and
$\varGamma$
we investigated, the particles did not follow straight paths; instead, the particles followed nearly sinusoidal paths, with a Strouhal number similar to what would be expected for inert particles. The primary effects of active buoyancy on the trajectory structure appears to be a stabilisation of this sinusoidal regime at higher than expected
${\textit{Re}}_p$
, which may also arise due to the additional wake momentum flux due to the non-Boussinesq density fluctuations. Again, further detailed study of the particle wake and scalar plume shed from the particle would help to illuminate the mechanisms responsible for wake stabilisation.
Finally, we note that, although firebrand transport is our primary motivation for this work, it is far from the only situation where active, dynamically modulated particle buoyancy is an important component of the physics. The thermal activity of particles in solar receivers, for example, has been shown to play a major role in the efficacy of their energy harvesting (Frankel et al. Reference Frankel, Pouransari, Coletti and Mani2016; Zamansky et al. Reference Zamansky, Coletti, Massot and Mani2016; Yang et al. Reference Yang, Wan, Zhou and Dong2022). Similarly, the heat transport in convective flows can be enhanced when they are laden with thermally expandable particles (Hu et al. Reference Hu, Wang, Jia, Zhong and Zhang2021). Large-scale convection can even be generated and maintained by thermally active particles (Agasthya et al. Reference Agasthya, Bartel, Biferale, Ehrhardt and Toschi2022). And in addition to these engineering examples, the active buoyancy modulation of biological ‘particles’ is essential for many marine organisms (Denton Reference Denton1971; Reynolds et al. Reference Reynolds, Oliver and Walsby1987; Gemmell et al. Reference Gemmell, Oh, Buskey and Villareal2016; Hamel et al. Reference Hamel, Sun, Gianasi, Montgomery, Kenchington, Burel, Rowe, Winger and Mercier2019; Arrieta et al. Reference Arrieta, Jeanneret, Roig and Tuval2020). Our results here indicate that the behaviour of actively buoyant particles is far from simple even for the case of settling in quiescent fluid, suggesting that there is a great deal of opportunity for further study in this area, particularly as the flow field becomes more complex.
Funding
This work was supported in part by the U.S. National Science Foundation under grant number CBET-2529151.
Declaration of interests
The authors report no conflict of interest.























































