Droplets move over viscoelastic substrates by surfing a ridge

Liquid drops on soft solids generate strong deformations below the contact line, resulting from a balance of capillary and elastic forces. The movement of these drops may cause strong, potentially singular dissipation in the soft solid. Here we show that a drop on a soft substrate moves by surfing a ridge: the initially flat solid surface is deformed into a sharp ridge whose orientation angle depends on the contact line velocity. We measure this angle for water on a silicone gel and develop a theory based on the substrate rheology. We quantitatively recover the dynamic contact angle and provide a mechanism for stick–slip motion when a drop is forced strongly: the contact line depins and slides down the wetting ridge, forming a new one after a transient. We anticipate that our theory will have implications in problems such as self-organization of cell tissues or the design of capillarity-based microrheometers.

C apillary interactions of soft materials are ubiquitous in nature and play a major role in self-organization of cell tissues 1 , in embryotic development 2,3 , in wound healing 4 , and in controlling the spreading of cancer cells 5,6 . Not least motivated by this, soft wetting [7][8][9] recently came to the attention of both physicists and biologists. Despite its potential for applications, such as the patterning of cells 10 or droplets 11 onto soft surfaces, or the optimization of condensation processes 12 , our fundamental understanding of the dynamics of soft wetting lags behind by far of what is known about rigid surfaces 13 .
Partial wetting of a liquid on a rigid (smooth) substrate is controlled by intermolecular interactions, whose strength is characterized by surface energies 13 . The motion of the threephase contact line is governed by the viscous dissipation in the liquid. A dissipation singularity arises at the moving contact line 14 and its regularization at the nanoscopic scale can result from various processes 13,15 . When the substrate is deformable, a sharp ridge forms below the contact line at the edge of the droplet [7][8][9] . The ridge geometry (Fig. 1a) originates from the coupling between elasticity and surface energy [16][17][18][19][20] . The problem is inherently multi-scale and non-local, even at equilibrium, due to the long range of elastic interactions [21][22][23] .
Pioneering experiments have shown that the softness drastically slows down the wetting dynamics 24,25 in comparison to rigid solids. This viscous braking in the wetting of soft solids has been attributed to a viscoelastic force, as discussed in several recent experimental articles 11,26 . The theoretical description of moving contact lines over soft solids is so far limited to global dissipation arguments 17 , which, at least for wetting of rigid solids, are known not to capture the entire physics behind the process.
In this paper, the physical mechanism that governs soft wetting dynamics is revealed. We measure the dynamical wetting of small water droplets on a rheologically characterized silicone gel and discover a saturation of the dynamical contact angle for large speeds, associated with a maximum friction force. Driving the contact line motion beyond this maximal force eventually leads to a dynamical depinning where the contact line surfs down the wetting ridge, providing a mechanism for recently observed stick-slip motion 26 . We develop a theoretical framework for dynamical soft wetting, suitable for any substrate rheology. The dynamic wetting angle is calculated from the velocity-dependent shape of the deformed solid (Fig. 1c). The experimental results are matched quantitatively, including saturation/dynamical depinning. The latter arises from an upper limit of the viscoelastic braking effect, which, by exploring different rheologies, is found to be a robust, universal feature of soft wetting and should thus be relevant far beyond droplets on silicone gels. In addition, the analysis captures recent X-ray measurements on the slow growth of wetting ridges when a drop is deposited on a substrate 27 (Fig. 1b).

