Parametric oscillations of the sessile drop

Abstract

Waves are formed on the surface of a sessile drop driven through substrate vibrations oriented at a slanting angle from the normal. A mathematical model is derived, which leads to an infinite system of coupled Mathieu equations governing the wave dynamics that are solved using Floquet theory. The spatial structure of the waves is described by the mode number pair $[\ell,m]$, where $\ell$ and $m$ are the polar and azimuthal mode numbers, respectively. Limiting cases corresponding to horizontal and vertical vibrations are discussed with predictions agreeing well with prior literature. We focus our results on three drop motions – (1) harmonic $[1,1]$ rocking mode, (2) harmonic $[2,0]$ pumping mode, and (3) subharmonic rocking $[1,1]$ mode – as they depend upon the slanting angle, static contact angle, and contact-line conditions, which we assume to be either pinned or freely moving with fixed contact angle. New theoretical predictions are tested through experiments over a range of parameters, showing good agreement.

1. Introduction

Droplet manipulation on substrates is relevant to technologies such as inkjet printing (Calvert 2001), spray cooling (James, Smith & Glezer 2003; Vukasinovic, Smith & Glezer 2007b; Tsai et al. 2012) and phase-change heat transfer (Daniel, Chaudhury & Chen 2001; Macner, Daniel & Steen 2014), and is often facilitated by external driving forces associated with acoustic fields (Marston & Apfel 1980; Chastrette et al. 2022), aerodynamic effects (Bouwhuis et al. 2013), magnetic fields (Hill & Eaves 2012), electric fields (Feng & Beard 1991; Oh, Ko & Kang 2008) and substrate vibrations (Brunet, Eggers & Deegan 2007; Chang et al. 2015). The temporal response of the droplet is often determined by the particular type of applied forcing. In this paper, we study the oscillations of a sessile drop on a mechanically vibrated substrate whose axis of vibration is oriented at an angle $\chi$ from the normal. This gives rise to parametric oscillations typically associated with Faraday waves, which contribute to the unidirectional drop motion (Costalonga & Brunet 2020; Guo, Shui & Chen 2020) and climbing drops (Brunet et al. 2007; Sartori et al. 2015).

Faraday (1831) observed the formation of surface waves on a vertically vibrated liquid bath that oscillated with frequency $f$ one-half that of the driving frequency $f_{d}$, i.e. $f=f_d/2$, above a critical threshold driving acceleration $a$. This is termed a subharmonic response. Benjamin & Ursell (1954) showed that Faraday waves were solutions of a generalized Mathieu equation,

(1.1)\begin{equation} \ddot{\boldsymbol{y}}+\boldsymbol{c} \dot{\boldsymbol{y}}+ (\boldsymbol{P}-2A\boldsymbol{Q}\cos f_{d}t)\boldsymbol{y} =0, \end{equation}

where $\boldsymbol {y}$ is the surface disturbance. Floquet theory is often invoked to solve (1.1) and gives rise to a set of instability tongues for a subharmonic, harmonic and superharmonic wave response in which the wave is observed. There is a large volume of literature on Faraday waves on planar interfaces (Miles & Henderson 1990). However, there is limited research on Faraday waves on drops, particularly theoretical modelling. Ebo-Adou & Tuckerman (2016) use Floquet theory to predict the instability tongues for a free drop subject to radial forcing. Maksymov & Pototsky (2019) develop a simple hydrodynamic model of a one-dimensional pancake-shaped drop, and analytically determine the response amplitudes for the first harmonic and first superharmonic. To our knowledge, there has been no theoretical investigation of Faraday waves on a sessile drop with arbitrary contact angle $\alpha$ and slanting angle $\chi$. The aim of this study is to fill this gap in the literature. Unlike the case of radial forcing, slanted vibrations lead to a set of coupled Mathieu equations, which we solve using the Floquet theory method proposed by Kidambi (2013) for a brimful circular cylinder.

The frequency spectrum for the sessile drop exhibits a rich structure that has been predicted theoretically by Bostwick & Steen (2014) and verified experimentally under both terrestrial (Chang et al. 2015) and microgravity (McCraney et al. 2022) conditions, and numerically by Sakakeeny & Ling (2020). Here, the drop will oscillate with frequency $\lambda _{\ell,m}$ and mode shape defined by the mode number pair $[\ell,m]$, where $\ell$ is the polar mode number, and $m$ is the azimuthal mode number, restricted such that $\ell +m=\textrm {even}$ due to symmetry considerations associated with enforcement of the no-penetration condition at the solid support. For reference, the $\ell +m=\textrm {odd}$ modes have been predicted recently by Ding & Bostwick (2022), and are associated with a flux condition on the solid support. Figure 1 shows typical mode shapes, which can be classified by their mode number pair $[\ell, m]$ into zonal $[\ell,0]$, lateral $[\ell,1]$, sectoral $[\ell,\ell ]$ and tesseral $[\ell,m\neq \ell ]$ shapes. The spectrum of an inviscid sessile drop depends upon the static contact angle $\alpha$ and contact-line boundary condition. The special case of a hemispherical $\alpha =90^\circ$ base state with fully mobile contact line exhibits a high degree of symmetry and inherits the Rayleigh–Lamb spectrum for the free drop, including the well-known spectral degeneracy with respect to azimuthal mode number $m$ (Rayleigh 1879; Lamb 1924). Deviations from this special base state lead to spectral splitting and spectral reordering, which have been organized using a symmetry-breaking perspective into a periodic table of droplet motions by Steen, Chang & Bostwick (2019). Recent work by Zhang, Zhou & Ding (2023) has revealed symmetry-breaking due to gravitational effects, demonstrating that the lowest-frequency mode is transformed from a zonal mode to a sectoral mode.

Figure 1. Sessile drop mode shape $[\ell,m]$ classification into (a) zonal $[\ell,0]$, (b) lateral $[\ell,1]$, and (c) sectoral $[\ell,\ell ]$ modes.

