Fully determined three-dimensional velocity �eld in a divergence-free convection experiment with rigid boundary conditions

Various velocity measurements techniques for 1 ﬂuid dynamics experiments have been developed over 2 the last four decades. These include Particle Image Ve-3 locimetry (PIV) and Particle Tracking Velocimetry (PTV), 4 both of which rely on imaging the ﬂow with passive 5 particles. While these techniques now allow for an ac-6 curate determination of velocity ﬁelds at high particle 7 densities, they still suﬀer from signiﬁcant inaccuracies 8 when it comes to measurements near boundaries, due 9 to smaller particle concentrations and imaging diﬃcul-10 ties in these regions. To improve such measurements, 11 we developed a correction technique for measured ve-12 locity ﬁelds. Using a PTV approach, we reconstructed 13 the 3D velocity ﬁeld from an experimental setup de-14 signed to study Rayleigh-B´enard convection in a tank. 15 We corrected the obtained velocity ﬁeld to ensure its 16 accurate, rigid boundary conditions corresponding to 17 the tank walls, while keeping the ﬂow divergence-free. 18 We then assessed the quality of the corrected velocity 19 ﬁelds by comparing them with the outputs of 3D numer-20 ical simulations that reproduced the exact experimen-21 tal conditions. The experimentally obtained divergence-22 free velocity ﬁelds with exact boundary conditions com-23 pare well with the results of the numerical simulations, 24 which validates our approach. This correction technique 25 could be useful for a wide range of ﬂuid dynamics ex-26 periments where accurate measurements near the walls 27 are needed. These include the advection of passive trac-28 ers, where the boundary conditions need to be exact, 29 or turbulence measurements, where the ﬂow gradient 30 at the walls must be precisely determined.


