Localized Faraday patterns under heterogeneous parametric excitation

TL;DR

This paper investigates how localized and heterogeneous parametric excitation influences Faraday wave patterns, revealing the formation of localized subharmonic patterns and providing a theoretical model that aligns well with experimental results.

Contribution

It introduces a combined experimental and theoretical study of Faraday waves under localized forcing, extending understanding beyond uniform excitation scenarios.

Findings

Localized forcing induces subharmonic wave patterns.

Theoretical model predicts pattern evolution accurately.

Onset of instability is altered by heterogeneity.

Abstract

Faraday waves are a classic example of a system in which an extended pattern emerges under spatially uniform forcing. Motivated by systems in which uniform excitation is not plausible, we study both experimentally and theoretically the effect of heterogeneous forcing on Faraday waves. Our experiments show that vibrations restricted to finite regions lead to the formation of localized subharmonic wave patterns and change the onset of the instability. The prototype model used for the theoretical calculations is the parametrically driven and damped nonlinear Schr\"odinger equation, which is known to describe well Faraday-instability regimes. For an energy injection with a Gaussian spatial profile, we show that the evolution of the envelope of the wave pattern can be reduced to a Weber-equation eigenvalue problem. Our theoretical results provide very good predictions of our experimental observations provided that the decay length scale of the Gaussian profile is much larger than the pattern wavelength.

Tables

Table 1. Table 1: The table shows our experiment parameters including the forcing amplitude Γ Γ \Gamma , the width of the injection region σ i subscript 𝜎 𝑖 \sigma_{i} and the width of the wave envelope σ w subscript 𝜎 𝑤 \sigma_{w} . In all cases, f = 14.86 ​ Hz 𝑓 14.86 Hz f=14.86\,\mathrm{Hz} . Errors estimated from the confidence intervals of the fitted parameters are negligible ( ≲ 0.5 ​ mm less-than-or-similar-to absent 0.5 mm \lesssim 0.5\,\mathrm{mm} ).

$n$	${\Gamma }_c$	$\Gamma$	${\sigma }_i$ (mm)	${\sigma }_w$ (mm)
6	0.370	0.443	48	42
7	0.345	0.453	54	46
8	0.345	0.392	61	51
9	0.331	0.382	68	57
10	0.331	0.382	69	66
11	0.331	0.382	72	74
12	0.331	0.382	73	85
13	0.331	0.382	80	92

Equations

$${\partial }_t\psi =-i\nu \psi -iA|\psi {|}^2\psi -iB{\partial }_{x^{\prime }{x}^{\prime }}\psi -\mu \psi +\gamma \bar{\psi },$$

$${\partial }_t\psi =-i\nu \psi -i|\psi {|}^2\psi -i{\partial }_{xx}\psi -\mu \psi +\gamma \bar{\psi }.$$

$$\gamma \left(x\right)={\gamma }_i\exp \left(-\frac{x^2}{2{\sigma }_i^2}\right),$$

$$\left[{\left({\partial }_t+\mu \right)}^2+{\left(\nu +{\partial }_{xx}\right)}^2\right]\psi -\gamma \left(x\right)\psi =0.$$

$$\left[{\left({\partial }_t+\mu \right)}^2+{\left(\nu +{ϵ}^2{\partial }_{XX}\right)}^2\right]\psi -{\gamma }^2\left(X\right)\psi =0,$$

$$\psi \left(X,t\right)=[A_0\left(X\right)+ϵA_1\left(X\right)+\mathcal{O}\left({ϵ}^2\right)]\text{mbox}{e}^{-i\omega t}{e}^{\left[\frac{i}{ϵ}\underset{X_s}{\overset{X}{\int }}k(X)dX\right]},$$

$${\left(-i\omega +\mu \right)}^2+{\left(\nu -k^2\right)}^2-{\gamma }^2\left(X\right)=0.$$

$${\partial }_y^2{A}_0+\left({\beta }^2-\alpha y^2\right)A_0=0,$$

$$A_0(x)=e^{-\frac{x^2}{2{\sigma }_w^2}},\mathrm{with}{\sigma }_w={\left(\frac{\sqrt{\nu }}{\mu }{\sigma }_i\right)}^{1/2}.$$

$${\psi }_c(x,t)=C_i(t)(X_k+i{Y}_k)e^{i\sqrt{\nu }x}+c.c.+h.o.h.,$$

$${\partial }_t{C}_i=\delta C_i-\frac{9}{2\mu }|C_i{|}^4{C}_i.$$

$$-i\omega =-\mu \pm \sqrt{{\gamma }^2(X)-(\nu -k^2{)}^2}+\mathcal{O}(ϵ).$$

$$-i\omega =-\mu +f(X,k)-\frac{9}{2\mu }|C{|}^4,$$

$${\partial }_{t}C=\frac{2\nu }{\mu }{\partial }_x^{2}C+\left(\delta -\frac{\mu }{2{\sigma }^2}{x}^2\right)C-\frac{9}{2\mu }|C{|}^{4}C,$$

$${\partial }_t{D}_0=\delta D_0-\frac{9}{2\sqrt{3}\mu }{D}_0^5,$$

$$A_0\propto \text{mbox}exp\left[\frac{-i\int F(k,\omega ,\mu ,{\gamma }_0,\nu ,X)\mathrm{d}X}{\partial \omega /\partial k}\right],$$

