2.1. The small-amplitude approximation

The assumption that the amplitude of the sheet remains small, or equivalently that $\varDelta \ll 1$, implies that the geometric relations (2.6) reduce to $\partial y_{sh}/\partial s\simeq \theta$ and $\partial x_{sh}/\partial s\simeq 1-\frac {1}{2}(\partial y_{sh}/\partial s)^2$. The nonlinearity in the derivative of $x_{sh}(s,t)$ is retained in the leading order of the theory so as to satisfy the constraint on the excess length (2.9a). Indeed, in the small-amplitude approximation, this constraint is given by

(2.11)\begin{equation} \varDelta=\frac{1}{2}\int_0^1\left(\frac{\partial y_{sh}}{\partial x}\right)^2{\rm d}\kern0.7pt x. \end{equation}

Here, we replace the arclength coordinate of the sheet with the Eulerian coordinate of the fluid, $s\simeq x$, according to our level of approximation. Correspondingly, the balance of forces and moments on the sheet (2.8) reduces to

(2.12)\begin{equation} \frac{\partial^2 y_{sh}}{\partial t^2}+\frac{\partial^4 y_{sh}}{\partial x^4}+F_x(t) \frac{\partial^2 y_{sh}}{\partial x^2}+\left[p_{u}(x,0,t)-p_{d}(x,0,t)\right]=0, \end{equation}

where the lateral compression, $F_x(t)$, is a function that depends solely of time. In addition, the pressure difference that the fluid exerts on the sheet, i.e. the last term in (2.12), is calculated at the sheet–fluid interface, $y=0$.

Thus far, we have approximated only the elastic part of the model. To further simplify the hydrodynamic part, we need to estimate the order of its corresponding fields. Given that the initial energy of the sheet scales linearly with the excess length, $E_{sh}(0)\propto \varDelta$, and that the total energy of the system is conserved, the energy of the fluid is, at most, proportional to $E_{f}\sim \varDelta$. Therefore, if we approximate the energy of the fluid by $E_{f}\sim |{\boldsymbol {v}}|^2 \ell$, where $|{\boldsymbol {v}}|$ is the typical velocity in the chamber and $\ell$ is the decay length of the disturbances in the flow, we obtain $|{\boldsymbol {v}}|\sim \sqrt {\varDelta /\ell }$. Furthermore, if we assume that the order of approximation of a derivative over the potential function, with respect to either a spatial dimension or time, does not change, then we can approximate Bernoulli's equation (2.2b), and the kinematic boundary conditions by

(2.13a)\begin{gather} \lambda p_i+\frac{\partial\phi_i}{\partial t}=c(t), \end{gather}

(2.13b)\begin{gather}\frac{\partial y_{sh}}{\partial t}=\left(\frac{\partial \phi_i}{\partial y}\right)_{y=0}. \end{gather}

These approximations will further be verified a posteriori in § 4, where we analyse the nonlinear dynamics of the system. In particular, we compare the results of this approximation with the numerical solution of the nonlinear model (2.2)–(2.9). Note that since the continuity equations (2.2a) are already linear in the potential functions, they remain unchanged in our approximated model.

This completes the reduction of our model to the small-amplitude limit. In summary, (2.2a) and (2.11)–(2.13), supplemented by the linearized form of the boundary conditions (2.5), (2.4), (2.9b) and (2.9c), form a closure and describe the coupled dynamics of the sheet and the fluid in the small-amplitude approximation. A comment is necessary regarding this simplified formulation. In accordance with the derivation in Appendix B, it can be shown that the reduced model emanates from the minimization of the action $\mathcal {S}=\int _0^{\rm T}\mathcal {L}\,{\rm d}t$, where

(2.14)\begin{align} \mathcal{L} &= \int_0^1\left[\frac{1}{2}\left(\frac{\partial y_{sh}}{\partial t}\right)^2-\frac{1}{2} \left(\frac{\partial^2y_{sh}}{\partial x^2}\right)^2+F_x(t) \left(\frac{1}{2}\left(\frac{\partial y_{sh}}{\partial x}\right)^2-\varDelta\right)\right. \nonumber\\ &\quad \left.\vphantom{\left(\frac{\partial y_{sh}}{\partial t}\right)^2} +\frac{1}{\lambda}\left[\phi_{d}(x,0,t)-\phi_{u}(x,0,t)\right] \frac{\partial y_{sh}}{\partial t}\right]{{\rm d}\kern0.7pt x} \nonumber\\ &\quad -\frac{1}{2\lambda}\int_{0}^{{L_y}/{2}} \int_0^1|\boldsymbol{\nabla}\phi_{u}|^2\,{\rm d}\kern0.7pt x\,{\rm d}y- \frac{1}{2\lambda}\int_{-{L_y}/{2}}^0 \int_0^1|\boldsymbol{\nabla}\phi_{d}|^2\,{\rm d}\kern0.7pt x\,{\rm d}y, \end{align}

with respect to the elastic fields $y_{sh}(x,t)$ and $F_x(t)$ and the hydrodynamic fields $\phi _{i}(x,y,t)$. In the next section, we employ a modal expansion of these fields and combine it with the Lagrangian formulation to derive a simplified set of equations that are dependent solely on time.

2.1.1. Modal expansion

The continuity equations (2.2a), and their corresponding boundary conditions on the fluid–chamber interfaces (2.4), are satisfied when the potential functions are given by

(2.15)\begin{align} \phi_i(x,y,t) &= a_0(t)\left(y\pm \frac{L_y}{2}\right)+\sum_{m=1}^{\infty}\cos ({\rm \pi} m x)\left[a_m(t)\exp\left({{\rm \pi} m\left(y\pm \frac{L_y}{2}\right)}\right) \right. \nonumber\\ &\quad \left.+\,c_m(t)\exp\left({-{\rm \pi} m \left(y\pm\frac{L_y}{2}\right)}\right)\right], \end{align}

where $a_m(t)$ ($m=0,1,2,\ldots$) and $c_m(t)$ ($m=1,2,3,\ldots$) are unknown time-dependent coefficients and the $\pm$ signs correspond to the solutions of the potential functions in the lower and the upper parts of the chamber, respectively. Similarly, we expand the solution of the sheet's height function in the following normal modes:

(2.16)\begin{equation} y_{sh}(x,t)=\sum_{n=1}^{\infty}A_n(t)\sin({\rm \pi} n x), \end{equation}

where the functions $\sin ({\rm \pi} n x)$ automatically satisfy the boundary conditions on the sheet edges (2.9b) and (2.9c), and $A_n(t)$ are as-yet unknown coefficients.

With these expansions, the solution to our problem reduces to finding the unknown coefficients, $a_m(t)$, $c_m(t)$ and $A_n(t)$, and the compression force $F_x(t)$, from the solution of the force balance equation (2.12), Bernoulli's equation (2.13a), the kinematic boundary conditions (2.13b) and the geometric constraint (2.11). Equation (2.16) involves infinite summation over the modes of the height function, but, in practice, we will truncate this series at $n=N$. A closed system of equations is then obtained when the coefficients of $a_m(t)$ and $c_m(t)$ are truncated at $N-1$.

However, instead of directly using these equations, we take a different – yet equivalent – approach, by utilizing the Lagrangian formulation (2.14). To this end, we follow the analysis in Appendix C and substitute the potential functions (2.15), and the height function (2.16), into the Lagrangian (2.14). We then integrate over the spatial coordinates. Thereafter, we minimize the Lagrangian with respect to $a_m(t)$ and $c_m(t)$ and express these coefficients in terms of $A_n(t)$. Substituting $a_m(t)$ and $c_m(t)$ back into the Lagrangian gives

