The 3-wave resonant interaction model: spectra and instabilities of plane waves

The three wave resonant interaction model (3WRI) is a non-dispersive system with quadratic coupling between the components that finds application in many areas, including nonlinear optics, fluids and plasma physics. Using its integrability, and in particular its Lax Pair representation, we carry out the linear stability analysis of the plane wave solutions interacting under resonant conditions when they are perturbed via localised perturbations. A topological classification of the so-called stability spectra is provided with respect to the physical parameters appearing both in the system itself and in its plane wave solution. Alongside the stability spectra, we compute the corresponding gain function, from which we deduce that this system is linearly unstable for any generic choice of the physical parameters. In addition to stability spectra of the same kind observed in the system of two coupled nonlinear Schrödinger equations, whose non-vanishing gain functions detect the occurrence of the modulational instability, the stability spectra of the 3WRI system possess new topological components, whose associated gain functions are different from those characterising the modulational instability. By drawing on a recent link between modulational instability and the occurrence of rogue waves, we speculate that linear instability of baseband-type can be a necessary condition for the onset of rogue wave types in the 3WRI system, thus providing a tool to predict the subsequent nonlinear evolution of the perturbation.


Introduction
The three wave resonant interaction (3WRI) system is a well known universal model which describes the non-dispersive propagation of waves in weakly nonlinear media [1,2].Despite its dispersionless nature, it possesses solitons (see, e.g., [3][4][5][6] and references therein), which occur as result of the balance between nonlinearity and the energy flow due to the mismatch of the velocities of the three interacting waves (see, e.g., [1,7,8]), and it is not due to the balance between nonlinearity and dispersion as it happens for the 'traditional' solitons [9].Such solitons are ubiquitous in many nonlinear systems of physical interest (e.g. in nonlinear optics [10]) and they have been observed experimentally in stimulated Raman scattering [11], in Brillouin finger ring laser [12] and in quadratic optical media [13].Among its soliton solutions, rational and semi-rational solutions [14][15][16][17] have enjoyed renewed interest as they can be used to model rogue wave-type phenomena [18,19], recently observed in oceans [20], in optics [21], in the atmosphere [22], in Bose-Einstein condensates [23], and in capillary-waves [24].They were also observed in laboratory experiments with superfluid helium [25] and generated in water tanks [26].In dispersive systems, modulational instability (MI) has been proposed as a possible cause of the onset of rogue waves [27,28]; however, not all kinds of MIs have been considered, but only the so-called baseband MI [27,28], that is, a fundamental feature of nonlinear dispersive systems for which a periodic perturbation on an unstable plane wave or continuous wave background (CW) grows when the perturbation frequency is around zero (zero not included) (see, for instance, [29]).In [14,17], novel families of rational and semi-rational solutions on non-vanishing background of the integrable 3WRI system, which is instead non-dispersive, have been obtained via Darboux Dressing transformations (see, e.g., [30,31]), and they have contributed to our understanding of extreme phenomena in multi-component resonant processes.
In this work, we carry out the linear stability analysis of plane wave solutions of the 3WRI model (1), motivated by the fact that this resonant process is observed in different physical contexts such as in plasma wave interactions [32,33], in nonlinear dielectric media [34], in planetary wave interactions [35], in shear flows [36] and internal waves [37].The non-dispersive nature of the 3WRI system makes the linear stability analysis of its solutions more intriguing.Indeed, we expect the mechanism leading to the linear instability in the 3WRI system to differ from the one associated with MI in dispersive nonlinear systems, the latter being due to the interplay between dispersion and a weak nonlinearity [38].In particular, in dispersive multi-component systems, the coupling among the components can suppress or enhance the MI of the single component [39][40][41] and the dispersion plays a key-role in the occurrence of MI.For instance, it has been shown that two coupled sine-Gordon equations in the non-dispersive regime are stable with respect to MI [39].The 3WRI system admits also dispersive shock waves [42], similar to those admitted by well known dispersive systems such as the nonlinear Schrödinger (NLS) and the Korteweg-de Vries (KdV) equation.Despite the relevance of the 3WRI system, the linear stability analysis of its plane wave solutions has not received much attention in the literature, to the best of our knowledge.This is possibly due to algebraic difficulties in dealing with a multi-component system and solutions with nonzero boundary conditions.In the paper, we follow [43] and adopt the approach presented in [44,45] exploiting the integrability of the system, in particular its Lax Pair formulation, to study the linear stability of its plane waves solution and classify the so-called stability spectra and associated gain functions.We provide a complete topological classification in the physical parameters space.Moreover, we observe that all topologies of the stability spectra of a system of two coupled nonlinear Schrödinger equations (CNLS) detailed in [44,46], as well as those of a long wave-short wave interaction system (YON) detailed in [45], are also present for the 3WRI system.Indeed, there exist open curves (branches) associated with an MI-baseband-type instability and closed curves (loops) associated with an MI-passband-type instability.However, there exist new topological features (twisted loops) in the stability spectra of the 3WRI model which do not exist for the CNLS system.They are associated with an isolated real point of the stability spectrum and, at this point, the solutions are always linearly unstable; on the contrary, the CNLS system is linearly stable for any real values of the spectral parameter outside the gaps [44,46].Moreover, the gain function associated with these twisted loops is neither baseband-type nor passband-type, being nonzero in the whole neighbourhood of the zero wave number.Like in the CNLS case, also for the 3WRI system, we show that for a generic choice of physical parameters, the plane waves solutions are linearly unstable.In Sect.2, we recall the Lax pair formulation of the 3WRI model and its plane wave solution.Here, we briefly describe also the relation between the linearised equation and the squared eigenfunctions, following [44], and we introduce the concept of stability spectrum associated with the plane wave solution.In Sect.3, we identify regions in the parameters space associated with different topologies of the stability spectrum and we construct the corresponding gain functions.We deduce that the plane wave solution of the 3WRI system is linearly unstable for any generic choice of the physical parameters.In the Conclusions, we summarise our results, and in particular we highlight the link between instabilities of baseband-type and the occurrence of rogue waves phenomena, proving a spectral characterisation of rogue waves type solutions, already known in the literature.

