Elastic spheres can walk on water

Incited by public fascination and engineering application, water-skipping of rigid stones and spheres has received considerable study. While these objects can be coaxed to ricochet, elastic spheres demonstrate superior water-skipping ability, but little is known about the effect of large material compliance on water impact physics. Here we show that upon water impact, very compliant spheres naturally assume a disk-like geometry and dynamic orientation that are favourable for water-skipping. Experiments and numerical modelling reveal that the initial spherical shape evolves as elastic waves propagate through the material. We find that the skipping dynamics are governed by the wave propagation speed and by the ratio of material shear modulus to hydrodynamic pressure. With these insights, we explain why softer spheres skip more easily than stiffer ones. Our results advance understanding of fluid-elastic body interaction during water impact, which could benefit inflatable craft modelling and, more playfully, design of elastic aquatic toys.

W ater-skipping has been studied for centuries with diverse motivations including: the ancient art of stone skipping [1][2][3][4] ; naval application [5][6][7][8][9][10] ; water-surface craft 11,12 ; and biological 13,14 and biomimetic 15 water-walking. While water ricochet of rigid objects has been well studied, the physics underlying the water impact of highly deformable elastic solids remains poorly understood 12,16,17 . Compliant bodies such as inflatable boats 12 and elastic aquatic toys 18,19 exhibit behaviour that is not readily explained within the traditional framework for rigid objects. For such elastic bodies, an understanding of the coupling between the material response and hydrodynamic loading is essential in unravelling the overall dynamics.
An object obliquely impacting a water surface with sufficient inertia will carve a cavity on the air-water interface 20 and experience a pressure-driven hydrodynamic force dependent on object velocity, geometry and orientation 4,8,[21][22][23] . Water-skipping occurs when the upward vertical component of this force is large enough to lift the object off the water surface 8 . Studies born from naval applications ranging from cannonball skipping tactics 5,8 to the dam-busting Wallis bomb 10 have revealed an upper bound on the impact angle b o (angle between the free-surface and object velocity vector) below which rigid spheres will skip on water 7,8,10,21 . Disk-shaped stones are more amenable to skipping, particularly if one orients the stone at just the right angle 2 . Further research has revealed more details regarding the oblique water impact of these and other canonical rigid body geometries [23][24][25][26] . The referenced ricochet events are dominated by inertia, with negligible contributions from viscous and surface tension forces 4,20 . In this regime, the physics of water-skipping also generalize to the water-walking ability of basilisk lizards 13,14 and some birds 14,27 , and to surface craft slamming 11 .
In this work, we investigate the skipping of deformable elastic solid spheres on water. We observe that elastic spheres can skip for impact angles nearly three times larger than predicted for rigid spheres. Experiments and numerical modelling show that the spheres deform throughout impact in response to elastic waves propagating in the material. In some cases these elastic waves actually interact with the air-water interface to create nested cavities. We determine how the deformed geometry scales with material properties and initial impact kinematics. Using an analytical model to relate deformation to the hydrodynamic lift force, we identify the mechanisms by which elastic spheres skip so readily on water. Furthermore, we compute the normal and tangential restitution coefficients and find, surprisingly, that they display analogous behaviour to liquid droplets bouncing on inclined liquid films 28 . Based on our findings about single impact events, we explain how elastic spheres are able to achieve multiple successive skips on water.

