Stable solution of induced modulation instability

In this paper,we discussed the nonlinear evolution of modulation instability and the steady-state process of induced modulation instability in sine-oscillatory response nonlocal nonlinear media. With plane wave plus perturbation as initial conditions, we simulated the long-term evolution of modulation instability in the nonlocal nonlinear Schrodinger equation with sine-oscillatory response numerically. For the input of modulated wave, the approximate analytical solution of the stable solution of the equation is obtained under the assumption that only the fundamental wave and the first harmonic wave are present. For the input of modulated wave with arbitrary harmonic waves, we obtained the exact numerical solution of the stable solution of the induced modulation instability.

Modulation instability (MI) is considered as an important problem in nonlinear physics. A series of studies on MI have been carried out since 1960s, such as fluids 1 , plasma 2 , nonlinear optics 3-6 , among others. MI signifies the exponential growth of a weak perturbation of the amplitude of the wave as it propagates. The gain leads to amplification of sidebands, which breaks up the otherwise uniform wave front and generates fine localized structures. MI means the energy exchange between waves at different frequencies in essence, especially referring to the energy exchange between plane wave and perturbation at high-order harmonics 6 . The energy of plane wave is transformed into that of perturbation at high-order harmonics. Therefore, perturbation is amplified to lead to the generation of localized structures on plane wave during transmission. Theories have validated that perturbation can only trigger MI in a specific frequency range 3,4 . Reported was the first experimental observation of the modulational instability using single-mode fibers 5 . In nonlinear Kerr media, bright solitons can survive within the frequency range in which MI is generated; while dark soliton exist in the frequency range in which MI does not occur 4,7,8 . On the other hand, MI of a nonlinear system leads to the growth of amplitude of modulated wave when the modulated frequency of the input signals of the modulated wave within the frequency range in which MI is generated. The process is called induced modulation instability 3,9,10 . Hasegawa first proposed a theoretical scheme of generating ultra-high-rate pulse trains by utilizing induced modulation instability in fibers 3 , which has been experimentally verified by Tai et al. 10 .
A majority of researches on MI focus on discussing the linear process of MI during short-term transmission, that is processing the condition under small perturbation by using linear approximation method while ignoring the perturbation of high-order harmonics. However, the increment of perturbation is not much lower than the amplitude of plane wave after optical field evolves for a long term, so linear approximation method is not applicable. It is necessary to come up with new methods to discuss the nonlinear process when the increment of perturbation is not small after perturbation evolving for a long term. Beeckman first explored the nonlinear process of induced modulation instability in nonlocal nonlinear Kerr media during long-term evolution by using analytical and numerical methods 11 . Beeckman mainly discussed the nonlinear process of induced modulation instability in Kerr media with Gaussian response function [12][13][14][15][16] and exponent response function [17][18][19] and found that the evolution of optical fields exhibited quasiperiodicity in the process.
Sine-oscillatory response function is also a common one apart from Gaussian and exponent response functions of nonlocal nonlinear Kerr media, which was first put forward by Nikolov during the discussion of quadratic soliton 20 . Recently, researches have shown that it is feasible to attain a sine-oscillatory response function when the boundary of the system satisfies specific conditions in a bounded system 21 . Linear analysis has been carried out on MI in nonlocal nonlinear Kerr media with sine-oscillatory response 6 . In this paper, we discuss the long-term (nonlinear) evolution of the optical field under MI in nonlocal nonlinear Schrodinger equation with sine-oscillatory response through numerical method, and the characteristics of the evolution was also discussed. The stable solutions (approximate analytical solution and exact numerical solution) under induced modulation instability of nonlocal nonlinear Schrodinger equation with sine-oscillatory response was also found and the relationship between the beamwidth in the stable solution and two degrees of freedom (DOF) in the system was surveyed.