Stability spectra associated with the plane wave solutions
The propagation in 1 + 1 dimensions of three resonating waves is modelled by the following system of three coupled partial differential equations (PDEs): where the asterisk denotes the complex conjugate, u j = u j (x, t) are complex valued functions of space x and time t, s j are signs such that s 2 j = 1 for j = 1, 2, and c 1 = c 2 are the nonzero velocities associated to u 1 and u 2 , respectively.Indeed, without loss of generality, we have set the velocity associated with u 3 to zero, that is equivalent to choosing a reference frame co-moving with the solution u 3 .The ordering of the velocities is otherwise generic at this stage.Hereinafter, subscripted variables stand for partial differentiation with respect to x or t.It is worth highlighting that (1) might look different from the one usually given in the literature (see e.g.[3] and [42]), however it contains all the possible three waves resonant interaction systems, irrespectively of the choice of signs and ordering of the velocities.This allows us to study a single Lax Pair (see (7)), instead of dealing with at least three different Lax Pairs (e.g.[3], [42]) depending on the choice of the ordering of the velocities, this being advantageous from the point of view of the classification of the stability properties of the plane waves solutions.System (1) models the interaction of three resonating waves and, as such, it has been the subject of extensive investigation (see, for instance, [3]).It can be derived via multiscale analysis of a large class of nonlinear wave equations with weak dispersion and nonlinearity (see, for instance, [4]) and therefore it represents a universal model able to capture the lowest order corrections, due to nonlinear effects, to linear propagation.We have already pointed out that, contrary to other integrable nonlinear wave equations such as the nonlinear Schrödinger equation, the Korteweg-de Vries equation or the sine-Gordon equation, this system has no (linear) dispersion.The second peculiar characteristic is that no self-interaction occurs, namely, the forcing of each wave (i.e. the non-homogeneous term in the propagation Eq. ( 1)) is just the product of the (complex conjugated) amplitudes of the other two waves.The plane wave solutions of (1) can be written as with a j constant amplitudes, η j frequencies and ν j wave numbers.Frequencies and wave numbers are related one another via the resonant conditions Thus, the solution u 3 reads and the nonlinear dispersion relations are System (1) is integrable, as it admits a Lax Pair (e.g., see, [3]) namely, it can be written as the compatibility condition of two linear differential equations for an unknown auxiliary matrix-valued function ψ: where X and T are matrix-valued functions depending on u(x, t) and on a complex variable κ, called spectral parameter, as follows: Here Note that, U and V are related as follows: where the notation [M 1 , M 2 ] stands for the commutator In this setup, for all κ ∈ C, Eq. ( 1) is retrieved by imposing the compatibility condition ψ xt = ψ tx for system (6): A perturbation of u j (x, t) in the form u j + δu j corresponds to a perturbation of X and T in the form X + δX and T + δT where, δX = δX(x, t) and δT = δT (x, t) are solutions of the linearised equation Equation ( 13) is the evolution equation for the perturbation δU of the matrix U , that is If δU is a localised perturbation, then, as proven in [44], a solution of the linearised Eq. ( 14) reads with where M = M (κ) is an arbitrary, nonzero constant matrix.The matrix Ψ, which is essentially a squared eigenfuction in the spirit of [47], for any given arbitrary M , satisfies the linear ordinary differential equations (ODEs): which are compatible with one another because of (12).The solution F plays the same role as the exponential solution of linear equations with constant coefficients; namely, by varying κ, it provides the set of "Fourier-like" modes of the linear PDE (13).As a consequence, any sum and/or integral of F over the spectral parameter κ is a solution δU of the linearised Eq. ( 13).As in [44], here we assume that the perturbation δU has the integral representation In order to explicitly compute F , we need first to compute the solution of the Lax Eq. ( 6).At this stage, it is convenient to set this can be done without loss of generality 1 .In this setup, the fundamental matrix solution of the Lax pair reads where the matrix R(x, t) is given by with I the 3 × 3 identity matrix and and where W = W (κ) and Z = Z(κ) are two x, t-independent matrices and with the property [W, Z] = 0, so that W (κ) and Z(κ) can be diagonalised by the same matrix G(κ), namely In terms of W , Z, R and M , the solution F of the linearised (14) reads for each value of κ.Then, as shown in [46], in order to make ( 27) explicit, we write M (κ) as a linear combination of the nine matrices with where δ jk is the Kronecker delta (δ jk = 1 if j = k and δ jk = 0 otherwise), in the form: where μ jm are (arbitrary) scalar functions.Plugging this decomposition into (27), we obtain the following representation of where w j = w j (κ) and z j = z j (κ), j = 1, 2, 3, are the eigenvalues of W and Z, respectively, and they are assumed to be simple, for a generic value of κ, and where F (jm) (κ) are the x, t-independent matrices defined as Requiring the solution to be localised in the variable x implies (see [44]) the necessary condition that the functions μ jm (κ) be vanishing for j = m, μ jj = 0, j = 1, 2, 3.
Hereafter, it is convenient to introduce the following re-scaled wave amplitudes α j and the re-scaled spectral parameter λ along with a re-scaling of w-and z-roots, respectively, as w → q 2 w and z → q z, and setting q = 1 with effect of removing the dependence on q from the roots of the characteristic polynomials of W and Z.
Note that we are assuming here both c 1 and c 2 different from zero; the cases where at least one of the velocities is zero is a degenerate case (see [48] for details).Moreover, let p 1 , p 2 and p 3 be real parameters defined as where s j , c 1 and c 2 are as in (1), entailing p 3 = 0, ±1.In the following, we will assume c 1 + c 2 = 0, and postpone the study of the counter-propagating case p 3 → ∞ to future investigation.In the new parameters, the characteristic polynomial P W (w; λ) of W reads whereas the characteristic polynomial P Z (z; λ) of Z reads The differences in (31), play the roles of "eigen-wave numbers" and "eigen-frequencies", respectively.We observe that the matrices W and Z are related to one another as follows: This induces a relation between the eigenvalues w j and z j so that Given p 1 , p 2 and p 3 , we are interested in finding values of the complex spectral parameter λ such that the locally perturbed plane waves u j are bounded in space, and see if they are linearly stable or unstable in time.In other words, we look for those values of λ corresponding to real eigen-wave numbers k j ; to such eigen-wave numbers there correspond in general complex eigen-frequencies ω j : if the imaginary part of ω j is nonzero, linear instability develops in time.
The stability spectrum S for the plane wave solutions of the 3WRI system (1) is defined as the locus of the λ-plane identified with C such that, for fixed values of the physical parameters p 1 , p 2 and p 3 , there exist at least two eigenvalues w and w m for the matrix W , for some and m, for which (w − w m ) ∈ R (namely, for which one of the k j is real).
The key-concept in understanding the linear stability is the difference of the eigenvalues w j , which is related to the difference of the eigenvalues z j via (39).Indeed, it is possible to show [48] that Proposition 1.The plane wave solution of the 3WRI system is stable against the perturbation δU integrated over strictly real values of λ (see 18), with the exception of a finite number of intervals on the real line, called gaps, and possibly of isolated, real points within such gaps.
Following [44], we construct the cubic polynomial P W (ξ; λ) ≡ P W (ξ; λ; p 1 , p 2 , p 3 ) in the variable ξ, whose ξ-roots are the squares of all the possible differences of the roots of the characteristic polynomial P W (w; λ) (35), namely A ξ-root (resp.a λ-root) of P W (ξ; λ) is a polynomial root of the equation P W (ξ; λ) = 0, solved with respect to ξ (resp.λ).For the sake of simplicity, we will refer to P W (ξ; λ) as the polynomial of the squares of the differences; it is a 2-variate polynomial in ξ and λ, which can be explicitly obtained in terms of the coefficients of P W (w; λ) (see [44,48]), and takes the expression: Therefore, for any fixed parameters p 1 , p 2 and p 3 , the spectrum S is the locus of the λ-roots of P W (ξ; λ) for all ξ ∈ R, ξ ≥ 0, that is The stability spectrum S of the solution U coincides with the subset of the complex λ-plane over which the integration (18) runs.This subset guarantees the boundedness of the solution F , and hence that of the perturbation δU .We are now in the position to redefine the perturbation in (18) as As in standard linear stability analysis, the given solution U is linearly stable if the initial perturbation δU (x, t 0 ) remains small as time grows, say for t > t 0 .
In the following section, we classify the stability spectra S in the parameter space.

