Complex Solitary Wave Dynamics, Pattern Formation, and Chaos in the GainLoss Nonlinear Schrödinger Equation
Abstract
A numerical exploration of a gainloss nonlinear Schrödinger equation was carried out utilizing over 180000 core hours to conduct more than 10000 unique simulations in an effort to characterize the model’s six dimensional parameter space. The study treated the problem in full generality, spanning a minimum of eight orders of magnitude for each of three linear and nonlinear gain terms and five orders of magnitude for higher order nonlinearities. The gainloss nonlinear Schrödinger equation is of interest as a model for spin wave envelopes in magnetic thin film active feedback rings and analogous driven damped nonlinear physical systems. Bright soliton trains were spontaneously driven out of equilibrium and behaviors stable for tens of thousands of round trip times were numerically identified. Nine distinct complex dynamical behaviors with lifetimes on the order of ms were isolated as part of six identified solution classes. Numerically located dynamical behaviors include: (1) Low dimensional chaotic modulations of bright soliton trains; (2) spatially symmetric/asymmetric interactions of solitary wave peaks; (3) dynamical pattern formation and recurrence; (4) steady state solutions and (5) intermittency. Simulations exhibiting chaotically modulating bright soliton trains were found to qualitatively match previous experimental observations. Ten new dynamical behaviors, eight demonstrating long lifetimes, are predicted to be observable in future experiments.
1 Introduction
Related forms of the nonlinear Schrödinger equation are used to explore nonlinear phenomena in many distinct physical systems: GinzburgLandau equations describe the envelope evolution of modelocked lasers, and superconductivity [1]; the cubic nonlinear Schrödinger equation treats deep water waves [2] and the dynamics of spinwave envelopes in magnetic thin films [3, 4]; a driven damped nonlinear Schrödinger equation models excitonpolariton and magnon BoseEinstein condensates (BECs) [5]; and the GrossPitaevskii equation models the mean field of atomic and molecular BECs [6, 7]. With the increasing supply of cheap computing power these systems have become the subject of extensive and sometimes rigorous numerical study.
Magnetic thin films have demonstrated potential as a versatile toy system for experiments on fundamental nonlinear dynamics [8]. Over the past two decades yttriumirongarnet (, YIG) magnetic thin films have been fruitfully studied by numerous experimental groups and have demonstrated a rich variety of nonlinear phenomena. These include bright and dark envelope solitons [9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 15, 20, 21, 22], soliton trains [23, 24, 25], möbius solitons [26], FermiPastaUlam and spatial recurrence [27, 28], soliton fractals [29], random solitons [30], chaotic spin waves [31, 32, 33], multiple solitons [34], and chaotic solitons [35, 36].
A majority of these phenomena were observed on active feedback rings; such feedback structures are ubiquitous within science and physics in general. Rings in particular are commonly used to study wave dynamics when one seeks a quantized wavenumber, periodic pumping, selfgeneration, or other resonant phenomena. Active rings, so called for the presence of periodic linear amplification, allow for the direct compensation of the major loss mechanisms present within a system. This permits one to drive the system into quasiconservative regimes, enabling the observation of dynamics on scales of several to tens of thousands of round trip times. Dynamics with lifetimes of this order would otherwise be prohibited by the presence of dissipation.
Yet, within the context of spin waves in magnetic thin films, little work has been carried out to develop an adequate theory for describing the rich range of behaviors evident within these recent experimental works. The integrable cubic nonlinear Schrödinger equation (NLS), while successful in quantitatively describing both dark and bright soliton trains, is unable to reproduce more complex phenomena such as the chaotic oscillation of soliton envelopes. However, there has been significant effort within the modelocked laser community to study analogous driven and damped systems. Works on dissipation terms and saturation [37, 38, 39, 40], the study of dissipative solitons dynamics [41, 42, 43, 44, 45] and other numerical studies of the cubic quintic complex GinzburgLandau equation [46, 47] are highly relevant to the development of a driven damped model for spin waves in magnetic thin films. These works explore the dynamics of solitary waves and their associated wave equations under the influences of gain and loss. For example, the dynamics of near steadystate dissipative solitons have been considered in detail; such studies include rigorous mappings of stable and unstable regions of parameter space [43, 44, 45]. Similarly, initial transient behaviors have been the subject of significant research efforts by the modelocked laser community. Transients are of interest for potential applications in signal processing and communication. To date there have been no efforts towards the characterization of long lifetime dynamics of soliton trains driven from equilibrium within active feedback rings. The work presented here demonstrates such an effort for a generalized nonlinear Schrödinger equation with a focus on applications to nonlinear spinwaves propagating in an active feedback ring.
Our paper is outlined as follows. Section 2 introduces the model to be studied along with the associated operating limits; here the methodology and scope of the simulations are explicitly defined. Experimental contexts for the work are also considered in section 2. Results in the form of eleven unique complex dynamical behaviors are presented and categorized in sections 37. Section 3 contains simulations of chaotically modulating soliton trains. Spatially symmetric and asymmetric solitary wave interactions are presented in section 4. Four examples of dynamical pattern formation are given in section 5. Two cases of steady state evolution are reported in section 6. Intermittent solutions are discussed in section 7. Finally, a discussion of numerical convergence and solution robustness is given in section 8. The work is summarized in section 9. Animations of each dynamical behavior discussed in sections 37 are available at http://mines.edu/~jusander/GLNLS_dynamics/.
2 Model Overview
Motivated by the works discussed in section 1, we propose a generalized governing equation for spin waves in magnetic thin film active feedback rings: the gain loss nonlinear Schrödinger equation (GLNLS),
(1) 
where is a dimensionless complex magnetization amplitude defined as ; here is the dynamic magnetization while is the saturation magnetization; is the dispersion coefficient; and are the cubic and quintic nonlinearity coefficients, respectively; is the ‘temporal’ evolution coordinate; is the ‘spatial’ coordinate of propagation boosted to the group velocity of the envelope; and , , and are the linear, cubic, and quintic gains (losses) if positive (negative). All parameters are taken to be real as the complex nature of the coefficients is explicitly accounted for in (1). The local intensity of the magnetization amplitude is given by . The norm and energy at a given time, , are defined as
(2) 
and
(3) 
respectively, where the integrals are taken over the length or circumference, , of the feedback ring. All norms, intensities and energies given within figures and animations are scaled by , and respectively where corresponds to the initial condition used during numerical simulation. Numerical values given within the text are not scaled. The specific choice of initial condition is discussed later in this section.
The gain and loss present within the GLNLS may be viewed as an expansion of saturable loss expressions studied separately by Ablowitz and Akhmediev [38, 39, 43, 44]. The higher order nonlinearity, , may be used either as a saturation of cubic nonlinearity or an additional selfsteepening; both cases are studied in the literature [2]. The GLNLS omits other terms commonly included in cubic quintic complex GinzburgLandau equations such as spectral filtering, periodic pumping, and integral mean terms, as they are not needed in this physical context. We are likewise compelled by Occam’s Razor to choose the simplest possible model which nevertheless reproduces measurements in magnetic thin film active feedback rings. NLSlike equations may be derived in magnetic thin films by use of a slowly varying envelope approximation [48], more rigorously through conservation considerations and a Hamiltonian formalism [9, 49], or directly from Maxwell’s equations using multiscale methods [50].
The operating limits of the GLNLS are motivated principally by experimental work on the excitation of chaotic solitons in YIG stripbased active feedback rings [35]. A block diagram of the active feedback ring experiment is shown in figure 1. The ring is comprised of a nonlinear propagation medium, in this case a magnetically saturated crystalline YIG thin film, connected via two transducers to an electronic feedback loop. The electronics loop is constructed of a directional coupler, allowing real time observation at an oscilloscope and/or spectrum analyzer, and an amplifier/attenuator pair for real time adjustment of ring gain. The GLNLS demonstrated qualitative agreement with the low dimensional chaotic modulation of a bright soliton train, as will be discussed in section 3. Detailed experimental results are discussed in Wang et al. [35]. These experiments indicate that nonlinearity and dispersion are the dominant sources of envelope shaping for chaotic spin wave solitons and that the losses present in the ring are fully compensated for by the amplifier. This imposed two constraints on modeling. (1) The coefficients and must be orders of magnitude larger than , , and . (2) The linear amplifier must compensate both the linear and nonlinear losses present in the film, requiring a net averaged (over many round trips) linear gain, . Likewise, the dissipative terms represent the net gain and loss processes occurring in the ring averaged over several round trip times. One expects the use of this approximation, and therefore the model, to be valid when the time scale of envelope modulation is much greater than the soliton round trip time. All simulations were performed using adaptive time step CashKarp RungeKutta for temporal evolution and pseudospectral techniques for spatial propagation [51, 52]. Periodic boundary conditions modeled propagation around a ring.
Parameter space explorations were explicitly chosen to encompass the GLNLS operating regime for magnetic thin film systems, while extending into other limits that could be of interest to other systems where the GLNLS is a useful model. Along with the previously mentioned restriction on the sign of only cases with cubic losses, , were considered. Both instances of saturating