Results
Elastic sphere skipping phenomena. Prior research has reported an upper bound on the impact angle of b omax ¼ 18°= ffiffiffiffiffiffiffiffiffiffiffiffi r s =r w p below which rigid spheres (density r s ) will skip on water (density r w ) 7,8,10,21 . We have found that elastic spheres skip at much larger values of b o , raising the question of how the elastic response enables this enhanced skipping behaviour. To investigate the mechanisms underlying elastic sphere skipping, we film the water impact of custom-made elastomeric spheres with a high-speed camera viewing from the side. Rigid and elastic spheres having nearly the same radius R and density r s are shown experimentally impacting the water in Fig. 1a,b. Each sphere strikes with approximately the same speed U o and impact angle b o , but the elastic sphere has a shear modulus G that is four orders of magnitude smaller than that of the rigid sphere material. Within a few milliseconds after impact, the elastic sphere deforms dramatically and rides along the front of a cavity on the air-water interface before lifting off the surface. By the time the elastic sphere is two diameters above the surface (tE25 ms), the rigid sphere has plunged nearly the same distance below it. The elastic sphere evidently experiences a larger upward vertical force from the water.
The extreme sphere deformation is more evident in Fig. 1c, which shows that the water-contacting surface assumes the shape of a disk with a larger radius than that of the undeformed sphere. The disk is oriented at an attack angle a that, unlike for skipping stones 4 , changes in time throughout the impact. Large amplitude oscillations excited by the impact can persist in the sphere after lifting off the surface ( Fig. 1b and Supplementary Movie 1) or even while the sphere is still in contact with the water (Fig. 1d and Supplementary Movie 2). In the latter case, the sphere vibrations form a group of nested cavities, or a so-called matryoshka cavity 29 , named after Russian nesting dolls (Fig. 1d). This phenomenon is in contrast to rigid sphere skipping, for which the cavity is asymmetric, but smooth (Fig. 1e).
Sphere deformation modes. To examine the sphere deformations more thoroughly, we implement a fully coupled numerical finite-element model in Abaqus (ref. 30) (see Methods section). Figure 2a shows the results of a numerical simulation carried out with the same sphere material properties (R, G and r s ) and impact conditions (U o and b o ) as the experimental test shown in Fig. 2b. The numerics reveal an elastic wave propagating around the circumference of the sphere in a counter-clockwise direction. We classify this type of wave propagation, depicted in the line drawing of Fig. 2c, as vibration mode 1 À . In some cases this elastic wave impacts the air-water interface at time t w , thus initiating a matryoshka cavity (as seen numerically and experimentally in Fig. 2a,b). While tempting to attribute these kinematics to rigid body rotation, we find that the elastic wave propagation time t w is typically much smaller than the measured period of rigid rotation (see Methods and Supplementary Movie 2). The Abaqus numerical model predicts two additional vibration modes, generally occurring with increasing impact speed and/or decreasing shear modulus. In mode 2 (Fig. 2d, Supplementary Movie 1), the sphere assumes an ellipsoidal shape with oscillating major and minor axes. In mode 1 þ (Fig. 2e) an elastic wave again propagates around the circumference of the sphere, but in the clockwise direction. Finally, we observe that the attack angle a of the deformed water-contacting face evolves as a result of the elastic wave propagation (Fig. 2a,

f and Supplementary Movie 2).
Skip-enhancing mechanisms. Based on our experimental observations and numerical simulations (Figs 1 and 2), we hypothesize that the elastic response of the sphere enhances skipping through two mechanisms: (1) by taking on the shape of a flat disk with an increased wetted area; and (2) by acquiring a favourable attack angle.
To connect the suggested skip-enhancing mechanisms to the vertical force acting on the sphere, we propose an analytical model of the coupled fluid-structure interaction. The sphere is idealized as an incompressible, neo-Hookean hyperelastic solid 31 with shear modulus G, radius R and density r s . During water impact, the sphere deforms into a disk-shaped ellipsoid inclined at attack angle a(t) to the water surface (Fig. 3a). A set of equations can then be written for the motion of the centre of mass (COM) and sphere deformation in terms of general forces and tractions acting on the body (see Methods section).
To couple the sphere response to the fluid loading, we extend an existing hydrodynamic force model for circular disk-shaped  skipping stones 4 by approximating the water-contacting face as a circular disk with radius l eq R, where l eq approximates the sphere deformation (see Fig. 3a and Methods section). The force is modelled as where U B and b B describe the velocity of the water-contacting face and the wetted area S w is proportional to (l eq R) 2 .
The equation for the vertical motion of the COM is then where d 2 is the vertical coordinate of the COM and g is gravity.
We can now predict how the hypothesized mechanisms relate to skipping. First, a larger area S w increases the magnitude of the hydrodynamic force term in equation 2. For a neo-Hookean material, S w increases with decreasing shear modulus G for a given applied compressive stress 31 . Second, a smaller attack angle increases cos a thereby increasing the vertical force component that lifts the sphere off the surface. We note that the governing equations for the sphere deformation predict a steady-state solution in which the attack angle evolves as a t ð Þ ¼ 1 response to a circumferential wave propagating in mode 1 À (see Supplementary Note 1). While the sphere deformation during impact is not steady-state, we nonetheless expect a to be governed by the speed of elastic waves propagating in the sphere through the characteristic distance R such that a t ð Þ / 1 we predict that a more compliant sphere (smaller G) will assume a smaller rate of change of a, thereby increasing the upward vertical force that enables skipping. Numerical simulations verify the expected dependence of these two mechanisms on G. Measurements from the Abaqus results show that the rate of change of the attack angle scales as for mode 1 À deformations (Fig. 3b). Additionally, we find that the maximum value of l eq achieved during impact, which we call l max , increases with a decrease in the dimensionless term G=r w U 2 o (Fig. 3c), which is the ratio of material stiffness to hydrodynamic pressure. Therefore, a smaller G yields a larger stretch and larger wetted area, as well as a smaller rate of change of a, as predicted.
To confirm that these mechanisms indeed enhance skipping, we perform experiments and simulations over a range of impact conditions and sphere properties and measure the minimum impact speed required to skip U min , as a function of G (Fig. 4). Above a certain value of G (E10 3 À 10 4 kPa, depending on b o ), we recover the rigid skipping regime, in which U min is independent of shear modulus but is very sensitive to b o . For rigid spheres impacting above r w , prior research suggests spheres may broach (that is, become completely submerged before exiting), but not skip 8 (Fig. 4f). For stiffness values below the rigid regime, the elasticity of the sphere becomes important and U min decreases monotonically with decreasing shear modulus. Our analytical model accurately predicts the experimental and numerical results in this regime. The minimum speed is also much less sensitive to b o in the elastic skipping regime and as a result we observe skipping at impact angles nearly three times larger than predicted for rigid spheres ( Fig. 4 and Supplementary Fig. 1). While our results show that reducing shear modulus has the predicted effect on wetted area and attack angle (Fig. 3b,c), it is unclear whether one or both of these mechanisms are responsible for the observed improved skipping performance. To isolate the mechanisms we consider the limiting case of small G, for which ! 0 and thus cos aE1 over typical impact timescales.
The expected dependence of U min on G in this limit can be rationalized by scaling analysis of equation 2 (see Methods and Supplementary Note 2), which gives 2R where t c is the collision time (that is, time in which the sphere is in contact with the water). For threshold skipping cases, we expect the characteristic acceleration 2R=t 2 c to be small compared with gravity and thus, to first order, 3U 2 min l 2 max =8pRr ? % g. Furthermore, in the small G limit our numerical modelling shows that l max / G=r w U 2 min À Á À 1=3 (Fig. 3c). Applying this dependence and solving for U min , we find U min / gr ? R ð Þ 3=10 G=r w ð Þ 1=5 . Figure 4a shows that U min approaches the G 1/5 relation in the limit of small G, indicating that the only mechanism by which reducing shear modulus enhances skipping in this limit is through the increased wetted area. However, for larger G (4E10 kPa), U min deviates from the G 1/5 relation as the coupling between shear modulus and a becomes important ( Fig. 4a-d). As stiffness continues to increase, ultimately the amplitude of the deformations (affecting both S w and a) become negligible and the sphere is effectively rigid (Fig. 4a,e). Consequently, we conclude that decreasing shear modulus below this rigid boundary causes an increase in the upward vertical force that promotes skipping through both of the hypothesized mechanisms, save for the limit of small G where decreasing shear modulus only affects lift by increasing S w .
Skipping regimes. As the impact events devolve from clear skipping to water entry, we observe a transitional regime characterized by a matryoshka cavity, in which the sphere still skips ( Fig. 1d; Supplementary Movie 2). The matryoshka phenomenon occurs when the total contact time of the sphere with water t c is longer than the wave time t w associated with mode 1 À elastic wave propagation, such that t c /t w 41. We define t w as the time from impact until the circumferential elastic wave strikes the air-water interface (Fig. 2a,b). Experiments over a range of sphere properties and impact conditions reveal that t c /t w is governed by the impact angle b o and the ratio G=r w U 2 o (Fig. 5). The dependence on these terms can be rationalized by scaling analysis with U o replacing U min in equation 3. When t c Et w , the characteristic sphere acceleration is much greater than g such that 2R=t 2 c % 3U 2 o l 2 max =8pRr ? and we find t c / R ffiffiffiffiffi r ? p =U o l max (see , the elastic response of the sphere affects skipping. In this elastic regime, the minimum speed required to skip U min decreases monotonically with decreasing shear modulus G, as shown by both experimental and numerical results (triangle and star markers, respectively). In the limit of small G we expect U min pG 1/5 , which is confirmed by the numerics. As stiffness increases above GE10 kPa, U min deviates from this relation as larger G augments the rate of change of a, thereby reducing the upward vertical force component. Our analytical model also captures this change in behaviour (solid coloured lines). For shallow impact angles, U min becomes insensitive to shear modulus for G4E10 3 kPa, indicating that rigid sphere skipping behaviour is recovered. The transition between the elastic and rigid skipping regime occurs at larger G as b o increases. In the rigid skipping regime, U min is very sensitive to b o (also evident in Supplementary Fig. 1). The lower bound of the rigid skipping regime is inferred from the dark triangle symbols, which occur for experiments with b o ffiffiffiffiffi r ? p o14:5°. The coloured triangle markers at Methods section). Furthermore, based on the propagation speed of mode 1 À elastic waves, we expect t w / 1= _ a / R= ffiffiffiffiffiffiffiffiffiffi G=r s p , which is confirmed experimentally (Supplementary Fig. 2). Combining the scaling for each time gives t c =t w / G=r w U 2 o À Á 1=2 =l max . We can now examine the evolution of the timescale ratio in the vicinity of the transitional regime in the limit of shallow (b o -0) and steep impact angles, making use of the relationship between l max and G=r w U 2 o shown in Fig. 3c. Here, we consider steep impact angles to be b o ! b o EM , where b o EM is the maximum impact angle at which we have observed elastic sphere skipping b oEM % 47°À Á . In the shallow b o limit, t c / t w E1 occurs at values of G=r w U 2 o 40:1, for which l max -1 (Fig. 3c); thus, we anticipate t c =t w / G=r w U 2 o À Á 1=2 for shallow angles. For steep b o , the transitional regime occurs for G=r w U 2 o o0:1, for which l max / ðG=r w U 2 o Þ À 5=12 (Fig. 3c), and we expect t c =t w / ðG=r w U 2 o Þ 11=12 . These limiting relations capture the evolution of t c /t w observed experimentally (Fig. 5) and provide insight into the differences observed at different impact angles. We see that for steeper b o , the sphere deformation has a larger effect on the collision time, which gets manifested as a higher sensitivity of t c /t w to G=r w U 2 o . Based on our findings regarding the transitional skipping regime, we hypothesize that the same dimensionless parameters (that is, b o and G=r w U 2 o ) govern all deformation modes and associated skipping behaviour. An empirical regime diagram indeed shows that these parameters classify all observed skip types (Fig. 6). Mode 1 þ is promoted by shallow b o , large U o and/or small G. As stiffness becomes large relative to hydrodynamic pressure, the vibration type traverses the mode 2 and mode 1 À regimes. Our analytical model correctly predicts the boundary between the mode 1 À skipping (t c /t w o1) and mode 1 À transitional (t c /t w 41) regimes (marked by red line on Fig. 6).
A small, radius-dependent overlap region exists between the transitional regime and water entry. Finally, we have observed that the vibration mode associated with the transitional and water entry regimes is exclusively mode 1 À .
Skipping sphere rebound. We quantify the rebound characteristics of skipping spheres by computing the normal and tangential coefficients of restitution, e n ¼ U e n =U o n and e t ¼ U e t =U o t , where U e is the exit velocity and n and t refer to components normal and tangential to the free-surface, respectively. Here, we find some interesting similarities between skipping elastic spheres and liquid droplets bouncing on inclined liquid layers 28 . Following the work on bouncing liquid droplets, we plot e 2 n and e 2 t as a function of r w U 2 on =G (Fig. 7), which is equivalent to the normal Weber number ðWe o n ¼ r w U 2 on R d =sÞ with shear modulus G replacing the Laplace pressure s/R d , where R d is the droplet radius and s is surface tension of the liquid droplet. Interestingly, we find e 2 n / ðr w U 2 on =GÞ À 1=2 and e 2 t / ðr w U 2 on =GÞ À 1=3 , which are the same scaling relations found for liquid droplets bouncing on inclined liquid layers 28 . The restitution coefficients deviate from these trends when t c /t w 41. In these cases, sphere vibrations interact with the water via the matryoshka cavity causing a significant decrease in bouncing efficiency. For both elastic skipping spheres and bouncing droplets, the restitution coefficients are always less than one as part of the initial translational kinetic energy goes into post-impact vibrations in the sphere or droplet 28,32 .
We build on the bouncing droplet analogy to speculate on the lower bound of validity of the U min pG 1/5 scaling relation in the x-marks: water entry). The experimental data shown are for all three elastic sphere radii. The coloured regions provide a visual boundary of the skipping (blue), transitional (yellow) and entry (red) regimes. Our analytical model accurately predicts the boundary between mode 1 À skipping and transitional events; the red line is found from the intersection of the analytical results with t c /t w ¼ 1 in Fig. 5. Given initial conditions of the first skip, our analytical model also predicts successive impact types in a multiple skip trajectory (blue squares). A multiple skip experiment validates these predictions within the limits of our test facility (purple markers; R ¼ 26.2 mm, G ¼ 12.3 kPa, r* ¼ 0.937). Images from the multiple skip experiment are shown in Fig. 8. Characteristic error bars are shown in grey. limit G-0 (Fig. 4). When the relative magnitude of droplet surface tension becomes small for liquid droplets impacting liquid layers, bouncing does not occur and the droplet completely merges with the liquid layer 28 . We conjecture about a similar limit for elastic spheres with G-0 and impacting with We ¼ r w U 2 o R=s w % 1, where s w is the surface tension of water. In this limit, we expect the surface tension force from the water to act prominently on the sphere 33 and to inhibit sphere reformation, thus preventing recovery of translational kinetic energy from deformational potential energy during impact. As a result, we hypothesize that sphere skipping would ultimately cease in this limit. We anticipate these dynamics would become relevant when GEs w /R (see Supplementary Fig. 3 and Supplementary Note 3). To validate these predictions is beyond the scope of the present work.