(2.17)\begin{equation} \mathcal{L}\left[A_1,\ldots,A_N,F_x\right]= T_{nk}\frac{{\rm d}A_k}{{\rm d}t}\frac{{\rm d}A_n}{{\rm d}t}-V_{nk}A_k A_n+F_x(t)C_{nk}A_k A_n-\varDelta F_x(t), \end{equation}

where Einstein's summation rule is implied for repeated indices, and we define the following symmetric matrices:

(2.18a)\begin{gather} T_{nk}=\frac{1}{4}\delta_{nk}+\frac{L_y}{2\lambda}W(n,0)W(k,0)+\sum_{m=1}^{N-1} \frac{2}{{\rm \pi} m \lambda}\tanh\left(\frac{{\rm \pi} m L_y}{2}\right)W(k,m)W(n,m), \end{gather}

(2.18b)\begin{gather}V_{nk}=\frac{{\rm \pi}^4}{4}n^2k^2\delta_{nk},\quad C_{kn}=\frac{{\rm \pi}^2}{4}nk\delta_{nk}, \end{gather}

where $\delta _{nk}$ is the Kronecker delta, and $W(n,m)= ({n}/{{\rm \pi} })({(1-(-1)^{n+m})}/{(n^2-m^2)})$ for $n\neq m$ and zero otherwise.

Two comments are in order regarding this Lagrangian. First, since the matrix $T_{nk}$ is coupled to the kinetic terms in the Lagrangian, it takes on the role of a mass matrix in this formulation. This mass matrix has contributions from both the inertia of the sheet, i.e. the first term in $T_{nk}$, and the hydrodynamics of the fluid, i.e. the terms proportional to $1/\lambda$. The latter hydrodynamic terms are frequently referred to as added mass or virtual mass, since they describe an additional mass that the sheet appears to acquire when it accelerates in the fluid (Munk Reference Munk1924; Lighthill Reference Lighthill1960; Coene Reference Coene1992).

The second comment is related to the potential functions $\phi _i(x,y,t)$ that result from the minimization. Using (C2b), we substitute $c_m(t)=-a_m(t)$ in the potential functions (2.15) to obtain

(2.19)\begin{equation} \phi_i(x,y,t)=a_0(t)\left(y\pm \frac{L_y}{2}\right)+2\sum_{m=1}^{N-1}a_m(t)\cos({\rm \pi} m x)\sinh\left[{\rm \pi} m\left(y\pm \frac{L_y}{2}\right)\right]. \end{equation}

This solution implies that at $y=\pm L_y/2$ the velocity of the fluid is oriented only in the $y$ direction, i.e. $\partial \phi _i/\partial x(x,\pm L_y/2,t)=0$. It also implies that $\phi _i(x,\pm L_y/2,t)=0$, which, given (2.5) and (2.13a), yields a constant zero pressure along the inlet and the outlet walls of the chamber, $p_i(x,\pm L_y/2,t)=0$. We anticipate that these conditions will occur only when the disturbances that the sheet induces in the flow decay to zero. Therefore, the small-amplitude model holds strictly when $L_y\gg \ell$, where $\ell \simeq 1/({\rm \pi} m)$, for the smallest non-zero mode, is now explicitly identified as the decay length of the hydrodynamic disturbances.

Keeping in mind these limitations of the small-amplitude model, we go back to derive the equations for the coefficients $A_n(t)$. Given an initial volume difference, which, in turn, corresponds to an initial configuration of the sheet, i.e. a set of initial conditions for the coefficients $A_n(0)$, and keeping in mind that the system starts from rest, i.e. $({\rm d}A_n/{\rm d}t)(0)=0$, we can determine the dynamic evolution of the system from the minimization of (2.17) with respect to $A_n(t)$ and $F_x(t)$. This minimization yields $N+1$ algebraic differential equations that, in our matrix notation, read

(2.20a)\begin{gather} T_{nk}\frac{{\rm d}^2A_k}{{\rm d} t^2}+(V_{nk}-F_x(t)C_{nk}) A_k=0, \end{gather}

(2.20b)\begin{gather}C_{nk}A_kA_n=\varDelta. \end{gather}

Once $A_n(t)$ are determined from the solution of (2.20), the position of the sheet in time and in space is given by (2.16), and the hydrodynamic potentials and the pressure fields are determined from (2.15), (C2) and (2.13a). In the next section, we utilize this formulation to investigate the early time evolution. Then, in § 4 we use it to analyse the dynamics at later times.

3.1. Recap of the quasi-static solution

Following the analysis in the study of Oshri (Reference Oshri2021), the quasi-static evolution of the system is governed by two different branches of solutions, which we call ‘asymmetric’ and ‘symmetric’. Here, we recall the height functions in these branches.

On the one hand, when the initial volume difference is set as $0\leq v_{du}(0)\leq v_{du}^{cr}$, where $v_{du}^{cr}= ({2(3+{\rm \pi} ^2)}/{{\rm \pi} \sqrt {3(15+2{\rm \pi} ^2)}})\varDelta ^{1/2}$, the system is governed by the asymmetric branch. In this branch, the lateral compression is constant, $F_x(0)=4{\rm \pi} ^2$, and the height functions are given by

(3.1a)\begin{align} y_{sh}(x,0) &= \frac{p_{ud}(x,y,0)}{16{\rm \pi}^4}\left[2{\rm \pi}^2(1-x)x+1-\cos\left(2{\rm \pi} x\right)\right] \nonumber\\ &\quad +\frac{1}{\rm \pi}\sqrt{\varDelta-\frac{15+2{\rm \pi}^2}{768{\rm \pi}^6}p_{ud}(x,y,0)^2}\sin\left(2{\rm \pi} x\right), \end{align}

(3.1b)\begin{gather} p_{du}(x,y,0)=\frac{24{\rm \pi}^4}{3+{\rm \pi}^2}v_{du}(0), \end{gather}

where $p_{ud}(x,y,0)=p_{u}(x,y,0)-p_{d}(x,y,0)$ is the pressure difference between the upper and lower parts of the chamber. The potential energy of the sheet in this branch is given by

(3.2)\begin{equation} E_{as}=4{\rm \pi}^2\varDelta-\frac{6{\rm \pi}^4}{3+{\rm \pi}^2}v_{du}(0)^2. \end{equation}

Note that when $v_{du}(0)\rightarrow 0$ we know from (3.1b) that the pressure difference vanishes, $p_{ud}(x,y,0)\rightarrow 0$, and the elastic configuration converges to the second, asymmetric, mode of buckling, $y_{sh}(x,0)\rightarrow \sqrt {\varDelta /{\rm \pi} ^2}\sin (2{\rm \pi} x)$. The total energy of the system in this configuration is given by $E_{as}=4{\rm \pi} ^2\varDelta$. Note also that we considered solutions with an initial volume difference that is greater than zero. This is because the static solution has mirror symmetry around the $x$ axis. Solutions with $v_{du}(0)<0$ (and $p_{ud}(x,y,t)<0$) are obtained by a reflection of the height functions (3.1) around the horizontal axis.