Classification of stability spectra S
As we will illustrate in this section, the topological features of the spectrum S can be classified as follows.
Real values of the spectral parameter λ not belonging to S constitute a gap (G).A gap, containing an isolated point belonging to S, is called a split gap (SG).Non-real values of the spectral parameter λ, belonging to S, correspond to branches (B) and loops (L), which are open and closed curves, respectively.Figure-eight shaped loops, self-intersecting on the real axis, will be referred to as twisted loops (TL).
In order to classify the topological features of S, we need to investigate the algebraic properties of the polynomial of the squares of the differences P W (ξ; λ), defined in (40), which is a sixth degree polynomial in λ and a third degree polynomial in ξ.For fixed p 1 , p 2 and p 3 , we analyse the spectrum S following the approach introduced in [44,46] for the CNLS equations, and further investigated in [48] for the case of the 3WRI system.As explained in [44,46], the different components of the stability spectrum S result from the dynamics of the λ-roots of P W (ξ; λ) as ξ varies from 0 to ∞; in particular, the components are generated by collisions on the real axis of pairs of λ-roots, whenever after such collisions the colliding roots either leave or go into the real axis.In other words, collisions corresponding to real roots passing each other on the real axis are disregarded, as they do not originate additional components of the curve S.
To analyse these collisions, we consider the discriminant of P W (ξ; λ) with respect to λ, that is where ) are two real polynomials in the variable ξ, whose degrees are four and six, respectively, and Q 1,2 (0) = 0.The explicit forms of Q 1 (ξ) and Q 2 (ξ) are too long for the present paper and can be found in [48].We will refer to Q 1 and Q 2 as the even and odd parts of the discriminant Q(ξ), respectively.We are interested in changes of sign of Q(ξ; p 1 , p 2 , p 3 ).
If ξ > 0, the factor ξQ 2 1 (ξ) is strictly positive and so, to the purpose of the study of the sign of the discriminant Q(ξ), the only relevant part is Q 2 (ξ).As ξ varies towards ∞, Q 2 (ξ) will change sign each time ξ goes through a root of Q 2 (ξ): therefore, to each choice of p 1 , p 2 and p 3 , there corresponds a sequence of intervals of ξ, where Q 2 (ξ) is either positive or negative.In particular, the following proposition holds [44,48]: Proposition 2. Let ξ be the largest positive root of Q 2 (ξ).Then, for ξ > ξ, all the λ-roots of P W (ξ; λ) are real; therefore, the stability spectra always contains part of the real axis and never features a gap containing the point at infinity.
As p 1 , p 2 and p 3 vary, transitions from a sequence of intervals in ξ to another correspond to one of the intervals shrinking to a point, namely to a collision of a pair of roots of Q 2 (ξ).We therefore consider the discriminant of Q 2 (ξ) with respect to ξ, Δ ξ Q 2 (ξ), which factorises as where E j are polynomial factors and e j are their multiplicities.Then, the real-analytic varieties E j implicitly defined as are curves in the parameter space, defining the boundaries of regions corresponding to different topological structures of the spectra S.These are and the explicitly expression of E 5 can be obtained from ( 45) knowing E j , j = 1, ..., 4 (see [48]).All the multiplicities e j are equal to 1, with the only exception of e 3 = 3.For ξ = 0, Q(0) = 0 and therefore it requires a different approach; we will study this case in the next subsection.

