1. Introduction
The motion of a contact line governs a variety of phenomena, such as drops moving on surfaces, drop spreading and coating applications. This problem has been widely explored for Newtonian fluids over the last decades (de Gennes, Brochard-Wyart & Quéré Reference de Gennes, Brochard-Wyart and Quéré2004; Bonn et al. Reference Bonn, Eggers, Indekeu, Meunier and Rolley2009; Snoeijer & Andreotti Reference Snoeijer and Andreotti2013). Within the hydrodynamic framework of wetting, the macroscopic flow is often described using the Cox–Voinov theory (Voinov Reference Voinov1976; Cox Reference Cox1986), which connects the macroscopic contact angle to the microscopic angle. The relation between these angles is intricate, since the balance of the viscous shear forces (viscosity 
$\eta$) and capillary forces (surface tension 
$\gamma$) leads to a ‘bending’ of the interface in a region near the contact line. Macroscopically, the result can be expressed as a relation for the apparent outer (macroscopic) angle 
$\theta _{{app},o}$, defined in figure 1(a,b) for droplet spreading and dip coating, respectively. An explicit treatment of the problem is often given using the thin film (lubrication) framework, see e.g. Eggers & Fontelos (Reference Eggers and Fontelos2015), which leads to a dependence of the macroscopic angle on the contact-line speed 
$U$ through the form (Voinov Reference Voinov1976; Cox Reference Cox1986)
\begin{equation} \theta_{{app},o}^3 = \theta_e^3 \pm 9 {Ca} \ln(\ell_o/\ell_i), \quad {Ca} = \frac{\eta U}{\gamma}. \end{equation}
Here, 
$\theta _e$ is the (microscopic) equilibrium angle, and 
${Ca}$ the capillary number. The ‘
$\pm$’ depends on whether the contact line is advancing (
$+$) or receding (
$-$). The formula contains an inner length scale 
$\ell _i$, necessary to regularise the moving contact-line singularity (Huh & Scriven Reference Huh and Scriven1971), as well as an outer length scale 
$\ell _o$ that reflects the macroscopic geometry of the problem (e.g. the radius of a spreading drop). The strict validity of (1.1) is restricted to small contact angles and small capillary numbers, due to assumptions underlying the analysis (Voinov Reference Voinov1976; Cox Reference Cox1986; Eggers & Fontelos Reference Eggers and Fontelos2015). Yet, it has been successfully applied to a large class of problems (Snoeijer & Andreotti Reference Snoeijer and Andreotti2013).
Moving contact lines for viscoelastic liquids encountered in drop spreading (a) and dip coating (b). At a large scale, the flow is characterised by a macroscopic apparent contact angle 
$\theta _{{app},o}$, seen at an outer scale 
$\ell _o$. (c) Microscopic view of the interface. 
$({c}_i)$ The highest shear rates are encountered close to the contact line (slip length scale 
$\lambda$), where polymers become highly stretched. This gives rise to large normal stresses that ‘bend’ the interface to angles below the equilibrium contact angle 
$\theta _e$. In this paper we introduce 
$\theta _{{app},i}$ as the resulting apparent microscopic angle. 
$({c}_{ii})$ In a scale-free region 
$\lambda \ll x \ll \ell _o$, the interface slope follows the Voinov solution reflecting the visco-capillary balance; 
$\theta _{{app},i}$ serves as the apparent inner boundary condition at small scale. 
$({c}_{iii})$ At large scales, the slope is given by the apparent outer angle 
$\theta _{{app},o}$.
Many applications like the deposition of pesticides on plant leaves or ink-jet printing involve polymer solutions. Such fluids exhibit non-Newtonian rheology such as shear-thinning viscosity and viscoelasticity, which can affect the contact-line dynamics, as has been demonstrated through several experimental studies (Seevaratnam et al. Reference Seevaratnam, Suo, Ramé, Walker and Garoff2007; Wang et al. Reference Wang, Zhang, Lee and Peng2007; Wei, Garoff & Walker Reference Wei, Garoff and Walker2009a; Wei et al. Reference Wei, Rame, Walker and Garoff2009b; Han & Kim Reference Han and Kim2013; Kim & Rothstein Reference Kim and Rothstein2015). For the canonical problem of liquid spreading on a substrate, it was found that non-Newtonian drops spreading on completely wetting surfaces still followed the long-time power-law dynamics 
$R \sim t^{n}$, where 
$R(t)$ is the drop radius, and the exponent 
$n$ is slightly smaller than the Cox–Voinov–Tanner law 
$n=1/10$ (Rafaï, Bonn & Boudaoud Reference Rafaï, Bonn and Boudaoud2004). This comparatively weak effect could be attributed to both a shear-thinning viscosity and the consequence of normal stress effects. During the very early stages of drop spreading, ensuing after a droplet is gently brought into contact with a substrate, the spreading exhibits much larger exponents, with 
$n=1/2$ for low-viscosity Newtonian drops (Biance, Clanet & Quéré Reference Biance, Clanet and Quéré2004). Interestingly, the addition of polymers has no measurable effect on this rapid spreading (Bouillant et al. Reference Bouillant, Dekker, Hack and Snoeijer2022; Yada et al. Reference Yada, Bazesefidpar, Tammisola, Amberg and Bagheri2023). Similarly, during drop impact, the initial spreading phase is nearly unaffected by the presence of polymers (Bartolo et al. Reference Bartolo, Boudaoud, Narcy and Bonn2007; Gorin et al. Reference Gorin, Di Mauro, Bonn and Kellay2022; Sen, Lohse & Jalaal Reference Sen, Lohse and Jalaal2022).
By contrast, the retraction of a droplet, e.g. after an impact, is strongly influenced by the presence of polymers (Bergeron et al. Reference Bergeron, Bonn, Martin and Vovelle2000). The typical retraction velocity was found to decrease with the (first) normal stress coefficient 
$\psi$, following a scaling 
$U \sim 1/\sqrt {\psi }$ (Bartolo et al. Reference Bartolo, Boudaoud, Narcy and Bonn2007). This result can be interpreted from the balance between the normal stress near the contact line 
$\sim \psi (U/\ell )^2$ and the Laplace pressure 
$\gamma h'' \sim \gamma \theta _e^2/\ell$, giving a typical velocity
\begin{equation} U_\psi = \theta_e \sqrt{\frac{\gamma \ell}{\psi}}, \end{equation}
where 
$\ell$ is to be interpreted as a microscopic cutoff length (Bartolo et al. Reference Bartolo, Boudaoud, Narcy and Bonn2007), very much like 
$\ell _i$ in (1.1). The best fit of experimental data gave 
$\ell$ to be of the order of tens of microns (Bartolo et al. Reference Bartolo, Boudaoud, Narcy and Bonn2007), which is larger than the expected nanometric length used in Newtonian fluids. Alternatively, this result can be obtained from a balance of the capillary force 
$\gamma ( \cos {\theta } - \cos {\theta _e} ) \sim \gamma \theta _e^2/2$ and viscoelastic force 
$\psi U^2/ \ell$ induced by the extra tension along the streamlines (Bartolo et al. Reference Bartolo, Boudaoud, Narcy and Bonn2007). The addition of polymers also changes the drop shedding behind a sliding droplet (Xu et al. Reference Xu, Clarke, Rothstein and Poole2018).
In fact, polymers have been observed to be strongly stretched at a contact line (during retraction) (Smith & Bertola Reference Smith and Bertola2010) and stretched and depleted at the contact line (during spreading) (Shin, Bong & Kim Reference Shin, Bong and Kim2016) and are sometimes left on the surface after the retraction (Zang et al. Reference Zang, Wang, Geng, Zhang and Chen2013). This could be induced by either non-hydrodynamic polymer–surface interactions (Guyard et al. Reference Guyard, Vilquin, Sanson, Jouenne, Restagno and McGraw2021), or shear-induced migration (Han & Kim Reference Han and Kim2013; Shin & Kim Reference Shin and Kim2015). However, it has remained a challenge to disentangle the possible roles of all these mechanisms, which in part is due to the lack of a hydrodynamic framework for viscoelastic wetting.
The aim of the present study is to explore systematically the influence of the viscoelastic normal stress on the motion of contact lines. There have been several studies dedicated to quantifying normal stress effects in wetting flows. Numerical schemes have been developed using the phase-field method and polymer rheological models like Oldroyd-B or Giesekus (Yue & Feng Reference Yue and Feng2012; Wang, Minh & Amberg Reference Wang, Minh and Amberg2015; Wang, Do-Quang & Amberg Reference Wang, Do-Quang and Amberg2017), and results using these were found to be consistent with prior experimental observations (Wei et al. Reference Wei, Rame, Walker and Garoff2009b; Han & Kim Reference Han and Kim2013) in that the viscoelastic stresses, confined to a small region near the contact line, had only a mild effect on the contact-line motion of spreading drops. Boudaoud (Reference Boudaoud2007) introduced a phenomenological thin-film model that includes normal stress effects. Solving the resulting equation in the specific case of complete wetting, the advancing contact angle is found to follow
\begin{equation} \theta_{{app},o}^3 = 9 {Ca} \ln{\left( \frac{\ell_o}{\ell_{VE}}\right)}, \quad \ell_{VE} = \frac{\psi U}{2 \eta}. \end{equation}
This expression can be seen as a form of the Cox–Voinov law (1.1), with 
$\theta _e=0$ and the microscopic regularising length provided by a viscoelastic length 
$\ell _{VE}$. The latter corresponds to the product between the contact-line velocity 
$U$ and the typical relaxation time of the liquid 
$\psi /\eta$. However, it has remained unclear whether this theory can actually be derived from any tensorial constitutive relation, and how these results carry over to partial wetting conditions, for which the equilibrium contact angle 
$\theta _e$ must play a role. Also, as the equation derived by Boudaoud (Reference Boudaoud2007) (given as (1.3)) does not admit any solutions with partial wetting conditions, it cannot be applied to receding contact lines – which is crucial to understanding the observed strong effects of normal stresses on the retraction velocity (Bartolo et al. Reference Bartolo, Boudaoud, Narcy and Bonn2007).
In this paper, we study how normal stress affects both advancing and receding contact-line motion of partially wetting fluids. Through our analysis, we demonstrate the strong effect of the viscoelastic normal stress on receding contact lines which was previously observed in the experiments of Bartolo et al. (Reference Bartolo, Boudaoud, Narcy and Bonn2007). Our paper is based on an expansion of the stress to second order in the shear rate, an approximation which exhibits a normal stress effect. The so-called second-order fluid is appropriate for modelling slow and steady flows in the weak viscoelasticity limit (Tanner Reference Tanner2000; Morozov & Spagnolie Reference Morozov and Spagnolie2015; De Corato, Greco & Maffettone Reference De Corato, Greco and Maffettone2016). This constitutive relation was recently given a long-wave expansion (Datt, Kansal & Snoeijer Reference Datt, Kansal and Snoeijer2022), yielding a thin-film equation that can be used to describe contact-line motion in the presence of normal stress. The approach resembles somewhat the analysis by Han & Kim (Reference Han and Kim2014), who used corner solutions of the second-order fluid to estimate how normal stress affects the pressure inside the liquid. However, the thin-film description by Datt et al. (Reference Datt, Kansal and Snoeijer2022) results from a systematic long-wave expansion (naturally involving the gradient of pressure, rather than the pressure), and has the additional benefit of being applicable beyond the no-slip condition to relieve the moving contact-line singularity. We focus on the case of partial wetting using a Navier-slip condition as a microscopic regularisation, so that the slip length 
$\lambda _s$ appears as the inner scale. The central result of this paper consists of the modified form of the Cox–Voinov law
\begin{equation} \theta_{{app},o}^3 = \theta_e^3 - \frac{3}{4} \frac{\psi U^2\theta_e}{\gamma \lambda_s} \pm 9 {Ca} \ln{\left( \frac{\ell_o}{\lambda_s}\right)}, \end{equation}
which now exhibits an explicit dependence on the normal stress coefficient 
$\psi$. A dependence of the Cox–Voinov relation on 
$\sim U^2$ was also suggested by Kim & Rothstein (Reference Kim and Rothstein2015), but in that study the sign was proposed to be different for advancing (
$+$) and receding (
$-$) motions. Here, we find that the normal stress always leads to a smaller macroscopic contact angle. Comparing the first two terms on the right-hand side of (1.4), we identify a dimensionless normal stress parameter
\begin{equation} N = \frac{\psi U^2}{\gamma \lambda_s \theta_e^2}. \end{equation}
This dimensionless number 
$N \sim (U/U_\psi )^2$ can be seen as the viscoelastic equivalent of the capillary number, comparing the contact-line speed with the intrinsic velocity 
$U_\psi$ defined in (1.2). Indeed, as we will show, the derivation of (1.4) requires 
$N \ll 1$ in a way similar to 
${Ca} \ll 1$. The physical picture emerging from the viscoelastic Cox–Voinov theory is sketched in figure 1(c). The modification with respect to the Newtonian case appears only at small scales, where the bending of the interface due to normal stress dominates over the viscous bending. The form of the contact angle law (1.4) can then be interpreted in terms of an apparent microscopic angle 
$\theta _{{app},i}$, as shown in figure 1(c), of the form
\begin{equation} \theta_{{app},i}^3 = \theta_e^3 (1 - \tfrac{3}{4} N ). \end{equation}
This angle decreases with viscoelasticity, both for advancing and receding contact lines, as the viscoelastic normal stress parameter depends on 
$U^2$. We recall that (1.4) is derived under the assumption of weak viscoelasticity, i.e. 
$N \ll 1$. Importantly, however, when 
$N \gtrsim 1$, another regime of strong viscoelasticity can emerge. As a result, the complete wetting solution of (1.3) (Boudaoud Reference Boudaoud2007) is recovered, even when the substrate is partially wetting. Section 2 provides the detailed derivations of these results. Subsequently, § 3 explores the consequences of viscoelasticity for drop spreading, drop retraction and dip coating. Finally, the paper closes with a discussion in § 4.
2. Modified Cox–Voinov law
2.1. Second-order fluid thin-film equation
We investigate a moving three-phase contact line, by considering a liquid film on a rigid surface. We choose a reference frame in which the contact line is at a fixed position 
$x=0$ and the solid surface is moving at a speed 
$U$. We treat the problem with the lubrication approximation, which assumes relatively small angles. To include the effect of normal stress differences, we use the simplest viscoelastic fluid model containing normal stresses, i.e. the second-order fluid model (Tanner Reference Tanner2000). This model offers a good description of the polymer solution rheology for steady flows in the weakly viscoelastic limit.
The lubrication equation for steady contact-line motion with Navier-slip boundary condition has been derived recently for the second-order fluid (Datt et al. Reference Datt, Kansal and Snoeijer2022). Denoting the interface profile as 
$h(x)$ (see figure 1), the lubrication equation reads
\begin{equation} \gamma h''' ={-}\frac{3 \eta U}{h (h +3 \lambda_s )} - \frac{3}{4} \left( \psi \left(\frac{U}{h +3 \lambda_s } \right)^2 \right)' , \end{equation}
where 
$\lambda _s$ is the Navier-slip length. The introduction of a molecular mechanism, here the slip length, is necessary to regularise the well-known moving contact-line singularity (Huh & Scriven Reference Huh and Scriven1971; Bonn et al. Reference Bonn, Eggers, Indekeu, Meunier and Rolley2009; Snoeijer & Andreotti Reference Snoeijer and Andreotti2013). Setting 
$\psi =0$, we recover the standard Newtonian lubrication equation (Duffy & Wilson Reference Duffy and Wilson1997; Oron, Davis & Bankoff Reference Oron, Davis and Bankoff1997). The extra viscoelastic term in (2.1) is of the form of the gradient of normal stress, which is quadratic in shear rate. A very similar lubrication equation was obtained previously in Boudaoud (Reference Boudaoud2007), using scaling estimates for the components of the stress tensor. The resulting equation differs from (2.1) in two ways: slip is ignored and the prefactor in front of the viscoelastic term is not correct. We emphasise here that the second-order fluid thin-film equation is asymptotically equivalent to the Oldroyd-B model for weakly viscoelastic steady flows (see the Appendix in Datt et al. Reference Datt, Kansal and Snoeijer2022), ensuring the validity of (2.1). Importantly, the viscoelastic term scales as 
$1/h^3$ whereas the viscous term scales as 
$1/h^2$. Hence, we anticipate the normal stress to be negligible in the large 
$h$ limit with respect to the viscous term. The viscoelastic term is subdominant in the region far from the contact line, where viscous shear and capillary forces balance (see region (
$c_{ii}$) in figure 1). A consequence is that the structure of the Cox–Voinov law still holds with normal stress, and normal stress effects are important only at small scales.
We thus focus on the inner region near the contact line, in which the typical film thickness scale is 
$\lambda _s$. We therefore rescale 
$h$ and 
$x$ as
\begin{equation} X = \frac{x \theta_e}{3 \lambda_s }, \quad H(X) = \frac{h}{3 \lambda_s}, \quad \lambda = \frac{\lambda_s}{\theta_e}, \quad \delta = \frac{3 \eta U}{\gamma \theta_e^3} = \frac{3 {Ca}}{\theta_e^3}, \end{equation}
where the lateral scale has been rescaled by 
$3\lambda _s/\theta _e$ to enforce slopes of order 1 and to scale out the equilibrium angle. We also introduced a rescaled slip length 
$\lambda$ and capillary number 
$\delta$ following Eggers (Reference Eggers2005). The lubrication equation then becomes
\begin{equation} H''' ={-}\delta \left( \frac{1}{H^2 +H} - \frac{\ell_{VE}}{3 \lambda} \frac{ H'}{( H+1)^3} \right), \end{equation}
where viscoelastic effects are quantified by the dimensionless ratio 
$\ell _{VE}/\lambda$.
Equation (2.3) holds both for advancing (
$\delta >0$) and receding (
$\delta <0$) contact lines. In this section, we focus on advancing contact lines, for which the boundary conditions are (Eggers Reference Eggers2005)
\begin{equation} H(0)=0, \quad H'(0)=1, \quad H''(\infty) =0, \end{equation}
fixing the position of the contact line and the equilibrium angle condition, respectively. The third condition 
$H''(\infty )=0$, which can only be imposed for advancing contact lines (Duffy & Wilson Reference Duffy and Wilson1997), matches any outer solution that has an asymptotically small curvature toward the contact line.
We perform direct numerical integration of the system (2.3) and (2.4a–c), by using the continuation code AUTO-07P (Doedel et al. Reference Doedel, Champneys, Dercole, Fairgrieve, Kuznetsov, Oldeman, Paffenroth, Sandstede, Wang and Zhang2007). Typical profiles 
$H(X)$ are shown in figure 2(a) for 
$\delta = 10^{-2}$ and changing the ratio 
$\ell _{VE}/\lambda$ to vary the strength of normal stresses. The normal stress effects tend to lower the slope of the interface profile, although each of the profiles exhibits the same microscopic contact angle 
$H'(0)=1$ (see figure 2b). Asymptotically, in the large-
$X$ limit, the interface slopes cubed all follow the same behaviour, as
\begin{equation} H^{\prime 3} \simeq a + 3 \delta \ln{( e X )}, \end{equation}
and differ by a vertical offset position denoted 
$a$, which depends on 
$\ell _{VE}/\lambda$, as shown in figure 2(c). This asymptotic behaviour relates to the Voinov solution, which is classically found in moving contact-line problems, reflecting the visco-capillary balance (see region (
$c_{ii}$) in figure 1). We introduce in (2.5) a factor 
$e= \exp ({1})$ inside the logarithm for convenience, so that the Newtonian solution corresponds to 
$a=1$ (as shown below). Hence, as anticipated, the normal stress effect does not modify the large-scale asymptotic behaviour of the inner solution, but induces a non-trivial offset in the Voinov solution. In the next two subsections, we determine the offset 
$a$ by using asymptotic expansions, in the weak and strong viscoelastic limits.
Advancing contact line. (a) Interface profile 
$H(X)$ obtained from numerical integration of (2.3) and (2.4a–c), with reduced capillary number 
$\delta = 10^{-2}$. The inset shows a zoom near the contact line, illustrating that all profiles have the same microscopic contact angle. (b) Cube of the interface slope 
$H'(X)$ as function of distance from the contact line. Different colours correspond to varying normal stress effect quantified by 
$\ell _{VE}/\lambda$. (c) The parameter 
$a$, obtained from fitting (2.5) to the large-
$X$ limit of 
$H'(X)$, is plotted as a function of 
$\ell _{VE}/\lambda$. The dots indicate the cases of (a) and (b) using the same colour code. For 
$a>0$, the parameter 
$a=(\theta _{{app},i}/\theta _e)^3$ can be interpreted in terms of the apparent microscopic contact angle, and is well described by (2.10) (red dashed line). A regime of ‘apparent complete wetting’ emerges when 
$a<0$, captured by (2.21) (green dashed line).
2.2. Weak viscoelasticity: apparent microscopic angle
The standard approach to determining the interface profile near the contact line analytically is via a perturbation expansion of (2.3) in the small parameter 
$\delta$ (Hocking Reference Hocking1983; Eggers Reference Eggers2005)
\begin{equation} H(X) = H_0(X) + \delta H_1(X) + O(\delta^2). \end{equation}
Here, supposing weak viscoelasticity, we assume explicitly the ratio 
$\ell _{VE}/\lambda$ to be of order 
$\sim 1$. The zeroth-order solution satisfies 
$H_0'''=0$, which from the boundary conditions (2.4a–c) gives 
$H_0 = X$. At the first order in 
$\delta$, the interface profile follows
\begin{equation} H_1''' ={-} \frac{1}{X^2 + X} + \frac{\ell_{VE}}{3 \lambda} \frac{1}{( 1 + X )^3}, \end{equation}
with boundary conditions 
$H_1'(0)=0$ and 
$H_1''(\infty )=0$. Integration of (2.7) then gives
\begin{equation} H_1' ={-} X \ln{X} + (1+X) \ln{(1+X)} + \frac{\ell_{VE}}{6 \lambda} \left( \frac{1}{1+X} - 1\right) \simeq \ln{X} + 1 - \frac{\ell_{VE}}{6 \lambda} , \end{equation}
where the last step represents the 
$X \to \infty$ asymptotic behaviour. Combined with the zeroth-order solution, the large-
$X$ asymptotic behaviour of the interface slope in the small-
$\delta$ limit is
\begin{equation} H' \simeq 1 + \delta \left( \ln{X} + 1 - \frac{\ell_{VE}}{6 \lambda} \right) + O(\delta^2). \end{equation}Raising this result to the third power and comparing with (2.5), we obtain the sought-after offset constant
\begin{equation} a = 1 - \frac{1}{2} \frac{\delta \ell_{VE}}{\lambda}. \end{equation}
Compared with the Newtonian case where 
$a=1$, the interface slope indeed decreases due to the presence of normal stress effects. The dimensional form of (2.5) reads
\begin{equation} h^{\prime 3} = (\theta_{{app},i} )^3 + 9 {Ca} \ln{\left( \frac{e x \theta_e}{3 \lambda_s} \right)}. \end{equation}
Here, we introduce the apparent microscopic angle 
$\theta _{{app},i}$ (see figure 1), that depends on the normal stress parameter 
$N$ as
\begin{equation} (\theta_{{app},i} )^3= a \theta_e^3 = \theta_e^3 \left( 1 - \frac{3}{4} N \right), \end{equation}
where we used 
$\ell _{VE}/\lambda = 3 N/( 2 \delta )$. We identify (2.11) as a modified Cox–Voinov solution where the equilibrium Young contact angle 
$\theta _e$ has been modified to an apparent microscopic angle 
$\theta _{{app},i}$ (see (2.12)). This result corresponds to (1.4) in § 1.
The numerical data agree very well with this prediction for small-to-intermediate values of 
$\ell _{VE}/\lambda$ (see red dashed lines in figure 2c). However, a significant deviation is observed at large 
$\ell _{VE}/\lambda$, which coincides with 
$a < 0$ (grey zone in figure 2c). To understand the breakdown of (2.10), we recall that the underlying perturbation expansion assumed both 
$H'$ and 
$\ell _{VE}/\lambda$ to remain of order unity, which is clearly violated at large 
$\ell _{VE}/\lambda$. Instead, one observes the formation of a flat film region where 
$H' \simeq 0$ in the large viscoelastic limit (see figure 2a,b), indicating the onset of a new regime. The transition between the two regimes occurs when the apparent microscopic angle goes to zero 
$\theta _{{app},i} \to 0$, which occurs when
\begin{equation} \theta_e^3 \sim \frac{\psi U^2\theta_e}{\gamma \lambda_s} \implies \frac{\ell_{VE}}{\lambda} \sim \frac{1}{\delta}. \end{equation}
This estimate is consistent with figure 2(c), where 
$\delta = 10^{-2}$, for which the perturbation expansion (2.10) breaks down around 
$\ell _{VE}/\lambda \sim 1/\delta \sim 10^2$. The next section is dedicated to a description of the new regime at large 
$\ell _{VE}/\lambda$, named strong viscoelasticity.
2.3. Strong viscoelasticity: (apparent) complete wetting
We interpret the behaviour at large 
$\ell _{VE}/\lambda$ as a regime of ‘apparent complete wetting’. Indeed, the purple and brown solutions in figure 2(b) exhibit an extended flat zone with 
$H' \approx 0$, resembling a situation of the complete wetting condition, even though 
$H'(0)=1$. In the flat film region, the film curvature is very small, such that the capillary term is subdominant. Hence, the governing balance of (2.1) involves the viscous and viscoelastic terms, from which one identifies 
$\ell _{VE}$ as the natural lateral length scale. Physically, this implies that the information of the partial wetting nature of the problem, encoded in 
$\theta _e$, is lost across the flat film as it only modifies the profiles at a length scale 
$\sim \lambda$ (see figure 2b), but this is ‘screened’ by the flat wetting film on a length scale 
$\ell _{VE}$.
We proceed in the asymptotic analysis using 
$\ell _{VE}$ as the relevant length on the flat film region. We assume, and validate a posteriori, that we can drop the slip length 
$\lambda _s$ from (2.1), in the large 
$\ell _{VE}/\lambda$ limit. Hence, we introduce the rescaling
\begin{equation} \xi = \frac{x}{\ell_{VE}}, \quad \mathcal{H} = \frac{h}{\ell_{VE} (3 {Ca})^{1/3}}, \end{equation}for which (2.1) becomes
\begin{equation} \mathcal{H}''' ={-}\frac{1}{\mathcal{H}^2} + \frac{\mathcal{H}'}{\mathcal{H}^3}. \end{equation}
Here, we dropped the slip terms, which is consistent as long as 
$\ell _{VE}/\lambda \gg \delta ^{-1/3}$, which is indeed satisfied in the strong viscoelastic limit. We notice that (2.15) resembles the thin-film equation with a disjoining pressure (see p. 398 in Eggers & Fontelos Reference Eggers and Fontelos2015), which is used to model contact lines in complete wetting. Within this analogy, the effective ‘disjoining pressure’ is 
$-3\psi U^2 / (4h^2)$, similar to the standard van der Waals disjoining pressure 
$-A/(6{\rm \pi} h^3)$, where 
$A$ is the Hamaker constant. We employ boundary conditions
\begin{equation} \mathcal{H}(-\infty)=0, \quad \mathcal{H}'(-\infty)=0, \quad \mathcal{H}''(+\infty) =0, \end{equation}
where the first two conditions correspond to the formation of the film, replacing the partial wetting condition at 
$\xi =0$. Similarly to § 2.2, since viscoelasticity only acts at small scales, the large distance asymptotic behaviour is given by the Voinov solution (see (2.5)), which using the new scales reads
\begin{equation} \mathcal{H}'^3 \simeq 3 \ln \left(\frac{\xi}{\xi_0} \right). \end{equation}
The offset 
$a$ has been replaced by an unknown constant 
$\xi _0$, which remains to be found.
The system (2.15) and (2.16a–c) is exactly the problem analysed in Boudaoud (Reference Boudaoud2007) up to a trivial rescaling, which we now understand to emerge as the limit of large 
$\ell _{VE}/\lambda$. In the flat film region, the curvature is nearly zero and the capillary term drops out, such that the dominant balance involves the viscous stress and normal stress. Towards the flat film, (2.15) reduces to
\begin{equation} \mathcal{H}' \simeq \mathcal{H}, \end{equation}
to leading order. The viscous-viscoelastic solution takes the form of an exponentially growing solution 
$\mathcal {H} \simeq e^{\xi }$, with a prefactor that can be set to unity, owing to the translational invariance of (2.15) and (2.16a–c). A more refined analysis of the asymptotic behaviour of 
$\mathcal {H}$ in the 
$\xi \to -\infty$ limit leads to the asymptotic expression (Boudaoud Reference Boudaoud2007)
\begin{equation} \mathcal{H} \simeq {\rm e}^{\xi} + c_{+} \,{\rm e}^{ \frac{2}{3} \,{\rm e}^{{-}3\xi/2}} + c_{-} \,{\rm e}^{ - \frac{2}{3} \,{\rm e}^{{-}3\xi/2}}, \end{equation}
where 
$c_\pm$ are constants; using the boundary condition for 
$\xi \to -\infty$ leads to 
$c_{+} =0$. We solve (2.15) numerically, using (2.19) as an initial condition, where 
$c_-$ is treated as a shooting parameter to ensure 
$\mathcal {H}''(\xi \to \infty ) =0$. Fitting the numerically obtained slope 
$\mathcal {H}'(\xi )$ in the large-
$\xi$ limit, typically 
$\xi = 10^8$, with the Voinov solution (Bender & Orszag Reference Bender and Orszag1999; Eggers & Fontelos Reference Eggers and Fontelos2015), we find a universal value of 
$\xi _0 \approx 1.16$ (see Snoeijer & Eggers (Reference Snoeijer and Eggers2010) for details of the numerical calculation).
In dimensional units, the contact-line solution at large 
$x$ gives
\begin{equation} h'^3 = 9 {Ca} \ln{\left( \frac{x}{\ell_{VE} \xi_0}\right)}. \end{equation}
This has the same form as the Cox–Voinov solution for complete wetting, i.e. when 
$\theta _e=0$, where 
$\ell _{VE}$ now offers the regularisation length, replacing the microscopic scale. This was also demonstrated in Boudaoud (Reference Boudaoud2007), but we now provide the prefactor. In terms of the constant 
$a$, as defined by (2.5), we find
\begin{equation} a = 3 \delta \ln \left( \frac{3 \lambda}{e \ell_{VE} \xi_0} \right), \end{equation}
which agrees well with the numerical solution at large 
$\ell _{VE}/\lambda$ (figure 2(c), green dashed line). It justifies a posteriori that the equilibrium angle 
$\theta _e$ does not play a role in this regime.
3. Applications
We now apply the viscoelastic contact-line theory from the previous section to two paradigmatic wetting problems, namely drop spreading and dip coating. Classically, the interface profile in these problems can be described using matched asymptotics, where the Voinov contact-line solution corresponds to the inner solution, matched to an outer solution at macroscopic scale. In this section, we briefly recall the outer solutions of the two problems, which are unaltered by viscoelasticity, and discuss the matching to the newly derived viscoelastic inner solutions.
In applications such as drop spreading, the contact-line velocity 
$U$ is not a control parameter. Instead, the contact-line velocity is determined through the matching procedure. Using the Cox–Voinov framework, an intrinsic velocity scale 
$U^*$ naturally arises by comparing the viscous and capillary terms as
\begin{equation} U^* = \frac{\gamma \theta_e^3}{\eta}, \end{equation}which we will use as a rescaling velocity below. In these scaled units, we will see the appearance of a new dimensionless parameter
\begin{equation} N_0 = \left( \frac{U^*}{U_\psi} \right)^2 = \frac{\psi \gamma \theta_e^4}{\eta^2 \lambda_s }. \end{equation}
This parameter governs the relative importance of viscosity and normal stress in the context of partially wetting flows. Alternatively, 
$N_0$ can be interpreted as a scaled normal stress coefficient, made dimensionless only with material parameters.
3.1. Drop spreading and retraction
We study the spreading (and retraction) of a drop, of volume 
$\varOmega$, over a solid surface (figure 3). We consider drops smaller than the capillary length so that gravity can be ignored. The outer solution of the interface profile is obtained away from the contact line, where viscous and viscoelastic forces are not dominant. Therefore, at a given time, the shape adopted by the drop is close to the equilibrium shape of a static drop 
$h_0(r)$ as
\begin{equation} h_0(r)= \frac{2 \varOmega}{{\rm \pi} r^2} \left[ 1- \left( \frac{r}{R}\right)^2 \right], \end{equation}
where the drop meets the solid at the apparent angle 
$\theta _{{app},o} = -h_0'(r=R)$, which reads
\begin{equation} \theta_{{app},o}= \frac{4 \varOmega}{ {\rm \pi}R^3}. \end{equation}
As discussed before, close to the contact line (to be precise, the outer asymptotics of the inner region), the slope of an advancing contact-line solution exhibits a logarithmic dependence with the distance. Hence, the static equilibrium profile (3.3) cannot be matched directly to the Cox–Voinov theory. To resolve this issue, Hocking (Reference Hocking1983) showed that the first-order correction to the outer solution in capillary number, defined as 
$\eta \dot {R} /\gamma$, where 
$\dot {R}$ is the radius rate of change, takes the asymptotic form
\begin{equation} h'^3(x)= \theta_{{app},o}^3 + 9 \frac{\eta \dot{R}}{\gamma}\ln{\left( \frac{x}{R b}\right)}, \end{equation}
for 
$x=R-r \ll R$ (see figure 3a), and where 
$b=1/(2e^2 )$. This solution is to be matched to the newly derived viscoelastic inner solutions (2.11) and (2.20).
Drop spreading (a) and retraction (b,c). The schematic indicates the case where the initial drop is drawn in black and the equilibrium shape in grey. The drop contact radius rescaled by the equilibrium contact radius is plotted as a function of dimensionless time for various 
$N_0$, corresponding to numerical solutions of (3.9), where 
$c_1 = 10^7$. The black dashed lines show the Newtonian solution for 
$N_0=0$, and the grey dotted line displays the viscoelastic retraction velocity, corresponding to (3.12) evaluated at the initial time. (c) The initial retraction velocity scaled by 
$U^*$ as a function of 
$N_0$, computed from (3.8) using 
$R_{{initial}}/R_e=1.5$ (as in panel (b)). The orange line in large 
$N_0$ regime corresponds to the 
$U \sim \psi ^{-1/2}$ scaling as expected from (3.12). (a) Spreading: contact radius. (b) Retraction: contact radius. (c) Initial retraction velocity.
3.1.1. Spreading with weak viscoelasticity: partial wetting
For the partial wetting case, the drop spreads until the contact angle reaches the equilibrium contact angle 
$\theta _e$, or equivalently an equilibrium radius 
$R_{e}= ( 4 \varOmega / ({\rm \pi} \theta _e ) )^{1/3}$. For small capillary number, we can match the logarithmic dependencies of the outer solution (3.5) to the modified Cox–Voinov solution (2.11), which gives
\begin{equation} (\theta_{{app},o} )^3 = (\theta_{{app},i} )^3 + 9 \frac{\eta \dot{R}}{\gamma} \ln{\left( \frac{R \theta_{e}}{6 e \lambda_s} \right)}. \end{equation}
As usual, the matching procedure connects the microscopic angle to the macroscopic angle, but now the microscopic angle contains viscoelastic effects. Inserting 
$\theta _{{app},o}$ from the static solution (3.4) and 
$\theta _{{app},i}$ from (2.12), we obtain an ordinary differential equation (ODE) for the drop radius 
$R(t)$, as
\begin{equation} (\theta_{{app},o} )^3 = \left(\frac{4 \varOmega}{ {\rm \pi}R^3} \right)^3 = \theta_e^3 - \frac{3}{4}\frac{\psi \dot{R}^2 \theta_e}{\gamma \lambda_s} + 9 \frac{\eta \dot{R}}{\gamma} \ln{\left( \frac{R \theta_{e}}{6 e \lambda_s} \right)}. \end{equation}
Finally, using dimensionless variables with bar as 
$\bar {R} = R/R_e$ and 
$\bar {t} = t / (R_e/U^*)$, we can rewrite (3.7) as
\begin{equation} \frac{3}{4} N_0 \dot{\bar{R}}^2 - 9 \ln{( c_1\bar{R} )} \dot{\bar{R}} + \left( \frac{1}{\bar{R}^9} - 1 \right) = 0, \end{equation}
with 
$c_1 = (R_e \theta _e)/(6 e \lambda _s)$. Here, we recover the explicit dependency on 
$N_0$, which is the material parameter that, for a given fluid, quantifies the normal stress effect.
Equation (3.8) is a second-order polynomial equation for 
$\dot {\bar {R}}$, which we first solve to isolate 
$\dot {\bar {R}}$, leading to
\begin{equation} \dot{\bar{R}} = \frac{6\ln( c_1 \bar{R} )}{N_0} \left[1 - \sqrt{1 - \frac{N_0}{27 \ln^2\left( c_1 \bar{R} \right)}\left(\frac{1}{\bar{R}^9}-1\right)}\,\right], \end{equation}
and then integrate numerically. We notice that the second-order polynomial equation does not have a solution in real space for negative discriminant, which imposes a condition 
$27\ln ^2( c_1\bar {R} ) > N_0 (1/\bar {R}^9-1)$. Hence, in the spreading case, i.e. when 
$\bar {R}(\bar {t}=0)<1$, we cannot find solutions for arbitrarily large values of 
$N_0$. Mathematically speaking, the reason is that the apparent macroscopic angle vs 
$U$ has a maximum value when using the modified Cox–Voinov law. Therefore, if the outer geometry imposes an apparent macroscopic angle that is larger than this maximum, no contact-line speed can be found. The lack of solution is an artefact of the weak viscoelastic solution, which as also mentioned in the Discussion, is resolved by admitting the strong viscoelastic solution (2.20).
In figure 3(a) we show a typical solution of (3.9), with 
$\bar {R}(\bar {t}=0) = 0.95$ as an initial condition. Qualitatively, the normal stress effects speed up the spreading dynamics, although the quantitative change is very weak. This corroborates the conclusions of drop impact with polymer solutions, showing a weak dependence on the addition of polymers (Bartolo et al. Reference Bartolo, Boudaoud, Narcy and Bonn2007; Gorin et al. Reference Gorin, Di Mauro, Bonn and Kellay2022; Sen et al. Reference Sen, Lohse and Jalaal2022).
3.1.2. Spreading with strong viscoelasticity: (apparent) complete wetting
As shown in § 2.3, the advancing contact-line dynamics at large velocities deviates from the modified Cox–Voinov law used in the previous subsection. In the drop spreading problem, if the initial radius is small, or equivalently the apparent angle is much larger than the equilibrium angle, then the contact-line speed is large and lies within the strongly viscoelastic regime. In this case, the inner solution for partially wetting substrates becomes completely independent of the equilibrium angle, and takes the form (2.20). Repeating the matching procedure with this strongly viscoelastic inner solution, the drop radius follows from
\begin{equation} \theta_{{app},o}^3 = 9 \frac{\eta \dot{R}}{\gamma} \ln{\left( \frac{R b}{\ell_{VE} \xi_0}\right)} = 9 \frac{\eta \dot{R}}{\gamma} \ln{\left( \frac{2\eta R b}{\psi \dot{R} \xi_0}\right)}. \end{equation}This result is identical to that of a completely wetting viscoelastic drop; and it is of the same form as derived in Boudaoud (Reference Boudaoud2007).
During the spreading, which we initially assumed to be rapid, the contact-line velocity will decrease in time. For partially wetting drops, this means that, upon approaching the equilibrium angle, the strongly viscoelastic regime will give way to the weak regime discussed in the preceding subsection.
For completely wetting surfaces, however, the drop radius follows a scaling law 
$R\propto t^{1/10}$ with logarithmic corrections, which is analogous to the Newtonian dynamics known as Tanner's law (Bonn et al. Reference Bonn, Eggers, Indekeu, Meunier and Rolley2009; Eggers & Fontelos Reference Eggers and Fontelos2015). Using 
$R \propto t^{1/10}$ as an ansatz, the ratio between contact-line speed and drop radius inside the logarithm in (3.10) can be approximated as 
$\dot {R}/R \approx 1/(10t)$. Inserting (3.4) into (3.10), and integrating, we obtain a modified Tanner's law as
\begin{equation} R \approx \left( \frac{10}{9} \left( \frac{4 \varOmega}{ {\rm \pi}} \right)^3 \frac{\gamma}{\eta \ln{\left( c_2\dfrac{\eta t}{\psi}\right)}} t \right)^{1/10}, \end{equation}
where 
$c_2= 10/(e^2\xi _0)$. We note that the same equation has been derived using scaling arguments in Rafaï et al. (Reference Rafaï, Bonn and Boudaoud2004) and Boudaoud (Reference Boudaoud2007), but we now provide a prefactor as derived for the second-order fluid. Again, the spreading of drops is accelerated by normal stress effects with respect to the Newtonian dynamics via the logarithmic term in (3.11).
3.1.3. Drop retraction
If the initial drop radius is above the equilibrium radius, the drop is retracting (see figure 3b), for instance after an impact on a surface. Even though the contact line is now receding, the same framework as § 3.1.1 can be used as long as the change in contact angle is sufficiently small (cf. Eggers (Reference Eggers2005) and the discussion on dip coating). By consequence, (3.8) also governs the retraction dynamics of a partially wetting drop, but 
$\dot {\bar R}$ is now negative. In the retraction case, the discriminant of the second-order polynomial equation is always positive such that one can find a solution for any 
$N_0$. Interestingly, viscoelastic effects now oppose the motion and slow down the retraction (see figure 3b). In the large viscoelastic limit, for 
$N_0 \to \infty$, the viscous stress becomes negligible and the retraction velocity is set by balancing the normal stress and the capillary term of (3.7), which gives in dimensional form
\begin{equation} \dot{R}^2 = \frac{4 \gamma (\theta_e^3-\theta_{{app},o}^3)\lambda_s}{3\psi \theta_e}. \end{equation}
The dotted line in figure 3(b) indicates the retraction as estimated using the viscoelastic receding speed (3.12), with the apparent angle taken from the initial condition. This provides a very good description of the retraction velocity for large 
$N_0$ (the contact line only slowing down upon approaching the equilibrium radius). The initial retraction velocity is plotted as a function of 
$N_0$ in figure 3(c). In the small viscoelastic limit, quantified by small 
$N_0$, the viscoelastic term of (3.8) is subdominant and the retraction velocity typically scales with 
$U^*$. Conversely, for large viscoelasticity, the retraction velocity is set by 
$U_\psi \sim 1/\sqrt {\psi }$, leading to the scaling 
$U/U^* \sim N_0^{-1/2}$, observed at large 
$N_0$ in figure 3(c). Using the values of relevant parameters from the experiments by Bartolo et al. (Reference Bartolo, Boudaoud, Narcy and Bonn2007), we estimate 
$N_0$ to be in the range 
$\sim 10^3$ to 
$10^6$. The 
$N_0$ values chosen in figure 3(a,b) fall within this range.
In the large viscoelastic limit, this retraction velocity scales with the inverse square root normal stress coefficient, as observed experimentally in Bartolo et al. (Reference Bartolo, Boudaoud, Narcy and Bonn2007). These authors rationalised their experimental findings by using a phenomenological normal stress thin-film model and derived a retraction velocity of similar form as the one discussed here, namely 
$\dot {R}^2 = \gamma (\theta _e^2 - \theta _{{app},o}^2) \ell /(8 \psi )$, where 
$\ell$ is a microscopic cutoff length. Up to a prefactor, in the limit of 
$\theta _{{app},o} \to 0$, this phenomenological equation agrees with the rigorous calculation provided here.
3.2. Forced wetting transition
We now turn to dip coating, sketched in figure 4(a), in which a plate is pulled out of a liquid bath at a constant speed 
$U$. At small speed, the liquid–air interface is modified and the position of the contact line is shifted upward with respect to the equilibrium position. Above a critical speed 
$U_c$, however, there is no static contact-line solution anymore and a Landau–Levich liquid film is entrained. This transition is usually called forced dewetting (Eggers Reference Eggers2005; Snoeijer et al. Reference Snoeijer, Andreotti, Delon and Fermigier2007; Chan, Snoeijer & Eggers Reference Chan, Snoeijer and Eggers2012; Galvagno et al. Reference Galvagno, Tseluiko, Lopez and Thiele2014).
Forced wetting transition. (a) Schematic of the dip-coating problem. Below the critical speed 
$U_c$ (top), the contact line is displaced and the apparent contact angle is 
$\theta _{{app},o}$. Above the critical speed, a Landau–Levich film of thickness 
$h_\infty$ is entrained. (b) Critical dimensionless speed as function of dimensionless normal stress coefficient. The dots represent full-scale numerical integration of (3.13), with 
$\theta _p/\theta _e = 1$ and 
$\lambda _s/(\ell _\gamma \theta _e) = 10^{-4}$. The dashed line shows the asymptotic matching solution (3.17) based on the modified Cox–Voinov law.
The question we wish to address is how 
$U_c$ is affected by viscoelasticity. We restrict ourselves to the case of small plate angle 
$\theta _p$, so that lubrication theory can be described for both the contact line and the bath. Assuming a stationary state, the interface position can then be computed from the lubrication equation as
\begin{equation} \gamma h''' = \frac{3 \eta U}{h (h +3 \lambda_s )} - \frac{3}{4} \left(\psi\left(\frac{U}{h +3 \lambda_s } \right)^2\right)' + \rho g ( h' - \theta_p). \end{equation}
With respect to (2.1), the sign of the viscous term is flipped as the contact line recedes, and a gravitational term has been added, 
$\rho$ and 
$g$ being the liquid density and gravitational acceleration, respectively. The sign of the viscoelastic term is unchanged, as it is quadratic in contact-line speed. The boundary conditions are
\begin{equation} h(0)=0, \quad h'(0)=\theta_e, \quad h''(\infty) =0, \end{equation}
and are identical to those in § 2. In the Newtonian case, solutions only exist up to a critical velocity. Likewise, we have numerically determined 
$U_c$ for the viscoelastic case by finding the maximum velocity for which (3.13) and (3.14a–c) admit a solution. A typical numerical result is shown as the circles in figure 4(b), where we plot the normalised critical speed 
$U_c/U^*$ as a function of the viscoelastic parameter 
$N_0$. We thus observe that the normal stress effect lowers the critical speed of entrainment with respect to the Newtonian case (
$N_0=0$). As for the droplet retraction dynamics, any receding contact-line motion is inhibited by viscoelastic effects, now leading to an earlier entrainment.
Once again, these results can be described quantitatively from the modified inner solution for receding contact lines. The outer solution is again unaffected by viscoelasticity: far from the contact-line position, the interface profile is found by balancing the hydrostatic pressure with the Laplace pressure. This leads to the equation of a static meniscus, which is entirely characterised by an apparent angle 
$\theta _{{app},o}$ (see figure 4a) and which corresponds to the outer solution of the full-scale problem. The inner solution close to the contact line, at scales much smaller than the capillary length 
$\ell _\gamma = \sqrt {\gamma /(\rho g)}$, does not involve gravity. Importantly, and in contrast to the advancing case, the general solution of the inner problem does not have a vanishing curvature at infinity. Therefore, the matching procedure is more intricate and analytical progress can be made by invoking an intermediate zone, for 
$\lambda \ll x\ll \ell _\gamma$, where slip and gravity are neglected, such that (3.13) reduces to the visco-capillary balance 
$\gamma h''' = 3\eta U/h^2$. The latter equation has an exact solution (Duffy & Wilson Reference Duffy and Wilson1997), which has a Voinov asymptotic behaviour 
$h'^3 = -9 {Ca} \ln (x/\ell )$ at small 
$x$ and a constant curvature 
$h'' \to {cte}$ at large 
$x$. As shown in Eggers (Reference Eggers2005), this visco-capillary solution can be matched both to the outer static meniscus and to the slip region.
Hence, to understand the normal stress effects on dip coating, which appear only on the slip scale, the remaining task is to derive the slip-scale viscoelastic receding solution. Modifying the steps of § 2.2 along the lines of Eggers (Reference Eggers2005), but now including the normal stress term, we find that the large-
$x$ asymptotic behaviour of the viscoelastic inner solution follows
\begin{equation} h^{\prime 3} = (\theta_{{app},i} )^3 - 9 {Ca} \ln{\left( \frac{e x \theta_e}{3 \lambda} \right)}. \end{equation}
This result is identical to the advancing counterpart (2.11), with only a change of sign in the viscous term. Importantly, the only change with respect to the Newtonian case is a modification of the microscopic angle from 
$\theta _e$ to the apparent microscopic angle defined in (2.12). Following Eggers (Reference Eggers2005), we find the equation for the critical speed
\begin{equation} \theta_{{app},i}^3 = 9 \frac{\eta U_c}{\gamma} \ln\left( \frac{(\eta U_c/\gamma)^{1/3} \ell_\gamma \theta_e}{18^{1/3} {\rm \pi}[{Ai}(s_{max})]^2 \lambda_s \theta_p} \right), \end{equation}
where 
${Ai}(s_{max}) = 0.53566\ldots$ is the global maximum of the Airy function. The Newtonian result is recovered by replacing the left-hand side by 
$\theta _e^3$. Using the characteristic speed 
$U^*$ as a scale, the dimensionless critical velocity 
$\bar {U}_c = U_c/U^*$ is the solution of
\begin{equation} \bar{U}_c = \frac{1-\dfrac{3}{4} N_0 \bar{U}_c^2}{9 \ln (c_3 \bar{U}^{1/3}_c)}, \end{equation}
where 
$c_3 = 1/(18^{1/3} {\rm \pi}[{Ai}(s_{max})]^2)(\ell _\gamma \theta _e^2 )/( \lambda _s \theta _p ) \approx 0.423 (\ell _\gamma \theta _e^2 )/( \lambda _s \theta _p )$. Here again, we recover the material parameter 
$N_0$ that quantifies viscoelasticity, as discussed earlier.
Figure 4(b) shows excellent agreement between (3.17), shown as the dashed line, and the numerical solution. Hence, the decrease of the entrainment velocity with increasing viscoelasticity can be understood from the decrease of the apparent inner angle 
$\theta _{{app},i}$. In the limit of large 
$N_0$, one finds 
$\bar U_c^2 \simeq 4/(3N_0)$, which in dimensional units reads
\begin{equation} U^2_c \simeq \frac{4 \gamma \theta_e^2\lambda_s}{3\psi}. \end{equation}Note the similarity with (3.12) for drop retraction. This result suggests the emergence of a universal receding velocity for strong normal stress effect, independently of the outer geometry.
4. Discussion
In this paper we have studied the effect of viscoelastic normal stresses on contact-line motion. Using the lubrication equation for the second-order fluid, we have derived a modified form of the classical Cox–Voinov theory: the macroscopic apparent contact angle (
$\theta _{{app},o}$) is modified, due to a small-scale ‘apparent inner contact angle’ (
$\theta _{{app},i}$) that is lowered due to viscoelasticity both for advancing and receding contact-line motion. Figure 5(a) shows the force balance for a fluid volume in the inner region close to the contact line, with the two dominant forces: the capillary force and normal stress. The contribution of the normal stress on both lateral sides can be calculated by integrating the stress across the film thickness, and scales as 
$\sim \psi (U/h )^2 h= \psi U^2/h$. Hence, the left side near the contact line is dominant as 
$h \to 0$, and the global effect of normal stresses is to pull towards the contact line, regardless of whether the contact line is advancing or receding. This explains why the normal stress facilitates advancing and inhibits receding motions. The singularity of the viscoelastic force as 
$h \to 0$ is regularised by the slip length 
$\lambda _s$, such that the effective resulting viscoelastic force is 
$\psi U^2/\lambda _s$. This viscoelastic force in the inner region bends the interface and lowers the microscopic angle from 
$\theta _e$ to 
$\theta _{{app},i}$. We have worked out quantitatively the effect on drop spreading, drop retraction and dip coating, and found fundamental differences between advancing and receding contact lines. To summarise the findings, we turn to the schematics for the macroscopic apparent angle 
$\theta _{{app},o}$ vs speed provided in figure 5(b,c).
(a) Force balance on a control volume, highlighted by the blue area, in the inner region near the contact line, where capillary force 
$\gamma$ and the normal stress 
$\psi (U/h)^2$ are balanced. Macroscopic apparent contact angle vs speed for receding (b) and advancing (c) contact lines. In the advancing case, we need to distinguish between weak and strong viscoelastic behaviour. The cross-over between these regimes is indicated schematically by the black solid line, and occurs around 
$U_\psi$ defined in (1.2).
For receding contact lines, viscosity and normal stress provide forces that act in the same direction, both opposing the contact-line motion. Both viscosity and normal stress give rise to a ‘bending’ of the interface that tends to lower the macroscopic angle, as can be seen in figure 5(b). As a result, the retraction velocity of a droplet and the critical velocity for dip coating are lowered with respect to the Newtonian case. The relative importance of viscous and viscoelastic bending is governed by the parameter 
$N_0$ defined in (3.2). In the case of strong viscoelasticity, the second-order fluid model predicts a universal receding velocity (3.18), which is proportional to 
$U_\psi$ defined in (1.2), that gives both the critical speed for dip coating and droplet retraction. This prediction is consistent with experiments by Bartolo et al. (Reference Bartolo, Boudaoud, Narcy and Bonn2007), who in fact already introduced 
$U_\psi$ to explain experimental measurements of droplet retraction after impact. Specifically, they found the retraction velocity to be proportional to 
$1/\sqrt {\psi }$. However, they found that the microscopic length appearing in 
$U_\psi$ was tens of microns, which is much larger than a nanoscopic slip length. Bartolo et al. (Reference Bartolo, Boudaoud, Narcy and Bonn2007) remarked that the fitted lengths are comparable to the length of the fully extended polymer chains. Such large values of cutoff length may be due a shear-induced migration of the polymers which can lead to depletion layers (apparent slip lengths) of micrometres (Müller-Mohnssen, Weiss & Tippe Reference Müller-Mohnssen, Weiss and Tippe1990; Ma & Graham Reference Ma and Graham2005; Degré et al. Reference Degré, Joseph, Tabeling, Lerouge, Cloitre and Ajdari2006; Fang, Hsieh & Larson Reference Fang, Hsieh and Larson2007; Shin et al. Reference Shin, Bong and Kim2016). In general, non-hydrodynamic interaction between polymers and the surface can also contribute to large slip lengths (Mhetar & Archer Reference Mhetar and Archer1998; Gupta & Vanapalli Reference Gupta and Vanapalli2020). Regardless of the molecular mechanisms at play, we expect that such effects do not alter the phenomenology discussed here, including the interpretation in terms of an apparent inner angle 
$\theta _{{app},i} \to 0$ due to viscoelasticity. The presented theory provides a benchmark for the influence of the normal stress on receding contact-line motion, but molecular mechanisms might be required to provide an effective scale of regularisation.
For advancing contact lines, viscosity and normal stress provide forces that act in opposite directions: the normal stress actually facilitates contact-line motion. We have seen that we need to distinguish between weak and strong viscoelasticity. For weak viscoelasticity, there is a non-monotonic dependency of the macroscopic angle on contact-line velocity, as indicated by the orange dashed line in figure 5(c). The maximum arises due to a competition between viscous bending (increasing the angle) and viscoelastic bending (decreasing the angle). The cross-over to the strongly viscoelastic regime is given by 
$U \sim U_\psi$ (see (2.13)), which occurs when the normal stress effect is so strong that it bends the interface to form a precursor film of typical size 
$\ell _{VE}$. This strongly viscoelastic regime is indicated by the dashed blue line, and can be interpreted as an ‘apparent complete wetting’, since the contact-line motion becomes completely insensitive to 
$\theta _e$. The cross-over is schematically indicated by the black solid line. The exact details of the cross-over from weak to strong viscoelasticity are not universal, but depend on the material parameter 
$N_0$ defined in (3.2). For example, the speed at which the maximum appears in the weakly viscoelastic regime scales as 
$U_\psi /\sqrt {N_0}$; when 
$N_0$ is sufficiently large the maximum will thus emerge well before it is erased by the apparent complete wetting regime.
The presented analysis provides concrete predictions for experiments on wetting of viscoelastic liquids, rooted in a systematic long-wave expansion of the second-order fluid. Our findings suggest the possibility of a non-monotonic dependence of the contact angle with advancing velocity, which could potentially lead to instabilities. Such effects remain to be explored experimentally, but also deserve some further theoretical analysis. Namely, the calculations are based on the second-order fluid: while this lubrication model is equivalent to that obtained with Oldroyd-B fluid when viscoelasticity is introduced perturbatively (Ro & Homsy Reference Ro and Homsy1995; Datt et al. Reference Datt, Kansal and Snoeijer2022), its quantitative validity at strong viscoelasticity remains to be explored, more so when polymers in moving contact-line problems have been observed to be strongly stretched (Smith & Bertola Reference Smith and Bertola2010; Shin et al. Reference Shin, Bong and Kim2016). It would therefore be of interest to study contact-line motion for a variety of constitutive models.
Acknowledgements
The authors gratefully acknowledge discussions with D. Bonn and S. Kumar.
Funding
This work was supported by NWO through VICI grant no. 680-47-632.
Declaration of interests
The authors report no conflict of interest.
Open access