On the other hand, when the volume difference is set as $v_{du}^{cr}\leq v_{du}(0)< \sqrt {2\varDelta /3}$, the system is governed by the symmetric branch. In this case, the inextensibility of the sheet implies an upper limit on the volume difference. In the case of a hinged sheet, this limit is given by $\sqrt {2\varDelta /3}$. The height functions in this branch are given by the parametric solution:

(3.3a)\begin{gather} y_{sh}(x,0)=\frac{p_{ud}(x,y,0)}{8u^2}(1-x)x+ \frac{p_{ud}(x,y,0)}{16u^4}\left[1-\frac{\cos\left[2u(x-1/2)\right]}{\cos u}\right], \end{gather}

(3.3b)\begin{gather}p_{ud}(x,y,0)={-}\frac{16\sqrt{6}u^{7/2}\cos u}{\sqrt{6u+4u(6+u^2)\cos^2 u-15\sin (2u)}}\varDelta^{1/2}, \end{gather}

(3.3c)\begin{gather}v_{du}(0)={-}\frac{2\sqrt{2}\left[u(3+u^2)-3\tan u\right] \cos u}{\sqrt{3}u^{3/2}\sqrt{6u+4u(6+u^2)\cos^2 u-15\sin(2u)}}\varDelta^{1/2}, \end{gather}

where $u=\sqrt {F_x(0)}/2$ is a function of the lateral compression. Given the initial volume difference, $v_{du}(0)$, and the excess length, $\varDelta$, we can determine the lateral compression, $F_x(0)$, from (3.3c), and then substitute this solution into (3.3a) and (3.3b) to obtain the height profile. When $p_{ud}(x,y,0)\rightarrow 0$, the height function converges to the first, symmetric, mode of buckling, which is given by $y_{sh}(x,0)\rightarrow \sqrt {4\varDelta /{\rm \pi} ^2}\sin ({\rm \pi} x)$ and $v_{du}(0)=8\varDelta ^{1/2}/{\rm \pi} ^2$. The total elastic energy of this shape is given by $E_{s}={\rm \pi} ^2\varDelta$. An example of the evolution of the sheet and the $p_{ud}(x,y,0)$–$v_{du}(0)$ relation in this static solution, where $0\leq v_{du}(0)< \sqrt {2\varDelta /3}$, is plotted in figure 2. In the following sections, we use these height functions (3.1a) and (3.3a) as the base solutions for our perturbative time-dependent expansion.

Figure 2. Evolution of the static solution in the small-amplitude approximation, where $\varDelta =0.01$. (a) The volume difference, $v_{du}(0)$, as a function of the pressure difference, $p_{ud}(x,y,0)$. In the asymmetric branch the $p_{ud}(x,y,0)$–$v_{du}(0)$ relation is given by (3.1b), while in the symmetric branch it is given by (3.3b) and (3.3c). The volume difference at the asymmetric-to-symmetric transition, $v_{du}(0)=v_{du}^{cr}$, is labelled by ③. The pressure difference in the chamber vanishes when the sheet accommodates either the second or the first mode of buckling, labels ① and ④. As the volume difference approaches its limiting value $v_{du}(0)\rightarrow (2\varDelta /3)^{1/2}$, the pressure difference diverges. (b) Evolution of the sheet's profile as the volume difference increases; see the corresponding labelled numbers in (a). Note, that despite the relatively large change in the pressure difference in the symmetric branch, the elastic configurations remain almost unchanged.

Before we proceed, we emphasize that while the sheet's configuration evolves continuously from the moment that we open the valve, the static pressure difference, given in (3.1b) and (3.3b), changes instantaneously at $t=0$. This is because we assumed an incompressible fluid, in which the speed of sound is infinite. Nonetheless, the asymmetric and the symmetric modes of buckling, obtained respectively from (3.1a) and (3.3a) in the limit $p_{ud}(x,y,0)\rightarrow 0$, are exceptions. These configurations remain in static equilibrium, which can nevertheless be unstable, when the valve is opened. For this reason, in the next section, we investigate the linear stability of the system around these two limiting initial states.

3.2. Linear stability

To derive the linear stability around the two limiting scenarios, i.e. the second and first buckling modes, we assume that the sheet's height function is given by the static solution, up to a small perturbation that grows exponentially with time. Correspondingly, we first perturb the amplitudes of the normal modes and the lateral compression around the base solutions, i.e. $A_n(t)=A_n(0)+\epsilon \bar {A}_n\,{\rm e}^{\sigma t}$ and $F_x(t)=F_x(0)+\epsilon \bar {F}_x\,{\rm e}^{\sigma t}$, where $\bar {A}_n$ and $\bar {F}_x$ are unknown constants, $\sigma$ is the growth rate and $\epsilon \ll 1$ is an arbitrarily small parameter. Then, we substitute these perturbed functions into the equations of motion (2.20), and expand them up to linear order in $\epsilon$. The leading order of this expansion, order $\epsilon ^0$, is given by

(3.4a)\begin{gather} V_{nk}A_k(0)-F_x(0)C_{nk}A_k(0)=0, \end{gather}

(3.4b)\begin{gather}C_{nk}A_k(0)A_n(0)=\varDelta, \end{gather}

and the subleading order, order $\epsilon$, is given by

(3.5a)\begin{gather} \left[\sigma^2 T_{nk}+V_{nk}-F_x(0)C_{nk}\right]\bar{A}_k-C_{nk}A_k(0)\bar{F}_x=0, \end{gather}

(3.5b)\begin{gather}C_{nk}A_k(0)\bar{A}_n=0. \end{gather}

The above equations in the subleading order always have the trivial solution $\bar {A}_n=0$ and $\bar {F}_x=0$, unless their corresponding determinant vanishes. This condition gives the growth rate, $\sigma$. Once $\sigma$ is determined, its corresponding eigenfunction is obtained from the solution of (3.5). The hydrodynamic fields related to this eigenfunction are determined from (2.15) and (C2).

3.2.1. Linear stability around the second mode of buckling

When the initial configuration of the sheet is given by the second mode of buckling, the base solution is derived from (3.1a) in the limit $p_{ud}(x,y,0)\rightarrow 0$. This solution reads $y_{sh}(x,0)=\sqrt {\varDelta /{\rm \pi} ^2}\sin (2{\rm \pi} x)$ and $F_x(0)=4{\rm \pi} ^2$. A projection of this configuration on the normal mode expansion (2.16) gives $A_n(0)=\sqrt {\varDelta /{\rm \pi} ^2}\delta _{2n}$. As expected, this initial state exactly satisfies the equilibrium equations at order $\epsilon ^0$ (3.4). At the next order, i.e. order $\epsilon$, we find that $\bar {F}_x=0$ and $\bar {A}_n=0$ for all the even perturbations, i.e. $n=2,4,6,\dots$. Consequently, (3.5) yields linear and homogeneous equations that involve only the odd perturbations. These equations always have the trivial solution $\bar {A}_n=0$, except when the corresponding determinant vanishes. A tractable solution to this condition, which also gives a good approximation to the highest growth rate, is obtained at the lowest order when $N=2$. This solution reads

(3.6)\begin{equation} \sigma=\frac{\sqrt{3}{\rm \pi}^2}{\sqrt{1+\dfrac{8 L_y}{{\rm \pi}^2\lambda}}}. \end{equation}

In figure 3, we plot this analytical approximation for the growth rate as a function of $\lambda$, and compare it with the numerical solution of (3.5) for the case where $N=8$. In addition, we compare this analytical solution with the growth rate obtained from the linearization of (2.2)–(2.9), i.e. where $\varDelta$ is assumed to be finite. See Appendix D for the details of this solution.