Sessile drop mobility is achieved once the pinning forces at the three-phase contact line are overcome. Daniel & Chaudhury (2002) and Daniel et al. (2004) showed that a drop placed on a chemically anisotropic substrate will undergo translational motion once the contact-angle hysteresis has been overcome so as to mobilize the contact line. Here, substrate vibration is shown to help mobilize the contact line and induce such motions through resonant shape oscillations. This effect can become so strong that a drop may be driven up an inclined substrate against the weight of gravity (Brunet et al. 2007; Sartori et al. 2015). Noblin, Kofman & Celestini (2009) have shown that such motion can be decomposed into components associated with normal and tangential substrate vibrations. To explain and quantitatively predict directional drop motion, direct numerical simulations (Ding et al. 2018) and simplified flow models, e.g. lubrication theory (Savva & Kalliadasis 2014) and shallow-water models (Benilov & Billingham 2011; Bradshaw & Billingham 2016), have been put forth to simulate a drop climbing up a vibrating inclined plate. Costalonga & Brunet (2020) conducted quantitative experiments of directional motion of sessile drops induced by slanted vibrations to compare with existing models. It is generally accepted that to induce directional drop motion requires the simultaneous excitation of the asymmetric $[1,1]$ rocking mode and the symmetric $[2,0]$ pumping mode as the minimal ingredients (Noblin et al. 2009; John & Thiele 2010). It is possible to do so in an experiment with a single driving frequency. This is illustrated in figure 2, which plots the subharmonic and harmonic tongues for the rocking and pumping modes in the acceleration–frequency ($\varLambda \unicode{x2013}\omega$) parameter space for a typical set of experimental conditions. Note the overlap of the instability tongues for the subharmonic rocking mode and harmonic pumping mode, which gives rise to drop motion that we will refer to hereafter as a mixed mode. Thus the mixed mode is defined as the superposition of a subharmonic rocking mode and harmonic pumping mode. For reference, other mixed mode pairs have been observed experimentally that consist of a harmonic zonal mode and subharmonic tesseral or sectoral mode (Chang et al. 2015, figure 12). What is also interesting from figure 2 is that the harmonic $[1,1]$ rocking mode is isolated, which makes it easy to observe in experiment. Given the interest in such motions, we focus on presenting theoretical predictions for the $[1,1]$ and $[2,0]$ modes, and conduct preliminary experiments to verify these theoretical predictions.

Figure 2. Typical instability tongues for the $[1,1]$ and $[2,0]$ modes plotted in the acceleration–frequency ($\varLambda \unicode{x2013}\omega$) space. Here, the contact line is pinned, the slanting angle is $\chi =30^\circ$, and the contact angle is $\alpha =90^\circ$.

We begin this paper by describing our mathematical model of surface waves on a sessile drop subjected to slanted substrate vibrations. This gives rise to a set of coupled Mathieu equations, which are solved using Floquet theory to give the shape of the instability tongues, as they depend upon static contact angle $\alpha$, slanting angle $\chi$, and contact-line conditions. Experiments are performed to test theoretical predictions and are focused on the $[1,1]$ rocking and $[2,0]$ pumping modes, as relevant to drop transport on substrates. Theoretical predictions show good agreement with experimental observations for most cases, but there are exceptions, and we provide a rationalization for these differences. We finish by offering some concluding remarks.

2. Mathematical formulation

Consider a sessile drop defined by a base radius $r$ and contact angle $\alpha$ on a horizontal substrate, as shown in figure 3, and defined as

(2.1a–c)\begin{align} X(s,\varphi;\alpha)=\frac{\sin (s)}{\sin (\alpha)} \cos (\varphi),\quad Y(s,\varphi;\alpha)=\frac{\sin (s)}{\sin (\alpha)} \sin (\varphi),\quad Z(s;\alpha)= \frac{\cos (s) - \cos (\alpha)}{\sin (\alpha)},\end{align}

using arc-length-like $s\in [0,\alpha ]$ and azimuthal angle $\varphi \in [0,2{\rm \pi} ]$ as generalized surface coordinates (Bostwick & Steen 2014). The liquid is immersed in a passive gas, has density $\varrho$, and is assumed to be incompressible and inviscid. Gravitational effects are neglected. The plate is driven by an applied acceleration $a_{tt}$ that is oriented at an angle $\chi$ away from the vertical. This induces an interface disturbance $\eta$ and corresponding free surface flow that is described by a velocity field $\boldsymbol {U}$ and pressure field $P$ that satisfy the continuity and Euler equations, i.e.

(2.2a)\begin{gather} \boldsymbol{\nabla}\boldsymbol{\cdot} \boldsymbol{U}=0, \end{gather}

(2.2b)\begin{gather}\varrho\left(\frac{\partial}{\partial t}+\boldsymbol{U}\boldsymbol{\cdot} \boldsymbol{\nabla}\right) \boldsymbol{U}={-}\boldsymbol{\nabla}P+\varrho a_{tt}\boldsymbol{k}, \end{gather}

respectively.

Figure 3. Definition sketch of a sessile drop on a substrate that is mechanically vibrated at angle $\chi$ defined from the normal in (a) two-dimensional planar and (b) three-dimensional perspective views.

At the free interface, the interface shape $\eta$ is coupled to the flow field through the kinematic condition

(2.3)\begin{equation} \frac{\partial \eta}{\partial t}+\boldsymbol{U}\boldsymbol{\cdot}\boldsymbol{\nabla}\eta= \boldsymbol{U}\boldsymbol{\cdot}\boldsymbol{n},\end{equation}

and the Young–Laplace equation

(2.4)\begin{equation} P/\sigma=2 H,\end{equation}

which relates the pressure jump across the interface to the mean curvature $H$ there, with $\sigma$ the surface tension. Finally, the no-penetration condition is applied on the solid boundary:

(2.5)\begin{equation} \boldsymbol{U}\boldsymbol{\cdot}\hat{z}=0\lvert_{z=0}.\end{equation}

The orientation of the applied acceleration gives rise to both horizontal and vertical forcing, and can be expressed in component form as $a_{tt}\boldsymbol {k}=\cos \chi \,a_{vtt}\boldsymbol {z}+\sin \chi \, a_{htt}\boldsymbol {x}$. Here, we assume $a_{vtt}=a_{htt}=({\sigma }/{\varrho r^2})\varLambda \cos \omega t$, with $\varLambda \equiv A/r$, where $A$ is the driving amplitude, and $\omega$ is the driving frequency.

2.1. Linearized equations

Equation (2.2) admits the base state

(2.6a,b)\begin{equation} \hat{\boldsymbol{U}}=0,\quad \boldsymbol{\nabla}\hat{P}= \frac{\sigma}{r^2}\,\varLambda(\cos\chi\,\boldsymbol{z}+\sin \chi\,\boldsymbol{x})\cos \omega t . \end{equation}

We perform a linear stability analysis about this base state by decomposing the velocity and pressure fields as

(2.7a)\begin{gather} \boldsymbol{U}=\hat{\boldsymbol{U}}+ \epsilon\,\boldsymbol{u}= \epsilon \boldsymbol{\nabla}\phi, \end{gather}

(2.7b)\begin{gather}P=\hat{P}+ \epsilon p, \end{gather}

where $\epsilon$ is a small parameter. Here, we have assumed an irrotational velocity field that is described by the velocity potential $\phi$, such that $\boldsymbol {u}=\boldsymbol {\nabla } \phi$. Substituting (2.7) into (2.2) leads to the linearized equations

