# A model to predict the oscillation frequency for drops pinned on a vertical
planar surface

J. Sakakeeny1  C. Deshpande2  S. Deb2  J. L. Alvarado2,3  Y. Ling1
<EMAIL_ADDRESS>1Department of Mechanical Engineering, Baylor
University, Waco, TX 76643, USA 2Department of Mechanical Engineering, Texas
A&M University, College Station, TX 77843, USA 3Department of Engineering
Technology and Industrial Distribution, Texas A&M University, College Station,
TX 77843, USA

###### Abstract

Accurate prediction of the natural frequency for the lateral oscillation of a
liquid drop pinned on a vertical planar surface is important to many drop
applications. The natural oscillation frequency, normalized by the capillary
frequency, is mainly a function of the equilibrium contact angle and the Bond
number (Bo), when the contact lines remain pinned. Parametric numerical and
experimental studies have been performed to establish a comprehensive
understanding of oscillation dynamics. An inviscid model has been developed to
predict the oscillation frequency for wide ranges of Bo and contact angles.
The model reveals the scaling relation between the normalized frequency and
Bo, which is validated by the numerical simulation results. For a given
equilibrium contact angle, the lateral oscillation frequency decreases with
Bo, implying that resonance frequencies will be magnified if the drop
oscillations occur in a reduced gravity environment.

###### keywords:

Drops, contact lines

##  1 Introduction

Effective drop removal from a condenser surface is critical to highly-
efficient dropwise condensation, which is in turn important to many
applications such as condensation heat transfer (Yao et al., 2017) and water
harvesting (Dai et al., 2018). The recent advancement on superhydrophobic
surfaces (Boreyko & Chen, 2009b; Yao et al., 2017) have shown promising
enhancement of droplet mobility. Nevertheless, high mobility of drops on
superhydrophobic surfaces can be lost when the condensate nucleates within the
texture and the wetting falls into the Wenzel state. To restore the high drop
mobility, external forcing such as surface vibration is typically required to
excite the oscillation of the drops and to revert the wetting to the Cassie
state (Boreyko & Chen, 2009a; Yao et al., 2017). When the excitation frequency
matches the natural frequency, the oscillation amplitude of the drop is
maximized for a given energy input due to the resonance effect (Noblin et al.,
2004; Boreyko & Chen, 2009a; Yao et al., 2017). Therefore, it is highly
advantageous to predict the natural oscillation frequencies for drops on
surfaces of different wettability.

Similar to a free drop, the oscillation of a drop supported by a planar
surface can be characterized by spherical harmonic modes, with $n$ and $m$
representing the longitudinal and azimuthal wavenumbers, respectively (Strani
& Sabetta, 1984; Bostwick & Steen, 2014). Different from free drops, the two
$n=1$ modes, corresponding to the centroid translation in the directions
normal $(n=1,m=0)$ and tangential $(n=1,m=1)$ to the surface, trigger shape
oscillations for supported drops. Since the natural frequency and damping rate
increase with $n$ (Lamb, 1932), the $n=1$ modes typically dominate the
oscillations of supported drops (Strani & Sabetta, 1984; Noblin et al., 2004;
Sakakeeny & Ling, 2020, 2021). While the longitudinal mode $(n=1,m=0)$ is
axisymmetric, see figure 1(a), and has been extensively studied (Strani &
Sabetta, 1984, 1988; Noblin et al., 2004; Sakakeeny & Ling, 2020, 2021), less
attention has been paid to the lateral mode $(n=1,m=1)$, see figure 1(b),
which is the focus of the present study.

In the limit of zero-gravity, different inviscid theoretical models (Strani &
Sabetta, 1984; Celestini & Kofman, 2006; Bostwick & Steen, 2014) have been
developed to predict the frequency of different spherical-harmonic modes,
including the lateral mode of interest here. As the effect of gravity is
characterized by the Bond number (Bo), the aforementioned inviscid models are
strictly valid only in the limit of zero Bo, though the model predictions were
observed to agree reasonably well with experiment for drops with small Bo
(Chang et al., 2015). A simple capillary-gravity wave model has been suggested
by Noblin et al. (2004) to include the effect of gravity on the oscillation
frequency, but significant discrepancies with experiment was found for
surfaces with large contact angles (Yao et al., 2017). Numerical simulations
have been recently performed by Sakakeeny & Ling (2020, 2021) to investigate
the effect of Bond number on the axisymmetric zonal oscillation modes for
supported drops. Their results indicated that the drop oscillation frequency
varies significantly with Bo, even when the drop radius is smaller than the
capillary length ($\text{Bo}<1$). Furthermore, they also found that the
orientation of the supporting surface, with respect to the gravity
acceleration, also matters: while the drop oscillation frequency increases
with Bo for sessile drops, the frequency decreases with Bo for pendant drops.
To the knowlege of the authors, there is no former experimental or numerical
studies in the literature that systematically address the effects of Bo and
equilibrium contact angle on the oscillation of drops pinned on a vertical
planar surface. In the present study, a combined numerical, experimental and
theoretical analysis is performed for the lateral oscillation of a water drop
on a vertical planar surface. The goal of the study is to establish a model to
predict the oscillation frequency for drop at small but finite Bo. We focus
only on the drops with pinned contact lines, since they are the ones that
require oscillation excitation to restore high mobility.