Figure 3. Log–log plot of the highest growth rate as a function of $\lambda$ where $L_y=2$. Symbols correspond to the linear stability analysis of (2.2)–(2.9), and solid and dashed lines correspond to the growth rates obtained from the small-amplitude approximation. When $\lambda \gg 1$, the growth rate converges to the constant $\sigma \simeq \sqrt {3}{\rm \pi} ^2$, whereas when $\lambda \ll 1$ the growth rate is given by $\sigma \simeq \sqrt {3{\rm \pi} ^6/8}(\lambda /L_y)^{1/2}$. Note that the differences between the solid ($N=2$) and dashed ($N=8$) lines are almost not visible in the figure. While the growth rate is independent of the excess length in the small-amplitude approximation, the more general solution (symbols) shows that the growth rate increases with $\varDelta$.

Equation (3.6) is one of the central results in this paper. Several comments are in order regarding this solution. Firstly, note that the highest growth rate is always real and positive, i.e. the second mode of buckling is always an unstable state of the system.

Secondly, while (3.6) depends on the parameter $\lambda /L_y=\tilde {\rho }_{sh} \tilde {h}/(\tilde {\rho }_{\ell }\tilde {L}_y)$, this result is not general but depends on the order of approximation. Had we solved (3.5) with $N\geq 3$, the two parameters $\lambda$ and $L_y$ would have appeared independently in the solution. Nonetheless, comparing the lowest order solution, i.e. the solution with $N=2$, with that of, say, $N=8$, we find that the growth rate remains almost unchanged (compare the solid and dashed lines in figure 3). For this reason, we conclude that (3.6) well describes the growth rate in the limit $\varDelta \ll 1$.

Thirdly, while the solution of the small-amplitude approximation is independent of the excess length, $\varDelta$, the more general solution for finite values of $\varDelta$ does depend on this parameter. See figure 3 for a comparison. In particular, for a fixed value of $\lambda$, the growth rate increases with an increase in $\varDelta$.

Fourthly, the analytical solution of the growth rate (3.6), exhibits two different regions as a function of $\lambda$. When $\lambda /L_y\gg 1$, the growth rate is a constant, $\sigma \simeq \sqrt {3}{\rm \pi} ^2$, that coincides with the growth rate of a sheet that is uncoupled from an external fluid. However, when $\lambda /L_y\ll 1$, the growth rate exhibits the scaling $\sigma \simeq \sqrt {3{\rm \pi} ^6/8}(\lambda /L_y)^{1/2}$. These asymptotic solutions define two limiting behaviours of the system. The former scenario represents a ‘solid-dominated’ region, in which the pressure difference exerted by the fluid on the sheet is negligible compared with the inertia of the sheet. The latter scenario, where $\sigma \propto (\lambda /L_y)^{1/2}$, represents the opposite limit of a ‘fluid-dominated’ region, in which the inertia of the sheet is negligible compared with the pressure difference exerted by the fluid on the sheet. Similarly, these two regions are also reflected in the added mass term, $8L_y/({\rm \pi} ^2\lambda )$, in the denominator of the growth rate. When $\lambda$ is large, the added mass approaches zero and the sheet's inertia is not affected by the fluid's motion. In contrast, when $\lambda$ is small, the effective mass of the sheet increases, resulting in slower dynamics. It is worth noting that since the dynamics of the system becomes very slow in the fluid-dominated region, we would expect some aspects of this solution to align with our earlier, quasi-static solution in the asymmetric branch, as shown in (3.1). This is because the quasi-static solution describes a slow spontaneous relaxation of the system, where the inertia of the sheet is negligible. This convergence to the quasi-static solution is further demonstrated in the following analysis.

Fifthly, note that if we fix $L_y$, the matrices in (3.5) and (2.18) become diagonal in the solid-dominated region. Therefore, up to small corrections of the order $1/\lambda$, $\bar {A}_1$ alone is excited at the instability. Indeed, in figure 4(a), in which we plot the eigenfunction for the case where $\lambda =100$ for both $N=2$ and $N=8$, we find that the two solutions are almost identical. In contrast, in the fluid-dominated region, $\lambda \ll 1$ and the matrices in (3.5) have non-zero off-diagonal terms. Therefore, all the odd modes become coupled and are excited at the instability. Nonetheless, our numerical investigation indicates that $\bar {A}_n/\bar {A}_1\ll 1$ for all $n\geq 3$. Therefore, while we would expect the leading-order solution, i.e. $N=2$, to approximate the eigenfunction well, it will not coincide with the higher-order solution. Indeed, in figure 4(b), we compare the eigenfunctions obtained from the lowest order and the high orders of approximations and find finite differences between them. We note that the eigenfunction in the fluid-dominated region emerges from the quasi-static configuration (3.1a), when we expand (3.1a) in powers of $p_{ud}(x,y,0)$, extract the linear order of this expansion and normalize it in accordance with our convention. The agreement between this quasi-static profile and the eigenfunction obtained from the linear stability analysis is shown in figure 4(b) (open squares).

Figure 4. Eigenfunctions of the sheet's height function in both the solid- and fluid-dominated regions. In both panels, $L_y=2$ and open blue circles correspond to the linear stability analysis at finite $\varDelta$, i.e. derived from (2.2)–(2.9). The eigenfunctions are normalized such that $[y_{sh}(1/2,t)-y_{sh}(1/2,0)]/{\rm e}^{\sigma t}=1$ (numerically we choose $\hat {y}_{sh}(1/2)=1$; see Appendix D). (a) In the solid-dominated region $(\lambda =100)$, only one mode is excited, i.e. the low ($N=2$) and the high ($N=8$) mode approximations coincide. (b) In the fluid-dominated region ($\lambda =0.1$), all odd modes are excited. Therefore, the lowest approximation, $N=2$, does not coincide with that obtained with higher modes, $N=8$. Nonetheless, since $\bar {A}_n/\bar {A}_1\ll 1$, the differences between the low and the high orders of the approximations are still small. Open squares represent the quasi-static approximation obtained from (3.1).

Sixthly, in figure 5(a), we plot the flow field obtained from the linear stability analysis at finite $\varDelta$. Note that the maximum velocity of the fluid is obtained at the sheet's centre, where the eigenfunction of the sheet's height is maximized (see figure 4). Unlike the growth rate, for which the small-amplitude approximation provided us with a good estimation already at $N=2$, the flow field converges at a slower pace. Convergence to the spatially dependent solution, presented in figure 5(b), is obtained only when higher modes, say $N\geq 3$, are included. This is because each coefficient $a_m^i(t)$ in the fluid's potential functions (2.15) depends on all the excited modes; see (C2). In addition, in figure 5(b), we plot the eigenfunctions of the pressure fields. Note that the sheet moves upwards towards the higher-pressure field. This is because the sheet's motion drives the flow, and the pressure drop in the fluid acts to slow down the onset of the elastic instability. We note that there are no qualitative differences in the flow fields and the pressure distributions between the solid- and fluid-dominated regions. For this reason, the plots in figure 5 refer only to the fluid-dominated region.