$$\frac{\partial\omega\left(X\right)}{\partial k}\big{|}_{k=k_{0}}=0.$$

$$\psi \left(y,t\right)=[A_0\left(y\right)+{ϵ}^{1/2}{A}_1\left(y\right)+\mathcal{O}\left(ϵ\right)]e^{-i\omega t}{e}^{\left[\frac{i}{{ϵ}^{1/2}}\underset{y_s}{\overset{y}{\int }}k(y)\mathrm{d}y\right]}.$$

$$\omega (k_c)$$

$$\gamma$$

$$\omega$$

$$\mu \left(A_1-{\bar{A}}_1\right)=2\sqrt{\nu }{\partial }_y{A}_0.$$

$${\partial }_y^2{A}_0+\left({\beta }^2-\alpha y^2\right)A_0=0$$

$$A_0\left(x\right)=H_{\text{ν}}\left({\alpha }^{\frac{1}{4}}{\sigma }_i^{-\frac{1}{2}}x\right)e^{-x^2/2{\sigma }_w^2}$$

$${\sigma }_w^2\equiv \left(\frac{1}{\sqrt{\alpha }}\right){\sigma }_i=\left(\frac{2\sqrt{\nu }}{\mu }\right){\sigma }_i.$$

$$C(x,t):=\sum _{i=1}^{\infty }{\epsilon }^{i/4}{A}_i(x,T_1,T_2,\dots ),$$

$$\delta :={\delta }_0+\epsilon {\delta }_1+{\epsilon }^2{\delta }_2+{\epsilon }^3{\delta }_3+\mathcal{O}({\epsilon }^4).$$

$${\partial }_x^2{A}_1+\frac{\mu }{2\nu }\left({\delta }_m^c-\frac{\mu }{2{\sigma }^2}{x}^2\right)A_1=0,$$

$$A_1=D_0(T_1,\dots )A_1^{(0)}+\sum _{n=1}^{\infty }{D}_n(T_1,\dots )A_1^{(n)},$$

$$A_1^{(0)}(x)=e^{-x^2/2{\sigma }_w^2},$$

Full text

I Introduction

Pattern formation is a major area of nonlinear dynamics Nicolis and Prigogine (1977); Cross and Hohenberg (1993). During the last decades, a major progress has been achieved in understanding how an extended system with homogeneous conditions can spontaneously undergo from a basic homogeneous state to a self-organized pattern Cross and Greenside (2009); Pismen (2006). However, a renewed interest has come from the observation of spatially localized states in uniform and non-uniform systems Dewel and Borckmans (1989); Mercader et al. (2011); Moriarty and Holt (2011). In uniform systems, localized patterns arise in bi-stable regions. An extended pattern solution and an homogenous one coexist, setting up a family of solutions via snaking bifurcations Burke and Knobloch (2007). In heterogeneous media, a local spatial pattern can develop from the non-uniformity of system parameters, such as forcing or dissipation. In this latter scenario, the dynamical behavior of the system suffers modifications as corrections on the instability domains and threshold discretization Coulibaly (2006); Huerre and Monkewitz (1990); Ouarzazi et al. (1996); Monkewitz et al. (1993). In particular, the concept of global mode has been introduced to characterize the synchronized response of the system to the localization of the forcing.

Alligator’s water dance is a striking example in nature of heterogenous forcing. Crocodiles and alligators are able to create spectacular local Faraday waves—spatial stationary subharmonic responses—on the water surface through the infrasonic resonance of their lungs Powell (2011); Kofron and Farris (2015). The water dance is used as an advertisement call for mating purposes of male individuals and have shown to be crucial for reproduction. Direct observations of this phenomenon, both in animals in nature or captivity, suggest that the infrasonic radiating waves spread several kilometers under water. Once females approach, localized Faraday waves on the surface of the water will provide a visual signature of the size of the animal Kofron and Farris (2015); Dinets (2013).

In the literature, there have been many efforts to study parametrically forced systems. Some of them have focused on the dynamics of localized structures such as solitons Wu et al. (1984); Gordillo (2012), local defects Lü et al. (2004); Alexeeva et al. (2000), finite-size effects Clerc et al. (2011); Gordillo and García-Ñustes (2014), linear-depth gradients in a water trough Gordillo et al. (2011); Jiaren et al. (1992) and Gaussian parametric injection in optical systems Ye et al. (2013). However, very scarce studies explore the dynamics of localized Faraday waves induced by heterogenous forcing in laboratory conditions Moriarty and Holt (2011). A thorough understanding of these systems may help us to understand how female alligators can decode the male size from Faraday waves signals. Proper tuning of heterogeneous parameters may also allow us to engineer the outcome of subharmonic out-of-equilibrium systems that could be used for technological applications Jia et al. (2013); Foster et al. (2006).

In this article, we study both experimentally and theoretically the Faraday instability generated by localized forcing. Our experimental setup consists of a water channel with a deformable bottom. The system is theoretically modeled by the parametrically driven and damped non-linear Schrödinger (PDNLS) equation Miles (1984); Miles et al. (1990) with a spatial varying forcing parameter. Assuming a Gaussian profile for the injection, we use a WKBJ scaling technique Huerre and Monkewitz (1990) to derive an eigenvalue Weber equation that governs the pattern envelope. Consequently the response of the system is discrete and the solutions are shown to be Hermite polynomials with Gaussian modulation. From the weakly nonlinear analysis, we successfully describe the nonlinear saturation of the patterns close to the threshold of the instability of the fundamental Gauss-Hermite mode. The theoretical results are in very good agreement with experiments.