(2.8a)\begin{gather} \nabla^2\phi=0, \end{gather}

(2.8b)\begin{gather}p={-}\varrho\,\frac{\partial \phi}{\partial t}, \end{gather}

on the drop domain. Similarly, the kinematic condition (2.3) reduces to

(2.9)\begin{equation} \frac{\partial \phi}{\partial n}=\frac{\partial \eta}{\partial t},\end{equation}

and the no-penetration condition (2.5) reduces to

(2.10)\begin{equation} \frac{\partial \phi}{\partial z}=0.\end{equation}

The Young–Laplace equation (2.4) can be written as

(2.11)\begin{equation} P|_{\rho=R+ \eta}\approx\hat{P}|_{\rho=R}+ p|_{\rho=R}+\eta \left.\left(\frac{\partial \hat{P}}{\partial \rho}\right)\right|_{\rho=R}= \sigma(\boldsymbol{\nabla}\boldsymbol{\cdot}\boldsymbol{n})_{\rho=R+ \eta},\end{equation}

using a spherical coordinate system $(\rho,\theta, \varphi )$ with the static drop surface defined at $\rho =R$. These linearized equations can be simplified by defining $\hat {P}|_{\rho =R}=\sigma (\boldsymbol {\nabla }\boldsymbol {\cdot }\boldsymbol {n})_{\rho =R}$ and $p|_{\rho =R}=-\varrho ({\partial \phi }/{\partial t})$, and performing some algebraic manipulations to get

(2.12)\begin{align} & \varrho\,\frac{\partial \phi}{\partial t}+\sigma(\boldsymbol{\nabla}\boldsymbol{\cdot} \boldsymbol{n})_{\rho=R+ \eta}-\sigma(\boldsymbol{\nabla}\boldsymbol{\cdot}\boldsymbol{n})_{\rho=R} \nonumber\\ &\quad - \frac{\sigma}{ r^2}\,\varLambda(\cos\chi \cos\theta+\sin\chi \sin\theta \cos\varphi)\cos\omega t\eta=0. \end{align}

2.2. Dimensionless equations

Dimensionless variables are introduced,

(2.13a–e)\begin{equation} \rho^{*}=\rho/r,\quad \eta^{*}=\eta/r,\quad t^{*}=t\,\sqrt{\frac{\sigma}{\varrho r^3}},\quad \phi^{*}=\phi\,\sqrt{\frac{\varrho}{\sigma r}},\quad p^{*}=p\,\frac{r}{\sigma}, \end{equation}

and applied to the governing equations. Herein, we drop the $*$ and refer to dimensionless quantities. The resulting equations to be solved are Laplace's equation on the domain

(2.14)\begin{equation} \nabla^2\phi=0, \end{equation}

with no-penetration condition on the solid support $z=0$,

(2.15)\begin{equation} \frac{\partial \phi}{\partial z}=0, \end{equation}

and kinematic condition

(2.16)\begin{equation} \frac{\partial \phi}{\partial n}=\frac{\partial \eta}{\partial t} \end{equation}

and Young–Laplace equation

(2.17)\begin{align} & \frac{\partial \phi}{\partial t}-\sin^2 (\alpha)\left(\frac{\partial^2\eta}{\partial s^2}+\cot (s)\, \frac{\partial\eta}{\partial s}+2\eta+\frac{1}{\sin^2(s)}\,\frac{\partial^2\eta}{\partial \psi^2}\right) \nonumber\\ &\quad -\varLambda \eta (\cos\chi \cos\theta+\sin\chi \sin\theta \cos\varphi) \cos \omega t=0 \end{align}

on the free surface $\rho =R$. The governing equations are augmented with a contact-line boundary condition applied at $s=\alpha$, which we assume to be either pinned or free to move, i.e.

(2.18a)\begin{gather} \eta=0,\quad \text{pinned condition}, \end{gather}

(2.18b)\begin{gather}\frac{\partial \eta}{\partial s}-(\cos \alpha) \eta=0,\quad \text{free condition}, \end{gather}

and an integral condition

(2.19)\begin{equation} \int_{\varGamma}\frac{\partial \phi}{\partial n}\,{\rm d}\varGamma=0,\end{equation}

necessary to ensure volume conservation.

2.3. Derivation of the coupled Mathieu equations

We seek a solution $(\eta, \phi )$ to the governing equations (2.14), (2.15), (2.16), (2.17), (2.18), (2.19) in the form

(2.20a,b)\begin{align} \eta(x,t)=\cos m\varphi\sum_{\ell =1}^{\infty}a_{\ell}(t)\,V^{m}_{\ell}(x),\quad \phi(\rho,\theta,t)=\cos m\varphi\sum_{\ell =1}^{\infty}b_{\ell}(t)\, \rho^{\ell}\,P^{m}_{\ell}(\cos\theta),\end{align}

where $x=\cos s$. Here, $V^{m}_{\ell }$ are a set of orthonormal basis functions chosen to satisfy the contact line (2.18) and integral conditions (2.19). The idea is to choose basis functions from the natural oscillations problem, as described by Bostwick & Steen (2014, (4.10)), as

(2.21)\begin{equation} V^{m}_{\ell}(x)=\int_{b}^{1}G(x,y;m)\,(\rho^{\ell}\,P^{m}_{\ell}(y)+ \delta_{m,0}C)\,{{\rm d}y},\quad \ell=1,2,\ldots, N, \end{equation}

where the Green's function is defined in

(2.22) \begin{equation} G(x,y;m) = \left\{\begin{array}{@{}ll} \displaystyle\dfrac{1}{1-y^2}\,\dfrac{U(x;m)\,X(y;m)}{W(y;m)}, & b< y< x<1, \\ \displaystyle\dfrac{1}{1-y^2}\,\dfrac{U(y;m)\,X(x;m)}{W(y;m)}, & b< x< y<1 . \end{array}\right. \end{equation}

Here, $b \equiv \cos \alpha$, and $U$ and $X$ are the homogeneous solutions of

(2.23)\begin{equation} \left(2-\frac{m^2}{\sin^2(s)}\right)\frac{\partial \phi}{\partial n}+ \cot(s) \left(\frac{\partial \phi}{\partial n}\right)' + \left(\frac{\partial \phi}{\partial n} \right)''=0 \end{equation}

that satisfy the associated boundary conditions

(2.24a,b)\begin{equation} U=y_{1}(x;m),\quad X=y_{2}(x;m)-\frac{\tau_{2}(m)}{\tau_{1}(m)}\,y_{1}(x;m), \end{equation}

and $W$ is the Wronskian of the solutions $U$ and $X$.

The Green's function is parametrized by the azimuthal mode number $m$ and the transformed contact angle $b$. The functions $y_{1}$ and $y_{2}$ in the solutions $U$ and $X$ are given by

(2.25) \begin{align} \left.\begin{gathered} y_{1}(x;0)=P_{1}(x),\quad y_{2}(x;0)=Q_{1}(x),\quad y_{1}(x;1)=P_{1}^1(x),\quad y_{2}(x;1)=Q_{1}^1(x),\\ y_{1}(x;m\geq2)=(x+m)\left(\frac{1-x}{1+x}\right)^{m/2},\quad y_{2}(x;m\geq2)=\frac{x-m}{2m(m^2-1)}\left(\frac{1+x}{1-x}\right)^{m/2}, \end{gathered}\right\} \end{align}