Figure 5. (a) The flow field and (b) the hydrodynamic pressures obtained from the linear stability analysis of (2.2)–(2.9) at the highest growth rate. In both panels, $\varDelta =0.01$, $\lambda =0.1$ and $L_y=2$. The eigenfunctions are normalized as indicated in figure 4. This gives $[\text {Low}, \text {High}]=[0.98,2.94]$ in the colour bar of the flow field and $[\text {Low}, \text {High}]=[-78,78]$ in the colour bar of the hydrodynamic pressures. The solid black line corresponds to the initial configuration of the sheet, i.e. the asymmetric second mode of buckling. In (a), arrows represent the streamlines and colours represent the relative magnitudes of the velocity.

Lastly, in addition to the growth rate, the flow field and the hydrodynamic pressures, another experimentally measurable quantity is the system's ‘compressibility’, i.e. the change in the volume difference relative to the change in the pressure difference. Since the pressure difference varies in space, we define the compressibility as $\beta \equiv {\rm d}v_{du}/{\rm d}\bar {p}_{ud}$, where $\bar {p}_{ud}(t)$ is the average pressure drop on the sheet. In the small-amplitude approximation the average pressure drop on the sheet is given by $\bar {p}_{ud}(t)=\int _0^1 [p_{u}(x,0,t)-p_{d}(x,0,t)]\,{{\rm d}\kern0.7pt x}$ and the volume difference in the chamber is given by $v_{du}(t)=2\int _0^1 y_{sh}(x,t)\,{{\rm d}\kern0.7pt x}$ (Oshri Reference Oshri2021). Keeping in mind that the base solution is asymmetric, and therefore does not contribute to the above integral, we find that the leading order ($N=2$) is $v_{du}(t)/\bar {A}_1=(4/{\rm \pi} )\,{\rm e}^{\sigma t}$. In addition, the average pressure difference at this order is given by $\bar {p}_{ud}(t)/\bar {A}_1=2L_y\sigma ^2 \,{\rm e}^{\sigma t}/({\rm \pi} \lambda )$. Consequently, at the onset of the instability, the compressibility is a constant, which is given by

(3.7)\begin{equation} \beta=\frac{2\lambda}{\sigma^2 L_y}. \end{equation}

In figure 6, we plot this result as a function of $\lambda$ and compare the lowest-order approximation with the solution at finite values of $\varDelta$. On the one hand, since the growth rate is constant in the solid-dominated region, we find that the compressibility scales linearly with $\lambda$, i.e. $\beta \rightarrow 2\lambda /(3{\rm \pi} ^4 L_y)$. On the other hand, since in the fluid-dominated region the growth rate scales as $(\lambda /L_y)^{1/2}$, the compressibility converges to the constant $\beta \rightarrow 16/(3{\rm \pi} ^6)\simeq 5.54\times 10^{-3}$. As expected, this result for the fluid-dominated region is very close to the compressibility obtained in the quasi-static solution (3.1b), which gives $\beta =(3+{\rm \pi} ^2)/(24{\rm \pi} ^4)\simeq 5.50\times 10^{-3}$.

Figure 6. Log–log plot of the compressibility as a function of $\lambda$ close to the onset of the instability. Symbols correspond to the compressibility at finite values of $\varDelta$ and the solid black line corresponds to the analytical solution obtained from the small-slope approximation. While in the solid-dominated region $\beta \propto \lambda$, in the fluid-dominated region $\beta$ converges to a constant.

3.2.2. Linear stability around the first mode of buckling

This section analyses the system's stability around the first mode of buckling. In contrast to the second mode, which is always unstable, the first mode represents the minimum of the elastic potential energy. Therefore, the first mode is expected to remain stable and to yield oscillatory motion under a small dynamical perturbation. Since a purely imaginary growth rate represents this oscillatory motion, we set $\sigma \equiv i\omega$, where $i=\sqrt {-1}$, and look for the lowest frequency of oscillation at the onset of the instability.

When the initial state of the sheet is given by the first mode of buckling, we have from (3.3) that the base solution is given by $y(x,0)=(2\sqrt {\varDelta }/{\rm \pi} )\sin ({\rm \pi} x)$ and $F_x(0)={\rm \pi} ^2$. A projection of this configuration on the normal mode expansion (2.16) gives $A_n(0)=(2\sqrt {\varDelta }/{\rm \pi} )\delta _{1n}$. This initial state, as expected, satisfies the leading order of our perturbative expansion (3.4). Substituting this leading order in (3.5) and solving for the unknown constants, we find that $\bar {F}_x=0$ and $\bar {A}_n=0$ for all the odd modes ($n=1,3,5,\dots$). Consequently, (3.5) yields linear and homogeneous equations that involve only the even modes of the height function. As in the previous case that we considered, a tractable solution to the vanishing determinant condition is obtained when we cut the normal mode expansion at the smallest value, $N=2$. (We note that when $L_y\gg 1$, the solution obtained from the two-mode approximation ($N=2$) is pre-empted by a different branch of solutions. The new branch emanates from a higher-order correction in the modal expansion, and its details are beyond the scope of the present study.) This gives

(3.8)\begin{equation} \omega= \frac{2\sqrt{3}{\rm \pi}^2}{\sqrt{1+\dfrac{128}{9{\rm \pi}^3}\dfrac{\tanh({\rm \pi} L_y/2)}{\lambda}}}. \end{equation}

In figure 7(a), we plot this solution and compare it with the solution of (3.5) for the case of $N=8$. Note that our solution, (3.8) with $N=2$, approximates well the small-amplitude limit, $\varDelta \ll 1$. Higher-order corrections, e.g. with $N=8$, almost do not alter this solution; compare the solid and dashed lines in figure 7(a). In addition, we compare (3.8) with the eigenvalues obtained from the linearization of (2.2)–(2.9), where the small-amplitude approximation is relaxed. We observe that as $\varDelta$ increases, the oscillation frequency decreases with increasing excess length. However, the overall trend of the dependence of $\omega$ on $\lambda$ remains consistent with the small-amplitude solution. We also note that the relatively large decrease in $\omega$ when $\lambda \gg 1$ is a result of our chosen hinged boundary conditions. In contrast, systems with clamped boundary conditions display a much milder dependence on $\varDelta$, as reported in previous studies (Neukirch et al. Reference Neukirch, Frelat, Goriely and Maurini2012; Pandey et al. Reference Pandey, Moulton, Vella and Holmes2014).

Figure 7. The lowest oscillation frequency and the sheet's eigenfunction obtained from the linear stability analysis around the first buckling mode. In both panels $L_y=2$. (a) Log–log plot of the oscillation frequency as a function of the parameter $\lambda$. Solid and dashed lines correspond to the solution of (3.8) when $N=2$ and $N=8$, respectively. All symbols correspond to the linear stability analysis of (2.2)–(2.9) at finite $\varDelta$. (b) The sheet's eigenfunctions in the solid- and fluid-dominated regions. All eigenfunctions are normalized such that their height at $x=1/4$ is equal to one (numerically we choose $\hat {y}_{sh}(1/4)=1$; see Appendix D). Although only one mode is excited when $\lambda \gg 1$ and infinite modes are excited when $\lambda \ll 1$, the eigenfunctions of the two regions are almost identical. Inset: an example of the sheet's oscillations around the base solution. Dashed lines correspond to an illustration of the dynamic oscillations, and the solid line to the base solution.