The article is organized as follows. In section II we show our experimental setup and describe our measurement protocols. The theoretical description of the localized Faraday patterns are given in section III and numerical simulations, in section IV. We provide final remarks and conclusions in section V.

II Experimental setup and measurement protocol

Our experimental setup consists of a transparent rectangular water channel $15$-$\mathrm{mm}$ long, $490$-$\mathrm{mm}$ wide and $100$-$\mathrm{mm}$ deep, whose bottom has a central soft region $240$-mm wide manufactured in a soft silicone-elastomer (Shore hardness OO). The assembly rests over a system of 13 pistons evenly spaced ($\Delta x=16\,\mathrm{mm}$, each of them constrained to vertical motion by two fixed axial bearings. At the bottom of each piston, we assembled a tiny roller to be used as a follower. Using compressed springs, the pistons are pushed towards a set of rotary cams placed in a common horizontal axis. The axis is respectively coupled to a brushless motor with feedback (Model BLM-N23-50-1000-B).

The setup resembles the mechanical transmission system of a music box as shown in Fig. 1 and allows to deform the bottom of the channel with a spatial distribution. Cams are shaped in such a way that an oscillatory angular motion on the axis creates a vertical oscillatory motion on the piston. In this way, both the acceleration amplitude $\Gamma$ (normalized by the acceleration of gravity $g$) and frequency $f$ of oscillations can be easily programmed through the motor controller. The motion of any piston can be easily switched off by changing our special cams to circular ones.

The channel trough was filled with a Photoflo-water solution (concentration: 2%) up to $20$ $\mathrm{mm}$ deep. Under uniform forcing (e.g. frequency $f=14.86\;\mathrm{Hz}$), Faraday waves emerge above an acceleration threshold of $\Gamma_{c}\approx$ 0.3. The waves display a central node in the cross-wise direction. The emerging waves were visualized using a high-speed camera. The channel was front illuminated so a clear vertical cut of the flow at the wall can be observed. A small amount of white dye was added to the solution to improve visualization.

In Fig. 2, we display typical images of the observed waves for different numbers of excited pistons. The snapshots include also the reconstructed free-surface deformation at the wall $\eta(x,t)$ and the profile of the effective acceleration at the surface $\Gamma(x,t)$ at $t=t_{\mathrm{max}}$, *i.e. *when maximal deformations are observed through a cycle. The curves were calculated using standard edge-detection algorithms on the sequence of images with a sensitivity of $0.1\,\mathrm{mm}$.

Examples of the full spatiotemporal evolution of the free-surface $\eta\left(x,t\right)$ are plotted in Fig. 3. While obtaining $\eta\left(x,t\right)$ is straightforward from the wavy free surface, the effective acceleration at the surface, $\Gamma\left(x,t\right)$, was obtained by imaging the surface of the liquid slightly below the instability threshold $\Gamma_{c}$ where the response is still linear, and then rescaling by the driving amplitude $\Gamma$ used in the experiments. It is a known issue, that the fluid layer acts as a longpass filter of the bottom deformation so the effective deformation/acceleration at the free surface is a smoothened version of the bottom driving Hammack (1973); Jamin et al. (2015).

As shown in Figs. 2 and 3, the excitation of a reduced bottom region generates wave patterns that are spatially localized. The patterns oscillate at half the forcing frequency (parametric instability), which is the Faraday-waves signature. The observed patterns, that we will refer to as localized Faraday waves, have a standing-wave core that emits evanescent waves toward the unperturbed regions.

II.1 Pattern vs. injection length

To compare how the wave pattern length depends on the injection length, our measurement protocol was the following: First, we chose the frequency in such a way that the wavelength of the Faraday waves matches the inter-piston distance. Starting from the maximum number of excited pistons, i.e. $n=13$, we acquired a video sequence of highly resolved images in both space and time. Then, we sequentially decreased the number of excited pistons, hence reducing the length of the injection region.

To characterize the wave localization, we perform the following analysis on $\eta\left(x,t\right)$: First, we apply the temporal Fourier transform and extract the phase and amplitude for the dominant frequency $f/2$. In the standing-wave region, the wave displays a constant phase along $x$ that decreases linearly as we enter into the wave-emission region. The amplitude on the other hand, displays a smooth decay in the wave-emission region but a serrated shape on the standing-wave. The spatial envelope of the amplitude for the whole domain was obtained by fitting a Gaussian curve on the amplitude local maxima of the standing-wave region and the remaining tails of the wave-emission one. The width of the envelope is defined as $\sigma_{w}\equiv\text{HWHM}$, i.e. the half width at half maximum. Likewise, to characterize the injection localization, we obtained the envelope of $\Gamma\left(x,t\right)$ and straightforwardly obtain the $\mathrm{HWHM}$, $\sigma_{i}$. The two quantities, the wave-envelope width $\sigma_{w}$ and the injection-envelope width $\sigma_{i}$, were measured for runs with different number of pistons. Results are summarized in Table 1 and displayed in Fig. 4.

II.2 Onset of localized Faraday waves

To study how localized Faraday waves emerge, we designed two-fold protocol for a given injection length ($n=13$). First, we started from the flat state and increased the amplitude of oscillation with fine steps (0.1 mm, $\Delta\Gamma=0.086$), starting from $\Gamma_{i}=0.265$ up to $\Gamma_{f}=0.308$. For each given $\Gamma$, we waited $\approx 45$ min. and checked that the wave pattern was stationary throughout several cycles before making the measurements. The second protocol was the same but we started from $\Gamma_{f}$ and decreased sequentially down to $\Gamma_{i}$. The results are shown in Fig. 5. The system does not display hysteretic behavior, which is the signature of supercritical bifurcations.

III Theoretical description of localized Faraday waves

It has been shown that the hydrodynamical problem of the free surface of a fluid which is oscillated vertically in the vicinity of Faraday instability can be reduced to an amplitude equation for the envelope of the surface: the PDNLS equation Miles (1984); Miles et al. (1990). Subsequently, the equation has been derived in different context as nonlinear lattices Denardo et al. (1992), optical fibers Kutz et al. (1993), Kerr-type optical parametric oscillators Longhi (1996), easy-plane ferromagnetic materials exposed to oscillatory magnetic fields Barashenkov et al. (1991); Clerc et al. (2012) and parametrically driven damped chains of pendula Alexeeva et al. (2000). The governing equation for the envelope of the water surface displacement of the transversal mode is

[TABLE]

where $\psi(x^{\prime},t)$ stands for the complex envelope of the standing waves and $\bar{\psi}$, its complex conjugate; $t$ is the dimensionless time. Besides, $\nu$ is the detuning parameter which measures the frequency offset to the parametric resonance, $\mu$ is the damping parameter, and $\gamma$ stands for the amplitude of the parametric forcing. The parameters $A$ and $B$ are functions of the wavenumber $k$. In particular, $A\sim k^{2}$ and $B\sim 1/k^{2}$. The relation between the experiment quantities and the dimensionless parameters $\nu,\mu,\gamma$, as well as the envelope $\psi$ and the time, are given in detail in Refs. Miles (1984); Gordillo and Mujica (2014); Périnet et al. (2017). For our experiments, it can be shown that $\left|\psi\right|^{2}\sim 10^{-3}$, $\mu\sim 10^{-2},\nu\sim 10^{-1}$ and $\gamma\sim 10^{-1}$. Notice that the PDNLS equation (2) applies only in the limit $\nu\sim\mu\sim\gamma\ll 1$. For $\gamma>\mu$, $\nu>0$ and $\gamma^{2}<\nu^{2}+\mu^{2}$, the system exhibits subharmonic patterns with critical wavelength $k_{c}=\pm\sqrt{\nu}$, i.e., Faraday waves. Setting a dimensionless variable $x^{\prime}\equiv\sqrt{B}x$ where is $x$ is the dimensionless space variable, we can rewrite the equation in the following dimensionless form,

[TABLE]

III.1 Linear stability analysis

We extend Eq. 2 to heterogeneous systems by assuming that $\gamma\equiv\gamma\left(x\right)$ is a function describing the spatial profile of the forcing. Following our experimental results, see e.g blue dashed lines in Fig. 2, we assume that the injection profile $\gamma$ is a localized function satisfying three key features: i) $\gamma\left(x\right)$ is symmetric with respect to a given $x_{0}$, ii) $\gamma\left(x\right)$ has a single extremum at $x_{0}$ with non-vanishing second derivative, iii) decays to [math] as $x\rightarrow\pm\infty$. For the sake of simplicity, we choose $\gamma\left(x\right)$ to be a Gaussian function:

[TABLE]

where $\gamma_{i}$ is the forcing amplitude and $\sigma_{i}$ is a dimensionless standard deviation, $\sigma_{i}=\sigma^{\prime}_{i}/\sqrt{B}$, where $\sigma^{\prime}_{i}$ is a space parameter.

The degree of heterogeneity of the system can be modeled through the parameter $\epsilon\equiv 1/\sigma{}_{i}$, which should be small to satisfy the condition of slow spatial dependence in $\gamma\left(x\right)$. In the original physical variables, the condition $\epsilon\ll 1$ is equivalent to require that the pattern wavelength $\lambda$ is much smaller than the variation length-scale of the forcing. Indeed, $\epsilon=1/\sigma_{i}=\sqrt{B}/\sigma^{\prime}_{i}\sim\lambda/\sigma^{\prime}_{i}\ll 1$, thus $\lambda\ll\sigma^{\prime}_{i}$. Notice that this assumption is in agreement with our experimental observations.

We linearize Eq. 2 around the trivial homogeneous steady state and analyze the result and its complex conjugate (for details, see appendix A). Considering $\gamma(x)$ as a slowly varying function of space, we obtain after some algebra,

[TABLE]

Introducing the slowly varying variable $X\equiv\epsilon x$, Eq.(4) becomes

[TABLE]

where $\gamma\left(X\right)=\gamma_{i}\left(1-X^{2}/2\right)$ is the Taylor series expansion of $\gamma$ up to second order.

