Capillary surfers: wave-driven particles at a fluid interface
CCapillary surfers:wave-driven particles at a fluid interface
Ian Ho , † , Giuseppe Pucci , , † , Anand U. Oza & Daniel M. Harris , ∗ February 24, 2021 School of Engineering, Brown University, 184 Hope Street, Providence, Rhode Island 02912, USA Univ Rennes, CNRS, IPR (Institut de Physique de Rennes)—UMR 6251, F-35000 Rennes, France Department of Mathematical Sciences & Center for Applied Mathematics and Statistics, New Jersey Insti-tute of Technology, Newark, New Jersey 07102, USA ∗ corresponding author, † co-first author Active systems have recently attracted considerable interest for the possibility of extendingstatistical mechanics to incorporate non-equilibrium phenomena [1, 2]. Among active systems,self-propelled particles are natural or artificial objects that convert energy from the environ-ment into directed motion [3]. Vibrating platforms are suitable sources of diffuse energy toobserve macroscopic self-propelled particles [4, 5] and their collective behaviour [6–10] withtable-top experiments. Here we introduce capillary surfers: solid particles propelled by theirself-generated wave field on the surface of a vibrating liquid bath. The surfer speed and its in-teraction with the environment can be tuned through the particle shape, mass, fluid propertiesand the vibration parameters. The wave nature of the interactions among surfers allows formultistability of states and promises a number of novel collective behaviours. Capillary surfersbenefit from the interplay between dissipation and inertia manifesting as a field of unsteadysurface waves, and thus have the potential to fill the gap [11] between dissipation- [2, 12–15]and inertia-dominated active systems [16–18].
Self-propulsion of artificial particles has been achieved from nanometric to centimetric scale with a numberof different strategies [3]. At the macroscopic scale, artificial self-propelled particles can be obtained on table-top experiments via the rectification of mechanical vibrations. Solid asymmetric particles can self-propel withan internal source of energy [19] or by converting energy from the vibration of the platform on which theymove [4,20]. On the other hand, fluid interfaces are convenient platforms for the self-propulsion of natural [21]and artificial bodies, including solid particles [22–25], drops [5, 26, 27] and small-scale robots [28]. When the1 a r X i v : . [ phy s i c s . f l u - dyn ] F e b igure 1: A capillary surfer self-propels on a fluid interface due to its self-generated waves (SupplementaryVideo 1). (a) Surfer geometry. (b) Side view schematic of the experimental setup (not to scale). (c) Obliquewave field visualization, in which colours are obtained from the distorted reflection of a yellow and bluebackground on the fluid surface. (d) Dependence of the surfer speed U on the forcing frequency f and forcingacceleration γ , as collapsed by Eq. (1), with L = 4 . ± .
03 mm, w = 2 . ± .
03 mm, a = L/
2, theparameters being the same in (c). Inset: dependence of U on the asymmetry a/L for f = 100 Hz, γ = 3 . L = 6 . ± .
02 mm, w = 4 . ± .
02 mm. We use h + = 1 .
22 mm and h − = 0 .
82 mm throughout allexperiments.interface is vibrating, millimetric droplets may self-propel by bouncing on the sloped surface wave theygenerate in certain parameter regimes [5], and floating drops may be deformed and driven by Faraday surfacewaves [26, 27]. For floating drops, propulsion mechanisms based on asymmetric vortex generation [27] andwave radiation pressure [26] have been proposed, but the drops are difficult to manipulate directly due totheir self-adaptive nature. The capillary surfers we introduce here are highly tunable and propel as a resultof the asymmetric radiation pressure of the waves they generate on the vibrating liquid surface.In our experiments, capillary surfers had the rectangular shape represented in Fig. 1a. The centre ofmass was offset with respect to their in-plane geometrical centre. Surfers were deposited on a bath of water-glycerol mixture and were supported at the liquid-air interface by virtue of the equilibrium between theirweight, hydrostatic forces and surface tension. The contact line of the bath was pinned to the surfer’s baseperimeter. As a result of their mass asymmetry, surfers were slightly tilted in equilibrium (Fig. 1b) andthe deformation of the interface varied along their perimeter. The liquid bath was vertically driven by anelectromagnetic shaker with acceleration Γ( t ) = γ cos ωt , with ω = 2 πf and f the forcing frequency in therange 40–100 Hz. (See Methods for more details on experimental materials and procedures.)As soon as the bath was set into vibration, a surfer generated propagating surface waves as a result of the2elative vertical motion between the surfer and bath, a consequence of the surfer’s inertia. Correspondingly,the surfer moved along its long axis in the direction of its thinner half (Fig. 1b–d). In the absence of externalperturbations or manufacturing imperfections, the surfer moved with constant velocity U along a rectilineartrajectory. Taking inspiration from marine terminology, in the following we refer to the front and back of thesurfer as the “bow” and “stern”, respectively. We proceeded by analysing the motion of a single surfer as afunction of the forcing parameters, the surfer asymmetry, and its size (Fig. 1d and Supplementary FiguresS1 and S2). For a given surfer size and asymmetry and a fixed forcing frequency f , the surfer speed increaseswith the forcing acceleration γ . For a fixed acceleration amplitude, the speed decreases with frequency.In order to rationalise this behaviour, we consider the wave radiation stress generated by the surfer. Theradiation stress S of surface waves may be defined as “the excess flow of momentum due to the presence ofthe waves” [29]. For capillary waves of amplitude A and wavenumber k in the deep-water limit kH (cid:29) S = (3 / σA k , where H is the bath depth and σ the surface tension. Experimental measurements of thedependence of the surfer speed on its asymmetry show that the speed is maximised when the displacement ofthe centre of mass is maximised, which corresponds to the highest difference in the equilibrium deformationof the liquid surface between bow and stern (inset of Fig. 1d). The radiation stress generated by the surferthus exhibits a fore-aft asymmetry: waves of larger amplitude A + are generated at the stern, where theeffective mass is larger, while waves of smaller amplitude A − are generated at the bow. As a result, thesurfer experiences a net propulsive force F p = (3 / σk w ( A − A − ), where w is the surfer width. Thispropulsive force is resisted by viscous shear beneath the surfer. By approximating this flow as a locallyfully-developed Couette flow, the viscous force on a surfer with speed U is F D ≈ ηwLU/H . The surfer speedis ultimately prescribed by the balance of F p and F D , which yields U = 3 Hσk ηL (cid:0) A − A − (cid:1) = 3 Hσ ηL (cid:18) kγω (cid:19) (cid:0) q − q − (cid:1) , (1)where we assume that the wave amplitudes A ± are proportional to the bath forcing amplitude γ/ω withunknown constants of proportionality q ± . Fig. 1d shows that the experimental data may be collapsed usingthis scaling relationship (see raw data in Supplementary Figure S1), confirming the hypothesis of radiationstress as the mechanism underlying surfer propulsion. The discrepancy observed at large forcing amplitudemay be ascribed to the influence of nonlinear effects on the wave profile that arise for large wave amplitudes,since the expression for the radiation stress we used is derived for linear waves [29].Surfers may interact with each other through the propagating waves they generate on the fluid interface.We observed that two surfers of equal size and speed can arrange into up to seven qualitatively distinctmodes of interaction (Fig. 2 and Supplementary Video 2): head-to-head (c), back-to-back (d), orbiting (e),promenading (f), tailgating (g), t-bone (h) and jackknife (i). In the head-to-head (back-to-back) mode thetwo surfer sterns (bows) face each other (Fig. 2c,d). While these modes are static, the remaining five modesare dynamic. In the orbiting mode, the two surfers orbit around the system’s centre of mass (Fig. 2e). Inthe promenade mode, they proceed side by side with constant speed along a rectilinear trajectory (Fig. 2f).In the tailgating mode, they are aligned along their major axis, with the stern of one surfer pointing toward3igure 2: The seven interaction modes of surfer pairs (see also Supplementary Video 2). (a,b) Side view (a)and top view (b) schematics of the theoretical model, in which a surfer is represented by two unequal pointmasses m + and m − . (c–i) Modes observed in experiment (top) and theory (bottom), for f = 100 Hz and γ = 3 . L = 4 . ± .
03 mm, w = 2 . ± .
03 mm, a = L/
2. Parameters used for the theoreticalpredictions are detailed in the Methods. 4he bow of the other, and they move with constant speed along a rectilinear trajectory (Fig. 2g). In thet-bone mode, the two major axes are perpendicular to each other and the stern of one surfer points towardthe bow of the other, while they both follow a circular trajectory (Fig. 2h). The jackknife mode has a similarconfiguration except the two bows face each other (Fig. 2i). Note that should the vibration be eliminated,the surfers cease to propel and immediately collapse in towards each other under the influence of capillaryattraction.These interaction modes can be captured by a conceptually simple theoretical model for the surfers’horizontal dynamics. We model a surfer of mass m as a pair of point masses m + and m − , where m + > m − and m + + m − = m , connected by a massless rod of length l = L/ i th surfer by its centerof mass xxx i ( t ) ∈ R and orientation (unit) vector nnn i ( t ), which points from m + to m − . The masses m ± arechosen to represent the surfer’s asymmetric mass distribution in experiments (Fig. 1a), and are located at xxx i, ± = xxx i ∓ µ ∓ lnnn i where µ ± = m ± /m . The masses are assumed to oscillate vertically at the forcing frequency ω and amplitude γ/ω of the bath, and thus generate capillary waves of wavelength λ = 2 π/k given by thedispersion relation in the deep-water limit, ω = σk /ρ . Since the waves are of small amplitude, γ/ω (cid:28) λ ,their form may be deduced by solving the linearized hydrodynamic problem of a periodically oscillating pointforce acting on a fluid interface [30] (see Methods). Each surfer is a pair of capillary wave sources, whichmediate the interaction between surfers.Each surfer moves in response to three forces: a propulsive force due to radiation pressure F p nnn i , a dragforce proportional to its velocity ˙ xxx i due to the viscous shear stress under the surfer and a lateral capillarywave force proportional to the local gradient of the interfacial deformation. Specifically, each mass on a givensurfer experiences a force due to the linear superposition of capillary waves generated by all other surfers,as shown schematically in Fig. 2b. Both the surfers and capillary waves oscillate at the forcing frequency ω ,and the capillary interaction force is nonzero when time-averaged over the forcing period [30] (see Methods).Crucially, this force has an oscillatory spatial dependence Φ( kr ), oscillating between attractive and repulsiveon the capillary wavelength λ . The equations of motion for a collection of identical interacting surfers arethus m ¨ xxx i + mτ v ˙ xxx i = F p nnn i + F c (cid:88) α,β = ± µ α µ β (cid:88) j (cid:54) = i Φ( k | xxx j,β − xxx i,α | ) xxx j,β − xxx i,α | xxx j,β − xxx i,α | ,I ¨ θ i + Iτ v ˙ θ i = − lF c (cid:88) α,β = ± µ α µ β (cid:88) j (cid:54) = i αµ − α Φ( k | xxx j,β − xxx i,α | ) nnn i × xxx j,β − xxx i,α | xxx j,β − xxx i,α | · zzz, (2)where xxx i,α = xxx i − αµ − α lnnn i , I is the moment of inertia, F c = ( mg ) k/σ is the capillary force coefficient and τ v = mH/ηwL is the viscous timescale. The first and second equations account for the lateral force andtorque balances on each surfer, respectively, as detailed in the Methods. The trajectory equations contain asingle unknown parameter F p , whose value F p = mU/τ v is directly inferred from the experimentally measuredfree speed U of a single surfer in isolation. Numerical simulations of the model (2) successfully recover theinteraction modes of surfer pairs (see Methods for details), and the theoretically predicted wave forms are5ualitatively similar to those observed in experiments (Fig. 2c–i, lower panels).In each mode the two surfers exhibit discrete equilibrium spacings d ≈ nλ where n = 1 , , . . . . Themaximum order n max achievable per mode depends on the properties of the liquid, vibration, and surfers.We focused on the promenade mode as in our experiments it exhibited the largest order observed: n max = 4(Fig. 3). First, we fixed the forcing frequency and measured the promenade spacing for increasing forcingacceleration. The spacing d increases very slightly with forcing acceleration γ (Fig. 3e), but is significantlymore sensitive to the forcing frequency f . The normalized spacing d/λ is approximately independent of f ,showing that the equilibrium spacings are quantized by the capillary wavelength in the range of forcing fre-quencies explored (Fig. 3f). Numerical simulations of the theoretical model (2) exhibit excellent quantitativeagreement with the experimental results (Fig. 3e,f).Many-body experiments and simulations show that capillary surfers have potential as constituents of anovel active system (Fig. 4 and Supplementary Video 3). For instance, a many-body promenade mode hasbeen observed in experiments (Fig. 4a) and simulations (Fig. 4d). Experiments reveal a super-orbiting state ofeight surfers (Fig. 4b), and simulations show that a collection of nine surfers may execute a flocking state in anorderly lattice formation (Fig. 4c). While there have been extensive studies on overdamped active systems (i.e.bacterial suspensions) [2], mediated by viscous hydrodynamic forces that decay monotonically with distance,a collection of surfers has the peculiar feature of being characterised by wave-mediated interactions, whichresults in long-range spatially-oscillatory forces defined by alternating regions of attraction and repulsion.This feature is a consequence of fluid inertia, and responsible for the multistability of a discrete set ofinteraction states, as documented here for the promenade mode (Fig. 3). The influence of multistabilityis particularly evident in the exotic many-body mode shown in Fig. 4d: the spacing between surfers isapproximately either one or two capillary wavelengths, so the mode may be interpreted as an aggregate of n = 1 (Fig. 3a) and n = 2 (Fig. 3b) promenade modes.Novel collective behaviours at a fluid interface have been observed with camphor boats [22] and walk-ing droplets [10] previously, but only within a relatively narrow parameter space defined by the physicalconstraints of these systems. Furthermore, the richness of the two-body interaction landscape (Fig. 2) doc-umented here for the surfer system already far exceeds those documented for the aforementioned systems.Ultimately, the surfer system is a uniquely tunable and accessible experimental platform that has the poten-tial to fill the gap between active systems at low and high Reynolds numbers: the surfers used in the presentwork self-propel at intermediate Reynolds number Re = ρU L/η ≈ γ = 3 . f = 100 Hz. (e) Dependence of the dimensionlessinter-surfer spacing d/λ on the forcing acceleration γ for fixed forcing frequency f = 100 Hz. (f) Dependenceof d/λ on f just below the Faraday instability threshold. In experiments, f ranges from 50–100 Hz inincrements of 10 Hz, and the corresponding values of γ/g are 1.1, 1.5, 2.0, 2.3, 3.0 and 3.5. In the model, d is defined as the distance between the surfers’ centres of mass, while in experiments the surfers’ finite widthis accounted for by subtracting w . In experiments, L = 4 . ± .
03 mm, w = 2 . ± .
03 mm and a = L/ f = 100 Hz, γ = 3 . L = 2 . ± .
01 mm, w = 1 . ± .
01 mm, a = L/
2. (c, d) Parameters used for the theoretical predictions are detailed in theMethods.and can serve to help address outstanding questions while inspiring previously unforeseen directions.AUO acknowledges support from the Simons Foundation (Collaboration Grant for Mathematicians, AwardNo. 587006). GP acknowledges the CNRS Momentum program for its support. IH and DMH acknowledgethe support of the Brown Undergraduate Teaching and Research Award. Additionally, the authors thankProf. Roberto Zenit for use of the conical rheometer.IH, GP, and DMH designed, constructed, and performed the experiments, and analyzed the experimentaldata. AUO developed the theoretical model and numerical simulation, and analyzed the simulation data.GP, AUO, and DMH interpreted the experimental and theoretical results and wrote the paper.The authors declare that they have no competing financial interests.Correspondence and requests for materials should be addressed to Daniel M. Harris(email: daniel [email protected]).
References [1] Ramaswamy, S. The mechanics and statistics of active matter.
Annu. Rev. Condens. Matter Phys. ,323–345 (2010).[2] Marchetti, M. C. et al. Hydrodynamics of soft active matter.
Rev. Mod. Phys. , 1143 (2013).[3] Bechinger, C. et al. Active particles in complex and crowded environments.
Rev. Mod. Phys. , 45006(2016). 84] Dorbolo, S., Volfson, D., Tsimring, L. & Kudrolli, A. Dynamics of a bouncing dimer. Phys. Rev. Lett. , 044101 (2005).[5] Couder, Y., Protiere, S., Fort, E. & Boudaoud, A. Walking and orbiting droplets. Nature , 208–208(2005).[6] Narayan, V., Ramaswamy, S. & Menon, N. Long-lived giant number fluctuations in a swarming granularnematic.
Science , 105–108 (2007).[7] Aranson, I. S., Volfson, D. & Tsimring, L. S. Swirling motion in a system of vibrated elongated particles.
Phys. Rev. E , 051301 (2007).[8] Kudrolli, A., Lumay, G., Volfson, D. & Tsimring, L. S. Swarming and swirling in self-propelled polargranular rods. Phys. Rev. Lett. , 058001 (2008).[9] Deseigne, J., Dauchot, O. & Chat´e, H. Collective motion of vibrated polar disks.
Phys. Rev. Lett. ,098001 (2010).[10] S´aenz, P. J. et al.
Spin lattices of walking droplets.
Phys. Rev. Fluids , 100508 (2018).[11] Klotsa, D. As above, so below, and also in between: mesoscale active matter in fluids. Soft Matter ,8946–8950 (2019).[12] Sanchez, T., Chen, D. T. N., DeCamp, S. J., Heymann, M. & Dogic, Z. Spontaneous motion inhierarchically assembled active matter. Nature , 431–434 (2012).[13] Palacci, J., Sacanna, S., Steinberg, A. P., Pine, D. J. & Chaikin, P. M. Living crystals of light-activatedcolloidal surfers.
Science , 936–940 (2013).[14] Doostmohammadi, A., Ign´es-Mullol, J., Yeomans, J. M. & Sagu´es, F. Active nematics.
Nat. Commun. , 3246 (2018).[15] Soni, V. et al. The odd free surface flows of a colloidal chiral fluid.
Nat. Phys. , 1188–1194 (2019).[16] Attanasi, A. et al. Information transfer and behavioural inertia in starling flocks.
Nat. Phys. , 691–696(2014).[17] Oza, A. U., Ristroph, L. & Shelley, M. J. Lattices of Hydrodynamically Interacting Flapping Swimmers. Phys. Rev. X , 41024 (2019).[18] L¨owen, H. Inertial effects of self-propelled particles: From active Brownian to active Langevin motion. J. Chem. Phys. , 040901 (2020).[19] Giomi, L., Hawley-Weld, N. & Mahadevan, L. Swarming, swirling and stasis in sequestered bristle-bots.
Proc. R. Soc. A , 20120637 (2013). 920] Koumakis, N., Gnoli, A., Maggi, C., Puglisi, A. & Di Leonardo, R. Mechanism of self-propulsion in3D-printed active granular particles.
New J. Phys. , 113046 (2016).[21] Bush, J. W. M. & Hu, D. L. Walking on water: biolocomotion at the interface. Annu. Rev. Fluid Mech. , 339–369 (2006).[22] Suematsu, N. J., Nakata, S., Awazu, A. & Nishimori, H. Collective behavior of inanimate boats. Phys.Rev. E , 056210 (2010).[23] Snezhko, A., Belkin, M., Aranson, I. & Kwok, W.-K. Self-assembled magnetic surface swimmers. Phys.Rev. Lett. , 118103 (2009).[24] Chung, S. K., Ryu, K. & Cho, S. K. Electrowetting propulsion of water-floating objects.
Appl. Phys.Lett. , 014107 (2009).[25] Grosjean, G. et al. Remote control of self-assembled microswimmers.
Sci. Rep. , 16035 (2015).[26] Pucci, G., Fort, E., Amar, M. B. & Couder, Y. Mutual adaptation of a Faraday instability pattern withits flexible boundaries in floating fluid drops. Phys. Rev. Lett. , 024503 (2011).[27] Ebata, H. & Sano, M. Swimming droplets driven by a surface wave.
Sci. Rep. , 8546 (2015).[28] Chen, Y., Doshi, N., Goldberg, B., Wang, H. & Wood, R. J. Controllable water surface to underwatertransition through electrowetting in a hybrid terrestrial-aquatic microrobot. Nat. Commun. , 2495(2018).[29] Longuet-Higgins, M. S. & Stewart, R. W. Radiation stresses in water waves; a physical discussion, withapplications. Deep-Sea Research , 529–562 (1964).[30] De Corato, M. & Garbin, V. Capillary interactions between dynamically forced particles adsorbed at aplanar interface and on a bubble. J. Fluid Mech. , 71–92 (2018).10 ethods
Liquid bath
The bath was composed of a glycerol-water mixture and had diameter 10 cm. In order to precisely infer thevolume fraction of glycerol, we measured the mixture density ρ = 1 . ± . at T = 21 . ◦ C witha densitometer (DM 35 Basic, Anton Paar). Using an empirical formula [M1] the glycerol volume fractionwas determined to be 63.2%. The dynamic viscosity η = 0 . ± . · s was measured at T = 22 . ◦ Cusing a conical rheometer (ARES-G2, TA Instruments), a result that is in close agreement with the viscosity η = 0 . · s obtained at the same temperature from an empirical formula [M2].All experiments were conducted at the temperature T = 21 . ± . ◦ C, at which ρ = 1 . ± .
003 g/cm and η = 0 . ± . · s, where uncertainties are due to temperature fluctuations and estimated fromempirical formulae [M1,M2]. The value η = 0 . ± . · s was used for the scaling in equation (1) andFig. 1d. The surface tension σ = 66 . ± . T = 21 . ◦ C by using a pendant dropmethod [M3]. A graduated cylinder with resolution 1 mL was used to measure and pour 45 mL of liquid intothe bath. The resulting depth of the bath was H = 5 . ± .
06 mm. The Faraday instability threshold γ F of the vibrated bath was measured before, in between, and after sets of experiments in order to ensure thatthe physical properties of the liquid remained unchanged, and ranged in the interval γ F = 3 . Surfer manufacturing
Centimeter-sized surfers with width w = 1 . L = 2 . , which are chemically hydrophobic.An adhesive PTFE sheet with thickness 0 .
40 mm was attached to the upper surface of a PTFE sheet withthickness 0 .
82 mm. A laser cutter (Universal Laser Systems, VLS 4.60) was then used to cut out rectangularprofiles with rounded corners (with radius equal to 0 . L ) to avoid sharp corners along the contact line,which were observed to reduce the reproducibility of the surfer motion. In order to introduce an asymmetryin the surfer distribution of mass, a laser engraved line was used as a guide to cut out a portion of the upperadhered layer using a fine razor. The height of the surfer’s stern and bow were thus h + = 1 .
22 mm and h − = 0 .
82 mm, respectively.
Vibration setup
The bath container was made of acrylic and directly mounted on an electrodynamic shaker (The Modal Shop,2025E) connected to an amplifier (The Modal Shop, 2100E21). Two accelerometres (PCB Piezoelectronics,352C65) were attached at diametrically opposed ends of the container and their signals were acquired by acomputer through a data acquisition device (National Instruments, USE-6343). The reported accelerationwas the mean of the measurements of the accelerometres. A closed feedback loop was used to maintain themean acceleration at a specified target value. The acceleration difference between the accelerometres was11lso recorded at each frequency. The shaker was mounted on a isolated optical bench (Thorlabs, SDA75120)to reduce the influence of external vibrations. A monochrome USB camera (Allied Vision, Mako) with amacro lens for video acquisition was mounted above the bath and normal to its surface. In order to increasethe contrast of the video recordings, the bath’s base was constructed out of a black acrylic plate. The systemwas placed inside an acrylic box to isolate the bath from ambient air currents and contaminants.
Visualization
Tracking was performed with the surfers appearing as white objects on a dark bath on the recordings ofa camera placed vertically above the bath. The wave field was visualized using a semi-transparent mirrorplaced at 45 degrees relative to the vertical, and directly above the fluid bath. A light source with diffuserwas directed horizontally toward the mirror. A black backdrop was placed behind the mirror to minimizeexterior sources of light, and the camera was placed vertically above the mirror. An in-sync strobing effectwas obtained with the camera and triggered at a frequency that was an integer divisor of f . Out-of-syncvideos were recorded with the camera frame rate slightly detuned from f . High speed videos were recordedwith a Phantom Miro R311 color high-speed camera with oblique orientation with respect to the bath usinga reflection technique to visualize the waves [M4]. Preliminary tests
The speed of 15 surfers with L = 4 . w = 2 . γ = 3 . f = 100 Hz. The mean speed was U = 2 .
27 mm/s with standard deviation 0 .
25 mm/s.Among the 15 surfers, 6 were selected for the experiments as those exhibiting the straightest trajectories.Imperfections in the pinning of the contact line were observed to lead to curved trajectories.
Experimental procedures
Prior to each experiment, the liquid container was cleaned with ethanol, rinsed with deionized water anddried with clean compressed air. In order to avoid the effects of wave reflection or waves due to the oscillatingmeniscus, data were only recorded when the surfers were within 26.5 mm of the bath centre, where the surfervelocity was measured to be statistically independent of its proximity to the border. During single-surferspeed measurements, the fluid bath was changed for every 3 surfers tested. A linear least-squares fit wasused to calculate the time-averaged speed of the surfer from the distance traveled by the surfer as a functionof time. At f = 100 Hz the error on the time-averaged speed ranged from 0.0001–0.001 mm/s, with a 95%confidence bound for a single surfer. The standard deviation among 6 surfers ranged from 0.1–0.01 mm/s.The interaction modes of surfer pairs (Fig. 2c–i, top panels) were explored with 6 pairs, chosen among the6 surfers according to the best similarity in propulsion speed. Here the fluid bath was changed and refreshedfor every new pair. For the promenade spacing measurements at f = 100 Hz (Fig. 3e), the standard deviation12f the spacing among pairs of surfers ranged from 0.01–0.1 mm.For Fig. 1d (inset), the position of the surfer’s centre of mass was varied by varying the length of thePTFE upper layer in a surfer with L = 6 . ± .
02 mm and w = 4 . ± .
02 mm. These surfers were slightlylarger than those used elsewhere in the paper in order to ease the manufacturing process. Five differentpositions of the centre of mass were explored and six surfers were used per position.Surfers with five different sizes were used to investigate the effect of magnification. Six surfers per sizewere used. Each surfer was a scaled-up version of a sample surfer (so that the width, length and corner radiusof curvature were magnified by a constant) but with the overall thickness kept constant. Each surfer hadasymmetry parameter a = L/
2. Here the fluid bath was changed every time the experiments on two surferswere completed. These results are reported in the Supplementary Materials (Fig. S2).
Uncertainties and error propagation
The uncertainty on the scaling of the surfer speed was calculated by substituting the formula for the wavenum-ber of capillary waves in the deep-water limit, k = ( ρω /σ ) / , and ω = 2 πf into Eq. (1) and computing thesquare root of the variance formula for independent variables ∆ U = (cid:115)(cid:18) ∂U∂H (cid:19) ∆ H + (cid:18) ∂U∂σ (cid:19) ∆ σ + (cid:18) ∂U∂γ (cid:19) ∆ γ + (cid:18) ∂U∂f (cid:19) ∆ f + (cid:18) ∂U∂η (cid:19) ∆ η + (cid:18) ∂U∂ρ (cid:19) ∆ ρ + (cid:18) ∂U∂L (cid:19) ∆ L , (3) where ∆ H = 0 .
06 mm, ∆ σ = 0 . f = 0 . η = 0 . · s, ∆ ρ = 0 . and∆ L = 0.01–0.05 mm was taken as twice the standard deviation of measurements on six surfers per L . Theuncertainty in the forcing acceleration ∆ γ/γ = 0.01–0.06 was measured in the range of frequency 40–100 Hzas the largest discrepancy between the measurements of the two accelerometres. The horizontal error barsin Fig. 1d represent the resulting uncertainty ∆ U . Theoretical model
A single surfer of length L , width w , asymmetry a , stern (bow) heights h + ( h − ), mass density ρ s , and mass m (Fig. 1a) is modeled as a pair of point masses m + = aρ s wh + and m − = ( L − a ) ρ s wh − connected by arigid massless rod of length l = L/ xxx i and orientation (unit) vector nnn i of the i th surfer. The masses are located at xxx i, ± = xxx i ∓ µ ∓ lnnn i ,where µ ± = m ± /m . They move under the influence of wave forces FFF ± , time-averaged over the forcing period ω = 2 πf of the bath, and drag forces − D ± ˙ xxx i, ± , yielding the equations of motion m + (¨ xxx i − µ − l ¨ nnn i ) + D + ( ˙ xxx i − µ − l ˙ nnn i ) = FFF + , m − (¨ xxx i + µ + l ¨ nnn i ) + D − ( ˙ xxx i + µ + l ˙ nnn i ) = FFF − . (4)We now assume that D ± = m ± /τ v , where τ v = mH/ηwL is the viscous timescale obtained by computingthe shear stress due to a locally fully-developed Couette flow on the underside of the surfer. Adding the two13quations in Eq. (4), we obtain the trajectory equation for the center of mass, m ¨ xxx i + D ˙ xxx i = FFF + + FFF − , (5)where D = D + + D − . To model the rotational dynamics, we take the cross product of the first equation inEq. (4) with − µ − lnnn i , the second equation with µ + lnnn i , and add the two resulting equations: (cid:2) m − ( µ + l ) + m + ( µ − l ) (cid:3) nnn i × ¨ nnn i + 1 τ v (cid:2) m − ( µ + l ) + m + ( µ − l ) (cid:3) nnn i × ˙ nnn i = l [ µ + nnn i × FFF − − µ − nnn i × FFF + ] . (6)Writing nnn i = (cos θ i , sin θ i ), Eq. (6) reduces to I ¨ θ i + Iτ v ˙ θ i = l ( µ + nnn i × FFF − − µ − nnn i × FFF + ) · zzz, (7)where zzz is the unit vector in the (vertical) z -direction and I = m + ( µ − l ) + m − ( µ + l ) = µ + µ − ml is themoment of inertia.We now deduce an expression for the capillary interaction forces FFF ± . Consider a point mass m j withposition ( xxx j , z j ) oscillating periodically with frequency ω on the surface of an inviscid and incompressiblefluid bath with density ρ and surface tension σ . We neglect gravitational effects on the waves and thus theparametric forcing on them due to the bath vibration. We also assume the waves to be of small amplitude,so that the governing equations for the velocity potential φ j ( xxx, z, t ) and interfacial deformation h j ( xxx, t ) maybe linearized: (∆ + ∂ zz ) φ j = 0 , z ≤ ,∂ t h j = ∂ z φ j , z = 0 , − ρ∂ t φ j + σ ∆ h j + m j ¨ z j δ ( xxx − xxx j ) = 0 , z = 0 , (8)where ∆ = ∂ xx + ∂ yy is the 2D Laplacian, and z = 0 denotes the position of the undisturbed fluid interface.The first equation in Eq. (8) enforces incompressibility of the fluid, while the second and third equationsare the kinematic and dynamic free surface boundary conditions, respectively. We let ¨ z j = − g − ζ j ω cos ωt ,where ζ j is the oscillation amplitude of particle j . By generalizing the derivation presented by De Corato &Garbin [30] to account for the static (non-oscillatory) deformation of the interface, we find that the solutionto Eq. (8) is h j ( xxx, t ) = m j g πσ log | xxx − xxx j | − m j ζ j ω σ G ( k | xxx − xxx j | ) cos ωt, (9a)where G ( r ) = − H ( r ) − Y ( r ) + 2 Re [ H ( ςr ) − Y ( ςr )] , (9b) ς = e i π/ , H is the zeroth-order Struve function and Y is the zeroth-order Bessel function of the secondkind [M5]. The first term in Eq. (9a) is the static deformation of the interface induced by the particle’sweight, while the second is the dynamic deformation whose form was derived by De Corato & Garbin [30].The wave fields presented in Fig. 2c–i (lower panels), Fig. 3a–d (lower panels) and Fig. 4c,d correspond to14he linear superposition of these wave profiles over all point masses. In these plots, the first term in Eq. (9a)is neglected for the sake of facilitating visualization of the wave field’s oscillations.We proceed by deriving an expression for the lateral force FFF ij exerted by the mass m j on a mass m i with position ( xxx i , z i ), where ¨ z i = − g − ζ i ω cos ωt . The mass m i exerts a force m ¨ z i ddd on the interface, where ddd = ( −∇∇∇ h j , / (cid:112) |∇∇∇ h j | is the unit normal vector to the interface. The interface thus exerts an equaland opposite force − m ¨ z i ddd on the mass, the horizontal component of which is m ¨ z i ∇∇∇ h j in the limit of smalldeformations, |∇∇∇ h | (cid:28)
1. We thus obtain the expression
FFF ij = m i ¨ z i ∇∇∇ h j ( xxx i , t ) = − m i ( g + ζ i ω cos ωt ) ∇∇∇ h j ( xxx i , t ) . (10)Using Eq. (9) and the facts that H (cid:48) = H − and Y (cid:48) = − Y , we find that the time-average of Eq. (10) overthe oscillation period 2 π/ω is [30] (cid:104) FFF ij (cid:105) = − m i m j g kσ (cid:20) πk | xxx i − xxx j | + ξ i ξ j
24 ˜Φ( k | xxx i − xxx j | ) (cid:21) xxx i − xxx j | xxx i − xxx j | , (11)where ξ i = ζ i ω /g and ˜Φ( r ) = − G (cid:48) ( r ) = H − ( r ) − Y ( r ) − ς ( H − ( ςr ) + Y ( ςr ))] . (12)The function ˜Φ( r ) changes sign as r is progressively increased [30], so the capillary force between the massesoscillates between attractive and repulsive as the distance between them increases.We may now state the equations of motion for a collection of surfers, assuming that all of the surfers andtheir associated point masses have the same oscillation amplitude, ξ i = ξ . Combining Equations (5), (7) and(11) yields the trajectory equations (2), whereΦ( r ) = 12 πr + ξ
24 ˜Φ( r ) and F c = ( mg ) kσ . (13) Numerical Simulations
We proceed by non-dimensionalizing the trajectory equations (2) using xxx → kxxx and t → tkF p τ v /m :˜ m ¨ xxx i + ˙ xxx i = nnn i + ˜ F c (cid:88) α,β = ± µ α µ β (cid:88) j (cid:54) = i Φ( | xxx j,β − xxx i,α | ) xxx j,β − xxx i,α | xxx j,β − xxx i,α | , ˜ m ˜ l ¨ θ i + ˜ l ˙ θ i = − ˜ F c / − µ (cid:88) α,β = ± µ α µ β (cid:88) j (cid:54) = i αµ − α Φ( | xxx j,β − xxx i,α | ) nnn i × xxx j,β − xxx i,α | xxx j,β − xxx i,α | · zzz, (14)where µ is given by µ ± = 1 / ± µ , ˜ l = lk , xxx i,α = xxx i − αµ − α ˜ lnnn i , and the dimensionless parameters˜ m = kF p τ v m = kU τ v and ˜ F c = F c F p = F c τ v mU (15)are defined through the free speed U = F p τ v /m of a single surfer in isolation. We further assume that thesurfers’ oscillation amplitudes are equal to the forcing amplitude of the bath, ξ = γ/g , so the model (14)contains the single unknown parameter U that is inferred from experimental measurements of a single surfer15n isolation. Equation (14) is solved using a fourth-order explicit Runge-Kutta method in MATLAB, and theStruve functions in the expression for Φ are evaluated using the toolbox “Struve functions” developed by T.P. Theodoulidis.In accordance with the experiments, the simulations were conducted assuming a surfer with length L = 4 . w = 2 . a = L/
2, stern (bow) heights h + = 1 . h − = 0 . ρ s = 2 . and mass m = 0 .
026 g. The fluid is assumed to have dynamic viscosity η = 0 . · s, density ρ = 1 .
175 g/cm , surface tension σ = 66 mN/m and depth H = 5 mm, for which the viscoustimescale is τ v = 0 .
61 s.The interaction modes shown in Fig. 2c–i (and Supplementary Video 2) and Fig. 3a–d (lower panels) wereconducted using a forcing frequency f = 100 Hz and forcing acceleration γ/g = 3 .
3, for which the surfer freespeed is U = 1 . k = (cid:0) ρω /σ (cid:1) / = 1 . − . Theassociated dimensionless parameters are ˜ m = 2 .
2, ˜ l = 4 . µ = 0 . F c = 2 . × . For Fig. 3e, the forcingfrequency was fixed, f = 100 Hz, while the forcing acceleration γ/g = ξ was varied. The surfer free speed U was inferred from the experimental data in Supplementary Figure S1 using linear interpolation/extrapolation.For Fig. 3f both f and γ/g were varied in experiments, so linear interpolation/extrapolation was used todetermine both U and ξ in the simulations. For Fig. 4c,d (and Supplementary Video 3), the forcing frequencywas f = 100 Hz and the forcing acceleration was γ/g = 6 .
6. The corresponding values of U were extrapolatedfrom the experimental data in SI Fig. 1.The different modes presented in Fig. 2 and Fig. 4 were obtained by solving Eq. (14) with differentinitial conditions. Each curve in Fig. 3e was obtained by conducting a numerical sweep in the parameter γ ;specifically, we found a promenade mode of a given order n at the largest value of the forcing accelerationconsidered, γ = 3 . γ in steps of 0.05 g, using the final state of theprior simulation as the initial condition for the next. The curves terminate when a promenade mode couldnot be found. The same was done for Fig. 3f, except the sweep was done in the parameter f in steps of 1 Hz. References: Methods [M1] Volk, A. & K¨ahler, C. J. Density model for aqueous glycerol solutions.