Figure 1: Schematics of the (a) longitudinal ($n=1$, $m=0$) and (b) lateral
($n=1$, $m=1$) oscillation modes for a drop on a flat surface with pinned
contact lines. The red and green open circles indicate the longitudinal and
azimuthal nodal lines, while the red filled circle indicates the centroid of
the deformed drop. (c) Simulation results for the frequency spectra for the
maximum interfacial locations $l_{x}$ and $l_{y}$, for $\theta_{0}=$$$ and
$\text{Bo}=0$. The B-S model predictions (Bostwick & Steen, 2014) for modes
$(n,m)$ are shown (see vertical lines) for comparison.

##  2 Methods

###  2.1 Key parameters

The shape oscillation of a liquid drop is mainly controlled by surface
tension. The contact with the solid surface has a significant impact on the
oscillation dynamics. The natural frequency of the drop oscillation depends on
the equilibrium contact angle ($\theta_{0}$), contact line mobility, and Bo
(Bostwick & Steen, 2014; Sakakeeny & Ling, 2020, 2021). The key dimensionless
parameters and the ranges of values considered are summarized in Table 1. Due
to the small gas-to-liquid density and viscosity ratios, $r=\rho_{g}/\rho_{l}$
and $m=\mu_{g}/\mu_{l}$, the effect of the surrounding gas is generally small.
The effects of gravity and liquid viscosity, compared to surface tension, are
characterized by the Bond and Ohnesorge numbers,
$\text{Bo}=\rho_{l}gR_{d}^{2}/\sigma$ and
$\text{Oh}=\mu_{l}/\sqrt{\rho_{l}\sigma R_{d}}$, where $g$ and $\sigma$ are
gravity, and surface tension, while $R_{d}=(3V_{d}/4\pi)^{1/3}$ is the radius
based on the drop volume $V_{d}$. For millimeter-size water drops, Oh and Bo
are smaller than one. As a result, the oscillation frequency $\omega$
generally scales with the capillary frequency (Lamb, 1932),
$\omega_{c}=\sqrt{\sigma/(\rho_{l}R_{0}^{3})}$, where $R_{0}$ is the radius of
the equilibrium spherical-cap shape when gravity is absent, and
$R_{0}=R_{d}\left((2+cos\theta_{0})(1-cos\theta_{0})^{2}/4\right)^{-1/3}$. The
reason for using $R_{0}$, instead of $R_{d}$ is because $R_{0}$ better
represents the radius of curvature.

While $r$, $m$, and Oh are fixed in the simulation, the equilibrium contact
angle $\theta_{0}$ is varied from $None$ to $None$ to cover a wide range of
wettability, varying from hydrophilic to superhydrophobic surfaces. For a
given liquid, Bo can be varied by changing $V_{d}$ or $g$. In the present
study, we keep $V_{d}$ fixed and vary $g$, yielding $0\leq\text{Bo}\leq 0.78$.
Except for very small Bo, the simulation results also well represent cases for
drops with different $V_{d}$ under full gravity. Though Bo is smaller than
unity for all cases considered, namely $R_{0}$ is smaller than the capillary
length $\sqrt{\sigma/(\rho_{l}g)}$, it will be shown later that the effect of
gravity on the oscillation frequency can be significant.

r | $m$ | Oh | $\theta_{0}$ | Bo
---|---|---|---|---
$\rho_{g}/\rho_{l}$ | $\mu_{g}/\mu_{l}$ | $\mu_{l}/\sqrt{\rho_{l}\sigma R_{d}}$ | (°) | $\rho_{l}gR_{d}^{2}/\sigma$
$0.0012$ | $0.01$ | 0.00239 |  ${50}$ \- $150$ | 0 - 0.78
Table 1: Key parameters for simulations.

###  2.2 Simulation methods

The numerical simulations were conducted using the _Basilisk_ solver, in which
the Navier-Stokes equations are solved by a finite-volume approach. The sharp
interface is captured by the volume-of-fluid (VOF) method (Scardovelli &
Zaleski, 1999). The balanced-force method is used to discretize the surface
tension (Popinet, 2009). The height-function method is used to calculate the
interface curvature (Popinet, 2009) and to specify the contact angle at the
surface (Afkhami & Bussmann, 2009). Details of the numerical methods and
validation studies for the _Basilisk_ solver can be found in previous studies
(Sakakeeny & Ling, 2020, 2021) and the code website.

The computational domain is a cube with the edge length $L=8R_{d}$, see figure
2(a). The drop is symmetric with respect to the plane $z=0$ for the lateral
mode, so only half of the drop is simulated. The surface $x=0$ is in contact
with the drop and is taken as a no-slip wall. Other boundaries are slip walls.
Though modeling of mobile contact lines is challenging (Afkhami et al., 2018;
Snoeijer & Andreotti, 2013), here we consider only the asymptotic limit that
the contact lines are pinned. Physically this limit occurs when the hysteresis
effect is very strong and the oscillation amplitude is small. The pinning of
contact lines is achieved by imposing very small and large contact angles,
i.e., $\theta_{0}=$$$ and $None$, inside and outside of the contact line on
the surface $x=0$, see figure 2(b), to mimic the strong hysteresis effect. As
long as the contact angle $\theta$ remains bounded these two values when the
drop oscillates, the contact line cannot recede or advance. Validation for of
this treatment of contact lines can be found in our previous study (Sakakeeny
& Ling, 2021). The drop oscillations are initiated by specifying an initial
shape different from the equilibrium shape or by increasing the gravity for a
short period. The oscillation frequency will not be affected if the amplitudes
of excited oscillations are small.

