2.1. Modelling approach

Represented in figures 1(a) and 1(b) are simplified schematics of the mechanical configuration involved in recent experiments where strong attractive forces were generated with an SFL system (Colasante Reference Colasante2016; Weston-Dawkes et al. Reference Weston-Dawkes, Adibnazari, Hu, Everman, Gravish and Tolley2021; Omodia, Das & Mahadevan Reference Omodia, Das and Mahadevan2022). A flexible disk of radius $a$ and uniform thickness $t_d$ is located near a parallel wall. Glued coaxially to the opposite surface of the disk is a rigid, circular cylinder of radius $r_c < a$, to which a vibration actuator is affixed. The actuator exerts on the cylinder an oscillatory axial force that varies sinusoidally in time with frequency $\omega /(2{\rm \pi} )$ (where $\omega$ denotes the associated angular frequency), causing the disk to undergo dynamic bending. The nonlinear dynamics of the oscillatory airflow excited in the air layer separating the disk and the wall (Ramanarayanan et al. Reference Ramanarayanan, Coenen and Sánchez2022), the ‘squeeze film’, gives rise to a cycle-averaged pressure field in the film that provides a steady force attracting the disk toward the wall. This attractive force is able to support an applied pulling load of magnitude ${\mathcal {F}}_\ell$, steadied at the point of application by a vibration damper.

Application of Newton's second law for the rigid central assembly on which the load is exerted and the central portion of the disk to which it is affixed (see the inset in figure 1(a) for clarification) provides

(2.1)\begin{equation} {\mathcal{F}}_\ell + {\mathcal{F}}_m \cos(\omega t + \phi) + 2{\rm \pi} \int_0^{r_c} (p-p_a) r \, {\rm d} r + 2{\rm \pi} r_c \mathcal{V}_c = ({\rm \pi} r_c^2 \rho_d t_d + m_c) {\dfrac{\partial ^2h_c}{\partial t^2}} , \end{equation}

where $t$ and $r$ denote, respectively, the time and the radial distance from the axis of symmetry, ${\mathcal {F}}_m$ denotes the amplitude of the sinusoidal excitation force, $p(r,t)-p_a$ denotes the distributed gauge pressure acting on the disk surface (with $p_a$ denoting the ambient pressure), $\mathcal {V}_c$ represents the axial structural stress resultant within the disk at the critical radius $r_c$, $m_c$ quantifies the collective mass of the cylinder and other structures involved in the transfer of the vibrational force to the disk, $\rho _d$ denotes the uniform density of the disk and $h_c(t)$ represents the time-varying distance between the wall and the central portion of the disk. The phase shift $\phi$ included in the argument of the excitation force is introduced to facilitate the analysis, and will be discussed later. In writing (2.1), the gravitational acceleration is assumed to be negligibly small in comparison with that induced by the high-frequency excitation, and the weight of the central assembly is assumed to be small relative to the steady pulling load, both assumptions being consistent with the operating conditions used by Weston-Dawkes et al. (Reference Weston-Dawkes, Adibnazari, Hu, Everman, Gravish and Tolley2021).

Experimental data suggest that the central separation distance $h_c$ fluctuates periodically about a constant mean value $h_o$ if the imposed load ${\mathcal {F}}_\ell$ is steady. The mean distance $h_o$ grows as the load is quasi-statically increased, until a critical value is reached beyond which the attractive force becomes insufficient and the disk detaches abruptly from the wall, as depicted in figure 2(a) (Weston-Dawkes et al. Reference Weston-Dawkes, Adibnazari, Hu, Everman, Gravish and Tolley2021). It is of great interest to devise a solution of (2.1) that allows prediction of the ‘load capacity’ of an SFL system, the maximal pulling force $max[{\mathcal {F}}_\ell (h_o)]$ that can be supported without failure. In practical systems, known quantities are typically limited to the structural properties of the disk, the thermodynamic properties of the operating fluid and the amplitude ${\mathcal {F}}_m$ and frequency of the excitation force (Burroughs Reference Burroughs2021). Solution thus requires determining simultaneously the oscillating position $h_c(t)$ of the central assembly and the reactive aerodynamic and structural forces affecting its motion. The nonlinear fluid–structure dynamics involved complicates the solution, in that the periodic cylinder motion may exhibit multiple harmonics in addition to the excitation frequency $\omega$.

As shown in the following derivation, the problem simplifies in the case of small axial displacements of the disk relative to the mean separation distance $h_o$, allowing approximate solution of (2.1) with use of classical perturbation methods. Under the limit of small relative amplitudes, the cylinder may be assumed in the first approximation to oscillate sinusoidally along its axis with a given amplitude $b \ll h_o$ – as indicated in figure 1(b) – and an angular frequency equal to that of the excitation force, so that

(2.2)\begin{equation} h_c=h_o+b \cos(\omega t) .\end{equation}

Secondary frequencies arising from weak nonlinear dynamics will enter as asymptotic corrections to $h_c/h_o$ of order $(b/h_o)^2 \ll 1$ and smaller, to be discussed later. As indicated in (2.1), the cylinder displacement (2.2) is expected to exhibit a phase lag $\phi$ with respect to the excitation force, which is affected non-negligibly by the dissipative effects of viscous stresses in the squeeze film and material friction within the disk (Bettini Reference Bettini2016). In general, the evolving width $h(r,t)$ of the film may be expressed with use of the function

(2.3)\begin{equation} w(r,t)=h - h_o , \end{equation}

which denotes the position of the disk surface bounding the air film. Since the central portion of the disk is glued to the cylinder, it follows from (2.2) that, in the first approximation,

(2.4)\begin{equation} w = b \cos(\omega t) \quad{\rm for}\ 0 \le r \le r_c .\end{equation}

In pursuit of determining simultaneously the central oscillation amplitude $b$ and the reactive forces modulating its value, we introduce below a reduced theoretical description of the relevant fluid–structure dynamics, which couples the Reynolds equation governing the variation of air pressure $p$ in the squeeze film with the Kirchhoff–Love equation governing the dynamic displacement $w$ of the flexible annular portion of the disk $r_c \le r \le a$. We show that, when $b/h_o \ll 1$, the coupled equations can be solved with use of (2.4) as the driving kinematic boundary condition to provide analytical expressions for the reactive forces in terms of $b$. Substituting these expressions in (2.1) allows straightforward computation of the value of $b$ with relative errors of order $b/h_o \ll 1$. Solving the nonlinear system of equations that emerges at the following order in the asymptotic description for $b/h_o \ll 1$ allows analytical determination of the levitation force ${\mathcal {F}}_\ell$ with the same level of accuracy.

The asymptotic formulation will be shown below to provide promising agreement with experimental measurements of the pull-off curve ${\mathcal {F}}_\ell (h_o)$, the time-averaged overpressure in the squeeze film $\langle p \rangle (r) - p_a$ and the accompanying steady deformation of the disk $\langle w \rangle (r)$ (see figures 1(b) and 2 for clarification). Here, the angled brackets are used to denote the cycle-averaging operator