Introduction
In analog fluid dynamics experiments, the velocity field is by definition one of the key measurements.Consequently, over the years, a variety of velocity measurement techniques have been developed.These include measurements through correlation time of two thermal probes (Chavanne et al, 1997), laser doppler velocimetry (LDV) (Ashkenazi and Steinberg, 1999;Verdoold et al, 2006), Particle Image Velocimetry (Liot et al, 2017;Valori et al, 2017), shadowgraph techniques (Méthivier et al, 2022) and Particle Tracking Velocimetry (Ni et al, 2012;Liot et al, 2016).The thermal probes is the most limited as it only provides information on local velocity.While LDV techniques offer a high resolution as well as 3D measurements (Lowe and Simpson, 2009), they do not allow for large measurement volumes.PIV techniques rely on images of passive particles and on the search for a correlation peak using sliding search windows between two consecutive picture frames.They can be either 2D or 3D.Because of the correlation process, PIV techniques tend to smooth out details of the flow field that are smaller than the search window size.The shadowgraph technique uses the same correlation technique as the PIV techniques, but relies on thermal plumes instead of particles, implying that it is not limited by the need for particles as passive tracers of the flow.However, shadowgraph visualization averages the velocity field across the depth of the tank.Instead, Particle Tracking Velocimetry more efficiently determines the entire velocity field, while keeping a high resolution.3D PTV triangulates imaged particles at a given time step and finds matching particles at the following time steps (Scarano, 2012;Schröder and Schanz, 2023).Because it does not rely on correlations of groups of particles, this technique allows for a higher resolution than traditional PIV techniques.Recently, a particular PTV technique called Shake-The-Box (STB) was introduced (Schanz et al, 2013(Schanz et al, , 2016)).It relies on particle tracks prediction and a matching step to reduce computational time and to remove additional artifacts (Elsinga et al, 2011;Schanz et al, 2013Schanz et al, , 2016)).
The velocity field can then be retrieved by interpolating the tracks onto a fine grid under a zero-divergence constraint.This technique has been used for various measurements, such as high speed jets (Manovski et al, 2021;Sellappan et al, 2020) and stirring experiments (Fitschen et al, 2021;Hofmann et al, 2022).For both PIV and PTV, low particle densities near the walls result in an increase of the uncertainties (Kähler et al, 2012).This proves to be a general and recurrent limitation for various applications.For example, in turbulent boundary layer theory, the friction velocity is defined through the wall shear stress, which is itself defined with the flow gradient near the walls.In this case, an accurate determination of the flow near the walls is mandatory because the presence of errors will make it difficult to compare between different experiments (Bitter et al, 2011;Marusic et al, 2010).While the STB method has already proven to be efficient for high velocity gradients measurements (Schröder et al, 2015), the velocity measurements at rigid boundaries are still not rigorously equal to zero.Another example is the study of Large Scale Circulations (LSCs) in thermal convection that are subject to reversals as well as torsional oscillations (Brown and Ahlers, 2006;Xi and Xia, 2007;Funfschilling et al, 2008).Since the reversals could be caused by boundary layer instabilities (Sreenivasan et al, 2002), it is important to determine accurately the velocity in the boundary layers and near the walls.
Velocity field measurements can also be used to estimate mixing times and stretching fields by numerically seeding the measured flow, as done in Voth et al (2002,2003).While this work provided great insight into mixing processes, it only focused on a central area of the flow.Inaccurate determination of the flow and its properties at the domain boundaries could lead to spurious accumulation of tracers along the walls.For a complete determination of the stretching field in an entire flow, the flow then needs to be accurately determined at the boundaries, and its properties such as the divergence must also be accurate.To circumvent these common limitations of PTV approaches, we present a new correction method ensuring the accurate reconstruction of the velocity field near the walls.We applied and validated our approach using STB velocity measurements in Rayleigh-Bénard convection experiments conducted in a tank.The following Section 2 of this paper first details the experimental setup we used, including the velocity measurement approach, and presents the novel 120 velocity correction technique we developed.Section 3 121 describes the results obtained for the convection ex-122 periments and shows how they compare with results of 123 numerical simulations conducted under the same condi-124 tions.Section 4 discusses the discrepancies between the 125 analog and numerical experiments and the strengths 126 and limitations of our proposed velocity correction and 127 of the convection experiment setup, prior to concluding 128 remarks in Section 5.The experimental setup displayed in Figure 1 aims at generating Rayleigh-Bénard convective motions in controlled conditions, and at measuring the velocity field in the entire domain.The dynamics of this system is controlled by two dimensionless parameters.The first one is the Rayleigh number: and second the Prandtl number: where α is the thermal expansivity, ∆T is the imposed 133 temperature difference between the two horizontal sur-134 faces, g is the gravity, H is the height of the cell, κ is the 135 thermal diffusivity and ν the kinematic viscosity.This 136 configuration yields two outputs, the dimensionless heat 137 flux, the Nusselt number, N u = φ/(λ∆T /H) and the 138 dimensionless mean velocity, the Reynolds number, Re = 139 U H/ν, where φ is the mean heat flux, evaluated at the 140 top or at the bottom of the tank, λ is the thermal con-141 ductivity and U is a representative velocity of the fluid.142 Depending on the definition, it can refer to a large scale 143 speed, to a vertical or horizontal mean speed, as well as 144 to a Root Mean Squared (RMS) speed.

145
The setup is composed of a plexiglas tank with width 146 x, depth y and height z equal to 30 cm × 10 cm × 15 cm 147 respectively, and of two metallic heat exchangers lo-148 cated at the bottom and at the top of the tank.Both 149 heat exchangers are set to the desired temperature with 150 the help of two thermostated baths.Their tempera-151 tures are recorded throughout the experiments.We per-152 formed experiments with two different fluids.The first 153 one is a mixture of 75 wt.%water and 25 wt.% ethy-154 lene glycol, and the second one is a mixture of 87 wt.% 155 water, 10 wt.% UCON™ oil, and 3 wt.%salt.During an 156   The viscosities of the fluids used differ on average by a factor four, which allows extending the range of P r number used.We compute the thermal expansivity from polynomial fittings of the density curve.The reference density we used is defined at the average temperature of each experiment.The thermal expansivity ranges from 2.7 10 −4 K −1 to 5.4 10 −4 K −1 .
We aim at extracting the 3D velocity field in the entire domain, including the vicinity of the tank boundaries with the best signal-to-noise ratio.This requires sufficiently high particle concentrations to maximize the accuracy of the PTV approach, together with particle concentrations low enough to allow the visualization of the domain away from the cameras.We found a good balance for a particle concentration of 8 10 −4 g L −1 .Given the volume of the tank (4.5 L), the area captured by the cameras (30 × 15 cm 2 ) and their resolution, this corresponds to 0.035 particle-per-pixel resolution.While we found this concentration to be an optimum for our setup, the exact value may be different for another configuration depending on lighting conditions, tank depth, and camera angles.We use four Complementary Metal Oxide Semi-conductor (CMOS) Imager M-lite 2M cameras with 1936×1216 pixels resolution and maximum acquisition rate of 100 frames per second.We use a white LED Flashlight 300 panel for which the pulse length and frequency can be varied for a wide range.The panel is composed of 72 white high power, high efficiency LEDs, ensuring the illumination of the entire volume of the tank.The cameras record frames at rates varying from 1 to 10 Hz, with higher frequencies for higher velocities.The fast illumination from the LED panel allows us to avoid blurred particles, and to get sharp pictures, independently of the recording rate.To ensure that no relative displacement of the tank and the cameras occurs, the whole apparatus is set