where $P_{1}$, $Q_{1}$ and $P_{1}^1$, $Q_{1}^1$ are the Legendre functions of orders 0 and 1 and index 1, respectively. The parameters $\tau _{1}$ and $\tau _{2}$ are related to the contact-line boundary conditions, with

(2.26a)\begin{equation} \tau_{1}^{p}=y_{1}(b;m),\quad \tau_{2}^{p}=y_{2}(b;m),\end{equation}

for the pinned contact-line disturbance (superscript $p$), and

(2.26b)\begin{equation} \tau_{1}^{n}=y_{1}'(b;m)+\frac{b}{\sqrt{1-b^2}}\,y_{1}(b;m),\quad \tau_{2}^{n}=y_{2}'(b;m)+\frac{b}{\sqrt{1-b^2}}\,y_{2}(b;m),\end{equation}

for the natural contact-line disturbance (superscript $n$). The constant $C$ is defined as

(2.27)\begin{equation} C={-}\frac{\displaystyle\int_{b}^{1}\displaystyle\int_{b}^{1}G(x,y;0)\, \phi(y)\,{\rm d} y\,{\rm d}{\kern0.8pt}x}{\displaystyle\int_{b}^{1}\int_{b}^{1}G(x,y;0)\,{\rm d} y\,{\rm d}{\kern0.8pt}x}. \end{equation}

Substituting (2.20a,b) into (2.16) and (2.17) gives

(2.28a)\begin{equation} \sum_{\ell =1}^{\infty}b_{\ell}(t)\,\frac{\partial (\rho^{\ell}\,P^{m}_{\ell}(\cos\theta))}{\partial n}= \frac{{\rm d} a_{\ell}(t)}{{\rm d} t}\,V^{m}_{\ell}(x) \end{equation}

and

(2.28b) \begin{align} &\cos m \varphi \sum_{\ell =1}^{\infty}\frac{{\rm d} b_{\ell}(t)}{{\rm d} t}\, (\rho^{\ell}\,P^{m}_{\ell}(\cos\theta)) \nonumber\\ &\quad {}\times\left(-\sin^2\alpha\left((1-x^2)\,\frac{\partial^2 }{\partial x^2}-2x\, \frac{\partial}{\partial x}+2-\frac{m^2}{1-x^2}\right)V^{m}_{\ell}(x)\right)a_{\ell}(t) \nonumber\\ &\quad {}- V^{m}_{\ell}(x)\,a_{\ell}(t)(\varLambda\cos\chi \cos\theta+\varLambda\sin\chi \sin\theta \cos\varphi)\cos\omega t=0, \end{align}

respectively. We can project these equations onto the function space $\cos m\varphi \,V^{m}_{k}(x), k=1,$ $2,\ldots, N,$ via the inner product $\langle{\kern1.8pt}f,g\rangle =\int _{b}^{1}f(x)\,g(x)\,\textrm {d}{\kern0.8pt}x$, noting the identity

(2.29)\begin{equation} \int_{0}^{2{\rm \pi}}\int_{b}^{1}\cos^2 m \varphi\,\varLambda_{0}\,V^{m}_{\ell}(x)\, V^{m}_{k}(x)\sin\chi \sin\theta \cos\varphi \cos\omega t=0. \end{equation}

We project (2.28a) onto the function space to give

(2.30)\begin{equation} b_{\ell}(t)=\sum^{\infty}_{k=1}\beta_{\ell k}\,\frac{{\rm d} a_{k}}{{\rm d} t},\end{equation}

with

(2.31)\begin{equation} \beta_{\ell k}=\left(\int_{b}^{1}\frac{\partial (\rho^{\ell}\,P^{m}_{\ell} (\cos\theta))}{\partial n}\,V^{m}_{k}(x)\,{{\rm d}{\kern0.8pt}x}\right)^{{-}1} \int_{b}^{1} V_{\ell}^{m}(x)\,V_{k}^{m}(x)\,{{\rm d}{\kern0.8pt}x}. \end{equation}

Similarly, projecting (2.28b) onto the function space $V^{m}_{\ell }$ yields

(2.32)\begin{align} & \sum_{\ell =1}^{\infty}\frac{{\rm d} b_{\ell}(t)}{{\rm d} t}\,(\rho^{\ell}\,P^{m}_{\ell}(\cos\theta)) \nonumber\\ &\quad \times\left(-\sin^2\alpha\left((1-x^2)\,\frac{\partial^2 }{\partial x^2}-2x\, \frac{\partial}{\partial x}+2-\frac{m^2}{1-x^2}\right)V^{m}_{\ell}(x)\right)a_{\ell}(t) \nonumber\\ &\quad -\varLambda\cos\chi \cos\theta \cos\omega t\,V^{m}_{\ell}(x)\,a_{\ell}(t)=0. \end{align}

Finally, we can project (2.32) onto the function space $V^{m}_{p}$ and use (2.30) to obtain a system of $N$ coupled second-order ordinary differential equations:

(2.33)\begin{equation} \sum_{\ell =1}^{N}\gamma_{p\ell}\sum^{N}_{k=1}\beta_{\ell k}\, \frac{{\rm d}^2 a_{k}}{{\rm d} t^2}+(\tau_{p\ell}-\varLambda\cos \chi\,H_{p\ell}\cos\omega t)a_{\ell}=0,\quad p=1,2,\ldots, N, \end{equation}

with

(2.34a)\begin{gather} \gamma_{p\ell}=\int_{b}^{1}\rho^{\ell}\,P^{m}_{\ell}(\cos\theta)\,V^{m}_{p}(x)\,{{\rm d}{\kern0.8pt}x}, \end{gather}

(2.34b)\begin{gather}\tau_{p\ell}=\sin^2\alpha\int_{b}^{1}-\left((1-x^2)\,\frac{\partial^2 }{\partial x^2}-2x\, \frac{\partial}{\partial x}+2-\frac{m^2}{1-x^2}\right)V^{m}_{\ell}(x)\,V^{m}_{p}(x)\,{{\rm d}{\kern0.8pt}x}, \end{gather}

