On the limits of validity of nonparaxial propagation equations in Kerr media

Amnon Yariv California Institute of Technology 128-95, Pasadena, California 91125 Abstract: Recent generalizations of the standard nonlinear Schroedinger equation (NLSE), aimed at describing nonparaxial propagation in Kerr media are examined. An analysis of their limitations, based on available exact results for transverse electric (TE) and transverse magnetic (TM) (1+1)-D spatial solitons, is presented. Numerical stability analysis reveals that nonparaxial TM soltions are unstable to perturbations and tend to catastrophically collapse while TE solitons are stable even in the extreme nonparaxial limit. c © 2006 Optical Society of America


Introduction
The recent development of nanotechnology, where sub-wavelength features are typically present, force on us the necessity of studying optical propagation beyond the standard paraxial approximation which has served us so well up to now.In particular, this is true in the nonlinear propagation regime, and, more specifically, in the presence of optical Kerr effect where the problems of stability and catastrophic collapse play an important role which can be properly understood only in terms of nonparaxial effects.This can be done by starting ab initio from Maxwell's equations and solving them numerically or, more conveniently, by trying to derive a single approximate propagation equation which takes nonparaxial effects into account.Many attempts aimed at generalizing, as a first approximation, the paraxial approach by adding on the RHS of the standard (paraxial) nonlinear Schroedinger equation (NLSE) higher-orders terms of the order of the square of the smallness parameter = λ/σ (where λ and σ are the wavelength and the beam typical transverse dimension, respectively) have resulted in different equations which do not agree among themselves.While some of these nonparaxial nonlinear Schroedinger equations (NNLSE) contain intrinsic scalar approximations, two general approaches, [1,2]- [3,4], preserve the fully vectorial nature of the problem.The first approach is based on the analysis of light propagation in Fourier space, where dealing with nonparaxial effects is both simpler and more natural.The second approach recasts Maxwell equations into operatorial identities which are handled by means of a suitable iterative scheme, a process giving the field at any desired order in .The two approaches lead to two slightly different equations, one of which should represent the correct generalization of the standard nonlinear Schroedinger equation.Actually, both equations can be shown to be incapable of reproducing some exact results available in the literature about (1+1)-D spatial solitons as ab initio solutions of Maxwell's equations, both for TE [5] and TM fields.[6] This seems to indicate that writing a unique nonlinear nonparaxial propagation equation capable of going beyond the standard paraxial approximation may still be an open problem.

