A new integrable model of long wave–short wave interaction and linear stability spectra

We consider the propagation of short waves which generate waves of much longer (infinite) wavelength. Model equations of such long wave–short wave (LS) resonant interaction, including integrable ones, are well known and have received much attention because of their appearance in various physical contexts, particularly fluid dynamics and plasma physics. Here we introduce a new LS integrable model which generalizes those first proposed by Yajima and Oikawa and by Newell. By means of its associated Lax pair, we carry out the linear stability analysis of its continuous wave solutions by introducing the stability spectrum as an algebraic curve in the complex plane. This is done starting from the construction of the eigenfunctions of the linearized LS model equations. The geometrical features of this spectrum are related to the stability/instability properties of the solution under scrutiny. Stability spectra for the plane wave solutions are fully classified in the parameter space together with types of modulational instabilities.

We consider the propagation of short waves which generate waves of much longer (infinite) wavelength. Model equations of such long wave-short wave (LS) resonant interaction, including integrable ones, are well known and have received much attention because of their appearance in various physical contexts, particularly fluid dynamics and plasma physics. Here we introduce a new LS integrable model which generalizes those first proposed by Yajima and Oikawa and by Newell. By means of its associated Lax pair, we carry out the linear stability analysis of its continuous wave solutions by introducing the stability spectrum as an algebraic curve in the complex plane. This is done starting from the construction of the eigenfunctions of the linearized LS model equations. The geometrical features of this spectrum are related to the stability/instability properties of the solution under scrutiny. Stability spectra for the plane wave solutions are fully classified in the parameter space together with types of modulational instabilities. Finding out the dynamical behaviour of waves usually poses mathematical problems, which are very rarely solvable by analytical methods. More frequently, approximate solutions may be found by treating nonlinear terms of the equations of motion as small perturbations of the linear part. In this framework, investigating how the linear superposition of two or more plane waves changes due to nonlinear effects leads to predictions on wave-wave interaction. A well-known approach of this kind, the multiscale method (e.g. [1,2] and references therein), requires that Fourier wave amplitudes be small, say O(ε), and be functions only of appropriately ε-rescaled space and time coordinates, where ε is a small dimensionless parameter. In this way, one obtains partial differential equations in the rescaled variables which are generally simpler, and sometimes even integrable. When they are, these approximate wave equations, even if nonlinear, are exactly solvable by spectral techniques.
The best-known example of such model equations is the nonlinear Schrödinger (NLS) equation in a one-dimensional space: where ξ = ε(x − vt) and t 2 = ε 2 t are the rescaled space and time coordinates, and the asterisk indicates complex conjugation. This equation follows via the multiscale method from almost any real propagation equation and provides the lowest-order effect of the nonlinear terms on the linear solution εS(ξ , t 2 ) e i(kx−ωt) + εS * (ξ , t 2 ) e −i(kx−ωt) , t 2 = ε 2 t, ξ = ε(x − vt), given by the sum of two plane waves with amplitudes S and S * , and with wavenumber k, frequency ω and group velocity v. If more plane waves are superimposed, cross-interaction takes place if the wavenumbers and frequencies satisfy certain resonance conditions. Among the many cases of interest and physical significance, we consider here wave equations which model the interaction of two smallamplitude plane waves: one with a very long wavelength (k = 0) and a real amplitude L, and a second one with a much shorter wave length and a complex amplitude S. Equations of this kind have been derived and proposed with various motivations, which mainly come from plasma physics [3] and fluid dynamics [4,5].
On the mathematical side, long wave-short wave (LS) equations have been differently introduced because they are, at the same time, close to physical equations and integrable. This second distinctive property allows the construction of a broad range of solutions with relevant and useful insight into experimental observations. The focus of the present work is on LS models which are integrable.
Two very well-known examples of integrable LS models are the Yajima-Oikawa (YO) equation [3] iS t + S xx − LS = 0, L t = 2(|S| 2 ) x , (1.2) and the alternative integrable Newell (N) wave system, introduced in [5], where, in addition to an LS coupling, the short wave has the same self-interaction as the NLS equation; see (1.1). The two systems (1.2) and (1.3) are related by a Miura transformation (see [6]).
In [2,7,8] the YO equation has also been derived via multiscale method applied to a generic real wave equation, with the resonance condition that the the group velocity v S = ω (k S ) of the short wave equals the group velocity v L = ω (0) of the long wave, say v S = v L = v, where ω(k) is the linear dispersion relation, ω (k) = dω(k)/dk, and k S is the non-vanishing wavenumber of the short wave. waves, either in vector form, as in [9], or in matrix form as in [10] requiring a higher rank matrix Lax pair. A third LS wave scalar equation, which reads (1.4) has been derived in [11] from a 3 × 3 matrix Lax pair, which is however different from the Lax pair investigated in the present paper (see § §2 and 3), and is not discussed here. In §3, we show that the YO and N equations, (1.2) and (1.3), do not need to be treated separately. Indeed, these two model equations can be combined in just one system, namely which we refer to as YON equations. This system coincides with the YO equation (1.2) for α = 0, β = 1, while it reads as the N equation (1.3) by setting α = σ , β = 0 and by substituting the field L with σ L (σ = ±1). In fact, as we will show in §3, the novel system (1.5) is integrable for any real value of α and β. Thus it is likely to be relevant to applications as it is more flexible in modelling the LS interaction. Moreover, this unifying result makes our present overall analysis and discussion of the two model equations (1.2) and (1.3) simpler and more compact.
Other LS integrable equations may be sorted out via transformations of the wave fields. Indeed, the two integrable LS equations (1.5) and (1.4) may be given an equivalent, even simpler, form by performing the gauge transformation [5,11] S( while the same transformation, as applied to the integrable equations (1.4), leads to the LS system The integrability of (1.2), (1.3) and (1.4), and of their transformed versions via (1.6), allows the construction of special exact solutions of these model equations, such as solitons and rogue waves (e.g. [9,[12][13][14][15]). Integrability is also a key property to address one of the basic issues in nonlinear wave dynamics, namely that of linear stability against small changes of the initial conditions (e.g. [16] and references therein). For instance, the modulational instability of a periodic wave train, as in fluid dynamics and optics, has been well exemplified via the NLS equation (1.1) with focusing self-interaction (σ = −1). Indeed, the NLS equation has deserved special interest as it proves to yield, for focusing interaction, a simple description of the instability of a regular wave train, a phenomenon first observed by Talanov [17] in optics and by Benjamin & Feir [18] in a water tank. For this particular equation, understanding modulational instability of continuous wave solutions may be achieved via standard Fourier analysis. However this analysis fails if the unperturbed solution is not just a plane wave [19], or if several nonlinear wave trains have a resonant interaction. In general, when dealing with integrable systems, the stability of interacting plane waves is better treated by means of integrability techniques rather than by a Fourier approach. Indeed, the Lax pair plays a key role as it allows one to extend the linear stability analysis to other wave solutions, such as those constructed from plane waves by Darboux dressing methods [20].
The integrability method has been used to unveil the stability properties of plane wave solutions of two coupled NLS equations [16,21]. For these model equations, instabilities have been fully classified in terms of coupling constants, amplitudes and wavenumbers, including instability effects due to defocusing self-and cross-interactions. As for the LS resonance interaction, in addition to orbital stability [22] and transverse stability [23], also linear stability has been considered for solutions of particular models [5].
On the wake of those results, we investigate here the linear stability of continuous wave solutions of the YON equations (1.5). Although the presentation is intended to be self-contained, we refer the reader to [16,21] for specific proofs and further details of the general method. The main target of our approach is the computation of the stability spectrum (S), associated with the given plane wave solution. For a multicomponent system such as (1.5), this spectrum does not generically coincide with the ordinary Lax spectrum. It turns out to be a piecewise continuous curve in the complex plane of the spectral variable λ. In each point of this curve, one or more eigenmode solutions of the linearized equation of motion are well defined (see §2). The relevant stability-instability properties of the wave solution are readily conveyed by the geometrical properties of the spectral curve S itself. Additional relevant targets are the linearized equation eigenfunction frequencies and their associated gain functions, which tell us whether the wavenumber instability bands are of passband or baseband type [24]. This paper is organized as follows. In §2, we review our approach to linear stability; in particular, we sketch the role that integrability plays in investigating the linear stability of a given solution. In §3, the integrability of the unified YON model is established by showing its corresponding Lax pair. Moreover, we consider the plane wave solution of the YON model equations, and study its linear stability by constructing and classifying the corresponding stability spectra. We compute the gain function, characterize the instabilities and show that the plane wave solution is unstable for a generic choice of the physical parameters. We conclude with some remarks and outlooks in §4.

Integrable equations and small variation of their solutions
Our approach to the linear stability of solutions of integrable partial differential equations starts from the associated Lax pair, a characteristic feature of integrable systems: In the present context, the solution Ψ (x, t, λ) and the coefficients X(x, t, λ), T(x, t, λ) are 3 × 3 matrices. The latter two matrices are assumed to have the following polynomial dependence on the complex spectral variable λ: where Σ is the constant, traceless, diagonal matrix and Q(x, t) is λ-independent and off-diagonal, namely Q jj = 0, while its off-diagonal entries Q jm (x, t), j = m, are generically six complex-valued functions of x and t; A is a constant matrix, while B and C can be written in terms of commutators and anti-commutators of Σ, Q(x, t) and their powers. Then, the compatibility condition leads to an evolution equation for the matrix Q in the form of a system of six coupled nonlinear wave equations for the entries of Q. Among other reductions, this system can be reduced to just two equations describing the resonant interaction of two waves of physical interest (see next section). Before looking at any particular reduction, let us first turn our attention to the way integrability plays an essential role in investigating the linear stability of a given solution Q(x, t). In order to make this paper self-contained, we first briefly report the main ingredients and features of the method (see [16] and references therein). The starting point is the following theorem (see [16] for a proof), which provides the quite remarkable and well-known property of the matrix solution Ψ (x, t, λ) of the two ordinary differential equations (2.1). (2.4) and let Ψ (x, t, λ) be a corresponding fundamental solution of (2.1), that is det Ψ = 0, then the matrix In other words, suppose that Q(x, t) is such a given solution, and that Q(x, t) + δQ(x, t) is a second solution of the same evolution equation. Suppose moreover that the entries of δQ(x, t) are sufficiently small (this is certainly so at small initial times t, if the initial condition δQ(x, 0) is small enough), then the linearized evolution equation for δQ(x, t) is obtained by neglecting all the nonlinear terms. This evolution equation takes the form

Theorem 2.1. Let Q(x, t) be a given solution of the compatibility equation
for any complex λ. This is a system of six linear partial differential equations; its solutions F(x, t, λ), which depend on the spectral variable λ, play the role of Fourier-like modes if the values of λ lie in an appropriate subset S of the complex plane, defined below as stability spectrum because of its key relevance to finding out whether a given solution Q(x, t) is linearly stable or unstable. This spectrum, which depends on the particular solution Q(x, t), is defined so as to guarantee the boundedness of each solution F(x, t, λ), at any fixed t, on the entire x-axis. By keeping in mind this analogy with the Fourier analysis of solutions of linear partial differential equations with constant coefficients, we assume that the set of solutions F(x, t, λ), see (2.5), of equation (2.7) allows the representation for an appropriate choice of M(λ), and in an appropriate functional space of initial data: for instance, we have in mind the physical assumption that the entries of the initial condition δQ(x, 0) be localized wave packets. Although the proof [16] of the theorem above requires only local properties (differentiability up to sufficiently high order in x), the integral representation (2.8) should be justified by spectral methods. To our present purpose, that is assessing the linear stability of Q(x, t), we do not need to enter into this matter, while we devote our next section to detail the properties of the stability spectrum S, for a particularly simple choice of Q(x, t).
Before going into this investigation, we give an explicit and general construction of the matrix eigenfunctions F(x, t, λ) by taking advantage of the properties of the wedge product of threedimensional vectors. Consider first the matrix solution Ψ of the Lax pair (2.1) and its three column vectors ψ (j) , j = 1, 2, 3: (2.9) Without loss of generality, we assume that this fundamental solution has unit determinant, say where Ψ A is the adjugate matrix (or classical adjoint) of Ψ . In turn, the matrix of co-factors Ψ A can be expressed through its rows, namely where the superscript T stands for transposition that takes a column vector ψ A(j) into a row vector, so that the vectors ψ A(j) are the columns of Ψ A T . The computation of Ψ A is then provided by the expression (2.11) with As we need six eigenfunctions for each value of λ if Q(x, t) is generic (non-reduced), we choose to compute the eigenfunctions (2.10) by setting M equal to one of the six basis matrices M (jm) , The convenience of this choice stems from the algebraic identity where jnm is the parity of the permutation {j, n, m} of {1, 2, 3}, that is jnm = 1 if {j, n, m} is a cyclic permutation of {1, 2, 3}, and jnm = −1 otherwise. Thus, the six eigenfunctions finally take the expression with the side remark that this explicit expression of the eigenfunction F (jm) cubically depends on just two vector solutions of the Lax pair. Finally, the Fourier-like representation (2.8) can be more precisely written as where the six functions μ (jm) (λ) play the role of Fourier-like transform of δQ.

Interaction of long and short continuous waves
In this section, we consider reductions of the matrix Q in (2.2) which lead, via the compatibility condition (2.4), to a system of two evolution equations able to model the LS resonant interaction. To this purpose, the matrix Q should depend only on two wave fields: the amplitudes L(x, t), modulating a continuous wave with a (infinitely) long wave length, and S(x, t), which modulates the amplitude of a periodic short wave. As a consequence, one should impose the condition that L(x, t) be real and dispersion-free, and that S(x, t) be instead complex with a non-vanishing dispersion coefficient (which, without any loss of generality, is rescaled to unity). As we pointed out in §1, here we introduce the new integrable model (1.5), which includes the YO and N systems (1.2) and (1.3), respectively, by setting This system is integrable for any real value of α and β; it reduces to the YO equation From this choice of the matrix Q, and from the compatibility condition (2.4), the coefficients A, B, C of the matrix T (see (2.2)) turn out to be In these formulae, the parameters α, β are arbitrary but constrained to be real. These parameters may be considered as independent constants which are responsible for the LS cross-interaction.
Here we are concerned with the linear stability of a given solution S(x, t), L(x, t) of the system (1.5), under the assumption that this solution is bounded for all values of x. According to standard practice, we proceed by investigating the time evolution of a small variation δQ(x, t), namely δS(x, t), δL(x, t), of this solution. At the early stage of the evolution, these variations are assumed to be sufficiently small so as to keep only their linear contributions. Thus, they satisfy the linearized approximate equation (2.6): Moreover, and for the sake of simplicity, we assume that the initial values δS(x, 0) and δL(x, 0) are localized and bounded functions of x. With these assumptions in mind, the solution S(x, t), remain small for all later times t > 0 or, equivalently because of the integral representation (2.16), if the eigenfunctions F (jm) (x, t, λ), see (2.15), do not grow exponentially in time.
As testing a recipe goes by tasting the cake, we apply our approach to the simplest physically relevant solution of the LS equation (1.5), namely the continuous wave solution which introduces three independent real parameters: the two amplitudes a and b, and the wavenumber q of the short wave. The starting point is the construction of a corresponding fundamental solutionΨ (x, t, λ) of the Lax equations (2.1). This is given by the expression and where the ( Because of the trace expressions tr(W) = q, tr(W 2 ) = ν + 3ρ, alternative solution Ψ (x, t, λ) whose column vectors, see (2.9), are Here the three constant vectors f (j) (λ) are the eigenvectors of the matrix W(λ) (3.7), which are generically linearly independent. Moreover, we normalize them so that with the implication that det(Ψ (x, t, λ)) = 1. This construction finally yields the expression of the column vectors of Ψ (x, t, λ), see (3.9), namely with We are now in the position to obtain the expression of the eigenfunctions F (jm) (x, t, λ), which correspond to the continuous wave solution (3.4). By inserting in (2.15) the vectors (3.12), and by taking into account the relation 14) and the matrix identity (for two arbitrary vectors u, v) we end up with the expression where jnm is the parity of the permutation {j, n, m} of {1, 2, 3}. This formula explicitly shows that, apart from the λ-independent phase θ , the eigenfunctions F (jm) (x, t, λ) depend on x and t only via the exponentials e i(η j -η m ) . By taking into account the expression of η j , see (3.13), these exponentials for j = m take the familiar expression e ±i(k n x−ω n t) , k n = w n+1 − w n+2 , ω n = w 2 n+1 − w 2 n+2 , n = 1, 2, 3 mod 3, (3.17) where the wavenumbers k j (λ) and their corresponding frequencies ω j (λ) are explicitly defined as in terms of the eigenvalues w j of the matrix W, see (3.10), or, equivalently, in terms of the roots of the characteristic polynomial where the two parameters p and r have the expression Because of the requirement that the basic solutions (3.16) be bounded functions of x, we conclude that the spectral representation of δL(x, t) and of δS(x, t) requires integrating with respect to the complex spectral variable λ, see (2.16), over the subset of the complex λ-plane where at least one of the wavenumbers k j (λ) is real. These considerations then lead to the following.
Definition 3.1. The stability spectrum S is defined as the set of the complex values of λ such that at least one of the three wavenumber functions k 1 (λ), k 2 (λ), k 3 (λ) is real. In this construction, we notice that the dispersion relation between the wavenumber k j and its corresponding frequency ω j is parametrically defined by the pair of functions k j (λ) and ω j (λ) by varying the parameter λ only over the stability spectrum S. It follows that the continuous wave solution (3.4) is linearly stable if ω j (λ) is real for those values of j such that k j (λ) is real for any λ over the entire spectrum S. On the contrary, if for some j and λ ∈ S, ω j (λ) is not real, then this solution is linearly unstable. In this latter case, the relevant physical information is provided by the computation of the gain function Γ j (λ) = |Im(ω j (λ))|, λ ∈ S, Im(k j (λ)) = 0. (3.21) Before proceeding to classify the stability spectra S and gain functions (3.21), it is worth drawing the reader's attention to the following fact. If the present approach to linear stability is applied to continuous wave solutions of the NLS equation (1.1), then the stability spectrum introduced above by definition 3.1 coincides with the Lax spectrum (which follows from the Lax differential operator in the space variable x via standard definition) [19,25]. However, this coincidence of the two spectra is a peculiar property of the 2 × 2 matrix Lax equations. In the case of N × N Lax equations with N ≥ 3, the two spectra, i.e. the Lax spectrum and the stability spectrum, as defined here by (3.1) for N = 3, are generically different from each other. For such larger matrices, only the spectrum S as defined above is relevant to stability.
We begin our analysis of the stability spectra S and of the corresponding gains by observing that the expression of the coefficients of the characteristic polynomial P(w, λ) (3.19) depends on the wavenumber q of the short wave solution (3.4), on the parameter r given in (3.20), and on the parameter p (see (3.20)) and the spectral variable λ only via the combination λ 2 − p. This shows that it is sufficient to fix the value of the parameters q and r in the polynomial (3.19) and to define the stability spectrum as a curve in the complex plane of the variable λ 2 − p. This simplifies our task as the parameter space is reduced to the (q, r) plane. Thus, hereafter, we introduce the alternative, and more convenient, complex variable Λ defined as Λ = λ 2 − p, (3.22) with the implication that the parameter p becomes irrelevant to our characterization and classification of spectra. In fact, by a minor abuse of notation, we refer to the characteristic polynomial (3.19) as In order to make this change of variable from λ to Λ explicit, we denote as S Λ the stability spectrum in the complex Λ-plane by adopting the following definition. It is obvious that this spectrum looks different from the spectrum S. Furthermore, we note that the parameter q, if non-zero, can be rescaled to the value q = 1 by rescaling w by q, Λ by q 2 and r by q 3 , or, equivalently, by the change of variables w → qw, Λ → q 2 Λ, r → q 3 r. (3.24) Nevertheless, we find it convenient to keep q in our formulae, to numerically set q = 1 whenever q = 0, and to separately treat the case q = 0. Let us consider now the part of the spectrum S Λ which lies on the real axis Im Λ = 0. In this case, all coefficients of the characteristic polynomial P(w, Λ) are real, and therefore either the three zeros w 1 (Λ), w 2 (Λ), w 3 (Λ) are real, or one is real and two are complex conjugate. In the first case, the three wavenumbers k j are real, while in the second case none of them is real, which leads to the following.
Thus, as shown by this expression (3.25), the large and positive real values of Λ do belong to S Λ , while the large and negative real values of Λ do not. More explicitly, we approximately solve the equation P(w, Λ) = 0 around the point at infinity of the complex Λ-plane and obtain the following asymptotic expressions of the three roots w j : where the labelling index j is arbitrary. Consequently, if Λ is real, large and positive then also the three wavenumbers k j (Λ), see (3.18), are real and large. If instead Λ is real and negative, and its modulus is large enough, no real wavenumber k j exists and Λ does not belong to the spectrum. If Λ is real, large and positive, the w-discriminant (3.25) is positive, and it remains so while moving Λ along the real axis towards the origin, until it reaches its first zero, say Λ + , of the discriminant w P(w, Λ); see (3.25). At this zero, this discriminant generically changes its sign and one of the three wavenumbers k j (Λ) vanishes. In turn, the w-discriminant (3.25) is a cubic polynomial of the variable Λ with real coefficients. Therefore, it has either one or three real Λ-zeros. In the former case, the discriminant w P(w, Λ) is negative for all values Λ < Λ + and hence, according to our proposition 3.3, the real part of the spectrum S Λ is the semi-axis Λ + ≤ Λ < +∞. If instead the discriminant w P(w, Λ) has three real zeros, say Λ 0 < Λ − < Λ + , then the interval Λ − < Λ < Λ + cannot belong to the spectrum since in this interval the discriminant w P(w, Λ) is negative (see proposition 3.3). Thus, this interval is a finite gap of the spectrum. Indeed, in this case the real part of the spectrum consists on the finite interval Λ 0 < Λ < Λ − where the discriminant w P(w, Λ) is positive, and the semi-axis Λ + ≤ Λ < +∞. This finite gap should be considered as a distinctive feature of the spectra in our classification. Indeed, the Λ-discriminant of the discriminant w P(w, Λ), that is Λ w P(w, Λ), depends only on the two parameters q, r and serves our classification purpose, as summarized by the following. Let us now turn our attention to the complex values of Λ that are not on the real axis and yet belong to the spectrum, Λ ∈ S Λ , Im(Λ) = 0. To this purpose, it is far more convenient to introduce the new polynomial

Proposition 3.4. Let Λ w P(w, Λ) be the Λ-discriminant of the discriminant (3.25), that is
whose roots ζ j (Λ), j = 1, 2, 3, are the squares of the differences of the Λ-dependent roots w j (Λ) of the polynomial P(w, Λ), see (3.18), The coefficients of the polynomial P(ζ , Λ) (3.28a) can be computed explicitly in terms of the coefficients of the polynomial P(w, Λ) (3.23) (see [21]), and read γ 2 = −2(3Λ + q 2 ), The definition (3.28a) shows that the w-discriminant (3.25) of the polynomial P(w, Λ) is simply related to this new polynomial (3.28a) as it reads w P(w, Λ) = −P(0, Λ). (3.29) This property allows one to rewrite our proposition 3.4 on gap characterization in terms of P(0, Λ). The spectrum S Λ itself can be redefined as the locus of the Λ-zeros of P(ζ , Λ) which correspond to non-negative values of the variable ζ , say ζ ≥ 0. Thus, to each ζ ≥ 0 there correspond three values of Λ which belong to S Λ . In particular, since this polynomial P(ζ , Λ) also takes the expression with real coefficients for real ζ , we conclude the following.

Proposition 3.5. The spectrum S Λ is symmetric with respect to the real axis.
Indeed, if Λ and Λ * are roots of P(ζ , Λ) for a given non-negative ζ , then both Λ and Λ * belong to the spectrum. Thus the spectrum S Λ is a symmetric piecewise smooth curve in the complex Λ-plane. The transition from a triplet of real Λ-roots of P(ζ , Λ) to a pair of complex conjugate Λ-roots and one real Λ-root (or vice versa) originates from a collision of two real (or two complex conjugate) Λ-roots, while changing the value of the variable ζ . This happens at a zero ζ j of the discriminant Q(ζ ) = Λ P(ζ , Λ), (3.31a) namely, if Q(ζ j ) = 0, provided this discriminant changes its sign. This discriminant Q(ζ ) turns out to be of fifth degree. However, it factorizes as After the collision where the discriminant Q(ζ ) changes its sign, the two Λ-roots scatter off the real axis. The change of sign makes the factor Q 2 1 (ζ ), see (3.31b), irrelevant to this analysis, so that hereafter we focus our attention on the factor Q 2 (ζ ) only. This polynomial has generically three simple roots which depend only on the parameters q, r, and may be either all three real, if the discriminant of Q 2 (ζ ) is positive, or two of them are complex conjugate and one is real, if this discriminant is instead negative, where μ, ρ and ζ 3 are real. Thus, our classification is mainly based on the function Q 2 (ζ ) and on its discriminant ζ Q 2 (ζ ). In this context, the necessary existence condition of a gap, which follows from our proposition 3.4, can also be established as Q 2 (0) = −ζ 1 ζ 2 ζ 3 = 4r(8q 3 − 27r) > 0. In the following classification, we find it more convenient to play with the sign of qr rather than with the sign of q and of r separately. Thus, for instance, the existence condition of a gap of the spectrum reads (qr)[8q 4 − 27(qr)] > 0 which is never satisfied if qr < 0 and it is satisfied only for 0 < qr < (8/27)q 4 . Consider first a spectrum S Λ whose parameters q, r satisfy the positive discriminant inequality (3.33a), namely (qr)[16q 4 − 27(qr)] > 0. In this case, for large enough values of positive ζ , namely for ζ 3 < ζ < +∞, as well as for ζ 1 < ζ < ζ 2 , the function Q 2 (ζ ) is positive, Q 2 (ζ ) > 0, while Q 2 (ζ ) < 0 for ζ 2 < ζ < ζ 3 and for −∞ < ζ < ζ 1 . In order to look for distinctive features of the spectrum, we have to analyse the following four different cases separately according to the positiveness of the three ζ -roots of Q 2 (ζ ).
In this case, two Λ-roots collide for ζ = ζ 3 and, for 0 ≤ ζ < ζ 3 , the spectrum acquires two complex conjugate curves in the Λ-plane. We refer to this complex part of the spectrum as branch (B), since the two endpoints at ζ = 0 of the two Λ-roots trajectories do not generically coincide with each other (figure 1a). No gap is possible and we term this spectrum type 0G 1B 0L, or of B-type. (iii) ζ 1 < 0 < ζ 2 < ζ 3 : the same Vieta relation used above necessarily requires 0 < qr < (8/27)q 4 , which is also compatible with the positive discriminant condition. Two Λ-roots collide for ζ = ζ 3 , get off the real axis and collide again for ζ = ζ 2 on the real Λ-axis thereby forming one complex closed curve, which we term loop (L) (figure 2a). In this case, there is no branch. However, a gap exists because its existence condition (qr)[8q 4 − 27(qr)] > 0, see our proposition 3.4, is satisfied. Thus, we call this spectrum of type 1G 0B 1L, or of LG-type. (iv) 0 < ζ 1 < ζ 2 < ζ 3 : this case is not allowed by two Vieta relations, which lead to the inequalities 6q 4  Let us consider now those spectra S Λ whose parameters q, r satisfy the negative discriminant inequality (3.33b) instead, namely (qr)[16q 4 − 27(qr)] < 0. Again we distinguish all possible cases according to the ζ -roots (3.33b).
These observations above conclude our characterization and classification of spectra, which are summarized in table 1. As stated after our definition 3.2 of the spectrum S Λ , see also (3.24), the classification as in table 1 has been derived with the assumption that q = 0. For the sake of completeness, we now remove this limitation by computing the spectrum for q = 0 with r = 0 (see below for r = 0). We first observe that in this special case the polynomial P(ζ , Λ) reads complex part of the spectrum, we note that Q 2 (ζ ) = ζ 3 − 108r 2 , see (3.31b), and therefore the only real zero of Q 2 (ζ ) is the positive number ζ 3 = (108r 2 ) 1/3 . For this value of ζ , the polynomial P(ζ , Λ) must have a double real Λ-root Λ B , which is found to be Λ B = (1/2)ζ 3 = (1/2)(108r 2 ) 1/3 = 2 1/3 Λ 0 .
As ζ is moved from ∞ towards the origin, two real Λ-roots collide with each other at Λ B for ζ = ζ 3 to give rise to a branch, whose off real parts terminate for ζ = 0 at the complex conjugated points Λ 0 exp(2iπ/3), Λ 0 exp(−2iπ/3). Thus, for q = 0, r = 0 the spectrum is of type 0G 1B 0L. A similar analysis can be carried out at any zero of discriminants of cubic polynomials, in particular to discuss threshold phenomena. Indeed, the inequalities shown in table 1 divide the plane (q, r) into regions with different types of spectra. At the boundaries of these regions, the inequalities in table 1 turn into equalities which characterize thresholds. One simple and interesting instance of such thresholds occurs at the special value r = 0. In this special case, the spectrum S Λ is readily found because the roots of the characteristic polynomial (3.23) P(w, Λ) = (w − q)(w 2 − Λ) are explicit. The corresponding spectrum is entirely real and we will call it of type 0G 0B 0L. Indeed, the wavenumber functions k j (Λ) have the simple expression 1L where however the gap is reduced to just one point.
Let us now turn our attention to the time variable t and to the stability of the plane wave solution (3.4). This point requires investigating whether the frequencies ω 1 (Λ), ω 2 (Λ), ω 3 (Λ), see (3.18), are real numbers for Λ ∈ S Λ . Our starting point is their relation to the wavenumbers k j (Λ) This relation (3.36) is clearly not a dispersion relation since each frequency ω j is written in terms of two wavenumbers, k j and k j+1 . In order to compute the dispersion relation, we need to eliminate the variable Λ among the algebraic relations P(w, Λ) = 0, P(ζ , Λ) = 0 and ω j = k j (q − w j ). This requires a tedious but straightforward computation aimed to factorize high degree polynomials. This is done with the help of the additional algebraic relation R(ξ , Λ) = 0 with the following specifications: where the coefficients δ n (Λ) can be explicitly computed as in [16], and have the following expression:  with the characteristic polynomial of the matrix W 2 , that is The upshot of this calculation is the three-branch dispersion relation H(±k j , ±ω j ) = 0, where H(k, ω) is the polynomial It is worth noticing that this dispersion relation H(k, ω) = 0 can be obtained by means of a standard Fourier approach to the linearized equations. In other words, the continuous wave solution (3.4) is so special that its stability properties can be investigated by both Fourier expansion and our Lax pair method. However, since the Fourier expansion approach works only for this particular solution, we have turned our focus on the integrability methods in order to explicitly apply the mathematical formalism (see §2) in a simple context. The next, and harder, task will be to extend the present analysis to other known solutions, such as solitary waves, say solitons, periodic and rogue waves. Let us now consider the linear stability of the plane wave solution (3.4) for generic values of q and r. Our analysis has shown that the Λ-spectra have a real component and a complex, off real one. The former lies on the real axis, and at each point all three wavenumbers k 1 (Λ), k 2 (Λ), k 3   do not cause instabilities. On the contrary, if Λ is off the real axis, namely on a branch or on a loop, which are the complex components of S Λ , only one of the three wavenumbers, for instance k 3 (Λ), is real, while the other two k 1 (Λ), k 2 (Λ) are non-real. On branches and loops, while k 3 (Λ) is real, the corresponding frequency ω 3 (Λ) is necessarily complex with a non-vanishing imaginary part. This follows from the expression (3.36). Precisely, since ω 3 = (1/3)k 3 (2q + k 3 + 2k 1 ), the real frequency is The corresponding instability is characterized by the gain function (3.21) which is non-vanishing as far as Λ remains on the non-real component of the spectrum. We conclude this analysis with the following. Examples of these two types of spectra have been numerically computed, by evaluating the Λ-roots of P(ζ , Λ) in (3.34), as ζ varies from 0 to the largest, positive, real ζ -root of Q 2 (ζ ) in (3.31b), as detailed in appendix C of [16]. In figure 1, we show the B-type spectrum (figure 1a); the corresponding gain function Γ on the branch (figure 1b), which proves that this instability is of baseband type, namely waves are unstable |k| = 0; the functions k j (Λ), for j = 1, 2, 3, if Λ is real, together with the function k 3 (Re(Λ)) on the branch (figure 1c); the real frequency on the branch as a function of k = k 3 (figure 1d). The same functions are plotted in figure 2b-d, for an LG-type spectrum, which is shown in figure 2a. In particular, figure 2b shows that in this case the instability is of passband type, namely, waves are stable for sufficiently small values of |k|.

Summary and conclusion
Research on instabilities of nonlinear waves has witnessed a renewed interest in recent years, following the study of the modulational instability of continuous wave solutions of the focusing NLS equation. We exploited integrability to construct a method to predict whether a nonlinear wave, described by an integrable, nonlinear system, is linearly stable against small perturbations. Here, we displayed how this approach works in the resonant coupling of long and short waves. As a side result, we first introduce a new model equation which combines two very well-known integrable models, namely the YO equation (1.2) and the Newell equation (1.3), into the more general YON equation (1.5), which features two arbitrary, real parameters, and thus it is likely to fit a broader range of physical contexts. This outcome is similar to the one that proves that the Korteweg-de Vries (KdV) and the modified-KdV equations are just two particular cases of the Gardner equation [26,27].
A common feature of all these LS equations is that a long wave always arises as generated by short waves. This process, which is guessed by observing the evolution equations (1.5) themselves, is revealed by the first two conservation laws ρ nt = g nx , n = 1, 2, where ρ n is the density and g n is the current. The first conservation law is the long wave equation itself, with ρ 1 = L and g 1 = 2|S| 2 . The second conservation law turns out to have the conserved density ρ 2 = |S| 2 − (α/2)L 2 and the current g 2 = 2 Im(SS * x ) − 2αL|S| 2 . In this work, we give the 3 × 3 Lax pair associated with the new YON model and display how this Lax pair allows us to construct the basic solutions, or eigenfunctions, of the linearized model equations. These computations are specialized to the continuous wave (or plane-wave) solutions (3.4) of the YON equations. In this case, the explicit construction of the eigenfunctions of the linearized wave equations allows us to define, and compute, their associated wavenumbers k and frequencies ω. Not only are these physical quantities related to each other by an explicit or implicit dispersion relation (3.40), as in standard Fourier analysis, but they are also parametrically related as k = k(Λ), ω = ω(Λ), where the variable Λ lives in the complex plane. In turn, the complex variable Λ is simply connected to the spectral variable λ, which appears in the Lax pair. The key point is the definition of the stability spectrum, an algebraic curve whose geometrical properties convey relevant physical information, and which generically differs from the Lax spectrum. Eigenfunctions, wavenumbers and frequencies are defined on the stability spectrum. In particular, on the real Λ part of the spectrum the eigenfunctions are bounded functions of time (stability), while on the non-real part of the spectrum the eigenfunctions exponentially grow with time (instability). This strictly complex part of the spectrum brings its own information: it displays one open branch, if the instability is of baseband type, or a closed loop, if the instability is instead of passband type.
All the spectra are classified by their dependence on essentially one real parameter r; see (3.20). This is a simple function of the amplitudes a and b of the short and, respectively, long waves of the wavenumber q of the short wave, and of the self-and cross-coupling constants α and β of the model. For a generic choice of all these parameters, only two types of spectra, and therefore of perturbed wave behaviours, emerge from our analysis. Plots illustrate these two typical stability spectra. Their properties show that all values of the parameters lead to instabilities, with the only exception of the case r = 0, which instead guarantees stability since the associated spectrum has