Discussion
Perhaps the most mesmerizing manifestation of elastic sphere water impact is continual skipping across water (Supplementary Movie 3). To confirm that our physical description of a single skip generalizes to multiple skip events, we predict the placement of successive impacts on the regime diagram (Fig. 6). An experimental investigation using isolated water tanks validates the predictions and shows the sphere traversing through the vibration modes with each impact (Fig. 8). As to how repeated skipping is sustained over very long skipping trajectories (Supplementary Movie 3), we gain insight from the behaviour of the restitution coefficients (Fig. 7). First, e t is consistently larger than e n , which causes b o to decrease and thus become more favourable with each skip. Second, the restitution coefficients actually become larger as U o n decreases, until the sphere enters the mode 1 À skipping and transitional regimes. Therefore, one could say it becomes easier to skip with every skip.
While toy elastic balls may bestow upon the casual sportsman the ability to break the world stone skipping record (88 skips by K. Steiner, Guinness World Records), we believe the physics underlying the elastic sphere impact are common to the large deformation hydroelastic response of surface-riding and skipping compliant bodies. Models of these structures, such as inflatable boats, typically ignore extreme elastic deformation even though it is known to affect drag, stability and slamming loads 12 . The mechanisms of form and force augmentation, as well as the secondary vibration-induced fluid interactions that we have revealed, can be exploited for functional advantage and incorporated into higher-fidelity hydroelastic models.