The nonlinear nonparaxial propagation equations
Let us recall briefly that, in the frame of the paraxial approximation, the propagating field can be described in terms of its two mutually orthogonal components, transverse to the propagation direction.The longitudinal component is assumed to be much smaller than the transverse one and can be evaluated once the transverse component is known.Accordingly, propagation is described by an equation for the two transverse components, the paraxial (or parabolic) equation, which is valid as far as the longitudinal component is small compared to the transverse one.In the specific case of nonlinear propagation in a Kerr medium, after writing the electric field E(x, y, z, t) propagating in the z-direction as E(x, y, z, t) = A(x, y, z) exp (ikz − iωt) , (1) where k = ωn 0 /c (n 0 being the refractive index of the background linear medium), the paraxial (2+1)-D equation for the transverse component reads (see, e.g., [2]) where 2 and n 2 is the nonlinear refractive-index coefficient.The longitudinal component A z can be approximately expressed by In the (1+1)-D case, it is possible to describe the evolution of a single linearlypolarized field component.For example, after introducing the normalized coordinates (ξ, η, ζ) = (kx, ky, kz) and the normalized amplitude U (ξ, ζ) = (|n 2 |/n 0 ) 1/2 A y (x, z) of the linearly y-polarized TE field component, Eq.( 2) yields which is the standard paraxial form of the (1+1)-D NLSE (an identical equation holds true for the x-polarized TM field component In order to generalize Eq.( 3) beyond the paraxial approximation one has to rely on the (2+1)-D nonparaxial propagation equations present in [1,2,3,4].The distinction between TE and TM field component is no longer meaningful since the propagation equation is a fully vectorial one.However, if only one of the two transverse components is initially different from zero and depends on a single transverse Cartesian coordinate (a geometry separating TE and TM fields), it is still possible to write a (1+1)-D equation for the single transverse polarization component.By inspecting these equations (see next Section), the natural generalization of Eq.( 3), that is the NNLSE which includes higher-order nonparaxial terms up to the second order in the smallness parameter, is of the form where a, b, c and d are suitable real coefficients.Note that the form of the above equation is the most general one compatible with the inclusion of terms up to the second order in , the invariance under the reflection ξ → −ξ and the description of forward propagating beams.Substituting into Eq.( 4) a spatial soliton solution of the form U (ξ, ζ) = exp (iβζ) u(ξ), with u and β real, yields where the prime stands for the derivative with respect to ξ and A = a + b, B = c + d.Equation ( 5) can be integrated by following a general procedure (see, e.g., [7]).After defining f = u 2 , so that df /dξ = 2u u = u df /du, Eq.( 5) becomes (6) Equation ( 6) can be readily integrated yielding where α = B/A and u 0 is a fixed value of u.By inspecting Eq.( 7), it is evident that the condition 1 + 2γAu 2 > 0 has to be fulfilled for the last term in the right hand side to be real for every value of α.This implies that, for γ = 1 (bright solitons), u 2 can take on any value as long as A > 0, while, for A < 0, one must have u 2 < 1/(2|A|).For γ = −1 (dark solitons), u 2 can take any value as long as A < 0, while, if A > 0, u 2 is restricted by the condition u 2 < 1/(2A).Taking the ξ-derivative of Eq.( 7) yields Equations ( 7) and ( 8) can be used to derive β.For bright solitons (γ = 1), Eq.( 8) is identically satisfied (0 = 0) by taking the limit ξ → ±∞; in the same limit, by choosing u 0 to be the peak amplitude of the soliton (so that f (u 0 ) = 0), Eq. ( 7) yields For dark solitons (γ = −1), by taking the limit ξ → ±∞ and u 0 = 0 in Eqs. ( 7) and (8), we obtain a set of equations in f (0) and β whose solution yields where ±u ∞ is the asymptotic value of u when ξ → ±∞.

Comparison of the available NNLS equations with the exact results
The propagation coefficients β obtained from the NNLSE Eq.( 4), can now be compared to the exact ones obtained, both for bright and dark solitons, by solving ab initio Maxwell's equations.For the TE polarization, the NNLSE obtained in [1,2] reads or, by using a prime to indicate differentiation with respect to ξ, The corresponding equation derived in [3,4] is For the TM polarization, the NNLSE derived in [1,2] is while the one derived in [3,4] reads By using the results of the previous Section, in particular Eqs. ( 9) and (10), we can now compare, both for bright (γ = 1) and dark (γ = −1) solitons, the values of the different propagation constants with the exact ones associated with the exact solutions.In particular, for exact the TE solitons the propagation coefficient is [5] and for bright and dark solitons respectively, while, for the exact TM solitons, [6] These expressions can be compared, to the lowest significant nonparaxial order in u 4 0 or u 4 ∞ , with those obtained for the most general NNLSE (see Eqs.( 9) and ( 10)), where the coefficients A and B are deduced by inspecting Eqs.(12,13,14,15) (note that A and B are both negative for the TE polarization and both positive for the TM one).
For dark solitons, for both TE and TM polarizations, the NNLSE's are unable to produce the nonparaxial correction terms of the order u 4  ∞ predicted by the exact results.The above considerations imply that, in the frame of bright soliton propagation, three out of the four possible NNLSE's (two for TE and one for TM configurations) are unable to yield the relation between β and u 0 (or u ∞ ) predicted by the exact theory.Only Eq.(15) provides the correct expression.For dark solitons, both TE and TM equations are unable to provide the correct relation between the propagation constant and the asymptotic amplitude.The differences between the approximated equations derived in [1,2] and [3,4] results from the different adopted approaches.Basically, the scheme of Ref. [1,2] relies on a detailed analysis of Maxwell equations in the transverse Fourier space.After the dominant contributions are recognized, the relevant propagation equations are derived by inverse fourier transform.On the other hand, the approach of Ref. [3,4], representing a generalization of Lax et.al [8], is based on a suitable iteration scheme which extracts from Maxwell equations a propagation equation describing the dynamics to the relevant order of asymptotic expansion.The main difference between the two schemes lies in the method in which the longitudinal component E z of the electromagnetic field is dealt with.In the former approach [1,2] the exact Fourier spectrum of E z is exploited to obtain exact equations for the transverse part (E x , E y ) that are subsequently approximated taking into account the paraxial expansion.In the approach of [3,4], the equations for the longitudinal and transverse components of the electromagnetic field are asymptotically expanded in parallel.

Stability analysis of nonparaxial, exact, TE and TM solitons
In this section we present a comparison between the spatial stability properties of the approximate nonparaxial solitons and of the temporal stability properties of the exact ones, leading to the conclusion that TE nonparaxial solitons are stable while TM are not.
Figures 1(a)-1(d) show respectively the numerical solutions of Eqs.( 12), ( 13), ( 14), (15) obtained by beam propagation method (BPM) simulation, the initial conditions for the numerical analysis being the pertinent analytical soliton profile of the corresponding equation with initial (peak) amplitude u 0 = 0.42.Note that while both equations for the TE case exhibit stable propagation of their soliton solutions, the TM solutions of Eqs.( 14) and ( 15) are unstable and undergo catastrophic collapse as they propagate.This difference stems from the fact that all the coefficients of the nonlinear nonparaxial terms have opposite sign in the cases of TE and TM (see Eqs.( 12), (13) and Eqs.( 14), ( 15)).
We have investigated the stability characteristics of the exact TE and TM solutions of the Maxwell equations using finite difference time domain (FDTD) simulations [9].The propagation of solitons with an amplitude of A 0 = 0.4 (i.e., well beyond the non-paraxial approximation) for the TE and TM cases is shown respectively in Figs.2(a) and 2(b).The spatial resolution of the FDTD simulations was λ/100 and the initial conditions were the exact soliton solutions for the TE and TM cases derived in [6].The FDTD simulations confirm that also the exact TM soliton solution of the Maxwell equation is unstable.As the TM beam propagates, temporal and longitudinal oscillations evolve, eventually leading to a catastrophic collapse of the beam (see Fig. 2(b)).On the contrary, the TE soliton solution (Fig. 2(a)) appears to be stable even for highly non-paraxial initial conditions.
It should be emphasized that temporal instabilities cannot be deduced from harmonic propagation equations such as Eqs.( 2) and (4) because they incorporate the assumptions of a single frequency and slowly varying envelope in z.To our knowledge, such instabilities in spatial Kerr solitons were not studied or observed previously.

Conclusions
Current versions of the nonparaxial nonlinear Schroedinger equation aimed at generalizing its standard paraxial counterpart have been tested by comparing some of their predictions about bright and dark soliton propagation with exact results available in the literature.Only one of the different versions, and only in the case of bright TM soliton, provide the correct relation between the propagation constant and the soliton peak amplitude.These results seems to indicate that the problem of generalizing the NLSE is still an open one which requires a more rigorous analysis of the asymptotic approaches that has led to the current versions.Moreover, the existence of temporal and longitudinal instabilities (at least in the TM case) indicate that generalized nonlinear propagation equation do not provide a complete description of the field evolution.However, the available versions of the NNLSE's provide results that are close enough to the exact ones to justify their use in a qualitative basis for describing nonparaxial propagation features which would otherwise require the full solution of Maxwell's equations.

Fig. 2 .
Fig. 2. (a) (2.04 MB) Movie of FDTD simulation showing the stable propagation of a TE polarized nonparaxial soliton with initial amplitude of A = 0.4.(b) (2.37 MB) Movie of FDTD simulation showing the unstable propagation of a TM polarized nonparaxial soliton with initial amplitude of A = 0.4