Gaps, branches and loops
As shown in [44,46], end points of gaps on the real line and of branches in the complex λ-plane correspond to setting ξ = 0 in P W (ξ; λ), or equivalently, to multiple w-roots of the characteristic polynomial P W (w; λ).When ξ = 0, it turns out that where P Z (ζ; λ) is the polynomial of the squares of the differences of the z-roots of P Z (z; λ), and where Gaps and branches appear or disappear at the multiple-zeros of the polynomial P Z (0; λ), thus at the zeros of the discriminant Δ λ P Z (0; λ), namely when Indeed, as we will see later in this section, R(λ) does not play a role in the classification of gaps and branches, but it provides conditions for the existence of split-gaps.The curves E j associated with the three polynomial factors E j in (51), with j = 1, 2, 3, bound the regions in the (p 1 , p 2 )-plane characterised by different numbers of gaps and branches.Following this observation, we will classify the spectra fixing the parameter p 3 and varying p 1 , p 2 (this is trivial as long as one is interested only in E 1 , E 2 and E 3 , which do not depend on p 3 , but it informs our approach in view of the full classification, when E 4 and E 5 too are taken into account).In particular, we observe that the curves E 1 , E 2 and E 3 are the same that appear in the classification of gaps and branches of the CNLS system [44,46].Consequently, as in [44,46], the following proposition holds: The structure of the gaps depends only on p 1 and p 2 , see (47a), (47b) and (47c).Nevertheless, because of the relation (48), P W (0; λ) can be zero also when R(λ) = 0, namely when This corresponds to a double-root of P W (0; λ).Since, R(λ) is a first degree polynomial and it has real coefficients, the equality R(λ) = 0 may be verified just in one point of the real spectrum S.Moreover, for values of p 1 and p 2 varying in an interval without changing the number of gaps, this point can move inside a gap or can coincide with the endpoint of a gap, depending on p 3 .In other words, we expect that the regions in the (p 1 , p 2 )-plane in which the zero of R(λ) is inside a gap move by varying p 3 .A tedious but easy calculation reveals that when R(λ) = 0, P Z (0; λ) can be either positive or zero (see 49), ( 50); we have only two possible scenarios: 1. R(λ) = 0 with P Z (0; λ) > 0: the double-zero (52) falls within a gap, namely, the spectrum features a split-gap; 2. R(λ) = 0 with P Z (0; λ) = 0: the double-zero( 52) is a triple-zero of the polynomial P W (0; λ) and it coincides with the endpoint of a gap.Thus, the condition for existence of a split-gap can be obtained by studying when the resultant with respect to λ of P Z (0; λ) and R(λ) is zero, namely Res λ (P Z (0; λ), R(λ)) = 0.The expression of this resultant contains two polynomial factors: a squared polynomial, which does not change sign and does not correspond to any transition, while the other coincides with −E 4 (see (47d)).Consequently, we have Proposition 4. Fixing the parameter p 3 , the curve E 4 (47d), identifies in the (p 1 , p 2 )-plane the transition curve for the existence of split gap.In particular, fixing the values of p 3 , the values of the parameters p 1 and p 2 for which the polynomial E 4 is negative correspond to regions where there is a split gap.
Proof.Recall the expressions of P Z (0; λ) and R(λ) in ( 49) and ( 50), respectively, and let τ j and θ be their roots, so that the resultant with respect to λ can be written as Let us suppose, say, τ j are all real and distinct, such that we have two gaps; without loss of generality we assume τ 1 < τ 2 < τ 3 < τ 4 .If θ coincides with the end point of a gap, so that θ = τ j for some j, then the resultant (53) vanishes.If, instead, θ falls within a gap, that is , namely if we have a split-gap, then the resultant (53) is positive.Since the resultant can be factorised as a squared polynomial times −E 4 , this implies E 4 < 0. Finally, if θ falls between two gaps, that is τ 1 < τ 2 < θ < τ 3 < τ 4 , then the resultant in ( 53) is negative, and E 4 > 0.
The same argument can be straightforwardly extended to the case of two distinct real roots and of two complex conjugate roots.
By taking into account the regions associated with split-gaps, one has the following Proposition.
Proposition 5. Fixing the parameter p 3 , split-gaps exist in the regions defined as follows: where To classify gaps as p 3 changes, we take advantage of the symmetries in the (p 1 , p 2 )-plane.Indeed, since Δ λ P Z (0; λ) in ( 51) is invariant under transformations p 1 → −p 1 and p 2 → p 2 , the curves E 1 , E 2 and E 3 are symmetric with respect to the p 2 -axis.The curve E 4 depends also on p 3 , and it has symmetries p 1 → −p 1 , p 2 → p 2 and p 3 → −p 3 .In particular, by changing p 3 into −p 3 , the curve E 4 is reflected with respect to the p 2 -axis.Therefore, following Proposition 5, in order to understand the regions in the (p 1 , p 2 )-plane associated with different numbers of gaps, branches and split gaps, it is sufficient to consider the two cases, |p 3 | > 1 and 0 < |p 3 | < 1 (see Fig. 1).So far we have considered the case ξ = 0 in which at least two w-roots of the characteristic polynomial P (w; λ) coincide.Considering the λ-roots of the polynomial P Z (0; λ), this case can be summarised as follows: 4 distinct real roots, or 2 distinct real roots and 1 double real root (2G 0B or 1G 1SG 0B); 2 distinct real roots and 2 complex conjugate roots, or 1 double root and 2 complex conjugate roots (1G 1B or 0G 1SG 1B); 2 pairs of complex conjugate roots (0G 2B), see Proposition 3. Let us consider the polynomial P W (ξ; λ) for ξ > 0. As we have observed at the beginning of Sect.3, this case can be studied by analysing Q 2 (ξ), requiring also the inclusion of the curve E 5 in the (p 1 , p 2 ), for a fixed p 3 , when defining the bounding regions of different topology.In particular, the only additional features of the spectrum S that form for ξ > 0 are closed continues curves that start and end on the real axis.They come in two different kinds: loops and twisted loops.Loops are generated by a pair of real λ-roots, which collide for some value of ξ > 0 and scatter into the complex plane as a pair of complex conjugate roots, as ξ increases, to collide again on the real axis for another, larger, value of ξ.Twisted loops are figure-eight shaped loops, self-intersecting on the real axis.The self intersection point of the twisted loop corresponds to the double λ-root (52) and thus corresponds to the isolated point within a split-gap (subfigure on the left of Fig. 4, Fig. 5 and subfigure on the left of Fig. 6).
Adopting the same approach as in [46], the following relations linking the number of branches, loops and twisted loops can be derived: where ξ + , ξ − and ξ cc stand for the positive, negative and complex conjugate ξ-roots of the polynomial Q 2 (ξ), and where the symbol # stands for 'number of'.Moreover, since Q 2 (ξ) is a sixth degree polynomial, the relation #ξ + + #ξ − + #ξ cc = 6 must be verified.This implies that, because of (54) and keeping in mind that #T L = #SG, we have also

