1. Introduction
Convection driven by a temperature difference is omnipresent in nature, for example in oceans (Marshall & Schott Reference Marshall and Schott1999), in the atmosphere (Hartmann, Moy & Fu Reference Hartmann, Moy and Fu2001), inside the Earth (Cardin & Olson Reference Cardin and Olson1994; Jones Reference Jones2011) and in the outer layer of the Sun (Cattaneo, Emonet & Weiss Reference Cattaneo, Emonet and Weiss2003). It is also relevant to numerous technological applications, ranging from the cooling of electronic devices to the ventilation in buildings (Incropera Reference Incropera1999; Linden Reference Linden1999). In many cases, thermal convection is coupled to other physical processes, such as phase changes (Kim Reference Kim2009), transport and dispersion of inertial particles (Shaw Reference Shaw2003; Ackerman et al. Reference Ackerman, Kirkpatrick, Stevens and Toon2004), leading to complex and interesting dynamics.
Rayleigh–Bénard convection (RBC) has been extensively studied as a simplified paradigm for thermal convection (Castaing et al. Reference Castaing, Gunaratne, Heslot, Kadanoff, Libchaber, Thomae, Wu, Zaleski and Zanetti1989; Ahlers, Grossmann & Lohse Reference Ahlers, Grossmann and Lohse2009; Chillà & Schumacher Reference Chillà and Schumacher2012; Xia Reference Xia2013). Buoyancy-driven flow in RBC is induced by heating a fluid layer from below and cooling it from above. The system is characterized by two control parameters: the Rayleigh number 
$\textit {Ra}$ and the Prandtl number 
$\textit {Pr}$, defined as
\begin{equation} \textit{Ra} = \alpha gH^3 \Delta T/\nu \kappa \quad \text{and} \quad \textit{Pr} = \nu/\kappa, \end{equation}
respectively. Here, 
$\Delta T$ is the temperature difference applied across the fluid layer, 
$H$ is the height of the fluid layer and 
$g$ is the gravitational acceleration. The parameters of the fluid, 
$\alpha$, 
$\nu$ and 
$\kappa$, are the thermal expansion coefficient, kinematic viscosity and thermal diffusivity, respectively. The heat transport by the convective flow is characterized by the dimensionless Nusselt number,
\begin{equation} \textit{Nu} = QH/\chi \Delta T, \end{equation}
where 
$Q$ is the heat flux (power per unit area) through the convection cell and 
$\chi$ is the thermal conductivity of the fluid. The dependence of 
$\textit {Nu}$ on 
$\textit {Ra}$ and 
$\textit {Pr}$, 
$\textit {Nu} = f(\textit {Ra}, \textit {Pr})$, has been the focus of numerous experimental and theoretical studies (for review see e.g. Ahlers et al. Reference Ahlers, Grossmann and Lohse2009).
An important question in thermal convection is how to modify the fluid structures, change the heating-passing mechanisms and perhaps enhance the heat transport efficiency. Various approaches have been proposed and tested, such as creating roughness on the top and bottom plates (Shen, Tong & Xia Reference Shen, Tong and Xia1996; Du & Tong Reference Du and Tong1998), exploiting multiphase working fluid (Zhong, Funfschilling & Ahlers Reference Zhong, Funfschilling and Ahlers2009; Lakkaraju et al. Reference Lakkaraju, Stevens, Oresta, Verzicco and Lohse2013; Dabiri & Tryggvason Reference Dabiri and Tryggvason2015; Guzman et al. Reference Guzman, Xie, Chen, Rivas, Sun, Lohse and Ahlers2016; Wang, Mathai & Sun Reference Wang, Mathai and Sun2019; Liu et al. Reference Liu, Chong, Ng, Verzicco and Lohse2021), imposing spatial confinement (Huang et al. Reference Huang, Kaczorowski, Ni and Xia2013), partitioning the convection cell (Bao et al. Reference Bao, Chen, Liu, She, Zhang and Zhou2015) and manipulating the coherent structure (Chong et al. Reference Chong, Yang, Huang, Zhong, Stevens, Verzicco, Lohse and Xia2017).
Another potential approach to tune heat transport in RBC is to couple thermal convection with inertial particles. In particle-laden flows, the particle inertia, i.e. lagging in responses to the changes in the fluid flows, can arise from the density mismatch between the particle and the local fluid or from the finite size of the particle (Ouellette, O'Malley & Gollub Reference Ouellette, O'Malley and Gollub2008; Cartwright et al. Reference Cartwright, Feudel, Károlyi, de Moura, Piro and Tél2010). Previous studies have focused on micron- or millimetre-scale particles. Particles may settle out of the fluid when there is a mismatch between the densities of the particles and the background fluid. A turnover of the system dynamics and resuspension of settling particles have been identified beyond a critical particle concentration (Koyaguchi et al. Reference Koyaguchi, Hallworth, Huppert, Sparks and Stephen1990). The settling particles may form porous layers that reduce the heat transport (Joshi et al. Reference Joshi, Rajaei, Kunnen and Clercx2016). The distribution of thermally expandable point-like particles in RBC has also been studied numerically in the ‘soft turbulent’ regime, with the assumption that the presence of the particles does not modify the flow (Alards et al. Reference Alards, Kunnen, Clercx and Toschi2019). However, to the best of our knowledge, no enhancement of 
$\textit {Nu}$ has been reported in turbulent RBC with solid inertial particles. When the particle size is much larger than the Kolmogorov length scale (Kolmogorov Reference Kolmogorov1941), the effect of the motion and geometry of particles on fluid flows must be considered, and whether 
$\textit {Nu}$ can be enhanced is still an open question.
The study of thermal convection with large inertial particles has been challenging both numerically and experimentally. In numerical simulations, one has to properly resolve the fluid–particle coupling, as well as particle–particle and particle–wall collisions (Balachandar & Eaton Reference Balachandar and Eaton2010; Mathai, Lohse & Sun Reference Mathai, Lohse and Sun2020). In experiments, as particle size increases, gravitational settling due to density mismatch becomes more pronounced. Therefore, uniformity of the particle density and a close density match to the background fluid are essential to ensure and maintain stable suspensions. In contrast to infinitesimal passive tracer particles commonly used to map the velocity field, inertial particles may depart from the local flow and show various behaviours, such as mixing, separation and aggregation (Fung & Vassilicos Reference Fung and Vassilicos2003; Saw et al. Reference Saw, Shaw, Ayyalasomayajula, Chuang and Gylfason2008; Sudharsan, Brunton & Riley Reference Sudharsan, Brunton and Riley2016).
Here, we experimentally study the effect of inertial rod-like particles suspended in convecting flows on the overall heat transport. The particle lengths are centimetre-scale, much larger than the Kolmogorov length scale and the thickness of the thermal boundary layers (BLs). The particles are made of polydimethylsiloxane (PDMS). Besides their large size, the particles’ inertia is further enhanced by their large thermal expansion coefficient (see table 1): the particles are lighter than the fluid in hot regions and heavier than the fluid in cold regions. For the first time, we observe an increase in 
$\textit {Nu}$ in turbulent RBC with suspensions of solid inertial particles due to the interactions between particles and the background flows, especially within the BLs.
Parameters of the working fluid and PDMS particles. Densities and the particle length are measured in the experiments. Other parameters are from empirical equations or interpolated using known experimental values (Mark Reference Mark1999; Cheng Reference Cheng2008; Volk & Kähler Reference Volk and Kähler2018).
2. Experimental procedures
Our experiments have been conducted in a cubic cell of size 
$20\times 20\times 19.6 \ \textrm {cm}^{3}$ with the cell height 
$H=19.6$ cm. The sidewalls consist of four glass plates of thickness 0.5 cm. The top cooling plate and the bottom heating plate are made of surface-anodized aluminium. The top plate is cooled by passing temperature-regulated water through its internal grooves. The bottom plate is heated by a film heater at constant power. The temperatures of the top plate and bottom plate are monitored by thermistors embedded within them roughly 1 mm away from the fluid contact surfaces, which leads to short thermal response time and allows the thermistors to perceive the temperature fluctuations in the BLs (Verzicco Reference Verzicco2004). The bulk temperature 
$T_c$ is measured by two thermistors extended into the centre of the convection cell.
The suspended particles are made of PDMS, which has been widely used in microfluidics (for review, see Sackmann, Fulton & Beebe Reference Sackmann, Fulton and Beebe2014) due to its advantageous properties, such as chemical stability and mechanical flexibility after curing. However, the properties relevant to our purposes are different. First, its density after curing is within about 3 % of that of water. A stable suspension of particles, which is crucial to our experiments, is achieved by tuning the background fluid: a dilution of glycerol in water. The density of the working fluid is carefully adjusted so that the particles are neutrally buoyant at 
$27\,^\circ \textrm {C}$: particles can suspend for hours even in quiescent fluid. The working fluid is then used for all the subsequent experiments, which are conducted in the range 
$2\,^{\circ }\mathrm {C}\lesssim \Delta T \lesssim 32\,^{\circ }\mathrm {C}$ and 
$3\times 10^{8}\lesssim \textit {Ra} \lesssim 4\times 10^{9}$. The 
$\textit {Ra}$ range is about 2–3 orders of magnitude larger than the 
$\textit {Ra}$ values in Alards et al. (Reference Alards, Kunnen, Clercx and Toschi2019), and about one order of magnitude smaller than that in Wang et al. (Reference Wang, Mathai and Sun2019) and Joshi et al. (Reference Joshi, Rajaei, Kunnen and Clercx2016). More importantly, PDMS has a large thermal expansion coefficient that is more than three times that of the working fluid (table 1), making the particle thermally responsive and more active than the background fluid.
All PDMS particles are handmade. After many trials, we found that a rod-like geometry yields high precision and high quality. To make these particles, the PDMS fluid and the curing agent (SYLGARD 184 Silicone Elastomer Kit, Dow Corning) are first mixed at a 10 : 1 volumetric ratio and degassed for 30 min by a vacuum pump. The mixture is then poured into moulding tubes of desired diameters, and cured at room temperature for 48 h. The cured PDMS rods are then cut into desired lengths. Making spherical particles (a simpler geometry than rod) using a similar procedure turns out to be impractical, especially when making thousands of identical particles.
Figure 1 presents photos of the suspended particles in turbulent convection with particle volume fraction 
$\phi = 1\,\%$, where 
$\phi$ is defined as the ratio of the total volume of particles to the volume of the convection cell. The heating power and cooling temperature are pre-adjusted to produce a bulk temperature of 
$T_{{c}} = 27 \pm 0.10\,^\circ \textrm {C}$ for no-particle convection experiments. For experiments with particles, two control modes are used. The first one is constant-
$T_{c}$ mode: for different particle volume fractions 
$\phi$, the heat flux 
$Q$ is regulated through a feedback control, and once the system reaches dynamic equilibrium with 
$T_{{c}} = 27\pm 0.15\,^{\circ }\textrm {C}$, 
$Q$ is fixed and no longer changes over time. The second mode is constant-
$Q$, in which the same 
$Q$ is used for different values of 
$\phi$. When measuring 
$\textit {Nu}$ in response to the effects of added particles, constant-
$T_{c}$ mode is used to rigorously maintain the bottom–top symmetry of the particle distributions in the convection cell. Constant-
$Q$ mode is mainly used to measure temperature fluctuations.
The suspended PDMS particles of length 2.5 cm and diameter 0.45 cm with volume fraction 
$\phi = 1\,\%$ interact with the background turbulent flows in RBC at 
$\textit {Ra} \approx 2.1\times 10^{9}$. (a) A 1/100 s exposure photo shows how particles are randomly distributed in the convection cell, and reach the top and bottom plates occasionally. (b) A 1 s long exposure photo shows different degrees of motion blurs that reflect the turbulent background flows.
3. Results and discussions
In figure 2(a), we show the compensated Nusselt number, 
$\textit {Nu}/\textit {Ra}^{0.314}$, as a function of 
$\textit {Ra}$ in a semi-log scale. The measured values of Nusselt number 
$\textit {Nu}_0$ from no-particle experiments follow a power law with 
$\textit {Ra}$, and the scaling exponent agrees with previous experiments (see, e.g. Niemela et al. Reference Niemela, Skrbek, Sreenivasan and Donnelly2000). We then measure 
$\textit {Nu}$ with the suspensions of particles of different diameters and compute the relative changes of 
$\textit {Nu}$, 
$\Delta \textit {Nu}/\textit {Nu}_0 = (\textit {Nu}-\textit {Nu}_0)/\textit {Nu}_0$. As shown in figures 2(a) and 2(b), 
$\Delta \textit {Nu}/\textit {Nu}_0$ depends strongly on 
$\textit {Ra}$, and an optimal value of 
$\textit {Ra}$ exists, around which the 
$\textit {Nu}$ enhancement is most obvious. For particles of diameter 0.30 cm, the maximum increase of 
$\textit {Nu}$ is achieved around 
$\textit {Ra}\approx 2\times 10^{9}$; for particles of a larger diameter 0.45 cm, the maximum is located at a slightly lower 
$\textit {Ra}$. As 
$\textit {Ra}$ becomes larger than the optimal values, 
$\Delta \textit {Nu}/\textit {Nu}_0$ decreases slowly. However, lowering the values of 
$\textit {Ra}$, 
$\Delta \textit {Nu}/\textit {Nu}_0$ decreases steeply and even turns negative at sufficiently small 
$\textit {Ra}$. Figure 2(c) shows that 
$\Delta \textit {Nu}(\phi )$ increases faster as 
$\phi$ is increased from 0 but slows down for larger values of 
$\phi$. It is evident that the addition of the passive inertial particles increases the heat-transport efficiency in RBC. For higher values of 
$\phi$, the particles are likely to get trapped at the corners of the convection cell and form porous layers covering the top or bottom plate, making the mixture no longer a mobile suspension.
With constant-
$T_{{c}}$ mode, enhancement of heat transport with PDMS particles. (a) 
$\textit {Nu}/\textit {Ra}^{0.314}$ as a function of 
$\textit {Ra}$ in semi-log scale. Grey filled circles show data measured without particles; other symbols represent measurements with particles of different diameter 
$d$. Error bars represent uncertainties due to temperature fluctuations and heat loss to the ambient air. (b) Data from (a), plotted as relative increase of 
$\textit {Nu}$, 
$\Delta \textit {Nu}/\textit {Nu}_0$ (left axis), as a function of 
$\textit {Ra}$. Upper 
$x$ axis shows 
$\Delta T$ without particles. (c) 
$\Delta \textit {Nu}/\textit {Nu}_0$ as a function of 
$\phi$ at 
$\textit {Ra} \approx 3.4\times 10^{9}$ for 
$d=0.45$ cm and 
$L=2.50$ cm (red circles). Blue triangles show data for 
$d=0.30$ cm and 
$L=2.50$ cm at 
$\phi =5.0\,\%$; green squares show data for 
$d=0.30$ cm and 
$L=1.25$ cm at 
$\phi =2.5\,\%$ and 
$5.0\,\%$. Particle aspect ratio plays a less important role for 
$Nu$ than 
$\phi$, as can be seen by comparing rods of different geometries.
The Stokes number 
$St$ of the rod-like particles can be estimated as 
$St \sim \bar {\rho }d^2 H^{-1} h(\lambda ) v_f/\nu$, where 
$\bar {\rho }$ is the density ratio of particle to fluid, 
$\bar {\rho } = \rho _p/\rho _f$, 
$v_f$ is the characteristic velocity of the large-scale circulation and 
$h(\lambda )$ is a function of the particle aspect ratio 
$\lambda$ with 
$\lambda =L/d$ (Voth & Soldati Reference Voth and Soldati2017). We estimate 
$v_f$ as 
$v_f\approx (0.2\nu /H) \sqrt {\textit {Ra}/\textit {Pr}}$, where the prefactor 0.2 accounts for the correction at large Pr (Silano, Sreenivasan & Verzicco Reference Silano, Sreenivasan and Verzicco2010). For the 
$Ra$ range in our experiments, 
$v_f \approx 0.58\text {--}2.6\ \textrm {cm}\ \textrm {s}^{-1}$, and 
$St \approx 0.04\text {--}0.16$ for particles of 
$L = 2.5$ cm and 
$d=0.3$ cm, which is orders of magnitude larger than the values reported in previous experiments and simulations (Ouellette et al. Reference Ouellette, O'Malley and Gollub2008; Lopez & Guazzelli Reference Lopez and Guazzelli2017), indicating a finite inertia effect afforded by the large particle size and deviations of particle trajectories from fluid flows. To see the effect of particle aspect ratio, we manually cut each particle of 
$d=0.30$ cm in half. As shown in figure 2(c), for the same volume fraction at 
$\phi =5\,\%$, 
$\Delta Nu$ of particles of 
$L=1.25$ cm (right green square) is nearly the same as that of particles of double length (blue triangle); for the same number of particles, the 
$Nu$ of particles of 
$L=1.25$ cm at 
$\phi =2.5\,\%$ (left green square) is slightly smaller than that of particles of double length at 
$\phi =5.0\,\%$ (blue triangle). Given the current large particle sizes at centimetre scale, our data (within experimental accuracy) show that the particle aspect ratio has only a small effect on 
$Nu$. We speculate that the aspect ratio is more important for thinner particles of 
$St \ll 1$, and more experiments are needed in the future to further explore the effect of elongated particle shapes.
Figure 3(a) shows the probability distribution function of the bottom plate temperature 
$T_{{b}}$ around 
$\textit {Ra}$ values where particles enhance 
$\textit {Nu}$. With constant-
$Q$ mode, i.e. the heating power of the bottom plate is fixed to be the same as that used in the no-particle experiment, 
$T_{{b}}$ decreases as 
$\phi$ is increased, indicating a decrease of 
$\Delta T$ and thus an increase of 
$\textit {Nu}$. This is in line with the increase of 
$\textit {Nu}$ shown in figure 2(c) under constant-
$T_{{c}}$ mode, in which 
$\Delta T$ remains relatively fixed within 
$\pm 0.4\,^{\circ }\textrm {C}$ and the heating power of the bottom plate increases as 
$\phi$ is increased. Meanwhile, as shown in figure 3(b), the standard deviation 
$\sigma _T$, which characterizes the magnitude of the temperature fluctuations, decreases as 
$\phi$ is increased for both top and bottom plates. For different values of 
$\phi$, the normalized bottom temperatures, 
$(T_{{b}}-\langle T_{{b}} \rangle )/\sigma _{T_{{b}}}$, collapse onto each other and closely resemble a Gaussian distribution (figure 3c). It is well accepted that thermal BLs impose the most resistance on heat transport and dominantly determine 
$\textit {Nu}$ (see e.g. Castaing et al. Reference Castaing, Gunaratne, Heslot, Kadanoff, Libchaber, Thomae, Wu, Zaleski and Zanetti1989). Much of the heat flux is carried by the thermal plumes, which are emitted from the thermal BLs and contrast with the background flow in temperature and momentum (Zocchi, Moses & Libchaber Reference Zocchi, Moses and Libchaber1990; Kadanoff Reference Kadanoff2001). The decrease in 
$\sigma _T$ observed in figure 3(b) suggests that particles act as mixers of the flow and temperature fields near the BLs. The large temperature fluctuations caused by the emission and arrival of thermal plumes are partially smoothed out by the motion of the particles.
With constant-
$Q$ mode, 
$T_{{t}}$ remains within 
$10.5 \pm 0.1\,^{\circ }\textrm {C}$. (a) Probability distribution function (p.d.f.) of 
$T_{b}$ for 
$\textit {Ra}\approx 3.9\text {--}4.2\times 10^9$ for no-particle experiments (
$\phi =0$, grey solid circles) and experiments with particles of diameter 0.30 cm and length 2.50 cm (red symbols) at different values of 
$\phi$; blue symbols show particles of diameter 0.30 cm and length 1.25 cm. (b) Standard deviation 
$\sigma _T$ of the bottom temperature 
$T_{b}$ (red circles) and top temperature 
$T_{t}$ (blue triangles) as a function of 
$\phi$. (c) P.d.f. of normalized 
$T_{b}$ for different values of 
$\phi$ in semi-log scale. The green curve shows a standard Gaussian distribution. Statistics are obtained from time series of 2400 s at a sampling rate of 0.5 Hz.
The path and orientation of each particle seem to be complex, as depicted in figure 1; more so are the interactions between particles and BLs. Several forms of interaction between particles and BLs are observed (see figure 4 and the supplementary movie available at https://doi.org/10.1017/jfm.2021.809). Particles may collide with the top and bottom plates and then return to the bulk (figure 4a); they may also slide along the top and bottom plates, within or close to the BLs (figure 4b,c).
Three typical examples of how particles of length 2.5 cm and diameter 0.45 cm interact with bottom BL at 
$\textit {Ra} \approx 2.1\times 10^{9}$. Successive snapshots (white cylinders) showing particles moving from left to right, as extracted from a video with a time interval of 1.0 s for (a,c) and 0.7 s for (b). The background pictures are taken from the last frame. (a) A particle briefly collides twice with the BL and returns to the bulk. (b) A particle slides along the BL. (c) A particle hovers above the BL.
The estimated thickness of thermal BL, 
$\delta = H/2 \textit {Nu} \approx 0.1 \text {--} 0.2\ \text {cm}$, is close to the particle radius, and the particles are exposed to large temperature variations when they get close to the thermal BLs. Due to the large thermal expansion coefficient of PDMS, the particles become lighter than the local fluid when in contact with the hot BL, and heavier than the local fluid when in contact with the cold BL. At the largest value of 
$Ra$, the maximum and minimum density ratio 
$\bar {\rho }$ are about 
$1.01$ and 
$0.99$, respectively, taking into account the thermal expansion coefficients of both the particles and the fluid. Although there is only 
$1\,\%$ density variation, the particles may significantly deviate from fluid flows due to their large sizes. We estimate the buoyancy-based settling/rising velocity in quiescent fluid as 
$v_g \approx [|\bar {\rho }-1|\times g L]^{1/2}$ (Auguste & Magnaudet Reference Auguste and Magnaudet2018). When the particle length is 
$L = 2.50$ cm, 
$v_g \approx 1.5\text {--}5.0\ \textrm {cm}\ \textrm {s}^{-1}$ for the 
$Ra$ range in our experiments. If the particle diameter is taken as the characteristic scale, 
$v_g \approx 0.55\text {--}1.7\ \textrm {cm}\ \textrm {s}^{-1}$. Both estimates are comparable to 
$v_f$, and therefore, the particles become more ‘active’ than the local fluid within the BL, making them more efficient in carrying fluid and heat into the bulk. As a result, the standard deviation of the bulk temperature 
$\sigma _{T_{{c}}}$ increases with 
$\phi$ under both constant-
$T_{{c}}$ and constant-
$Q$ modes (figure 5a). As shown in figure 5(b), the normalized 
$T_{{c}}$ displays exponential distributions with relatively large deviations from the means (see, e.g. Castaing et al. Reference Castaing, Gunaratne, Heslot, Kadanoff, Libchaber, Thomae, Wu, Zaleski and Zanetti1989). The enhancement of fluctuations in incident turbulence has been observed due to the agitation of bubbles (Alméras et al. Reference Alméras, Mathai, Lohse and Sun2017). In RBC with vapour-bubble nucleation, the temperature fluctuations are decreased, probably owing to vapour bubbles absorbing and releasing heat (Lakkaraju, Toschi & Lohse Reference Lakkaraju, Toschi and Lohse2014; Guzman et al. Reference Guzman, Xie, Chen, Rivas, Sun, Lohse and Ahlers2016). This effect is negligible in our experiments since no latent heat (no phase transition) is involved and the specific heat capacity of PDMS is much smaller than that of the background fluid. It has also been observed that the bulk temperature fluctuations are higher in convection cells with rough surfaces than in those with smooth surfaces, as is the heat transport efficiency (Du & Tong Reference Du and Tong2001).
Statistics of the bulk temperature 
$T_{{c}}$ obtained from time series of 2400 s at a sampling rate of 0.5 Hz. (a) 
$\sigma _{T_{c}}$ as a function of 
$\phi$ with constant-
$T_{{c}}$ (red circles) and constant-
$Q$ modes (blue triangles). (b) Normalized distributions of 
$T_\mathrm {c}$ for different values of 
$\phi$ in semi-log scale with the best-fit exponential distribution (green curve), which has better agreement with the experimental data than the best-fit normal distribution (blue curve).
To make particles more ‘active’ than the background fluid, it is necessary to have large 
$\Delta T$ and sufficient time for particles to thermally adapt to 
$\Delta T$. The latter is determined by two time scales. The first is the thermal response time of the particles, 
$\tau _p \sim d^{2}/4\kappa _p$, where 
$\kappa _p$ is the thermal diffusivity of PDMS. The second is the characteristic time for which the particles remain close to the BLs, 
$\tau _f$, which decreases as 
$\textit {Ra}$ is increased. We estimate 
$\tau _f$ as 
$\tau _f \approx H/v_f$. For the 
$\textit {Ra}$ range in our experiments, the ratio of the two characteristic times is 
$\bar {\tau }=\tau _p/\tau _f\approx 0.5\text {--}2.4$. At large 
$\textit {Ra}$, 
$\bar {\tau }>1$, although 
$\Delta T$ is large, 
$\tau _f<\tau _p$ and the particles have insufficient time to warm up or cool down in response to the high/low temperature of the bottom/top BL before returning to the bulk. At small 
$\textit {Ra}$, 
$\bar {\tau }<1$ and 
$\tau _f>\tau _p$, but 
$\Delta T$ is small, leading to relatively small variations in the densities of the particles, i.e. the particles are not ‘active’ enough even after a full thermal relaxation to the temperatures of the BLs. Indeed, we observe at small 
$\textit {Ra}$ that the particles occasionally accumulate and stay near the bottom and top plates for relatively long times, blocking the thermal plumes and thus decreasing 
$\textit {Nu}$. At intermediate 
$\textit {Ra}$, where 
$\Delta T$ is relatively large and 
$\tau _f \approx \tau _p$, 
$\Delta \textit {Nu}/\textit {Nu}_0$ reaches the largest value. In this regime, a particle gains enough density difference and returns to the bulk dragging some fluid with it. It may be able to enter the BL of the other side, switch to another temperature and density extreme, and repeat this process. The intrusion of the particles into the thermal BLs is probably the mechanism most responsible for enhancing the overall heat-transport efficiency, reflected in the increase of 
$\textit {Nu}$.
To directly verify that it is the particle–BL interaction that leads to the enhancement of heat transport efficiency, as illustrated in figure 6(a), we mount two polythene mesh sheets outside the thermal BLs, 1 cm away from the top and bottom plates, respectively. Each mesh has about 
$56\,\%$ opening area, and its thickness or wire diameter is 0.1 cm. The mesh opening has a diamond shape of size 
$0.3 \times 0.3$ cm, smaller than the 
$d=0.45$ cm particles, and thus prevents them from passing through. The mesh sheets have little effect on the heat transport at large 
$\textit {Ra}$ (figure 6(b), blue triangles). We then add particles, which are kept within the two meshes and circulating in the bulk without interacting with the BLs. As evidenced in figure 6(b) (red squares), the difference between the measured 
$\textit {Nu}$ values with and without particles is negligible. Indeed, particle–BL interaction is crucial in the enhancement of heat transport.
(a) Schematic of experiment with particles and meshes, which prevent particles from getting close to the BLs. (b) 
$\textit {Nu}/\textit {Ra}^{0.314}$ as a function of 
$\textit {Ra}$ without mesh or particles (grey filled circles), with mesh and no particles (blue triangles), and with mesh and particles (red solid squares) of diameter 0.45 cm and length 2.50 cm with 
$\phi = 6\,\%$. Enhancement of 
$\textit {Nu}$ is no longer observed.
4. Conclusion
We have experimentally studied a cubic RBC system with suspensions of inertial and expandable rod-like PDMS particles. The heat transport may be enhanced or reduced depending on the Rayleigh number. Our measurements demonstrate that the particle's large thermal expansibility has strong effects on the heat transport. Our results may shed light on a new approach to control the heat transport in thermal convection with suspensions of thermally responsive particles without modification to the classical convection system. Our results may also have implications for the transport and mixing of particles in complex flows and in confined environments.
Engineered materials with voids filled with air, such as silicon foams, may further enhance the heat transport, since the thermal expansion coefficient of air is more than ten times that of water. Particles made of these materials are expected to gain large density variations across temperature differences in thermal convection, but density match may be difficult to achieve. In this work, we have not quantified the three-dimensional dynamics of the particles, which could be constructed using two-dimensional images taken by multiple cameras (Parsa et al. Reference Parsa, Calzavarini, Toschi and Voth2012). In convection cells with extreme aspect ratios where multiple rolls exist, the dynamics is expected to be more complex since particles may move across different convection rolls (Solomon & Gollub Reference Solomon and Gollub1988). We speculate that in a confined rectangular cell (with a narrower width), the interactions between particles and BLs may become more frequent, and 
$\textit {Nu}$ may be further enhanced.
Supplementary movie
Supplementary movie is available at https://doi.org/10.1017/jfm.2021.809.
Acknowledgements
We thank M. Huang for helpful discussions and C. Sun for useful references. We also thank the anonymous reviewers for their insightful criticisms and constructive suggestions.
Funding
S.Y.H. and K.Z.W. gratefully acknowledge support from the MacCracken Fellowship provided by New York University. J.Z. acknowledges support from NYU Shanghai and partial support by Tamkeen under the NYU Abu Dhabi Research Institute grant CG002.
Declaration of interests
The authors report no conflict of interest.