,

,

,

,
for a total of eight possible values of and , five choices for , and 17 unique choices for . Ignoring cases with solely gains present we performed 5,470 unique simulations. An additional 1,530 simulations were undertaken with random parameters. The value for any single parameter in these simulations was generated by multiplying a pseudo random number between zero and one, from the uniform distribution, by an order of magnitude and sign chosen at random, again with uniform weight, from a parameter’s allowed values, as defined above. To avoid ambiguity all statements in this paper concerning the relative size of GLNLS parameters refer to the order of magnitude and not the sign.
Simulations began as a bright soliton initial condition obtained via imaginary time relaxation [53], the ground state solution to the GLNLS with , , , and set as zero, with . Experimentally this corresponds to a stable bright soliton circling within a YIG stripbased active feedback ring, a solution analogous to a soliton train. Physical GLNLS parameter values are obtained by fitting this initial condition to experimentally observed bright soliton train conditions. This choice of units also fixes the ratio of used in simulations, while the amplitude of and dictate the simulation timescale. We assumed that the dimensionless spinwave intensity is directly proportional to the spinwave power, , since experimental measurements of voltage are taken across a diode with quadratic behavior and are generally taken to be proportional to power. Values typical for a chaotic soliton experiment are , the round trip time; , the transducer separation; , the electronic loop propagation time; , the group velocity; , the cubic nonlinearity; and , the dispersion. Using these parameters one finds where is the scaled temporal unit used in simulations. This relation may be used to immediately transform code values for , , and , which share units of inverse time, to physical values. For example the largest studied linear gain is which matches the order of experimentally approximated linear losses for magnetostatic backward volume spin waves in YIG thin films [54]. Experimentally a time series is recorded at the detection transducer with the full waveform being captured once a round trip after the signal has propagated a length between the transducers and passed through the electronics loop. The length of the ring, , is taken to be the transducer separation, , as the propagation delay is orders of magnitude smaller than the round trip time, . Simulations explicitly model the entire feedback loop at the group velocity of the waveform. A time series may be reconstructed from numerical data by concatenating the simulated waveform after a temporal evolution of or a spatial evolution of . In this work we adopt the former convention to ease the direct comparison of simulations to the power vs. time data often observed experimentally for spin waves in magnetic thin films. Such a reconstructed time series is labeled throughout the paper. A time series of solitary wave peak intensity at successive round trips is useful in studying modulating single solitary wave trains and is defined by
(4) 
where is the round trip time and is the total number of round trips.
Over 180000 core hours were utilized to conduct more than 10000 unique simulations and convergence studies. An initial study of 3500 simulations was undertaken to explore the extent of transient effects and the numerical convergence behaviors of the GLNLS. A summary and analysis of the subsequent 7000 simulations, corresponding to over 3 TB of data, are presented in sections 37. Approximately 1500 simulations were evaluated in detail; the remaining simulations were spot checked for consistency. Dozens of complex dynamical behaviors were identified during the course of simulation. We call this system complex because it displays a rich variety of dynamical behaviors, including chaos, robust emergent solitarywave features, and generally multiple scales in both space and time. Solution types were divided into three stability cases, with each case corresponding to roughly of observed dynamics. The three cases are stable, intermittent and unstable. Temporally stable solutions demonstrated substantial observable lifetimes, greater than or 7000+ round trips, and robustness to variations in initial conditions of at least . Evolution was found to be least sensitive to changes in and and most sensitive to perturbations in . In general the effect of changes in initial conditions tended to degrade the lifetime of dynamical behaviors and push solutions towards the intermittent case. Nine temporally stable distinct dynamics and two separate cases of intermittency are discussed below.
3 Chaotic Modulation
The chaotic modulation of stable solitary wave trains was observed for solutions containing strongly saturated cubic nonlinearity, , and the lowest studied ring gains, , with matching orders of cubic and/or quintic losses. A single bright soliton is observed to circulate within an active feedback ring while exhibiting complex modulations in peak intensity. Low ring losses are anticipated for this solution type, as experimentally observed chaotically modulating soliton trains have lifetimes measured in seconds. The presence of a single stable bright soliton suggests that nonlinearity and dispersion are the dominate forces in peak shaping. These are two conditions used during the derivation of the GLNLS, (1), as discussed previously in section 2.
The chaotic nature of measured time series was verified by using standard phase space reconstruction techniques available in the open source Nonlinear Time Series (TISEAN) package to arrive at a stable correlation dimension, [55]. The correlation dimension, a phasespace invariant, was estimated via computation of the correlation sum for increasing embedding dimensions of the time series [56]. The standard embedding procedure of Taken and Sauer was followed using timedelayed reconstruction of the time series [57, 58]. The time delay was chosen as the first minimum of autocorrelation to maximize the linear independence of the time delayed vectors. As the phase space was reconstructed from a single time series, a Theiler window of ten times the single round trip time was used to avoid the misinterpretation of temporal correlation as geometrical structure on the attractor [59]. If the correlation dimension was observed to saturate with increasing embedding dimension the time series was said to have a stable correlation dimension. If the stable correlation dimension was not an integer then the system was said to be chaotic.
We further required the correlation dimension to be stable across a wide range of embedding parameters as one expects the reconstructed attractor to be invariant under smooth transformations. This requirement was extremely conservative as it was computationally onerous and sensitive to noise. However, such a requirement forbids the optimization of phase space invariants by the tuning of embedding parameters, and the requirement of saturation across embedding dimension eliminates any assumptions required to study a single reconstruction. Additional indicators of chaos include broadband spectra and positive Lyapunov exponents [56]; note both these properties are shared with noise so a finite correlation dimension is necessary to demonstrate chaotic, rather than random, motion. The principle challenge to finding a stable correlation dimension was isolating a stationary solution.
Two examples, with peak variations of and about their mean, are shown in figure 2 and figure 3, respectively. Percent peak variation is defined as
(5) 
where was previously defined in Equation (4) and is the sample variance. These values are chosen to match the peak variation of two low ring gain chaotic solitary wave trains observed experimentally by Wang et al. [35]. In both figures panel (a) shows the intensity of the single soliton peak for 2500 consecutive round trips while panel (b) shows 100 round trips as would be observed experimentally, as in Figure 4(d) below, and each vertical line is in fact a bright soliton of finite width. The single soliton peak intensity is immensely complex on inspection and is at worst random and at best chaotic or quasiperiodic. Figure 2(c) and 3(c) show the phase space reconstruction for an embedding dimension of 2 and a time delay of 1, also known as a return plot, for 100000 round trips of the full time series. The finite width and structure of the reconstructed attractor is one indication of chaotic, as opposed to random, motion. In figure 2(d) and 3(d) is shown correlation dimension versus embedding dimension for each variation case. Both cases saturate above an embedding dimension of 15 to a correlation dimension of and , respectively. Error estimates are confidence intervals given by two times the standard deviation for values of for embedding dimensions above saturation. This low dimensional chaos closely matches the low ring gain experimental observations by Wang et alwhere variation yields a correlation dimension of . However the numerically generated peak variation does not reproduce the high dimension chaos, , observed experimentally at matching variations [35]. The cause of collapse at embedding dimensions 6, 16 and 26 for the modulation case has not been rigorously determined but is robust against reasonable perturbations in embedding parameters. The periodicity of the effect suggests the cause is related to sensitivities in the correlation sum to temporal correlations and finite time series. The embedding procedure is also sensitive to time series periodicity, which is present in these low dimensional examples [57, 58]. Low dimensional chaos often presents as widened Fourier peaks rather than pure broadband spectra. The oscillation of for low embedding dimension is a common phenomenon as the embedding procedure is not an accurate reconstruction of phase space unless the embedding dimension is at least twice the box counting dimension of the system’s attractor [56].
We find numerically that amplitude of peak modulation and the dimensionality of the chaos are principally dependent on the magnitude of the saturating quintic nonlinearity, . The presence of both a linear gain and nonlinear loss term is necessary for a stable correlation dimension to be determined. Chaotic modulations of the train envelope are the most complex examples of a more general modulation behavior. Parameter space explorations yielded examples of bright soliton trains with no, periodic, multiperiodic or quasiperiodic modulations. We note these types of deepening modulations were experimentally observed as the first generations of soliton fractals [29].
4 Symmetric and Asymmetric Interacting Solitary Waves
When more than one solitary wave propagate with differing group velocities, enabling dynamics such as collisions, we say these waves interact. Two distinct cases where the spatial features of solitary wave interactions are symmetric or asymmetric under rotation are discussed below.
4.1 Symmetric interaction
Symmetric interaction solutions are highly complex, but ordered, gain driven interactions between a number of intensity peaks varying from two to more than twenty. These solutions evolve in intricate and complicated patterns but maintain symmetry in space under a rotation of rads. The solution intensity exhibits a constrained modulation about a stable mean, but is energetically unstable. The energy of the system grows approximately linearly in time and is closely correlated, with a correlation coefficient of , to the timeaveraged number of peaks present in the system. The sample correlation coefficient is a measure of the linear correlation between two variables and is defined as
(6) 
where is the standard deviation of ; and are the number of peaks and system energy at the round trip. This relationship suggests every intensity peak present in the system has similar energy. Peaks undergoing symmetric interactions also demonstrate persistence in time under collisions and have linear or constant phases, both characteristics of bright solitons. Further, individual intensity peaks may also be fit to a profile when they are spatially isolated from other peaks circulating the ring. A typical example is illustrated in figure 4(a) by a spatiotemporal plot of intensity across 800 round trips, each vertical slice shows the waveform on the ring at a specific round trip. There exists a stark symmetry in dynamics with respect to a rotation by rads. Figure 4(b) and (c) show the scaled norm and energy, respectively, for the same time frame. Over these 800 round trips we note the norm varies about a stable mean by while the system energy increases by . A reconstructed time series of data presented in panel (a) is shown in figure 4(d) to indicate what the behavior would look like if measured experimentally at a single observation point and discretely in time. We note that the symmetry demonstrated by the spatiotemporal intensity plot, figure 4(a), is not evident in what appears to be a highly noisy time series. Whether the symmetry observed numerically persists when the iterative nature of amplification and transmission delays in an electronic feedback loop are considered remains an open question.
Symmetric interactions are observed to evolve in systems with linear gains between and and cubic losses of the same, or , orders of magnitude. This long time stable evolution requires a nearbalanced system where linear gain is the dominant force and peak growth is meaningfully restricted by the presence of nonlinear losses. No solutions had initial peak growth above prior to the initial splitting event.
The dynamics within this regime demonstrate a characteristic splitting process, diagrammed in figure 5(a)(d). The initial bright soliton modulates and grows until the domination of linear gain over nonlinear loss in lowintensity regions yields a nonzero intensity floor. Energy enters the system until these low lying excitations reach intensities where attractive nonlinearity and dispersion may shape the excitation into a stable solitary wave close in form to the wellknown hyperbolic secant. The new peak then begins to interact with its neighbors. This same procedure results in the generation of a second, then third, and so on, intensity peak. Thus, in contrast to more typical nonlinear partial differential equations which give rise to fixed soliton dynamics for all times, the GLNLS here displays a particular soliton dynamics on long but not infinitely long timescales. This gives rise to the possibility of a new form of integrability which is relevant on long but not infinite times, and may require the development of new mathematical formalisms, in particular a multiscale approach in time. The timing of the initial splitting event varies from to where is defined as the moment gain and losses are turned on. The effect of quintic loss/gain is superficial to this solution category until orders above when it begins to dominate the dynamics. Quintic losses (gains) result in slower (faster) rates of initial splitting, but do not have any meaningful impact on the rate of energy gain. This splitting process is stabilized (weakened) by the addition of an attractive (repulsive) quintic nonlinearity term of the same order as the cubic present in the system. Higher orders of quintic nonlinearity destroy the stability, driving the dynamics into the intermittent regimes described later in section 7. This solution type demonstrates a high sensitivity to initial conditions, which is discussed in section 8. A single round trip of a symmetric interaction solution closely resembles the multipeaked solitons previously reported by Wu et al. [34].
4.2 Asymmetric Interaction
Asymmetric interaction solutions are loss driven solutions which behave similarly to the symmetrical case discussed in section 4.1 but do not maintain a spatial symmetry with respect to rotations around the ring. The number of interacting peaks was observed to vary from five to twenty depending on the parameters of the simulation. The total number of peaks is conserved, in an average sense, after spatial symmetry about the feedback ring center breaks and is closely correlated, , to the system’s energy. Here is the sample correlation coefficient defined by (6). An example is shown in figure 6(a) by a spatiotemporal plot of intensity over 10000 round trips. Symmetry about the feedback ring center, , can be seen breaking near round trip 4300. An animation of this symmetry breaking is available online. A scaled intensity plot of a single symmetric (asymmetric) round trip is shown in figure 6(b) (figure 6(c)). The interacting peaks are seen to be nodeless, and evolve about a nonzero, , intensity floor. The stability of the asymmetric interaction solution type is demonstrated in figure 6(d) and (e) showing scaled norm and energy, respectively, over the same 10000 round trips. Normalization varies about a stable mean by while energy modulates by ; this stands in contrast to the energetic instability inherent to symmetric interactions, section 4.1.
Asymmetric interactions are observed in systems with linear gains, , between and 1 and cubic loss, , of the same order of magnitude. Quintic gains and losses, , are stable up to this same order of magnitude. The number of peaks and peak height increased (decreased) with the presence of attractive (repulsive) quintic nonlinearity of the same order as the cubic. Higher orders of quintic nonlinearity push the solution into intermittency, a temporally unstable class of solutions discussed in section 7. The solution intensity floor varies with parameter choice, including nonlinearity, but trends towards the constant intensity which satisfies the energy balance of the GLNLS. The balance is given explicitly by the expression
(7) 
where , , and are the same terms as in (1) discussed in section 2 and we choose the smallest positive solution. For the simulation shown in 6 we have , , and corresponding to an average solution intensity of which closely matches the numerically observed value of . Error bounds are given by two times the standard deviation of intensity across all available round trip data.
This regime demonstrates a characteristic splitting process, diagrammed in figure 5(e)(h). An initial bright soliton initial condition grows and flattens into a plateau under the influence of a strong linear gain and saturating nonlinear losses. Once the nonzero plateau expands to fill the feedback ring the central peak undergoes a splitting procedure similar to that observed for symmetric interactions, diagrammed in figure 5(a)(b). The domination of linear gain over nonlinear losses in low amplitude regions produces small peaks. These smaller excitations grow until the system’s attractive nonlinearity and dispersion shape them into solitary wave intensity peaks. Unlike the process for symmetric interactions, section 4.1, this splitting process occurs within the first of evolution, where is defined as the moment gains and losses are turned on, and saturates within the first yielding an energetically stable excitation. The amplitude of intensity peaks relative to the plateau intensity varies from to , but the peak heights measured from the plateau mean are of the same order as those observed in symmetric interactions.
This solution type demonstrates a high sensitivity to initial conditions, which is discussed in 8. This sensitivity and the highly complex nature of the evolution are hallmarks of chaotic dynamics. However, attempts to arrive at a converged correlation dimension using the methods discussed in section 3 were inconclusive. Such sensitivity is typically characterized by a positive Lyapunov exponent [56]. While a careful determination of the largest Lyapunov exponent requires a rigorous reconstruction of phase space we may estimate the exponent numerically by evolving nearby trajectories in time. Direct measurement suggests a Lyapunov exponent between and . This rate of trajectory separation is of the same order as that observed experimentally ( [35]) for the modulating soliton train discussed previously in section 3.
5 Dynamical Pattern Formation
Four distinct robust dynamical patterns which demonstrate lifetimes of at least or 7000 round trips were located during GLNLS parameter space exploration. Solutions of this group differ from previously discussed solution behaviors in that they exhibit a periodic recurrence of their characteristic dynamic. Self organization of this kind is common in open nonlinear systems [60]. These examples are discussed to demonstrate the breadth of pattern formation supported by the GLNLS under fixed choice of and . The regions of parameter space supporting dynamical pattern formation violates the assumptions underlying the derivation of the GLNLS in the context of magnetic spin waves, as discussed in section 2, owing to the high order of quintic nonlinearity and losses which drive evolution. However, the GLNLS is a useful model in a variety of systems including laser cavities, as discussed in section 1, and these dynamics may appear in such contexts.
5.1 Central Peak Recombination
Central peak recombinations exhibit a complex 5 peak solitary wave recombination pattern with a periodicity of 180250 round trips, depending on parameter choices. This behavior is driven by a strongly attractive quintic loss, and a linear gain of with quintic loss of . The presence of quintic nonlinearity has a severely negative impact on the behavior lifetime. The median wave height of central peak recombination solutions satisfies the energy balance equation, (7). For the example shown in figure 7 we predict an average intensity of , corresponding to the parameters have , and , which closely matches the observed numerical average intensity, . The error estimate is defined as in section 4.2.
An example of central peak recombination is shown in figure 7. Panel (a) shows a spatiotemporal plot of scaled intensity over 3500 round trips and (b) shows the scaled energy over the same round trips. The periodicity of the recombination is evident in the spatiotemporal plot, which contains 16 periods. The dynamics are most readily described starting when the central peak collapses. Immediately following the collapse, the peaks on either side of the ring center propagate towards the middle of the ring and recombine into a new central peak matching the original peak’s amplitude. At the same time the outlying peaks split into two. The innermost of these new peaks grows until one finds three central peaks of equal amplitude. At this point the central peak undergoes collapse and the process repeats. A single bright solitary wave propagates unperturbed along the edge of the ring. This process is animated in an attached movie, available online.
5.2 Complex Copropagation
The complex copropagation solution was so named as it resembles the steady state copropagation solution (see section 6.2 below) and is likewise energetically stable. It differs primarily in that the waveform undergoes complex, but periodic, modulation. The dynamics also occur on a nonzero density floor satisfying the GLNLS energy balance equation, (7). The example shown in figure 7(c)(d) was simulated with the parameters , , , resulting in an anticipated average intensity of . This prediction closely matches the numerically observed intensity , where the error is defined as in section 4.2. Like central peak recombination the complex copropagation behavior is driven by a large quintic loss. The dynamical patterns demonstrated by these solutions also require a large attractive quintic nonlinearity. The parameter space region which supports these behaviors is characterized primarily by large, negative quintic terms: . The smallest linear gain which compensates these driving nonlinear losses is . These solutions are in general insensitive to the choice of cubic loss, with any value smaller than supporting the observed dynamical pattern.
Figure 7(c)(d) illustrates this behavior. Panel (c) shows a spatiotemporal plot of scaled intensity over 12000 round trips and (d) shows the scaled energy over the same round trips. In panel (d) the energy has been offset for clarity, . The behavior is characterized by the spatiotemporal plot which shows two spatially stable bright solitary waves occupying the center and edges of the ring. The central solitary wave is flanked on each side by a set of two periodically oscillating solitary waves for a total of six large peaks being equispaced around the ring. Six additional small amplitude peaks occupy the space between each larger wave. The entire waveform breathes between two distinct energy states with a period of 750 round trip times. The frequency of oscillation matches that predicted by a simple twolevel quantum system where . For the GLNLS we have amd , as defined in section 2. Taking the average energy difference between states, see figure 7, one predicts an angular frequency of compared to the observed oscillation frequency of . An animation of the breathing is available online.
5.3 Spatial Shifting
Spatial shifting solutions are simulations which exhibit energetically stable evolution with a welldefined and periodic shifting of the spatial location of the dynamical behaviors. In all observed cases the underlying energetically stable dynamics are evenly distributed bright solitary waves copropagating on an intensity floor which satisfies the GLNLS energy balance given by equation (7). For the example shown in figure 8(a)(b) we have , and corresponding to which closely matches the numerically observed average intensity of . The solitary waves spontaneously split at a constant periodicity and reform into an identical set of copropagating peaks with a spatial shift defined as where is the feedback ring length and is the number of peaks present in the simulation. All peak properties as well as the splitting dynamics remain consistent through multiple periods. A strong attractive quintic nonlinearity is required to support this dynamical behavior, as seen previously with central peak recombinations and complex copropagation in sections 5.2 and 5.1. Spatial shifting is seen in simulations with quintic nonlinearities of and moderate linear gains of . Cubic losses near support this behavior, while quintic losses were found to be unimportant until above values of where they dominated the dynamics.
An example of temporal shifting is illustrated in figure 8(a)(b). Panel (a) shows two bright solitary waves copropagating while undergoing a spatial shift of every 22000 round trips. The shifting event occurs over 1500 round trips. Panel (b) shows the solution’s scaled energy over these same round trips; the energetic stability of the copropagation regimes is demonstrated. The energy profile of each shifting event was found to be identical.
5.4 Breathers
Solitary wave breathers on a ring are characterized by a single solitary wave which undergoes a periodic disappearance of the peak and reappearance at the other side of the ring. The frequency of breathing increases with system energy. The solution is not energetically stable and breathing frequency increases until the system reaches a new dynamical behavior. Numerically observed lifetimes were never less than 20000 round trips, or 3 milliseconds. The wave breathing is driven by a strong quintic loss, , with comparatively weak linear gain, , and cubic loss, , terms. The solution type is sensitive to the presence of quintic nonlinearity with any magnitude above , whether attractive or repulsive, pushing the dynamics into the intermittent regime, discussed later in section 7. Linear gain dominates during low intensity periods of the breathing behavior, resulting in a nonzero intensity floor which ultimately drives the collapse of stable breathing.
Figure 8(c)(d) contains a typical example of solitary wave breathing. The periodic spatial shifting of the bright solitary wave is seen as a relocation from the center of the ring to the other side in panel (c). An average breathing period of 1200 round trips is observed in this example. A positive linear trend in energy, see panel (d), is the result of linear gain causing growth in lowintensity regions. A periodic high rate of energy growth matches the low intensity period following the collapse of bright solitary waves. An animation of the breathing behavior is available online.
6 Steady State Solutions
Simulations which evolved into energetically stable static wave forms were named steady state solutions. Two distinct steady state solutions were isolated from the parameter space exploration: multipeaked solitary waves and copropagating solitons.
6.1 Multipeaked Solitary Waves
Multipeaked solitary waves were characterized by energetically stable, to machine precision, nodeless complex waveforms that evolve without exhibiting any time dependence in their intensity. The shape of the wave and the number of principle peaks varies from two to eight in studied cases, depending on parameter choice. Symmetric and asymmetric waveforms were observed. Multipeaked solitary waves were observed for any linear gain, , below 1 and cubic losses, , of orders of magnitude. The impact of quintic losses and gains principally affected the median wave height according to the GLNLS energy balance equation, (7). For the multipeaked solitary wave shown in figure 9 we have , and corresponding to an estimated average intensity of which closely matches the observed value, . The error was previously defined in section 4.2. As with previous examples, high values of relative to lead to the term dominating dynamics and the solution leaving the steady state solution class. Positive, or saturating, values of quintic nonlinearity lead to reductions in secondary peak heights while attractive values leads to the presence of additional principal peaks via further shaping of secondary peaks. The overall shape of the multipeaked solitary wave, including the number of principal and secondary peaks, is dependent on the choice of parameters.
A typical example is shown in figure 9 of a symmetric multipeaked solitary wave with two principle and two secondary peaks. Figure 9(a) shows a spatiotemporal plot of scaled intensity over 12000 round trips with each vertical slice showing the intensity across a single round trip. Panel (b) is the scaled intensity plot of the final round trip and panel (c) shows the static solution energy over the same evolution period. Not all multipeaked solitary waves travel at the group velocity as the example in figure 9. This solution type is the most commonly observed long time behavior in studied simulations and was one of the behaviors present in a majority of the intermittent cases, discussed further in section 7.
6.2 Copropagating Solitary Waves
Copropagating solitary waves are the second steady state isolated during parameter space exploration. Copropagating solitary waves are time independent solutions where identical bright solitary waves propagate alongside one another without interacting. Similar soliton solutions of the simple cubic NLS are well studied and the number of solitons is found to be proportional to the power of the initial condition relative to the value of [2]. Periodic boundary conditions require solutions have an even number of nodes. Within studied solutions was observed to vary from two to eight. Figure 10 shows the same physical quantities as plotted in figure 9, with panel (c) again demonstrating the energetic stability of the solution type. The peak shape, shown in panel (b), is not consistent with either bright or dark solitons. The example plotted in panel (a) exhibits a modulation in peak heights with a variance of about the mean. While the variation is not visible in panel (a) it can be observed in an animation of the evolution, available online.
Stable copropagation was observed only in an isolated region of parameter solutions with and . Quintic gains, , of orders higher than the cubic present or quintic nonlinearities, , with magnitude higher than (the lowest order studied) drove the solution out of the steady state and in general pushed solutions into the intermittent class, discussed in section 7. Lower orders of quintic gain did not have any meaningful effect on stability, the number of peaks or peak height.
7 Intermittent Solutions
Intermittent solutions demonstrate numerous distinct dynamical behaviors as the waveform evolves in time. The lifetime of these behaviors ranges from hundreds of round trips to hundreds of thousands. This corresponds to up to before the dynamics transitions from one behavior to another. These solutions are robust to at least variation of initial conditions in the sense that they do not degrade to noise or experience blowup. Such perturbations do have significant effects on the relative lifetime of each dynamical behavior and even the types of behaviors a simulation exhibits. Quantitative matching of the intermittent dynamics to experiment will offer a challenge due to their highly transient nature; however, qualitative behaviors should be observable experimentally. In general, intermittent solutions spend a majority of their time in aperiodic evolution between distinct dynamical behaviors. Intermittent solutions can exhibit all of the behaviors previously described as temporally stable for a finite numbers of round trips. Intermittency is the typical dynamic exhibited when stable solutions are perturbed and is therefore not observed only in isolated regions of parameter space. Stable solutions are robust to variations in initial conditions, as previously stated in section 2. Intermittency is observed when perturbations exceed , however it bears mention that the necessary value is ultimately highly dependent on both solution type and the parameter being perturbed. Hundreds of intermittent simulations were identified during the study.
Two illustrative examples are shown in figure 11. The same physical quantities are shown as in figure 8. Panel (a) shows a typical simulation with three distinct multipeaked solitary wave regimes separated by two aperiodic regimes exhibiting splitting, modulation and copropagation behaviors. The energy is shown in plot (b) and was relatively constant during each of the multisolitary wave regimes. The aperiodic regimes demonstrate significantly lower energy than the finite lifetime multisolitary wave excitations. Panel (c) shows a simulation which exhibits periods of complex four solitary wave copropagation interspersed with periods of aperiodic dynamics. The lengths of successive periods of dynamical behavior are highly variable and sensitive to both parameter choice and initial condition. This sensitivity makes numerical convergence difficult to demonstrate, as discussed below in section 8.
8 Numerical Convergence and Quantitative vs. Qualitative Robustness
Simulations used well established algorithms, fifth order adaptive CashKarp RungeKutta in time and pseudospectral methods in space [51, 52]. Simulations were run with a spatial grid of 256 and a single step truncation error of . The maximum number of time steps performed in a single simulation was . Initial conditions were generated using imaginary time propagation and a single step truncation error of . Stability of initial conditions was confirmed via real time propagation. To machine precision all initial states exhibit zero change in energy when propagated in real time for steps. Initial conditions for different spatial resolutions have fixed differences in energy owning to discretization. This discretization error decreases exponentially for increasing spatial resolution: , and . Results for each solution class were compared across grid sizes of 512 and 1024 and a single step truncation error of to verify numerical convergence for algorithms used for both space and time propagation. We present convergence data for two distinct groupings of solution classes, those which exhibit extreme sensitivity to initial conditions and those which do not. The former demonstrate a qualitative robustness, while the latter are quantitatively robust.
Convergence can be demonstrated by relative difference. Given two sets, and , of data consisting of directly comparable observations then the relative difference at the entry is defined as
(8) 
This quantity offers a simple, unitless measure of the relative difference between two quantities.
8.1 Quantitative Robustness
Solution classes which did not demonstrate a marked sensitivity to initial conditions were numerically converged in a traditional manner. A distinct measurable, in this case energy, is quantitatively compared across successive time steps under different spatial and temporal resolutions. Convergence data is graphically displayed in figure 12 for the initial condition used during simulations, figure 12(a), as well as four categorical behaviors, panels (b)(e). In each case the solid blue line compares spatial resolutions of 256 and 512 grid points, while the dashed red line compares the spatial resolutions of 512 and 1024 grid points. Figure 12(a) shows the fixed discretization error discussed previously in section 8 while figure 12(b)(e) demonstrate the spatial convergence of each dynamical behavior over the entire evolution period is as good or better than that of the initial conditions. The greatest observed single time step relative spatial resolution error was . The greatest observed single time step relative temporal resolution difference was . The solution types listed here were quantitatively converged:

Complex copropagation

Spatial shifting

Breathers

Multipeaked solitary waves

Copropagating solitary waves
8.2 Qualitative Robustness
A subset of observed dynamical behaviors, from both the temporally stable and intermittent categories discussed in sections 37, demonstrate an extreme sensitivity to initial conditions. These solutions were robust to variations in initial conditions and parameters of at least in the sense that such perturbations did not yield a shift in their categorization. However, changes in initial energy or in loss parameters of the order and lower resulted in distinct dynamics within that categorical behavior and shifts in the starting and ending times. We note shifts in spatial resolution introduce variations of this order to the relaxed initial condition. Therefore solutions exhibiting this sensitivity may not be converged numerically in the traditional sense.
An illustrative example is given by the asymmetrical interaction behavior, discussed in section 4.2, after spatial symmetry about the feedback ring center has broken. Figure 13(a)(c) depicts scaled intensity across the same 2500 round trips for three different choices of spatial resolution: 256, 512 and 1024 spatial grid points respectively. The behaviors across the three spatial resolutions are qualitatively similar with each exhibiting an asymmetric solitary wave interaction which is characteristic of the solution type. However, the detailed dynamical behaviors of each case are quantitatively different. Further, as shown graphically in Figure 13(d), the three solutions have energies which vary by less than . The behaviors listed below all demonstrated sensitivity similar to that discussed here with a relative energy difference no greater than :