A particularly well-suited approach to find solutions in the limit $\epsilon\ll 1$ is the WKBJ approximation Bender and Orszag (1978); Dewel and Borckmans (1989); Huerre and Monkewitz (1990); Monkewitz et al. (1993); Ouarzazi et al. (1996); Coulibaly (2006). Thus, we propose the following expansion

[TABLE]

where $A\left(X\to\infty,t,\epsilon\right)=0$, $k\left(X\right)$ is the pattern wavenumber, and $X_{s}$ is the source point where the heterogeneous profile is centered. Substituting the expansion (6) in (5), we obtain at zeroth-order ($\epsilon^{0}$) the dispersion relation

[TABLE]

Calculating higher orders of $\epsilon$ (see appendix A), the solution of (5) can be expressed in terms of a carrier wave (wavelength $\sqrt{\nu}$ at dominant order) and an envelope $A_{0}\left(y\equiv\sqrt{\epsilon}x\right)$, which obeys

[TABLE]

where $\alpha\equiv\mu^{2}/4\nu$ and $\beta^{2}\equiv\mu\left(\gamma_{1}/2\nu+i\omega_{1}/2\nu\right)$. Here $\gamma_{1}$ and $\omega_{1}$ have been introduced as $\gamma\simeq\mu+\epsilon\gamma^{\left(1\right)},\quad\omega\simeq\epsilon\omega^{\left(1\right)}.$ Equation (8) is a linear eigenvalue problem with a discrete set of solutions $A_{0,m}\left(y\right)$, each given in terms of the $m$-th Hermite polynomial $H_{m}$ modulated by a Gaussian function; i.e. $A_{0,m}\left(y\right)=H_{m}\left(\alpha^{1/4}y\right)e^{-\left(\sqrt{\alpha}/2\right)y^{2}}$. The eigenvalue problem also requires $\beta^{2}/\sqrt{\alpha}=2m+1$. Further calculations show that the related quantities $\gamma_{m}^{\left(1\right)}=\left(2m+1\right)\sqrt{\nu}$ and $\omega_{m}^{\left(1\right)}=0$ are now discrete. Likewise, $\gamma_{m}=\mu+\left(2m+1\right)\left(\frac{\sqrt{\nu}}{\sigma_{i}}\right),\omega_{m}=0.$ Notice that the corrections of the mode thresholds $\gamma_{m}^{c}$ are inversely proportional to the parameter $\sigma_{i}$. This is a counterintuitive result: smaller volumes of water require higher forcing to generate patterns compared to larger ones.

We can infer that the first emerging mode in our experiments is the fundamental one ($m=0$), which in terms of $x$ reads,

[TABLE]

This means that all the solutions derived from (5), including the fundamental mode (9), display localization, which is the key qualitative feature of localized Faraday waves.

Equation (9) also shows that the width $\sigma_{w}$ of the envelope $A_{0}$ scales as the square root of the injection-region width $\sigma_{i}$, which agrees our experimental observation [see Fig. (4)]. Indeed, to make quantitative comparisons, we first determine the dimensionless quantities of Eq. (2) in terms of experimental parameters. The formulas provided in Miles (1984); Gordillo and Mujica (2014) can be used to directly compute $\gamma_{i}$ and $\nu$. The formulas available in the literature for the damping coefficient $\mu$ however seem to be not appropriate for our experimental conditions since in our setup an extra shear occurs on the fixed wall. Hence we considered $\mu$ as a phenomenological parameter, which we estimated by fitting (9) on the experimental data of Fig. 4. To test the validity of our results we fitted a power law $\sigma_{w}=a\sigma_{i}^{p}$ and found $\left(a,p\right)=\left(2.27,0.480\right)$. The exponent $p$ is remarkably consistent with the predicted square-root dependence in (9).

III.2 Weakly nonlinear analysis

To describe the nonlinear saturation of the unstable global modes, we have done a weakly nonlinear analysis of the system close to the spatial instability. We introduce a bifurcation parameter $\delta\equiv\gamma_{0}-\gamma_{0}^{(1)}$ and a slowly varying amplitude $C_{i}(t)$ on the oscillations of the critical mode, i.e.

[TABLE]

where h.o.h. denotes the higher order harmonics. At the first order of nonlinearity, one obtains that $C_{i}$ is governed by the well known normal form previously derived by Coullet et al. Coullet et al. (1994),

[TABLE]

Here, $C_{i}$ corresponds to the amplitude of oscillations of a single oscillator. To introduce the spatial dependence of the amplitude $C$, we first consider the dispersion relation (7), obtained from the WKBJ formalism at order $\epsilon^{0}$,

[TABLE]

Assuming that the nonlinear part of Eq. (11) has slow variations in space and time compared to the linear terms, then Eq. (11) in the Fourier space reads $-i\omega=-\mu+\gamma_{0}-9|C_{i}|^{2}/2\mu$. Thus, from Eq. (12) we notice that to introduce the spatial dependence in our system we must map $\gamma_{0}\to[\gamma^{2}(X)-(\nu-k^{2})^{2}]^{1/2}$. One obtains for the fully heterogeneous system in the Fourier space the following expression:

[TABLE]