Figure 2: (a) The simulation domain. (b) Closeup of the field of liquid volume
fraction $f$ on the contact plane ($x=0$), showing the pinned contact line
circle and the small and large $\theta_{0}$ specified inside and outside of
the contact line. (c) Closeup of the symmetric plane ($z=0$), showing the
adaptive mesh and VOF field.

An adaptive octree mesh is used for spatial discretization, see figure 2(c),
and the maximum refinement level is $\mathcal{L}=10$, corresponding to
$R_{d}/\Delta_{\min}=128$, where $\Delta_{\min}$ is the minimum cell size. The
adaptation thresholds for the liquid volume fraction ($f$) and the fluid
velocity components ($u,v,w$) are 0.0001 and 0.01, respectively. It is
verified that the refinement level and error thresholds are sufficient to
yield mesh-independent results. The simulations were performed on the cluster
_Kodiak_ using 144 CPU cores (Intel E5-2695 V4). Each simulation case takes
about 15 days to reach the time $t\omega_{c}\approx 110$. The simulation time
has been verified to be sufficiently long for accurate measurement of
frequencies from the spectra. To vary both $\theta_{0}$ and Bo, there are in
total about 30 cases simulated.

Figure 1(c) shows the frequency spectra of the maximum interfacial locations
in $x$ and $y$ directions, $l_{x}$ and $l_{y}$, for $\theta_{0}=$$$ and
$\text{Bo}=0$. The different oscillation modes are clearly identified: while
the spectrum for $l_{y}$ shows the $m=1$ modes, that for $l_{x}$ indicates the
$m=0$ modes. The frequencies predicted by the inviscid model of Bostwick &
Steen (2014) (B-S) for $(n,m$) modes are shown for comparison, see the
vertical lines. The simulation results generally agree very well with the
inviscid model, though small discrepancies were observed for the high-order
(3,0) and (3,1) modes. It can be seen that the lateral mode considered in the
present study exhibits the lowest frequency.

###  2.3 Experimental methods