(2.34c)\begin{gather}H_{p\ell}=\int_{b}^{1} \cos\theta\,V^{m}_{\ell}(x)\,V^{m}_{p}(x)\,{{\rm d}{\kern0.8pt}x}. \end{gather}

As is customary, we have kept only a finite number of modes for $N$ for numerical computation. Equation (2.33) can be written in matrix form as

(2.35) \begin{equation} \boldsymbol{\mathsf{A}}\,\frac{{\rm d}^2\boldsymbol{y}}{{\rm d} t^2}+ (\hat{\boldsymbol{\mathsf{P}}}-2\varLambda\cos\chi\,\hat{\boldsymbol{\mathsf{Q}}}\cos\omega t)\boldsymbol{y}=0, \end{equation}

where $\boldsymbol {y}=(a_{1},a_{2},\ldots,a_{N})$ is a coefficient vector, and the components of the matrices $\boldsymbol{\mathsf{A}}$ and $\hat {\boldsymbol{\mathsf{P}}}$ are given by

(2.36a–c) \begin{equation} {\mathsf{A}}_{ij}=\sum_{\ell=1}^{N}\gamma_{i\ell}\beta_{\ell j},\quad \hat{{\mathsf{P}}}_{ij}=\tau_{ij},\quad \hat{{\mathsf{Q}}}_{ij}={\mathsf{H}}_{ij}/2. \end{equation}

To facilitate a solution, we rescale (2.35) by defining $T=\omega t/2$, $\boldsymbol{\mathsf{P}}=4\boldsymbol{\mathsf{A}}^{-1}\hat {\boldsymbol {\boldsymbol{\mathsf{P}}}}$, $\boldsymbol{\mathsf{Q}}=4\boldsymbol{\mathsf{A}}^{-1}\hat {\boldsymbol{\mathsf{Q}}}$, which gives

(2.37) \begin{equation} \omega^2\,\frac{{\rm d}^2\boldsymbol{y}}{{\rm d} T^2}+ (\boldsymbol{\mathsf{P}}-2\varLambda\cos\chi\,\boldsymbol{\mathsf{Q}}\cos 2T)\boldsymbol{y}=0. \end{equation}

2.4. Solution method

Floquet theory is used to determine the stability of (2.37). We follow the method suggested by Kumar & Tuckerman (1994) and Kidambi (2013) by fixing the values of the Floquet exponent $\mu$ for both harmonic and subharmonic tongues to trace out the shape of the instability tongues. In general, this is a computationally efficient method.

We begin by seeking a solution of the form

(2.38)\begin{equation} \boldsymbol{y}(T)={\rm e}^{\mu T}\left(\sum_{l={-}L}^{L}\boldsymbol{\xi}_{l}\exp({{\rm i}2lT})\right),\end{equation}

with $\boldsymbol {\xi }_{l}$ a vector with components $\xi _{l}^{n}$, $n=1,\ldots, N$. Substituting (2.38) into (2.37) gives

(2.39) \begin{equation} \sum_{l={-}L}^{L}\left(\omega^2(\mu+i2l)^2 \xi_{l}^{n}+ \sum_{j=1}^{N}({\mathsf{P}}_{nj}-2\varLambda \cos\chi\,{\mathsf{Q}}_{nj}\cos2T)\xi_{l}^{j}\right)=0. \end{equation}

For a given $\varLambda$, we can compute the stability boundary by setting $\mu =0$ for the harmonic case and $\mu =i$ for the subharmonic case, and solving for the unknown $\omega$ from the generalized eigenvalue problem $(\boldsymbol {M}\omega ^2+\boldsymbol {W})\boldsymbol {v}=0$, where $\boldsymbol {v}=(\xi _{-L}^{1},\xi _{-L+1}^{1},\ldots,\xi _{L-1}^{1},\xi _{L}^{1}, \xi _{-L}^{2},\ldots,\xi _{L-1}^{N},\xi _{L}^{N})^{T}$.

3. Theoretical predictions

For the results presented here, we use a truncation $N=10$ and $L=20$ to produce sufficient convergence for the Floquet theory calculations. We begin by verifying our model against prior results in the literature. We then examine the limiting cases of horizontal and vertical vibrations, before discussing the more general case of slanted vibrations. Our focus is on the low-frequency $[1,1]$ and $[2,0]$ modes as the main contributors to directional drop transport, as discussed in the Introduction.

3.1. Radial vibration verification of solution method

To our knowledge, only Ebo-Adou & Tuckerman (2016) have considered the parametric oscillations of a free spherical drop, and they focus on the case where motion is driven by a time-harmonic radial acceleration. A sessile drop with hemispherical base state $\alpha =90^\circ$ and free contact line has a subset of motions that may be extended smoothly by symmetry operation to the full drop (Steen et al. 2019). This allows us to verify our Floquet theory. It is straightforward to adapt our analysis to this limiting case, where (2.17) becomes

(3.1)\begin{equation} \frac{\partial \phi}{\partial t}-\sin^2 (\alpha)\left(\frac{\partial^2\eta}{\partial s^2}+\cot (s)\, \frac{\partial\eta}{\partial s}+2\eta+\frac{1}{\sin^2(s)}\, \frac{\partial^2\eta}{\partial \psi^2}\right) -\varLambda \eta \cos \omega t=0. \end{equation}

The steps from § 2.3 may then be applied to generate a set of uncoupled Mathieu equations from which the subharmonic and harmonic tongues are computed by Floquet theory. Figure 4 compares the shape of the subharmonic instability tongues for the $[2,0]$, $[4,0]$, $[6,0]$ hemispherical sessile drop modes to those for the free drop, which reveals excellent agreement, verifying our Floquet theory for sessile drops.

Figure 4. Verification of Floquet theory for sessile drops with radial forcing. Subharmonic instability tongues for the $[2,0]$, $[4,0]$, $[6,0]$ hemispherical sessile drop modes with free contact line compared with those predicted by Ebo-Adou & Tuckerman (2016) for the full spherical drop.

3.2. Horizontal vibration

Horizontal vibration $\chi =90^\circ$ is a limiting case where (2.37) becomes

(3.2)\begin{equation} \omega^2\,\frac{{\rm d}^2\boldsymbol{y}}{{\rm d} T^2}+\boldsymbol{\mathsf{P}}\boldsymbol{y}=0. \end{equation}

Here, there is no parametric term implying that the temporal response is harmonic, which has been observed experimentally by Khan & Eslamian (2019) for surface waves on a water layer with pinned contact line, and argued mathematically by Or (1997) for a liquid film with infinite lateral extent. We similarly report a harmonic drop response in experiments for slanting angle close to $\chi =90^\circ$, which we will discuss in the next subsection. For reference, the frequency spectrum for the hemispherical sessile drop subject to horizontal forcing has been reported by Lyubimov, Lyubimova & Shklyaev (2004).

