Drop impact on a solid surface: short time self-similarity
UUnder consideration for publication in J. Fluid Mech. Drop impact on a solid surface : short timeself-similarity
J U L I E N P H I L I P P I,P I E R R E - Y V E S L A G R É E and
A R N A U D A N T K O W I A K
Sorbonne Universités, UPMC Univ Paris 06, CNRS, UMR 7190 Institut Jean Le Rondd’Alembert, F-75005 Paris, France.(Received April 23, 2015)
The early stages of drop impact onto a solid surface are considered. Detailed numericalsimulations and detailed asymptotic analysis of the process reveal a self-similar structureboth for the velocity field and the pressure field. The latter is shown to exhibit a maxi-mum not near the impact point, but rather at the contact line. The motion of the contactline is furthermore shown to exhibit a ’tank treading’ motion. These observations are ap-prehended at the light of a variant of Wagner theory for liquid impact. This frameworkoffers a simple analogy where the fluid motion within the impacting drop may be viewedas the flow induced by a flat rising expanding disk. The theoretical predictions are foundto be in very close agreement both qualitatively and quantitatively with the numeri-cal observations for about three decades in time. Interestingly the inviscid self-similarimpact pressure and velocities are shown to depend solely on the self-similar variables( r/ √ t, z/ √ t ). The structure of the boundary layer developing along the wet substrateis investigated as well, and is proven to be formally analogous to that of the boundarylayer growing in the trail of a shockwave. Interestingly, the boundary layer structure onlydepends on the impact self-similar variables. This allows to construct a seamless uniformanalytical solution encompassing both impact and viscous effects. The depiction of thedifferent dynamical fields enables to quantitatively predict observables of interest, suchas the evolution of the integral viscous shearing force and of the net normal force.
1. Introduction
The impact of a liquid drop onto a rigid surface results in a rapid sequence of eventsending, in the inertial limit, in spreading (Eggers et al. et al. (2002), that evidenced a “kinematicphase” where the drop merely resembles a truncated sphere and spreads as the square-root of time. This phase precedes the apparition of the liquid lamella.Probably one of the first depiction of the very first instants of drop impact datesback to Engel (1955). With the help of high-speed cinematography, Engel captured thechronology of events triggered by drop impact. He noted in particular the unvarying shapeof the drop apex during the earliest moments of impact, which might be surprising due tothe incompressible character of the liquid. Engel put forward the possible roles of inertia,viscosity or surface tension to explain this observation. Actually, the physical mechanismunderpinning this behaviour is best illustrated with Figure 1a. There, the numerically a r X i v : . [ phy s i c s . f l u - dyn ] A p r J. Philippi, P.-Y. Lagrée and A. Antkowiak
Figure 1.
Close-ups of increasing magnitude on the pressure field developing inside an impactingdrop in the inertial limit. The pressure field is extracted from Navier-Stokes
Gerris computationsof a drop impacting a solid surface at early times (note that the surrounding gas dynamics iscomputed as well, but not represented here). Noticeably the motion is essentially pressureless(and therefore corresponds to a free fall) except in a concentrated region in the contact zone.The successive close-ups on pressure field structure in the contact region reveal a pressure peaknear the contact line (the physical parameters are here Re = 5000, We = 250, tU/R =4 × − .The total size of the numerical axisymmetric domain is 2 R × R , and the adaptive mesh haslocally a mesh density corresponding to 32768 × computed pressure field within an impacting drop is represented shortly after impact(details to follow in the paper). It is readily seen that the structure of the pressure field isextremely concentrated near the contact zone, as in Hertz’ classic elastic contact problem.Conversely the outer region is essentially pressureless. This strong inhomegeneity in thepressure distribution therefore explains why, in the absence of any pressure hindrance,the upper part of the drop freely falls even after impact while remaining undeformed.The pressure concentration in the early stages of impact was first identified by Josserand& Zaleski (2003). From the key remark that the extent of the pressure concentration zonescales with the contact radius, these authors conjectured a self-similar structure for thepressure field and evolution with time as 1 / √ t – an hypothesis comforted by numericalresults. Though sufficient to detect hints of self-similarity, numerical simulations were rop impact on a solid surface : short time self-similarity Figure 2.
Top: Time sequence of the pressure field developing inside an impacting drop(Navier-Stokes
Gerris computation, fixed spatial magnification). Bottom: Corresponding traceof the pressure exerted by the drop on the solid substrate. The physical parameters for thissimulation are Re = 5000 and We = 250. The snapshots correspond respectively to times tU/R = 10 − , − and 10 − . nonetheless unable to reveal the inner structure of the contact zone until recently, es-sentially because of the very large scale ratio between this zone and the drop size. Theincrease of computational performance along with the development of adaptive numericaltechniques for two-phase flows (Popinet 2009) now allow to unravel the intimate structureof the contact zone, see Figures 1b,c and 2. These snapshots reveal a quite complex struc-ture for the pressure, which counter-intuitively exhibits sharp maxima near the contactline, and not on the axis as in steady stagnation point flows. Interestingly this structureis reminiscent of typical pressure field structures observed in the water entry of solidobjects, and evidenced by Wagner in the context of alighting seaplanes (Wagner 1932).In such problems a solid object impacts a flat liquid surface at a given velocity. Dropimpact may be viewed as water entry’s opposite, for a liquid object impacts a rigid flatsurface at a given velocity (see Fig. 3). It is therefore likely that the analytical tech-niques developed since the thirties to describe with great precision the flow generatedwith the impact of an object, and proven to be in close agreement with experimental data(Howison et al. et al. (2005)proposed a theoretical investigation of two-dimensional drop impact on a thin fluid layerand described the different regions and scalings of importance for the flow dynamics.In particular, they reveal the radius of contact between the two liquids as a key lengthscaling the problem, analogously to the problem of water entry where the wet lengthof the solid is also determining and consistently with the observations of Josserand &Zaleski (2003).The central motivation of the present study is to revisit the problem of a single spher-ical drop impacting a smooth flat solid surface at early times at the light of Wagner’stheory of impact, understand the dynamic fields structure and elucidate the short-timeself-similar behaviour discerned in earlier studies. To develop a consistent theory, theapproach followed throughout the manuscript will be to confront and cross-test sys-tematically the theoretical predictions with detailed and accurate numerical simulationsperformed with a Navier-Stokes multiphase flow solver. As a side note, we essay to makethe paper self-contained whenever possible. In §2 we formulate the hypotheses and theo-retical framework of the problem, and describe the short-time drop impact dynamics inthe context of Wagner’s theory. We put forward in particular a so-called “Lamb analogy”mirroring the flow within the impacting drop with the one induced with a flat rising ex-panding disk. In §3 we demonstrate that the Wagner flow can be recast as a self-similarsolution for the drop impact problem. The nature of the near-axis stagnation flow and ofthe near contact line flow and pressure maxima are also discussed. Numerical results ob-tained with Gerris (Navier-Stokes solver, VOF, adaptive mesh) taking into account surfacetension, surrounding gas and viscous effects are compared with the theoretical prediction.The structure of the boundary layer is examined in §4 and is found to be reminiscentof the viscous boundary layer leaved in the trail of a shockwave (Mirels analogy). Theinviscid Wagner flow and this viscous boundary layer are found to depend on the same
J. Philippi, P.-Y. Lagrée and A. Antkowiak
Figure 3.
Sketch of the impacting drop before contact (left), and shortly after impact (middle).The shape the drop would assume in absence of wall is outlined with a dashed line, and thecontact line position is marked with red dots. This problem may be viewed as the dual of theclassic water entry of a solid object (right). self-similar variable. This allows us to build a uniform, seamless solution encompassingboth impact and viscous effects. We conclude in §5 by discussing the obtained resultsand the limits of the present investigation (such as the role of air), and by remindingobservables of interest, such as the net impacting force of the total viscous shearing forceexerted by an impacting drop.
2. Model
Theoretical framework & hypotheses
We consider throughout this study an idealized drop impact where a perfectly sphericalliquid drop collides with a flat rigid surface. Though classic, this model situation relies ona number of physical hypotheses detailed in the following. Starting with the initial perfectspherical shape assumption, we may identify several typical causes for deviations fromsphericity such as capillary drop oscillations during free fall (Engel 1955; Thoroddsen et al.
U R/c during timescales of order
U R/c , where U stands forthe impact velocity, R for the drop radius and c for the celerity of compressive sound wavesin the liquid (( Uc ) (cid:28) UtR (cid:28) − representative of e.g. a raindrop falling at terminal velocity, wereckon that acoustic effects matter only in a micron-sized region over a few nanoseconds(see Weiss & Yarin 1999, for a discussion and references). The following discussion willtherefore be limited to those cases where the impact velocity is much lower than the speedof sound, as the falling raindrop, where acoustic effects can harmlessly be neglected andan incompressible description remains accurate. The high pressure and stresses generatedupon impact can result in marked erosion or yielding (Rein 1993). Furthermore substratedeformation has recently been shown to significantly alter drop impact in the limit ofvery soft (Mangili et al. et al. rop impact on a solid surface : short time self-similarity et al. i.e. Froude Fr = U /gR , Weber We = ρU R/σ and Reynolds Re = ρU R/µ numbers are all large with respect to unity. Here g denotes the gravity, σ the liquid-gassurface tension, ρ the liquid density and µ its viscosity. These assumptions underpinthe choice a purely inertial description free of these effects in the following. However,locally these phenomena might become more important or even dominant, examplesbeing viscosity near the boundaries or capillarity in high-curvature region. In section §4we will address viscous effects and develop a boundary-layer correction to the inviscidsolution, and capillary effects will eventually be discussed in §5.2.2. Governing equations and analogy with the water entry problem
Problem statement
We consider a perfectly spherical liquid drop of radius R and density ρ impacting nor-mally a flat rigid ground with velocity U , see Fig. 3. Neglecting for now the developmentof viscous rotational boundary layers, we assume that the fluid motion following impactis irrotational, axisymmetric and can be described with the scalar potential φ ( r, z ), i.e. the fluid velocity u ( r, z ) satisfies u ( r, z ) = ∇ φ ( r, z ). Incompressibility requires φ to be anharmonic potential satisfying Laplace’s equation, here written in cylindrical coordinates:1 r ∂∂r (cid:18) r ∂φ∂r (cid:19) + ∂ φ∂z = 0 . (2.1)The liquid dynamics obeys the unsteady form of Bernoulli’s conservation equation: ∂φ∂t + 12 |∇ φ | + pρ = const . (2.2)This set of equations is completed by appropriate boundary equations. At the wall z = 0,the condition of impermeability reads ∂φ∂z = 0 for 0 ≤ r ≤ d ( t ) , (2.3)where d ( t ) stands for the contact line position, an unknown of the problem. The positionof the free surface is tracked with the kinematic condition:d S d t = 0 , (2.4)where S ( r, z, t ) is a function vanishing on the free surface. Expressing normal stresscontinuity at this interface yields the following dynamic boundary condition: p = 0 at the free surface . (2.5)Note that atmospheric pressure as here been taken as the reference pressure.Anticipating the forthcoming analysis of the contact region, we further note that thefree fall behaviour outside the contact region can be recast into the following far-fieldcondition: φ → − U z far from the contact region . (2.6) J. Philippi, P.-Y. Lagrée and A. Antkowiak
Figure 4.
Flow in an impacting drop in the fixed frame and in the drop frame computedwith
Gerris . Left, top to bottom : Streamlines and pressure map within an impacting drop forRe = 5000 and We = 250 in the laboratory frame at different post-impact times (¯ t = 5 × − ,5 × − and 10 − ). The overall velocity field resembles a stagnation point flow in a near-wallregion whose extent is scaling with the wet area, and a uniform downwards flow outside. Right:Same velocity field in a reference frame moving with the drop initial velocity, evidencing a bypassmotion near the contact line and an overacceleration of the free surface towards the wall. This condition allows to identify the constant in (2.2) as U .Now nondimensionalising the problem using the inertial scales R , ρ and U , we introducethe following quantities: r = R ¯ r, z = R ¯ z, t = RU ¯ t, φ = U R ¯ φ, p = ρU ¯ p, (2.7)and rewrite the equations into their dimensionless counterparts:1¯ r ∂∂ ¯ r (cid:18) ¯ r ∂ ¯ φ∂ ¯ r (cid:19) + ∂ ¯ φ∂ ¯ z = 0 in the liquid, (2.8) ∂ ¯ φ∂ ¯ t + 12 | ¯ ∇ ¯ φ | + ¯ p = 12 in the liquid, (2.9) ∂ ¯ φ∂ ¯ z (¯ r, ¯ z = 0 , ¯ t ) = 0 over the wet area ¯ r < ¯ d (¯ t ) , (2.10)¯ p = 0 on the free surface, (2.11)d ¯ S d¯ t = 0 on the free surface. (2.12)Finally the nondimensional far-field condition reads:¯ φ = − ¯ z far from the contact region. (2.13)As posed, the problem entirely depends on the wet area extent ¯ d (¯ t ), whose dynamicshas still to be determined. In the following, we investigate the near-contact line flow toclarify this wetting dynamics.2.2.2. Contact line motion: numerical observations
To shed light on the contact line dynamics, detailed numerical simulations of impact-ing drops were carried out with
Gerris (see §6). Figure 4 represents typical streamlines rop impact on a solid surface : short time self-similarity Figure 5.
Advancing contact line shortly after impact. In the earliest moments following impact,the motion of the free surface near the contact zone is essentially directed downwards. Thesketches show the position of the contact line for two successive instants, and illustrate thefact that the horizontal extension of the wet area is governed by the vertical movement of theinterface. extracted from the simulations, shortly after impact. On the left panel it can be seenthat the motion within the impacting drop far from the contact zone is vertical, uniformand pointing downwards, corresponding merely to the free flight behaviour − U e z . Nearthe wall though, the flow is deflected and exhibits a stagnation point-like structure, ina region whose extent scales with the wet area. To investigate further the nature of thiscorrective flow, we represent on the right panel of Fig. 4 the streamlines in a referenceframe moving with the initial velocity of the drop. There it appears that the flow windsaround the contact line, revealing that (i) the liquid near the contact line falls fasterthan free-flight and (ii) rather than being pushed by a sweeping motion, the contact lineprogresses via a tank-treading movement, analogous to the rolling motion evidenced inprevious studies of advancing contact lines (Dussan V. & Davis 1974; Chen et al. d ( t ) can be obtained fromthe knowledge of velocity field at the free surface. Such a kinematic condition expressingthe contact between a liquid surface and a solid object has actually been used in thecontext of the water entry of solid objects for about 80 years, and is currently referred toas Wagner condition . In the following we depict the analogies between these two liquidimpact problems, and use them to derive a simple fluid mechanical model for the dropimpact.2.2.3.
Analogy with the water entry problem
The modern understanding of the liquid motion and forces generated by an impactingobject in water originates in the pioneering work of Wagner in the early thirties (Wagner1932). The primary motivation of Wagner was to provide a detailed characterization ofthe impulsive forces generated with impact – already known to be of sufficient amplitudeto induce bouncing (ricochet), and even possibly structural failure of alighting seaplanesor slammed ships (Nethercote et al. “gleiche
J. Philippi, P.-Y. Lagrée and A. Antkowiak
Figure 6.
In the reference frame of the falling drop, the flow induced by impact may be seen asthe one induced by a flat rising disk (Lamb disk analogy). The winding motion is here representedwith orange arrows, and the radial expansion of the disk with the wet area is indicated withpurple arrows. The motion of the disk itself is given by the red arrow.
Tragflügelbewegung” – equivalent aerofoil motion) is typically found to wind around theplate and therefore to promote jetting or splashing. The knowledge of this flow fieldthen allows to determine the motion of the free surface, and finally provides the neededcondition in the determination of the wet length d ( t ).Analogously, for the drop impact problem, our numerical simulations evidence similarflow features and winding motion. These observations advocate for the use of a waterentry-analogue description, where the flow induced by drop impact would correspond tothe one induced by a flat expanding disk in normal incidence, which extent is given thewet area (see figure 6). Following this vision of drop impact as a dual version of the waterentry problem, we adopt from now on the corresponding formalism to describe the fluidmechanics of impact.2.3. Leading-order description for the drop impact problem
Interested in the early-time behaviour of the impact-induced flow, we set out by exam-ining time-dependent solutions of system (2.8–2.12) in the vicinity of the contact zone.To this end, we introduce ε as a measure of the wet region: d ( t ) /R = O ( ε ) (see Fig. 7).This ε is the fundamental small parameter of our problem.As typical in two-phase phenomena, the lengthscales for the dynamical fields and forthe geometry of the free surface differ in this problem. Starting by considering the spacevariables ¯ r and ¯ z on which depend the dynamical fields (such as the velocity potential φ or the pressure p ), we introduce the following rescaling: ¯ r = ε r ˜ r and ¯ z = ε z ˜ z , where˜ r and ˜ z are O (1) quantities and ε r and ε z are gauge functions. From the structure ofLaplace operator, we expect the dynamical fields to display identical length scales in eachdirection, so that ε r = ε z = ε .Insights into the relevant lengthscales for the description of the free surface geometrycan be gained by decomposing the position of the surface into that of a translatingsphere ¯ z S (¯ r, ¯ t ) plus a surface disturbance ¯ h (¯ r, ¯ t ) (see Fig. 3b). Assuming the drop fallswith constant velocity, the shape of the unperturbed translating sphere obeys ¯ r + (¯ z S − (1 − ¯ t )) = 1. Sufficiently close to the contact area, we introduce gauge functions forthe vertical position of the moving sphere ¯ z S and the time ¯ t : ¯ z S = ε z S ˜ z S and ¯ t = ε t ˜ t .The equation for the sphere surface can be approximated by ε z S ˜ z S = ε ˜ r − ε t ˜ t . Aspreviously the determination of these scaling functions is obtained by dominant balance rop impact on a solid surface : short time self-similarity Figure 7.
Scalings in the contact zone. At the earliest times only a very small portion (of order ε ) of the drop touches the wall. The fluid sets into motion with impact is in a region of extent ε in every direction. The air wedge confined between the wall and the drop presents an angle oforder ε as well. The colormap illustrates the pressure distribution. The physical parameters forthis simulation are Re = 5000 and We = 250. This snapshot corresponds to a nondimensionaltime ¯ t = 10 − . The position of the contact line is here ¯ d = 1 .
73 10 − / . arguments: ε z S = ε t = ε . Note that at short times the intersection radius between thesphere and the impacting plane is given by ˜ r intersect = √ t .We remark that as in the original study of Wagner, a scale separation between ¯ z S and ¯ z exists ( small deadrise angle hypothesis , see e.g. Oliver 2002). This scale separation arisesbecause the drop typical radius of curvature ( O (1)) is very large in front of the otherlengthscales of the problem, see Fig. 7.We now turn on to the surface perturbation ¯ h (¯ r, ¯ t ), that embodies the impact-inducedflow. Recalling that ¯ h represents a perturbation around a falling sphere, we can expressthe position of the free surface with the following implicit equation: ¯ S (¯ r, ¯ z, ¯ t ) = ¯ z − (cid:0) ¯ z S − ¯ h (¯ r, ¯ t ) (cid:1) . Introducing an appropriate gauge function ε h such that ¯ h = ε h ˜ h we obtainby dominant balance analysis that ε h = ε . It follows that:˜ S (˜ r, ˜ z, ˜ t ) = ˜ z −
12 ˜ r + ˜ t + ˜ h (˜ r, ˜ t )= 0 on the free surface . (2.14)The kinematic boundary condition derives from the previous equation. At the free surface,we have: d ¯ S d¯ t = 1 + ∂ ˜ h∂ ˜ t |{z} O (1) − ε φ ˜ r ∂ ˜ φ∂ ˜ r | {z } O ( ε φ ) + ε φ ∂ ˜ h∂ ˜ r ∂ ˜ φ∂ ˜ r | {z } O ( ε φ ) + ε φ ε ∂ ˜ φ∂ ˜ z | {z } O ( ε φ /ε ) = 0 , (2.15)where ¯ φ = ε φ ˜ φ . It is impossible here to keep all terms at the same order; the dominantbalance between the vertical velocities ∂ ˜ h/∂ ˜ t and ∂ ˜ φ/∂ ˜ z implies that ε φ = ε . At leadingorder, the kinematic boundary condition is therefore reduced to :1 + ∂ ˜ φ∂ ˜ z + ∂ ˜ h∂ ˜ t = 0 on the free surface. (2.16)It proves convenient to introduce a translation of the velocity potential such that ˜ φ = − ˜ z + ˇ φ . This translation merely accounts to analyse the problem in the falling-dropreference frame. The kinematic boundary condition is then simply rewritten as: ∂ ˜ h∂ ˜ t = − ∂ ˇ φ∂ ˜ z on the free surface. (2.17)0 J. Philippi, P.-Y. Lagrée and A. Antkowiak ¯ r = ε ˜ r ¯ z = ε ˜ z ¯ t = ε ˜ t ¯ p = ε − ˜ p ¯ φ = ε ˜ φ (¯ u, ¯ v ) = (˜ u, ˜ v ) . Table 1.
Résumé of the most important asymptotic scales of the problem.
Inserting these different scaled variables into Bernoulli’s equation, we obtain: ε p ˜ p + 1 ε ∂ ˇ φ∂ ˜ t + 12 ∂ ˇ φ∂ ˜ r ! + − ∂ ˇ φ∂ ˜ z ! = 12 in the liquid , (2.18)where ¯ p = ε p ˜ p . The scale of the pressure, ε p = ε , is here seen to be as large as the contactzone is small – as expected in an impact problem. At leading order, Bernoulli’s equationis therefore reduced to: ˜ p = − ∂ ˇ φ∂ ˜ t in the liquid . (2.19)It follows from this equation that the constant pressure Dirichlet boundary conditionon the free surface p = 0 can be recast as a condition for the potential at the free surface: φ = const, where the constant is arbitrary. Without loss of generality, we set from nowon this constant to zero.Finally, as classic in water wave theory, we exploit the shallowness of the gap betweenthe free surface and the plane to transfer the boundary condition at the free surface ontothe plane (see e.g. Van Dyke 1975, §3.8).Summarizing, the model problem takes the following expression:1˜ r ∂∂ ˜ r ˜ r ∂ ˇ φ∂ ˜ r ! + ∂ ˇ φ∂ ˜ z = 0 in the liquid, (2.20) − ∂ ˇ φ∂ ˜ t = ˜ p in the liquid, (2.21)the locus ˜ d (˜ t ) of the contact line is determined with the Wagner condition:˜ h (˜ r, ˜ t ) = 12 ˜ r − ˜ t for ˜ r = ˜ d (˜ t ) , (2.22)so that the boundary conditions at ˜ z = 0 read:ˇ φ = 0 for ˜ r > ˜ d ( t ) , (2.23) ∂ ˜ h∂ ˜ t = − ∂ ˇ φ∂ ˜ z for ˜ r > ˜ d ( t ) , (2.24) ∂ ˇ φ∂ ˜ z = 1 for ˜ r < ˜ d (˜ t ) , (2.25)and the far-field behaviour is given by: ˇ φ → r, ˜ z → ∞ (2.26)˜ h → r → ∞ . (2.27)Finally the corresponding model geometry is sketched Fig. 8. We remark that the previousset of equations resembles to that of the classic water entry problem, and can be solved rop impact on a solid surface : short time self-similarity Figure 8.
Leading order outer problem for times of order ε . using the methodology described in e.g. Oliver (2002). In the next section though we willpresent an alternate method based on self-similar solutions.
3. Self-similar solutions and numerical simulations
A self-similar problem
To reveal the self-similar nature of our problem, we classically look in the following forscale invariance (Darrozès & François (1982)). We start by expressing the fact that anyvariable ˜ q in (˜ r, ˜ z, ˜ t, ˇ φ, ˜ h, ˜ d, ˜ p ) can be rewritten as ˜ q = λ q ˆ q , where ˆ q is a rescaled variableand λ q a numerical stretching coefficient embodying the change of scale. Inserting thesevariables into the governing equations, it is straightforward to see that invariance ofLaplace equation through this stretching requires λ r = λ z . Similarly, expressing theinvariance of Wagner condition yields λ h = λ t , λ r = √ λ t and λ d = √ λ t . The sameoperation performed on the additional boundary conditions finally imposes λ φ = √ λ t and λ p = 1 / √ λ t . Note that λ t remains here as the sole stretching parameter.The pressure field can be written as an implicit function of time and space as fol-lows: F (˜ p, ˜ r, ˜ z, ˜ t ) = 0. Upon using the previous scale invariance arguments, this rela-tion may be rewritten as F (ˆ p/ √ λ t , √ λ t ˆ r, √ λ t ˆ z, λ t ˆ t ) = 0. A simple algebraic manip-ulation allows to remove the λ t dependence for all but one variables, so that finally G ( √ ˆ t ˆ p, ˆ r/ √ ˆ t, ˆ z/ √ ˆ t, λ t ˆ t ) = 0, for any λ t . Remarking that for a given ˆ t , this function hasto cancel whatever the choice of the scale λ t , it readily appears that the last variable issuperfluous. In other words, a relation linking √ ˆ t ˆ p to ˆ r/ √ ˆ t and ˆ z/ √ ˆ t only must exist.The pressure field may therefore be rewritten explicitly as:˜ p = 1 √ ˜ t P (cid:18) ˜ r √ ˜ t , ˜ z √ ˜ t (cid:19) . (3.1)With a similar reasoning, and upon introducing the self-similar variables ξ = ˜ r/ √ ˜ t and η = ˜ z/ √ ˜ t , we readily obtain :ˇ φ (˜ r, ˜ z, ˜ t ) = p ˜ t Φ( ξ, η ) , ˜ h (˜ r, ˜ t ) = ˜ t H ( ξ ) and ˜ d (˜ x, ˜ t ) = p ˜ t δ, (3.2)where Φ and H are unknown functions of the self-similar variables and δ a constantrepresenting the (fixed) position of the contact line in self-similar space. This allows us2 J. Philippi, P.-Y. Lagrée and A. Antkowiak to formulate the self-similar version of the drop impact problem :1 ξ ∂∂ξ (cid:18) ξ ∂ Φ ∂ξ (cid:19) + ∂ Φ ∂η = 0 in the liquid, (3.3) P ( ξ, η ) = 12 (cid:18) − Φ( ξ, η ) + ξ ∂ Φ ∂ξ + η ∂ Φ ∂η (cid:19) in the liquid, (3.4)the boundary conditions at η = 0 take the following form: H − ξ ∂ H ∂ξ = − ∂ Φ ∂η for ξ > δ, (3.5) ∂ Φ ∂η = 1 for ξ < δ, (3.6)Φ = 0 for ξ > δ, (3.7)the far-field behaviour is: Φ → ξ, η → ∞ (3.8) H → ξ → ∞ , (3.9)and the self-similar version of Wagner condition is finally given by: H ( ξ ) = 12 ξ − ξ = δ. (3.10)This problem can now be solved in several steps.3.2. Self-similar potential
In this geometry, Laplace equation can be solved with variable separation, leading to afamily of elementary cylindrical harmonic solutions with an exponential behaviour in η and an oscillatory one in ξ . We recompose by summation and obtain :Φ( ξ, η ) = Z ∞ C ( k ) J ( kξ ) e − kη d k. (3.11)The weight function C ( k ) is determined with boundary conditions (3.6) and (3.7), leadingto the following pair of dual integral equations: Z ∞ k C ( k ) J ( kξ ) d k = − ξ < δ , (3.12 a ) Z ∞ C ( k ) J ( kξ ) d k = 0 for ξ > δ . (3.12 b )Solving these dual integral equations using the technique described in Sneddon (1960),we obtain a closed-form expression for the weight function: C ( k ) = 2 π δk cos( kδ ) − sin( kδ ) k = 2 π dd k (cid:18) sin( kδ ) k (cid:19) . (3.13)Anticipating the description of the contact line dynamics, we now derive ∂ Φ /∂η at thesubstrate level η = 0: ∂ Φ ∂η = − π Z ∞ kδ cos( kδ ) − sin( kδ ) k J ( kξ ) k d k, (3.14) rop impact on a solid surface : short time self-similarity e.g. Sneddon 1995, table IV,page 528). This allows us to obtain the following explicit expression for ∂ Φ /∂η for η = 0: ∂ Φ ∂η = 1 for ξ < δ and ∂ Φ ∂η = − π δ p ξ − δ − arcsin (cid:18) δξ (cid:19)! for ξ > δ . (3.15)3.3. Wagner condition and contact line dynamics
With the help of the vertical velocity expression in the near wall region just derived, wecan rewrite the kinematic boundary condition (3.5) as: H ( ξ ) − ξ ∂ H ∂ξ ( ξ ) = 2 π δ p ξ − δ − arcsin (cid:18) δξ (cid:19)! for ξ > δ . (3.16)This inhomogeneous differential equation can be solved using variation of parameters, i.e. looking for a solution of the form H ( ξ ) = ξ f ( ξ ). This gives: (cid:20) f ( ξ ) (cid:21) + ∞ δ = − π Z ∞ δ ξ δ p ξ − δ − arcsin (cid:18) δξ (cid:19)! d ξ. (3.17)Upon using the far-field decaying behaviour of H (see equation (3.9)), this last equationreduces to f ( δ ) = δ − so that at the contact line the drop deformation is: H ( δ ) = δ f ( δ ) = 12 . (3.18)In the self-similar space, the Wagner condition therefore takes the following remarkablysimple form: 12 = 12 δ − , (3.19)from which we finally derive the position of the contact line: δ = √ . (3.20)It is interesting to remark that the contact line motion ˜ d (˜ t ) = √ t just predictedwithin the framework of Wagner theory is quite close from the rough truncated sphereapproximation ˜ r intersect = √ t (Rioboo et al. Gerris along with ourtheoretical prediction. Noticeably the superposition between theory and numerical resultsis excellent, at least until the moment of formation of a liquid corolla (here indicated witha red dot).3.4.
Analogy with the normal motion of an expanding disk in an infinite mass of liquid
In §2.2.3 we proposed to visualize the flow in an impacting drop as the one induced by aflat rising disk expanding radially as the wet area (see also Fig. 6). We are now in a posi-tion to formally justify this water-entry analogy. The axisymmetric flow induced by ‘themotion of a thin circular disk with velocity U normal to its plane, in a infinite mass of liq-uid’ is for example analysed in Lamb’s classic textbook §101 (Lamb 1932). After derivingsome elementary axisymmetric solutions of Laplace equation of the form exp( ± kz ) J ( kr )in §100, Lamb examined a variety of axisymmetric potential flows. Among those was theone (later connected to the flow around a flat circular disk in normal incidence) whereat the symmetry plane z = 0 the potential takes the value φ = C √ a − r for r < a and4 J. Philippi, P.-Y. Lagrée and A. Antkowiak
Figure 9.
Comparison between the theoretical position of the contact line as a function oftime deduced from Wagner theory ˜ d (˜ t ) = √ t (dashed line) and the position of the contact lineextracted from Gerris computations of an impacting drop at Re = 5000 and We = 250 (greydots). The red dot marks the birth of the corolla. φ = 0 for r > a , with a the disk radius. The solution for this problem was stated underthe following integral representation: φ ( r, z ) = − C Z ∞ e − kz J ( kr ) dd k (cid:18) sin kak (cid:19) d k. (3.21)And from ‘a known theorem in Electrostatics’ , Lamb obtained the expression for thevertical velocity in the symmetry plane: − (cid:18) ∂φ∂z (cid:19) z =0 = πC for r < a , (3.22 a ) C (cid:18) arcsin (cid:16) ar (cid:17) − a √ r − a (cid:19) for r > a . (3.22 b )This corresponds precisely to the flow within the impacting drop, after posing C = − /π and a = δ , thereby justifying formally our initial analogy between the impact-inducedflow with the one associated with a flat rising disk rapidly expanding with the wet area.Setting C = 2 U/π , Lamb remarked that the above potential indeed describes the flowwinding around a flat disk moving at velocity U . He further noted that a simple expressionfor the fluid half-space kinetic energy could be derived from the previous relation: T disk = 43 ρa U . (3.23)This expression can immediately be transposed into the (nondimensional) kinetic energyof the impact-induced flow within the drop:˜ T = 4 √ t / , (3.24)or, equivalently, into its dimensioned counterpart: T = 4 √ ρU / R / t / . (3.25)We emphasize that this expression is derived within the frame of the falling drop and, assuch, represents the kinetic energy of the defect flow associated with impact. Although adirect physical interpretation of this quantity is not straightforward, we will see in §3.7 rop impact on a solid surface : short time self-similarity Structure of the velocity field
We now investigate the structure of the velocity field in the contact region and searchfor exact closed-form expressions and convenient approximations for this field.3.5.1.
Integral representation of the velocity field
In the fixed frame, the velocity field ˜ u (˜ r, ˜ z, ˜ t ) inside the impacting drop can formally bederived from the (untranslated) potential − η + Φ. Following the arguments developed in§3.1, this velocity field is simply related to the self-similar velocity field U ( ξ, η ) via therelation: ˜ u (˜ r, ˜ z, ˜ t ) = U ( ξ, η ) . (3.26)where the components of the self-similar vector field U = ( U ξ , U η ) are: U ξ ( ξ, η ) = ∂ Φ ∂ξ , (3.27 a ) U η ( ξ, η ) = − ∂ Φ ∂η . (3.27 b )Inserting the expression of the self-similar potential determined previously yields thefollowing integral representation for the vector field components: U ξ ( ξ, η ) = − π Z ∞ √ k cos( √ k ) − sin( √ k ) k e − kη J ( kξ ) d k, (3.28 a ) U η ( ξ, η ) = − − π Z ∞ √ k cos( √ k ) − sin( √ k ) k e − kη J ( kξ ) d k. (3.28 b )A closed-form expression is unfortunately not accessible in the general case. In the fol-lowing however we calculate the value of these integrals at some particular places.3.5.2. Closed-formed expressions for the velocity field along the axis and the substrate
Simple analytical solutions for the velocity field can be obtained from (3.28) at preciselocations. Along the symmetry axis for example, where ξ = 0, the properties of integralsof exponentials allow to write: U η ( ξ = 0 , η ) = − π (cid:18) arctan (cid:18) √ η (cid:19) − √ η η (cid:19) for η ≥ . (3.29)This last result is confronted with numerical velocity profiles extracted from Gerris compu-tations in Fig. 10. The nice agreement between the theoretical solution and the numericalprofiles seen in the self-similar space (Fig. 10b) here holds over more than a decade intime.Analogously, analytical forms for (3.28) can also be obtained along the substrate plane η = 0 by exploiting the properties of Hankel transforms (Sneddon 1995). An expressionfor the vertical velocity is already provided with equation (3.15), after inserting δ = √ U ξ ( ξ, η = 0) = 2 π ξ p − ξ for 0 ≤ ξ < √ . (3.30)6 J. Philippi, P.-Y. Lagrée and A. Antkowiak
Figure 10.
Left : Axial velocity profiles along the axis extracted from
Gerris computations attimes ¯ t = 3 × − , × − , − , × − and 10 − . Right : Comparison between the analyt-ical prediction for the axial velocity given by equation (3.29) (red dashed line) and numericalsolutions obtained with Gerris , rescaled in the self-similar space (blue solid lines). The physicalparameters for this simulation are Re = 5000 and We = 250.
This unphysical inviscid slip velocity ˜ u e ( ξ ) cannot be observed in our simulations encom-passing viscous effects. But this quantity is nonetheless relevant for it corresponds to theedge velocity of the viscous boundary layer (studied in detail in §4).3.5.3. An unusual stagnation point flow
In the very vicinity of the origin, the first order power series of the velocity field (3.28)reads: U ξ ( ξ, η ) ’ π √ ξ, (3.31 a ) U η ( ξ, η ) ’ − π √ η, (3.31 b )or, equivalently, in dimensioned variables: u r ( r, z, t ) ’ π √ r UR r √ t , (3.32 a ) u z ( r, z, t ) ’ − π √ r UR z √ t . (3.32 b )Though simple, this peculiar structure for the impact-induced unsteady stagnation pointflow is nonetheless counter-intuitive and could not have been inferred from simple dimen-sional analysis. Noteworthy enough, this result is at variance with the typical structureof the later intermediate flow associated with spreading ¯ u ’ (¯ r/ ¯ t, − z/ ¯ t ) (see e.g. Eggers et al. et al.
Beyond the stagnation point: a remark on the overall velocity field structure
The previous approximation for the impact flow is valid in a small region near the origin.To further investigate the limits of this representation we show Fig. 11 different radialvelocity profiles corresponding to various locations ξ . The collapse of the numerical pro-files taken at various ¯ r and ¯ t (but such that ¯ r/ √ ¯ t = ˜ r/ √ ˜ t is constant in each figure)onto the theoretical profiles is again an illustration of the relevance of the self-similarrepresentation. But it is also to be noted that while the stagnation point ansatz disre-gards any radial velocity variation in η , the profiles exhibit a sensible variation along thevertical coordinate η . This variation is best depicted with Fig. 12 where theoretical radial rop impact on a solid surface : short time self-similarity Figure 11.
Self-similar radial velocity as a function of η for ξ = 0.125, 0.25, 0.5, 1 and 1.5. Thesevelocities are rescaled by the outer solution of the boundary layer ˜ u e ( ξ ) given by equation (3.30).Blue solid lines represent the numerical solutions extracted from Gerris computations in theself-similar space for ¯ t = 5 × − , − , × − and 10 − (Re = 5000 and We = 250). The reddashed line represents the theoretical solution U ξ ( ξ, η ) given by equation (3.28 a ). Note that theboundary layer is so thin that it is almost indistinguishable (see also Fig. 19). velocity profiles taken at different values ξ have been represented. Noteworthy enough,profiles corresponding to ξ (cid:46) U ξ ( ξ, η ) = u e ( ξ ) f ( η ) in this region. Forlarger values of ξ though, significant deviations from this behaviour arise and variableseparation cease to hold: U ξ ( ξ, η ) = ˜ u e ( ξ ) g ( ξ, η ) with g ( ξ, η = 0) = 1 for ξ (cid:38) J. Philippi, P.-Y. Lagrée and A. Antkowiak
Figure 12.
Evolution of the analytical self-similar radial velocity given by equation (3.28 a ) asa function of η for ξ = 0.125, 0.25, 0.5, 1 and 1.5. These velocities are rescaled by the outersolution of the boundary layer U ξ ( ξ, η = 0) = ˜ u e ( ξ ), equation (3.30). Figure 13.
Comparison between the flow pattern within an impacting drop (left) and arounda rapidly expanding disk (Lamb analogy, right) in the self-similar space. In both cases, thestreamlines are represented in the moving frame. The red dots represent the theoretical positionof the contact line ξ = √
3. The numerical streamlines represented on the left are derived fromthe velocity field computed with
Gerris at ¯ t = 10 − (Re = 5000 and We = 250). The theoreticalstreamlines shown on the right correspond to isovalues of Ψ( ξ, η ) defined in equation (3.34) (notethe correspondence with Lamb’s figure page 145). Flow pattern, contact line bypass and Lamb analogy
We now define the self-similar stream function Ψ( ξ, η ) in the drop reference frame fromthe potential: ξ ∂ Ψ ∂η = ∂ Φ ∂ξ , (3.33 a ) − ξ ∂ Ψ ∂ξ = ∂ Φ ∂η . (3.33 b )By integration of the previous relations we deduce the following expression for Ψ( ξ, η ):Ψ( ξ, η ) = 2 π Z ∞ √ k cos( √ k ) − sin( √ k ) k e − kη ξJ ( kξ ) d k, (3.34)up to a constant. Formally, Ψ( ξ, η ) is the stream function describing the winding flowaround a flat rising disk (Lamb 1932, §108). Figure 13 offers a comparison betweenthe streamlines of this Lamb analogy and the ones computed with Gerris for the dropimpact problem in the self-similar space. A good qualitative agreement between theanalytical and the numerical streamlines is noticeable, comforting the expanding diskanalogy followed here. Interestingly the winding motion around the contact line, as wellas the falling velocity overshoot near this region, are both captured with this analogy rop impact on a solid surface : short time self-similarity Figure 14.
Comparison between the pressure field developing inside an impacting drop (left)and around a rapidly expanding disk (Lamb analogy, right). The pressure field represented onthe left is extracted from
Gerris computations and represented in the self-similar space (¯ t = 10 − ,Re = 5000, We = 250; The isovalues are: 0.12, 0.24, 0.36, 0.48, 0.6, 0.72, 0.84). The self-similartheoretical pressure field represented on the right is given by equation (3.4) (isovalues: 0.13,0.28, 0.445, 0.57, 0.73, 0.9, 1.2). Though isovalues have been slightly changed between the twopanels, theoretical and numerical results are in a good overall agreement. and can be correlated with the peculiarities of the winding flow near the edge of a risingdisk. 3.6. Self-similar pressure
From the knowledge of the velocity potential we are now in a position to derive thepressure field as the time derivative of the potential. In the self-similar space, the pressurefield is given by equation (3.4). Figure 14 proposes a comparison between the structureof the self-similar pressure extracted from numerical computations performed with
Gerris and the theoretical prediction. There it can be seen that the overall structure of thepressure field developing in the impacting drop, and in particular the pressure peak inthe vicinity of the contact line already pinpointed out in Fig. 1, nicely matches withthe theory. Interestingly, the structure just described is at variance with the pressuredistribution around a flat disk rising steadily (Lamb’s original problem). Indeed in sucha configuration the pressure is expected to be maximal in the stagnation point area,whereas in our model problem the pressure peaks near the contact line/disk edge. This isa consequence of the motion unsteadiness: the pressure is here dominated by the ∂ ˇ φ/∂ ˜ t contribution rather than the steady ˜ ∇ ˇ φ term.As in §3.5.2, closed form expressions for the pressure can be obtained along the axisand the substrate plane. The radial structure of the self-similar pressure across the wetarea reads P ( ξ, η = 0) = π √ − ξ for 0 ≤ ξ < √
3. This analytical prediction is con-fronted with
Gerris numerical results in Fig. 15. After a transient numerical initializationphase (corresponding to the red curves), the pressure profiles collapse on the self-similaranalytical solutions (blue curves). In accordance with the overall pressure field structuredepicted earlier, the pressure radial profile presents a local minimum at ξ = 0 and amaximum in the vicinity of the contact line, that is for ξ = √ J. Philippi, P.-Y. Lagrée and A. Antkowiak
Figure 15.
Left: Pressure trace on the substrate ¯ z = 0 obtained from the numerical simulationsbetween ¯ t = 5 × − and 10 − for Re = 5000 and We = 250. The color code for each decade isthe same as in Fig. 17. Note that the curves are equally distributed within each decade. Right:Same as in the left but in the self-similar space. The black dashed curve represent the analyticalsolution for P ( ξ, η = 0). Figure 16.
Pressure along the axis ¯ r = 0 obtained with Gerris for ¯ t = 5 × − , − , × − and10 − (same as in Fig. 10 right) represented in the self-similar space (Re = 5000 and We = 250).The red dashed line is the analytical solution for P ( ξ = 0 , η ). Fluctuations of the pressure aroundthe theoretical prediction is to be related with numerical projections errors. Similarly the expression for the self-similar pressure along the symmetry axis can alsobe obtained analytically: P ( ξ = 0 , η ) = √ π (3+ η ) for η ≥
0. Figure 16 compares this lastresult with rescaled axial pressure profiles extracted from numerical simulations. Thereagain the agreement between the computations and the theory is seen to hold for a largetime span.The structure of the pressure field in the vicinity of the origin can be inferred fromthese last results. Reexpressing the pressure cuts determined in terms of ˜ r , ˜ z and ˜ t , weget: ˜ p (˜ r, ˜ z = 0 , ˜ t ) = 3 π √ t − ˜ r , (3.35 a )˜ p (˜ r = 0 , ˜ z, ˜ t ) = 3 √ tπ (3˜ t + ˜ z ) . (3.35 b )From these two relations, we deduce the following expansion for the pressure near theorigin: ˜ p (˜ r, ˜ z, ˜ t ) = √ π ˜ t − (cid:18) r t − ˜ z t (cid:19) + . . . (3.36)This expression provides with a local approximation for ∂ ˇ φ/∂ ˜ t from which, after time in-tegration and space differentiation, we readily recover the stagnation point flow structure rop impact on a solid surface : short time self-similarity Figure 17.
Time evolution of the pressure ¯ p (0 , , ¯ t ) measured at the origin in the numericalsimulations with Gerris (Re = 5000 and We = 250). Note that each decade is represented witha different colour. The theoretical prediction ¯ p (0 , , ¯ t ) = √ π ¯ t − is superimposed with a blackdashed line. found earlier: (˜ u ˜ r , ˜ u ˜ z ) = π √ (˜ r/ √ ˜ t, − z/ √ ˜ t ) . This near-axis behaviour emphasizes againthat simple intuitive dimensional analysis suggestion r/t and − z/t is here not relevant.The leading order term for the pressure at the origin follows:˜ p (0 , , ˜ t ) = √ π ˜ t − , or, with dimensions p (0 , , t ) = ρU / π r Rt . (3.37)This result extends the t − scaling law proposed by Josserand & Zaleski (2003) on thebasis of scaling arguments. A comparison between this theoretical prediction and Gerris numerical simulations is proposed Fig. 17, using the color code of Fig. 15. After a numeri-cal transient phase, the pressure rapidly reaches the self-similar regime. Remarkably, thisshort-time similarity regime holds for almost 4 decades in time and is nicely capturedwith Wagner theory. Eventually a departure from similarity occurs when ¯ t becomes oforder 1, and the pressure promptly drops to 0.3.7. Normal force induced by drop impact
Building on the last set of results, we deduce the total net normal force imparted by animpacting drop on the underlying substrate at early times. Integrating the pressure onthe wet surface, we have:˜ F (˜ t ) = 1 √ ˜ t Z Z S P ( ξ, η = 0) d S = 2 π p ˜ t Z √ ξπ p − ξ d ξ = 6 p t. (3.38)The dimensional counterpart of this net total force induced by the drop on the substratetherefore reads: F ( t ) = 6 √ ρU / R / √ t, (3.39)where the force is seen to increase as t for short times. Interestingly F ( t ) could have beeninferred directly from energy arguments, with no knowledge of the pressure distribution.Indeed, writing the global kinetic energy conservation for the upper semi-infinite space,we have: dd t T = − I p u · n d S, where T = (cid:18)Z Z Z ρu V (cid:19) . (3.40)In the context of a flat rising disk, the kinetic energy reduces to T disk = ρa U (Lamb§102). This expression can immediately be transposed to the impacting drop problem sothat T = 4 √ ρU / R / t / (see equation (3.25) in §3.4). The power of pressure forcesthen follows as dd t T = 6 √ ρU R (cid:0) UtR (cid:1) / . Dividing this power by U , we recover exactly2 J. Philippi, P.-Y. Lagrée and A. Antkowiak the previously obtained result for the net normal total force. This alternate derivationof the normal force provides with yet an other illustration of the relevance of Lamb’sanalogy for the drop impact problem.
4. Matching with the viscous solution
The inertial limit (large Reynolds number hypothesis) investigated so far has allowedus to model the flow within an impacting drop as the winding motion of an inviscid fluidaround an expanding disk, appropriately described by an harmonic potential obeyingthe unsteady Bernoulli equation (§2.2). Actually the agreement between the correspond-ing theoretical results and numerical Navier-Stokes computations carried out with
Gerris (encompassing viscous effects) comforted this approximation, see e.g.
Figs 10 for veloc-ity, 16 for pressure or 9 for contact line motion comparisons. Most presumably, viscouseffects are here dominating only in very thin boundary layers developing along the wetsubstrate. And indeed, even if the overall agreement between the radial velocity profilesand the inviscid solution is evident, a careful examination of Figure 11 reveals the pres-ence of these thin layers in the very vicinity of the solid wall. Even if spatially confined,these boundary layers nonetheless play a key role when comes e.g. the question of theerosion potential of an impacting drop. Consequently we now set out to describe the in-ner structure of these viscous layers and to match it to the previously determined outerinviscid solution. Viscous shear stresses and total erosion potential are eventually brieflydiscussed. 4.1.
A simple boundary layer problem?
Typically, the (inviscid) slip velocity ˜ u e (˜ r, ˜ t ) = π ˜ r/ √ t − ˜ r , here first introduced equa-tion (3.30), and the no-slip condition at the substrate, trademark of real fluids, arereconciled through the introduction of a viscous boundary layer. According to the classicboundary layer theory ( e.g. Schlichting 1968), the transverse scale of this layer is Re − / ,so that an appropriate inner coordinate ˜ Z can be defined via ˜ z = Re − / ˜ Z . The mostsimple idea at this point is to think that the outer variables scales defined in §2.3 implythat the non-linear terms of the boundary layer equation are negligible when comparedto unsteady and viscous terms, so that this equation would simply read: ∂ ˜ U r ∂ ˜ t = − ∂ ˜ p∂ ˜ r + ∂ ˜ U r ∂ ˜ Z , (4.1)where capitalized variables refer to boundary layer quantities. Considering that Eulerequation in the inviscid outer domain reduces to ∂ ˜ u e ∂ ˜ t = − ∂ ˜ p∂ ˜ r , the boundary layer equationcan be recast as the following diffusion equation for the defect velocity: ∂∂ ˜ t ( ˜ U r − ˜ u e ) = ∂ ∂ ˜ Z ( ˜ U r − ˜ u e ) . (4.2)The corresponding solution can then be expressed as a convolution between the forcingterm and the Green function of the heat equation:˜ U r = ˜ u e (˜ r, ˜ t ) − ˜ Z √ π Z ˜ t exp (cid:18) − ˜ Z t − τ ) (cid:19) ˜ u e (˜ r, τ )(˜ t − τ ) d τ. (4.3)Unfortunately this solution leads to a paradoxical cancelling of shear stresses at the wall.We conjecture that this unreasonable result stems from the fact that the sharp longi-tudinal variations associated with the contact line motion have here been disregarded.Specifically non linear terms do balance unsteady terms, at least near the contact line rop impact on a solid surface : short time self-similarity r = √ t . As a result, the boundary layer actually grows from this moving pointboth in space and time. While a comprehensive analysis of this problem demands a care-ful balance of each term likely resulting in a non linear boundary layer problem, beyondthe scope of the present study, we nonetheless propose in the following an approximationbased on an analogy with boundary layers developing behind shockwaves.4.2. Approximation of the drop impact boundary layer via an analogy withshock-induced boundary layers
We now depict qualitatively the inner viscous structure of the velocity field by using asimple analogy. First remembering the tank-treading movement in the vicinity of thecontact line observed and discussed in §2.2.2, we point out the violent change in radialvelocity when passing through the contact line. In other words, the contact line embodiesa neat discontinuity where the slip velocity sees its value suddenly change from 0 to ˜ u e .Building on this observation, we consider in the following the contact line as a kind ofshock wave sweeping the substrate, and seeding a boundary layer in its trail (see figure18). This problem is classic in compressible flows and was solved by Mirels (1955) inthe context of a shock tube (see Schlichting 1968, for more details). In this study, afluid initially at rest is swept by a shockwave travelling at celerity U s in the direction x and instantly acquires an impulse of velocity U ∞ in the process. Behind the normalshockwave is left a growing viscous boundary layer.The Ansatz for Mirel’s solution is to introduce η m = z/ p t − x/U s as the self-similarvariable. This variable not only takes into account time variations but also longitudinaleffects from the shock backwards in x . Disregarding any pressure gradient but consideringboth unsteady and nonlinear effects, the momentum equation may be rewritten in termsof η m and of the velocity U ∞ f ( η m ): f ( η m ) + 12 ( η m − U ∞ U s f ( η m )) f ( η m ) = 0 , with f (0) = f (0) = 0 , and f ( ∞ ) = 1 . (4.4)Note that compressible effects have here been absorbed via an appropriate Lees-Dorod-nitsyn’s transformation (see Stewartson 1964). Two limiting cases clearly emerge fromthe picture. For large U ∞ /U s (and after a rescaling and a change of sign due to the choiceof origin), the velocity profile tends to a Blasius profile. Conversely, for small values ofthe velocity ratio, the velocity rather adopts an error function profile. Note that profilescorresponding to intermediate values of this ratio can be found in Schlichting’s textbook.From this sound result we may by analogy transpose this approach to the drop impactproblem (see Fig. 18). Obviously the outer solution for the drop impact problem is quitemore complex for neither U s nor U ∞ are constant. The core idea consists in drawing aparallel between the shock (at position U s t ) and the contact line (at position √ t ) onthe one hand, and between the steady slip velocity U ∞ and ˜ u e (˜ r, ˜ t ) on the other hand.Following this simple analogy the longitudinal velocity is approximated with:˜ U r (˜ r, ˜ z, ˜ t ) = 2˜ rπ √ t − ˜ r f ˜ z p ˜ t − ˜ r / √ Re ! . (4.5)where f is solution of an equation which is analogous to Eq. (4.4). The so-called com-posite solution (Van Dyke 1975), which is an expansion valid in the ideal fluid and in the4 J. Philippi, P.-Y. Lagrée and A. Antkowiak
Figure 18.
Left: Sketch of the contact line during its motion and of the growing boundarylayer in its trail, analogous to that developing behind a shockwave. Right: Shockwave-inducedboundary layer, reproduced from the german edition of Schlichting textbook (Schlichting 1968).Notations are from Schlichting, with a correspondence between x and ˜ r . Note that in the shock-wave case, U ∞ and U s are both constant. boundary layer, then follows:˜ u comp r = − π Z ∞ √ k cos( √ k ) − sin( √ k ) k e − k ˜ z √ ˜ t J ( k ˜ r √ ˜ t ) d k ++ 2˜ rπ √ t − ˜ r f ˜ z p ˜ t − ˜ r / √ Re ! − ! . (4.6)In practice we approximated f with erf function. Figure 19 proposes a comparison be-tween the numerical velocity profiles extracted from Gerris computations and this approx-imation, which proves to provide a fairly good description for the flow. As a side note,we remark that replacing the error function with Blasius profile yields a slightly moremarked deviation between theory and numerical results. That said we chose not to tunethe velocity ratio appearing in equation (4.4) as (i) this is too speculative and (ii) suchadjustment is certainly beyond the limits of our analogy.It is interesting to note that for ξ smaller than √
3, Mirel’s self-similar variable η m tends to Re / η = ˜ Z/ √ ˜ t , so that the vertical structure for solution (4.5) now simplyinvolves erf (cid:16) Re / η (cid:17) . Actually, this solution is merely the purely diffusive solution ofequation (4.2) for constant forcing (˜ u e constant in time).From the previous results we may extract several quantities, such as the displacementthickness or the locus of iso-velocities. The displacement thickness δ can readily beestimated with f = erf as: δ = 1 √ Re Z + ∞ (cid:18) − ˜ U r ˜ u e (cid:19) d ˜ Z = 2 √ π √ Re q ˜ t − ˜ r / . (4.7)Similarly, isolines for the velocity can be extracted both for Gerris computations and forboundary layer theory. Figure 20 provides with a qualitative comparison between theoryand numerical results, and it can be remarked that the overall prediction is more thanjust qualitative.Though the velocity ratio (2 r/π/ √ t − ˜ r ) / ( √ t/ / ˜ t ) (counterpart of U ∞ /U s in equa-tion (4.4)) is infinite near the shock, we note that the agreement between numerical andtheoretical solutions is actually surprisingly good. We conjecture that whenever this ratiodecreases to a value lower than one, i.e. near the centre and for large times where thisratio behaves as (4˜ r ) / (3 π ) so tends to 0, the error function approximation emerges as rop impact on a solid surface : short time self-similarity Figure 19.
Inner-boundary layer radial velocity profiles at different locations ξ : 0.125, 0.25,0.5, 1 and 1.5. Blue solid lines correspond to numerical solutions obtained with Gerris at¯ t = 5 × − , − , × − and 10 − for Re = 5000 and We = 250 and represented in theself-similar space. Note that velocities are rescaled by their maximum value. The red dashedlines stand for the theoretical composite solution (equation (4.6)) blending the self-similar vis-cous boundary layer solution with the self-similar Wagner inviscid solution for impact. Thecomposite solution is also rescaled by the edge velocity ˜ u e ( ξ ) given by equation (3.30). the solution of equation (4.4). Eventually we remark that an in-depth analysis of thesephenomena demands a more involved description for the boundary layer (such as theInteractive Boundary Layer theory, see e.g. Lagrée 2010) and/or a deeper analysis of theWagner region (Oliver 2002; Korobkin 2007; Oliver 2007).4.3.
Estimation of the shear stress and the total drag
With this boundary layer solution, we are now in a position to provide with an estimationof the wall shear stress ˜ τ = ∂ ˜ U r /∂ ˜ Z (cid:12)(cid:12) ˜ Z =0 , i.e. the viscous component of the stress whichhas been disregarded so far. And indeed this quantity is of paramount importance as faras raindrop-induced erosion of erodible beds is concerned (Ellison 1945; Rein 1993; Lagrée2003; Leguédois et al. J. Philippi, P.-Y. Lagrée and A. Antkowiak
Figure 20.
Left: Isolines for the radial velocity ˜ u r extracted from the numerical simulations.Right: Theoretical isolines for the radial velocity given by the composite expansion ˜ u comp r . Notethat the transverse scale has here been stretched to visualize the boundary layer. for function f as before), we readily obtain:˜ τ (˜ r, ˜ t ) = 2 √ rπ Re / (3˜ t − ˜ r ) . (4.8)This theoretical prediction is confronted Fig. 21 with numerical profiles for the shearstress extracted from Gerris computations, and is shown to nicely agree with observations.From this local distribution for the stress we may infer the total drag induced with adrop impact, by integration over the wet area:˜ D (˜ t ) = Z π Z √ t ˜ τ (˜ t ) ˜ r d˜ r d θ. (4.9)Unfortunately this integral diverges due of the 1 /x singularity developing in the nearcontact line region, and visible from Fig. 21 left. Such singularities are usually a signatureof an additional physics in the diverging region, not taken into account in the model.And indeed, Fig. 21 right reveals that the calculated shear stress significantly deviatesfrom the theoretical prediction at some small distance ∆ from the contact line positionto reach a maximum value. Now integrating the local shear stress up to ˜ r = √ t − ∆,where ∆ is this small cut-off length, we can provide an estimation for the drag at leadingorder in log(∆): ˜ D (˜ t ) = 3 r ˜ tπ Re (cid:18) − (cid:18) ∆ √ ˜ t (cid:19) − (cid:19) . (4.10)Upon noting that this quantity can be dimensionalised with ρU R , the expression forthe total drag in dimensioned variables follows: D ( t ) = 3 √ π µ ρ U R √ t − ∆ q UtR − . (4.11)Noticeably, the departure from the theoretical prediction pinpointed out in Fig. 21 rightseems to occur at a precise location in self-similar variables, therefore suggesting a √ ˜ t timedependence for ∆. From the numerical computations the value of ∆ / √ ˜ t can be estimatedto be around 0 .
03. Note that this is obviously a crude estimation, which nonetheless allowsto propose the following estimate for the impact-induced drag: D ( t ) ’ . µ ρ U R √ t. (4.12)To further refine this prediction, the true nature of the cut-off length ∆ needs to beclearly identified. Several candidates for governing this quantity naturally emerge, withfor example the viscous 1 / Re regularisation length in the vicinity of the contact lineregion or the inertial matching with the Wagner inner layer of typical size ( d ( t ) /R ) ).This requires further investigation. rop impact on a solid surface : short time self-similarity Figure 21.
Left: Numerical and theoretical shear stress distribution underneath the drop, repre-sented in the self-similar space (the numerical data are taken at times ¯ t = 5 × − , 10 − , 5 × − and 10 − ). Right: Same data represented as a function of the distance from the contact line (logplot). This representation reveals a cut-off distance ∆ from which the 1 /x singularity is screened.Importantly the numerical mesh size has been chosen to be small enough (∆ x = 5 × − ) toensure the resolution of the fine-scale motion in the vicinity of the contact line.
5. Further comments and conclusion
Capillary phenomena as well as possible aerodynamic effects from the surrounding gashave been disregarded so far. In this last part we shall estimate their influence on impactand discuss natural extensions of the present work. A summary and general conclusionthen follow in §5.5. 5.1.
Influence of capillary phenomena
The Weber number provides with a global measure of the ratio of available kinetic en-ergy to surface energy. For low values of this number (with respect to unity), dropstypically bounce (Richard & Quéré 2000) while preserving their shape or gently spread(Pasandideh-Fard et al. ∼ et al. (cid:29) e.g. Oliver (2002), the typical extent of this intermediateWagner region associated with highly curved interfaces is found to be O ( ε ), to be com-pared both with the O ( ε ) size of the main impact region considered throughout this paper(see §2.3) and with the O ( ε ) thickness of the lamella. In the framework of our first-ordertheory, we therefore do not anticipate appreciable deviations stemming from this zone.Conversely, for a correct description of the ejected liquid sheet feeding conditions and ofthe pressure fall-off near the lamella root reported Fig. 15, an accurate representation ofthis matching region appears mandatory.The contact line is an other region where marked effects from capillarity are to be ex-pected. Drop impact is characterized with fast motions near the contact line. This violentdynamic wetting phenomenon can arguably bring about issues in our description of im-pact. Actually Blake et al. (1999) demonstrated that nonlocal hydrodynamics could playa significant role in the dynamic contact angle selection. Based on experimental data,8 J. Philippi, P.-Y. Lagrée and A. Antkowiak
Blake et al. further put forward the possible ‘mutual interdependence’ between the phe-nomena in the near contact line region and the far-field hydrodynamics. This complexinterplay was further confirmed in the context of drop impact by Šikalo et al. (2005),but especially for the late receding phase. Interestingly, these authors demonstrated thatthe early evolution of the contact angle was quite insensitive to the experimental con-ditions and fairly well captured by the contact angle of a truncated sphere. This niceagreement certainly advocates for a predominance of inertial effects over capillary cor-rections emanating from the dynamic contact line, at least in the early stage of impact.And indeed, remembering that shortly after impact the fluid motion in the contact areais essentially vertical, it appears likely that the point of contact can be determined withmere inertial arguments. In our simulations, dynamic effects have been disregarded inthe description of the contact angle, which has been set to the constant value π/
2. Theagreement between our simulations and the purely inertial theory is again an indicationof the unimportance of dynamic wetting. It might further be interesting to note that thesurface energy gained by wetting the solid is of the order of 1 / We when rescaled by theinitial kinetic energy. Again this heuristically rules out any leading effect from wettingin the short-term dynamics. This ratio evolves with time though, and ultimately wet-ting phenomena become dominant, as evidenced by the late t / spreading behaviourconsistent with Tanner’s law in the experiments of Rioboo et al. (2002).5.2. Influence of ambient air
For about a decade or so, there has been an increasing realization of the role played bysurrounding air in liquid impact in general, and drop impact in particular. Following keyexperiments performed by Xu et al. (2005) on air-induced splash triggering, a number ofstudies have focused on the events preluding liquid sheet ejection. The first significanteffect of surrounding air is to impart a dimple-like deformation in the bottommost regionof the drop (see experimental observations of the dimple obtained by Thoroddsen et al. et al. et al. (2003) first depicted theoretically this process by couplinglubrication in the squeezed air film and potential flow inside the drop. These authorsnotably evidenced the presence of off-axis pressure peaks. While more recent studiesraised doubt about the link between this dimple formation and splash triggering per se –that might merely be a secondary independent consequence of the presence of surroundinggas (Duchemin & Josserand 2011), this gas pocket is nonetheless formed over timescalesand lengthscales overlapping that of the phenomenon reported in the present paper (Mani et al. et al. t impact corresponding to the moment wherethe drop and the solid are only a grid cell apart. Interestingly our results reveal thatthe pressure at the origin (measuring now the entrapped bubble pressure) is fairly wellcaptured by relation (3.37) after replacing ˜ t with the true time from impact ˜ t − ˜ t impact ,that is: ˜ p (0 , , ˜ t ) = √ π (˜ t − ˜ t impact ) − . (5.1)This agreement between our prediction and a simulation incorporating air entrapmenteffects not only validates and extends our results beyond the initial scope of Wagner rop impact on a solid surface : short time self-similarity Figure 22.
Left: Time evolution of the pressure as measured under an impacting drop withair-induced dimple formation ( i.e. bubble entrapment) taken into account. The red trace moni-tors the pressure at the origin. The physical parameters of this
Gerris simulation are Re = 5000and We = 250. The superimposed black dashed line corresponds to the theoretical solution¯ p (0 , , ¯ t ) = √ π (¯ t − ¯ t impact ) − delayed from ¯ t impact , where ¯ t impact is the real impact time (in thenumerical simulations, ¯ t impact corresponds to the time at which liquid and solid are just one gridcell away). Right: Close-up of the bottommost point of the drop in the numerical simulation (at¯ t = 1 . × − ). The position of the interface is materialised with a blue line. The colormapillustrates the distribution of the pressure field within the drop and in the gas layer. Noteworthyenough the isopressure lines seamlessly cross the interface, revealing the transparency of thedimple to pressure. Note that for the sake of clarity, the vertical scale has here been magnifiedby a factor 22. impact theory (disregarding air effects), but also suggests that the results of the presentmanuscript correspond to the far-field behaviour of an impacting drop in presence ofsurrounding gas. This observation outlines the appealing prospect of describing both thedimple geometry and associated dynamical fields by analytical means.5.3. Main results
In this paper, the short-term dynamics of a drop impacting a rigid substrate has beenelucidated. A self-similar solution for the impact-induced flow has in particular beenunraveled and matched to a self-similar viscous boundary layer. This solution has beenintensively validated with numerical
Gerris computations, and this constant cross-testingbetween asymptotic theory and multiphase adaptive flow simulations is one of the keyfeature of the present approach. In the course of this investigation, several importantresults have been substantiated. These results allow both for a simple yet accurate qual-itative depiction of drop impact along with an in-depth quantitative understanding ofthis phenomenon. These key results are summarised in the following: • A fundamental analogy between the water entry of a solid object (Wagner’s originalproblem) and drop impact exists, • During the earliest moments post-impact, the contact line follows a tank-treadingmotion. There is in particular no contact line sweeping motion, • The impact-induced flow is concentrated in the contact zone, and the far-field merelycorresponds to an undisturbed rigid-body motion reducing to a global free-flight at ve-locity U . There is no global or large-scale drop deformation during impact, • The position of the contact line is given by the simple relation d ( t ) = √ RU t .Though simple, this locus does not correspond to the cut radius of a truncated sphere, • The wet footprint extent of the drop dictates the size of the impact-induced per-turbed flow,0
J. Philippi, P.-Y. Lagrée and A. Antkowiak ( r = 0 , z = O ( d ( t )) > + ) (0 ≤ r < d ( t ) = √ RUt, z = 0 + ) p (0 , z, t ) = √ ρU R π √ t URt + z p ( r, + , t ) = ρU Rπ √ URt − r u r (0 , z, t ) = 0 u r ( r, + , t ) = Urπ √ URt − r u z (0 , z, t ) = Uπ arctan (cid:16) √ URtz (cid:17) − √ U Rt z URt + z − U u z ( r, z, t ) = 0( d ( t ) (cid:29) r, d ( t ) (cid:29) z > + ) (0 ≤ r < d ( t ) = √ RUt, z = O (Re − / d ( t )) u r ( r, z, t ) = π √ U R − r √ t U r ( r, z, t ) = Urπ √ URt − r erf (cid:18)q ρURµ ( URt − r ) z (cid:19) u z ( r, z, t ) = − π √ U R − z √ t τ ( r, t ) = √ ρ U R µ π r URt − r p (0 , , t ) = ρU R π p t F ( t ) = 6 √ ρU R √ t D ( t ) ’ . µ ρ U R √ t Table 2.
Summary of the main results of the paper in dimensioned form. The top part of thetable refers to ideal fluid results (left: closed-form results along the axis of symmetry, right:along the substrate). The left middle part sums up inviscid stagnation point results, and theright middle part summarises the viscous boundary layer results. Observables such as the netnormal force F ( t ) and tangential force D ( t ) are reminded as well in the bottom part of the table. • There is a consistent analogy between the impact-induced flow within the drop andthe flow induced by a flat rising expanding disk (Lamb’s analogy), • The impact pressure is to be associated with the unsteady Bernoulli contribution − ∂ t φ . It cannot be inferred from usual inertial steady contribution − ρU , • As a corollary to the previous point, the impact pressure is extremal at the contactline. It is not maximal at the stagnation point, • A full three-dimensional self-similar solution for the impact-induced flow of an in-viscid drop exists and matches quantitatively realistic numerical data on drop impact, • Analytical solutions for this flow have been presented in integral forms (with someexplicit closed-form expressions along some particular locations), see table 5.3 (in dimen-sional form), • An original inviscid stagnation point structure with an unexpected r/ √ t slip velocitydevelops in the vicinity of the origin. The velocity field structure markedly differs fromthe classic r/t prediction occurring for later times, • An approximate self-similar solution for the viscous boundary layer seamlessly matchwith the inviscid impact-flow (analogy with Mirels shockwave problem), • Self-similar variables have the same structure z/ √ t both in the outer region and inthe boundary layer, • From the knowledge of the distribution of the dynamical fields across the wet area,the expressions for the normal and tangential total force on the substrate are provided, • The asymptotic solution is found to be numerically valid over several decades intime. This solution was found to be insensitive to air-induced dimple formation, • For times of order one, the present results remain at least qualitative. rop impact on a solid surface : short time self-similarity
Perspectives
The results obtained in the present manuscript offer several appealing prospects. On therole of surrounding air first, the last results of §5.2 support the idea that the dimplegeometry and characteristics could be derived analytically, with a far-field correspondingto the here presented flow. The question of the role of air in suppressing the pressuredivergence (Josserand et al. et al. (2013) or granular media, (Lagrée et al.
To conclude
Within the numerous limits carefully drawn along this paper, a consistent asymptoticdescription of the dynamics and geometry of drop impacting a solid surface has beenproposed. The results may simply be summarised through three analogies: Wagner waterentry (drop impact being the dual of this problem), Lamb’s disk winding flow (thataccurately represents the flow induced with the impact) and Mirels shockwave-inducedboundary layer (remarkably capturing the boundary layer developing in the contactline’s trail). The original strategy developed throughout this paper has been to validatethose three analogies through a constant confrontation between numerical simulationsand asymptotic analysis. Our study revealed that very powerful state-of-the-art adaptivecodes now allow to probe all the dynamic features of realistic violent events such asdrop impact, but in the meantime, it also emphasised again how powerful and usefulasymptotic analysis is in providing an in-depth understanding of such phenomena and inuncloaking the raw data delivered by the code. Finally our study brought to light someinteresting features and observables (such as the particular stagnation point structure,pressure distribution, contact line motion, viscous total drag force) never observed to dateneither in simulations nor in experiments. This certainly arouses the exciting prospectof their unveiling in future experimental studies.2
J. Philippi, P.-Y. Lagrée and A. Antkowiak
Figure 23.
Typical mesh structure refined adaptively by
Gerris flow solver during a simulation.
6. Appendix
Gerris flow solver
All the numerical simulations were performed with the open-source code
Gerris (freelydownloadable at http://gfs.sourceforge.net – see also Popinet 2003, 2009; Lagrée et al.
Gerris is a solver of the incompressible Navier-Stokes equationstaking into account multiple phases and surface tension. The code makes use of a finite-volume approach and of a Volume-of-Fluid (VoF) method for an accurate descriptionof the transport of the interfaces between two-phase flows. It also features an adaptivemesh refinement procedure allowing for both a precise description of flows with largescale separation and a reduction of computational costs. Typically in our simulations thefinest grid is chosen to be concentrated along free surfaces and within the contact zone tofully capture the features of the pressure field and of the boundary layers (see Fig. 23). Inthese areas the corresponding local resolution usually corresponds to 4096 × × Gerris , but, to be consistent with the post-impact theorydeveloped in this paper the simulations disregarded air cushioning and dimple formation(except explicitly specified, see §5.2). To avoid dimple formation in this multiphase flowsimulation, the initial configuration is set to a slightly truncated liquid sphere alreadytouching the solid surface. The liquid is initialised with a constant downward velocity.The initial sphere penetration ¯ r = 10 − is at most one grid cell deep (for example thegrid spacing is ∆ x ’ × − for 4096 × × × p (0 , , ¯ t ) = √ π ¯ t − (see equation (3.37)) around ¯ t = 5 × − and leaves the self-similar regime at around¯ t = 6 × − . We remark that after a transient period, both numerical solutions give con-sistent information and collapse onto the theoretical solution over almost three decades.Note that the occurrence of sporadic glitches in the numerical solution (see e.g. Figs rop impact on a solid surface : short time self-similarity Figure 24.
Comparison between the pressure measured in the
Gerris simulations at the originfor two different maximum level of mesh refinement (red dots) and the theoretical prediction(black dashed line). The left panel corresponds to a simulation where the maximal grid densityis 2048 × × p (0 , , ¯ t ) = √ π ¯ t − , see equation (3.37). After a short transient, bothsimulations quickly reach the same self-similar asymptotic regime. Figure 25.
Close-up of Fig. 24 right. Note that the numerical evolution for the pressure isslightly above (about 7 %) the analytical solution.
17, 24 or 25) are to be related with the classic difficulty of computing the pressure inprojection methods, such as the one implemented in
Gerris (Brown et al.
Gerris parameter file
A minimal
Gerris parameter file allowing to reproduce the results presented in this manuscriptis provided here for the reader’s convenience (Note for referees and Editors: this sectionwill likely be suppressed in the final version for this script will be on a dedicated examplepage in
Gerris website):
Define t0 1e-4Define Re 5000Define ReAir 277778Define We 250Define VAR(T,min,max) (min + CLAMP(T,0,1)*(max - min))Define RHO_EAU 1000.Define RHO_AIR 1.Define RHO(T) VAR(T, RHO_AIR/RHO_EAU, RHO_EAU/RHO_EAU) J. Philippi, P.-Y. Lagrée and A. Antkowiak rop impact on a solid surface : short time self-similarity GfsBox { right = Boundary{BcDirichlet P 0BcNeumann U 0BcNeumann V 0 }top = Boundary{BcDirichlet P 0BcNeumann U 0BcNeumann V 0 }left = Boundary{BcDirichlet U 0BcDirichlet V 0 }}1 2 right1 3 top3 4 top
Thanks
We warmly thank Stéphane Popinet, Pascal Ray, Christophe Josserand and StéphaneZaleski from Institut Jean le Rond ∂ ’Alembert for enlightening discussions on splashphenomena, and Neil Balmforth from the Department of Mathematics of the Universityof British Columbia for interesting discussions on boundary layers. REFERENCESAntkowiak, A., Audoly, B., Josserand, C., Neukirch, S. & Rivetti, M.
Proc. Natl Acad. Sci.U.S.A. (26), 10400–10404.
Blake, T. D., Bracke, M. & Shikhmurzaev, Y. D.
Physics of Fluids (1994-present) (8), 1995–2007. Brown, David L, Cortez, Ricardo & Minion, Michael L
Journal of Computational Physics (2),464–499.
Chen, Q., Ramé, E. & Garoff, S.
J. FluidMech. , 49–66.
Darrozès, J.-S. & François, C.
Mécanique Des Fluides Incompressibles (Lecture Notesin Physics) . Springer-Verlag.
Duchemin, Laurent & Josserand, Christophe
Physics of Fluids (1994-present) (9), 091701. Dussan V., Elizabeth B. & Davis, Stephen H.
Journal of Fluid Mechanics , 71–95. Eggers, Jens, Fontelos, Marco A., Josserand, Christophe & Zaleski, Stéphane
Physics of Fluids(1994-present) (6), –. Ellison, W. D.
Eos, Transactions American Geophysical Union (3), 415–429. Engel, Olive G
J. Res. Nat. Bur. Stand. (5),281–298. Howison, S. D., Ockendon, J. R., Oliver, J. M., Purvis, R. & Smith, F. T.
J. Fluid Mech. (-1), 1–23.
Howison, S. D., Ockendon, J. R. & Wilson, S. K.
J. Fluid Mech. , 215–230.
Josserand, Christophe, Ray, Pascal & Zaleski, Stephane
Manuscript in preparation . Josserand, Christophe & Zaleski, Stephane
Physics of Fluids (6), 1650–1657. J. Philippi, P.-Y. Lagrée and A. Antkowiak
Korobkin, A.A.
Journal of EngineeringMathematics (1-4), 121–139. Lagrée, Pierre-Yves
Physics ofFluids (1994-present) (8), 2355–2368. Lagrée, Pierre-Yves
Interactive Boundary Layers in: Asymptotic Methods in FluidMechanics: Survey and Recent Advances (CISM International Centre for Mechanical Sci-ences) . Springer.
Lagrée, P.-Y., Staron, L. & Popinet, S. µ ( i )-rheology. J. Fluid Mech. ,378–408.
Lagubeau, Guillaume, Fontelos, Marco A., Josserand, Christophe, Maurel, Agnès,Pagneux, Vincent & Petitjeans, Philippe
Journal of Fluid Mechanics , 50–60.
Lamb, Horace
Hydrodynamics . Cambridge University Press.
Lee, Ji San, Weon, Byung Mook, Je, Jung Ho & Fezzaa, Kamel
Phys. Rev. Lett. , 204501.
Leguédois, Sophie, Planchon, Olivier, Legout, Cédric & Le Bissonnais, Yves
Soil Sci. Soc. Am. J. pp. 30–37.
Mangili, S., Antonini, C., Marengo, M. & Amirfazli, A.
Soft Matter , 10045–10054. Mani, Madhav, Mandre, Shreyas & Brenner, Michael P.
Journal of Fluid Mechanics (-1), 163–185.
Mirels, Harold
NACA
TN 3401 , 25p.
Nethercote, WCE, Mackay, M & Menon, BC
Some warship slamming investigations .Defence Research Establishment Atlantic.
Oliver, James Mark
Oliver, J. M.
J. Fluid Mech. , 59–86.
Pasandideh-Fard, M., Qiao, Y. M., Chandra, S. & Mostaghimi, J.
Physics of Fluids (1994-present) (3), 650–659. Popinet, Stéphane
Journal of Computational Physics (2), 572 – 600.
Popinet, Stéphane
Journal of Computational Physics (16), 5838–5866.
Pruppacher, H. R. & Beard, K. V.
Quarterly Journal of theRoyal Meteorological Society (408), 247–256. Rein, M.
Fluid DynamicsResearch (2), 61–93. Renardy, Y., Popinet, S., Duchemin, L., Renardy, M., Zaleski, S., Josserand, C.,Drumright-Clarke, M. A., Richard, D. & Quéré, D.
J. Fluid Mech. , 69–83.
Richard, D. & Quéré, D.
EPL (Europhysics Letters) (6),769–775. Rioboo, R., Marengo, M. & Tropea, C.
Experiments in Fluids (1), 112–124. Schlichting, Hermann
Boundary-layer theory . New York: McGraw-Hill.
Šikalo, Š., Wilhelm, H.-D., Roisman, I. V., Jakirlić, S. & Tropea, C.
Physics of Fluids (1994-present) (6), –. Smith, F. T., Li, L. & Wu, G. X.
J.Fluid Mech. , 291–318.
Sneddon, Ian N.
Glasgow Mathemat-ical Journal , 108–110. Sneddon, Ian N.
Fourier Transforms . Dover. rop impact on a solid surface : short time self-similarity Staron, Lydie, Lagrée, Pierre-Yves, Ray, Pascal & Popinet, Stéphane
Journal of Rheology (1978-present) (4). Stebnovskii, S.V.
Journal of Applied Mechanics and Technical Physics (1), 66–69. Stewartson, Keith
The Theory of Laminar Boundary Layers in Compressible Fluids .(Oxford Mathematical Monographs), Clarendon Press.
Stow, C. D. & Hadfield, M. G.
Proceedings of the RoyalSociety of London. Series A, Mathematical and Physical Sciences (1755), 419–441.
Stow, C. D. & Stainer, R. D.
Journalof the Meteorological Society of Japan. Ser. II (5), 518–532. Thoroddsen, S. T., Etoh, T. G., Takehara, K., Ootsuka, N. & Hatsuki, Y.
J. Fluid Mech. , 203–212.
Van Dyke, Milton D
Perturbation methods in fluid dynamics . Parabolic Press.
Villermaux, E. & Bossa, B.
J. Fluid Mech. , 412–435.
Wagner, Herbert
ZAMM - Journal of Applied Mathematics and Mechanics / Zeitschrift für AngewandteMathematik und Mechanik (4), 193–215. Weiss, Daniel A. & Yarin, Alexander L.
J. Fluid Mech. ,229–254.
Xu, Lei, Zhang, Wendy W. & Nagel, Sidney R.
Physical Review Letters (18), 184505–4. Yarin, A. L. & Weiss, D. A.
J. Fluid Mech.283