(2.5)\begin{equation} \langle * \rangle = \frac{\omega}{2{\rm \pi}} \int_{t}^{t+2{\rm \pi}/\omega} * \, {\rm d} t.\end{equation}

2.2. Description of the squeeze-film gas dynamics

In modelling the gas flow, it will be assumed that the squeeze film is slender, i.e. $h_o \ll \lambda$, with the characteristic flexural wavelength $\lambda$ of the disk being comparable to its radius $a$. Additionally, the viscous time scale $h_o^2/(\mu _a/\rho _a)$, where $\mu _a$ and $\rho _a$ denote respectively the ambient dynamic viscosity and density of the gas, is taken to be much smaller than the oscillation period $2{\rm \pi} /\omega$, so that the momentum balance is dominated by viscous forces and the energy balance by transverse heat conduction (since the Prandtl number of air is of order unity (Rohsenow, Hartnett & Cho Reference Rohsenow, Hartnett and Cho1998, p. 2.4)). In this approximation, the gas temperature across the film equals that of the solid boundaries, which is assumed here to equal the ambient gas temperature. In the associated limit of isothermal lubrication theory, the gas flow is determined, with relative errors of order $(h_o/a)^2 \!\ll\! 1$ and $h_o^2\omega /(\mu _a/\rho _a) \!\ll\! 1$, by the familiar equations (Langlois Reference Langlois1962; Ramanarayanan et al. Reference Ramanarayanan, Coenen and Sánchez2022)

(2.6a–d)\begin{equation} {\dfrac{\partial \rho}{\partial t}} + \frac{1}{r}{\dfrac{\partial (\rho r u)}{\partial r}} + {\dfrac{\partial (\rho v)}{\partial y}} = 0 ,\quad {\dfrac{\partial p}{\partial r}} = \mu_a {\dfrac{\partial ^2 u}{\partial y^2}} , \quad {\dfrac{\partial p}{\partial y}} = 0 \quad{\rm and}\quad \frac{p}{\rho} = \frac{p_a}{\rho_a} , \end{equation}

where $y$ is the normal distance to the wall, $u$ and $v$ denote respectively the radial and axial components of the flow velocity and $\rho$ and $p$ denote respectively the variable gas density and pressure, the subscript ‘$a$’ denoting ambient values found in the unperturbed surroundings.

The flow velocity in the film must satisfy the no-slip condition $u=0$ on both bounding surfaces $y=0$ and $y=h$ as well as the corresponding no-penetration conditions given by $v=0$ at $y=0$ and $v=\partial w/\partial t$ at $y=h$. Note that the condition $u(y=h)=0$ ignores negligibly small radial displacements of material points on the disk, of the order of $(b/h_o)^2(h_o/a)^2a \ll a$. The associated gas pressure must satisfy the conditions of regularity at the central axis of the film and relaxation at its edge

(2.7a,b)\begin{equation} \partial{p}/\partial{r} = 0 \ {\rm at} \ r = 0 \quad {\rm and} \quad p = p_a \ {\rm at} \ r = a , \end{equation}

respectively, the latter of which neglects variations of pressure existing beyond the edge, which are negligibly small in the lubrication limit, as explained in Appendix A.

The problem defined above affords straightforward reduction to the single equation

(2.8)\begin{equation} \frac{\partial }{\partial t}[(h_o+w) p]-\frac{1}{12 \mu_a r} \frac{\partial }{\partial r} \left[(h_o+w)^3 p r \frac{\partial p}{\partial r} \right]=0 , \end{equation}

the relevant Reynolds lubrication equation for compressible flow (Harrison Reference Harrison1913; Langlois Reference Langlois1962; Ramanarayanan et al. Reference Ramanarayanan, Coenen and Sánchez2022), which relates the disk deformation $w(r,t)=h-h_o$ with the coupled film pressure $p(r,t)$, the latter independent of $y$, as follows from (2.6c). The evolving pressure is to be determined by integration of (2.8) subject to the boundary conditions (2.7). Solution requires knowledge of the disk position $w$, whose value for $0 \le r \le r_c$ is determined in the first approximation by the presumed cylinder motion, as stated in (2.4), while its unknown value in the annular region $r_c < r \le a$ must be obtained from an analysis of the disk structural dynamics, as described below.

2.3. Description of the disk structural dynamics

In analysing the oscillations of the deformable portion of the disk extending for $r_c \le r \le a$, it will be assumed that the disk is made of a homogeneous, isotropic material, and that its axial deflections, comparable in magnitude to the central amplitude $b$, are much smaller than its uniform thickness $t_d$, which in turn is much smaller than its radius $a$, i.e. $b \ll t_d \ll a$. The excitation frequency is assumed to be sufficiently low for the characteristic flexural wavelength $\lambda$ of the disk to be comparable to the disk radius $a$, consistent with the slender-flow assumption drawn in § 2.2. Under these conditions the disk can be assumed to undergo pure bending (Timoshenko & Woinowsky-Krieger Reference Timoshenko and Woinowsky-Krieger1959, pp. 47–49), whence its dynamic deformations, influenced non-negligibly by the squeeze-film overpressure $p(r,t)-p_a$, can be described with use of the Kirchhoff–Love equation (Kelly Reference Kelly2013; Ducceschi Reference Ducceschi2014)

(2.9a,b)\begin{equation} \mathcal{D} \bar\nabla^4 w + \rho_d t_d {\dfrac{\partial ^2w}{\partial t^2}} = p-p_a , \quad{\rm where}\ \mathcal{D} = \frac{E_d t_d^3}{12(1-\nu_d^2)},\end{equation}

involving the axisymmetric Laplacian operator $\bar \nabla ^2 = [\partial ^2/\partial r^2 + (1/r)\partial /\partial r]$ and the flexural rigidity $\mathcal {D}$ of the disk, the latter defined in terms of the disk Young's modulus $E_d$, mass density $\rho _d$ and Poisson's ratio $\nu _d$. Correspondingly, the axial stress resultant $\mathcal {V}_c$ affecting the dynamics of the central assembly (2.1) can be expressed as (Yang Reference Yang2005)

(2.10)\begin{equation} {\mathcal{V}}_c ={-} {\mathcal{D}} [ \partial ( \bar{\nabla}^2 w )/\partial r ]_{r=r_c} .\end{equation}

Note that the displacement of the neutral plane of the disk has been substituted in (2.9a) and (2.10) with that of its film-adjacent surface $w$ (2.3), introducing negligible relative errors of order $(b/h_o)(t_d h_o/a^2) \ll 1$. Also, in-plane stresses induced due to fluid shear acting on this surface (Timoshenko & Woinowsky-Krieger Reference Timoshenko and Woinowsky-Krieger1959, pp. 378–380) are smaller than the expected overpressure by a factor of order $(b/h_o)(h_o/a)^2 \ll 1$ (Watanabe & Hara Reference Watanabe and Hara1997; Tulchinsky & Gat Reference Tulchinsky and Gat2016) and, hence, neglected when writing (2.9a).