Gain function
The solution U is linearly stable if the initial perturbation δU (x, t 0 ) remains small as time grows, say for t > t 0 .The expression of the perturbation δU is given via the matrix F in (43), where F is a linear This plot is a magnification of the 0 G 2B 2 L region displayed in the previous figure, which would otherwise be too small to be properly visualised and labelled combination of the "Fourier-like modes e i(xkj −tωj ) .Let w 1 and w 2 denote the two eigenvalues of W such that their difference k 3 = w 1 − w 2 is real.For λ ∈ S, let k 3 be the real eigen-wave-number and ω 3 is the corresponding eigen-frequency.Since ω 3 (k 3 ) is in general complex, its imaginary part represents the growth rate of the perturbation, which we will denote as Γ(k 3 ).Then, if ω 3 (k 3 ) is strictly real, that is Γ(k 3 ) = 0 for all k 3 real, the system is linearly stable.On the contrary, if ω 3 (k 3 ) is non-real, that is Γ(k 3 ) = 0 for some k 3 , then the solution F , and so the perturbation δU , grows in time, leading to the linear instability of the system.We define the gain function |Γ(k 3 )| as the absolute value of the growth rate, and we say that when the gain function is positive for some values of k, then the system is linearly unstable for those wave numbers.
The gain function can be computed by constructing an implicit relation, linking the eigen-frequencies ω to the eigen-wave-number k.This can be obtained as follows: we start by constructing the ideal generated by P W (k 2 ; λ), P Z (ω 2 ; λ), C(ω; σ), S Z (σ; λ) where P W (k 2 ; λ) ≡ P W (ξ; λ), P Z (ω 2 ; λ) ≡ P Z (ζ; λ) are the polynomials of the square of the differences of the eigenvalues of W and Z, respectively, C Z (ω; σ) = p 3 ω σ+(1 + p 3 ) ω +p 2 p 3 −k is the Cayley-Hamilton relation obtained from (38), and the polynomial whose roots are the sums of the pairs of eigenvalues z j of Z, so that the three σ-roots are We then apply the Buchberger's algorithm to compute the Gröbner basis for the ideal (see, e.g., [49]), so to obtain a polynomial in the new basis featuring only k and ω 2 .This polynomial reads as follows: For a fixed eigen-wave-number k 3 , the three ω-roots of H(ω, k 3 ) = 0 are the three possible values of the eigen-frequency ω 3 .The multi-valuedness of ω 3 (k 3 ) results from the polynomial of the square of the differences P W (ξ; λ) being of degree six, and thus entailing six values of λ for a given k 3 .By writing the equations for the real and imaginary part of ω in H(ω, k) = 0, and of λ in P(k 2 , λ), one obtains a system of four real equations for four real unknowns.By analysing this latter system, it is then possible to show that, if p 1 , p 2 and p 3 are such that E j = 0, j = 1, 2, 3, then the imaginary part of ω j is zero, for j = 1, 2, 3 and for all k > 0, only if the imaginary part of λ is zero for all k (see [48] for details).This leads to the following Proposition: Proposition 7.For a generic choice of p 1 , p 2 and p 3 (i.e.such that E j = 0), the plane wave solution of the 3WRI system ( 1) is linearly unstable.