231
The main parameters for the experiments we conducted 232 are listed in Table 1.picture, where all non-zero intensity pixels are considered to be particles.The result is shown in Figure 4b.
The cleaning processing allows one to see more clearly all particles in the tank and shows that the tank is seeded homogeneously.We use the retrieved pictures to perform a volume self-calibration (Wieneke, 2008).This second calibration is performed by computing 3D positions of particles and using the error as a correction to the first calibration.This second calibration process reduces our average spatial error to about 0.1 pixel.
The STB method consists of the following steps.A fit is performed to the last known positions of tracked particles.This fit is used to predict the next position of the particles.The predicted particle is moved in every direction of space until its position matches the particle position identified on the images.This step is hereafter referred to as the shaking step.Once the known tracks have been updated, a triangulation is used to find new particles.Particles whose intensity falls below a certain threshold are removed.Newly found particles are then added to the list of tracks.
In tomographic-PIV or 3D PTV, particles are found at the intersection of lines of sight, reconstructed from images.This reconstruction can also lead to the appearance of spurious particles, the so-called ghost particles.Ghost particles induce measurement errors, and will increase the computational time of tomographic-PIV and PTV (Elsinga et al, 2011).With the STB method, the shaking step is of particular interest because it removes efficiently ghost particles.The STB method proves to be computationally more efficient than the tomographic PIV method (Schanz et al, 2013).Moreover, because it does not rely on correlation of groups of particles, it has a higher spatial resolution than PIV methods (Schanz et al, 2013(Schanz et al, , 2016)).We conducted several series of tests to determine the optimal set of pre-processing and STB parameters, our goal being to extract large numbers of long tracks.
For each experiment, we record 20 pictures at different recording rates and run the previously described STB processing steps on each image series.For each recording rate, we check both the number of detected particles and the ratio of continuous tracks to new tracks detected.This allows us to choose the best recording frequency for each experiment.We then record longer movies that we analyze following the procedure described above.
We show the result for a typical experiment in Fig- ure 4c.The tracks are colored according to their corresponding velocity, and exhibit a LSC, as well as low velocity regions, located close to the boundaries.All the experiments we conducted presented such a LSC,  the one observed in (Schanz et al, 2013) and in (Schanz et al, 2016).This implies that we cannot follow parti-316 cles during the entire duration of the recorded experi-317 ment, making it impossible to use these tracks directly 318 to study long-term fluid trajectories.However, because 319 of their high density and good spatial sampling of the 320 domain, we can perform an interpolation of the track 321 velocities onto a fine grid to retrieve the 3D velocity 322 field in the entire tank.To do so, we use the method 323 implemented in DaVis (Schneiders and Scarano, 2016; 324 Jeon et al, 2018).This step allows one to project track 325 data (velocity and acceleration) onto a regular grid,

326
while enforcing additional physical constraints on the 327 flow.The method relies on a cost function minimiza-328 tion J expressed as: The first two terms on the right hand side of Equa-330 tion (3) J u and J a aim at minimizing the difference be- tween the velocity and acceleration on the chosen grid 332 and the velocity and acceleration of the tracks, and are 333 expressed as: where u is the velocity on the chosen grid, u PTV is the 335 velocity on the tracks, a is the acceleration on the cho-336 sen grid and a PTV is the acceleration on the tracks.345 The weighing coefficient G that appears in Equation ( 3) is defined as the ratio of particle tracks to the number of estimate the the viscous boundary layer to be thicker 360 than 6.7 mm (Grossmann and Lohse, 2000).We there-361 fore choose a resolution of 2.28 mm, corresponding to 362 12 pixels, which ensures that at least four points lie in the viscous boundary layer.
The result of this reconstruction is shown in Fig- ure 4d.This approach reconstructs well the velocity field already seen with the tracks, in particular the LSC mentioned above.The magnitude of the velocity field decreases close to zero in the vicinity of the rigid boundaries of the tank (marked by the black box in Figure 4d), but increases again outside of the tank.The velocity field is not exactly zero at the boundaries.The presence of this unphysical, non-zero velocity field outside the tank is only due to the extrapolation mentioned above, and is an artifact of the reconstruction method used.The velocity field is reliable wherever the tracks are present, i.e., within the tank, but becomes less accurate when closer to the boundaries of the tank, because the track particle density decreases near the walls.These non-physical velocities that do no conform to the rigid boundary conditions along the tank walls is problematic when one is interested in reconstructing accurately the long-term trajectories, or simply to characterize the velocity structure near domain boundaries (e.g., viscous or thermal boundary layers).Consequently, additional treatment is required to enforce the relevant boundary conditions to the reconstructed velocity field, which we discuss below.

