Energy shedding during nonlinear self-focusing of optical beams

: Self-focusing of intense laser beams and pulses of light in real nonlinear media is in general accompanied by material losses that require corrections to the conservative Nonlinear Schr¨odinger equations describing their propagation. Here we examine loss mechanisms that exist even in lossless media and are caused by shedding of energy away from the self-trapping beam making it to relax to an exact solution of lower energy. Using the conservative NLS equations with absorbing boundary conditions we show that energy shedding not only occurs during the initial reshaping process but also during oscillatory propagation induced by saturation of the nonlinear effect. For pulsed input we also show that, depending on the sign and magnitude of dispersion, pulse splitting, energy shedding, collapse or stable self-focusing may result.


Introduction
Self-focusing of light in nonlinear media is a well studied phenomenon first predicted in [1] and modelled during the 1960s [2]. Since Abbe's diffraction theory limits the minimum spot size of a focused beam, self-focusing offers potential means to generate beams that are sub-diffraction limited. This has potential uses in fields such as optical imaging and micromachining. In general there are mechanisms of loss during beam propagation in real media such as multiphoton excitation [3], and in many media the critical power for self-focusing is above the ionisation threshold of the medium leading to optical damage or plasma defocusing, thus arresting the collapse of the beam. Some media, in particular organic compounds, have significantly large third-order optical nonlinearities [4] and present the possibility of stable, self-focused propagation below the diffraction limit with limited loss due to material processes. Even in idealized lossless media, however, there will still be the shedding of energy during entry to and propagation in the medium, occurring in different forms depending on the particular balance of diffraction, dispersion, nonlinearity and saturation of the nonlinear effect.
In this paper we examine forms of energy shedding in lossless media using the Nonlinear Schrödinger equation (NLS) and show that, even though the NLS is conservative, in real media some energy will always be moved away from any self-trapped solution except in the case where the input beam matches the soliton parameters that correspond to the exact solutions of the equation. The energy shedding occurs on entry to the nonlinear medium due to the well known process of reshaping of the beam to the characteristic sech 2 profile [2]. Depending on the level of saturation of the nonlinear effect there may also be energy shedding during propagation or even symmetry breaking of the beam depending on the choice of parameters [5]. With pulsed input, dispersion also plays a role in determining the forms of energy shedding, with pulse splitting along the temporal axis and subsequent collapse occurring during propagation in normally dispersive media [6] compared to potentially stable spatio-temporal self-focusing occurring in anomalously dispersive media [7]. We quantitatively examine the levels of energy shedding in these different regimes to examine the most appropriate choice of conditions for stable self-focused propagation of cw and pulsed input beams.