Conclusions
We provide a classification of the stability spectra for the 3WRI system in the parameters space, and for each spectrum, we give the associated gain function.We find that for a generic choice of the parameters p 1 , p 2 and p 3 , the plane wave solution of the 3WRI system is linearly unstable (see Proposition 7).The parameters p j are related to the physical parameters (amplitudes a j , velocities c j and signs s j , j = 1, 2, 3) via the formulae (34) (Figs.4,5,6).
In analogy to the CNLS system where baseband MI is believed to be a sufficient and necessary condition for the onset of rogue wave-type localised peaks [27,28,50], we speculate that linear instability of baseband type can be a necessary condition for the same phenomenology in the 3WRI system, that is, the linear stability analysis of the system provides information about its subsequent nonlinear evolution.For instance, the plane wave solution used as a seed in the Darboux transformation generating the bright-dark-bright rogue wave triplet in the middle column of Figure 15 in [51], corresponds to a stability spectrum of 2B 0G 0L type (see Fig. 7), p 1 = 0.0, p 2 = −8/3 and p 3 = 3/5, i.e., c 1 = 4, c 2 = 1, a 1 = 1/2 √ 3, a 2 = 2/ √ 3 and s 1 = s 3 = −s 2 .Along this line, we have carried out preliminary numerical tests by numerically integrating, via the method of lines, the initial value problem from perturbed plane waves with different kinds of perturbations at time t = 0. We observe that, when the plane waves solution is associated with a spectrum that does not feature branches or twisted loops, namely when the gain function has only passband components, a random perturbation develops as a wave train.On the contrary, if the corresponding spectrum features at least one branch and no twisted loops, namely if the gain function has a baseband component, then a random perturbation develops as a series of space-time localised coherent structures.Finally, if the corresponding spectrum features a twisted loop, we observe that the gain function has neither a baseband nor a passband component, being nonzero in the neighbourhood of the zero wave number.In this case, we consistently observe blow up at finite time in our numerical explorations.All this phenomenology, in particular the connection between twisted loops in the spectrum, blow-ups and explosive instability will be matter of future investigations.