Equation (2.9a) must be integrated for $r_c \le r \le a$ with the four boundary conditions

(2.11) \begin{equation} \begin{cases} w - b \cos(\omega t)={\dfrac{\partial w}{\partial r}}=0 & {\rm at}\ r = r_c \\ \left( {\dfrac{\partial ^2}{\partial r^2}} + \dfrac{\nu_d}{r}{\dfrac{\partial }{\partial r}} \right) w = {\dfrac{\partial }{\partial r}} ( \bar{\nabla}^2 w) = 0 & {\rm at}\ r = a \end{cases}, \end{equation}

which state that the annular portion of the disk is clamped at its inner edge $r=r_c$, where it follows the driving motion of the cylinder (2.4), and that neither bending moments nor axisymmetric axial stresses are supported at its free outer edge $r=a$ (see Timoshenko & Woinowsky-Krieger Reference Timoshenko and Woinowsky-Krieger1959, pp. 83–84 and Yang Reference Yang2005, pp. 840–841).

2.4. Dimensionless formulation and governing parameters

It is useful to reformulate the fluid–structure problem defined in §§ 2.2 and 2.3 in terms of appropriate dimensionless variables, beginning with the temporal and spatial variables $\tau =\omega t$ and $\xi =r/a$. The disk displacement is scaled with the central amplitude $b$ to give $W(\xi,\tau )=w/b$. We may then use (2.3) to write $h/h_o=1+\varepsilon W$, where

(2.12)\begin{equation} \varepsilon = {b}/{h_o} \end{equation}

defines the relative central oscillation amplitude. It follows from (2.4) that, in the first approximation, $W=\cos \tau$ for $\xi \le {\xi _c}$, where ${\xi _c} = r_c/a$ defines the dimensionless clamping radius. The characteristic scale of film overpressure, deduced readily from (2.6a,b) (see Appendix A), is used to define the variable $P=12 (p-p_a)/(\varepsilon \sigma p_a)$, where

(2.13)\begin{equation} \sigma = \frac{12 \mu_a \omega a^2}{p_a h_o^2} , \end{equation}

the classical dimensionless ‘squeeze number’ (Langlois Reference Langlois1962), quantifies the degree of gaseous compressibility in the film induced by viscous retardation of radial airflow.

The Reynolds lubrication equation (2.8) can now be written in the dimensionless form

(2.14)\begin{equation} \sigma {\dfrac{\partial }{\partial \tau}} [(1+\varepsilon W)P] - \frac{1}{\xi}{\dfrac{\partial }{\partial \xi}} \left[ (1+\varepsilon W)^3 \left(1+\frac{\varepsilon\sigma}{12}P \right) \xi {\dfrac{\partial P}{\partial \xi}} \right] + 12{\dfrac{\partial W}{\partial \tau}} = 0 ,\end{equation}

while the Kirchhoff–Love equation (2.9a) takes the form

(2.15)\begin{equation} \nabla^4 W + K^4 {\dfrac{\partial ^2W}{\partial \tau^2}} = {\mathcal{C}}^6 P ,\end{equation}

where $\nabla ^2 = [\partial ^2/\partial \xi ^2 + (1/\xi )\partial /\partial \xi ]$. As seen from (2.15), consideration of the disk dynamics introduces two additional dimensionless parameters defined by

(2.16a,b)\begin{equation} K^4 = a^4 \omega^2 \frac{12\rho_d (1-\nu_d^2)}{E_d t_d^2} \quad{\rm and}\quad {\mathcal{C}}^6 = a^6\frac{ 12 \mu_a \omega (1-\nu_d^2) }{ E_d (t_d h_o)^3 } , \end{equation}

where $K$ represents the well-known elastic wavenumber characterizing the bending of an isolated disk, i.e. one for which ${\mathcal {C}}=0$ (Warren Reference Warren1930; Yang Reference Yang2005; Ramanarayanan & Sánchez Reference Ramanarayanan and Sánchez2022). In the proceeding analysis, the ratio of the scale of film overpressure to that of disk inertia, ${\mathcal {C}}^6/K^4$, is found to be an important quantity that characterizes the ‘degree’ of flow–structure coupling. Interestingly, in the extreme case of strong coupling ${\mathcal {C}}^6/K^4 \gg 1$, the parameter $K$ is found to disappear from the formulation and is replaced by ${\mathcal {C}}$ as the relevant flexural wavenumber.

As follows from (2.7), (2.14) must be integrated with the boundary conditions

(2.17)\begin{equation} \partial{P}/\partial{\xi} = 0 \ {\rm at} \ \xi = 0 \quad {\rm and} \quad P = 0 \ {\rm at} \ \xi = 1 , \end{equation}

while, based on (2.11), (2.15) must be integrated with

(2.18) \begin{equation} \begin{cases} W - \cos \tau=\partial{W}/\partial{\xi} = 0 & {\rm at}\ \xi = {\xi_c} \\ \left( {\dfrac{\partial ^2}{\partial \xi^2}} + \dfrac{\nu_d}{\xi}{\dfrac{\partial }{\partial \xi}} \right) W = {\dfrac{\partial }{\partial \xi}} ( \nabla^2 W ) = 0 & {\rm at}\ \xi = 1 \end{cases}. \end{equation}

The problem defined above is to be addressed below for order-unity values of the distinguished dimensionless parameters $\sigma$, $K$ and ${\mathcal {C}}$, in order to develop a unifying theoretical description that reconciles the remarkable recent discoveries with the well-understood physics of conventional SFL systems which, together, span a wide parametric range (see, for example, Salbu Reference Salbu1964; Hatanaka et al. Reference Hatanaka, Yoshikazu, Kentaro, Sadayuki and Yoshiki1999; Sadayuki Reference Sadayuki2002; Zhao Reference Zhao2010; Yoshimoto et al. Reference Yoshimoto, Shou and Somaya2013; Colasante Reference Colasante2015; Weston-Dawkes et al. Reference Weston-Dawkes, Adibnazari, Hu, Everman, Gravish and Tolley2021). Analytical progress will be made using classical asymptotic methods that exploit the limit of small relative amplitudes $\varepsilon = b/h_o \ll 1$, as similar historical approaches have afforded clear insights regarding the simplified physics of SFL configurations that feature extreme values of $\sigma$, $K$ and/or ${\mathcal {C}}$ (Taylor & Saffman Reference Taylor and Saffman1957; Langlois Reference Langlois1962; Ramanarayanan et al. Reference Ramanarayanan, Coenen and Sánchez2022).

3.1. Solving for the presumed displacement amplitude $b$

The reactive aerodynamic and structural forces affecting the oscillation of the central assembly can now be expressed, with errors of order $\varepsilon \ll 1$, in the dimensionless forms