Results
Experiments. Experiments were performed using water drops on a silicone gel (cf. Methods section for details). This system was previously used in static 28 and transient 27 experiments. Figure 2a shows the rheology of this gel; similar data were reported in ref. 29. The storage G 0 and loss G 00 moduli are related by Kramers-Kronig relation: they originate from the same relaxation function C(t). More precisely, the complex shear modulus obey the relation m o ð Þ G 0 þ iG 00 ¼ io A silicone gel is a reticulated polymer formed by polymerizing small multifunctional prepolymers: contrarily to other types of gels, there is no liquid phase trapped inside. Such cross-linked polymer networks exhibit scale invariance that yields power-law response of the form 17,30,31 : where G is a static shear modulus and G is the gamma function. The associated complex modulus reads m The value of n is not universal but depends on the stoichiometric ratio r between reticulant and prepolymer, with n typically between 2/3 and 1/2 (ref. 32). The best fit in Fig. 2a gives an  Lines are best fit of m ¼ G(1 þ (ito) n ), giving n ¼ 0.55, G ¼ 1.2 kPa and t ¼ 0.13 s. (b) Dynamic angle j ¼ y À y eq for water on the silicone gel (symbols). Data are averaged over 10 independent experiments, error bars represent the s.d. The small-v behaviour exhibits the same power-law as G 00 . Dashed line is the best fit of the asymptote equation (8). Solid line corresponds to equation (7), describing the full range of velocities. The critical angle of depinning j crit ¼ 39 AE 3 , measured separately, is plotted at an arbitrary velocity. exponent n ¼ 0.55. Note that the associated effective viscosity G 00 /oBG(t/o) n is very large, beyond 10 Pa.s, over the entire frequency domain. This will imply that the dissipation mainly occurs in the solid, not in the liquid. Figure 2b shows the dynamical angle j as a function of velocity, both of which are measured while the droplet slowly relaxes over time towards its equilibrium shape (after the injection phase). The resulting j versus v is not sensitive to the history of the relaxation process, and the dynamics can thus be considered quasi-steady. The log-log plot reveals a power-law relation between j and v at small velocities, with an exponent equal to the rheological value of n ¼ 0.55, within error bars. Such a power-law dependence is similar to previously obtained results 24,25 . For large velocities, we report a striking saturation of the dynamical contact angle. Neither the small-v power law nor the saturation can be explained by dissipation mechanisms in the liquid, and one needs to account for the dissipation within the solid. Long et al. 17 addressed this using a global dissipation argument based on the solid rheology, but this argument fails to capture key features such as the saturation.
When the drop is kept inflating with a large over-pressure, we observe a 'depinning' of the wetting front, with a sudden increase of its velocity, as the dynamical angle reaches a value j crit % 39 AE 3 . This compares well to the saturation of 37°o bserved during relaxation (after injection), as indicated in Fig. 2.
When forcing the contact angle beyond this angle, the contact line dynamically depins from the wetting ridge, surfing it, until a new ridge forms. Note that such a depinning, leading to stick-slip motions, had already been reported in refs 26,33. Our current measurements show that this is a direct consequence of the saturation of the dynamic contact angle.
Theoretical framework. A liquid drop deposited on a substrate exerts a capillary traction on the surface 7,34-37 . While the resulting elastic deformation has been computed and measured for static situations 8,9,18,20,23,28 , the traction becomes time dependent in the case of dynamical wetting. Here we consider a single straight contact line, for which the elastic problem is two-dimensional. Our goal is to find the deformation of the solid, h(x, t), resulting from the time-dependent capillary traction, T(x, t). For simplicity, we consider the same surface energy g s for the wet and the dry parts of the substrate, and assume that there is no Shuttleworth effect, i.e. that g s does not depend on strain 16,19 . The resulting traction on the solid is then purely normal, and reads T(x, t) þ g s q xx h, the latter term being the solid Laplace pressure 9,23,38,39 . The theory is rigorous for small slopes (q x h) 2 o o1, but can be extended in a semi-quantitative way to finite slopes. The shape of the deformed substrate h(x, t) follows from the normal substrate displacements. Inside a purely elastic material, the displacements adapt instantaneously to changes in the capillary traction; the problem is therefore essentially static. For realistic soft materials, however, the displacements are delayed with respect to the imposed forcing. For small deformations, this is captured by a linear stress-strain rate relation where C is the relaxation function previously introduced (see equation (1)) and p is the pressure. Like ref. 17, we apply a Fourier transform with respect to time (noted by '^'): The mathematical problem defined by mechanical equilibrium, r Áŝ ¼ 0, the constitutive relation (equation (3)) and the traction at the free surface is identical to the static problem, but features dependences on the frequency o. The time-dependent traction can therefore be solved analogously to (refs 9,22,23), by an additional spatial Fourier transformation (noted by 'B'): where q is the wavenumber. The Green's functionĜðq; oÞ is the product of the time kernel m(o) À 1 by the space kernel K(q). For an incompressible layer of thickness Left-right symmetry and volume conservation are reflected by K(q) ¼ K( À q) and K(0) ¼ 0. Sharp features in the solid profile, like the solid contact angle, are found in the large-q asymptotics for which K(q)C(2|q|) À 1 .
The moving contact line. We now apply our theory to a contact line moving at a constant velocity v, which induces a traction This reflects the normal force per unit contact line that is exerted by the liquid on the solid, while y is the liquid angle at the location of the cusp. For simplicity we consider that the drop size is much larger than the substrate thickness, in which case the Laplace pressure inside the liquid can be neglected 28 . We indeed verified that the finite drop size has a negligible influence on the resulting motion: the relevant scale for the dynamics is the size of the ridge g s /G, which is much smaller than the drop size. This also justifies a two-dimensional model. Another important simplification comes from the quasi-steady nature of the droplet relaxation: temporal changes of contact angle and contact line velocity are small in our experiments (dy/dto ot À 1 ), so that the process can be modelled by a constant velocity. According to equation (4), the capillary traction induces a wetting ridge moving at a velocity v (see Methods section): in the co-moving frame. In real space this gives profiles such as shown in Fig. 1c. The motion induces a left-right symmetry breaking: the asymmetric deformation of the solid results in a tilt angle jðvÞ of the cusp. Since the liquid is close to equilibrium, the change in the solid angle j induced by the dissipation taking place inside the solid directly yields a change in the liquid angle y (Fig. 1d). Hence y ¼ y eq þ j, where y eq is the equilibrium liquid angle by Neumann's law. The liquid contact angle y gets deviated away from y eq due to the viscoelastic forces in the substrate. Figure 2b shows the calculated tilt curve for the gel used in experiments. It predicts not only the power-law behaviour at low velocity but also presents a saturation of the tilt angle. The tilt angle quantifies the velocity-dependent viscoelastic force between the solid and liquid phases. For a well-established moving ridge, it behaves as a resistive force increasing with the velocity. When the drop is forced to inflate with a driving force larger than the maximal braking force, the contact line can no longer remain 'pinned' to the steadily moving solid ridge and surfs the gel wave. To investigate further the relation between the tilt j and the substrate constitutive relation, we use the gel rheology (equation (1)), and expand (equation (7)) in the small-v asymptotics (and hence small j, that is, sinyEsiny eq ), which gives (see Methods section): where the characteristic velocity scale emerges as v* ¼ g s /(Gt). Note that the outer length scale (thickness of substrate) does not appear. This expression can be simply interpreted. At vanishing response time t, a deformation matching the static ridge would propagate at a velocity v, pushing the substrate material up and down at a characteristic frequency o equal to the velocity v divided by the characteristic width of the ridge Bg s /G. The perturbation introduced by a finite t is encoded by the loss modulus G 00 (o). As the characteristic strain is set by the slope of the ridge Bgsiny/g s , one obtains dimensionally equation (8).
The scaling law j / ðv = v Ã Þ n thus simply carries over the low frequency behaviour G 00 (o)po n , which is a robust mechanism valid beyond the small-slope approximation of our theory. At small v (small o), dissipation will dominantly occur in the solid because no1, while the loss modulus of a Newtonian liquid G 00 liq / o. In contrast to the tilt angle j, we find that that the solid angle y s does not depend on v. This result can be derived analytically from the large-q asymptotics of equation (7),h ' g siny g s q À 2 , valid for all v and arbitrary m(o). In real space, this implies a slope discontinuity which is the small-slope limit (g/g s o o1) of Neumann's contact angle law. Physically, equation (9) reflects that y s is determined by surface tensions only 9,18,41 : bulk viscoelastic stresses are not singular enough to contribute to the contact angle selection. This feature remains true for arbitrary angles 23 . To match Neumann's law quantitatively, our theory must be corrected at large slopes to take geometric nonlinearities into account. This can be achieved phenomenologically in equation (4) by replacing g s ! ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ffi g 2 s À g 2 4 q . Indeed, the Neumann condition for y ¼ 90°r eads 2g s sin a ¼ g, where a is the angle of the solid interface with the horizontal. Small-slope theory gives 2 h 0 j j¼ 2tana ¼ g g s , and hence lacks a factor cos a ¼ ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi . For the first time, we reveal that the exponent of the dynamical contact angle directly originates from the gel rheology. The rheological parameters being calibrated independently, the dynamic contact angle can be fitted to the model to extract the solid surface tension. Using equation (8), which is valid for small slopes, we find g s ¼ 16 mN m À 1 . This is a reasonable value, though a bit lower than the value previously derived from Neumann's law 28 . We think this difference can be attributed to the small-slope nature of our theory: condering the phenomenological correction for large slopes gives a value g s ¼ 39 mN m À 1 , in close agreement with ref. 28. The solid line in Fig. 2b shows the prediction from equation (7), providing an excellent description over the full range of velocities. The model captures also the saturation, though the value for j crit is slightly overestimated.
Depinning and growth of a new wetting ridge. How can the contact line escape pinning, without dragging the capillary wedge along with it? To answer this question, let us consider the recent experiments investigating the growth of a wetting ridge after depositing a droplet on a silicone gel 27 . The substrate was observed to only very slowly establish the final shape of the wetting ridge-such a delay in growth (or decay) of wetting ridges would explain how a sufficiently rapid contact line could escape from the ridge. However, the solid angle y s (cf. Fig. 1a) appeared very quickly and remained constant during the entire growth of the ridge 27 .
These features of ridge growth can all be explained by considering our theory for a traction that is suddenly imposed at the time of deposition (t ¼ 0), so that where y is the liquid contact angle and Y(t) is the Heaviside step function. Combining this traction with equation (4), one can compute the resulting h(x, t) for any rheology m(o). An example of the evolution of the wetting ridge is shown in Fig. 1b. A movie is given in Supplementary Movie 1. First, the theory recovers the experimental finding that y s is constant at all times. As for the moving contact line, this result can be derived analytically from the large-q asymptotics of equation (4), which again results in equation (9): the asymptotics are independent of the rheology and the history of the traction, but entirely governed by the surface tensions.
Second, the theory explains why, contrarily to the rapid appearance of y s , the global shape of the ridge evolves much more slowly. Figure 3 shows the evolution of the central height of the ridge, h(x ¼ 0), towards its static value h N , for the two idealized rheological models. The relaxation towards the equilibrium height is algebraic for the gel model, with an exponent directly following that of rheological relaxations (Fig. 3a, as t À 1/2 for n ¼ 1/2). This clarifies the complex evolution of the wetting ridge of the silicone gel in ref. 27: small-scale characteristics such as y s are dominated by surface tension and relax quickly, while largescale features inherit the relaxation dynamics from the bulk rheology. This means that immediately after depinning, where the contact line exhibits a rapid motion, the solid cusp cannot adapt quickly. The liquid will slide down the wetting ridge, which appears frozen on the timescale of the depinning. During this phase it is clear that the liquid dynamics will be important-still, the onset of the depinning can be explained quantitatively without invoking the fluid dynamics inside the liquid, because the saturation angle j max coincides with the observed depinning angle j crit (Fig. 2b).
Robustness and interpretation. The theory can be applied to any viscoelastic substrate, assuming that it is probed in the linear regime. Generic reticulated polymer networks possess a long-time entropic elasticity 42,43 , that is characterized by a static shear modulus G. Such networks become viscoelastic when excited over timescales shorter than a certain response time t. To investigate the robustness of the phenomenology that was observed experimentally and reproduced quantitatively by our model, we will consider a different rheological limit. When cross-linking long polymer chains, one forms an elastomer rather than a gel. Assuming a single (Rouse) timescale t to characterize the onset of entanglements, the rheology can be idealized as 17,31,43 : This single timescale response is also referred to as standard linear model and has frequently been used to describe the transition from rubber to glass behaviours. In general, several relaxation times must be introduced to capture quantitatively the rheology of actual elastomers. The wetting ridge relaxation, which follows the rheological relaxation, becomes exponential for the standard linear model (Fig. 3b). The case of a moving contact line with the rheology (equation (11)) is given by the solid lines in Fig. 4, showing the tilt curves for various parameters b. As for the gel case, the tilt has an upper bound and this therefore appears a robust feature of soft wetting. Note that the maximum depends on g/g s and b (or n for the gel case), illustrating that the value of j max will in general depend on the details of the rheology (for example, j max in ref. 26 is much smaller than that in the presented experiments). The dynamic contact angle for a standard linear solid can actually be captured in a simple analytical form. For this we consider the limit b-N while keeping t 0 bt constant-this corresponds to the Kelvin-Voigt model with a frequency-independent effective viscosity Z ¼ Gt 0 . Intriguingly, this limit turns out to be singular: the high-frequency behaviour of equation (11) becomes purely viscous and gives a non-integrable singularity of the dissipation. This singularity could already be anticipated from equation (8), since the Kelvin-Voigt rheology has G 00 Bo 1 , while equation (8) presents a divergence for n ¼ 1.
In fact, this viscoelastic singularity is the soft-solid analogue of the classical Huh and Scriven paradox for viscous contact line motion 14 . This is demonstrated from the large-q asymptotics of equation (7) in the Kelvin-Voigt limit, giving a slope close to the contact line (see the Methods section): with gE ¼ Euler's constant and v* ¼ g s /(Gt 0 ). This expression reveals a logarithmic divergence of the slope, in perfect analogy to the Cox-Voinov result for liquid contact lines 44,45 .
Contrarily to the viscous-liquid singularity, the presence of an instantaneous elastic response, that is, a finite value of b, is sufficient to regularize the divergence. This is illustrated in Fig. 5, which shows the slope ahead of the moving contact line as a function of the distance to the contact line. The dashed line corresponds to the Kelvin-Voigt limit equation (12), showing the logarithmic steepening of the slope. For finite b, the slope saturates on approaching the moving contact line. The saturation wavenumber is found q 0 B(vt) À 1 , which corresponds to a length where ' is a dynamical regularization length that depends linearly on the velocity of the contact line; for v $ Oð1Þ this scale is still much smaller than the substrate thickness, by a factor b À 1 . The physical origin of the regularization lies in the instantaneous elasticity in the high-frequency limit, which applies at frequencies beyond Bbt.
Inserting the regularization length into the Kelvin-Voigt limit equation (12), we identify the tilt and get an analytical expression of the dynamic liquid contact angle (the strict validity of the analysis requires small slopes, that is, small j; we therefore replaced siny ¼ sinðp = 2 þ jÞ ' 1): This result is analogous to the Cox-Voinov law 44,45 in partial wetting of viscous fluids. In that case, a similar logarithmic factor linking microscopic and macroscopic scales appears, for arbitrary contact angles 45 , and the resulting expression for small Ca is of the form equation (14). Interestingly, the analysis reveals that the relevant dimensionless velocity for soft wetting is not the classical liquid capillary number Ca ¼ vZ ' = g, based on the liquid viscosity Z ' , but the 'solid capillary number' Equation (14) closely follows the numerical results (Fig. 5, dashed and solid lines, respectively). The moving contact line singularity is avoided altogether when G 00 has an exponent no1, as was the case for the power-law gel, in perfect analogy to shear-thinning fluids moving on a rigid substrate.

Discussion
We have shown how contact lines can surf on a wetting ridge, and that this governs the remarkable spreading of drops on viscoelastic substrates. We have quantified this dynamics by  (14)). The upper bound for the viscous braking force is robust with respect to the details of the rheology. measuring the dynamic contact angle of water on a silicone gel for a wide range of velocities and described, for the first time, a saturation of the dynamical angle for large velocities. This saturation is in harsh contrast to wetting dynamics of rigid solids and leads to depinning, where the contact line slides down the ridge until a new wetting ridge has had time to grow and sustain a steady motion-hence, explaining the remarkable stick-slip motion 26 found recently on soft solids. We develop a theory that identifies a robust maximum in viscous braking force that correctly predicts the onset of dynamical depinning. In addition, our theory captures the unsteady growth of a wetting ridge 27 . This work provides a framework for viscoelasto-capillary dynamics valid beyond droplets, and should be applicable, for example, within a biological context. It also opens a new perspective where droplets can be used as microrheometers, since the length scales probed by the droplets are given by the elastocapillary length (that is, a few microns), and the tilt saturation occurs at velocities that are directly related to the relaxation timescale.