Boundary conditions enforcement for the retrieved velocity field
To properly account for the relevant boundary conditions along the tank walls, the velocity field, u, resulting from the reconstruction has to be corrected.In the frame of our fluid dynamics experiments we aim at implementing rigid boundary conditions, i.e., velocity field equal to zero on the boundaries, without affecting the velocity field inside the tank.Additionally, the experimental fluids we consider have very small density variations (Figure 2), hence the velocity field to be extracted is nearly solenoidal.Under these conditions, one cannot only set the velocity at the boundaries to zero, which would propagate a non-zero divergence inside the tank.Instead, we derived a velocity correction based on the Helmholtz-Hodge decomposition whereby the raw velocity field, u, is expressed as the sum of a divergence-free field, v, and a curl-free field, ∇P : Taking the divergence of Equation ( 6) above and using the fact that ∇•v = 0 yields the following Poisson equation for the scalar field P : The correction method described above modifies the 456 velocity field along the boundaries, but may also prop-457 agate into the volume.In the following, we evaluate the extent of such propagation, showing that the amplitude of the correction is overall small, and remains essentially confined to the boundaries.
A straightforward way to evaluate the magnitude and the extent of the correction is to compute the difference between the raw and the corrected velocity fields.Figure 5 shows the RMS of this difference for experiment 2 in Table 1.We found these results to be similar for other velocity fields from the same experiment taken at different times, or from other experiments we performed with distinct values of the governing parameters.
Figure 5 shows that the correction is overall small.It is close to zero in the bulk of the domain and the largest differences lie along the boundaries.As mentioned above this is due to (i) the smaller density of particle tracks close to the boundaries, and (ii) the absence of constraints on rigid boundaries during the interpolation process of the Davis STB algorithm used to extract the raw velocity field.Consequently, we do not expect the velocities to be accurate in the vicinity of domain boundaries prior to applying the correction.Therefore, the observed differences between the raw and corrected fields remain confined to the boundaries, as expected.
In greater details, Figure 6 shows the histogram distributions for each component of the velocity field for the same case displayed in Figure 5.All velocity field components have a symmetrical distribution, centered around about zero, with a peak at v x , v y , and v z =0.While the raw distributions are smooth, the corrected distributions exhibit much sharper peaks at 0 velocity, and smaller count values around zero compared to the raw distributions.This stems from the fact that upon the velocity correction, we replace raw (extrapolated) values, that are close to zero, by velocity values that are exactly equal to zero to enforce the proper rigid boundary conditions.
Since the STB extraction induces velocity extrapolations in regions that lack particle tracks, the resulting magnitudes of velocity components in the vicinity of the domain boundaries are proportional to those located in the bulk.This would imply that the subsequent correction scales with the magnitude of velocity components in the bulk.To verify this hypothesis, we computed the RMS velocity for the entire tank in both cases.For all our experiments, we found that the raw velocities are systematically higher than the corrected velocities with a RMS difference ranging between 2% and 6%.Fig. 6 Histograms of the velocity components for experiment 2 in Table 1 along the (a) x, (b) y, and (c) z directions, prior to and after performing the correction using the Helmoltz-Hodge decomposition.The histograms are truncated for visibility.For V = 0, the counts are 63242, 89328 and 76302 for V x , V y and V z respectively.
Furthermore, our velocity correction also ensures a divergence-free (or a weakly divergent) velocity field to comply with the fact that our experimental fluids have negligible density variations.
Overall, these tests demonstrate that the velocity correction we apply to enforce the correct boundary conditions is only important in the vicinity of the domain boundaries and weak away from them.It allows obtaining a velocity field that is consistent with our rigid boundaries, which represents a considerable improvement.