Symmetric interaction

Asymmetric interaction

Chaotic Modulation

Central Peak Recombination
This sensitivity is a typical property of evolution towards and around a strange attractor. Major attempts were made to quantify the dimensionality of the attractor as discussed in section 3 by the authors. No stable correlation dimensions were located for any of the complex dynamical behaviors which demonstrated sensitivity to initial conditions, excluding chaotic envelope modulation. A Lyapunov exponent was estimated for the asymmetric interaction case, see section 4.2. Further quantification and exploration of the phase space properties of GLNLS solution types and the GLNLS itself are warranted.
8.3 Unstable
Any simulation which evolved into a trivial result, degraded to noise, blew up or decayed to zero is considered to be unstable. Approximately of studied simulations are unstable. This is not surprising as the explored parameter space includes cases when either gain (or loss) dominate the dynamics by orders of magnitude. A subset of simulations are observed to degrade into noise due to spatial resolution issues. Pseudospectral methods rely on discrete fast Fourier transforms (dFFT) which provide excellent convergence for wellbehaved curves. However, if the spatial features of , the complex spin wave amplitude, approach the length of the numerical lattice spacing singularities may appear and the dFFT algorithms will no longer converge locally. These errors will continue to grow and propagate over time in an evolution. Further exploration of these cases with finer spatial resolutions is prohibited by computational resource constraints. In contrast, all results previously presented in sections 36 are converged in both space and time, as demonstrated in section 8.1 and 8.2.
9 Conclusions
We report the numerical identification of nine distinct long lifetime complex dynamical behaviors as part of six broad solution classes of the gainloss nonlinear Schrödinger equation (GLNLS), (1). Behaviors were located during an extensive numerical exploration of six dimensional parameter space. A minimum of eight decades were examined for each gain term while five decades of higher order nonlinearities were considered at fixed dispersion and cubic nonlinearity. The GLNLS served as a driven damped model of long lifetime spin wave dynamics in magnetic thin film active feedback rings and analogous driven damped nonlinear physical systems. Agreement of GLNLS low dimensional chaotic modulating bright soliton trains with experimental measurements [35] was discussed in detail. We predicted additional GLNLS dynamical behaviors including two distinct steady state solutions, four unique examples of stable dynamical pattern formation and the intricate spatially symmetric/asymmetric interactions of solitary wave peaks. Finally we reported the existence of intermittent regimes within GLNLS parameter space, a phenomena typical of chaotic dynamical systems. Two unique examples of intermittency were presented which demonstrated finite periods of two distinct dynamical behaviors. The variety of presented GLNLS solution types matches the scope of dynamical behaviors observed experimentally in YIG film spin wave systems, as well as predicting new behaviors that can be tested in present experiments. The GLNLS thus presents a simple yet viable and fundamental model for driven, damped nonlinear waves propagating in dispersive mediums.
We neglected the periodic effect of amplification within the feedback ring, so the gain and loss terms presented in this work represented averaged quantities. Highly variable solution types such as the symmetric and asymmetric interactions potentially violate the GLNLS operating regime, with gain and loss driven dynamics occurring on the scale of a single round trip. Future study of this limit and adiabatically driven soliton trains is warranted. A fine grained exploration of parameter space may also be justified to identify distinct domains of stability for each observed behavior. In the future a rigorous study of GLNLS phase space would be useful to determine the cause of intermittency and potentially locate chaotic attractors of higher dimension.
This material is based upon work supported under grants number NSF PHY0547845, NSF DMR0906489, NSF PHY1067973, and NSF PHY1207881.
References
Footnotes
 A typical expression for saturable gain is given by where and are control parameters. Expanding denominator to third order yields , hence positive quintic gain being named saturating.