Experiments have also been performed to measure the oscillation frequency of
the lateral mode. The drop oscillations were excited by the lateral vibration
of the supported substrate. Distilled water drops were placed on the vertical
substrate, which was fixed to a drop stand, see figure 3(a). The stand was
attached to a sound speaker (Infinity Reference 860 W, Infinity Inc.) with an
amplifier system (Russound P75-2 Channel Dual Source 75 W). The speaker
provided sine wave vibrations of different frequencies to the substrate using
a connected online tone generator. An accelerometer (352C04, PCB
Piezoelectronics) was used to measure vibration frequencies and acceleration
of the substrate. A National Instruments Data Acquisition System ((NI
USB-5132) was used to interface the accelerometer with the computer. A high-
speed camera (Photron SA3) with a halogen lamp (Fiber-Lite MI series, Dolan-
Jenner Inc.) was used to capture the drop deformation.

Figure 3: (a) Experiment setup. (b) Superimposition of drop edges on the
central plane from 100 snapshots for the polystyrene surface and drop volume
4.78 µL. The red dashed line indicates the equilibrium shape. The contact
angles at the top and bottom of the drop, $\theta_{t}$ and $\theta_{b}$, are
measured. The maximum contact angle hysteresis is defined as
$\Delta\theta=\max(\theta_{b})-\min(\theta_{t})$. (c) Variation of
$\Delta\theta$ as a function of the surface vibration frequency $f_{v}$. The
resonance angular frequency $\omega=2\pi f$ is identified as the peak of the
function.

The different substrates and drop volumes tested are summarized in Table 2.
Five different commercially available substrates with increasing $\theta_{0}$
were used. For each substrate, four different droplet volumes were tested. The
wax substrate was prepared in-house by firmly attaching wax paper on
polystyrene surface. All substrates, except for the wax one were cleaned with
isopropyl alcohol, distilled water, and air dried before testing. The ranges
of parameters considered in the simulation are wide enough to cover all the
experimental cases.

Substrate |  $\theta_{0}$ (°) |  $V_{d}$ (µL)
---|---|---
Al 6061 | 74.2 $\pm$ 1.80 | 3.30, 3.79, 4.29, 4.78
Cu 110 | 85.5 $\pm$ 1.30 | 3.30, 3.79, 4.29, 4.78
Polystyrene | 93.0 $\pm$ 1.06 | 3.30, 3.79, 4.29, 4.78
PTFE | 100.6 $\pm$ 0.70 | 3.30, 3.79, 4.29, 4.78
Wax | 105.0 $\pm$ 0.60 | 3.30, 3.79, 4.29, 4.78
Table 2: Experimental cases.

For a given acceleration of the substrate, the excited drop oscillation
amplitude is maximized when the substrate vibration frequency $\omega_{v}$
matches with the natural frequency $\omega$. The acceleration amplitude is
controlled to be small to avoid drop sliding on the surface. Different
$\omega_{v}$ were imposed, and the range of $\omega_{v}$ tested was estimated
based on the theoretical model of Celestini & Kofman (2006) (C-K). For each
$\omega_{v}$, 100 drop images were recorded. The drop contours on the central
plane were identified using the Canny method in MATLAB and superimposed as
shown in figure 3(b). The contact angles at the top and bottom of the drops
$\theta_{t}$ and $\theta_{b}$ are measured. The maximum contact angle
difference, $\Delta\theta=\max(\theta_{b})-\min(\theta_{t})$, is used to
characterize the excited drop oscillation amplitude. The resonance frequency
$\omega=2\pi f$ is then identified by the local maximum of $\Delta\theta$ as a
function of $\omega_{v}=2\pi f_{v}$, see figure 3(c). A small $f_{v}$
increment of 0.5 Hz were used near the local maximum for accurate measurement.

##  3 Results and discussion

###  3.1 Oscillation frequency

The simulation and experimental results for the oscillation frequency of the
lateral mode $\omega$ are shown in figure 4. The results of previous
experiments by C-K and Sharp et al. (2011), and theoretical models by C-K and
B-S are also included for comparison. It is observed that the experimental
data for different drop volumes collapse approximately, when $\omega$ is
normalized by $\omega_{c}$. This scaling relation $\omega\sim\omega_{c}$
indicates $\omega\sim V_{d}^{3/2}$, which is consistent with previous
observations for minor gravity effect (Noblin et al., 2004). The normalized
frequency $\omega/\omega_{c}$ for different $V_{d}$ still varies slightly due
to the effect of Bo. For the experiments of C-K
($0.04\lesssim\text{Bo}\lesssim 0.19$) and Sharp
($0.16\lesssim\text{Bo}\lesssim 0.66$), the surface is horizontal,
$\omega/\omega_{c}$ increases with Bo (Sakakeeny & Ling, 2021), while for the
present experiment ($0.12\lesssim\text{Bo}\lesssim 0.16$) that uses vertical
surfaces, $\omega/\omega_{c}$ decreases with Bo (as shown later in figure 7).
When Bo decreases, the experimental results for different surface orientations
approach to the simulation results for $\text{Bo}=0$.

The B-S model assumes Bo=0, so the equilibrium drop is a spherical cap. The
inviscid linear stability theory is used to predict the frequency of the
oscillation. The B-S predictions generally agree well with simulation results,
see figure 4. Small discrepancy are observed for $\theta_{0}\gg$$$, which may
be due to the simplified contact-line boundary conditions used in the B-S
model (Sakakeeny & Ling, 2020). The C-K model predicts the frequency by a
different approach. The shape oscillation of the drop is modeled as a harmonic
oscillator, so the frequency is $\omega=\sqrt{k/m_{d}}$, where $m_{d}$ is the
drop mass. The effective spring constant $k$ is related to surface tension,
which serves as a restoring force to bring the drop back to its equilibrium
state. The C-K model is more convenient to use since an explicit expression of
$\omega$ is given, but it significantly under-predicts the frequencies, as
shown in figure 4 and also in the original paper of C-K.

Figure 4: Present simulation, model, and experiment results for the frequency
of the lateral oscillation mode as a function of the equilibrium contact angle
$\theta_{0}$ for $\text{Bo}=0$, compared with previous models and experiments
by Celestini & Kofman (2006) (C-K), Sharp et al. (2011), and Bostwick & Steen
(2014) (B-S). The experiments of C-K and Sharp are for sessile drops on a
horizontal surface, for which the normalized frequency increases with Bo,
while the present experiments are for a vertical surface, where the frequency
decreases with Bo. The arrows indicates the variation of the normalized
frequency as Bo increases. When Bo decreases, the experimental results for
different surface orientations approach to the simulation and model results
for $\text{Bo}=0$.

###  3.2 An inviscid model to predict oscillation frequency

A new inviscid model is developed in this study. The goal is to achieve an
explicit expression to accurately predict the oscillation frequency.
Furthermore, the new model incorporates the effect of finite Bo, which is
ignored in B-S and C-K models. Taking the equilibrium state for Bo=0 as the
reference state, where the y-coordinate of the centroid and the drop surface
area are $y_{c}=y_{c0}=0$ and $S=S_{0}$, respectively. The surface energy of
the drop is $E_{s}=\sigma(S-S_{0})$. The surface area $S-S_{0}$ can be
expanded as a function of $y_{c}-y_{c0}$,

$\displaystyle\frac{S-S_{0}}{S_{0}}=\eta\left(\frac{y_{c}-y_{c0}}{R_{0}}\right)^{2}+\xi\left(\frac{y_{c}-y_{c0}}{R_{0}}\right)^{4}.$
(1)

Due to symmetry, only the even order terms remain and the terms higher than
fourth order are truncated, assuming small oscillation amplitudes. The
relation between $S$ and $y_{c}$ during oscillation is assumed to be similar
to that for a static equilibrium drop under different body forces (Celestini &
Kofman, 2006). As a result, the parameters $\eta$ and $\xi$ depend only on
$\theta_{0}$ and can be obtained from the solution for static drops. The
_Surface Evolver_ code was used to obtain the static solutions for different
body forces for a given $\theta_{0}$, which are then used to measure $\eta$
and $\xi$ following Eq. (1). The procedures are similar to our previous study
(Sakakeeny & Ling, 2020) and thus are not repeated here. The results for
$\eta$ and $\xi$ are shown in figure 5(a). For convenient use of the results,
fitted correlations are provided,

$\displaystyle\eta(\theta_{0})$
$\displaystyle=a_{0}\left[\exp\left(a_{1}(\cos\theta_{0}+1)\right)-1\right]\,,$
(2) $\displaystyle-\xi(\theta_{0})$
$\displaystyle=b_{0}\left[\exp\left(b_{1}(\cos\theta_{0}+1)\right)-1\right]\,,$
(3)

where the fitting coefficients $[a_{0},a_{1}]=[0.52.0.99]$ and
$[b_{0},b_{1}]=[0.20,2.11]$. In the asymptotic limit $\theta_{0}\to$$$, both
$\eta$ and $\xi$ approach zero, since in that limit the lateral mode does not
induce deformation. The static solutions also indicate that $\eta>0$ and
$\xi<0$ for all $\theta_{0}$.

Figure 5: Model parameters (a) $\eta$ and $\xi$ and (b) $\zeta$ as functions
of $\theta_{0}$. (c) $E_{k}$ as a function of $E_{kc}$ and (d) temporal
evolution of $E_{k}$ and the approximation $\zeta E_{kc}$ for $\theta_{0}=$$$
and $\text{Bo}=0$.

The gravitational potential energy is expressed as
$E_{g}=m_{d}g(y_{c}-y_{c0})$, and the overall potential energy is
$E_{p}=E_{s}+E_{g}$. The equilibrium state is located at the minimum of
$E_{p}$, namely $dE_{p}/dy_{c}=0$, which yields

$\displaystyle 2\eta(y_{c1}/R_{0})+4\xi(y_{c1}/R_{0})^{3}+\text{Bo}_{1}=0\,,$
(4)

where $\text{Bo}_{1}=mgR_{0}/\sigma
S_{0}=0.42(2+\cos\theta_{0})^{1/3}(1-\cos\theta_{0})^{-1/3}\text{Bo}$. The
cubic equation can be solved for the equilibrium centroid position $y_{c1}$,
near which $E_{p}$ can be expanded as

$\displaystyle
E_{p}=\frac{1}{2}\frac{d^{2}E_{p}}{dy_{c}^{2}}(y_{c}-y_{c1})^{2}=\sigma
S_{0}\eta_{1}\left(\frac{y_{c}-y_{c1}}{R_{0}}\right)^{2}\,,$ (5)

where $\eta_{1}=\eta+6\xi(y_{c1}/R_{0})^{2}$. Though the exact solution exists
for Eq. (4), we use the approximation $y_{c1}\approx-{\text{Bo}_{1}}/{2\eta}$
by neglecting the cubic term. Then we get the approximation

$\eta_{1}\approx\eta+\frac{3\xi}{2\eta^{2}}Bo_{1}^{2}\,.$ (6)

The kinetic energy of the drop fluid is
$E_{k}=\rho_{l}\int_{V_{d}}|\bm{u}|^{2}/2\ dV$, which is evaluated in the
simulation by the integral over the domain $E_{k}=\rho_{l}\int
f(u^{2}+v^{2}+w^{2})dV$, where $f$ is the liquid volume fraction, and $u$,
$v$, and $w$ are the velocity components. If the drop moves as a rigid body,
the kinetic energy would be $E_{kc}=m_{d}v_{c}^{2}/2$, where $v_{c}=dy_{c}/dt$
is the y-component of the centroid velocity. Due to the internal flow with
respect to the centroid induced by the shape oscillation, $E_{k}>E_{kc}$
(Sakakeeny & Ling, 2021). Furthermore, as shown in figure 5(c) that, $E_{k}$
is approximately a linear function of $E_{kc}$ when the drop oscillates. This
is due to the fact that the temporal evolutions of energy of the internal flow
and the magnitude of $v_{c}$ are in phase. As a result, the kinetic energy
correction factor, defined as $\zeta=E_{k}/E_{kc}$, is approximately a
constant in time and can be measured from the simulation results as shown in
figure 5(c). Then we can approximate

$\displaystyle E_{k}\approx\zeta E_{kc}=\zeta m_{d}v_{c}^{2}/2\,.$ (7)

The temporal evolutions of $E_{k}$ and $\zeta E_{kc}$ are plotted in figure
5(b) and the latter is affirmed to be generally a good approximation for all
time. Further tests show that $\zeta$ varies little over Bo and is only a
function of $\theta_{0}$, see figure 5(b). The fitted function for $\zeta$ is
given as

$\displaystyle\zeta(\theta_{0})$
$\displaystyle=c_{1}\mathrm{erf}\left[(\cos\theta_{0}+1)/c_{2}\right]+1$ (8)

where the fitting coefficients $[c_{1},c_{2}]=[0.14,0.93]$. In the asymptotic
limit $\theta_{0}\to$$$, $\zeta\to 1$ since there is no shape oscillation and
the lateral mode corresponds to rigid-body translation of the drop.

The total energy $E=E_{p}+E_{k}$ is constant in time, so $dE/dt=0$. With Eqs.
(5) and (7) we obtain

$\displaystyle k(y_{c}-y_{c1})+m_{d}\frac{d^{2}(y_{c}-y_{c1})}{dt^{2}}=0\,,$
(9)

where $k=2\sigma S_{0}\eta_{1}/(\zeta R_{0}^{2})$. The oscillation frequency
$\omega=\sqrt{k/m_{d}}$ can be expressed as

$\displaystyle\frac{\omega^{2}}{\omega_{c}^{2}}=\frac{12\left(\eta+\frac{3\xi}{2\eta^{2}}\text{Bo}_{1}^{2}\right)}{\zeta(2+\cos\theta_{0})(1-\cos\theta_{0})}\,.$
(10)

In the limit of $\text{Bo}=0$, Eq. (10) reduces to

$\displaystyle\frac{\omega^{2}}{\omega_{c}^{2}}=\frac{12\eta}{\zeta(2+\cos\theta_{0})(1-\cos\theta_{0})}\,.$
(11)

The results of Eq. (11) are plotted in figure 4, which are shown to agree
remarkably well with the simulation results. The C-K model for Bo=0 can be
written as

$\displaystyle\frac{\omega_{CK}^{2}}{\omega_{c}^{2}}=\frac{6\eta}{(2+\cos\theta_{0})(1-\cos\theta_{0})}\,.$
(12)

Comparing the Eqs. (11) and (12) for Bo=0, it is shown that the present model
yields better predictions due to the incorporation of contribution of the
kinetic energy of the internal flow through $\zeta$.

###  3.3 Model predictions for finite Bond numbers

Since numerical solutions are required for the model parameters $\eta$, $\xi$
and $\zeta$, the present model, similar to the C-K model, is semi-analytical.
Nevertheless, once these parameters are established, as shown in Eqs.
(2),(3),(8), the explicit expression Eq. (10) is ready to predict the natural
frequencies for the lateral oscillation of drops pinned on a vertical surface
with different $\theta_{0}$ and Bo.

The predictions of $\omega/\omega_{c}$ by the present, B-S, and CK models are
compared with the present experimental results for finite Bo in Table 3. The
results are also plotted in figure 6. Since the C-K model neglects the effect
of kinetic energy of the fluid motion within the drop, it significantly under-
predicts the frequency, up to 21%. In contrast, the B-S model, which neglects
the effect of Bo, over-predicts the frequency, up to 17%. As shown in figure
4, the present model yields very similar results as the B-S model for
$\text{Bo}=0$ in the range of $\theta_{0}$ considered in the present
experiments. By incorporating the contribution of Bo, the present model yields
more accurate predictions than the other two models for all cases. The
relative errors are less than 12% for all cases. Excluding the cases for
copper, the errors are less than 9%.

The present model over-predicts the frequency for all cases, compared to the
experimental results. There are a couple potential reasons that contribute to
the discrepancy between the model and experimental results. Similar to the B-S
and C-K models, the contact area is taken to be a circle with radius as
$R_{c}=R_{0}\sin{\theta_{0}}$ in the present model. This ideal contact area
corresponds to the spherical-cap shape of the drop when gravity is absent. The
contact area is then fixed during oscillation, since the contact line is
pinned. In the experiment, the contact area may be slightly smaller than ideal
contact area considered in the model, since gravity is present when the drop
is deposited on the surface. The smaller contact area will result in weaker
constraint on the drop and a lower oscillation frequency. Furthermore, though
in the experiment the oscillation amplitude is controlled to be very small,
the contact line may still not be fully pinned and may have oscillated in a
very small amplitude near the equilibrium position. In contrast, the contact
line is perfectly pinned in the model and simulation. The mobility of the
contact line will result in a reduction of the oscillation frequency for the
same $\theta_{0}$ and Bo, according to our previous study (Sakakeeny & Ling,
2021).

Substrate |  $V_{d}$ (µL) |  $\theta_{0}$ (°) | Bo | $(\frac{\omega}{\omega_{c}})_{\text{exp}}$ | $(\frac{\omega}{\omega_{c}})_{\text{mod}}$ |  $\varepsilon_{\text{mod}}$ (%) | $(\frac{\omega}{\omega_{c}})_{\text{BS}}$ |  $\varepsilon_{\text{BS}}$ (%) | $(\frac{\omega}{\omega_{c}})_{\text{CK}}$ |  $\varepsilon_{\text{CK}}$ (%)
---|---|---|---|---|---|---|---|---|---|---
Al 6061 | 3.30 | 74.60 | 0.116 | 2.726 | 2.878 | 5.3 | 2.924 | 7.3 | 2.175 | -20.2
| 3.79 | 74.50 | 0.127 | 2.586 | 2.773 | 6.7 | 2.827 | 9.3 | 2.097 | -18.9
| 4.29 | 72.15 | 0.138 | 2.745 | 3.013 | 8.9 | 3.054 | 11.3 | 2.282 | -16.9
| 4.78 | 73.60 | 0.149 | 2.699 | 2.926 | 7.8 | 2.976 | 10.2 | 2.218 | -17.8
Cu 110 | 3.30 | 85.85 | 0.116 | 2.163 | 2.333 | 7.3 | 2.396 | 10.8 | 1.760 | -18.7
| 3.79 | 84.80 | 0.127 | 2.134 | 2.376 | 10.2 | 2.441 | 14.4 | 1.794 | -15.9
| 4.29 | 87.10 | 0.138 | 2.022 | 2.276 | 11.2 | 2.342 | 15.8 | 1.720 | -14.9
| 4.78 | 84.10 | 0.149 | 2.115 | 2.401 | 11.9 | 2.472 | 16.9 | 1.817 | -14.1
Polystyrene | 3.30 | 91.55 | 0.116 | 1.959 | 2.106 | 7.0 | 2.156 | 10.0 | 1.586 | -19.0
| 3.79 | 92.20 | 0.127 | 1.890 | 2.080 | 9.1 | 2.130 | 12.7 | 1.568 | -17.0
| 4.29 | 93.50 | 0.138 | 1.867 | 2.029 | 8.0 | 2.082 | 11.5 | 1.531 | -18.0
| 4.78 | 91.75 | 0.149 | 1.903 | 2.091 | 9.0 | 2.148 | 12.9 | 1.580 | -16.9
PTFE | 3.30 | 100.10 | 0.116 | 1.702 | 1.812 | 6.1 | 1.834 | 7.8 | 1.361 | -20.2
| 3.79 | 99.85 | 0.127 | 1.702 | 1.817 | 6.4 | 1.846 | 8.5 | 1.367 | -19.7
| 4.29 | 101.20 | 0.138 | 1.683 | 1.772 | 5.0 | 1.802 | 7.0 | 1.334 | -20.7
| 4.78 | 101.05 | 0.149 | 1.686 | 1.773 | 4.9 | 1.807 | 7.2 | 1.338 | -20.6
Wax | 3.30 | 105.10 | 0.116 | 1.550 | 1.660 | 6.8 | 1.668 | 7.6 | 1.244 | -19.7
| 3.79 | 103.95 | 0.127 | 1.543 | 1.691 | 8.8 | 1.708 | 10.7 | 1.270 | -17.7
| 4.29 | 104.65 | 0.138 | 1.574 | 1.667 | 5.6 | 1.683 | 6.9 | 1.254 | -20.3
| 4.78 | 105.30 | 0.149 | 1.562 | 1.644 | 5.0 | 1.661 | 6.3 | 1.240 | -20.6
Table 3: Comparison between experimental results of the normalized frequency
$\omega/\omega_{c}$ and model predictions with the present, Bostwick & Steen
(2014) (B-S) and Celestini & Kofman (2006) (C-K) models.  Figure 6: Comparison
between experimental results of the normalized frequency $\omega/\omega_{c}$
and model predictions with the present, B-S and C-K models.

In the experiment, Bo is varied only in a very small range
($0.12<\text{Bo}<0.15$). Drops with larger Bo are very unstable and may slide
even when the surface vibration amplitude is very small. To better pin the
contact line and to enable experiment for larger Bo, a more sophisticated
surface treatment, such as that by Chang et al. (2015), is required, which
will be relegated to future work. In order to validate the model for higher
Bo, the model predictions are compared with the fully-resolved simulation
results.

The model and simulation results of the normalized frequency
$\omega/\omega_{c}$ as a function of Bo for $\theta_{0}=$$$ are shown in
figure 7(a). It is seen that $\omega/\omega_{c}$ decreases with Bo and
excellent agreement between the present model and the simulation results is
observed. It is conventionally considered that the effect of gravity on the
the supported drop oscillation is negligible if $\text{Bo}<1$, namely when the
drop radius is smaller than the capillary length (Chang et al., 2015).
However, it can be seen that the frequency decreases for more than 25% when Bo
increases from 0 to 0.62.

Figure 7: Comparison between the model and simulation results for oscillation
frequency as a function of Bo and $\theta_{0}$.

The variation of $\omega/\omega_{c}$ over Bo can be better depicted if
$\omega$ is normalized by its $\text{Bo}=0$ counterpart
$\omega_{\text{Bo}=0}$, From Eq. (10), it can be shown that

$\displaystyle{\omega^{2}}/{\omega_{\text{Bo}=0}^{2}}=1+({3\xi}/{2\eta^{3}})\text{Bo}_{1}^{2}\,.$
(13)

The equilibrium drop results show that $\xi<0$ for all $\theta_{0}$ (see
figure 5(a)), indicating that that the curvature of the $S$-$y_{c}$ curve (Eq.
(1)) decreases as $y_{c}$ increases. As a result, $\omega/\omega_{c}$
decreases with Bo, see figure 7(a). Furthermore, the leading correction term
for ${\omega^{2}/\omega_{\text{Bo}=0}^{2}}$ scales with $\text{Bo}^{2}$.
Therefore, for finite and small Bo, $\omega/\omega_{c}$ decreases
quadratically with Bo. It is difficult to vary Bo for a wide range in
experiment to examine the scaling relation, but that can be done in the
simulation. We have change $g$ from 0 to 0.9 m/s2, yielding
$0\leq\text{Bo}\leq 0.78$. The model and simulation results for
$\theta_{0}=$$,$$$, and $None$ are plotted in figure 7(b). Since drops with
large $\theta_{0}$ become unstable for large Bo, the range of Bo for
$\theta_{0}=$$$ is smaller than other cases. The model predictions (Eq. (13))
agree very well with the simulation results, validating the quadratic scaling
law revealed in the model.

In the present study, Bo is varied by changing the gravity $g$. Therefore, the
results imply that the natural frequency of the lateral oscillation mode will
increase if the drops oscillations occur in a reduced gravity environment.
This knowledge is important to the surface-vibration dropwise condensation
applications in the space station.

The present model can also be used to estimate the frequencies for the drops
of different volumes. For a given gravity, Bo increases with $V_{d}$. As a
result, $\omega/\omega_{c}$ decreases as $V_{d}$ increases. Since the model is
based on an inviscid assumption and the viscous effect on the oscillation
frequency (Lamb, 1932; Prosperetti, 1980) is not accounted, it is valid only
for $\text{Oh}\ll 1$. For low-viscosity liquids like water, Oh for millimeter-
drops is generally small. Yet, as Oh increases as the drop size decreases,
caution is required if the model is applied to very small drops.

## Acknowledgement

This work was supported by the startup fund at Baylor University and the
National Science Foundation (1942324). The Baylor High Performance and
Research Computing Services (HPRCS) have provided the computational resources
that have contributed to the research results reported. The authors also thank
Prof. Stéphane Zaleski for the helpful discussions on numerical model of
pinned contact lines.

## Declaration of Interests

The authors report no conflict of interest.

## References

  * Afkhami et al. (2018) Afkhami, S., Buongiorno, J., Guion, A., Popinet, S., Saade, Y., Scardovelli, R. & Zaleski, S. 2018 Transition in a numerical model of contact line dynamics and forced dewetting. J. Comput. Phys. 374, 1061–1093.
  * Afkhami & Bussmann (2009) Afkhami, S. & Bussmann, M. 2009 Height functions for applying contact angles to 3d VOF simulations. Int. J. Numer. Meth. Fluids 61 (8), 827–847.
  * Boreyko & Chen (2009a) Boreyko, J. B. & Chen, C.-H. 2009a Restoring superhydrophobicity of lotus leaves with vibration-induced dewetting. Phys. Rev. Lett. 103, 174502.
  * Boreyko & Chen (2009b) Boreyko, J. B. & Chen, C.-H. 2009b Self-propelled dropwise condensate on superhydrophobic surfaces. Phys. Rev. Lett. 103, 184501.
  * Bostwick & Steen (2014) Bostwick, J. B. & Steen, P. H. 2014 Dynamics of sessile drops. part 1. inviscid theory. J. Fluid Mech. 760, 5–38.
  * Celestini & Kofman (2006) Celestini, F. & Kofman, R. 2006 Vibration of submillimeter-size supported droplets. Phys. Rev. E 73, 041602\.
  * Chang et al. (2015) Chang, C.-T., Bostwick, J. B., Daniel, S. & Steen, P. H. 2015 Dynamics of sessile drops. part 2. experiment. J. Fluid Mech. 768, 442–467.
  * Dai et al. (2018) Dai, X., Sun, N., Nielsen, S. O., Stogin, B. B., Wang, J., Yang, S. & Wong, T.-S. 2018 Hydrophilic directional slippery rough surfaces for water harvesting. Sci. Adv. 4, eaaq0919.
  * Lamb (1932) Lamb, H. 1932 Hydrodynamics. Cambridge University Press.
  * Noblin et al. (2004) Noblin, X., Buguin, A. & Brochard-Wyart, F. 2004 Vibrated sessile drops: Transition between pinned and mobile contact line oscillations. Eur. Phys. J. E 14, 395–404.
  * Popinet (2009) Popinet, S. 2009 An accurate adaptive solver for surface-tension-driven interfacial flows. J. Comput. Phys. 228 (16), 5838–5866.
  * Prosperetti (1980) Prosperetti, A. 1980 Normal-mode analysis for the oscillations of a viscous-liquid drop in an immiscible liquid. J. Mecanique 19, 149–182.
  * Sakakeeny & Ling (2020) Sakakeeny, J. & Ling, Y. 2020 Natural oscillations of a sessile drop on flat surfaces with mobile contact lines. Phys. Rev. Fluids 5, 123604.
  * Sakakeeny & Ling (2021) Sakakeeny, J. & Ling, Y. 2021 Numerical study of natural oscillations of supported drops with free and pinned contact lines. Phys. Fluids 33, 062109.
  * Scardovelli & Zaleski (1999) Scardovelli, R. & Zaleski, S. 1999 Direct numerical simulation of free-surface and interfacial flow. Annu. Rev. Fluid Mech. 31, 567–603.
  * Sharp et al. (2011) Sharp, J. S., Farmer, D. J. & Kelly, J. 2011 Contact angle dependence of the resonant frequency of sessile water droplets. Langmuir 27, 9367–9371.
  * Snoeijer & Andreotti (2013) Snoeijer, J. H. & Andreotti, B. 2013 Moving contact lines: scales, regimes, and dynamical transitions. Annu. Rev. Fluid Mech. 45, 269–292.
  * Strani & Sabetta (1984) Strani, M. & Sabetta, F. 1984 Free vibrations of a drop in partial contact with a solid support. J. Fluid Mech. 141, 233–247.
  * Strani & Sabetta (1988) Strani, M. & Sabetta, F. 1988 Viscous oscillations of a supported drop in an immiscible fluid. J. Fluid Mech. 189, 397–421.
  * Yao et al. (2012) Yao, C. W., Garvin, T. P., Alvarado, J. L.., Jacobi, A. M., Jones, B. G. & Marsh, C. P. 2012 Droplet contact angle behavior on a hybrid surface with hydrophobic and hydrophilic properties. Appl. Phys. Lett. 101, 111605.
  * Yao et al. (2017) Yao, C.-W., Lai, C.-L., Alvarado, J. L., Zhou, J., Aung, K. T. & Mejia, J. E. 2017 Experimental study on effect of surface vibration on micro textured surfaces with hydrophobic and hydrophilic materials. Appl. Surf. Sci. 412, 45–51.