3.3. Vertical vibration

Vertical vibration $\chi =0^\circ$ is another limiting case that is relevant to most experimental studies. Here, (2.37) becomes

(3.3)\begin{equation} \omega^2\,\frac{d^2\boldsymbol{y}}{dT^2}+(\boldsymbol{\mathsf{P}}-2\varLambda\boldsymbol{\mathsf{Q}}\cos2 T)\boldsymbol{y}=0. \end{equation}

Figure 5 plots the instability tongues for the three modes $[2,0]$, $[1,1]$, $[2,2]$ with smallest natural frequencies for the case of vertical vibration. For both harmonic and subharmonic tongues, note that the pinned contact-line modes are shifted to higher frequency relative to the free contact-line modes. This is due to the fact that the pinned contact line is more constrained than the free contact line, which is known to lead to higher frequencies according to variational principles (Bostwick & Steen 2015). Also, the instability tongues for the free contact line are notably wider than those for the pinned contact line, a feature that is more prominent for the harmonic tongues (cf. figure 5b). A comparison between modes shows that the $[2,0]$ mode has the largest instability tongue, which we attribute to the fact that zonal modes are dominated by vertical motions, making them easier to be excited by vertical vibration. Figure 6 shows how the shape of the instability tongues depends upon the contact angle $\alpha$. Here, the bandwidth increases with $\alpha$, suggesting that larger drop inertia makes it easier to excite waves.

Figure 5. (a) Subharmonic and (b) harmonic instability tongues for the $[2,0]$, $[1,1]$ and $[2,2]$ modes plotted in the acceleration–frequency ($\varLambda - \omega$) space, contrasting pinned and free contact-line conditions for $\alpha =60^\circ$.

Figure 6. (a,c) Subharmonic and (b,d) harmonic instability tongues for the (a,b) $[2,0]$ and (c,d) $[1,1]$ modes plotted in the acceleration–frequency ($\varLambda - \omega$) space, depending upon contact angle ($\alpha =60^\circ, 90^\circ, 120^\circ$).

3.4. Slanted vibration

Equation (2.37) describes the most general case for slanted vibration, as defined by the slanting angle $\chi$. Here, the $\cos \chi$ coefficient on the parametric forcing term strongly affects the shape of the instability tongues. Our focus is on the $[1,1]$ harmonic mode, $[2,0]$ harmonic mode and $[1,1]$ subharmonic mode, as these are the easiest to excite in experiment for small driving frequency and related to translational drop motion. Figure 7 plots the instability tongues for a hemispherical $\alpha =90^\circ$ drop showing an overlap of the $[1,1]$ subharmonic mode with the $[2,0]$ harmonic mode. This occurs because the natural frequencies for these spatial modes are related by $2\omega _{1,1}\approx \omega _{2,0}$ for $\alpha =90^\circ$. This gives rise to the mixing mode, which we will discuss in detail in the next section when we describe our experimental results. We also note that the bandwidth for the instability tongues decreases with increasing slanting angle $\chi$, which makes sense given that the coefficient $\cos \chi$ on the parametric term also decreases with $\chi$.

Figure 7. Instability tongues for $[2,0]$ harmonic and $[1,1]$ subharmonic modes plotted in the acceleration–frequency ($\varLambda - \omega$) space, depending upon the slanting angle $\chi =0^\circ$ and $60^\circ$, with pinned contact line and $\alpha =90^\circ$.

4. Experiment

Experiments were performed and compared with the theoretical predictions discussed previously. Given the recent interest in directional drop motion on substrates, we focus on characterizing the motions for the $[2,0]$ harmonic pumping mode and the $[1,1]$ subharmonic rocking mode for pinned contact-line conditions, as they depend upon the slanting angle $\chi$ and drop contact angle $\alpha$.

4.1. Experimental procedure

Experiments are performed using the set-up shown in figure 8. A droplet of distilled water is deposited onto a substrate through a needle using an AL-1000 syringe pump at low flow rate to induce pinch-off of a droplet of known volume $V$. The liquid properties of distilled water include the density $\rho =986.2\pm 2.4\,\textrm {kg}\,\textrm {m}^{-3}$, surface tension $\sigma = 71.06\pm 0.01\,\textrm {mN}\,\textrm {m}^{-1}$ and dynamic viscosity $\mu = 1\,\textrm {mPa}\,\textrm {s}$. The droplet contact line is pinned by a circular depression with depth $0.1\,\textrm {mm}$ and diameter $D=3.3\,\textrm {mm}$, which is machined into a solid piece of either aluminium or copper. The substrate is secured to a tilting angle table (TAT) (Sherline 3750-LAZ), which consists of two plates separated by angle $\chi$, the bottom of which is mounted on a Labworks ET-139 mechanical shaker. The shaker is driven by a Labworks power amplifier (PA-138) with feedback from an accelerometer (PCB J352C33) that is controlled by VibeLab Pro software (VL144x-7.09) to provide sinusoidal vibrations with acceleration $a$ and frequency $f$ at an angle $\chi \in [0^\circ, 90^\circ ]$ with limiting cases of vertical $\chi =0^\circ$ and horizontal $\chi =90^\circ$ vibration.

Figure 8. Experimental schematic shown in (a) front and (b) side views.

Droplet dynamics are captured by a Phantom VEO 410L high-speed camera at over 1000 fps. The camera is positioned on a three-axis stage and equipped with a LAOWA Ultra Macro $25\,\textrm {mm}$ lens. An optical diffuser (Edmund optics, $127\times 178\,\textrm {mm}$ sandblasted glass) is used to diffuse the backlight (MultiLED LT-V9-15). Both the shaker and the camera are placed on separate optic tables to minimize interference. The shaker rests on a $2.5\,\textrm {cm}$ thick vibration isolation pad for bench isolation. Image calibration is carried out through ImageJ using the outer diameter of the needle as the known length scale. An in-house MATLAB code is developed in which the drop image undergoes Gaussian filtering, a morphological closing operation using disk-shaped structuring elements, and finally binarization. The contact angle $\alpha$ is measured by fitting the interface shape to a second-order polynomial for the 50 closest points adjacent to the contact line, and computing the slope of the tangent line to give $\alpha$. Given $\alpha$ and the fixed wetted radius $r=D/2$, the drop volume can be determined by assuming a spherical cap shape to give $V$. The shape change dynamics are tracked using the binarized image from which the interface shape and centre of mass are readily computed. The resulting signal is subjected to a fast Fourier transform to extract the resonant frequency. Following each experiment, the droplet is removed using compressed air, and periodic cleaning of the substrate is performed using an ultrasonic cleaner (MGUC500).