The frequency $\omega$ exemplifies the two limiting scenarios that we encountered for the growth rate in (3.6). On the one hand, in the solid-dominated region, where $\lambda \gg 1$, the frequency converges to the constant $\omega \simeq 2\sqrt {3}{\rm \pi} ^2$ that coincides with the frequency of oscillation of a sheet that is uncoupled from a fluid flow. On the other hand, in the fluid-dominated region, where $\lambda \ll 1$, we find the scaling $\omega \simeq \sqrt {27{\rm \pi} ^7/32}[\lambda /\tanh (L_y{\rm \pi} /2)]^{1/2}$, i.e. $\omega \propto \lambda ^{1/2}$. Alternatively, since the added mass is given by $({128}/{9{\rm \pi} ^3})({\tanh ({\rm \pi} L_y/2)}/{\lambda })$, when $\lambda \gg 1$ the sheet's inertia is almost unaffected by the fluid motion, but when $\lambda \ll 1$, the added mass increases and slows down the dynamics.

These two scenarios are also manifested in the eigenfunction of the sheet's height function. In the solid-dominated region, (3.5) becomes diagonal, up to corrections of the order of $1/\lambda$, and essentially only the second mode, i.e. $\bar {A}_2$, is excited at the instability, while in the fluid-dominated region (3.5) has non-zero off-diagonal terms, and all the even modes are excited. Nonetheless, our numerical investigation of the solution of (3.5) in the fluid-dominated region indicates that the ratio $\bar {A}_n/\bar {A}_2$ remains small for all $n$. Therefore, in both regions, deviations of the eigenfunction from the second mode are almost not visible in figure 7(b).

The oscillations of the sheet induce rotational flow in the chamber, whose magnitude decreases monotonically as $y$ approaches the far distant walls; see figure 8(a). Indeed, the solution for $\omega$ (3.8) becomes independent of $L_y$ as the vertical dimension of the chamber increases. Furthermore, the fluid's maximum velocity is obtained close to the centre point, $x=1/2$, where the sheet's velocity equals zero. This is because the sheet moves up and down around this centre line and drives a net flux across it. The velocity of this flux is maximized at the centre of the sheet. It is important to note that there is a slight asymmetry in the flow patterns observed in the upper and lower regions of the chamber, which may be attributed to the initial non-zero volume difference in the base solution.

Figure 8. (a) The flow fields and (b) the hydrodynamic pressure fields obtained from the linear stability analysis of (2.2)–(2.9) when $\varDelta =0.01$, $\lambda =0.1$ and $L_y=2$. The normalization of the eigenfunctions is as in figure 7(b). This normalization implies $[\text {High}, \text {Low}]=[5.3,26.5]$ in the flow fields and $[\text {High}, \text {Low}]=[-560,560]$ in the pressure fields. In (a), arrows represent the streamlines and colours represent the relative magnitudes of the velocity.

The pressure fields induced by the elastic oscillations are plotted in figure 8(b) and show an asymmetric profile in correlation with the eigenfunction of the sheet (figure 7b). Since only even modes are excited at the instability, the average pressure difference on the sheet, $\bar {p}_{ud}(t)$, vanishes. Therefore, in this case, there is no analogue to the compressibility calculated in § 3.2.1.

4.1. The elastohydrodynamic energetic interplay

Since the initial configuration of the sheet is close in shape to the second mode of buckling and the total energy of the system is conserved, the total energy at any $t>0$ is given by (3.2). In the limit $v_{du}(0)\ll v_{du}^{cr}$, this energy is approximated as $E_{as}\simeq 4{\rm \pi} ^2 \varDelta$, i.e. the energy of the second mode of buckling. In addition, since the first mode of buckling is the global minimizer of the elastic sheet's potential energy, the potential energy of the sheet cannot fall below $E_{s}={\rm \pi} ^2\varDelta$. Therefore, at most, our system can convert $\delta E=E_{as}-E_{s}\simeq 3{\rm \pi} ^2 \varDelta$ of the initial potential energy either into kinetic energy of the sheet or into hydrodynamic energy of the fluid. We remind the reader that the kinetic and potential energies of the sheet, $E_{sh}^{k}(t)$ and $E_{sh}^{p}(t)$, and the energy of the fluid, $E_{f}(t)$, are given by the first, second and third terms, respectively, on the right-hand side of (2.10). In the small-amplitude approximation, these energies reduce to $E_{sh}^{k}(t)=(\lim _{\lambda \rightarrow \infty } T_{kn})({{\rm d}A_k}/{{\rm d}t})({{\rm d}A_n}/{{\rm d}t})$, $E_{sh}^{p}(t)=V_{kn}A_k A_n$ and $E_{f}(t)=T_{kn}({{\rm d}A_k}/{{\rm d}t})({{\rm d}A_n}/{{\rm d}t})-E_{sh}^{k}$.

The typical evolution of the three components of the energy, i.e. the kinetic energy of the sheet, the potential energy of the sheet and the energy of the fluid, is plotted in figures 10(a) and 10(b). These plots are obtained from the solution of (2.20) with $N=2$ in the solid- and the fluid-dominated regions, $\lambda =10$ and $\lambda =0.1$ respectively, where the initial configuration is given by (3.1) with $v_{du}(0)=0.01v_{du}^{cr}$. In both cases, we find that, after some initial delay, the potential energy of the sheet drops from $E\simeq 4{\rm \pi} ^2 \varDelta$, i.e. the energy of the second mode of buckling (3.2), to $E\simeq {\rm \pi}^2\varDelta$, i.e. the energy of the first mode of buckling. The energy released in this process, $\delta E\simeq 3{\rm \pi} ^2\varDelta$, is converted to the kinetic energy of the sheet and the energy of the fluid, while the total energy remains fixed. In the solid-dominated region (figure 10a) the kinetic energy of the sheet becomes much larger than the energy of the fluid, while as $\lambda$ decreases, the opposite picture emerges, i.e. the energy of the fluid becomes much larger than the kinetic energy of the sheet (figure 10b). In both scenarios, shortly after the initial peak, a portion of the kinetic energy is converted back into potential energy of the sheet. As a result, the sheet tends to return to its configuration of the second buckling mode, thereby decreasing the height of the sheet's midpoint, as shown in figure 9(c).

Figure 10. The energetic interplay between the three components of the total energy in (a) the solid-dominated ($\lambda =10$) and (b) the fluid-dominated ($\lambda =0.1$) regions of the system. In these two panels, we solve (2.20) for the case where $N=2$, $\varDelta =0.01$, $L_y=2$, and $\sigma$ is given by (3.6). The initial configuration of the sheet is given by (3.1), where $v_{du}(0)=0.01v_{du}^{cr}$. After some initial delay, the potential energy of the sheet drops from $E_{sh}^{p}(t\ll 1)\simeq 4{\rm \pi} ^2 \varDelta$ to $E_{sh}^{p}({t}_p)\simeq {\rm \pi}^2\varDelta$. Open blue circles represent the potential energy obtained from the solution of the nonlinear equations (2.2)–(2.9). The energy released from the sheet is divided between the kinetic energy of the sheet, $E_{sh}^{k}(t)$, and the hydrodynamic energy of the fluid, $E_{f}(t)$, such that the total energy, $E$, remains constant. In the solid-dominated region, $E_{sh}^{k}(t_{p})\gg E_{f}(t_{p})$, while in the fluid-dominated region we find the opposite relation, $E_{sh}^{k}(t_{p})\ll E_{f}(t_{p})$. (c) The kinetic energy of the sheet and the fluid at $t=t_{p}$ as a function of $\lambda$, i.e. $E_i(t_{p})=E_{sh}^{k}(t_{p}), E_{f}(t_{p})$. A logarithmic scale is used on the $x$ axis. The dotted and the dashed-dotted lines correspond to our analytical solution from the two-mode approximation, and the colour symbols correspond to the numerical solution of (2.2)–(2.9), where the initial conditions are similar to those used in (a,b).