(3.15)\begin{equation} \frac{\displaystyle \int_0^{r_c} (p-p_a) r\,{\rm d} r}{\varepsilon p_a a^2} ={-}\mathrm{Re}\{ F_c \,\mathrm{e}^{\mathrm{i}\tau} \} , \quad{\rm where}\ F_c = {\xi_c} \left[ \frac{{\xi_c}}{2} + \frac{\bar A}{\varsigma} J_1(\varsigma{\xi_c}) \right] ,\end{equation}

and, based on (2.10),

(3.16)\begin{equation} \frac{ r_c \mathcal{V}_c }{ \mathcal{D} b/a^2 } ={-}\mathrm{Re}\{ V_c \,\mathrm{e}^{\mathrm{i}\tau} \} , \quad{\rm where}\ V_c = \xi_c \sum_{n=1}^3 \gamma_n^3 [ A_n J_1(\gamma_n\xi_c) + B_n Y_1(\gamma_n\xi_c) ] ,\end{equation}

respectively. Substituting into the governing equation of motion (2.1) the presumed central displacement (2.4), the central pressure force (3.15) and the structural impedance (3.16) reveals a linear relationship between the unknown displacement amplitude $b$ and the (generally) known excitation amplitude ${\mathcal {F}}_m$, given by

(3.17)\begin{equation} \frac{b}{{\mathcal{F}}_m} = | \mathcal{X} |^{{-}1} , \quad{\rm where}\ \mathcal{X} = 2{\rm \pi}\left( \frac{a^2}{h_o} p_a F_c + \frac{\mathcal{D}}{a^2} V_c \right) - \omega^2( m_c+ {\rm \pi}r_c^2 \rho_d t_d )\end{equation}

and an absolute value denotes the modulus of a complex quantity. The associated phase difference $\phi$ is given simply by $\phi = {\rm arg}(\mathcal {X})$. Note that this estimate for the value of $b$, accurate with relative errors of order $\varepsilon =b/h_o \ll 1$, does not require knowledge of the load ${\mathcal {F}}_\ell$ which enters in (2.1) as an asymptotic correction of order $\varepsilon$.

4.1. Time-averaged squeeze-film overpressure and disk deformation

The steady pressure distribution $\langle P_1 \rangle (\xi )$, independent of the steady disk deformation $\langle W_1 \rangle (\xi )$ under the present perturbative formulation, can be determined by straightforward integration of (4.1) to give

(4.3)\begin{equation} \langle P_1 \rangle(\xi) = 3\int_\xi^1 \left\langle W_0 {\dfrac{\partial P_0}{\partial \xi}} \right\rangle (x) \,{\rm d}\kern0.7pt x - \frac{\sigma}{24} \langle P_0^2 \rangle , \end{equation}

where $x$ serves as a dummy integration variable. Substitution of the expressions given in (3.3) provides