Simulation Model
In dielectric media self-focusing occurs due to a local change in the refractive index (n) caused by an applied electric field far from resonance inducing a dipole moment on the molecules or atoms and resulting in an electrostrictive effect that causes a localized increase in the photon density. For a medium with a sufficiently large, positive Kerr coefficient (n 2 ), diffraction of the light is balanced and the beam becomes self-trapped when the power reaches a critical value, P cr = αλ 2 /(4πn 0 n 2 ) [8]. Here λ is the wavelength of the light and the parameter α depends on the beam geometry (for a Gaussian profile α ≈ 1.8962). n 0 is the linear refractive index, with n = n 0 in linear media, and n 2 acting as a small intensity dependent correction to the refractive index in nonlinear media: n = n 0 + n 2 I. Notably the critical power does not depend on the intensity and self-focusing occurs when the peak power exceeds P cr . A widely accepted model for the evolution of the complex envelope of the electric field under the slowly varying envelope approximation and for 1D cw and short pulse regimes (> 100 f s) is the cubic NLS: where ζ is the retarded time coordinate normalized by the Rayleigh length (Z R ) and η is the normalized spatial variable. E is the electric field normalised by the square root of the input peak intensity, I 0 , such that the input power is normalized to unity, i.e.
for input field E 0 where x = η in the 1D case, x = (η, ξ ) in the 2D case. The nonlinear coefficient is γ = k 0 n 2 Z R I 0 /n 0 , where k 0 is the free space wavenumber. The first term on the RHS of (1) describes diffraction and the second term describes an instantaneous third-order nonlinearity.
In 2D, however, the NLS is not suitable for simulation without some changes as it allows unbounded collapse within a finite distance [9]. In order to avoid this problem the saturable Nonlinear Schrödinger equation is used, with saturation of the nonlinearity derived from Maxwell-Bloch equations included in the model to prevent collapse. This makes the model more physical and results into: where ξ is the second transverse spatial coordinate and σ = I 0 /I sat is the saturation parameter for saturation intensity I sat . We can see that as σ tends towards zero, the cubic nonlinearity is regained. In order to model spatio-temporal input propagating through nonlinear media, the generalized Nonlinear Schrödinger equation is required [10] since it includes dispersion phenomenologically derived from the Taylor expansion of the wavenumber in Fourier space: where T is the normalized temporal variable and β = k 2 Z R /T 2 0 is the normalized group velocity dispersion (GVD) coefficient for input pulse length T 0 . This new term accounts for dispersion in the medium with positive β giving normal dispersion and negative β providing anomalous dispersion.
In all cases an initial Gaussian profile is used corresponding to TEM 00 laser modes focused to a diffraction limited spot at the edge of the nonlinear medium before propagation. The equation is propagated forward in time by use of a split-step method, with the linear operators propagated in the frequency domain and the nonlinearity evaluated on the alternate steps by use of a Runge-Kutta method. Transverse absorbing boundary conditions are included at the edge of the computational domain in order to remove shed energy from the simulation and avoid unphysical reflections back towards the propagating beam.

Cubic Nonlinearity
On entry to a linear medium a 1D cw beam diffracts as expected. Diffraction, however, is increasingly slowed as the nonlinear coefficient is increased up to a value of γ = 5.8 at which point the beam has sufficient power to self-trap. This is illustrated in Fig. 1(a) where it can be seen that diffraction is balanced with self-focusing leaving a stable, self-trapped beam after an initial period of reshaping. Figure 1(b) shows the self-trapping on a log intensity scale and shows that there are low levels of energy shed away from the self-trapped beam during the reshaping. This initial shedding always occurs for Gaussian input, accounting for approximately 2% energy loss at the point of self-trapping. As Fig. 2 shows, the amount of shed energy increases as the nonlinear coefficient increases, up to 25% at γ = 19.2. At this point the first higher-order oscillating soliton solution and a second self-trapped state are reached. Figure 2 also shows that increasing the nonlinear coefficient leads to a decreasing beam width between the first self-trapped state and the second, reducing the width to half the initial, diffraction limited, size at γ = 19.2.

Saturating Nonlinearity
Significantly different behavior occurs when the saturating NLS [Eq. (3)] is used, even with only one spatial variable. Although the behavior is very similar to the cubic 1D case up to the point of self-trapping, self-focusing in the saturating case does not result in a constant soliton solution by its strictest definition since the equation is no longer integrable [11,12]. Rather than reshaping and relaxing to a constant profile, the beam has a long lived oscillation in amplitude [ Fig. 3] and width [ Fig. 4(a)], undergoing periodic focusing and defocusing. Figure 3 shows the evolution of the amplitude of the beam during propagation as well as the total power. It can be seen that rather than just shedding energy during the initial reshaping process, there is now a constant shedding of energy during propagation that occurs during the defocusing cycles of the oscillation [see Fig. 4(b)]. As predicted in [13,14], depending on σ and the input power, there are two possible soliton solution branches that have the same width but different amplitudes. The lower branch exists for input power corresponding to P cr , with solutions below the branch untrapped and solutions above and close to the branch radiating energy during defocusing cycles until the power has been reduced to that of the lower branch. If the solution is far above the lower branch (input power sufficient to excite the first higher-order solution) the evolution displays symmetry breaking instabilities [see Fig. 5]. The solutions of the upper branch are not of as much interest in terms of stable sub diffraction limited focusing, as they require significantly more power, thus increasing the impact of modulation instability [12]. The effect of an increase of the saturation parameter is to increase the magnitude of nonlinear coefficient required to get the same level of self-focusing. Similar behavior occurs for each level of saturation though on different length scales, with the oscillations taking longer to damp during propagation in higher saturation regimes [see Fig. 6]. Notably it also has an effect on the critical power required for self-focusing to dominate over diffraction, with the point of self-trapping requiring an increased nonlinear coefficient for increased saturation.   order solutions. The effect of varying σ is the same in 2D as in 1D -higher saturation requiring a larger nonlinear coefficient to get the same results. There are also distinct differences. The magnitude of γ required to cause a particular level of self-focusing is larger in 2D and the increase in intensity during focusing is significantly larger than for 1D [see Fig. 7]. The beam focusing is also more significant and when the oscillations have decreased to negligible levels after propagation over 100s of Rayleigh lengths, the beam width may be as small as (1/10)w 0 [see Fig. 8].