Our two-mode approximation allows us to quantify these findings and to further estimate the maximum values of $E_{sh}^{k}( t_{p})$ and $E_{f}(t_{p})$ as a function of $\lambda$, where $t_{p}$ denotes the time at which the potential energy of the sheet reaches its minimum value. Indeed, at $E_{sh}^{p}(t_{p})={\rm \pi} ^2\varDelta$, the sheet is close in shape to the first mode of buckling, i.e. $A_1(t_{p})\simeq 2\sqrt {\varDelta }/{\rm \pi}$. The constraint on the excess length (2.20b) then implies that $A_2(t_{p})\simeq 0$ and that $({\rm d}A_1/{\rm d}t)(t_{p})\simeq 0$. In addition, the derivative of the second mode at that moment, $({\rm d}A_2/{\rm d}t)(t_{p})$, is obtained from (4.1) when we substitute $E=4{\rm \pi} ^2\varDelta$ for the total energy. Taken together, these approximations yield the following maximum energies at $t=t_{p}$:

(4.2a)\begin{gather} \frac{E_{sh}^{k}(t_{p})}{3{\rm \pi}^2\varDelta}= \frac{1}{1+\dfrac{128 \tanh({\rm \pi} L_y/2)}{9{\rm \pi}^3\lambda}}, \end{gather}

(4.2b)\begin{gather}\frac{E_{f}(t_{p})}{3{\rm \pi}^2\varDelta}=\frac{1}{1+\dfrac{9{\rm \pi}^3\lambda}{128 \tanh({\rm \pi} L_y/2)}}. \end{gather}

In figure 10(c), we compare these maximum energies with the numerical solution of (2.2)–(2.9). Overall, we find a good fit between the analytical approximation and the numerical solution over the entire range of $\lambda$.

Using (4.2), we find that in the solid-dominated region, $\lambda \gg 1$, most of the initial energy is converted into the kinetic energy of the sheet, i.e. $E_{sh}^{k}(t_{p})\simeq 3{\rm \pi} ^2 \varDelta -\delta _1$ and $E_{f}(t_{p})\simeq \delta _1$, where $\delta _1=128\varDelta \tanh ({\rm \pi} L_y/2)/(3{\rm \pi} \lambda )$, while in the fluid-dominated region, $\lambda \ll 1$, most of the energy is converted into the energy of the fluid, i.e. $E_{sh}^{k}(t_{p})\simeq \delta _2$ and $E_{f}(t_{p})\simeq 3{\rm \pi} ^2\varDelta -\delta _2$, where $\delta _2=27{\rm \pi} ^5\varDelta \lambda /[128\tanh ({\rm \pi} L_y/2)]$. Furthermore, if we estimate the average instantaneous velocity of the fluid, $\bar {v}(t_{p})$, by using $E_{f}(t_{p})\propto \bar {v}(t_{p})^2 \ell /\lambda$, we find that in the solid-dominated region the average velocity is proportional to $\bar {v}(t_{p})\propto (\varDelta /\ell )^{1/2}$, while in the fluid-dominated region it is proportional to $\bar {v}(t_{p})\propto (\lambda \varDelta /\ell )^{1/2}$; namely, while $E_{f}(t_{p},\lambda \ll 1)\gg E_{f}(t_{p},\lambda \gg 1)$, their corresponding velocities present the opposite relationship, i.e. $\bar {v}(t_{p},\lambda \ll 1)\ll \bar {v}(t_{p},\lambda \gg 1)$. This is because the momentum of the fluid $p_{f}\propto \bar {v}(t_{p})/\lambda$ – but not its velocity – increases in the fluid-dominated region.

4.2. The peak time

In the previous section, we demonstrated through energetic considerations that the sheet tends to release most of its stored potential energy. In this section, we investigate the time it takes for the system to release this energy. Given an initial configuration of the sheet that is close in shape to the second mode of buckling, i.e. (3.1) with $v_{du}(0)\ll v_{du}^{cr}$, we aim to find the time $t=t_{p}$ at which the potential energy of the sheet first drops to the minimum value, $E_{sh}^{p}(t_{p})\simeq {\rm \pi}^2 \varDelta$. To do so, we first use (2.20b) and (3.6) to eliminate $A_2(t)$ and $\lambda$ in favour of $A_1(t)$ and $\sigma$, respectively, in (4.1). Then, we substitute the energy of the initial configuration (3.2) into (4.1) and integrate it between $t\in [0,t_{p}]$. This gives

(4.3)\begin{equation} \sigma t_{p}=\int_{A_1(0)}^{A_1(t_{p})} \frac{\sqrt{432{\rm \pi}^3\varDelta+\left[9{\rm \pi}\left({-}12{\rm \pi}^4+\sigma^2\right)+ \dfrac{16\left(3{\rm \pi}^4-\sigma^2\right)}{L_y\coth({\rm \pi} L_y/2)}\right]A_1^2}}{6\sqrt{3}{\rm \pi}^{3/2}\sqrt{\left(A_1^2-\dfrac{8v_{du}(0)^2}{3+{\rm \pi}^2}\right) \left(4\varDelta-{\rm \pi}^2A_1^2\right)}}{\rm d} A_1, \end{equation}

where $A_1(t_{p})\simeq 2\varDelta ^{1/2}/{\rm \pi}$ is approximately the amplitude of the sheet at time $t_{p}$ and $A_1(0)=2\int _0^1 y_{sh}(x,0)\sin ({\rm \pi} x)\,{{\rm d}\kern0.7pt x}=32v_{du}(0)/(3{\rm \pi} +{\rm \pi} ^3)$ is the projection of the initial configuration (3.1) on the first mode of the sheet.

In figure 11(a), we fix the excess length and the vertical dimension of the chamber at $\varDelta =0.01$ and $L_y=2$, respectively, and plot $\sigma t_{p}$ for $\lambda \in [10^{-3},10^3]$ ($\sigma \in [0.4,\sqrt {3}{\rm \pi} ^2]$). We find that at a given initial volume difference, $v_{du}(0)$, the peak time, $\sigma t_{p}$, changes by less than 5 % over more than six orders of magnitude in $\lambda$. While $\sigma t_{p}$ is almost independent of $\lambda$, it does depend strongly on the initial configuration of the sheet, $v_{du}(0)$, and the excess length, $\varDelta$. An analytical approximation of this dependence can be extracted from (4.3), if we assume that the system is in the solid-dominated region, where $\sigma \simeq \sqrt {3}{\rm \pi} ^2$. Under this assumption, we can integrate the right-hand side of this equation and take the limit $v_{du}(0)\ll 1$ of the resulting expression. (We use Mathematica (Wolfram Research 2018) for the symbolic integration.) This gives

(4.4)\begin{equation} \sigma t_{p}\simeq \ln \left(c\frac{\varDelta^{1/2}}{v_{du}(0)}\right), \end{equation}