(4.4) \begin{align} \langle P_1\rangle = \begin{cases} \displaystyle \frac{1}{2} \mathrm{Re} \left\{ 3\left[ \bar{\varPi}({\xi_c}) - \bar{\varPi}(\xi) + \int_{\xi_c}^1 \varPi'^*({x}) \varOmega({x}) {\rm d}\kern 0.06em {x} \right] - \dfrac{\sigma}{24} \bar{\varPi}\bar{\varPi}^* \right\} & 0 \le \xi \le {\xi_c} \\ \displaystyle \frac{1}{2} \mathrm{Re} \left\{ 3\int_\xi^1 \varPi'^*({x}) \varOmega({x}) {\rm d}\kern 0.06em {x} - \dfrac{\sigma}{24} \varPi\varPi^* \right\} & {\xi_c} \le \xi \le 1 \end{cases} ,\end{align}

having made use of the identity $\langle \mathrm {Re}\{ \mathcal {G} \mathrm {e}^{\mathrm {i}\tau }\} \mathrm {Re}\{ \mathcal {H} \mathrm {e}^{\mathrm {i}\tau }\} \rangle = \mathrm {Re} \{ \mathcal {G}^* \mathcal {H} \} / 2$, where $\mathcal {G}$ and $\mathcal {H}$ are complex spatial functions and an asterisk denotes a complex conjugate.

Integration of (4.2) provides the accompanying time-averaged deformation of the disk

(4.5) \begin{align} \langle W_1 \rangle &= {\mathcal{C}}^6 \left( {\mathcal{I}}_1({\xi_c}) \ln\left( \frac{{\xi_c}}{\xi} \right) + \frac{1}{64} \int_{\xi_c}^\xi x^4 [ 4\ln(x)-5 ] \langle P_1 \rangle ' (x) \, {\rm d}\kern0.7pt x \right. \nonumber\\ & \quad - \left.\frac{\ln(\xi)}{16} \int_{\xi_c}^\xi x^4 \langle P_1 \rangle ' (x) \, {\rm d}\kern0.7pt x \ +\ \right|_{\xi_c}^\xi {\mathcal{I}}_2(\xi) \nonumber\\ & \quad + \left.\frac{ {\xi_c}^2 [ 1+2\ln(\xi/{\xi_c}) ] - \xi^2 }{ 2[ {\xi_c}^2-1 + 2/(1-\nu_d) ]}\left[ {\mathcal{I}}_1({\xi_c})+ \frac{1}{16} \int_{\xi_c}^1 x^4 \langle P_1 \rangle ' (x) \, {\rm d}\kern0.7pt x \right] \right) \end{align}

expressed here in terms of the radial pressure gradient

(4.6) \begin{equation} \langle P_1 \rangle ' = {\dfrac{{\rm d} \langle P_1 \rangle }{{\rm d} \xi}} = \begin{cases} -\dfrac{1}{2} \mathrm{Re} \left\{ \bar{\varPi}'^* \left( 3 + \dfrac{\sigma}{12}\bar{\varPi} \right) \right\} & 0 \le \xi \le {\xi_c} \\ -\dfrac{1}{2} \mathrm{Re} \left\{ \varPi'^* \left( 3\varOmega + \dfrac{\sigma}{12}\varPi \right) \right\} & {\xi_c} \le \xi \le 1 \end{cases}, \end{equation}

and the auxiliary integral functions

(4.7) \begin{equation} \left.\begin{array}{c@{}} \displaystyle {\mathcal{I}}_1(\xi) = \dfrac{\xi^2}{4} \left[ \dfrac{\xi^2}{4}\langle P_1 \rangle + \ln(\xi)\int_\xi^1 x^2 \langle P_1 \rangle ' (x) \, {\rm d}\kern0.7pt x - \int_\xi^1 x^2 \ln(x) \langle P_1 \rangle ' (x) \, {\rm d}\kern0.7pt x \right]\\ \displaystyle {\mathcal{I}}_2(\xi) = \dfrac{\xi^2}{16} \left[ \dfrac{\xi^2}{4}\langle P_1 \rangle + [2\ln(\xi) - 1] \int_\xi^1 x^2 \langle P_1 \rangle ' (x) \, {\rm d}\kern0.7pt x - 2 \int_\xi^1 x^2 \ln(x) \langle P_1 \rangle ' (x) \, {\rm d}\kern0.7pt x \right] \end{array}\right\}.\end{equation}

4.2. An analytical expression for the levitation force

Consistent with the simplified description of gas dynamics presented in § 2.2, the steady, attractive squeeze-film levitation force can be expressed as ${\mathcal {F}}_\ell = -2{\rm \pi} \int _0^a \langle p-p_a \rangle \, r {\rm d}r$, which follows from a cycle-averaged statement of Newton's second law for the periodically oscillating system cartooned in figure 1(a), comprising both the central assembly and the annular portion of the disk. (This definition for ${\mathcal {F}}_\ell$ may be obtained also by integrating the Kirchhoff–Love equation (2.9a) across the annular portion of the disk $r_c \le r \le a$, substituting the resulting expression for the stress resultant $\mathcal {V}_c$ (2.10) into the equation of motion (2.1) and taking the time average thereof.)

Based on the scalings introduced in § 2.4, the force ${\mathcal {F}}_\ell$ can be normalized with the supplied excitation amplitude ${\mathcal {F}}_m$ and the resulting relative oscillation amplitude $\varepsilon = b/h_o$ to provide a dimensionless expression for the operating efficiency of the system

(4.8a,b)\begin{equation} \frac{{\mathcal{F}}_\ell}{\varepsilon {\mathcal{F}}_m} = \frac{{\rm \pi} \mu_a \omega a^4}{h_o^3} \frac{\langle F_L \rangle}{|\mathcal{X}|} , \quad{\rm where}\ \langle F_L \rangle = \frac{12{\mathcal{F}}_\ell}{\varepsilon^2 \sigma p_a {\rm \pi}a^2} ={-}2 \int_0^1 \langle P_1 \rangle \xi \, {\rm d} \xi \end{equation}

and the complex constant $\mathcal {X}$ is defined in (3.17). Integration by parts in (4.8b) yields $\langle F_L \rangle = \int _0^1 \xi ^2 \langle P_1 \rangle ' \, {\rm d} \xi$, wherein substitution of the steady pressure gradient (4.6) provides

(4.9)\begin{align} \langle F_L \rangle &={-} \frac{1}{2}\mathrm{Re}\left\{ \frac{12}{\sigma} {\bar A}^* \left( \frac{2{\xi_c}}{\varsigma^*} [ 2J_1(\varsigma^*{\xi_c}) - \varsigma^*{\xi_c} J_0(\varsigma^*{\xi_c}) ] - \varsigma^*{\bar A} {\mathcal{I}}_{JJ}({\xi_c} ; \varsigma,\varsigma^*) \right) \right. \nonumber\\ & \quad + \left.\left.\sum_{n=1}^3 \sum_{m=1}^3 \frac{ \gamma_n^* ( K^4 - {\gamma_n^*}^4 )}{{\mathcal{C}}^6} \left[ 3 + \frac{\sigma( \gamma_m^4 - K^4 )}{12{\mathcal{C}}^6} \right] \right|_{\xi_c}^1 \mathcal{I}_L(\xi) \right\} , \end{align}

involving the auxiliary function

(4.10)\begin{align} \mathcal{I}_L(\xi) & = A_mA_n^*{\mathcal{I}}_{JJ}(\xi ; \gamma_m,\gamma_n^*) + A_mB_n^*{\mathcal{I}}_{JY}(\xi ; \gamma_m,\gamma_n^*) \nonumber\\ &\quad + B_mA_n^*{\mathcal{I}}_{YJ}(\xi ; \gamma_m,\gamma_n^*) + B_mB_n^*{\mathcal{I}}_{YY}(\xi ; \gamma_m,\gamma_n^*) , \end{align}

and the four operators ${\mathcal {I}}_{JJ},{\mathcal {I}}_{JY},{\mathcal {I}}_{YJ}$ and ${\mathcal {I}}_{YY}$, defined collectively by

(4.11)\begin{align} {\mathcal{I}}_{FG}(\xi ; \beta_{F},\beta_{G}) &= 2\beta_{G}\xi [ \beta_{G} F_0(\beta_{F}\xi)G_1(\beta_{G}\xi) - \beta_{F} F_1(\beta_{F}\xi)G_0(\beta_{G}\xi) ] / ( \beta_{F}^2-\beta_{G}^2 )^2 \nonumber\\ &\quad + \xi^2 [ \beta_{G} F_0(\beta_{F}\xi)G_0(\beta_{G}\xi) + \beta_{F} F_1(\beta_{F}\xi)G_1(\beta_{G}\xi) ] / (\beta_{F}^2-\beta_{G}^2 ) , \end{align}

where each $F_m$ and $G_m$ represents a Bessel function of the first or second kind ($J_m$ or $Y_m$) of $m$th order (Rosenheinrich Reference Rosenheinrich2019, p. 301).

6.1. System dynamics under one-way coupling

The general formulation derived above simplifies in the limit ${\mathcal {C}} \to 0$ to describe one-way-coupled systems, as detailed in Appendix B. It is worth reviewing here the widely studied behaviour of such systems – those for which the oscillations of the disk are unaffected by the squeeze-film overpressure and, hence, the flow–structure physics is governed purely by the parameters $\sigma$ and $K$. An increase in $\sigma \propto \omega a^2/h_o^{2}$, realized commonly by reducing the separation distance $h_o$, is known to result in a magnification of the typically repulsive squeeze-film force ${\mathcal {F}}_\ell < 0$ in a manner proportional to the inverse square of $h_o$ (Taylor & Saffman Reference Taylor and Saffman1957; Zhao Reference Zhao2010). This magnification is found to occur regardless of the flexibility of the oscillator (Li et al. Reference Li, Liu and Ding2014), the latter quantified by the elastic wavenumber $K$. A transition to attractive forces ${\mathcal {F}}_\ell >0$ is found to occur for systems with sufficiently low frequency $\omega$ (Sadayuki Reference Sadayuki2002) and/or oscillator size $a$ (Yoshimoto et al. Reference Yoshimoto, Shou and Somaya2013; Andrade et al. Reference Andrade, Ramos, Adamowski and Marzo2020). For these systems, attraction emerges only beyond a critical separation distance $h_o$, i.e. below a critical value of $\sigma$.

The range of operating conditions $\omega$ and $a$ for which attraction can be found, as well as the resulting attractive load capacity, are found to depend strongly on the waveform exhibited by the oscillating disk (Ramanarayanan & Sánchez Reference Ramanarayanan and Sánchez2022). In the case of rigid-body oscillation, i.e. $K=0$, repulsive forces exist for nearly all values of $\sigma$ (Zhao Reference Zhao2010). Attractive forces can be produced under a very limited range of operating conditions, characterization of which requires modelling the effects of fluid inertia that are neglected in the lubrication approximation drawn in the present study, but the corresponding load capacities are several thousandfold smaller and, thus, of limited practical utility (Yoshimoto et al. Reference Yoshimoto, Shou and Somaya2013; Andrade et al. Reference Andrade, Ramos, Adamowski and Marzo2020). (Analysis of such conditions requires considering separately the contributions of fluid compressibility and inertia with use of the relevant acoustic wavenumber $K_a = a \omega / \sqrt {p_a/\rho _a}$ and the Womersley number $\alpha = [h_o^2\omega /(\mu _a/\rho _a)]^{1/2}$, which compare with the oscillation period the time scales of lateral pressure equilibration and transverse viscous diffusion, respectively. These two parameters enter implicitly in the present lubrication limit $\alpha ^2 \ll 1$ through the single squeeze number $\sigma =12 K_a^2 /\alpha ^2$ (Melikhov et al. Reference Melikhov, Chivilikhin, Amosov and Jeanson2016; Ramanarayanan et al. Reference Ramanarayanan, Coenen and Sánchez2022).) In contrast, for critical values of $K$ that correspond to the natural frequencies of the disk (see Appendix B), the disk performs resonant standing-wave oscillations, such as those depicted in figure 3(a), and produces powerful attractive forces for a substantial range of values of $\sigma$ (Ramanarayanan & Sánchez Reference Ramanarayanan and Sánchez2022). Note that excitation of the disk at a non-natural frequency may result in different flexural waveforms that, instead, only magnify the repulsive squeeze-film force. This is often seen in practical systems where the disk is affixed to a Piezoelectric transducer of much greater mass (Da Silva Reference Da Silva1980; Li et al. Reference Li, Liu and Ding2014), yielding a mechanical assembly whose optimal, resonant frequencies presumably differ substantially from the native resonant frequencies of the disk.

Figure 3. Behaviour of an SFL system operating with an elastic wavenumber of $K=6.2$ and a squeeze number of $\sigma =20$. The flexural disk oscillations that occur under (a) a state of one-way coupling (${\mathcal {C}}=0$) are juxtaposed with those that occur for states of (b) moderate two-way coupling (${\mathcal {C}}=3$) and (c) strong two-way coupling (${\mathcal {C}}>5$). Shown below in (d–f) are the associated radial distributions of time-averaged pressure and disk displacement. Dots and circles in (d–f) represent numerical solutions computed for $\varepsilon =0.1$.

Motion of a centrally excited disk at its second resonant mode, i.e. $K \approx 6.2$, is visualized in figure 3(a) for a value of the squeeze number $\sigma =20$ that gives rise to an attractive force. As seen from the blue curve in figure 3(d), the distribution of steady pressure generated along the gas layer features two local minima whose radial locations correspond closely with those of the two nodal circles in the standing flexural wave, $r/a \approx 0.4$ and $0.8$. In the absence of aerodynamic forcing, the disk oscillates sinusoidally about an undeformed mean position, as indicated by the red line in figure 3(d). While this one-way-coupled description has been shown (Ramanarayanan & Sánchez Reference Ramanarayanan and Sánchez2022) to successfully reproduce the magnification of attractive load capacity found in flexural SFL systems (Colasante Reference Colasante2015; Weston-Dawkes et al. Reference Weston-Dawkes, Adibnazari, Hu, Everman, Gravish and Tolley2021), the predicted distribution of the causal time-averaged overpressure exhibits glaring disagreement with that measured in experiments, the latter exemplified in figure 2(b). Furthermore, setting ${\mathcal {C}}=0$ precludes characterization of the accompanying time-averaged deformation of the disk, a sample measurement of which is also reproduced in figure 2(b). Accounting for two-way fluid–structure coupling by considering non-zero values of ${\mathcal {C}}$ yields greatly improved agreement with these observations, as described below.

6.2. System dynamics under two-way coupling

As seen from figures 3(b) and 3(e), increasing the value of the coupling parameter to ${\mathcal {C}}=3$ yields but minor deviation from the one-way-coupled dynamics depicted in figures 3(a) and 3(d). The two nodal points are slightly disturbed and give way to localized regions of minimal amplitude, symptomatic of a damped oscillator (Geist & McLaughlin Reference Geist and McLaughlin1994). The corresponding steady pressure distribution $\langle P_1 \rangle (\xi )$ is visibly skewed but maintains a similar fundamental shape, and the disk oscillates about an approximately linearly deformed mean shape $\langle W_1 \rangle (\xi )$. Much more pronounced changes are observed when the strength of the aerodynamic forcing is further increased by setting ${\mathcal {C}}=5.25$, for which value the disk exhibits travelling-wave-type oscillations whose amplitude is suppressed severely beyond $r/a \approx 0.3$, as shown in figure 3(c). As seen from figure 3(f), the steady gauge pressure then assumes a negative value at the central axis $r=0$ and the steady disk deformation varies non-monotonically with radial distance, both in qualitative congruence with the experimental measurements reproduced in figure 2(b) (Weston-Dawkes et al. Reference Weston-Dawkes, Adibnazari, Hu, Everman, Gravish and Tolley2021). The disappearance of the first nodal region indicates that the collection of sand on the disk at $r/a \approx 0.3$ observed by Weston-Dawkes et al. (Reference Weston-Dawkes, Adibnazari, Hu, Everman, Gravish and Tolley2021) may not constitute a Chladni pattern, as the authors suggest, and owes instead to gravity-driven transport toward the valley formed by the mean disk deformation. Interestingly, in the limit of strong coupling ${\mathcal {C}} \gg 1$, considered separately in Appendix B, the flexural wavenumber of the disk is found to scale with ${\mathcal {C}}$, in place of the elastic wavenumber $K$ that applies for ${\mathcal {C}} \sim 1$, and the steady gauge pressure is seen to relax sharply from a negative value at the centre of the film toward zero across a small region $r/a \lesssim 5/{\mathcal {C}}$. For strongly coupled systems, reducing the degree of coupling ${\mathcal {C}} \propto 1/\sqrt {h_o}$ thus widens this central region, as observed by Weston-Dawkes et al. (Reference Weston-Dawkes, Adibnazari, Hu, Everman, Gravish and Tolley2021) during the pull-off process.

For each of the operating conditions depicted in figure 3, the analytical predictions of the time-averaged pressure distribution and disk deformation are confirmed with use of a finite-difference solution of the governing equations (2.14) and (2.15), results of which are visualized in the form of dots and circles in figures 3(d–f). As can be seen, for the moderately small relative amplitude $\varepsilon = 0.1$ used in the integration, the numerical solutions display satisfactory agreement with the asymptotic results. The governing equations were discretized with use of second-order central-space finite-difference approximations and marched in time asynchronously, the Reynolds equation with a forward-Euler scheme and the Kirchhoff–Love equation with use of a second-order central-difference scheme (Michael Reference Michael1963; Moin Reference Moin2010). The system was initialized with ambient pressure throughout the film $p(r)=p_a$, a uniform disk displacement $w(r)=b$ and zero disk velocity $\partial w/\partial t = 0$. Accurate characterization of the disk motion for ${\mathcal {C}} \neq 0$ required a sufficiently fine spatial discretization $\Delta \xi = \Delta r/a$ and, in turn, stable convergence to periodicity of the film pressure required a sufficiently small time step $\Delta \tau = \omega \Delta t$. For instance, the solution for ${\mathcal {C}}=3$ represented in figure 3(e) was obtained using $\Delta \xi = 0.01$ and $\Delta \tau = 6.3\times 10^{-4}$. The restriction on the time step for numerical stability generally loosened when increasing the degree of damping; for instance, the solution for ${\mathcal {C}}=5.25$ represented in figure 3(f) was obtained using $\Delta \tau = 8.4\times 10^{-4}$. Note finally that the precise time-averaged shape of the deforming disk $\langle W \rangle$ continued to fluctuate noticeably even after the associated levitation force had converged with a convergence ratio of less than $0.5\,\%$.

7.1. Reactive forces that modulate disk oscillation

Exemplified in figures 4(a) and 4(b) is the variation with the governing parameters of the dimensionless amplitude $|V_c|$ of the oscillatory reactive forces affecting a centrally excited disk, which is evaluated using (5.2). The black dot-dashed curve in each panel, corresponding to a state of one-way coupling ${\mathcal {C}}=0$, exhibits zeros for critical values of the wavenumber $K \approx 3.0, 6.2, 9.4, \ldots\, $, which represent the resonant modes of the flexurally oscillating disk (Warren Reference Warren1930). As shown in figure 4(a), increasing the strength ${\mathcal {C}}$ of the aerodynamic forcing tends to uniformize $|V_c|$ for a range of wavenumbers $K$, thereby eliminating the first few resonant modes. The affected range of $K$ and the number of suppressed modes grow as ${\mathcal {C}}$ is increased. For sufficiently large values of $K$ that satisfy $K^4 \sim {\mathcal {C}}^6$, the variation of $|V_c|$ gradually returns, in accordance with a restored three-term balance in the Kirchhoff–Love equation (2.15). Noting that $K^4/{\mathcal {C}}^6 \propto \omega$, this result agrees well with the observed decrease of the damping ratio of flexural squeeze-film bearings for increasing near-resonant operating frequencies $\omega$ (Pandey & Pratap Reference Pandey and Pratap2007). It must be emphasized here that, although $|V_c|$ exhibits local minima of decreasing value as $K$ grows, it never vanishes if ${\mathcal {C}} > 0$ – the dissipative effect of viscous aerodynamic damping precludes idealistic resonance. The curves in figure 4(b) reveal that an increase in the squeeze number $\sigma$ mitigates modal suppression. The reappearing minima of $|V_c|$ are displaced along the $K$-axis, with greater displacement occurring for larger values of ${\mathcal {C}}$, i.e. the shift of the optimal near-resonant frequency $\omega \propto K^2$ increases with the degree of damping, as anticipated for a harmonic oscillator (Zhang et al. Reference Zhang, Xu and Jiang2004; Bettini Reference Bettini2016).

Figure 4. Variation with the elastic wavenumber $K$ of the dimensionless amplitude of the forces modulating the oscillation of a centrally excited disk, as defined in (5.2). Curves are shown for two distinct values of the squeeze number $\sigma$ – (a) σ = 1 and (b) $\sigma = 500$ – and compared for states of one-way coupling (${\mathcal {C}}=0$), moderate two-way coupling (${\mathcal {C}}=3$) and strong two-way coupling (${\mathcal {C}}>5$).

7.2. Parametric dependences of the operating efficiency

Plotted in figure 5(a) is the variation with the squeeze number $\sigma$ and the coupling parameter ${\mathcal {C}}$ of the levitative efficiency provided by a centrally excited disk that is oscillating near its first natural frequency, i.e. $K = 3$. The efficiency is computed in the normalized form ${\mathcal {F}}_\ell / (\varepsilon {\mathcal {F}}_m)$, with use of (5.6). In practical systems, the value of the elastic wavenumber $K$ is determined by the structural properties of the oscillator and the selected operating frequency, while those of $\sigma \propto 1/h_o^{2}$ and ${\mathcal {C}} \propto 1/\sqrt {h_o}$ vary strongly with the evolving separation distance $h_o$ during the pull-off process (see § 2.1 for clarification). Thus, figure 5(a) may be interpreted qualitatively to represent the variation of the normalized efficiency with the inverse of the separation distance.

Figure 5. (a) Variation with the squeeze number $\sigma$ and coupling parameter ${\mathcal {C}}$ of the levitation efficiency (4.8a) for a selected value of the elastic wavenumber $K$. (b) Variation with $K$ of the local maximum in efficiency observed for ${\mathcal {C}} \sim 1$, juxtaposed with the universal value approached for ${\mathcal {C}} \gg 1$. Vertical lines in (b) represent resonant wavenumbers of the undamped disk.

Note first that the efficiency, which scales with ${\mathcal {C}}^6$ as per (5.6), vanishes in the limit of one-way coupling ${\mathcal {C}} = 0$ for all values of $\sigma$. On the other hand, when ${\mathcal {C}}=20$, the efficiency is finite and varies weakly with $\sigma$; further computations reveal that this behaviour is found to occur regardless of the value of $K$, consistent with the observation of modal suppression discussed in § 7.1. In fact, as detailed in Appendix B, in the limit ${\mathcal {C}} \gg 1$, the efficiency is found to converge to a constant value ${\approx }0.6495$, independent of both $\sigma$ and K.

It is apparent from figure 5(a) that the efficiency varies non-monotonically with both $\sigma$ and ${\mathcal {C}}$, and exhibits a local maximum for a critical combination of the two parameters. Represented in figure 5(b) is the value of this extremum, $max(\varepsilon ^{-1} {\mathcal {F}}_\ell / {\mathcal {F}}_m)$, for various wavenumbers $K$. The corresponding critical values of $\sigma$ and ${\mathcal {C}}$, omitted from the figure for clarity, increase nearly monotonically with $K$. For $K \gtrsim 8.3$, the local maximum exceeds the aforementioned universal value approached in the limit ${\mathcal {C}} \gg 1$, the latter indicated in the figure by a dashed horizontal line. The curvature of the function $[max(\varepsilon ^{-1} {\mathcal {F}}_\ell / {\mathcal {F}}_m)](K)$ diminishes as $K$ grows, subsequent peaks which appear for increasing critical values of $K$, marked by dots, being less accentuated. This suggests a favourable reduction in the sensitivity of the load capacity to the precise excitation frequency (Davis et al. Reference Davis, Gabai and Bucher2018) when operating near higher resonant modes. Note, in connection with the concluding discussion in § 7.1, that the critical values of $K$ corresponding to subsequent peaks are increasingly displaced from the respective ideal resonant wavenumbers, the latter represented in the figure by solid vertical lines.

7.3. Sample calculations for a practical system

While the squeeze number $\sigma$ and the coupling parameter ${\mathcal {C}}$ are formally independent parameters governing the fluid–structure dynamics, during the pull-off process they evolve together as dictated by their respective dependences on the growing separation distance, i.e. $\sigma \propto 1/h_o^{2}$ and ${\mathcal {C}} \propto 1/\sqrt {h_o}$. Practical determination of the load capacity $max[{\mathcal {F}}_\ell (h_o)]$ thus warrants a system-specific calculation, as exemplified below.

Shown in blue in figure 6(a) is the predicted pull-off curve ${\mathcal {F}}_\ell (h_o)$ for an SFL system operating in standard sea-level air under conditions similar to those used in the recent experiments by Weston-Dawkes et al. (Reference Weston-Dawkes, Adibnazari, Hu, Everman, Gravish and Tolley2021). A vibration actuator with excitation amplitude ${\mathcal {F}}_m = 74.5$ gf and frequency $\omega /(2{\rm \pi} ) = 200$ Hz is mounted to a cylinder of finite radius $r_c = 10$ mm. The assembly, having an effective mass $m_c = 10$ g, is affixed concentrically to a polyester plastic disk of radius $a = 30$ mm and thickness $t_d = 150\ \mathrm {\mu }$m, giving rise to an operating elastic wavenumber of $K \approx 3.75$. As $h_o$ is increased from approximately $15$ to $185\ \mathrm {\mu }$m, the squeeze number and coupling parameter vary respectively from $[\sigma,{\mathcal {C}}] \approx [10.77,12.43]$ to $[0.07,3.54]$. The red curve in the same figure shows the accompanying variation of the central oscillation amplitude $b(h_o)$.

Figure 6. Sample calculations for a disk with thickness $t_d = 150\ \mathrm {\mu }$m made of polyester, for which $E_d = 3.65$ GPa, $\rho _d = 1.38\ {\rm g}\ {\rm cm}^{-3}$ and $\nu _d = 0.48$ (Sasmita et al. Reference Sasmita, Candra, Judawisastra and Priambodo2019; MatWeb n.d.), in a gas where $p_a = 101325$ pa, $\mu _a = 1.81\times 10^{-5}\ {\rm kg}\ ({\rm m}\ {\rm s})^{-1}$ and $\rho _a = 1.225\ {\rm kg}\ {\rm m}^{-3}$. (a) Variation with the separation distance $h_o$ of the attractive force ${\mathcal {F}}_\ell$ and the central oscillation amplitude $b$, and (b) variation with the disk radius $a$ of the area-normalized load capacity, the latter compared with measurements by Weston-Dawkes et al. (Reference Weston-Dawkes, Adibnazari, Hu, Everman, Gravish and Tolley2021).

For small separation distances $h_o \approx 20\ \mathrm {\mu }$m, the amplitude $b$ vanishes due to augmented aerodynamic resistance (${\mathcal {C}} \gtrsim 10$) and the levitation force, which scales with the square of $b$ as seen from (4.8b), also decays. Note that this is in stark contrast with the transition to strong repulsion found in systems where effects of two-way coupling on the oscillator dynamics are negligible (Ramanarayanan & Sánchez Reference Ramanarayanan and Sánchez2022), i.e. those with heavier oscillators driven by powerful actuators (Hatanaka et al. Reference Hatanaka, Yoshikazu, Kentaro, Sadayuki and Yoshiki1999; Sadayuki Reference Sadayuki2002; Andrade et al. Reference Andrade, Ramos, Adamowski and Marzo2020). As the disk is pulled away from the wall, the amplitude rises due to a weakening of aerodynamic resistance and approaches asymptotically a value that is modulated only by structural bending forces within the oscillator. The force ${\mathcal {F}}_\ell$ also rises and reaches a maximal value $max({\mathcal {F}}_\ell ) \approx 56$ gf at a distance $h_o \approx 82\ \mathrm {\mu }$m, beyond which it decays gradually. In practical systems, displacing the disk beyond the critical distance requires the exertion of a pulling force that exceeds the load capacity $max({\mathcal {F}}_\ell )$ and yields, instead, an abrupt failure of levitation, as exemplified in figure 2(a). Points on the predicted pull-off curve for distances $h_o \gtrsim 82\ \mathrm {\mu }$m thus represent states of unstable equilibria that occur only in experiments involving artificial control of the separation distance (Sadayuki Reference Sadayuki2002; Yoshimoto et al. Reference Yoshimoto, Shou and Somaya2013; Andrade et al. Reference Andrade, Ramos, Adamowski and Marzo2020).

Furthermore, it is found that the predicted load capacity $max({\mathcal {F}}_\ell )$ varies strongly and non-monotonically with the disk radius $a$. As $a$ is increased from an initial value equal to $a=r_c=10$ mm, $max({\mathcal {F}}_\ell )$ is seen to grow and reach a value of approximately $81$ gf when $a \approx 18.5$ mm, before declining and approaching oscillatorily a value of approximately $45$ gf when $a \gtrsim 50$ mm. The fact that the variation of the load capacity subsides when $a$ is made sufficiently large is consistent with the aforementioned universal behaviour of the fluid–structure problem in the limit of strong coupling ${\mathcal {C}}^6 \gg K^4 \sim 1$. The value of the area-normalized load capacity $max({\mathcal {F}}_\ell )/({\rm \pi} a^2)$ thus collapses when $a$ is increased beyond ${\approx }20$ mm, as shown in figure 6(b). As seen from the included inset, this prediction agrees qualitatively with the observations of Weston-Dawkes et al. (Reference Weston-Dawkes, Adibnazari, Hu, Everman, Gravish and Tolley2021), who alluded to the possible benefit of increasing the number of oscillators rather than their size when upscaling an SFL system.

In principle, use of the asymptotic formulation introduces relative errors that scale with $b/h_o$, $h_o/a$ and $h_o^2/(\mu _a/\rho _a/\omega )$, the latter of which is of order unity for this sample SFL system when $h_o \gtrsim 75\ \mathrm {\mu }$m. Nevertheless, the scale of forces ${\mathcal {F}}_\ell$ seen in figure 6(a) in the rising portion of the pull-off curve $h_o \lesssim 75\ \mathrm {\mu }$m compares favourably with that reported in the experiment (Weston-Dawkes et al. Reference Weston-Dawkes, Adibnazari, Hu, Everman, Gravish and Tolley2021) for small distances of order $100\ \mathrm {\mu }$m, as seen from figure 2(a), in stark contrast with the theoretical predictions produced previously with use of a one-way-coupled flow–structure description (Ramanarayanan & Sánchez Reference Ramanarayanan and Sánchez2022). It can be anticipated that accurate description of the observed growth of ${\mathcal {F}}_\ell$ to larger values ${\approx }500$ gf at greater distances $h_o \approx 2500\ \mathrm {\mu }$m, seen in figure 2(a), as well as the associated area-normalized load capacities $max({\mathcal {F}}_\ell )/({\rm \pi} a^2) > 2$ seen in the inset of figure 6(b), requires generalizing the present formulation to account for effects of fluid inertia that are non-negligible when $h_o^2 \sim \mu _a/\rho _a/\omega$ (Ramanarayanan et al. Reference Ramanarayanan, Coenen and Sánchez2022).