4.2. Results

Two types of motion are observed in experiment: (1) low driving amplitude harmonic pure modes, and (2) high driving amplitude mixed modes. The latter involves the simultaneous excitation of a subharmonic mode and a harmonic mode at a single driving frequency. Here, we quantify the drop response for these motions, as it depends upon the experimental parameters.

4.2.1. Pumping and rocking modes

Pure $[2,0]$ pumping and $[1,1]$ rocking modes are excited at comparatively low driving acceleration ($a=0.25\,g$ in this case), as shown in figures 9 and 10, respectively. The pumping mode is dominated by vertical motion, and the rocking mode by horizontal motion.

Figure 9. Harmonic pumping mode $[2,0]$ over a complete cycle of oscillation with period $T$, which has been excited at acceleration $a=0.25\,g$ and frequency $f=128\,\textrm {Hz}$ using vertical $\chi =0^\circ$ vibration, with contact angle $\alpha =60.7^\circ$ and volume $V=4.6$ ml.

Figure 10. Harmonic rocking mode $[1,1]$ over a complete cycle of oscillation with period $T$, which has been excited at acceleration $a=0.25\,g$ and frequency $f=63.5\,\textrm {Hz}$ using horizontal $\chi =90^\circ$ vibration, with contact angle $\alpha =60.7^\circ$ and volume $V=4.6$ ml.

4.2.2. Harmonic instability tongues

The harmonic rocking mode is isolated, as there is no adjacent mode to interact with. Here, the wave emerges gradually, making it experimentally challenging to identify the boundary of these tongues. This is in stark contrast to the threshold instability typically associated with the emergence of subharmonic Faraday waves. As such, we employ an artificial criterion to detect the harmonic drop response, which we use to present an equivalent harmonic instability tongue. The procedure that we use is to construct a frequency response curve by tracking the dynamic contact angles $\alpha _L,\alpha _R$ on the left and right sides of the drop during a frequency sweep with fixed driving acceleration, as shown in figure 11. This metric is more suitable for our experiments than tracking the centre of mass (COM) for two main reasons: (1) any rocking motion corresponds to a change in dynamic contact angles, and (2) this method focuses on accurately capturing the region near the contact points, whereas the calculation of the COM relies on the entire drop shape, making it susceptible to the influence of any defects in the interface shape due to light reflection. Figure 11(a) plots the time trace of the dynamic contact angles with the contact-angle difference $|\alpha _R-\alpha _L|$, an effective metric to quantify the drop response. Figure 11(b) plots the time trace of the contact-angle difference and identifies the local maxima (red dots) with associated average value (red solid line), as well as the average value of the overall signal (blue dashed line). This gives the average value representing the maximum difference in dynamic contact angles for that specific frequency, which we repeat across multiple frequencies (with maximum interval $2\,\textrm {Hz}$) to generate the frequency response curves shown in figure 11(c) for different slanting angles $\chi$. Each response curve exhibits a peak near the resonant frequency and an associated bandwidth, consistent with the resonance plots of Celestini & Kofman (2006). The bandwidth can be defined by choosing some fraction of the maximum drop response, for example 65 % shown as the green dashed line. We note that this choice is arbitrary and will change the corresponding value of the bandwidth, but it will not change the interpretation of the drop response if the same metric is used for all response curves. This procedure is repeated for various driving accelerations.

Figure 11. Procedure for obtaining the harmonic frequency response curve by measuring the time trace of (a) the dynamic contact angles at two sides of the drop, $\alpha _L,\alpha _R$, for $a= 0.25\,g$, $f=63\,\textrm {Hz}$ and $\chi = 60^\circ$. (b) The difference in contact angles $|\alpha _R-\alpha _L|$ is used as a metric for the frequency response, where the red dots denote the peak points of the signal, and the red line is the average value of these points, while the dashed blue line represents the overall average. (c) Frequency response curve plotting $|\alpha _R-\alpha _L|$ against driving frequency $f$ for two different slanting angles, $\chi = 30^\circ$ (squares) and $\chi = 60^\circ$ (circles). The green dashed line corresponds to 65 % of the peak value of the response curve at $\chi = 60^\circ$.

Figure 12 plots the bandwidth, computed using the procedure described in the previous paragraph, in the acceleration–frequency domain for two different slanting angles ($\chi = 30^\circ$ and $60^\circ$) to produce an equivalent harmonic instability tongue. It is important to note that we have not observed a pronounced rocking mode below $a=0.05\,g$ and $a=0.1\,g$ for slanting angles $\chi = 60$ and $30^\circ$, respectively, for $\alpha = 60.7^\circ$ and $V=4.6$ ml. Below these critical accelerations, the size of the dynamic contact angle difference is comparable to the contact-angle hysteresis, and the periodic signature of the rocking mode is not observed in the corresponding time trace. Here, decreasing the slanting angle $\chi$ increases the threshold acceleration where the rocking mode is prominent, while also decreasing the bandwidth. This is because the horizontal component of the vibration force, which contributes to the rocking mode, decreases as the slanting angle decreases, shifting the response curves downwards (cf. figure 11c). Consequently, we observed no harmonic rocking mode at $\chi = 0^\circ$. Figure 6(d) shows that there exists a theoretical rocking mode at vertical forcing ($\chi = 0^\circ$), but its bandwidth is so small that is difficult to observe in our experiments. Increasing the static contact angle $\alpha$ decreases the threshold acceleration and increases the bandwidth of the instability tongue, consistent with our theoretical results (cf. Figure 6d).

Figure 12. Harmonic instability tongue for the $[1,1]$ rocking mode with $\alpha = 60.7^\circ$ and $V=4.6$ ml for $\chi = 60^\circ$ (red) and $\chi =30^\circ$ (blue), and $\alpha = 74.2^\circ$ and $V=6.5$ ml for $\chi =30^\circ$ (green).

It is important to note that hysteresis is a factor in our measurements. To be precise, the measurements obtained during the sweeping process (over either frequency or acceleration) in the increasing direction differ from those obtained in the reverse direction due to the influence of wave size and nonlinearities (Benjamin & Ursell 1954). To mitigate this, we conduct each measurement using a new drop, with the sweep performed in the increasing direction.

4.2.3. Mixed mode