Methods
Sphere fabrication and material properties. In order to control material properties, custom elastomeric spheres were fabricated from a high performance platinum-cure silicone rubber called Dragon Skin produced by Smooth-On, Inc., which consists of two liquid constituent parts. Once the two constituents are mixed, the material sets without requiring heat treatment. The shear modulus was varied by adding a silicone thinner to the mixture before setting, which reduces the material shear modulus by decreasing the polymer cross-linking density. Sphere  Experiments reveal the dependence of the square of the normal (a) and tangential (b) restitution coefficients on r w U 2 on =G (grey-filled symbols denote mode 1 À transitional skips; marker shapes same as for Fig. 6; r* ¼ 1±0.08). We find that e 2 n and e 2 t scale with r w U 2 on =G in the same way that for bouncing liquid drops 28 they scale with We on ¼ r w U 2 on R d =s. Thus, shear modulus G for skipping elastic solid spheres plays the role of the Laplace pressure s/R d for bouncing liquid droplets. Because rebound is more efficient in the tangential direction, b o decreases (thus becoming more favourable) with each skip. In addition, both restitution coefficients increase with decreasing impact speed until the sphere enters the mode 1 À skipping and transitional regimes. These two effects combine to enable multiple skip trajectories ( ARTICLE materials with three different shear moduli were fabricated by adding 0, 1/3 and 1/2 parts thinner by mass ratio. Before setting, the liquid mixture was placed in a vacuum chamber to remove any entrained air. For our experiments, spheres were fabricated by curing the liquid mixture in smooth, machined aluminium moulds to produce spheres with three different radii: 20.1 ± 0.8 mm; 26.2 ± 0.8 mm; and 48.8±0.9 mm. A rigid sphere with R ¼ 25.8±1 mm was fabricated from Nylon DuraForm PA using selective laser sintering. The uncertainty on each sphere radius represents the 95% confidence interval based on several independent measurements. A thin Lycra casing was loosely fitted around each sphere in order to prevent undesired particles from adhering to the surface and to reduce the friction between the sphere and the launching mechanism from which it was fired. The silicone rubber was so compliant that traditional uniaxial 'dogbone' testing on our Instron machine was not feasible as the forces generated were too small to be reliably measured. To overcome this, we performed a test in which the spheres were compressed on the Instron to generate a quasi-static force-displacement curve. This test set-up was then numerically modelled in Abaqus with the sphere material described by a neo-Hookean hyperelastic constitutive model, parameterized by the shear modulus G 31 . We then varied G to find the value that produced the best fit between the numerically simulated and experimentally measured force-displacement curves. The elastomeric spheres used in our experiments had shear moduli of 97.2, 28.5 and 12.3 kPa corresponding to 0, 1/3 and 1/2 parts thinner, respectively. According to the manufacturer of the rigid sphere (3D Systems-Quickparts Solutions), the elastic modulus of the selective laser sintering Nylon DuraForm PA material is 1.59 Â 10 6 kPa, which-assuming a Poisson's ratio of 0.4-gives a shear modulus of G ¼ 5.66 Â 10 5 kPa.
Sphere skipping experiments and data processing. Spheres were launched at the water surface from a variable-angle, pressure-driven cannon consisting of a pressure chamber for compressed air, a sliding cylindrical piston and interchangeable barrels. Sliding the cylindrical piston allowed air to flow from the pressure chamber into the barrel, thus forcing the sphere to accelerate out of the barrel and strike the water surface with impact speed U o and angle b o . Impact events were illuminated with diffuse white back lighting and filmed with either NAC GX-3 or Photron SA3 high-speed cameras acquiring at 1,000-2,000 frames per second (fps).
The impact speed U o and angle b o were measured from images of the sphere before water impact using a cross-correlation algorithm. The mean uncertainties on U o and b o are ±1.09 m s À 1 and ±1.75°, respectively (computed at 95% confidence level). The same algorithm was used to measure the exit speed and angle of the sphere after lifting off the surface. Also, the rigid body rotation of the sphere was estimated after skipping by tracking reference markers on the exterior of the Lycra casing. It was found for mode 1 À skip types that the rotational kinetic energy was typically o8% of the translational kinetic energy after water exit. Furthermore, for these impacts the wave time t w was typically o30% of the period of rigid body rotation measured after skipping. The collision time t c , wave time t w and vibration mode classification were all determined from manual inspection of the high-speed images. The minimum skipping speed U min was found experimentally by performing successive experiments with identical conditions but with increasing speed until the sphere skipped.
Abaqus numerical model. Details of the Abaqus numerical model of the elastic sphere impact are contained in reference 30 ; clarifications relevant for the present work are summarized here. The finite-element model uses the built in coupled euler-lagrange functionality of Abaqus/Explicit, which couples the contact region between the Lagrangian (sphere) and Eulerian (fluid) domains using a penalty method. Direct numerical simulation of the compressible Navier-Stokes equations is performed in the Eulerian domain. For the solid, conservation of momentum is solved with an incompressible neo-Hookean constitutive model describing the sphere. For all numerical model results presented herein, the mesh consists of eight-noded Eulerian hexahedral elements with spatial resolution of 3 mm. The sphere radius R, density r s and shear modulus G, as well as the impact speed U o and angle b o were set to match experimental values. The three-dimensional computational domain consists of a water tank (length ¼ 30 R, depth ¼ 6.5 R) with a symmetry plane coinciding with the plane of motion. For the numerical results presented in Fig. 4a, U min is the average of the impact speeds for a skipping case and the non-skipping case with the nearest U o . The numerical marker error bars reported in Fig. 4a represent ±1/2 of the speed difference between the two cases.
Analytical model of elastic sphere skipping. An approximate analytical approach to modelling the impact between a compliant elastomeric sphere and a fluid surface is outlined (for a complete derivation, see Supplementary Note 1). Here, we describe the sphere deformation and motion using a set of reduced, scalar generalized coordinates that are governed by a system of ordinary differential equations (ODEs). We begin by defining a fixed Cartesian coordinate system {e 1 , e 2 , e 3 } (Fig. 3a). We assume the sphere moves only in the e 1 À e 2 plane and undergoes no rigid body rotation. The sphere deformation is first described by a rigid displacement d that describes the motion of the COM in terms of the generalized coordinates d 1 , d 2 . This is followed by a volume preserving stretch V that deforms the sphere into an ellipsoid. The coordinates of a material particle before deformation (x) and after deformation (y) are related by y ¼ d þ Vx. The velocity and acceleration fields follow as The stretch V and its time derivatives _ V, € V can be written in terms of the principal stretches l 1 , l 2 and l 3 ¼ 1/l 1 l 2 , which are aligned with the body-fixed {m 1 , m 2 , m 3 } coordinates, respectively (see Fig. 3a where s is the Cauchy stress tensor, D is the stretch rate, b are body forces, t are traction forces and V and A denote integration over the volume and surface of the deformed solid, respectively. The term involving the external traction can be re-written as where F represents the resultant hydrodynamic force acting on the solid and the second term on the right-hand side represents the virtual power associated with a force dipole tending to distort the elastomer. Expressing equations 4 and 5 in terms of the generalized coordinates (see Supplementary Note 1) and then setting each of d _ d 1 , d _ d 2 , d _ l 1 , d _ l 2 , d _ a to be non-zero in turn yields a set of coupled second order nonlinear ODEs in terms of general forces and tractions: where F h and F v are the horizontal and vertical force components in {e 1 , e 2 , e 3 } coordinates, respectively, and t i and y i are the components of the traction vector and position vector on the surface of the ellipsoid in {m 1 , m 2 , m 3 } coordinates, respectively. Next, we propose an approximate model for the hydrodynamic forcing on the sphere. It is beyond the scope of this work to derive an analytical expression for the dynamic pressure distribution over the wetted sphere surface. Rather, our goal is to generate a simplified model that captures the first-order sphere motion and deformation during an oblique water impact. Thus, the hydrodynamic force F is computed using equation 1, which follows from work on skipping stones 4 by considering the deformed sphere to be a circular disk. The disk radius is described by an equivalent principal stretch l eq computed by equating the area of the equivalent circular disk with the cross-sectional area of the deformed sphere in the m 1 -m 3 plane, p(l eq R) 2 ¼ p(l 1 R)(R/l 1 l 2 ), which gives l eq ¼ 1= ffiffiffiffi ffi l 2 p . Furthermore, S w , U B and b B can be written in terms of the generalized coordinates describing the sphere deformation (see Supplementary Note 1). Without describing the pressure distribution, we cannot specify the centre of pressure and, thus, cannot define the y 1 coordinate at which the traction vector acts in equation 8, which governs the attack angle a. To overcome this, we determine a from our Abaqus numerical model (Fig. 3b). With a prescribed, the remaining ODEs for d 1 , d 2 , l 1 and l 2 (equations 6,7,9 and 10) can be solved without further simplification. Inserting the