Methods
Wetting experiments. The silicone gels (Dow Corning CY52-276) are prepared by curing the mixed components onto glass slides, yielding 0.8-mm thick substrates. The rheology was determined using a MCR 501 rheometer (Anton Paar). Dynamic contact angles were measured using droplets of MilliQ water dispensed from a clean Hamilton syringe. First, a small droplet (2-20 ml) was placed onto the substrate, leaving the syringe needle attached to the droplet. Then, the contact angle of the droplet was increased by quickly injecting water (3-20 ml with 2-8 ml s À 1 ) into it. After the injection phase the drop relaxes quasistatically, causing the contact line velocity to decay slowly. The advancing motion of the contact line and the relaxation of the contact angle were imaged at 50 Hz with a long-distance video microscope. The droplet contour was extracted with sub-pixel resolution, and velocities down to Bnm s À 1 could be detected. The measured contact angles were translated to tilt angles j by subtracting y eq E106°±1°.
The moving contact line. Fourier transforming equation (6) from x to q and from t to o preserves the d-shape of the traction: Inserting the above into equation (4), the inverse transform to the time domain yieldsh The only explicit time dependence appears in the phase factor that shifts the profile in x-direction linearly with time. The transformation to the co-moving frame is done by multiplication with e À iqvt , which cancels the only explicit time dependence, and one obtains equation (7). The slopes are evaluated by multiplication with À iq before the inverse transform (in the co-moving frame): h 0 (x) is a real function because <½ À iqh q ð Þ ¼<½iqh À q ð Þ and I½ À iqh q ð Þ ¼ À I½iqh À q ð Þ. h 0 (x) can be split into a symmetric and an antisymmetric part, where the symmetric part is given by the inverse transform of the real part of iqhðqÞ: The antisymmetric part is obtained form the imaginary part: The solid angle y s is given by the (antisymmetric) slope discontinuity at x ¼ 0 and is thus encoded in the backward transform of the imaginary part. The discontinuity is caused by the large-q asymptotics alone.
Þo O q ð Þ, which is the case for the exponential-and power-law (no1) relaxation (but not for the Kelvin-Voigt model), it is independent of rheology: lim x!0 þ h 0 ðxÞ À lim x!0 À h 0 ðxÞ The rotation of the wetting ridge is given by the symmetric part of h 0 (x) and thus obtained by the backward transform of the real part, evaluated at x ¼ 0: With the symmetry property K(q) ¼ K( À q) and with small, positive v, equation (22) simplifies to (primes omitted): j ¼ g siny G sin np = 2 p Z 1 0 qðqvtÞ n KðqÞ ðg s = GÞq 2 KðqÞ þ 1 ð Þ 2 dq: ð23Þ In the limit of thick elastic layers, K(q) ¼ (2|q|) À 1 . After non-dimensionalizing the integration variable as q 0 ¼ g s G q, one obtains (primes omitted): where v* ¼ g s /(Gt) is the characteristic velocity.
Growth of a wetting ridge. Here we give the full derivation of the time-dependent wetting ridge shape after the deposition of a droplet. We only discuss the result for the exponential relaxation model. An analogous calculation can be performed for the power-law relaxation.
In the following, we non-dimensionalize x with h 0 , q with h À 1 0 , t with bt, o with (bt) À 1 and h with g siny/G. With this scaling, the Fourier transform of the time kernel for exponential relaxation (equation (11)) reads The space kernel in scaled variables is The traction equation (10) is transformed tỗ equations (16)- (18) are inserted into the general expression equation (3), which yields:ĥ with the dimensionless parameter a s ¼ g s /(Gh 0 ). a s compares the elastocapillary length for the solid surface tension to the layer thickness h 0 . The inverse Fourier transform to the time domain yields Fourier transformation to real space is performed numerically.