where $f(X,k)\equiv[\gamma^{2}(X)-(\nu-k^{2})^{2}]^{1/2}$. To consider the growth of modes with wavenumber $k_{c}=\pm\sqrt{\nu}$ due to the parametrically extended Gaussian excitation, we consider a Taylor expansion of the function $f$ in (13) for $X\equiv\epsilon x$ and $k\sim k_{c}=\pm\sqrt{\nu}$. Neglecting terms of order $\mathcal{O}(\nu^{3/2},X^{4},k^{4})$ and taking the inverse Fourier transform of (13), one obtains that the amplitude $C$ of the critical mode is governed by

[TABLE]

which is a dynamical Weber-like equation with a quintic nonlinearity. In the homogeneous limit ($\sigma\to\infty$), Eq. (14) is similar to the normal form in León et al. (2014) for a spatial supercritical quintic bifurcation in a homogeneously driven magnetic system. Notice that after taking the linear limit of (14) one recovers the Weber equation (8) for $\omega=0$.

To describe the nonlinear saturation of the fundamental mode, we have used a multiple-scale expansion in Eq. (14) to derive an evolution equation for its amplitude. After some straightforward calculations—detailed in appendix B—one obtains that the amplitude $D_{0}$ of the fundamental mode is governed by

[TABLE]

from which follows the stationary solution $D_{0}^{s}=(2\sqrt{3}\mu\delta/9)^{1/4}$. This result agrees with the experimental scaling law $\eta_{\max}\propto(\Gamma-\Gamma_{0})^{1/4}$ of Fig. 5 and describes the evolution of a supercritical quintic bifurcation. Furthermore, our theory predicts a scaling coefficient $(2\sqrt{3}\mu/9)^{1/4}=0.482$ based on fundamental quantities ($\mu=0.14$), which is shown as a dashed line in Fig.5 with no fitting parameters.

IV Numerical simulations

As a final check, we also performed direct numerical simulations of Eq. (2) with no-flux boundary conditions and $\gamma\left(x\right)$ given by the expression (3). The first goal is to determine if the linear approximation we did to obtain Eq. (5) remained valid for the set of $\left(\mu,\nu,\gamma,\sigma_{i}\right)$ that we chose. Using a 400-point spatial grid with resolution $\mathrm{d}x=0.25$, we used finite-differences of second order of accuracy for the space derivatives of Eq. (2). For the time integration we ran a fourth-order Runge-Kutta scheme with a time step $\mathrm{d}t=0.0001$.

To compare the solutions of the PDNLS equation with experiments, it is important to remark that Eq. (2) gives only a stroboscopic evolution of the surface instability. Using $\phi\equiv\psi\exp(i\pi ft)$, where $\psi$ is the solution to the PDNLS equation (2), from $\Re e(\phi)$ we can recover the non-stroboscopic picture of the Faraday patterns. We plot the result in Fig. 6(a), showing that the numerical solutions of Eq. (2) not only successfully reproduce the envelope of the localized Faraday waves but also its evanescent waves in agreement with our experiments (shown in Fig. 3).

An interesting feature that we put under test is that according to our linear stability analysis, even for $\gamma_{0}>\mu$, there will be no pattern formation if $\gamma_{0}<\gamma_{0}^{(0)}$. In this case, the amplitude of the injection is greater than the dissipation but is too localized in space to sustain an instability. We have confirmed this prediction in our experiments and numerical simulations. Indeed, any initial perturbation on the system will be eventually dissipated and end up by decaying into the homogeneous stable solution $A_{h}=0$.

However, if we increase $\gamma_{0}$ until we reach the region $\gamma_{0}^{(0)}<\gamma_{0}<\gamma_{0}^{(1)}$, a localized pattern will appear due to the instability of the fundamental Gauss-Hermite mode. This case is shown in the numerical simulation of Fig. 6(a), where we used the homogeneous solution $A_{h}=0$ as the initial condition and added a small-amplitude additive noise of order $10^{-3}$. These small fluctuations are enough to trigger the Faraday instability. According to equation (15), the amplitude of the pattern begins to grow at an exponential rate $D_{0}\propto\exp(\delta t)$. Nonlinear contributions become important as the amplitude of the instability grows, and the maximum amplitude of the pattern, $|C|$, saturates due to the quintic nonlinearity of Eq. (15), as evidenced in Fig. 6(b). In Fig. 6(c), we show that the real part of the solution after saturation displays an envelope that has a nearly Gaussian profile.

To give a more detailed insight of the stability properties of the Gauss-Hermite modes, we have numerically computed the spectrum of the linearized system following a similar procedure as in Ref. Clerc et al. (2010). First, we took the numerical solution of the PDNLS equation (2) after the envelope of the Faraday pattern has become steady. Then, we calculate the set of eigenvalues and eigenfunctions of the linear operator that describes the dynamics of small perturbations around the Faraday-pattern solution. Typical spectra are shown in Fig. 7 as well as the eigenfunctions obtained numerically for the first two eigenvalues. We expect from Eq. 4 that the spectrum of the linear operator is degenerate. We have confirmed this fact in our numerical results. For $\gamma_{0}^{(0)}<\gamma_{0}<\gamma_{0}^{(1)}$ only the real part of the first eigenvalue is positive, as showed in Fig. 7(a). In this case, the Faraday pattern is formed due to the instability of the fundamental Gauss-Hermite mode. If we increase the value of $\gamma_{0}$ until we reach the region $\gamma_{0}^{(1)}<\gamma_{0}<\gamma_{0}^{(2)}$, the eigenvalue of the first antisymmetric mode crosses the imaginary axis, as shown in figure 7.b. In this case, the Faraday pattern is formed due to the contributions of modes $m=0$ and $m=1$, which are both unstable.