where $c\simeq 2.9$. While the scaling in (4.4) agrees well with our numerical solution of the nonlinear model, there is a small deviation in the numerical prefactor. The best fit to the nonlinear model gives $c\simeq 2.0$ (see figure 11b). We note that the independence of $L_y$ on the right-hand side of (4.4) is a result of our assumption that $\lambda \gg 1$. Had we derived this scaling using the fluid-dominated region, we would have found that the integral in (4.3) does depend on the vertical dimension of the chamber. However, this dependence diminishes to zero when $L_y\gtrsim 1$, as is required by our fourth assumption in § 2.

Figure 11. The peak time, $\sigma t_{p}$, as a function of the system's parameters. (a) The peak time as a function of $\lambda$, where $\varDelta =0.01$ and $L_y=2$, for two different values of the initial volume difference. Dashed lines correspond to the numerical integration of (4.3) and symbols with corresponding colours represent the numerical solution of (2.2)–(2.9). For a given $v_{du}(0)$, the time $\sigma t_{p}$ changes by less than $5\,\%$ over six orders of magnitude of the parameter $\lambda$. (b) Comparison between the analytical scaling (4.4) and the solution obtained from the nonlinear model. A logarithmic scale is used on the $x$ axis. Symbols correspond to the numerical solution of the nonlinear model. While the scaling of the analytical solution agrees well with the numerical data, the prefactor $c\simeq 2.9$ slightly overestimates the numerical prediction.

We also note that the logarithmic divergence of $\sigma t_{p}$ in the limit $v_{du}(0)\rightarrow 0$ is similar to the divergence of the period of a pendulum near the separatrix (Butikov Reference Butikov1999). Within this analogy between the two problems, the small initial deviation of the pendulum from the unstable vertical position is analogous to the small initial volume difference $v_{du}(0)$. The separatrix of the pendulum is analogous to the trajectory $v_{du}(0)=0$ in the phase space of our system.

4.3. The $\bar {p}_{ud}$–$v_{du}$ relation

In this section, we investigate the behaviour of the $\bar {p}_{ud}$–$v_{du}$ relation at moderate times. To this end, we first use the two-mode approximation to calculate the pressure difference in the chamber. From Bernoulli's equation (2.13a), and the normal mode expansion (2.15) and (C2), we have that the average pressure difference on the sheet is given by $\bar {p}_{ud}(t)=({2L_y}/{{\rm \pi} \lambda })({{\rm d}^2A_1}/{{\rm d}t^2})$. In addition, using (2.16), we find the volume difference as a function of time, $v_{du}(t)=2\int _0^1 y_{sh}(x,t)\,{{\rm d}\kern0.7pt x}=4A_1(t)/{\rm \pi}$. Thereafter, we solve (2.20) numerically with $N=2$ and plot the parametric solution ($\bar {p}_{ud}(t)$, $v_{du}(t)$) in the range $t\in [0,t_{p}]$.

The results of these solutions are plotted in figures 12(a) and 12(b) for the solid-dominated ($\lambda =100$) and fluid-dominated ($\lambda =0.01$) regions, respectively. Qualitatively, the two regions exhibit similar behaviour. The pressure difference in the chamber increases from almost zero up to a maximum positive value, from which it rapidly decreases and becomes negative. The backward pressure is maximized at the peak time $t_{p}$; see the time-dependent behaviour of $\bar {p}_{ud}(t)$ in the insets of these figures. Quantitatively, however, the two profiles are considerably different, because the maximum backward pressure is much larger in the fluid-dominated region (figure 12b) than in the solid-dominated region (figure 12a).

Figure 12. The $\bar {p}_{ud}$–$v_{du}$ relation and the maximum backward pressure. In all three panels, open green triangles correspond to the numerical solution of (2.2)–(2.9) and solid lines correspond to the solution of the two-mode approximation. The volume difference as a function of the average pressure difference on the sheet is plotted in (a) the solid-dominated ($\lambda =100$) and (b) the fluid-dominated ($\lambda =0.01$) regions. In both panels, the solid black line corresponds to the solution of (2.20) with $N=2$, $\varDelta =0.01$, $L_y=2$ and $v_{du}(0)=0.01v_{du}^{cr}$. The blue points correspond to times $t=0$ and $t=t_{p}$. The dashed grey lines correspond to the static solution (3.1b) and (3.3b). The insets show $\bar {p}_{ud}(t)$ as a function of time, where $\sigma$ is given by (3.6). In the fluid-dominated region, the maximum backward pressure ($\bar {p}_{ud}(t_{p})\simeq -200$) is much larger than that in the solid-dominated region ($\bar {p}_{ud}(t_{p})\simeq -1$). In addition, the evolution of the $\bar {p}_{ud}$–$v_{du}$ relation in the fluid-dominated region follows the static solution, except for some deviations close to the asymmetric-to-symmetric transition. (c) The absolute value of the maximum average pressure difference on the sheet.

The transition from positive to negative pressure differences can be explained as follows. Initially, the system exhibits a ‘negative feedback’ between the sheet and the fluid, meaning that the sheet's motion drives the fluid's dynamics, which, in turn applies a positive pressure difference and resists the sheet's motion. As the system evolves, the fluid continuously gains kinetic energy and reduces its resistance to the sheet's motion. Then, at some instant, the pressure difference vanishes, $\bar {p}_{ud}=0$, and the resistance of the fluid is almost eliminated. Beyond this moment, the pressure difference becomes negative, and the system exhibits a ‘positive feedback’, meaning that the sheet transfers energy to the fluid, which, in turn, enhances the sheet's motion. This positive feedback accelerates the system's dynamics and creates a spike of pressure drop in the chamber. The process terminates when the system meets the constraint on the maximum volume difference.

To estimate the magnitude of the pressure spike, $\bar {p}_{ud}(t_{p})$, we use the two-mode approximation. Recalling that near the peak time when the sheet is close in shape to the first mode of buckling, i.e. $A_1(t_{p})\simeq 2\varDelta ^{1/2}/{\rm \pi}$ and $A_2(t_{p})\simeq 0$, and that $E\simeq 4{\rm \pi} ^2\varDelta$, we find from (2.20b) and (4.1) an expression for $({{\rm d}^2A_1}/{{\rm d}t^2})(t_{p})$ as a function of the system's parameters. This gives the maximum backward pressure:

(4.5)\begin{equation} \bar{p}_{ud}(t_{p})\simeq{-}\frac{432{\rm \pi}^5L_y\varDelta^{1/2}}{9{\rm \pi}^3\lambda+128\tanh({\rm \pi} L_y/2)}. \end{equation}

This analytical approximation compares well with the numerical solution of (2.2)–(2.9) (see figure 12c). Therefore, in the solid-dominated region, the maximum backward pressure decays to zero as $\bar {p}_{ud}(t_{p},\lambda \gg 1)\simeq -48{\rm \pi} ^2 L_y\varDelta ^{1/2}/\lambda$, whereas in the fluid-dominated region it is independent of $\lambda$, i.e. $\bar {p}_{ud}(t_{p},\lambda \ll 1)\simeq -(27{\rm \pi} ^5/8)L_y\varDelta ^{1/2}/\tanh ({\rm \pi} L_y/2)$.

Note also that the $\bar {p}_{ud}$–$v_{du}$ relation approximately follows the static solution (3.1b) and (3.3b) in the fluid-dominated region (see the dashed lines in figure 12b). This is because the dynamics of the fluid in this region is much slower than that in the solid-dominated region. Nonetheless, deviations between the two solutions, static and dynamic, are observed close to the asymmetric-to-symmetric transition. These deviations may be attributed to inertial effects in the dynamics solution.