Comparison with numerical simulations
Reproducing the laboratory results with numerical sim-

535
We reproduced the exact experimental conditions 536 (temperature difference between top and bottom plates, 537 viscosity, thermal expansivity and aspect ratios) with 538 StreamV, a 3D Finite-Volume numerical code that solves 539 for the Navier-Stokes equations coupled to the conser-540 vation of energy for a fluid under the Boussinesq ap-541 proximation (Samuel and Evonuk, 2010;Samuel, 2012).542 These governing equations are written below in dimen-543 sionless form: where v is the velocity vector, t is the time, p is the 548 dynamic pressure, η is the dynamic viscosity, α is the 549 thermal expansivity, T is the temperature (all in dimensionless units).e z is the vertical unit vector point- via Ra = α(T =0.5) η(T =0.5) × Ra 0 and P r = η(T = 0.5) × P r 0 .
Their values are listed in Table 1.
The dimensions of the model domain, and the tem- has "forgotten" its initial conditions.

Scaling laws
To characterize the dynamical behavior for each exper- overturn times.Figure 7 shows the Reynolds number as 601 a function of the Rayleigh number.The Prandtl num-602 ber, which depends exclusively on fluid properties, can-603 not be considered as constant for each of the two experi-604 mental fluids we considered, and ranges between 10 and 605 100, due to the temperature dependence of fluid prop-606 erties.Based on the literature, we account for these P r 607 variations by compensating Re with P r b with b = 0.7 608 (Chavanne et al, 1997;Méthivier et al, 2022).

609
Because we used approximate fits to the thermal 610 expansivity and viscosity, laboratory experiments and 611 their numerical replications can exhibit slight differ-612 ences in Ra and P r.The differences are not systematic and do not exceed 5% of the desired Ra and P r values.Due to these Ra − P r discrepancies between analog and numerical experiments, we do not expect Re values to coincide exactly for each analog experiment and its numerical replication.Yet, Figure 7 shows that the obtained Re values from the numerical experiments reproduce well the analog values, even if the numerical results seem to be systematically higher than the experimental results.One should note that even when Ra and P r -estimated at the average temperature -are equal for laboratory and numerical experiments, there may still be discrepancies in viscosities and thermal expansivities.These discrepancies can account for the systematic shift we observed.For comparison, we plot ReP r 0.7 ∝ Ra 0.5 , in Figure 7.The 0.5 exponent value corresponds to various findings in laboratory experiments (Ashkenazi and Steinberg, 1999;Chavanne et al, 1997;He et al, 2015), numerical simulations (Silano et al, 2010) as well as theoretical considerations (Grossmann and Lohse, 2000;Stevens et al, 2013).While this power law agrees well with our data Figure 8, we plot our data points as well as findings from previous laboratory experiments (Chavanne et al, 1997;Méthivier et al, 2022;He et al, 2015).This plot includes a wider range of Ra number, allowing us to see the overall trend more clearly than with only our experiments.We can see that low Re data points are distant from the blue dashed curved and high Re data Fig. 9 Spatial power spectral density for numerical and laboratory experiments, along the x direction.This graph is computed from experiment 2 in Table 1.The shaded areas represent a 95% confidence interval, defined as two times the standard deviation over time.
points collapses on the same curve, as reported in pre-666 vious Rayleigh-Bénard convection studies.