Finally, we have verified with several numerical simulations that the Gauss-Hermite modes turns unstable at the values of $\gamma_{0}$ predicted by the linear stability analysis. Figure 8 shows the real part of the eigenvalues of the Gauss-Hermite modes from $m=0$ to $m=3$ as a function of the injection $\gamma_{0}$. The vertical dashed lines indicates the theoretical values of the thresholds of instability for each of the modes. It is clear that each of the Gauss-Hermite modes turns unstable at the predicted values of $\gamma_{0}$.

In summary, the results are very consistent. The theory does not only match well the full-numerical simulation but also the experimental results.

V Conclusions

In conclusion, we have designed an experimental setup with energy injection in a spatial region whose extent can be controlled. The setup consists of a quasi one-dimensional rectangular water channel with a soft deformable bottom that can be forced with a set of pistons. Above certain threshold of vibration amplitude, the water surface destabilizes into subharmonic Faraday waves that are localized in space and emit evanescent waves.

Assuming that the width of the injection region is larger than the characteristic pattern wavelength, we developed a WKBJ approximation and a weakly nonlinear analysis on the prototype model, i.e., the heterogenous PDNLS equation, to describe our experimental observations. Using this framework we have: i) derived the spatial profile of the observed patterns; ii) showed the emission of evanescent waves; iii) computed the dependence of the envelope width on the length of the injection region; iv) described how the parameter space of the onset of instability modifies and discretizes with injection localization; v) determined that localized Faraday waves emerge via a supercritical quintic bifurcation. The results presented here are helpful to understand the impact of parameter heterogeneities on pattern-formation processes in general extended physical systems. It is noteworthy that this work is the first complete study of localized Faraday waves in laboratory conditions and provide key physical and mathematical insights on alligator’s water dance. Further studies on localized injection with a shapeable bottom is in progress.

Acknowledgements.

M.A.G-N. thanks for the financial support of grant FONDECYT 11130450. L. G. was partially supported by Conicyt FCHA/Postdoctorado Becas Chile 74160007, Conicyt PAI/IAC 79160140 and FONDECYT/Iniciación 11170700. S.C., M.A.G-N., J.F.M. and M.T. thanks ECOS-Sud n° C15E06. J.F.M. thanks CONICYT/Doctorado Nacional 21150292. H. U. was partially supported by Conicyt FCHA/Beca de Doctorado en el Extranjero, Becas Chile N*∘* 72180269.

Appendix A Solutions of linearized PDNLS equation under heterogeneous forcing $\gamma(x)$

At the next order $\mathcal{O}\left(\epsilon\right)$, we obtain the equation, $2k\left(\nu-k^{2}\right)\partial_{X}A_{0}-A_{0}\left(\nu-3k^{2}\right)\partial_{X}k=0$. Using relation (7), we can now deduce an expression for $A_{0}$, i.e.

[TABLE]

where $F(k,\omega,\mu,\gamma_{0},\nu,X)\equiv\frac{\nu-3k^{2}}{-i\omega+\mu}\partial_{X}k$. It is clear that for $\partial\omega/\partial k=0$, expression (16) is singular. The points where the singularity $\partial\omega/\partial k=0$ take place are called turning points Bender and Orszag (1978). On these points, the WKBJ approach is not longer valid. Let be $\omega_{0}\equiv\omega\left(k_{0}\right)$, where $k_{0}$ is such that

[TABLE]

By definition, if the contour in the complex plane $\mbox{Im}[\omega_{0}\left(X\right)]$ exhibits a contact line or pinch between two branches $X^{\pm}(\omega_{0})$ in a turning point $X^{t}$ that also verifies $\partial\omega_{0}\left(X^{t}\right)/\partial_{X}=0$, then we deal with a double turning point Huerre and Monkewitz (1990); Monkewitz et al. (1993); Dewel and Borckmans (1989); Ouarzazi et al. (1996); Coulibaly (2006). In that case, the scaling law that rules the system dynamics in the region close to the double turning point is $X=\epsilon^{1/2}y$. In consequence, the expansion takes the following new form,

[TABLE]

As in the homogeneous case, we have two solutions for $\partial\omega/\partial k=0$, i) $k_{c}=0$ and ii) $k_{c}=\pm\sqrt{\nu}$ with $\nu>0$, which correspond respectively to the angular frequencies

[TABLE]

Imposing the second condition ($\partial\omega/\partial X=0$) for the double turning point, it follows that $X^{t}=0$ for both cases. On the contrary, for different $k_{c}$, we have different values of $\gamma^{t}$. For Faraday waves ($\nu>0$, $\gamma>\mu$), the critical wavelength is $k_{c}=\pm\sqrt{\nu}$ and $\gamma^{t}=\mu$.