Fig. 1 .
Fig. 1.Regions in the (p 1 , p 2 )-plane (p 1 on the horizontal axis, p 2 on the vertical axis) associated with different numbers of gaps (G), branches (B) and split gaps (SG), for p 3 = −0.6 (left) and p 3 = 2 (right).The gap-branch topology of the small, unlabelled region bounded by E 2 , E 3 and E 4 , in the plot on the right, is 2G 0B

Fig. 3 .
Fig.3.Regions of different topology in the (p 1 , p 2 )-plane (p 1 on the horizontal axis, p 2 on the vertical axis) for p 3 = −0.6,when −7 < p 1 < −5 and −7 < p 2 < −5.This plot is a magnification of the 0 G 2B 2 L region displayed in the previous figure, which would otherwise be too small to be properly visualised and labelled

Fig. 4 .Fig. 5 .
Fig. 4. Stability spectrum (left) with 1 G 1SG 0B 1 L 1TL, obtained for p 1 = −90, p 2 = 60, p 3 = −0.6 and its associated gain function (right).As it is always observed numerically in the case of stability spectra with twisted-loops and splitgaps, the gain function features a nonzero component in the neighbourhood of the zero wave-number.The passband-type instability is associated with the presence of a loop in the stability spectrum

Fig. 6 .Fig. 7 .
Fig. 6.Stability spectrum (left) with 0 G 1SG 1B 0 L 1TL obtained for p 1 = −4, p 2 = 2, p 3 = −0.6 and its associated gain function (right).As it is always observed numerically in the case of stability spectra with twisted-loops and split-gaps, the gain function features a nonzero component in the neighbourhood of the zero wave-number.The baseband-type instability is associated with the presence of a branch in the stability spectrum