667
Overall, the comparison between analog and numer-668 ical experiments shows a reasonable agreement between 669 the two datasets.Both datasets exhibit a power law de-670 pendence of Re with Ra close to the expected law with 671 a 0.5 exponent value suggested by previous analog, nu-672 merical and theoretical works.To further compare the flow characteristics obtained 675 with the numerical code with those of the laboratory 676 experiments, we computed the spatial spectrum of the 677 velocity field.To do so, we interpolated the data from 678 the numerical code onto the same grid as the exper-679 imental data.We then performed a 3D fast Fourier 680 transform on each component of the velocity field.Fig- 681 ure 9 displays the sum of the power spectral densities, 682 averaged over y and z directions.Here, we also per-683 formed a time average over 8 overturns.Since we are 684 interested in comparing relative differences, we used the 685 normalized Power Spectral Density (PSD).
686 Figure 9 shows that numerical and laboratory ex-687 hibit the same behavior over a large range of spatial 688 frequencies.Most of the energy lies at small frequencies 689 in k x , which correspond to the LSC developed inside the 690 domain (Figure 4c).The PSD values then decay with 691 increasing k x .Taking into account the confidence in-692 terval, we can consider that the numerical experiments 693 reasonably reproduces the spatial variations observed 694 in the PTV experiments.We still notice a dip in the 695 spectrum at frequencies higher than 40 m −1 , which cor-696 responds to a spatial scale of 2.5 cm.This dip of un- focus on measuring 2D and/or thin 3D domain, or mean velocity, this work stands among the first to extract the complete 3D velocity field in a Rayleigh-Bénard convection laboratory experiment (Paolillo et al, 2021;Godbersen et al, 2021).Reproducing our laboratory experiments with a numerical code has allowed us to validate our method.Because we enforce a zero divergence at every step, the resulting velocity field is adapted to advection problems, such as the estimation of outgassing times of magma ocean (Salvador and Samuel, 2023;Lebrun et al, 2013).Using other velocity measurements for the advection of tracers would not be possible with partial measurements as all the components of the velocity field are needed.Since the velocity boundary conditions are exact with our approach, it is also possible to use this technique for measurements of the flow gradient near the walls in turbulent problems (Bitter et al, 2011;Marusic et al, 2010).The correction we performed significantly improves the accuracy of near wall measurements.An accurate velocity field, with accurate physical properties could also prove useful in the study of LSCs reversals (Xi and Xia, 2007;Funfschilling et al, 2008).Here, measurements of the flow at the walls could allow one to successfully estimate the viscous damping, and assess its importance (Brown and Ahlers, 2006;Sreenivasan et al, 2002).Most of the STB computational time comes from the velocity reconstruction and depends of the grid size used to extract the velocities.We used a grid resolution of 2.28 mm along each direction.With this resolution, the total computing time for 2000 frames is about one week on a modest desktop computer (32 GB RAM, Intel(R) Xeon(R) W-2133 CPU 3.6 GHz, Nvidia quadra p400).This computational time does not depend on the Ra number of the experiment, making this technique suitable for studying high Ra cases (which could be obtained using a different experimental fluid).

Conclusions
We proposed an improvement of the PTV Shake-The-Box method to extract the 3D velocity field in the entire domain of a fluid dynamics experiment including in the vicinity of the domain boundaries.In these regions the velocity field tends to be less accurate due to the rarefaction of particle tracks in these regions.Our new method relies on a post-processing correction of the raw velocity field.We tested this new method on a laboratory experiment setup to study Rayleigh-Bénard convective motions in a tank.

Fig. 1 Fig. 2
Fig. 1 Experimental setup used embedded in a rigid structure.(A) high repetition rate LED panel.(B) Plexiglas tank, filled the working fluid, illuminated in blue.(C) CMOS cameras, mounted on the rigid frame and oriented with different angles. 172

Fig. 3
Fig. 3 Temperature variations of the viscosities for the two experimental fluids used, whose compositions are given in Section 2.1.

208
into a rigid mechanical structure.The entire apparatus 209 is shown in Figure 1.The cameras are all recording 210 the entire tank volume, but are tilted relative to one 211 another.Prior to performing the experiments, we in-212 troduce a regular grid into the tank, with round marks.213 To extract the in-depth information from the pictures 214 taken, we map the oblique views of the cameras onto 215 the grid, using a polynomial fit.The average deviation 216 to the marks after this calibration process is 0.3 pixels.
experiments as follows.First, we fill 219 the tank with the fluid, making sure there are no air 220 bubbles left in the tank.The convection starts soon 221 after the heat exchangers are connected to the ther-222 mostated baths.The temperatures in the heat exchang-223 ers are continuously monitored by thermocouples and 224 Platinum sensors.Local vertical profiles are recorded 225 using two linear probes equipped with thermocouples 226 placed at several different heights.Once all tempera-227 tures have reached a stable value, after about 1 h, we 228 remove the vertical probes from the tank and we record 229 long series of images of the convecting fluid.An exam-230 ple of the pictures we record is shown in Figure 4 (a).
particle trajectories and to extract 235 the velocity field, we use the STB method, which was 236 introduced in Schanz et al (2013) and in Schanz et al 237 (2016).We apply our method in two stages.The first 238 stage is performed using the LaVision™ commercial soft-239 ware DaVis, while the second stage is performed sepa-240 rately with a home-made post-processing freeware code 241 STB correction (Samuel, 2023).
not clearly visible in the raw im-244 ages (Figure 4a), therefore to process efficiently the 245 recorded images, we first perform a cleaning step, aim-246 ing at keeping only the seeded particles in our pic-247 tures.Starting from the recorded images, we compute 248 the minimum intensity over time (over 11-frame win-249 dows) and the local minimum intensity (averaged over 250 four pixels) and subtract those from each picture.We 251 reduce the background intensity by removing a con-252 stant intensity of 100 from each image and artificially 253 increase the signal by applying a scaling factor of 10 to 254 each pixel intensity.This results in a dark background 255