After making the corresponding replacements in the parameter expansions and provided that $\epsilon\ll 1$, the spatial forcing takes the form $\gamma\left(x\right)=\gamma_{i}(1-\epsilon\left(y-y^{t}\right){}^{2}/2+\mathcal{O}\left(\epsilon^{2}\right)$. We hence introduce small deviations from the turning point $\left\{X^{t},\gamma^{t}\right\}$ through

[TABLE]

Thus, we get the expression for the forcing in the turning point, $\gamma\left(y\right)=\gamma^{t}+\epsilon\left[\gamma^{\left(1\right)}-\gamma_{0}^{t}\left(y-y^{t}\right){}^{2}/2\right]+\mathcal{O}\left(\epsilon^{2}\right)$. At dominant order. $k=k_{c}=\sqrt{\nu}$ and Eq. (18) reduces to $\psi\left(y,t\right)\sim\mbox{{exp}}\left[i\int\sqrt{\nu}\,\mathrm{d}y/\epsilon^{1/2}\right]$ Dewel and Borckmans (1989). Next, we replace $\psi\left(y,t\right)$ in (4) and analyze the equations in orders of $\epsilon$. At order $\mathcal{O}\left(\epsilon^{0}\right)$, we obtain that $A_{0}=\bar{A_{0}}$. At order $\mathcal{O}\left(\epsilon^{1/2}\right)$, we get the relation

[TABLE]

Finally at order $\mathcal{O}\left(\epsilon\right)$, using $A_{0}=\bar{A_{0}}$ and relation Eq. (23), we obtain a Weber equation that describes the linear behavior of the signal envelope $A_{0}$,

[TABLE]

where \text{\alpha}\equiv\mu^{2}/4\nu and $\beta^{2}\equiv\mu\left(\gamma^{\left(1\right)}/2\nu+i\omega^{\left(1\right)}/2\nu\right)$. The solutions of (24) are Hermite polynomials with a Gaussian modulation, i.e. $A_{0}=H_{m}\left(\alpha^{1/4}y\right)e^{-\left(\sqrt{\alpha}/2\right)y^{2}}$. Due to the discrete spectrum of the linear operator in (24), we also require $\beta^{2}/\sqrt{\alpha}=2m+1$, which imposes conditions over $\gamma^{\left(1\right)}$ and $\omega^{\left(1\right)}$. In terms of the original variable $x$, the solution is,

[TABLE]

where $H_{m}$ are the Hermite polynomials and the standard deviation $\sigma_{w}$ is given as

[TABLE]

We have shown that a weak spatial dependence of $\gamma$ generates a modulation on the wave pattern given by an amplitude equation equivalent to a Weber-equation eigenvalue problem.

Appendix B Nonlinear saturation of the fundamental mode

In this appendix, we give the derivation of the evolution equation for the amplitude of the fundamental Gauss-Hermite mode, when the system is close to the threshold of instability $\gamma_{0}^{c}$. Let us define $\delta_{m}^{c}:=\gamma_{m}^{c}-\mu$, where $\gamma_{m}^{c}=\mu+\sqrt{\nu}(2m+1)/\sigma_{i}$. We perform a multiscale development Peyrard and Dauxois (2004) in eq. (14) introducing new variables with different time-scales according to $T_{i}:=\varepsilon^{i}t$ (with $i=1,2,\dots$). Here, $\varepsilon$ is a small adimensional parameter introduced only for the multiscale analysis. We search the field $C(x,t)$ in the form of a perturbative development of functions of the different time-scales according to

[TABLE]

with a similar expansion in the bifurcation parameter

[TABLE]

Considering all these developments in eq. (14), we proceed to analyse the system at each order of $\varepsilon$.

B.1 Order $\varepsilon^{1/4}$

At order $\varepsilon^{1/4}$, Eq. (14) reads

[TABLE]

which is the Weber equation. Solutions of Eq. (27) are given by the Gauss-Hermite polynomials, which are denoted here as $A_{1}^{(n)}(x)\equiv H_{n}(\alpha^{1/4}x/\sigma_{i}^{1/2})\exp(-x^{2}/2\sigma_{w}^{2})$. Thus, the general solution of Eq. (27) can be written as

[TABLE]

which is a linear combination of the Gauss-Hermite modes with time-dependent coefficients. The fundamental mode is given by

[TABLE]

where $\sigma_{w}^{2}\equiv\sigma_{i}\sqrt{\nu}/\mu$. Thus, the saturation of the fundamental mode (29) will be given by the time evolution of the coefficient $D_{0}(T_{1},\ldots)$ in eq. (28).

B.2 Order $\varepsilon^{5/4}$

At order $\varepsilon^{5/4}$, Eq. (14) can be written as the linear problem

[TABLE]

where

[TABLE]

The linear problem (30) can be solved for $A_{5}$ only if $b$ is in the image of the operator $L$. According to the Fredholm alternative Pismen (2006), at least one solution for $A_{5}$ exists if $\exists!|v\rangle\in\ker L^{\dagger}$ such that $\langle v|b\rangle=0$. Notice that $A_{1}^{(0)}\in\ker L^{\dagger}$, since $LA_{1}^{(0)}=0$ and $L=L^{\dagger}$. Thus, the Fredholm alternative gives

[TABLE]

The modes $A_{1}^{(n)}$ for $n\geq 2$ are all stable near the threshold of instability of the fundamental mode. Thus, the amplitudes $D_{n}(T_{1})$ for $n\geq 2$ decays exponentially in time. Once the pattern has completely evolved, one simply obtains

[TABLE]

which is an even function in space. Inserting Eq. (34) in equations (32) and (33), one obtains the solvability condition

[TABLE]

from which follows Eq. (15).