Proposition 6 .
(1) different typologies of stability spectra of the 3WRI system(1)are given in the Tablebelow

.
Finally, besides the symmetries of the curves on the (p 1 , p 2 )-plane, we observe that the polynomial Q 2 (ξ) is invariant under the simultaneous transformations p 1 → −p 1 and p 3 → −p 3 .This invariance, together with the symmetry of the curve E 4 , which is reflected with respect to the p 2 -axis when p 3 goes −p 3 , entails that the regions of different topologies listed in Proposition 6 are simply reflected about the p 2 -axis on ZAMP Fig.2.The figure on the left displays regions of different topology in the (p 1 , p 2 )-plane (p 1 on the horizontal axis, p 2 on the vertical axis) for p 3 = −0.6,when−4< p 1 < 4 and −4 < p 2 < 4. The topology of the unlabelled region bounded by E 3 and the convex component of E 5 bounding the 1G 1B 0L region is also 1G 1B 0L.If a region without a label is reached by crossing a component of the curve E 4 starting from a neighbouring region with a label, then the topology of the region without label coincides with the topology of the region with a label that one started from, but with the number of gaps increased by one and one twisted-loops replaced by a loop.The figure on the right displays regions of different topology in the (p 1 , p 2 )-plane (p 1 on the horizontal axis, p 2 on the vertical axis) for p 3 = −0.6,when−100< p 1 < 100 and −100 < p 2 < 100.The topology of the very narrow region on the left, bounded by E 3 and E 5 , is 0G 1B 1TL.In both figures, the labels do not include the counting of the split gaps, as this always coincides with the number of twisted loops the (p 1 , p 2 )-plane if we exchange p 3 with −p 3 .For this reason, in view of the classification of the loops, in Figs.2 and 3below, we choose one value of p 3 , namely p 3 = −0.6.