Zero GVD Regime
Equation (4) was first integrated with β and σ = 0. The resulting spatial profile evolution is similar to an equivalent 1D simulation. The temporal profile does not spread or focus as the pulse propagates, as expected in cases of no dispersion, but it does exhibit an oscillation in the intensity of the centre of the pulse as a consequence of the transient oscillation in the transverse spatial dimension after entering the medium [see Fig. 9]. Spatially, the entire pulse is not compressed. Only those transverse 'slices' that have power above the critical level experience self-focusing, so that the front and back of the pulse diffract and only the middle section is trapped [see Fig. 10].
The inclusion of saturation (σ > 0) does not qualitatively change the results as the focusing is in effect 1D and there is no danger of collapse. The only noticeable effect is that the intensity does not reach the same peak levels, resulting in a slower dynamic. After propagation for around  and resulting beam width (red) dependence on the nonlinear coefficient for a 2D beam in saturable media (σ = 0.1). As in the 1D case, increasing the nonlinear coefficient causes first self trapping (γ = 12.4) and then self-focusing. Unlike the 1D case the resulting beam width is significantly smaller, reducing to as much a 10% of the input, and the power lost due to energy shedding is increased, varying from 4% at the point of self-trapping up to 40%. Varying the saturation parameter again changes the length scales on which the same behavior occurs.
100 Rayleigh lengths the pulse has settled to a steady shape with both temporal and spatial sections displaying narrowing [see Fig. 11] with up to 33% of the energy having been shed. The temporal length is decreased to roughly 50% by the shedding of untrapped energy at the front and rear of the pulse and the spatial width is reduced to about 1/3 of the input by selffocusing. Fig. 9. Left -evolution of the pulse profile on the temporal axis. Right -evolution of the profile on the spatial axis on entering a nonlinear medium. Without GVD there is no focusing or defocusing on the temporal axis, the splitting and reforming of the peak is due to oscillation along the spatial axis. Spatial evolution is similar to that of an equivalent 1D spatial simulation. Fig. 10. Left -initial pulse shape. Right -pulse shape after propagation to a stable profile. The front and back of the pulse have diffracted away, while the transverse slices with enough power to self-focus have become trapped.

Normal GVD Regime
Equation (4) was next integrated with positive β to provide normal dispersion during propagation. The results show that with normal GVD acting on the pulse as it propagates there is no approach to a stable steady state [see Fig. 12]. In normally dispersive media, the spreading of the pulse in conjunction with the nonlinearity induced self-phase modulation causes the pulse Fig. 11. Left -initial (red dash) and final (blue) temporal cross section. Right -initial and final spatial cross section. The pulse has narrowed in both dimensions during propagation, due to self-focusing in the spatial dimension, resulting in the sech 2 distribution, but due to untrapped energy being shed in the temporal dimension, resulting in a parabolic profile.
to split along the temporal axis due to the new frequencies being created at the front of the pulse travelling faster, and new frequencies at the back of the pulse travelling slower. Each time the pulse splits it reduces the power in the central spatial regions until it is below the critical level for self-focusing everywhere. The effect of increasing the magnitude of γ is only to prolong the duration of the splitting behavior. Inclusion of a saturating nonlinearity did not change the results significantly as at no point does the intensity reach a level at which the saturation may be relevant.