The close proximity of the instability tongues for the subharmonic $[1,1]$ rocking mode and the harmonic $[2,0]$ pumping mode leads to the mixed mode shown in figure 13. The observed subharmonic response implies that it takes time $2T$ for the drop to complete one cycle of oscillation, whereas the substrate oscillation completes one cycle in time $T$. A characteristic feature of the mixing mode is that the COM has both horizontal and vertical components, as shown in the time trace of figure 14(a), with the subharmonic behaviour attributed to the horizontal motion (i.e. rocking mode) of the drop COM, and the harmonic behaviour to the vertical motion (i.e. pumping mode) of the drop COM. This is more clearly illustrated by plotting the COM motion in the phase space, as shown in figure 14(b), which demonstrates that two cycles of vertical motion occur for every one cycle of horizontal motion. Note that the trajectory is not symmetric because the drop motion is biased in the forcing direction. The shape of this trajectory is expected to vary with the slanting angle, among other parameters, but such an exploration is beyond the scope of this paper.

Figure 13. Mixed mode excited with $a=2.2\,g$, $f=128\,\textrm {Hz}$, $\chi =60^\circ$, $\alpha =60.7^\circ$ and $V=4.6$ ml. The period of the substrate oscillation is $T$, whereas the period of the drop oscillation is $2T$.

Figure 14. The COM motion for the mixing mode plotted (a) as a time trace and (b) in the phase space, showing the trajectory corresponding to the experiment shown in figure 13.

The subharmonic instability tongues are obtained by conducting an acceleration sweep at a constant driving frequency. As the driving acceleration surpasses the threshold value, the COM of the drop experiences a sudden and pronounced horizontal oscillation corresponding to the emergence of the subharmonic $[1,1]$ rocking mode. This process is repeated over a range of driving frequencies. Figure 15(a) plots the instability tongues comparing experiment (dashed line type) to theory (solid line type), depending upon the slanting angle $\chi$. These curves shift upwards with increasing $\chi$, indicating that higher driving accelerations are required to excite the mixing mode. This trend is further highlighted in figure 15(b), which plots the normalized threshold acceleration $a_{th}/g$ for the mixed mode against $\chi$, showing an increasing trend. For the special case of horizontal forcing ($\chi =90^\circ$), we did not observe a subharmonic $[1,1]$ rocking mode, which agrees well with our theoretical predictions (§ 3.2).

Figure 15. (a) Instability tongue for the mixed mode plotted in the acceleration–frequency ($\varLambda - \omega$) space. Experimental data are given by the symbols. Theoretical predictions are shown in solid lines for the subharmonic $[1,1]$ rocking mode and harmonic $[2,0]$ pumping mode for $\alpha =60.7^\circ$ and $V=4.6$ ml with pinned contact line. (b) Threshold acceleration $a/g$ against slanting angle $\chi$ at the resonance frequency $f=128$ Hz.

Comparing experimental results with theoretical predictions, we observe that the right-hand side of the tongue shows good agreement with the subharmonic $[1,1]$ mode, whereas the left-hand side has deviations. In addition, theory does not predict an overlap of the instability tongues for the subharmonic $[1,1]$ rocking and harmonic $[2,0]$ pumping modes for these experimental conditions. This suggests that the harmonic $[2,0]$ pumping mode is not a Faraday wave, but rather a meniscus wave excited by low-amplitude vibrations. The interactions between this harmonic meniscus wave with the subharmonic rocking mode may manifest themselves through the shape of the instability tongue. Our model does not account for such interactions, which could explain the discrepancy with experiment. The presence of the meniscus wave is known to affect the shape of the Faraday wave instability tongue, as has been reported for waves in a cylindrical geometry by Batson, Zoueshtiagh & Narayanan (2013, figure 14). Low-amplitude meniscus waves have also been reported in the theoretical model of Maksymov & Pototsky (2019) for pancake-shaped sessile drops. These prior observations support our conclusion that the harmonic $[2,0]$ mode is indeed a meniscus wave.

5. Concluding remarks

We have derived a mathematical model of the parametric oscillations of the sessile drop subject to slanted vibrations, which gives rise to an infinite system of coupled Mathieu equations that is solved by Floquet theory. The limiting case of horizontal ($\chi =90^\circ$) vibrations gives rise to a purely harmonic drop response, consistent with experimental observations (Sharp, Farmer & Kelly 2011; Khan & Eslamian 2019), whereas the more general case $\chi \neq 90^\circ$ produces parametric oscillations defined by subharmonic and harmonic instability tongues in the driving frequency–amplitude space. We show how the shape of the instability tongues depends upon the static contact angle $\alpha$, slanting angle $\chi$, and contact-line boundary conditions.

Experiments have been performed to validate theoretical predictions focusing on the $[1,1]$ rocking mode, the $[2,0]$ pumping mode, and the mixed mode that is the superposition of a subharmonic rocking mode and harmonic pumping mode, the latter of which is highly relevant to directional drop motion (Brunet et al. 2007; Noblin et al. 2009; Guo et al. 2020). The agreement between theoretical predictions and experimental observations is generally quite good, with the exception of the shape of the instability tongue for the mixed mode (cf. figure 15a). Here, the theoretical instability tongue for the harmonic $[2,0]$ pumping mode does not overlap with the instability tongue for the subharmonic $[1,1]$ rocking mode, as we both expect and observe experimentally for the mixing mode. We conjecture that this is because the harmonic pumping mode is not a true Faraday wave associated with a threshold instability, but rather an edge wave that can be excited at comparatively low driving acceleration. This has been observed experimentally for zonal modes in general (Vukasinovic, Smith & Glezer 2007a; Chang et al. 2015). Similar wave dynamics has been observed in Faraday waves in a cylindrical container, including axisymmetric harmonic edge waves excited at low driving amplitude (Shao et al. 2021b), and mode mixing of a harmonic axisymmetric wave with a subharmonic asymmetric waves at a single driving frequency (Shao et al. 2021c). Shao et al. (2021a) conduct a set of experiments that clearly demonstrate that the presence of a meniscus at the container sidewall is responsible for generating harmonic edge waves, which are suppressed whenever the interface is perfectly flat. This observation can help to explain the reason why sessile drop zonal modes exhibit a harmonic response, and is perhaps not surprising given that a sessile drop necessarily has a meniscus.

Dissipation gives rise to a non-trivial threshold acceleration, which we do not predict given that our analysis is focused on inviscid fluids with either pinned or free contact-line boundary conditions. In the future, extensions to our model could be made to account for sources of dissipation, including bulk viscous dissipation (Prosperetti 1980; Strani & Sabetta 1988; Bostwick & Steen 2013) and contact-line dissipation (Davis 1980; Bostwick & Steen 2015). To incorporate the former into a model would require one to resolve the rotational component to the velocity field, while the latter can be accommodated through the Hocking (1987) contact-line condition. Recent interest in complex fluids with a non-trivial rheology suggests model extensions to include viscoelastic effects (Tamim & Bostwick 2021), which are known to produce novel dynamics in classical Faraday wave systems (Müller & Zimmermann 1999).