Alfredo García-Cid

Pablo Gutiérrez

Claudio Falcón

Sébastien Aumaître

Eric Falcon

Statistics of injected power on a bouncing ball subjected to a randomly vibrating piston

Keywords: number(s): 05, 45, -a, 47, 20, Ky

ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

I. INTRODUCTION

Dissipative systems display striking dynamical behaviors when they are driven by an external forcing. A nonequilibrium steady state (NESS) can be achieved by balancing the effect of the external forcing and the inherent dissipative mechanisms of the system under study. In physics, the number of examples is vast. For instance, in fluid mechanics, a typical example of NESS is the so-called fully developed hydrodynamic turbulence [START_REF] Kolmogorov | Proc. R. Soc. Lond. A[END_REF][START_REF] Frisch | Turbulence: The Legacy of A. N. Kolmogorov[END_REF], where the energy pumped at a large scales is transported towards small scales where it is dissipated by viscosity. In granular materials, a granular gas is a NESS achieved by injecting energy through the boundaries of a container into a dilute sample of particles which will later dissipate it via inelastic collisions, converting kinetic energy into plastic deformations [START_REF] Mehta | Granular Physics[END_REF][START_REF] Aranson | Granular Patterns[END_REF][START_REF] Aumaître | [END_REF][6]. Similar examples of nonequilibrium steady states can be found in electronic circuits [7][8][9], turbulent thermal convection [10], two-dimensional turbulence [11], wave turbulence in elastic plates [12] and surface waves [13], molecular machines [14], and so forth. In these steady states, the injected power or other global quantities (i.e., quantities averaged over space and/or time such as energy, momenta, and number of particles) will fluctuate strongly and depend on both the injection mechanism and the impedance of the system to the driving. Thus, in a NESS these global quantities cannot be regarded as constant as they depend on the interaction between the system and the injection mechanism, even in the case of time-periodic (sinusoidal) forcing [15]. When the driving can be described statistically by a Gaussian random process with zero mean, their fluctuations can be several times larger than their mean values [START_REF] Aumaître | [END_REF]6,[11][12][13]. It prevents the use of a large set of theoretical tools which are of fundamental use in the thermodynamical limit [START_REF] Zwanzig | Nonequilibrium Statistical Mechanics[END_REF][START_REF] Kubo | Statistical Physics II: Nonequilibrium Statistical Mechanics[END_REF] which allow the prediction and control of systems with a large degree of components.

* cfalcon@ing.uchile.cl A simple example of how an out-of-equilibrium balance between injection and dissipation generates a NESS is an inelastic bouncing ball driven vertically by a randomly moving piston in the presence of gravity. In this simple system energy is pumped into the bouncing ball by the moving piston, which is dissipated mainly by plastic deformations of the bead and viscous drag with the air surrounding the ball. The dynamics of inelastic bouncing balls excited by vertical vibrations has been researched extensively theoretically, experimentally, and numerically as examples of permanent chaotic behavior when the piston moves periodically [18][19][20][21]. On the other hand, only a handful of theoretical studies have reported on the properties and features of the NESS dynamics of the inelastic bouncing ball system when the motion of the piston is quasiperiodic [22] or when the bead is subjected to a random acceleration [23]. It must be noticed that, experimentally, the injection of energy (or other quantities) into a dissipative system is usually accompanied by fluctuations, affecting the statistical properties of the NESS under study. Thus, an experimental study on the properties and features of the NESS achieved by a randomly excited inelastic bouncing ball is lacking.

In the present paper, we focus our work on the statistics of the injected power furnished by the randomly moving piston which is necessary to maintain an inelastic ball bouncing constantly. In Sec. II we present the experimental setup and measuring techniques used to quantify the injected power I n on each impact of the bouncing ball. Section III shows the experimentally measured probability distribution functions (PDF) for the force F n and velocity v n of the moving piston as the bead makes the nth impact onto its surface, and with this data we compute the PDF of I n . This distribution has its most probable value close to zero and displays two asymmetric exponential tails, depending on the restitution coefficient, the piston acceleration, and its frequency content. In Sec. IV we show that by using a simple modified bivariate distribution for F n and v n we can be compute theoretically the distribution of I n . Finally, in Sec. V we discuss the possible applications of the model and how it can be generalized to granular gases.

II. EXPERIMENTAL SETUP AND MEASUREMENT TECHNIQUES

The experimental setup is shown in Fig. 1. An 8 mm bead (made out of steel, brass, copper, polyurethane, or nylon, with ±1 μm tolerance in diameter, ±1 μm in sphericity with masses presented in Table I) is placed inside a vertical plexiglass enclosing tube (outer diameter 15 mm, inner diameter 9 mm), which is polished on its inside part to reduce sliding friction with the bead. The tube is screwed to a plexiglass piston driven randomly in time by a minishaker (B&K 4810). An accelerometer (B&K 4803) screwed to the oscillating piston measures the acceleration a(t) of the whole setup via a charge amplifier (B&K 2635). The random forcing of the plexiglass piston, generated by a desktop computer using LabView via an acquisition card (NI-6022) through a power amplifier (B&K 2706), is a zero-mean colored Gaussian noise with a controlled frequency bandwidth [f l ,f u ] and rms amplitude normalized with respect to gravity = a(t) 2 /g. In our experiments f u At the bottom end of the tube a piezoelectric impact sensor (PCB-200B02) is attached to the oscillating piston, making impact with the bead. The impact signal is conditioned via a signal conditioner (PCB-480C02) for data treatment. From the acceleration signal, the velocity of the piston is computed by analog integration of the signal through the signal conditioner. Both force and acceleration signals are digitally acquired at 750 kHz, numerically low-passed filtered at 50 kHz to avoid the electrical resonance of the force sensor at 70 kHz, and stored in the desktop computer.

A. Coefficient of restitution r

In order to understand the dissipative motion of an inelastic bead bouncing on a randomly moving piston, we first measure the experimental coefficient of restitution (COR) r for different beads of different materials following a scheme similar to the one proposed in Ref. [24]. To do this, the impact sensor is screwed on an aluminum block (10 × 20 × 4 cm 3 ), on top of which an aluminum tube (outer diameter 20 mm, inner diameter 15 mm) is placed. A micrometric screw holds the bead inside the tube at a height h o = 5 mm over the impact sensor, allowing a controlled release of the bead with almost negligible air entrainment. The initial impact velocity v o coming from the free falling of the bead can be estimated v o = √ 2h o g 30 cm/s. Horizontal and rotational velocities are also negligible for this experimental configuration.

r is defined as the ratio between the bead velocity before and after an impact and it is assumed to be constant [24]. Thus, the COR can be computed as r = v n+1 /v n . The measurement of v n and v n+1 is difficult, but we can compute its ratio using simple particle mechanics relating r with the time of flight of the particle.

B. Calculation of r via time of flight

The total time elapsed d n between the nth and n+1th impact can be measured from the series of force peaks (as the one shown in Fig. 2), which can be used to compute r. d n can be expressed as d n = τ n + t n , where τ n is the duration of the nth collision and t n is the time of flight in between shocks. These times are measured from the temporal trace of the impact sensor, which consists of a series of peaks showing the time events where the bead impacts the sensor bouncing repeatedly over the piston in shortening time intervals as shown in Fig. 3.

Assuming negligible dissipation in between these intervals the ball is performing parabolic flights lasting t n . d n is computed as the temporal difference between two consecutive maxima of the force signal (discriminating between the real maxima and those arising from the piezoelectric inertial oscillation with a temporal period ∼100 μs). τ n is computed as the temporal width of the impact curve (which is always smaller than ∼ 100 μs). Thus, τ n t n and d n t n . We compute the COR from the set of t n using point particle dynamics [24]:

t n = r × t n-1 = r n × 8h o g . ( 1
)
t (s ) F(t), v(t)(a .rb.units) d n τ n Γ τ n (μs) (a) (b) (c) t (

C. Dissipation effects

As stated above, in this approximation, no viscous drag, sliding friction, or velocity-dependent effects on r are taken into account. Taking the natural logarithm on both sides we compute r as the slope of log t n versus n, as is shown in Fig. 3. r is defined as the mean slope of five different trials. This procedure was used to compute the COR for five different materials, as presented in Table I. The COR was also measured following the same scheme but using the enclosing tube to constrain the bead horizontal dynamics. This shifts r values to lower ones, as air entrainment becomes important and even dominant in the bead dynamics. Values for the COR with (r c ) and without confinement (r) are shown also in Table I.

The difference between r and r c appears as a combination of the effect of air entrainment and wall friction due to confinement, affecting the measurement of the COR as the initial and terminal velocity after a bounce of the inelastic ball. In the literature, the problem of the terminal speed change due to confinement (which is responsible for the shift to lower values of r) has been extensively studied for different tube and particle diameters, terminal speeds, confinement fluids, and flow regimes (see Ref. [25] and references therein). Thus, we will work using the effective r c in what follows and leave this effect for future studies.

In Fig. 2 we show how d n and τ n are defined from the acquired temporal traces. The values of d n change as piston moves, but its typical scale is still much larger than τ n . The values of τ n change slightly depending on : τ n is of the order of 50 μs when the piston is in repose and 65 μs when the piston is moving randomly. These differences are taken into account when computing the discrete values of the exerted force and velocity of the piston, as is shown below.

III. EXPERIMENTAL RESULTS

As shown in Fig. 2, a typical trace of the force and velocity for a steel bead (r c = 0.770, = 1.09, f l = 0) displays peaks and oscillations. The force signal consists of a series of peaks separated by the randomly spaced time intervals t n (as explained above). Between these intervals, the ball performs parabolic flights inside the tube. On the other hand, the velocity signal of the piston randomly fluctuates around zero. Around each impact n, we compute the discrete force F n and velocity v n by integrating the force and velocity signals in the time interval τ n . Each experimental run lasts 100 s. Between 10 and 500 collisions are observed, depending on the control parameters. To obtain sufficient data points for statistical analysis of all values of the control parameters (∼ 5 × 10 4 ), 100 runs are done for each value of these parameters.

A. Statistics of F n and v n as a function of

Let us begin by studying the statistics of F n as a function of shown in Fig. 4. There is a transition from an exponential to a convexly shaped PDF with increasing . All the moments of F n increase with increasing . For large the PDF of F n can be fitted accurately by a Gaussian curve with only positive events almost centered at F n = 0. This effect can be better observed in the normalized PDFs for the rescaled variable F n /σ (F n ), displayed in the inset of Fig. 4 (top), where a Gaussian fit is plotted (although the fit is not a normalized Gaussian with standard deviation equal to unity). The bottom of Fig. 4 shows the statistics for v n as a function of . All the curves can be accurately fitted by a Gaussian curve centered at v n = 0 with a standard deviation growing linearly with . As for F n , this effect can be better observed in the normalized PDFs for the rescaled variable v n /σ (v n ) shown in the inset of Fig. 4 (bottom), where a normalized Gaussian curve is plotted to emphasize this dependence.

B. Statistics of F n and v n as a function of r c

Fixing now at 1.00 and f l = 0 Hz, we change the bead material. Five different materials were used with different values of r c (see Table I). The PDFs of F n for different values of r are shown in Fig. 5 (top). A large difference between metallic and plastic beads is found in the typical values of F n . The rms value of metallic beads is 10 times larger than the one for plastic ones. In the normalized variable F n /σ (F n ), the PDFs do not collapse on the same curve [see inset of Fig. 5 (top)]. The order of the curves is nonmonotonic with r c . The PDFs of v n for different values of r c are shown in Fig. 5 (bottom). The different curves are very similar in scale and shape. A slight asymmetry towards positive values of v n can be observed on the curves when plotted in the rescaled variable v n /σ (v n ) [see Fig. 5 (bottom)], as can be seen when compared to a normalized Gaussian curve.

C. Statistics of F n and v n as a function of f l

Fixing now at 1.00 and using a brass bead (r c = 0.742), we change the low-frequency cutoff f l from 0 to 150 Hz. The PDFs of F n for different values of f l are shown in Fig. 6

(top).

As f l increases, the typical rms fluctuations of F n decrease and new local maxima appear in the PDFs. These maxima are related to the appearance of a pseudoperiodic motion of the , brass (×), bronze ( ), nylon ( ), and polyurethane (+) for = 1.00 and f l = 0 Hz. Inset: PDF for the normalized variable F n /σ (F n ) for steel (•), brass (×), bronze ( ), nylon ( ), and polyurethane (+) for = 1.00 and f l = 0 Hz. (b) PDF for v n for steel (•), brass (×), bronze ( ), nylon ( ), and polyurethane (+) for = 1.00 and f l = 0 Hz. Inset: PDF for the normalized variable v n /σ (v n ) for steel (•), brass (×), bronze ( ), nylon ( ), and polyurethane (+) for = 1.00 and f l = 0 Hz. bead as the bandwidth of the forcing becomes narrower around 200 Hz (the upper cutoff frequency of the selected bandwidth). When plotted in the rescaled variable F n /σ (F n ), the curves present a great level of collapse, although not complete. The PDFs of v n for different values of f l are shown in Fig. 6 (bottom). They are similar in shape, and their respective rms fluctuations decrease with increasing f l . A slight asymmetry can be observed on the curves when plotted in the rescaled variable v n /σ (v n ) [see Fig. 6 (bottom)], as can be seen when compared to a normalized Gaussian curve.

D. Injected power I n statistics as a function of , r, and f l

Using the discrete force F n exerted by the piston on the bead and velocity v n of the piston, we can compute the discrete power I n = F n × v n injected on the bouncing bead system by the piston as a function of the nondimensional rms acceleration , the COR r c and the lower cutoff frequency f l , as is shown in Fig. 7. Following the scheme presented above, we first fixed r c = 0.770 and f l = 200 Hz, and we computed .60, 0.68, 0.80, 0.87, 1.00, and 1.09 for steel (r c = 0.770) and f l = 0 Hz. Inset: PDF of the normalized variable I n /σ (I n ) for = 0.50, 0.60, 0.68, 0.80, 0.87, 1.00, and 1.09 for steel (r c = 0.770) and f l = 0 Hz. Dashed lines are exponential fits for both positive and negative tails. (b) PDF for I n for steel (•), brass (×), bronze ( ), nylon ( ), and polyurethane (+) for = 1.00 and f l = 0 Hz. Inset: PDF for the normalized variable I n /σ (I n ) for steel (•), brass (×), bronze ( ), nylon ( ), and polyurethane (+) for = 1.00 and f l = 0 Hz. Dashed lines are exponential fits for both positive and negative tails. (c) PDF for I n for f l = 0 (•), 50 ( ), 100 ( ), and 150 ( * ) Hz for = 1.00 and r c = 0.742 (brass). Inset: PDF for the normalized variable I n /σ (I n ) for f l = 0 (•), 50 ( ), 100 ( ), and 150 ( * ) Hz for = 1.00 and r c = 0.742 (brass). Dashed lines are exponential fits for both positive and negative tails.

the statistics of I n by changing . The PDFs display a cusp center at I n = 0 and tails which are asymmetrical; i.e., the probability of finding positive injected power events (I n > 0) is larger than the probability of finding negative injected power events (I n < 0). This means that the mean injected power I n is positive for every experimental run, which is needed in out-of-equilibrium systems as power will later be dissipated by internal mechanisms. For very low values of , the tails decrease slower than an exponential (the convexity of the slopes is positive). As increases the slopes become exponential and the typical scale of the fluctuations of I n increases. In the rescaled variable I n /σ (I n ) the curves collapse into a master curve which displays exponential tails. When we fix = 1.00 and f l = 200 Hz and we change r c , the PDFs of I n display a large difference between metallic and plastic beads (as in the case for F n ). The rms value of I n for metallic beads is 10 times larger than the one for plastic ones. In the normalized variable I n /σ (I n ), the PDFs collapse on the same curve [see inset of Fig. 7 (top)], displaying asymmetric exponential tails and a cusp close to I n = 0. Fixing at 1.00 and r c at 0.742, we change f l from 0 to 150 Hz. As before, the PDFs of I n for different values of r c are similar in shape, and their respective rms fluctuations decrease with increasing f l . In the rescaled variable I n /σ (I n ) all curves collapse into the same master curve, which displays exponential asymmetric tails.

In computing I n we have used the discrete multiplication instead of the integral over the collision time τ n of the product of force and velocity. The statistics of either computation method are almost equal, except for < 0.6 where the cusp near I n is rounder and smaller. This means that the approximation of I n as the direct multiplication of F n and v n is a valid choice. In this framework, the features of this PDF can be understood, as is shown in the next section.

IV. FEATURES OF THE INJECTED POWER STATISTICS

A way to understand the qualitative features of the PDF of injected power is to adapt the simple model proposed in Ref. [13] following Craig's XY distribution [26]. From the acquired data, we can see that the fluctuations of v n can be described accurately by a Gaussian random variable with zero mean. The fluctuations of F n cannot be described as the fluctuations of v n , as no negative values of F n can be achieved, even though the positive F n events are well fitted by a Gaussian curve. Thus, from the data presented above, F n can be approximated a semi-Gaussian variable, i.e., a random variable whose PDF is the product of a Gaussian and a Heaviside function (suitably normalized). We can model the PDF of F n , P (F n ) as

P (F n ) = 2 π 1 α exp -F 2 n /2α 2 , for F n 0. 0, otherwise. ( 2
)
where α is the standard deviation of the Gaussian function if no restriction for negative values of F n is in place. Using the approximation that the exerted force per impact F n the inelastic bead into motion, while it changes the momentum of the piston, we calculate the injected power I n = F n × v n (even in this case where F n is always positive). Here the injected power into the bouncing ball system is computed by the multiplication of the velocity and force of the piston, as in Ref. [13]. The PDF of I n can be derived from direct integration of the joint probability distribution function for F n and v n with suitable integration domains.

Let us use a simple model to compute the PDF of I n . In the case of two Gaussian variables which are correlated (linear statistical dependence), the joint probability distribution function is a bivariate distribution

P (F n ,v n ) = exp -1 2(1-c 2 ) v 2 n σ (v n ) 2 -2cv n F n σ (F n )σ (v n ) + F 2 n σ (F n ) 2 2πσ (v n )σ (F n ) √ 1 -c 2 , ( 3
) where c = (F n -F n )(v n -v n ) σ (F n )σ (v n )
is the correlation coefficient for the pair of random variables (F n ,v n ). Although other joint distribution functions can be displayed for two Gaussian variables with a given correlation [START_REF] Hogg | Probability and Statistical Inference[END_REF], we model our system with this joint distribution because of its simplicity in the calculation of the PDF of I n . Changing variables to I n = F n × v n and u n = F n and integrating Eq. ( 3) over F n on the domain [-∞,∞] we obtain the PDF for I n :

P (I n ) = exp cI n (1-c 2 )σ (F n )σ (v n ) πσ (F n )σ (v n ) √ 1 -c 2 K o |I n | (1 -c 2 )σ (F n )σ (v n ) , (4)
where K o (•) is the zeroth order modified Bessel function of the second kind. A simple expression of the above equation can be achieved by the steepest descent method,

P ( X) ∼ exp (c X -| X|)/ | X|, (5)
with

X = I n /[(1 -c 2 )σ (F n )σ (v n )]. In addition I n = F n × v n = c × σ (F n )σ (v n )
, which connects the correlation coefficient of both variables, their rms fluctuations, and the injected power into the system. This expression displays two exponential tails controlled solely by the correlation coefficient c, which is equal to the mean injected power normalized when v n = 0. We can obtain an expression for the PDF of I n by assuming that F n is semi-Gaussian and integrating it out in Eq. ( 3). The new PDF for I n has the same form as Eq. ( 4) with two modifications: by direct computation F n = √ 2/π α and σ (F n ) = α 2 (1 -2π ), which are functions of α. With these substitutions in Eq. ( 3), the calculation of P (I n ) is the same as the one above. To wit, we plot in Fig. 8 the PDFs of I n in the reduced variable X/(1 -X 2 ) with

X = I n √ 1 -2/π /[σ (F n )σ (v n )]
and the theoretical approximation for its PDF following Eq. ( 5). There is a good agreement for every COR between experimental and theoretically predicted curves. As predicted, the approximation to the theoretical PDF is good, although the negative tails depart from the predicted slope slightly for polyurethane. Notice that Eq. ( 4) uses no fitting parameters: all cumulants of X can be computed by knowing c.

V. DISCUSSION AND CONCLUSIONS

The discrete injected power I n needed to drive an inelastic bead into a NESS is studied. Power is injected via a randomly accelerating piston that drives the motion of the bead. Measuring both the mean force exerted by the piston over the inelastic bead and the mean velocity of the piston at each impact we can access the complete statistics of I n . We showed that by changing the rms fluctuations and frequency band of the random acceleration, and the restitution coefficient of the bead, the typical scale of the fluctuations of I n can change drastically, although the shape of the PDF of I n does not, as it is controlled solely by the normalized injected power c =

(F n -F n )(v n -v n ) /σ (F n )σ (v n ).
Using a simple model stemming from a bivariate distribution function the PDF can be predicted only knowing c, σ (v n ), and σ (F n ). This is shown in Fig. 8: as r c is similar for four out of five materials used in the experiment, the PDFs of the normalized injected power display the same shape for a given and f l and only for the case of polyurethane (r c = 0.634) shows a large deviation from the other curves. Thus, in this simple out-of-equilibrium case, one can predict and compute the PDF of I n with only a very small set of external parameters which control completely the injected power statistics into an out-of-equilibrium system.

The PDFs exhibited in Fig. 8 are actually quite generic of the injected power furnished by a random force. The simplest system of this type is a single particle submitted to a viscous drag and a random external force. This system can be described as an Ornstein-Uhlenbeck random process (a noise with exponential correlation [START_REF] Zwanzig | Nonequilibrium Statistical Mechanics[END_REF][START_REF] Kubo | Statistical Physics II: Nonequilibrium Statistical Mechanics[END_REF]). In the case using the joint bivariate normal probability function (3), the PDF of the injected power can be computed exactly [START_REF] Aumaître | Fluctuations and Instabilities in Out-of-Equilibrium Systems (Habilitation à diriger les recherches (HDR)[END_REF]. Experimentally, such PDFs have been observed for a randomly vibrated thin elastic plate [12] as well as in randomly forced surface waves [13]. The shape of these PDFs can be also achieved in turbulent shell model submitted to random force [START_REF] Aumaître | Fluctuations and Instabilities in Out-of-Equilibrium Systems (Habilitation à diriger les recherches (HDR)[END_REF]. This kind of fluctuations can be observed for other quantities resulting from the product of two centered random variables. The components of the turbulent Reynold stress tensor [START_REF] Pumir | [END_REF] or the rate of transferred energy in turbulence [30] displays the same features. They also appear in turbulent convection for the local heat transfer, which is mainly the product of temperature and vertical velocity fluctuations. Asymmetric PDFs with a cusp near zero and exponential tails were observed for the temporal fluctuations of the locally measured heat flux [31], for the instantaneous fluctuations in space of the numerically computed heat flux [32], and for the measurement of the heat flux using a Lagrangian probe [33]. There may exist other situations where such fluctuations are displayed. Even if the two variables involve in the product are not necessary Gaussian, the expression for the PDF of the injected power (4) is a useful benchmark to study the fluctuations of such a product. Thus, this simple yet powerful approach to the global fluctuations of an out-of-equilibrium system used in other NESS [8,[11][12][13] serves as a tool to study the statistical properties of different systems.

In that sense, the use of this approach can be thought as an opportunity to rethink the pertinence of some of the tools of out-of-equilibrium fluctuations, for instance, fluctuation theorems [34], which are used to compute parameters associated to dissipative constants and effective energy scales (as in the case of the granular temperature [START_REF] Aumaître | [END_REF]6]). Our work shows that from the statistics of I n we can distinguish what would be an effective dissipation coefficient in a dissipative out-of-equilibrium system. In our experiments, the different masses of the beads do not affect the PDF of the reduced variable X, even when the difference is one order of magnitude in mass. Only the different values of r c (which in our case depends on geometrical factors, viscous drag, friction, and plasticity of the material) affect the measured statistics of X as is shown in Fig. 8. Thus, r c is a relevant control parameter of the dissipative system as it controls the amount of power that can be injected into the bouncing bead. The idea of the construction of an effective dissipation coefficient in dissipative out-of-equilibrium systems (such as the case of r c ) can thus be exported to systems where internal dissipative mechanisms are not yet measured or not well understood. Simple models for the statistics of the injected power, as the one presented here or elsewhere [START_REF] Aumaître | Fluctuations and Instabilities in Out-of-Equilibrium Systems (Habilitation à diriger les recherches (HDR)[END_REF], can be used to quantify, for instance, typical dissipative scale differences between systems as internal parameters are changed, even though no complete information of the internal microscopic dynamics is known. The application of this scheme on other dissipative out-of-equilibrium systems is then needed to confirm its potential use.

FIG. 1 .

 1 FIG. 1. (Color online) Experimental setup. Inset: Typical temporal trace of the force signal for a steel ball (ball mass m = 8 g, restitution coefficient r = 0.923) at a given control acceleration = 1.09.

  FIG. 2. (Color online) (a) Temporal traces for the force (black) and velocity (red) traces between two impacts, defining d n for = 1.09. (b) Temporal traces for the force (black) and velocity (red) traces, around an impact for = 1.1. The shaded region defines τ n . Arbitrary units are used to place both curves in the same scale. (c) Dependence of τ n as a function of . Error bars stand for the standard deviation for τ n . The bead is made of steel.

1 FIG. 3

 13 FIG. 3. (Color online)Typical trace of force impacts of a steel ball on a piston on repose. Symbols (×) stand for peak detections. Inset: Semilog plot of time lag between impacts t n vs n for a steel ball. Best fitting slope is ln(-0.93).

FIG. 4 .

 4 FIG.4. (Color online) (a) PDF for F n for = 0.50, 0.60, 0.68, 0.80, 0.87, 1.0, 1.09, and 1.28 for steel (r c = 0.770) and f l = 0 Hz. Inset: PDF of the normalized variable F n /σ (F n ) for = 0.50, 0.60, 0.68, 0.80, 0.87, 1.00, and 1.09 for steel (r c = 0.770) and f l = 0 Hz. (b) PDF of v n for = 0.50, 0.60, 0.68, 0.80, 0.87, 1.0, and 1.09 for steel (r c = 0.770) and f l = 0 Hz. Inset: PDF of the normalized variable v n /σ (v n ) for = 0.50, 0.60, 0.68, 0.80, 0.87, 1.0, 1.09, and 1.28 for steel (r c = 0.770) and f l = 0 Hz. Dashed line is a normalized Gaussian. Arrows indicate how increases.

FIG. 5 .

 5 FIG.5. (Color online) (a) PDF for F n for steel (•), brass (×), bronze ( ), nylon ( ), and polyurethane (+) for = 1.00 and f l = 0 Hz. Inset: PDF for the normalized variable F n /σ (F n ) for steel (•), brass (×), bronze ( ), nylon ( ), and polyurethane (+) for = 1.00 and f l = 0 Hz. (b) PDF for v n for steel (•), brass (×), bronze ( ), nylon ( ), and polyurethane (+) for = 1.00 and f l = 0 Hz. Inset: PDF for the normalized variable v n /σ (v n ) for steel (•), brass (×), bronze ( ), nylon ( ), and polyurethane (+) for = 1.00 and f l = 0 Hz.

FIG. 6 .

 6 FIG.6. (Color online) (a) PDF for F n for f l = 0 (•), 50 ( ), 100 ( ), and 150 ( * ) Hz for = 1.00, and r c = 0.742 (brass). Inset: PDF for the normalized variable F n /σ (F n ) for f l = 0 (•), 50 ( ), 100 ( ), and 150 ( * ) Hz for = 1.00 and r c = 0.742 (brass). (b) PDF for v n for f l = 0 (•), 50 ( ), 100 ( ), and 150 ( * ) Hz for = 1.00 and r c = 0.742 (brass). Inset: PDF for the normalized variable v n /σ (v n ) for f l = 0 (•), 50 ( ), 100 ( ), and 150 ( * ) Hz for = 1.00 and r c = 0.742 (brass).

FIG. 7 .

 7 FIG.7. (Color online) (a) PDF for I n for = 0.50, 0.60, 0.68, 0.80, 0.87, 1.00, and 1.09 for steel (r c = 0.770) and f l = 0 Hz. Inset: PDF of the normalized variable I n /σ (I n ) for = 0.50, 0.60, 0.68, 0.80, 0.87, 1.00, and 1.09 for steel (r c = 0.770) and f l = 0 Hz. Dashed lines are exponential fits for both positive and negative tails. (b) PDF for I n for steel (•), brass (×), bronze ( ), nylon ( ), and polyurethane (+) for = 1.00 and f l = 0 Hz. Inset: PDF for the normalized variable I n /σ (I n ) for steel (•), brass (×), bronze ( ), nylon ( ), and polyurethane (+) for = 1.00 and f l = 0 Hz. Dashed lines are exponential fits for both positive and negative tails. (c) PDF for I n for f l = 0 (•), 50 ( ), 100 ( ), and 150 ( * ) Hz for = 1.00 and r c = 0.742 (brass). Inset: PDF for the normalized variable I n /σ (I n ) for f l = 0 (•), 50 ( ), 100 ( ), and 150 ( * ) Hz for = 1.00 and r c = 0.742 (brass). Dashed lines are exponential fits for both positive and negative tails.

FIG. 8 .

 8 FIG. 8. (Color online) PDF of the reduced variable X/(1 -X 2 ) with X = I n √ 1 -2/π /[σ (F n )σ (v n )] for steel (•), brass (×), bronze ( ), nylon ( ), and polyurethane (+) for = 1.00 and f l = 0 Hz. Dots show the approximation P (X) ∼ exp (cX -|X|)/ √ |X| for each value of r c .

TABLE I

 I Hz. Thus, f l and serve as control parameters.

				is fixed at 200
	. Masses and COR for different materials using spheres
	of 8 mm in diameter. The COR with and without confinement are r c
	and r, respectively.
	Material	m (g)	r c	r
	Steel	2.053 ± 0.010	0.770 ± 0.011	0.923 ± 0.006
	Bronze	2.235 ± 0.033	0.746 ± 0.012	0.861 ± 0.028
	Brass	2.383 ± 0.005	0.742 ± 0.036	0.860 ± 0.015
	Nylon	0.295 ± 0.006	0.757 ± 0.011	0.740 ± 0.042
	Polyurethane	0.308 ± 0.008	0.634 ± 0.024	0.667 ± 0.029

ACKNOWLEDGMENTS C.F. thanks Univ Paris Diderot for their Invited Researcher funding in Paris. E.F. thanks Ecos-Chili ECOS C07E07 for their Invited Researcher funding in Santiago. This work has been supported by FONDECYT 1130354, CONICYT FONDECYT Postdoctorado 3140550, ACT127 and European Space Agency Topical Team on granular materials No. 4000103461.


domain: phys nlin

text:
Pablo Gutiérrez

Sébastien Aumaître
email: sebastien.aumaitre@cea.fr

Clustering of floaters on the free surface of a turbulent flow: An experimental study

Keywords:

We present an experimental study of the statistical properties of millimeter-size spheres floating on the surface of a turbulent flow. The flow is generated in a layer of liquid metal by an electromagnetic forcing. By using two magnet arrays, we are able to create one highly fluctuating flow and another, more stationary flow. In both cases, we follow the motion of hundreds of particles floating at the deformed interface of the liquid metal. We evidence the clustering of floaters by a statistical study of the local concentration of particles. Some dynamical properties of clusters are exposed. We perform spatial correlations between particle concentration and hydrodynamical quantities linked with inertial effects; with vortical motion, and with horizontal divergence (corresponding to compressibility in the surface). From comparing these correlations, we propose the so-called surface compressibility as the main clustering mechanism in our system. Hence, although floaters are not passive scalar and move on a deformed surface, the scenario is similar to the one reported for passive scalar on an almost flat free surface of a turbulent flow.

Introduction

The motion of tracers in turbulent flows has attracted a lot of attention because of its impact in pollutant dispersion in the atmosphere and in the oceans. For instance, floating garbage concentrate in large litters in the middle of the ocean [START_REF] Martinez | Floating marine debris surface drift: Convergence and accumulation toward the South Pacific subtropical gyre[END_REF][START_REF] Law | Plastic accumulation in the North Atlantic subtropical gyre[END_REF][START_REF] Maximenko | Pathways of marine debris derived from trajectories of Lagrangian drifters[END_REF][START_REF] Cozar | Plastic debris in the open ocean[END_REF]. It also plays a role in the oceanic ecosystem (plankton mixing); in cloud dynamics and rain formation; or in ocean-atmosphere mixing.

At a fundamental level, it has been conjectured that a passive scalar advected by a turbulent flow has a highly intermittent dynamics, even if the flow itself is not intermittent [START_REF] Kraichnan | Anomalous scaling of a randomly advected passive scalar[END_REF]. The dynamics of finite size particles with inertia is even more intricate because it involves memory effects and a particular time scale, the Stokes time [START_REF] Maxey | Equation of motion for a small rigid sphere in a nonuniform flow[END_REF][START_REF] Balkovsky | Intermittent distribution of inertial particles in turbulent flows[END_REF][START_REF] Qureshi | Turbulent transport of material particles: An experimental study of finite size effects[END_REF][START_REF] Xu | Motion of inertial particles with size larger than Kolmogorov scale in turbulent flows[END_REF], which is controlled by the particle size and by the mismatch between the densities of the particle and the fluid. As a consequence, inertial particles do not sample the fluid uniformly, thus the phenomenon of preferential concentration emerges (see, for instance, [START_REF] Squires | Preferential concentration of particles by turbulence[END_REF][START_REF] Monchaux | Preferential concentration of heavy particles: A Voronoï analysis[END_REF][START_REF] Bourgoin | Focus on dynamics of particles in turbulence[END_REF]).

Because of buoyancy, particles of intermediate density stay on an interface between a heavy fluid (e.g. water) and a light fluid (e.g. air): so, particles behave as floaters. Hence they experience a compressibility effect [START_REF] Cressman | Compressible Flow: Turbulence at the surface[END_REF][START_REF] Boffetta | Lagrangian tracers on a surface flow: The role of time correlations[END_REF] induced by the motion of the heavy fluid, with sources of upwelling fluid, and sinks of downwelling fluid: floaters are attracted to fluid sinks, and they are expelled from fluid sources. Clustering has been observed for fictive, pointlike, fully passive particles floating on an almost-flat free surface. This has been explained by the compressible nature of the flow of particles. In order to isolate the compressible effect, these works combine experimental measurements of the surface flow and digital tracking of point-like fictive particles [START_REF] Cressman | Eulerian and Lagrangian studies in surface flow turbulence[END_REF][START_REF] Larkin | Time evolution of a fractal distribution: particle concentrations in free-surface turbulence[END_REF][START_REF] Lovecchio | Time persistence of floating-particle clusters in free-surface turbulence[END_REF]. To some extent, the present work expands these studies (i) to stronger flows inducing surface deformation and (ii) to more realistic finite size floaters subject to inertia and capillarity.

Floaters are also transported by surface waves. Stokes drift is responsible in the case of traveling waves [START_REF] Herterich | The horizontal diffusion of tracers by surface waves[END_REF][START_REF] Van Den Broeck | Stokes' drift: An exact result[END_REF][START_REF] Santamaria | Stokes drift for inertial particles transported by water waves[END_REF]. For standing waves, on the other hand, gentle transport of finite size floaters is observed in small scale experiments, where floaters eventually agglomerate [START_REF] Falkovich | Floater clustering in a standing wave[END_REF]. It is argued that inertia and surface tension play a fundamental role [START_REF] Falkovich | Floater clustering in a standing wave[END_REF], and the covering fraction of the surface is also relevant to trigger more complex collective scenarios [START_REF] Sanlı | From antinode clusters to node clusters: The concentration-dependent transition of floaters on a standing Faraday wave[END_REF][START_REF] Sanlı | Collective motion of macroscopic spheres floating on capillary ripples: Dynamic heterogeneity and dynamic criticality[END_REF]. Periodic motion in bounded domains can also induce mass transport as a consequence of viscous boundary layers. In acoustics, this phenomenon is well-known as acoustic streaming [START_REF] Batchelor | Introduction to Fluid Dynamics[END_REF], and it appeared to be relevant also for parametrically excited (Faraday) waves [START_REF] Douady | Experimental study of the Faraday instability[END_REF][START_REF] Gordillo | Measurement of the velocity field in parametrically excited solitary waves[END_REF].

The dynamics of particles floating on turbulent flows is thus subjected to several physical constrains, some of them leading to clustering. The aim of our study, rather than isolate these mechanisms, is to compare their contributions to clustering, in turbulent flows with a free surface mimicking those in natural contexts (like rivers or oceans). To do so, we pursue three main objectives: (a) to create and describe the flow; (b) to study the dynamics of real nonwetting particles floating on its surface, in particular their ability to form clusters; and (c) to compare the contributions to clustering of the different physical processes involved in our setup.

To generate the flow, we use a magnetohydrodynamical (MHD) forcing [START_REF] Bondarenko | Laboratory and theoretical models of plane periodic flow[END_REF][START_REF] Sommeria | Experimental study of the two-dimensional inverse energy cascade in a square box[END_REF]. Although it is challenging to perform measurements, this type of forcing has the advantage of producing a turbulent flow with important surface deformation, as it occurs in natural flows. Strong velocity fluctuations and surface deformation come from using a thin layer of liquid metal, of strong electrical currents and of a strong inhomogeneous magnetic field. It could be noticed that such MHD flows are also relevant for industrial applications [START_REF] Davidson | Magnetohydrodynamics in materials processing[END_REF]. A second advantage of our MHD setup is the ease in changing the geometry of the imposed magnetic field, i.e. controlling the dynamical properties of the flow. While a highly fluctuating flow is created by a regular array of magnets, a random array of magnets generates a much less fluctuating flow. In Section 2, the details of the experiment and the characteristics of the flow are described, for both magnet arrays.

In Section 3, we focus on the dynamics and clustering of floaters. Some mixing properties of floaters are exposed. It appears that above a threshold in the forcing, both magnet arrays exhibit similar mixing properties. Then, clusters are clearly identified by studying the statistical properties of the area of Delaunay triangles linking nearest neighbors. We observe stronger correlations in velocity for particles belonging to a cluster.

In Section 4, we evaluate the correlation between the particle concentration and some properties of the flow. The correlations measure the contributions to clustering in our setup. We show that the most probable clustering mechanism comes from the horizontal divergence, linked to a compressible effect for particles at the surface. Effects of curvature and capillarity are also discussed. Finally, in Section 5 we give the concluding remarks.

Creation and characteristics of a free-surface turbulent flow

This section concerns the turbulent flow under study. First we present our setup and the measurement techniques. Then, we give several orders of magnitude -set by dimensional arguments -that complement the measurements. Finally we present experimental results that reveal the main features of the flow.

Setup and methods

An electric current, of density J, and a magnetic field, B, generate a Lorentz force F L inside a conducting body, with F L = J × B. This force has been used to stir conducting fluids first by Bondarenko et al. [START_REF] Bondarenko | Laboratory and theoretical models of plane periodic flow[END_REF] in order to induce two-dimensional (2D) turbulence with a well-defined forcing wavelength. To do so, they used a uniform current and magnetic strips with alternating polarity. Later, Sommeria [START_REF] Sommeria | Experimental study of the two-dimensional inverse energy cascade in a square box[END_REF] applied a strong uniform magnetic field and space-dependent distribution of current to generate an almost 2D flow and to study the transition between large-scale structures. Since then, the technique has become a common tool to study 2D turbulence [START_REF] Tabeling | Experimental study of freely decaying two-dimensional turbulence[END_REF][START_REF] Rivera | External dissipation in driven two-dimensional turbulence[END_REF], instabilities [START_REF] Tabeling | Instability of a linear array of forced vortices[END_REF][START_REF] Dauxois | Nonlinear stability of counter-rotating vortices[END_REF][START_REF] Thess | Instabilities in two-dimensional spatially periodic flows. Part I: Kolmogorov flow[END_REF], chaotic mixing in 2D flows [START_REF] Voth | Experimental measurements of stretching fields in fluid mixing[END_REF] and wave-vortex interaction [START_REF] Falcón | Wave-vortex interaction[END_REF][START_REF] Gutiérrez | Surface waves propagating on a turbulent flow[END_REF]. We adopted a similar forcing. However, we used a layer of liquid metal which allowed us to reach high density currents (up to 1.5 × 10 5 A/m 2 ) with no need of high power nor cooling (the applied voltage is less than 1 V). The use of a horizontal current and vertical dipolar magnets to force the fluid layer, generates horizontal shear, which in turn generates strong vertical vorticity. This strong forcing distorts the interface and induces a vertical velocity component.

A diagram of the experiment is shown in Fig. 1. It is performed in a plastic (isolating) rectangular container, with a maximal working area of 40 × 50 cm 2 . This container is filled with a layer of Gallinstan up to a depth of H = 1 cm. Gallinstan is a liquid alloy at room temperature, made of gallium, indium and tin. 1 It has a density ρ = 6440 kg/m 3 . A current up to 600 A is delivered with a Sorensen DHP Series Power supply, by two brass electrodes placed along the container walls. Beneath the container, we can choose between two types of inhomogeneous magnetic fields B, which are produced by two different arrays of strong permanent Neodymium-Iron magnets of 20 mm diameter, as shown in Fig. 1. One array is made with regular lines of alternating polarity and the second one with magnets placed randomly. 2 Both have a mean distance between magnets of roughly l = 40 mm. At the bottom of the container, the magnetic field just above each permanent magnet is around 1200 G. The oxidation of the Gallinstan-air interface creates a thin solid film. To prevent it, the Gallinstan surface is covered by a layer of chlorite acid solution (at concentrations lower than a percent). The acid layer is thick enough (about 10 cm) to make the Gallinstan-acid interface insensitive to the boundary condition at the top of the acid layer. The interfacial tension (or simply the surface tension) was measured 3 to be γ = 0.5 N/m. The particles we use along this study are spherical, non-wetting, of diameter d = 1 mm and of density ρ p ∼ 0.3 ρ. Therefore, they are constrained to stay on the interface, floating on the liquid metal. In order to limit the particle-particle interaction and the particle feedback on the properties of the interface, we put only around N = 200 particles (corresponding to a filling fraction of order of 8 × 10 -4 ).

We obtain first the floaters positions. To do so, we use the configuration 1 shown in Fig. 1: a camera of 2000 × 1700 pixels 2 resolution allows us to detect the 200 particles in the whole container (d = 4 pixels) at a frame rate of 50 Hz. In each image, we obtain the coordinates of the particle centers, despite the difficulties induced by the reflective nature of the liquid metal surface. From the particles coordinates, we are able both to study the dynamics of particles (see Section 3), and to get estimates of the velocity field at the interface [START_REF] Sanlı | Collective motion of macroscopic spheres floating on capillary ripples: Dynamic heterogeneity and dynamic criticality[END_REF]. We access the latter estimates by using standard particle tracking velocimetry (PTV) algorithms. The spatial resolution of PTV makes this technique suitable for our flow. In particular, other procedures as particle image velocimetry induce spatial average that smooths intense events. The PTV avoids this problem. We compute trajectories using a multi-frame predictive tracking algorithm [START_REF] Ouellette | Programs available[END_REF][START_REF] Ouellette | A quantitative study of threedimensional Lagrangian particle tracking algorithms[END_REF], which is better suited for fast motion as particle velocity is used to predict the subsequent position (see [START_REF] Ouellette | A quantitative study of threedimensional Lagrangian particle tracking algorithms[END_REF] for comparison with other methods). Fig. 2 shows 1 From the safety datasheet acc, Guideline 93/112/EC of Germatherm Medical AG, the Gallinstan is made of 68.5% of Gallium, 21.5% of indium, 10% of Tin. Its kinematic viscosity is ν = 3.73 × 10 -7 m 2 /s, its electrical conductivity σ = 3.46 × 10 6 S/m.

2 The random distribution of magnets was obtained by choosing randomly the coordinates (x, y) of each magnet, together with their polarity. However, a balance in polarity is respected.

3 The interfacial tension was determined in a complementary experiment: in a smaller container we excited Faraday waves. By measuring simultaneously the wave frequency and wavenumber, the value of the interfacial tension was obtained after fitting the dispersion relation. representative examples of tracked particles for both the random and regular array of magnets. Supplementary movies (Movie 1) and (Movie 2) complement this picture (see Appendix A). Our PTV technique has an inherent restriction: it gives the velocity of the floaters at the surface instead of the one of the fluid. In other words, finite size floaters act as a filter for very fast or very small velocity fluctuations [START_REF] Qureshi | Turbulent transport of material particles: An experimental study of finite size effects[END_REF][START_REF] Xu | Motion of inertial particles with size larger than Kolmogorov scale in turbulent flows[END_REF]. Nevertheless, the obtained velocity does provide general properties of the flow, and consequently, the technique is commonly used in experimental fluid dynamics [START_REF] Cressman | Compressible Flow: Turbulence at the surface[END_REF][START_REF] Cressman | Eulerian and Lagrangian studies in surface flow turbulence[END_REF][START_REF] Larkin | Time evolution of a fractal distribution: particle concentrations in free-surface turbulence[END_REF][START_REF] Falkovich | Floater clustering in a standing wave[END_REF][START_REF] Sanlı | From antinode clusters to node clusters: The concentration-dependent transition of floaters on a standing Faraday wave[END_REF][START_REF] Sanlı | Collective motion of macroscopic spheres floating on capillary ripples: Dynamic heterogeneity and dynamic criticality[END_REF][START_REF] Voth | Experimental measurements of stretching fields in fluid mixing[END_REF]. We measure the surface elevation along a line using a classical triangulation technique (see Fig. 1, configuration 2), i.e. tracking the displacement of diffused light spots. This is difficult since the liquid metal interface is poorly diffusive and highly reflective. Indeed, we have to use a very sensitive camera to follow the diffused light, and we had to deal with direct reflective spots suddenly saturating the camera sensor. To recover the information lost due to these spots, we record the line displacement under two opposite angles. Hence, the bright spot in one angle is not seen in the other. Then the whole line displacement can be reconstructed [START_REF] Gutiérrez | Surface waves propagating on a turbulent flow[END_REF].

Dimensionless parameters

The dimensionless Navier-Stokes Equation, driven by an electromagnetic Lorentz forcing, exhibits a natural velocity scale U o = √ JBl/ρ, which balances the advection term and the Lorentz force. Here we use the forcing length l as the characteristic length of the flow. Thus one gets the Reynolds number Re = √ JB/ρ • l 3/2 /ν. In our device we can expect U o ∼ 30 cm/s and Re = 3 ×10 4 . Such estimates give a Kolmogorov length, η = ν 3/4 /ϵ 1/4 ∼ 3 × 10 -2 mm with ϵ ∼ U 3 o /l the energy flux by unit of mass. Note that another choice for the characteristic length, e.g. the container size, would give an unrealistic velocity and Reynolds number. In a thin fluid layer, friction on the bottom plate induces velocity damping. This friction term acts at all scales and induces interaction between structures of different sizes [START_REF] Tsang | Forced-dissipative two-dimensional turbulence: A scaling regime controlled by drag[END_REF]. In a liquid metal subject to an electromagnetic forcing, the friction is concentrated at a thin magnetic boundary layer where induction phenomena focus the electric currents and the velocity gradients [START_REF] Sommeria | Experimental study of the two-dimensional inverse energy cascade in a square box[END_REF].

= Re/Ha • H/l = √ Jl/σ νB ∼ 200. Although
the system is highly nonlinear, dynamical features are far from being those of isotropic turbulence or those of 2D turbulence. This is expected since the hypotheses of these two frameworks (isotropy and bidimensionality) are broken in our system. In particular, we observe important deformation of the surface, despite the main component of the vorticity is vertical. Gravity and capillarity govern the deformation of the surface depending on the length scale. Gravity dominates at scales larger than the capillary length l c = √ γ /ρg. At smaller scales capillarity prevails. In our experiment, l c corresponds to 3 mm. Two other dimensionless numbers are relevant in free surface flows. One is the Froude number, which is the ratio of the flow velocity over the characteristic velocity of gravity waves,

Fr = U o / √ gl ∼ 0.

General features of the observed turbulent flows

As we already noticed, floaters velocity gives an estimate of the actual flow velocity at the surface, which is useful to infer some general properties of the flow. We can first evaluate the root mean square (rms) of particle velocity



⟨V 2 P ⟩. This rms velocity normalized by U o is shown in Fig. 2(a). The value of  ⟨V 2 P ⟩ evolves between 0.5U o and 0.3U o . It is expected that the actual value is smaller than U o since the estimate is built on the maximum value of the magnetic field. Moreover we measure here the velocity of the floaters that can be smaller than the one of the sustaining fluid. Above a current of 200 A, the ratio decays for both magnet arrays. This may illustrate the fact that, when the forcing is increased, a larger part of the injected power goes to the vertical velocity component, which is excluded from our measurement. Below 200 A, the ratio saturates to a constant value for the random array whereas for the regular one, the ratio decreases continuously. This can be interpreted as a stronger bidimensionalization of the random array at low driving.

To analyze the spatial statistical features of the flow, we arbitrarily define a grid: we divide the container in N p = N w × N w cells, with N w typically equal to 50 or 100 (of 0.8 × 0.8 cm 2 or 0.4×0.4 cm 2 , respectively). This averaging procedure, used in [START_REF] Sanlı | Collective motion of macroscopic spheres floating on capillary ripples: Dynamic heterogeneity and dynamic criticality[END_REF], compensates out the inhomogeneity of particles repartition: as the measurement is long enough, all cells are visited by a significant number of particles. Hence, the averaged velocity u ij reflects the time averaged velocity field of the surface, u(x, y), and it can be used to measure the energy sustained by the mean-flow over the grid. In Fig. 2(b) we compare the energy ⟨u 2 ⟩ of this timeaveraged Eulerian velocity field, to the total kinetic energy of the particles. This gives the ratio of the energy contained in the meanflow. Almost 80% of the energy is contained in the mean-flow of the random array, compared to less than 40% in the case of the regular array. Hence the former is significantly less fluctuating than the latter. We will conveniently take advantage of this difference between the fluctuation properties of both arrays. The difference between both magnet arrays is also observable on the surface elevation induced by the forcing. We study the spatial variance of the elevation h(x, y, t) along a line perpendicular to the imposed current:

⟨ h 2 ⟩(t) = ⟨(h -⟨h⟩) 2
⟩, where ⟨•⟩ stands here for an average along a line. After time averaging, h seems to follow a power law as a function of the imposed current: ⟨ h 2 ⟩ ∝ I ζ , for both magnet arrays. However, the exponent ζ is different in both cases. It is around 1.6 for the regular array and around 2.2 for the random one. Moreover, the stationary part of the elevation, ⟨ h 2 ⟩, induces around 50% of the total surface elevation with the random magnet array, whereas it is less than 20% for the regular array. Thus, the results obtained studying the elevation h are consistent with the picture obtained from the averaged velocity field.

Dynamical properties of the floaters

We now focus on the dynamics of our finite-size floaters. Despite it is mainly governed by the one of the underlying flow, it may also be influenced by inertia, buoyancy and capillarity. We consider first the diffusion and mixing properties. Then, we study the instantaneous spatial distribution of particles. To do so, we focus on the statistical properties of the Delaunay triangles linking the nearest neighbors. The distribution of triangles obtained experimentally is compared with the one obtained from a homogeneous distribution of points. We are able to quantify the level of clustering and to determine a criterion defining clusters from the discrepancy of both distributions. Finally, the properties of the particles velocity inside a cluster are explored.

Particles diffusion and mixing

The usual way to quantify the diffusion properties of particles is to study the statistical properties of the displacement, to finite size effect. The statistical properties of the single-particle displacement are insensitive from the magnet array under study, despite the difference in their fluctuations.

R(t) =  (X(t) -X (0)) 2 + (Y (t) -Y (0))
In order to get a quantitative measurement of the mixing properties of both magnets arrays, we pixelate the container (of total volume V ) in N p squares of volume v i , as before. Then we compute the relative concentration in each cell ρ i . It is estimated by counting the number of particles n i that can be found in each cell i during the experimental run. We normalized it by the mean concentration. Thus, ρ i = (n i /v i ) • (V /N). Then we use a usual tool to quantify mixing: the variance of the relative concentration

σ (ρ i ) 2 = ⟨ρ 2 i -⟨ρ i ⟩ 2 ⟩.
The smaller is the σ (ρ i ) 2 , the better the mixing. Other tools to quantify mixing include the relative Rényi entropies [START_REF] Gouillard | Etude de l'advection chaotique dans des mélangeurs à tiges, en écoulements ouverts et fermés[END_REF]:

S q = 1 1 -α log  N p  i=1 ρ α i  , (1)
which range from 0 to 1. It reaches the limit S α = 1 for the perfectly homogeneous mixing. α = 1 corresponds to the usual Shannon entropy, α = 2 is related to the correlation entropy, and higher values of α, stress higher fluctuations [START_REF] Gouillard | Etude de l'advection chaotique dans des mélangeurs à tiges, en écoulements ouverts et fermés[END_REF]. We compute σ (ρ i ) and S α for both magnet arrays (up to α = 5 for the Rényi entropy). Results are shown in Fig. 3. The upper panel shows the difference of the concentration variance between the random array and the regular array, σ (

ρ i ) 2 = σ rg (ρ i ) 2 -σ rd (ρ i ) 2 ,
at various driving current intensities. The bottom panel shows the differences of the Rényi entropies between both arrays at five successive values of α. These quantities are estimated during the 60 s of statistically stationary regimes of the experiment. Below 200 A, there is a discrepancy. It shows that the regular array performs a better mixing. This discrepancy is more important for higher values of α. This result underlines that the difference increases when higher fluctuations of the concentration emerge. However above 200 A, both magnet arrays have the same mixing properties, despite that the flow produced by the regular array fluctuates more. It should be recalled that a transition around 200 A has been already observed in the kinetic energy of the particles driven by the random magnet array. Therefore, we can conclude that the mixing properties of both flows are equivalent above 200 A, despite their different temporal fluctuations.

Clustering characterization

We are now going to focus on the instantaneous spatial distribution of the floaters. To do so, we use the Delaunay triangles linking three nearest neighbors. We borrow this tool from the study of granular packing [START_REF] Bonamy | Multiscale clustering in granular surface flows[END_REF][START_REF] Aste | Emergence of Gamma distributions in granular materials and packing models[END_REF], and from more recent studies on clustering of inertial particles in fully developed 3D turbulence [START_REF] Monchaux | Preferential concentration of heavy particles: A Voronoï analysis[END_REF]. In order to quantify floaters concentration at the surface, we compute the area of the Delaunay triangles. Such triangulations are shown in Fig. 4, for a uniform distribution of 254 points (a) and for the same number of particles tracked on a snapshot of our experiment (b). In solid state physics and granular matter, these tessellations are used to study amorphous states. In the case of a random set of points, the tessellation gives a gamma-distribution P(A) for the elementary triangles area A [START_REF] Aste | Emergence of Gamma distributions in granular materials and packing models[END_REF], with:

P(A) = b a

Γ (a)

A a-1 exp(-bA),

(2) a = ⟨A⟩ 2 /σ (A) 2 and b = ⟨A⟩/σ (A) 2 . For uncorrelated points uniformly distributed, one expects an exponential distribution with a = 1 and b = 1/⟨A⟩ [START_REF] Castiglione | Evolution of triangles in a two-dimensional turbulent flow[END_REF]. This is indeed the case shown in Fig. 4(a), excepting small deviations due to constraints imposed by the container boundaries. This uniform distribution will be used as a reference hereafter. All excess of smaller areas from this reference, can be considered as a trace of clusters of correlated particles.

More precisely, an exponent a < 1 will be the signature of this excess of smaller areas. Indeed the probability density function (PDF)

of A diverges at vanishing values. Hence the most probable value of A is 0. The exponential cut-off at large values of A is given by the parameter b. We study the areas obtained for each snapshot of the experiment. As the number of followed particles and Delaunay triangles can change slightly from time to time, we normalized the area of each triangle, A i by the mean area of triangles at each instantaneous tessellation:

A i = A i • L x • L y /N t
where N t is the number of triangles of the instantaneous tessellation. Hence ⟨A⟩ = 1. Fig. 5 shows the PDF of these normalized areas obtained from 3000 successive snapshots of an experiment performed with the random array and with a driving intensity of I = 300 A (blue asterisks). It also shows the PDF of 3000 realizations of independent successive synthetic tessellations for sets of nearly 200 points uniformly distributed. 4 For both magnet arrays and all applied currents, the distribution of triangle areas follows a Gamma distribution (without any fitting parameters once the mean and the standard deviation are given). Note that, even if in the case of the synthetic uniform distribution one gets a and b close to one (within 20% of error due to container boundary), in the case of the experimental PDF one gets a = 0.311. This value, smaller than one, is responsible of the cusp observed near 0. The smaller a, the stronger is the cusp; i.e. the larger is the excess of smaller areas. Thus, a is indeed a signature of clustering.

Fig. 6 shows the value of a as a function of the driving current I for both magnet arrays. In both cases, except for the smallest intensity, one gets a decay with 0.24 < a < 0.5. There is therefore always a strong clustering. In the case of the regular array, the decay is almost linear, and a goes from 0.5 to 0.35. For the random array we can observe different behaviors below and above 200 A, once again. Below 200 A, the decay of a, going from 0.5 to 0.35 in 150 A, is faster than for the regular array. Above 200 A, the decay rate becomes of the same order for both magnet arrays.

The areas of the tessellation follow a Gamma distribution (2) both for the experimental points and for the synthetic set of uniformly distributed points. This allows to find a simple criterion to define particles inside a cluster. We consider that a particle is in a cluster if it belongs to a triangle with an area A smaller than a where the indices r and e refer to the synthetic reference distribution and experimental distribution respectively. Using (2) and neglecting the exponential decay at large A one gets the following critical value:

critical value A c . A c is chosen such that P e (A ≤ A c ) ≥ P r (A ≤ A c ),
A c =  b a r r • Γ (a e ) b a e e • Γ (a r )  1 ae -ar (3)
(with a r and b r close but not exactly equal to 1). This criterion has been used to define the cluster shown in Fig. 4, in which A c = 0.14.

We now have the tools to study cluster properties. Indeed, after knowing which particles belong to a cluster, we can compare the velocity fluctuations inside clusters (A < A c ) with the wholesample fluctuations (for any A). Fig. 7(a) exhibits the PDF of the fluctuations of one component of the velocity of particles belonging to clusters δv y = v y -V cy around the mean velocity, V c , of the cluster. These fluctuations are compared with those of the whole sample. The PDF for particles within the clusters is narrower and its shape, with exponential tails, differs from the nearly Gaussian whole-sample PDF. The flatness, ⟨(δv y -⟨δv y ⟩) 4 ⟩/ σ (δv y ) 4 of the fluctuations around the cluster velocity is equal to 5.1 whereas the flatness of the whole sample is 3.2, close to the value expected for Gaussian variables. This shows that the floaters belonging to the same cluster have correlated velocities. The same result can be obtained from the other velocity component. One can also look for correlations in the direction of displacement.

We define θ c as the angle between the velocity of a floater belonging to a cluster and the velocity of the cluster itself. θ c can be compared with the angle distribution of particles velocity in the whole sample. The distribution of θ c for clusters is narrower than angle distribution of the whole sample (which is almost uniformly distributed), as shown in Fig. 7(b). Both results mean that, as expected, motions inside the clusters are much more coherent that the unconditioned global ones. This strengthens our definition of clusters.

Finally we check that this coherence is conserved when the forcing is increased. The ratio of the conditioned over the unconditioned rms velocities is about 0.55 in the regular array whereas it was slightly smaller for the random one (about 0.45). It means that the cluster is slightly more coherent in the second case although the difference is not significant. In both cases the ratio does not evolve more than 10% when the current is increased. The rms fluctuations of the angles θ c are of the same order for both magnet arrays and the standard deviation is about 40%-50% smaller than the unconditioned case.

Comparing contributions to clustering

Now we explore the origin of clustering. To do so, we correlate the floater concentration with the properties of the flow. As shown before, above 200 A, both magnet arrays have similar mixing and clustering properties (see Figs 3-bottom and6). However, the random array is much less fluctuating, and the time averaged quantities are more representative of the flow properties. We use hereafter only the random array to take advantage of this property. It allows us to correlate in space various quantities averaged in time, for which we have a good spatial resolution. We get the time averaged Eulerian velocity field by the procedure described in Section 2.3. In each cell, we compute the 2D velocity component u i and v i and their derivatives, together with the floaters concentration ρ i , where the index i stands for the ith cell.

First we check a clustering mechanism similar to the one responsible of the clustering of completely passive tracers at almost-flat surfaces in turbulent flows [START_REF] Cressman | Eulerian and Lagrangian studies in surface flow turbulence[END_REF][START_REF] Larkin | Time evolution of a fractal distribution: particle concentrations in free-surface turbulence[END_REF][START_REF] Lovecchio | Time persistence of floating-particle clusters in free-surface turbulence[END_REF]. Here, as the fluid is incompressible, one can write the condition ∇ • u = 0 in terms of a horizontal divergence ∇ ⊥ •u ⊥ = ∂ x u+∂ y v = -∂ z w. In a completely two-dimensional fluid the horizontal divergence vanish, whereas it is non-zero in most practical situations [START_REF] Cressman | Eulerian and Lagrangian studies in surface flow turbulence[END_REF][START_REF] Larkin | Time evolution of a fractal distribution: particle concentrations in free-surface turbulence[END_REF][START_REF] Lovecchio | Time persistence of floating-particle clusters in free-surface turbulence[END_REF]. The horizontal divergence quantifies the presence of sources or sinks of fluid (as they produce variations in the vertical velocity w, even near the free surface where w becomes small). As floaters are constrained to stay at the surface, they cumulate or disperse depending on ∇ ⊥ • u ⊥ , experiencing a compressibility effect. We evaluate the horizontal divergence -∂ z w(x, y)| z=h(t) at the surface of the flow.

We compute

q i = ∂ z w i (x i , y i )| z i =h(t) = -∂ x u i (x i , y i )| z i =h(t) -∂ y v i (x i , y i )| z i =h(t).
Then, we estimate the local correlator r i , between q i and the normalized time-averaged concentration ρ i , previously introduced:

r i = q i • ρ i σ (q i ) • σ (ρ i ) , (4)
where X means X -⟨X⟩. For a driving current of 400 A, the global correlator ⟨r i ⟩, averaged over all the cells, is equal to 0.35. The local correlator ( 4) is larger than one on 19.8% of the pixelated area. This underlines spatial coherency between r i and ρ i . Although they are not identical, q i and ρ i present similar patterns. The correlations are not large, but they are still significant, considering the noise introduced by the coarse-grained gradient and the decoherence induced by time averaging.

We now compare the previous correlations with those between ρ i and other hydrodynamic quantities. For instance, the clustering of inertial particles in 3D turbulence is governed by the divergence of the inertial particle velocities, given at first order by ∇ • V i ∼ -τ s β i [START_REF] Balkovsky | Intermittent distribution of inertial particles in turbulent flows[END_REF][START_REF] Maxey | The gravitational settling of aerosol particles in homogeneous turbulence and random flow fields[END_REF], with

β i = ∇ • (u i • ∇u i ) (5)
and τ s = d 2 (|ρ p -ρ|)/(18νρ), the Stokes time, which takes values near 0.1 s in our conditions.

The spatial distribution of β i allows to define another local correlator

s i = β i • ρ i σ (β i ) • σ (ρ i ) . (6)
The global correlator ⟨s i ⟩ is about 0.07 and the local correlator overcomes unity only on 9% of the surface. Hence, no common pattern emerges from the comparison between these quantities. Thus, in our setup it is questionable to link the clustering with the inertial effects described, for instance, in [START_REF] Squires | Preferential concentration of particles by turbulence[END_REF][START_REF] Monchaux | Preferential concentration of heavy particles: A Voronoï analysis[END_REF].

Particles concentration can also be correlated to the vertical component of the vorticity, which is the strongest component of the vorticity with our forcing. In 3D flows, it was shown that light particles (air bubbles for instance) migrate to zones of intense vorticity, as the pressure is lower there [START_REF] Douady | Direct observation of the intermittency of intense vorticity filaments in turbulence[END_REF][START_REF] Cadot | Characterization of the low-pressure filaments in a three-dimensional turbulent shear flow[END_REF]. Despite the free surface differs from the flow in the bulk, one may expect a similar scenario for floating particles.

From our measurements, the vertical vorticity Ω i , can be evaluated on the N p cells. The global correlator with the particle concentration ⟨ Ω i • ρ i ⟩/(σ (ρ i ) • σ (Ω i )) gives 0.30, which is comparable, but smaller than ⟨r i ⟩ = 0.35, the global correlation with the horizontal divergence. High correlation between Ω i and ρ i implies that particles concentrate at regions of high vorticity, in opposition with the scenario proposed in [START_REF] Douady | Direct observation of the intermittency of intense vorticity filaments in turbulence[END_REF][START_REF] Cadot | Characterization of the low-pressure filaments in a three-dimensional turbulent shear flow[END_REF]. However, the global correlator, ⟨t i ⟩ = ⟨ Ω i • q i ⟩/(σ (q i ) • σ (Ω i )) between the vorticity and the horizontal divergence is 0.83, which is large. This strong correlation (much larger than the one between Ω i and ρ i ) may be explained by secondary flows that are induced in a shallow fluid layer around vertical vortices [START_REF] Kamp | Strain-vorticity induced secondary motion in shallow flows[END_REF]: upwelling flows merge at the vortex core whereas downwelling flows dive at the vortex edge. Therefore, we propose that the apparent correlation between vertical vorticity and the particles concentration is only the result of the correlation between Ω i and q i . In a complementary experiment, where secondary flows are absent, light nonwetting floaters move toward the axis of rotation. This experiment, where particles float at the parabolic surface of a fluid in solid-body rotation, will be the subject of a future work.

To summarize, the hierarchy of correlations we have computed suggests that clustering in our experiment is driven by horizontal divergence. To further discuss these effects, one can introduce a dimensionless factor [START_REF] Boffetta | Lagrangian tracers on a surface flow: The role of time correlations[END_REF][START_REF] Cressman | Eulerian and Lagrangian studies in surface flow turbulence[END_REF]

C = ⟨(∂ x u + ∂ y v) 2 ⟩ ⟨(∂ x u) 2 ⟩ + ⟨(∂ y u) 2 ⟩ + ⟨(∂ x v) 2 ⟩ + ⟨(∂ y v) 2 ⟩
that quantifies the degree of compressibility. It is zero for incompressible 2D flows and it takes values close to 0.5 near the surface, for free surface flows [START_REF] Boffetta | Lagrangian tracers on a surface flow: The role of time correlations[END_REF][START_REF] Cressman | Eulerian and Lagrangian studies in surface flow turbulence[END_REF]. We compute C for the random array of magnets, used to study correlations. It gives a value close to 1/6 for every forcing current I in our experiment, although the values are slightly larger when I is less than 200 A. The value 1/6 is low compared with those obtained previously near the surface [START_REF] Cressman | Eulerian and Lagrangian studies in surface flow turbulence[END_REF][START_REF] Lovecchio | Upscale energy transfer and flow topology in free-surface turbulence[END_REF]. It is not still clear for us if the coarsegraining process used to get spatially resolved fields reduces the C coefficient. Another explanation comes from the fact that we are only measuring a projection of the three dimensional deformation of the free surface (i.e. neglecting the vertical velocity). We can sustain this idea by interpreting the 1/6 as the C of the 2D projection of a 3D (homogeneous and isotropic turbulent) flow [START_REF] Boffetta | Lagrangian tracers on a surface flow: The role of time correlations[END_REF]. Our measurement giving C < 0.5 does not allow us to see any of the extreme events suggested in [START_REF] Lovecchio | Upscale energy transfer and flow topology in free-surface turbulence[END_REF].

Apart from the aforementioned physical mechanisms, our particles are sensitive to the surface tension because they are smaller than the capillarity length l c . Capillarity makes attractive particles of similar wetting [START_REF] Vella | The ''Cheerios effect[END_REF]. The attractive capillarity force between particles decays exponentially with the distance [START_REF] Dalbe | Aggregation of frictional particles due to capillary attraction[END_REF], and thus, it is significant on a characteristic length of order l c . This interaction length is an order of magnitude smaller than the mean free path between floaters in our experiment. Hence, due to the low filling fraction of particles, we expect that capillarity will be initially inefficient to agglomerate floaters. However, once the clusters are formed, the attraction could play a stabilizing role and could affect the cluster cohesion. Capillarity force also makes particles sensitive to the local curvature [START_REF] Falkovich | Floater clustering in a standing wave[END_REF]. However, our experiments with particles floating on the parabolic surface of a rotating fluid (where capillarity seems to play a crucial role), show that the motion is much slower than the one observed in the present experiment. A systematic study of the size and wetting properties of floating objects is a very important one [START_REF] Qureshi | Turbulent transport of material particles: An experimental study of finite size effects[END_REF][START_REF] Xu | Motion of inertial particles with size larger than Kolmogorov scale in turbulent flows[END_REF][START_REF] Falkovich | Floater clustering in a standing wave[END_REF][START_REF] Sanlı | From antinode clusters to node clusters: The concentration-dependent transition of floaters on a standing Faraday wave[END_REF][START_REF] Sanlı | Collective motion of macroscopic spheres floating on capillary ripples: Dynamic heterogeneity and dynamic criticality[END_REF][START_REF] Vella | The ''Cheerios effect[END_REF][START_REF] Dalbe | Aggregation of frictional particles due to capillary attraction[END_REF], although it is out of the scope of this work 5 .

A last issue is how spatial inhomogeneity of the floaters reflects the intermittent properties of the underlying flow. The high intermittency of the passive scalar stretched and folded by the velocity gradient is revealed by the anomalous scaling of the structure function of the concentration field [START_REF] Kraichnan | Anomalous scaling of a randomly advected passive scalar[END_REF]. It has been related to the ramp and cliff inhomogeneous structure of the passive scalar concentration [START_REF] Celani | Universality and saturation of intermittency in passive scalar turbulence[END_REF]. We do not reach enough resolution to compute the structure function of the concentration field of floaters in our experiment. However, it has been shown that some properties of turbulent flow, or others complex stretched flows, are enclosed in the time evolution of the shape of triangles in 2D (or tetrads in 3D) delimited by Lagrangian points passively advected by the flow [START_REF] Castiglione | Evolution of triangles in a two-dimensional turbulent flow[END_REF][START_REF] Pumir | Geometry of Lagrangian dispersion in turbulence[END_REF][START_REF] Merrifield | Scale-dependent statistical geometry in two-dimensional flow[END_REF]. Therefore a forthcoming work will be devoted to the study of the time evolution of such distorted triangles, delimited by floaters, in order to underline discrepancies with the passive scalar case [START_REF] Castiglione | Evolution of triangles in a two-dimensional turbulent flow[END_REF][START_REF] Merrifield | Scale-dependent statistical geometry in two-dimensional flow[END_REF]. Moreover, with the tool introduced to define particles belonging to a cluster, we should be able to study the triangle distortion evolution in relation to the particles ability to enter or escape from the clusters and thus to relate spatial inhomogeneity and intermittent properties.

Conclusions

In the first part of this article, we presented an experimental setup allowing us to generate strongly fluctuating free surface flows of liquid metal, thanks to an electromagnetic forcing. We apply to the fluid an electrical current I going from 25 A to 600 A. We used two kinds of magnet arrays, one regular and another random, to drive the flow. The two magnet arrays produce a different level of flow stationarity. Indeed, the ratio between the energy contained in temporal fluctuations and in the time averaged mean flow, is much larger with the regular array. This result is revealed both by velocity-field estimates and by the surface deformation. The random array, which produces a larger mean flow, exhibits several clues of a transition around 200 A: changes appear in the kinetic energy of particles, in the characteristic correlation time and in mixing properties. However, besides fluctuation properties of flows, both magnet arrays behave similarly above this value.

5 Although all the results presented here correspond to particles with a diameter of 1 mm, we also used larger particles (3 mm) to verify the reproducibility of the velocity measurements.

After presenting the flows, we focused on the dynamics of floaters. We show, in particular, that they do not mix uniformly but tend to form clusters, independent of the magnet array. This was identified by the study of the distribution of the Delaunay triangles areas. The statistics of the area of the triangle linking neighbor particles, follow Gamma distributions. These distributions are singular near zero, illustrating the tendency to form clusters, where the areas are very small. By comparing the triangles-area singularity with a uniform reference, a criterion defining clustered particles is obtained. Particles belonging to the same cluster have a coherent displacement.

To identify the main clustering mechanism, we study the correlations between surface concentration of the floaters and hydrodynamical quantities. One is linked with inertial effects; other with vortical motion; and other with horizontal divergence, that corresponds to compressibility in the surface. Knowing that both magnet arrays have similar mixing and clustering behaviors, this analysis is performed with the random array. By doing so, we benefit from stationarity, allowing us to consider time averaged quantities. Correlations suggest that the main clustering mechanism comes from the horizontal divergence in the surface, which induces a compressible effect in the floaters: they are expelled from the upwelling secondary flow at the vortex core and stretched at the vortex edge.

Fig. 1 .

 1 Fig. 1. (Color online) The experimental device. Top: diagram of the experiment and measurement configurations. A 1 cm layer of Gallinstan (GaIn Ti) is placed between two electrodes (E), over a magnet array (M.A.) in a container of 50 × 40 cm 2 . In configuration 1, the camera (C) records the position of particles floating on the Gallinstan. In configuration 2, the beam of the Laser diode (La) is transformed into a laser sheet by the cylindrical lens (Le) and projected on the surface of a mirror (M1). The diffused line, made by the laser sheet on the Gallinstan surface, is tracked with 2 opposite angles by a single Camera (C) by the two mirrors (M2) and the prisms (P). Bottom: Sketches of the random and regular magnet arrays used in the experiment. Black and white indicate the polarity of the magnets.

  2. Asit is not too far from unity, gravity waves generation cannot be completely discarded. The second number is the Weber number We = ρU 2 o L/γ (L is a characteristic length). It compares the kinetic energy of the flow and the surface tension energy. At the scale of the forcing L = l > l c surface tension is negligible (We ∼ 50). At the floater scale L = d < l c , surface tension cannot be neglected anymore (We ∼ 1). Thus our millimeter-size floaters are sensitive to capillarity.

Fig. 2 .

 2 Fig. 2. General features of the flows. (a) shows the evolution of RMS velocities of the floaters normalized by U o , versus the driving intensity I. Blue dots correspond to the random magnet array and red squares to the regular magnet array. (b) shows the ratio of the kinetic energy contained in a time-averaged Eulerian flow, divided by the total kinetic energy as a function of I. Two spatial resolutions of the Eulerian field are given for each magnets array: the random magnets array with a resolution of [50 × 50] (blue asterisks) and [100 × 100] (blue dots) and the regular magnet array at a resolution of [50 × 50] (red diamonds) and [100 × 100] (red squares). Error bars are estimated from the fluctuations in time and redundant measurements. Figures (c) and (d) show examples of tracked particles for the random and regular arrays of magnets, respectively. In both panels, dots trace the positions of all the particles recorded during 5 s at intervals of 0.1 s, for a forcing current of I = 250 A. See the motion of particles in the supplementary movies, for the random (Movie 1) and regular (Movie 2) array of magnets (see Appendix A). (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

Fig. 3 .

 3 Fig. 3. (Color online) Mixing properties as a function of I, the applied driving current. Top: Difference between the particles concentration variance of the regular magnets array and the random magnets array, as a function of I. Bottom: The difference between the Rényi entropies S α of the regular magnet array and the random magnet array, as a function of I. α = 1, * ; α = 2, •; α = 3, ; α = 4, ; q = 5, △.

Fig. 4 .Fig. 5 .

 45 Fig. 4. Position of particles and their corresponding Delaunay tessellation. (a) is for 254 uniformly distributed particles. (b) is for a similar number of particles obtained from a single experimental snapshot, with the random array at I = 300 A. In this last case, color marked points correspond to particles found in different clusters, where triangles areas are smaller than a threshold A c = 0.14 (see text). (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

Fig. 6 .

 6 Fig. 6. (Color online) Exponent a. The values of a are obtained by adjusting the PDF of the Delaunay triangles area with a Gamma distribution. It is done for both magnet arrays, as a function of I. Error bars are estimated from the accuracy of the adjustment.

Fig. 7 .

 7 Fig. 7. Statistics of the clusters dynamics. (a) PDF of the fluctuations around the mean velocity of a cluster for the y-component of the velocity of particles belonging to the cluster (red) compared to the unconditional fluctuations (blue). (b) Angle distribution of the direction of the particles displacement around the direction of the cluster, for particles belonging to a cluster (red). It is compared to the unconditioned distribution of directions (blue) and to the uniform distribution of angles (dashed line). (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

  2 , for all single particles that we are able to follow during a time t. It can be compared with the well-known Brownian motion. Within this framework one has ⟨R(τ ) 2 ⟩ = Dτ . The diffusion coefficient D is given by the Einstein formula: D = 2⟨V 2 p ⟩d 2 /(18ν), where ν is the fluid viscosity and d is the particle diameter. Unfortunately in our experiment, the range of scales between the forcing and the container size is too small to prevent finite size effects. We observe mainly a ballistic

		 ⟨V 2 p ⟩ as characteristic time unit. Above
	τ •	 ⟨V 2

transport until a time τ ∼ 10 • l/  ⟨V 2 p ⟩. It corresponds mostly to the mean time necessary to reach the border after a ballistic flight. All the displacement ⟨R(τ ) 2 ⟩ collapse in a single curve in the ballistic regime if we use l/ p ⟩/l ∼ 10 the displacement curve bends, probably due

P. Gutiérrez, S. Aumaître / European Journal of Mechanics B/Fluids 60 (2016)[START_REF] Batchelor | Introduction to Fluid Dynamics[END_REF][START_REF] Douady | Experimental study of the Faraday instability[END_REF][START_REF] Gordillo | Measurement of the velocity field in parametrically excited solitary waves[END_REF][START_REF] Bondarenko | Laboratory and theoretical models of plane periodic flow[END_REF][START_REF] Sommeria | Experimental study of the two-dimensional inverse energy cascade in a square box[END_REF][START_REF] Davidson | Magnetohydrodynamics in materials processing[END_REF][START_REF] Tabeling | Experimental study of freely decaying two-dimensional turbulence[END_REF][START_REF] Rivera | External dissipation in driven two-dimensional turbulence[END_REF][START_REF] Tabeling | Instability of a linear array of forced vortices[END_REF]

For each snapshot of the experiment, we computed an independent set of n ud uniformly distributed points. If the number of particles found in the snapshot varies, n ud varies as well.

Acknowledgments

We would like to thank V. Padilla for building the setup, C. Wiertel-Gasquet for helping us with the automation of the experiment, R. Candelier for the algorithm used to define clusters, F. Daviaud, M. Bonetti and G. Zalczer for helpful discussions. We are indebted to A. Heyries and D. Heyries for their work on the first stage of the experiment. We thank M. Hjeltman, C. Arratia and L. Gordillo for their careful reading of the manuscript. This work is supported by the ANR Turbulon (12-BS04-0005). PG also received support from the Triangle de la Physique (2009-050T Magnétonde) and CONICYT/FONDECYT postdoctorado N o 3140550.

Appendix A. Supplementary data

Supplementary material related to this article can be found online at http://dx.doi.org/10.1016/j.euromechflu.2016.06.009.