Anomalous GVD Regime
In the anomalous GVD regime there are significantly different results to normal or zero GVD. In our numerical simulations we identify a regime featuring focusing nonlinearity and anomalous dispersion in which a narrow and intense pulse can be formed. Figure 13 shows pulse propagation in the regime where diffraction and dispersion lengths are equal (i.e. |β | = L di f /L dis = 1).
As the effects balance in this regime, the pulse envelope evolves stably and symmetrically. The results are similar to 2D spatial propagation of a cw beam with a spatial dimension replaced by a temporal one. Saturating nonlinearity is required to arrest collapse and avoid the formation of a singularity. As in the 2D spatial case there is energy shedding on entry to the medium due to reshaping and shedding of energy during propagation due to the inclusion of saturation. The resulting beam width and power loss for varying γ are presented in Fig. 8. As the nonlinear coefficient is increased, more focused oscillations occur with larger and larger frequencies. Once the oscillations have damped down, our results show significantly reduced beam widths/lengths compared to the input. Figure 14 shows the resulting pulse sizes for β = −1, γ = 20 and γ = 50, with the achievement of a width as small as 1/10 the input size, again equivalent to the 2D spatial case. Although these results are very encouraging about possible regimes in the anomalous GVD regime where spatial and temporal compressions are observed in the presence of energy shedding, we note that there are regimes where for both Kerr and saturating nonlinearity pulse collapse is observed. Increasing the number of transverse grid points does not remedy the situation as the pulse is being compressed beyond the capabilities of the current simulation model and further terms in the equation would be necessary to model such ultra-short pulses. The case is somewhat similar to the inclusion of higher order dispersive terms in the generation of super-continuum radiation [10].

Conclusions
Dynamical energy shedding is an important phenomenon in the analysis of propagating beams in media with optical nonlinearities. We have shown that in lossless media, modelled by conservative Nonlinear Schrödinger Equations, there exist mechanisms that cause shedding of energy away from the self-trapping beam.
Reshaping of the profile of the initially Gaussian beam or pulse profile was observed in all cases, causing shedding of energy upon entry to the medium. In the 1D cw regime this is the only observed mechanism of energy shedding and accounts for up to 25% of the input energy being shed as a stable soliton is formed.
The inclusion of saturation in the model is required in order to simulate propagation of beams in two or more dimensions and introduced another mechanism of energy shedding. A long lived oscillation is induced and shedding of small amounts of energy occurs during the defocusing cycles of these oscillations as the beam relaxes towards a stable solution branch. We observed that with larger values of the saturation parameter, σ , a larger value of γ is required to get the same self-focusing effect, including increasing the critical power required for self-focusing. Unlike the results for media with a cubic nonlinearity, we also observe that higher-order solitons are unstable in the saturating regime and result in a symmetry breaking collapse of the beam.
For the case of 2+1D spatio-temporal pulses we have observed that in anomalous media, when the diffraction and dispersion lengths have equal magnitude, there is a symmetrical focusing of the pulse followed by a periodic alternation of defocusing and focusing that is very similar to the cw 2D spatial propagation in a saturable medium. We have shown that the selffocusing is increased with increasing γ, indicating that a tightly self-focused (2+1)D pulse exists and can stably propagate over long distance although it will be accompanied by greater energy shedding during the reshaping process. When anomalous dispersion and diffraction are significantly unequal, we have also observed asymmetrical self-focusing regimes that may require additional dispersive terms in the modelling.
For normal dispersion our results match those of simulations [15] that have been verified by experiments [6]. At low power or low nonlinearity the pulse diffracts and above a critical level the pulse splits until the power drops below critical level. These results confirm that the normal GVD regime is not favourable for self-focusing.