References
 Igor S. Aranson and Lorenz Kramer. The world of the complex GinzburgLandau equation. Rev. Mod. Phys., 74:99–143, 2002.
 Catherine Sulem and PierreLouis Sulem. The nonlinear Schrödinger equation: Selffocusing and wave collapse, volume 139. Springer, 1999.
 B.A. Kalinikos and A.N. Slavin. Theory of dipoleexchange spinwave spectrum for ferromagneticfilms with mixed exchange boundaryconditions. J. Phys. C, 19(35):7013–7033, 1986.
 A.N. Slavin and B.A. Kalinikos. Nonlineartheory of spinwaves in ferromagnetic films. Sov. Phys. Tech. Phys., 57(12):2387–2389, 1987.
 Iacopo Carusotto and Cristiano Ciuti. Quantum fluids of light. Rev. Mod. Phys., 85:299–366, Feb 2013.
 Lev Petrovich Pitaevskii and Sandro Stringari. BoseEinstein condensation. Int. Ser. Mono. Phys. Clarendon Press, Oxford, 2003.
 Lincoln D Carr, David DeMille, Roman V Krems, and Jun Ye. Cold and ultracold molecules: science, technology and applications. New J. Phys., 11(5):055049, 2009.
 Robert E. Camley and Robert L. Stamps, editors. Nonlinear Spin Waves in Magnetic Film Feedback Rings, volume 62 of Sol. State Phys., pages 163–224. Academic Press, 2010.
 A.N. Slavin and I.V. Rojdestvenski. Bright and dark spinwave envelope solitons in magneticfilms. IEEE Tran. Magn., 30(1):37–45, 1994.
 M.M. Scott, M.P. Kostylev, B.A. Kalinikos, and C.E. Patton. Excitation of bright and dark envelope solitons for magnetostatic waves with attractive nonlinearity. Phys. Rev. B, 71(17), 2005.
 M.P. Kostylev. Selfgeneration of bright spinwave envelope solitons in active ferromagneticfilm rings. J. Commun. Tech. Elec., 50(3):313–320, 2005.
 B.A. Kalinikos, N.G. Kovshikov, and C.E. Patton. Excitation of bright and dark microwave magnetic envelope solitons in a resonant ring. Appl. Phys. Lett., 75(2):265–267, 1999.
 N.G. Kovshikov, B.A. Kalinikos, C.E. Patton, E.S. Wright, and J.M. Nash. Formation, propagation, reflection, and collision of microwave envelope solitons in Yttrium Iron Garnet films. Phys. Rev. B, 54(21):15210–15223, 1996.
 B.A. Kalinikos, M.M. Scott, and C.E. Patton. Selfgeneration of fundamental dark solitons in magnetic films. Phys. Rev. Lett., 84(20):4697–4700, 2000.
 B.A. Kalinikos, N.G. Kovshikov, and A.N. Slavin. Envelope solitons of highly dispersive and low dispersive spinwaves in magneticfilms. J. Appl. Phys., 69(8):5712–5717, 1991.
 B.A. Kalinikos, N.G. Kovshikov, and A.N. Slavin. Experimental observation of magnetostatic wave envelope solitons in YttriumIronGarnet films. Phys. Rev. B, 42(13):8658–8660, 1990.
 B.A. Kalinikos, N.G. Kovshikov, and A.N. Slavin. Spinwave envelope solitons in thin ferromagneticfilms (invited). J. Appl. Phys., 67(9):5633–5638, 1990.
 M. Chen, M.A. Tsankov, J.M. Nash, and C.E. Patton. Backwardvolumewave microwaveenvelope solitons in YttriumIronGarnet films. Phys. Rev. B, 49(18):12773–12790, 1994.
 B.A. Kalinikos, N.G. Kovshikov, and C.E. Patton. Observation of selfgeneration of dark envelope solitons for spin waves in ferromagnetic films. JETP Lett., 68(3):243–247, 1998.
 B.A. Kalinikos, N.G. Kovshikov, and A.N. Slavin. Effect of magnetic dissipation on propagation of dipole spinwave envelope solitons in YttriumIronGarnet films. IEEE Trans. Magn., 28(5, Part 2):3207–3209, 1992.
 H. Benner, B.A. Kalinikos, N.G. Kovshikov, and M.P. Kostylev. Observation of spinwave envelope dark solitons in ferromagnetic films. JETP Lett., 72(4):213–216, 2000.
 A.N. Slavin and H. Benner. Formation and propagation of spinwave envelope solitons in weakly dissipative ferrite waveguides. Phys. Rev. B, 67(17), 2003.
 M. Wu, B.A. Kalinikos, and C.E. Patton. Generation of dark and bright spin wave envelope soliton trains through selfmodulational instability in magnetic films. Phys. Rev. Lett., 93(15), 2004.
 B.A. Kalinikos, N.G. Kovshikov, M.P. Kostylev, and H. Benner. Selfgeneration of spinwave envelope soliton trains with different periods. JETP Lett., 76(5):253–257, 2002.
 B.A. Kalinikos, N.G. Kovshikov, and C.E. Patton. Selfgeneration of microwave magnetic envelope soliton trains in Yttrium Iron Garnet thin films. Phys. Rev. Lett., 80(19):4301–4304, 1998.
 S.O. Demokritov, A.A. Serga, V.E. Demidov, B. Hillebrands, M.P. Kostylev, and B.A. Kalinikos. Experimental observation of symmetrybreaking nonlinear modes in an active ring. Nature, 426(6963):159–162, 2003.
 M. Wu and C.E. Patton. Experimental Observation of FermiPastaUlam Recurrence in a Nonlinear Feedback Ring System. Phys. Rev. Lett., 98:047202, 2007.
 M.M. Scott, B.A. Kalinikos, and C.E. Patton. Spatial recurrence for nonlinear magnetostatic wave excitations. J. Appl. Phys., 94(9):5877–5880, 2003.
 M. Wu, B.A. Kalinikos, L.D. Carr, and C.E. Patton. Observation of spinwave soliton fractals in magnetic film active feedback rings. Phys. Rev. Lett., 96(18), 2006.
 M. Wu, P. Krivosik, B.A. Kalinikos, and C.E. Patton. Random generation of coherent solitary waves from incoherent waves. Phys. Rev. Lett., 96(22), 2006.
 M. Wu, B.A. Kalinikos, and C.E. Patton. Selfgeneration of chaotic solitary spin wave pulses in magnetic film active feedback rings. Phys. Rev. Lett., 95(23), 2005.
 Aaron M. Hagerstrom, Wei Tong, Mingzhong Wu, Boris A. Kalinikos, and Richard Eykholt. Excitation of chaotic spin waves in magnetic film feedback rings through threewave nonlinear interactions. Phys. Rev. Lett., 102(20), 2009.
 A. V. Kondrashov, A. B. Ustinov, B. A. Kalinikos, and H. Benner. Chaotic microwave selfgeneration in active rings based on ferromagnetic films. Tech. Phys. Lett., 34(6):492–494, 2008.
 M. Wu, M.A. Kraemer, M.M. Scott, C.E. Patton, and B.A. Kalinikos. Spatial evolution of multipeaked microwave magnetic envelope solitons in Yttrium Iron Garnet thin films. Phys. Rev. B, 70(5), 2004.
 Z Wang, A Hagerstrom, J.Q. Anderson, W. Tong, M. Wu, L.D. Carr, R. Eykholt, and B.A. Kalinikos. Chaotic spinwave solitons in magnetic film feedback rings. Phys. Rev. Lett., 107:114102, 2011.
 Alexey B. Ustinov, Vladislav E. Demidov, Alexander V. Kondrashov, Boris A. Kalinikos, and Sergej O. Demokritov. Observation of the chaotic spinwave soliton trains in magnetic films. Phys. Rev. Lett., 106:017201, Jan 2011.
 Mark J. Ablowitz, G Biondini, and S Blair. Nonlinear Schrödinger equations with mean terms in nonresonant multidimensional quadratic materials. Phys. Rev. E, 63(4, Part 2), 2001.
 Mark J. Ablowitz, S.A. Ablowitz, and N. Antar. Damping of periodic waves in physically significant wave systems. Stud. Appl. Math., 121(3):313–335, 2008.
 Mark J. Ablowitz and Theodoros P. Horikis. Pulse dynamics and solitons in modelocked lasers. Phys. Rev. A, 78(1), 2008.
 Mark J. Ablowitz, Theodoros P. Horikis, and Boaz Ilan. Solitons in dispersionmanaged modelocked lasers. Phys. Rev. A, 77(3), 2008.
 N. Akhmediev, J.M. SotoCrespo, and G. Town. Pulsating solitons, chaotic solitons, period doubling, and pulse coexistence in modelocked lasers: complex GinzburgLandau equation approach. Phys. Rev. E, 63(5):056602, 2001.
 N Akhmediev and Adrian Ankiewicz. Dissipative solitons in the complex ginzburglandau and swifthohenberg equations. In Dissipative solitons, pages 1–17. Springer, 2005.
 N. Akhmediev, J. M. SotoCrespo, and Ph. Grelu. Dissipative solitons and their interactions. PAMM, 7(1):1130301–1130302, 2007.
 N. Akhmediev and A. Ankiewicz. Dissipative Solitons: from Optics to Biology and Medicine, volume 751. Springer, 2008.
 N. Akhmediev, A. Ankiewicz, J.E.M.I.A. SotoCrespo, and P. Grelu. Dissipative solitons: present understanding, applications and new developments. Int. J. Bif. Chaos, 19(08):2621–2636, 2009.
 E.N. Tsoy and N. Akhmediev. Bifurcations from stationary to pulsating solitons in the cubic–quintic complex Ginzburg–Landau equation. Phys. Lett. A, 343(6):417–422, 2005.
 M.N. Zhuravlev and N.V. Ostrovskaya. Dynamics of NLS solitons described by the cubicquintic GinzburgLandau equation. J. Exp. Theo. Phys., 99(2):427–442, 2004.
 D.D. Stancil and A. Prabhakar. Spin Waves: Theory and Applications. Springer, 2009.
 Pavol Krivosik and Carl E. Patton. Hamiltonian formulation of nonlinear spinwave dynamics: Theory and applications. Phys. Rev. B, 82:184428, 2010.
 H. Leblond. Rigorous derivation of the NLS in magnetic films. J. Phys. A, 34(45):9687, 2001.
 Victor Snyder. Quantum Fluctutations in the Dynamics of BoseEinstein Condensates. PhD thesis, Colorado School of Mines, 2008.
 W.H. Press, S.A. Teukolsky, W.T. Vetterling, and B.P. Flannery. Numerical recipes 3rd edition: The art of scientific computing. Cambridge University Press, 2007.
 R Kosloff and H TalEzer. A direct relaxation method for calculating eigenfunctions and eigenvalues of the Schrödinger equation on a grid. Chem. Phys. Lett., 127(3):223–230, 1986.
 M.M. Scott, C.E. Patton, M.P. Kostylev, and B.A. Kalinikos. Nonlinear damping of highpower magnetostatic waves in YttriumIronGarnet films. J. Appl. Phys., 95(11, Part 1):6294–6301, JUN 2004.
 R. Hegger, H. Kantz, and T. Schreiber. Practical implementation of nonlinear time series methods: The TISEAN package. Chaos, 9(2):413–435, 1999.
 H. Kantz and T. Schreiber. Nonlinear time series analysis. Cambridge Univ Pr, 2004.
 Floris Takens. Detecting strange attractors in turbulence. In Dynamical systems and turbulence, pages 366–381. Springer, 1981.
 T. Sauer, J.A. Yorke, and M. Casdagli. Embedology. J. Stat. Phys., 65(3):579–616, 1991.
 J. Theiler. Estimating fractal dimension. J. Opt. Soc. Am. A, 7(6):1055–1073, 1990.
 Mark C Cross and Pierre C Hohenberg. Pattern formation outside of equilibrium. Rev. Mod. Phys,, 65(3):851, 1993.