nonlinear evolutionary process of Mi
The evolution of the linearly polarized electric field transmitted along z-axis in 1 + 1 dimensional media with sine-oscillatory response satisfies nonlocal nonlinear Schrödinger equation 6,[21][22][23] where u(x, z) is the dimensionless slowly-varying complex amplitude of the optical field, x and z are, respectively, the dimensionless transverse coordinate and the dimensionless evolution coordinate, and the sine-oscillatory response function is expressed as follows 20,21 with w m being the characteristic length of the nonlinear response. Mathematically, the z coordinate stand for time in equation above 24,25 , but it is, physically, the longitudinal space (beam-propagation direction) coordinate for the problem of the optical-beam-propagation [26][27][28] . Previous research has been conducted on MI in nonlocal nonlinear Kerr media with sine-oscillatory response and mainly discusses the linear process of MI 6 . However, linear approximation method is not suitable owing to MI has been nonlinear when the optical field undergoes long-term evolution and the growing increment of perturbation cannot be ignored. Therefore, the nonlinear process of MI during long-term evolution would be discussed by employing numerical method. By taking infinite plane wave with perturbation as the initial input of Eq. (1), where the first term, the second term and k x represent plane wave, perturbation and frequency of perturbation, respectively. According to the previous research result of linear analysis on MI, it can be seen that the MI in the system described by Eq. (1) appears in the range of w k I w 1 (4 where I 0 denotes the intensity of plane wave) 6 . The amplitude of the perturbation will grow in the exponential form e gz in the initial evolution stage when the frequency of the perturbation is within the range, where g refers to the linear gain coefficient of perturbation. At different frequencies, g is expressed as x m in the range in which MI occurs to the system, whose corresponding frequency is the singular point in the range; the corresponding frequency at = g 0 is the cutoff point in the range in MI happens. It can be seen from Eq. (4) that a perturbation whose frequency is closer to the singular point shows a larger linear gain coefficient in Kerr media with sine-oscillatory response; by contrast, a perturbation whose frequency is closer to the cutoff frequency exhibits a lower linear gain coefficient.
MI appears in the range of < < . k 1 1 61 x is the singular point in the range in which MI occurs, at which linear gain coefficient is infinite; and = . k 1 61 x is the cutoff point in the range. By taking Eq. (3) as initial condition, the evolutionary process of the optical field under MI in media was simulated numerically, as shown in Fig. 1. In Fig. 1(a), it can be observed that the amplitude of perturbation exponentially grew in the short-term evolutionary process of the optical field, which validated the result of linear analysis on MI. However, it was impossible for perturbation to infinitely constantly increase. Linear approximation method was not applicable any more when the increment of perturbation grew to a point that cannot be ignored with the increase of evolution distance. It can be found from Fig. 1(b) that after growing in the initial evolution stage, perturbation started to attenuate and then evolved in a periodic-like oscillation form. www.nature.com/scientificreports www.nature.com/scientificreports/ A comparison was made between the analytical result in the evolutionary process of the optical field treated through linear approximation method and numerically simulated result without linear processing, as shown in Fig. 2. The solid and imaginary lines in the figure separately represent the analytical results and simulated results when inputting perturbation of different frequency. Linear approximation method was still applicable when perturbation was small in the short-term evolution of the optical field. In this case, the curve of intensity obtained through linear processing was favorably consistent with that attained through numerical simulation; linear approximation method was not applicable after the increment of perturbation grew to a large scale in the long-term evolution. Therefore, the evolution trend of the optical field obtained through linear approximation method gradually deviated from the simulated result. Additionally, it can be seen from Fig. 2 that the perturbation whose frequency was closer to the singular point exhibited faster increasing amplitude, more quickly reached the first peak and showed a larger peak; by contrast, the perturbation whose frequency was farther away from the singular point presented slower growth of the amplitude, reached the first peak later and attained a lower peak. It can be also seen in the Fig. 2 that the optical field was more irregularly evolved under the perturbation whose frequency was closer to the singular point.
Previous research indicated that perturbations would increase when their frequency fell in the range in which MI appeared; however, perturbations would be steadily transmitted and did not grow when their frequency was not found in the range 6 . It can be attained through numerical calculation that perturbations did not rise in the initial stage of short-term evolution of optical field when the frequency of perturbations was beyond the above range ( Fig. 3(a,c)). In this case, optical field evolved in an oscillation form. However, high-harmonic generations (e.g. second-and third-harmonic generations) of perturbations were generated in the system after optical field was subjected to a long-term evolution. MI was also found when high-harmonic generation of perturbations appeared in the range of the system in which MI occurred ( Fig. 3(b,d)). Figure 3(b,d) show MI separately induced by second-and third-harmonic generations of perturbations. On this basis, it can be seen that MI was also generated when high-harmonic generations of perturbations were found within the range in which MI appeared.