Fig. 4
Fig.4Illustration of the processing steps for experiment 1 in Table1.(a) Raw recorded image, with the axis definition.Some small bubbles can be seen at the top of the tank.(b) Pre-processed image with the parameters described in Section 2.2.Each bright dot will be considered a particle in the STB treatment.(c) Tracks obtained with the STB method.We show 11 consecutive frames so that trajectories of the particles can be seen.(d) 2D cross-sections of the velocity field obtained after the fine scale reconstruction.The XZ slice is located in the middle of the tank and allows us to see the flow field inside the tank.Boundaries of the tank are highlighted in black.

337A
tricubic interpolation from the grid to the parti-338 cle positions is performed.The weighing coefficient A 339 is computed from the ratio of the standard deviation 340 (σ) of the velocity to that of the acceleration, A = 341 σ u,P T V /σ a,P T V .342 A divergence-free velocity field is also enforced by 343 adding the following term to the cost function in Equa-

) 411 with
Neumann boundary conditions ∂P/∂n = 0 along 412 the domain boundary ∂Ω, where n is a unit vector 413 normal to ∂Ω.Note that in the above equations the 414 components of the raw velocity field u along the tank 415 boundaries have been set to zero.The above equation 416 is solved using a geometric multigrid method (Brandt, 417 1982) on a staggered grid (Patankar, 1980) perform-418 ing V-cycles and relying on a Jacobi smoother.Prior to 419 solving Equation (7), we perform two series of tri-linear 420 interpolations of the velocity field obtained from the 421 previously described STB procedure.A first interpola-422 tion projects the initial velocity field onto a Cartesian 423 co-located grid whose boundaries coincide exactly with 424 the tank walls.The number of grid points along each 425 direction is chosen to facilitate the subdivision of the 426 domain into a hierarchy of grids required by the multi-427 grid Poisson solver, while ensuring a spatial resolution 428 as close as possible to that of the input field, to re-429 duce interpolation errors.Then, a second interpolation 430 of the velocity field is performed from this co-located 431 grid onto a corresponding staggered grid, where Equation (7) is solved.In this grid the scalar field P is defined 433 at cell centers, and the velocity components are located 434 at the center of each cell surface.Therefore, upon solv-435 ing for the Poisson equation we impose rigid boundary 436 conditions to the right hand side of Equation (7) us-437 ing ghost nodes for velocity components parallel to the 438 boundaries (Patankar, 1980).Once the Poisson Equa-439 tion (7) is solved, we use Equation (6) to obtain v.This 440 corrected velocity field is then projected back onto the 441 original co-located grid using tri-linear interpolation.442 Convergence of the numerical solution for Equation (7) 443 is achieved whenever the RMS residual is less than 444 10 −5 , and the arithmetic average of ∇ • v normalized 445 by the RMS value of v is less than 10 −6 .This usually 446 requires 4-10 V-cycles depending on the raw velocity 447 field.Tests performed using more restrictive conditions 448 did not result in significant differences.This correction 449 step is performed using a Fortran code STB correction 450 (Samuel, 2023) parallelized using OpenMP directives, al-451 lowing for the correction of a 3D raw velocity field typ-