Stable solution of induced modulation instability
Approximate analytical solution of stable solution. To discuss the stable solution under induced modulation instability of Schrodinger equation in 1 + 1-dimensional media with sine-oscillatory response, it was necessary to introduce initial input www.nature.com/scientificreports www.nature.com/scientificreports/ ∑ = .
Based on Eq. (6), various harmonic separately satisfy the following equation It is supposed that only plane wave and first harmonic are taken into account, i.e. only having a 0 , a 1 and − a 1 . According to Eq. (9), plane wave and first harmonic are separately obtained as ⁎ ⁎ a ak a a a a a a a a w k x m x   Moreover, the intensity of plane wave in the optical field is By substituting Eqs. (14) and (15) into Eq. (10), the propagation constant of the optical field is attained as When only considering plane wave and first harmonic in the evolutionary process of the optical field, the approximate analytical solution of stable solution can be obtained as Simulation was carried out by taking the result of approximate analytical solution obtained using Eq. (17) at = z 0 as the initial input of Eq. (1). The evolution of the optical field based on the equation under the input condition is displayed in Fig. 4. As shown in the figure, not all approximate analytical solutions at any modulation frequency can be steadily transmitted and only the approximate analytical solutions at the modulation frequency approximate to the cutoff frequency are able to be steadily transmitted ( Fig. 4(a)); by contrast, the approximate analytical solutions at the modulation frequency greatly departing from the cutoff frequency greatly fluctuated in the transmission process ( Fig. 4(b)). It indicated that the approximate analytical solution obtained on condition of only considering plane wave and first harmonic was applicable when modulation frequency approached to the cutoff frequency, while was not suitable when modulation frequency greatly departed from the cutoff frequency. exact numerical solution of stable solution. The above approximate analytical solution was obtained in the case of only considering plane wave and first harmonic, and the approximation was suitable when the modulation frequency approached to the cutoff frequency. However, it was necessary to consider more number of harmonics during calculation when modulation frequency significantly departed from the cutoff frequency. Whereas, it was impossible to calculate the analytical solution when considering more enough harmonics, so the calculation was performed through numerical method. Based on normalizing condition calculation was conducted on Eq. (9) by using Newton iteration method to calculate the values of a n and β 29 . By taking a 0 , a 1 , and β in approximate analytical solutions as initial values, calculation was performed by constantly increasing the number of harmonics (that is increasing the number of equations) beginning with two harmonics. Iterative process was ended when harmonics with amplitude lower than 10 −4 appeared. In this way, the numerical solution of stable solution can be obtained as  Table 1. Additionally, a 0 , a 1 and β in approximate analytical solutions obtained when only considering plane wave and first harmonic are also listed in Table 1.
It can be seen from Table 1 that the approximate analytical solution was favorably consistent with the numerical solution when modulation frequency approached to the cutoff frequency within the range in which MI appeared; however, there was a large difference between the two solutions when modulation frequency greatly away from the cutoff frequency. To validate whether the numerical solutions can be steadily transmitted, transmission was conducted by taking the numerical solutions as the initial input of Eq. (1), that is  Figure 5 shows the simulated results of numerical solutions at different modulation frequencies in the evolution of the optical field. The figure revealed that the numerical solutions all can be steadily transmitted, which indicated that exact stable solutions can be attained by using the numerical method. By comparing Fig. 4(a) with 5(a), the approximate analytical solutions and numerical solutions all can be stably transmitted when the modulation frequency was approached to the cutoff frequency. It implied that the approximation method only considering plane wave and first harmonic was applicable in the case that the modulation frequency approached to the cutoff frequency; however, it was not suitable to the condition that the modulation frequency away from the cutoff frequency by comparing Fig. 4(b) with 5(b). It was also found from Table 1 that the difference between approximate analytical solutions and numerical solutions became increasingly larger, and therefore there were more harmonics in stable solutions when the modulation frequency more significantly departed from the cutoff frequency. www.nature.com/scientificreports www.nature.com/scientificreports/ Beamwidth of stable solution. According to Eq. (9), it can be obtained that there were n equations for describing the system when considering n harmonics and the equation set contained + n 3 unknown quantities: amplitude … − a a a , , , n 0 1 1 , nonlinear characteristic length w m , propagation constant β and modulation frequency k x . Therefore, three DOFs were present in the system. Hence, three conditions {Relation (18) [that is, supposing the power of stable solution as = P 1( = ∑ =− P a n n n n 2 and k x } were required when calculating stable solution in the section above (Exact numerical solution of stable solution). In this section, the characteristics of the stable solution with three DOFs was discussed. The beamwidth w r of stable solution within a single period was defined as to characterize the stable solution and discuss the change relationships of w r with three DOFs (P, w m and k x ). At first, the change laws of w r in stable solution with the other two DOFs w m and P were first explored on condition of keeping k x unchanged. The stable solution under different conditions were calculated through numerical method and the curve of w r with varying P under different w m was attained, as shown in Fig. 6. As the range in which MI appeared was related to w m , the selected values of w m in Fig. 6 were relatively approximated, so that the ranges in which MI appeared obtained under different values of w m showed a larger intersection (for convenience of comparison). As shown in the figure, w r decreased with increasing P under same k x and w m . The reason was that the cutoff frequency within the range in which MI occurred also grew with the increase of the power P 6 . It meant that k x would gradually depart from the cutoff frequency within the range in which MI happened with increasing P. Therefore, the number of harmonic in stable solution which needed to be considered rose, that is the spectral width of the stable solution grew, so that the beamwidth declined correspondingly. Similarly, the cutoff frequency also increased with the reduction of w m in the case of having same k x and P. In a similar way, the beamwidth was lowered. Figure 6 displays the relationships of the beamwidth of the stable solution with power and nonlinear characteristic length at a given modulation frequency. The change law of beamwidth can be attained through scale transform at different modulation frequencies, instead of numerically calculating the equation again. Nonlocal nonlinear Schrodinger Eq. (1) exhibited transformation invariance, so the invariant transformation 23 that Eq. (9) satisfied can be obtained It meant that Eq. (9) remained invariant via the transform (22). Therefore, the transformation relation of Σ 1 and Σ 2 in stable solution before and after transformation of the modulation frequency k x can be attained If the stable solution at a frequency k x had been calculated, the stable solution at any frequency k x can be acquired through the transform (23).
The stable solution of Schrodinger equation in nonlinear Kerr media with local, Gaussian and exponent responses have been found in some researches. Compared with the stable solution obtained in the study, the stable solution in these systems exhibited same characteristics, as shown in Fig. 7. The beamwidth w r of the stable solution decreased with increasing P while rose with the growth of w m within a single period on the premise of keeping k x unchanged. Local response corresponded to special nonlocal response at nonlocal characteristic www.nature.com/scientificreports www.nature.com/scientificreports/ length = w 0 m . Therefore, − w P r curve corresponding to nonlinear Kerr media with local response was found at the lowest part in Fig. 7. conclusions Two parts in a nonlocal nonlinear system with sine-oscillatory response were explored. On the one hand, it is nonlinear evolution of MI. The long-term evolution of optical field under MI of nonlocal nonlinear Schrodinger equation with sine-oscillatory response was numerically simulated by taking plane wave with perturbation as initial condition. Some characteristics in this process were obtained by simulating long-term nonlinear evolution of the optical field: firstly, perturbation started to attenuate after its amplitude grew to a certain value if its frequency was within the range in which MI occurred, and showing a quasi-periodic evolution in the long-term nonlinear process; secondly, perturbation whose frequency more approached to the singular point showed faster increasing amplitude, more quickly reached the first peak of its increment and exhibited a larger peak, moreover, the evolution of the optical field was increasingly irregular; thirdly, MI was also happened after long-term evolution if high-harmonic of perturbation were in the range in which MI appeared.
On the other hand, it is the steady-state process of induced modulation instability. The approximate analytical solution of stable solution under the induced modulation instability of nonlocal nonlinear equation with sine-oscillatory response was attained on condition of assuming only inputting plane wave and first harmonic (modulated wave); the exact numerical solution of stable solution was obtained when inputting the modulated wave with innumerable harmonics. By further discussing the stable solution, it can be found that fewer harmonics needed to be taken into account when calculating a stable solution if the modulation frequency of the stable solution was more approximated to the cutoff frequency in the range in which MI appeared. Additionally, the relationship of the beamwidth in the stable solution with two DOFs (power and nonlinear characteristic length) in the system were attained. On same conditions, the beamwidth reduced with increasing power while rose with the growth of nonlinear characteristic length.

Methods
We used the split-step Fourier method 4 to simulate the propagation of the optical beams. For the inputs Eqs. (3) and (20) at = z 0, the evolution of Eq. (1) was obtained by the split-step Fourier method. We used the Newton iteration method 29 to obtain the exact numerical solution of stable solution, i.e. the solution of Eq. (9). By taking a 0 , a 1 , and β in approximate analytical solutions as initial values, calculation was performed by the increase of the number of harmonics n beginning with = n 0, and 1. Iterative process was ended when the amplitude of the last harmonics was lower than 10 −4 . www.nature.com/scientificreports www.nature.com/scientificreports/