Fig. 5
Fig. 5 Root Mean Square of the difference of the velocity fields, prior and after reconstruction, for three different positions.The differences are normalized by the V RM S after the correction, averaged on the whole domain.The data shown here is from experiment 2 in Table 1.(a) At the wall furthest from the cameras, (b) In the middle of the tank, (c) At the wall closest to the cameras.
ulations is useful for two different reasons.Up to now, there has been numerous Rayleigh-Bénard laboratory convection experiments, considering various domain aspect ratios, various types of fluids with different properties, and different measurement methods.These differences complicate the direct comparison between different experiments.Moreover, most of these convection experiments focus on heat flux measurements and the velocity is only estimated locally and/or indirectly.Nu-531 merical experiments reproducing our exact experimen-532 tal conditions then allow us to compute and compare 533 Reynolds numbers directly and to derive the spatial 534 characteristics of the flow.
ing upwards, and Ra 0 and P r 0 are the Rayleigh and Prandtl number, computed at the top boundary.The code was successfully benchmarked against various numerical and analytical solutions (Salvador andSamuel,     2023;Samuel and Evonuk, 2010;Samuel, 2012).Both viscosity and thermal expansion are temperaturedependent, and we use the density and viscosity measurements performed on the fluids to reproduce the experiments.For each experiment, we perform a secondorder polynomial fit of the density measurements between the top and the bottom temperatures of the tank.We thus model the non-dimensional thermal expansivity as α = 1 + c T , where the constant c is derived independently for each experiment from the secondorder fit.The non-dimensional viscosity at the top of the tank is equal to 1 and follows an exponential law η = exp(−γT ), where γ was computed for each experiment via an exponential fit of the viscosity measurements between the top and bottom temperatures of the tank.The Rayleigh and Prandtl numbers of the numerical simulations are thus linked to Ra 0 and P r 0 perature and velocity boundary conditions are identical to those of the laboratory experiments.Initially, the fluid is at rest and the temperature field consists in a homogeneous dimensionless temperature of 0.5 bounded by top and bottom thermal boundary layers of dimensionless thicknesses 0.1.A small temperature perturbation is imposed at the bottom boundary to initiate motions.To avoid bias from the initial condition, each experiment is first run until reaching a statistical steadystate stage at which the top and bottom heat fluxes converge to the same value.At this stage, the system

Fig. 7 Fig. 8
Fig. 7 Comparison of Re variations with Ra between our numerical and laboratory experiments.The error bars correspond to a 95% confidence interval.They are computed from time series corresponding to 7-8 overturn times.
points at high Ra values, it does not capture the behavior at low Ra.Even though the Rayleigh number values remain relatively large (Ra > 10 7 ), the high Prandtl number values prevent the flow from reaching high Re numbers.Therefore, we do not expect experiments with low Ra values to exhibit the 0.5 exponent.Limiting our analysis to Re numbers larger than 200, we performed a fit on ReP r b = C 0 Ra a with a fixed exponent b on P r, for both numerical and laboratory experiments.For the numerical experiments, we obtain C 0 = 0.021 ± 0.015 and a = 0.59 ± 0.05.For the laboratory experiments, we obtain C 0 = 0.036 ± 0.020 and a = 0.56 ± 0.07.These values are close to the expected 0.5 exponent value.While the numerical exponent is larger than 0.5 (even when accounting for uncertainties), one should note that this value is strongly sensitive to the chosen value for b.With b = 0.8, we get a = 0.51 ± 0.07 and a = 0.53 ± 0.1 for laboratory and numerical experiments, respectively.Because we limited the above analysis to cases where Re > 200, the corresponding P r number does not vary much, which makes it difficult to obtain an accurate estimate of the b exponent.We can still compare our results with previous experiments by setting the exponent b to the well-accepted value.In

697
certain origin was systematically observed for all the 698 experiments.However, since it corresponds to low val-699 ues of the PSD, this dip does not play a role in the consists of using two fluorescent dyes, emitting in dif-736 ferent spectral bands upon excitation by a laser light 737 sheet, and whose emission respectively increase and de-738 crease with temperature.Using adapted color filters on 739 two different cameras and after an intensity calibration 740 procedure, this setup would allow one to retrieve the 741 temperature field, and to compute the heat flow in ver-742 tical cross-sections of the tank.743 Limitations and possible improvements of the velocity 744 correction method As most convection experiments only 745 Using a mixture of water and different fluids, and with temperature differences, we explored Prandtl number values ranging Voth GA, Haller G, Gollub JP (2002) Experimental measurements of stretching fields in fluid mixing.Physical review letters 88(25):254,501 Voth GA, Saint T, Dobler G, Gollub JP (2003) Mixing rates and symmetry breaking in two-dimensional chaotic flow.Physics of Fluids 15(9):2560-2566 Wieneke B (2008) Volume self-calibration for 3d particle image velocimetry.Experiments in fluids 45(4):549-556 Xi HD, Xia KQ (2007) Cessations and reversals of the large-scale circulation in turbulent thermal convection.Physical Review E 75(6):066,307