Spontaneously formed autofocusing caustics in a confined self-defocusing medium

Self-accelerating Airy beams, which are nondiffracting waves in the form of an Airy function that propagate in free space with constant acceleration, have received considerable attention in recent years. They are typically generated by manipulation of the phase front of the wave by means of specially designed optical elements. Here we show that autofocusing, radially symmetric Airy waves can form spontaneously as a laser beam propagates in a defocusing, nonlocal thermal nonlinear medium, inside a cylindrical channel with a reflective boundary. The beam forms a ringshaped optical caustic, which, following reflection from the boundary, converges to a focal point. We demonstrate this new method experimentally and numerically, and present a semi-classical analytical model for the wave dynamics that shows that the self-generated, radially symmetric wave is indeed a caustic with an Airy-function profile. In the hydrodynamic representation of the nonlinear wave equation, the ring-shaped caustic that we describe can be interpreted as a shock wave that forms as the “photonic fluid” bounces off the reflective boundary. These results suggest a very simple and accessible, yet mathematically accurate, way to obtain autofocusing radially symmetric Airy waves for various applications. © 2015 Optical Society of America

Self-accelerating Airy beams, namely nondiffracting waves in the form of an Airy function, with probability density functions (but not "centers of mass") that propagate in free space with constant acceleration, have received considerable attention since they were first proposed and demonstrated [1][2][3], and were shown to be useful for applications such as optical manipulation of small particles, microscopy, laser machining, and generation of light bullets and plasma channels in air [4][5][6][7][8].The idea has been generalized to nonlinear and nonlocal media [9][10][11] and to arbitrary convex trajectories [12], and was also applied to electron beams [13] and surface plasmons [14].The majority of works in this field have focused on one-dimensional Airy beams that propagate in a plane, but radially symmetric Airy waves, which exhibit abrupt autofocusing, have also been proposed and demonstrated [15,16].
The generation of Airy waves usually relies on manipulation of the phase front of the wave by means of a phase mask, diffractive optical elements, or wave-mixing processes.Here we describe spontaneously generated, ring-shaped autofocusing Airy waves that form when an optical caustic in a defocusing, nonlocal, thermal nonlinear medium is reflected from a cylindrical boundary.Rings in defocusing thermal media have been known for many years [17], although their relation to optical caustics was not explicitly mentioned.Here we use them to construct a ring-shaped caustic that autofocuses to a point.We first demonstrate the idea experimentally and numerically, and then study the wave dynamics analytically using a semi-classical model.We also briefly discuss the hydrodynamic representation of the nonlinear wave equation, interpreting the ring-shaped caustic as a shock wave that forms as the "photonic fluid" bounces off the reflective boundary.
In the paraxial approximation, the propagation of a coherent light beam through a nonlinear medium is usually described analytically by the nonlinear Schrödinger equation: Here ψ is the weakly z-dependent complex amplitude of the light propagating in the z direction, β 0 is the propagation constant, ∇ 2 ∂ 2 x ∂ 2 y is the transverse Laplacian, U x; y is the equivalent potential created by static spatial variations of the refractive index, and g is the nonlinear Kerr coefficient, representing the instantaneous nonlinear response of the medium.The last term in the equation thus describes the nonlinear refractive index changes due to local spatial variations of the light intensity.Strictly speaking, Eq. ( 1) is not applicable in the case of a thermal nonlinearity [17,18], where the light-induced refractive index variations are nonlocal in both space and time (due to heat diffusion and accumulation, respectively).It is common to take the spatial nonlocality into account by considering a nonlocal response kernel [19], but in fact, since in many cases a steady state is obtained after sufficiently long exposure to laser light (typically ∼1 s for laser intensities of a few watts), it is appropriate to describe the delayed, nonlocal response as a static potential U x; y, determined by the laser beam characteristics, the heat capacity and transport coefficients of the nonlinear medium, and the boundary conditions [20,21].
Consider a laser beam propagating in a thermal, self-defocusing medium [meaning that the index of refraction decreases as a function of the temperature, and corresponding to positive g in Eq. ( 1)] along the axis of a cylindrical channel with a reflective boundary [see Fig. 1(a)].The material that forms the boundary has very high thermal conductivity and heat capacity.In the steady state the heat flux produced by the laser is balanced by heat diffusion to the boundary, which acts as a heat sink.The heat diffusion is sustained by a radial temperature gradient, and the temperature profile that forms inside the channel, and peaks on its axis, is best described by a static repulsive potential with radial symmetry.As the laser beam propagates along the z direction, in the presence of the repulsive potential, its phase fronts acquire an increasing concave curvature, which corresponds to defocusing.(Note that the same result would be obtained in the case of a local, defocusing Kerr nonlinearity, but that the result here depends only weakly on the details of the light intensity distribution in the channel.)Now if the laser beam is sufficiently broad, reflection off the boundary forms another, inwardpropagating wave, with convex phase fronts, the curvature of which increases as a function of z.As shown in Fig. 1(b), this is a recipe for formation of a caustic.
To illustrate the above scenario we present one-dimensional numerical simulations of Eq. ( 1), for cylindrical symmetry, assuming a static potential that is a solution of the heat diffusion equation for uniform light intensity (i.e., heat generation) in the channel, namely U r −β 0 n 0 Δnr −β 0 n 0 Δn1 − r∕R 2 , where n 0 is the linear refractive index, r is the radial coordinate, and R is the radius of the channel.(Although this is not a selfconsistent model, it greatly simplifies the numerical simulation; moreover, we have verified that the results of the simulation are not very sensitive to the exact form of repulsive potential that we assume.) Figure 2 shows simulation results for input parameters that correspond to the experiment described later.Figure 2(a) is a contour plot of the light intensity as a function of the radial and z coordinates, and Fig. 2(b) is a plot of the intensity as a function of the radial coordinate at the exit of the 8-cm-long channel.The figure clearly shows the formation of a ring, which propagates inward as z increases.The cross section of the ring resembles an Airy function, which is known to be the typical mathematical form of a caustic [1].Below we will describe a semi-analytical model that shows that this conjecture is in fact accurate.At this point, however, we just show numerically that when propagating in free space, behind the nonlinear channel, the self-generated ring is also self-accelerating, and autofocuses to a point in space [Fig.2(c)].Finally, we find that the instantaneous Kerr nonlinearity [the last term on the right-hand side of Eq. ( 1)] does not play any role in the formation of the ring.
In the experiment we launch the beam of a frequency-doubled YAG laser (532 nm, continuous wave) into a 3-mm-diameter, 71-mm-long, circular channel, with polished, reflective walls, drilled in an aluminum block.The latter is enclosed in a cell with glass windows, filled with iodine-doped ethanol (∼20 ppm)-a defocusing nonlinear medium.The nonlinear index variations result from absorption by the iodine, which in turn leads to nonlocal, thermally induced changes of the index of refraction [17,18].Typically, after ∼1 s exposure to the laser (1-4 W) a steady state is reached, where the heat flux produced by the laser is balanced by heat diffusion to the aluminum block, which acts as a heat sink.Figure 3(a) shows an image of the beam at the output of the channel, obtained with a CCD camera, for a 1 W input beam with a ∼2.5 mm waist.A bright ring is clearly observed.Similar rings are obtained for a wide range of input powers   and beam waist sizes.As the input power is increased, the ring moves inward.Reducing the duty cycle of the laser, by means of a mechanical shutter, has the opposite effect.Using fluorescent dye on a glass cover slip as a screen allows imaging of the beam beyond the channel output.Figures 3(b) and 3(c) show images of the fluorescent screen when it is positioned at the output plane of the channel (b) and ∼6 cm behind it (c), with a 3.5 W input beam and a 27% duty cycle [about the same average power as in Fig. 3(a)], clearly showing that the outgoing beam focuses to a point.The focusing distance from the output plane can be controlled by varying the intensity and duty cycle of the laser.These findings are in excellent agreement with the simulations.
We now turn to the semi-classical model.Since the numerical simulations of Eq. ( 1) show that the instantaneous Kerr nonlinearity is not important, we replace Eq. ( 1) with a dimensionless Schrödinger equation for a quantum particle, with implicit cylindrical symmetry (i.e., only radial momentum): Here r is the normalized radial coordinate 0 ≤ r ≤ 1; t zU r0 , the mass m β 0 R 2 U r0 ∼ 10 5 corresponds to the values of the experimental parameters, and f r U r∕U r0 is a normalized potential satisfying f 0 1, f r > 1 ∞, as well as f 0 1 − 0, to simplify the behavior at the boundary.Furthermore, since the ring forms far from the center of the channel, the cylindrical term r −1 ∂ r in the Laplacian is neglected (we have numerically verified the excellent accuracy of this approximation for the beam configuration that we consider here).Introducing the notation ψr; t Ar; te iSr;t ; (3) where A and S are real, we obtain We solve the problem in the semi-classical limit in which A is a slowly varying function of position.The first term on the RHS of Eq. (4a) can then be neglected, and one obtains the Hamilton-Jacobi equation: while Eq. ( 4b) is a current-conservation equation.
To solve the problem semi-classically, we evolve a swarm of classical test particles in phase space fq; pg with the Hamiltonian H p 2 ∕2 m f r and the following initial conditions: qt 0 q 0 ∈ 0;1; pt 0 p 0 ∂S∕∂rj q 0 ;t0 0. The density of initial points along the segment [0, 1] is given by the initial distribution jψr; t 0j 2 A 2 r; t 0.
The phase-space curves that describe the distribution of classical particles as a function of time and that evolve according to the Hamiltonian flow are Lagrangian manifolds.At t 0 the Lagrangian manifold corresponding to our initial configuration coincides with the line p 0; q ∈ 0; 1.For sufficiently short times the position qt of a classical particle is a unique function q qq 0 ; t of its initial position q 0 , and one can determine q 0 for each r and t by inverting qq 0 ; t r.The standard semiclassical approach then yields ψr; t ∂q 0 ∂q 1∕2 A 0 q 0 r; te iSq 0 r;t;t ; where A 0 q 0 Aq 0 ; t 0 is the initial amplitude at the initial position of the particle, and Sq 0 r; t; t is the action R t 0 Lq; _ q; t 0 dt 0 calculated along the classical trajectory, which leads from q 0 to q. Figure 4 shows the evolution of the classical test particles in phase space and the corresponding density profiles A 2 r; t calculated from Eq. ( 6).
The situation gets more complex after some of the classical test particles bounce from the hard wall at r 1 and propagate backward toward r 0, against the potential f r, which typically occurs for t ∼ 250 (corresponding to a propagation of 5 cm in the experiment).This induces a fold catastrophe in the Lagrangian manifold.As seen in Fig. 5, there are now either one or three momenta p associated with each position q, which we denote by p 1 , p 2 , and p 3 .The latter corresponds to a particle (i.e., wave component) that still propagates outward, while p 1 and p 2 are the momenta of particles propagating inward.As we show below, using standard mathematical analysis, the joint contribution of these particles forms an interference peak, i.e., a caustic in the density profile at the turning point r R c (to keep the notation simple we do not write the time dependence of R c explicitly).
Since each r > R c can be reached by three different classical trajectories [q 1 t, q 2 t, and q 3 t in Fig. 5(b)], a trivial extension of the semi-classical approach leading to Eq. ( 6) yields  6); and the black dots are brute-force numerical solutions of Eq. ( 2); without loss of generality we take f r 1 cos πr∕2.Fig. 5. (a) Classical evolution of the Lagrangian manifold in phase space for t 400, corresponding to a propagation distance of 8 cm in the experiment-a caustic appears at r R c t.(b) The three classical trajectories q i t correspond to the three momenta p i in (a); q 1 t and q 2 t merge exactly for r R c t. ψr; t X b∈f1;2;3g ψ b r; t X b∈f1;2;3g e iμ b π∕2 A 0 q b r; t ∂q b 0 ∂q 1∕2 e iSq b r;t ; (7) where q b 0 q b t 0, the μ b 's are Maslov indices [μ 3 0 as in Eq. ( 6), and μ 1 ; μ 2 2; 1], and ψ 1 and ψ 2 are zero for r < R c .The intensity profile obtained from this crude semi-classical approximation is depicted in Fig. 6(a), for the same parameters and time as in Fig. 5. Also shown is the intensity profile calculated by numerical integration of Eq. ( 2). Figure 6(a) shows that near the caustic (r R c ) approximation Eq. ( 7) diverges-the two classical trajectories q 1 t and q 2 t coincide, and the stationary phase method, which is at the heart of the semi-classical approximation, must be refined, since the factors ∂q b 0 ∕∂q in Eq. ( 7) diverge for b 1; 2. Yet, since there are no caustics in momentum space, one may obtain a better solution in the vicinity of the caustic by Fourier-transforming the wave function at an earlier time, before the caustic forms, evolving the Fourier components in time, and finally Fouriertransforming back to real space.The last step must be performed with care because of the presence of two nearby saddle points in the integral (corresponding to the classical momenta p 1 and p 2 ).The result for ψ 1 ψ 2 involves the Airy function Ai [the contribution of ψ 3 in Eq. ( 7) is not modified]: where q q 1 0 R c ; t q 2 0 R c ; t is the initial position leading to q 1;2 t R c , and the quantities l and S c are defined by l 3  1 2 ; S c S q; t: This result, combining an isolated classical path and two others forming a caustic, is identical to rainbow scattering from a Lennard-Jones-type central potential [22].The corresponding intensity profile depicted in Fig. 6(b) validates the interpretation of the experimental and numerical results in terms of a caustic, with a typical Airy-function form.Note, however, that the divergence of Eq. ( 7) at r R c has been eliminated in Eq. ( 8) at the price of loss of accuracy away from the caustic, in the region 0.8 ≲ r < 1 [cf.Fig. 6(b)].Yet the discrepancy can be resolved by more uniform treatment, following [23].
The final result is A 0 q 2 0 r;t − A 0 q 1 0 r;t ∂q Figure 7(a) compares the results of this uniform approximation with a numerical solution of Eq. ( 2). Figure 7(b) is a sketch of the semi-classical trajectories and the position of the (timedependent) caustic obtained in this approximation.Considering all the approximations made, and in particular the impossibility of measuring the exact shape of the static potential, it is quite remarkable that the position of the ring in the experimental data is in excellent agreement with the prediction of the semi-classical model, for similar laser powers and propagation lengths.
In the hydrodynamic representation of the nonlinear wave equation, the Madelung transformation [24] maps Eq. ( 1) to coupled equations for the amplitude and phase, which have the form of continuity and Euler equations for an equivalent "photon fluid."This approach has been used successfully to model dispersive shock waves that appear around laser beams that propagate in isotropic self-defocusing media [25].Similarly, when the nonlinear term in Eq. ( 1) is predominant, the caustics can be interpreted as shock waves in the photon fluid.The increasing curvature of the reflected wave fronts as a function of z [see Fig. 1(b)] corresponds to increasing velocity of the partial waves as a function of time, leading to steepening of the wave and formation of a shock.We checked numerically and qualitatively that the caustic associated with the fold singularity of the Lagrangian manifold indeed turns into a dispersive shock wave when the nonlinear term is sufficiently large.
Finally, the simple method that we suggest for generating autofocusing, radially symmetric Airy waves has several advantages compared to established methods.Unlike a typical spatial light modulator, it can handle high average powers.Moreover, as the autofocusing distance is power-dependent, the autofocusing element can be combined with a pinhole in a nonlinear switching device, the pinhole position defining the threshold power above which transmission is high.The abrupt autofocusing of the non-Gaussian beam will result in a sharp transition.

Fig. 1 .
Fig. 1.(a) Setup described in the text.(b) Sketch explaining the formation of a caustic from partial waves reflected off the boundary.

Fig. 2 .
Fig. 2. (a) Numerical simulation of Eq. (1) with a 3-mm-diameter, 8-cm-long channel, a centered, 1-mm-diameter Gaussian input beam, the parabolic potential described in the text, U 0 9000 m −1 , β 0 1.181 × 10 7 m −1 , and gjA0j 2 100 m −1 ; (a) is a contour plot of the beam intensity as a function of the radial coordinate and the propagation distance inside the channel; (b) is the intensity profile at the output of the channel; (c) self-acceleration and focusing to a point of the output beam as it propagates in free space behind the channel.See also Supplement 1.

Fig. 3 .
Fig. 3. Experimental: (a) image of the output plane of the 3 mm channel for a 1 W input beam (note: the cylindrical channel has an opening, i.e., a groove on one side, which results in a cutoff piece that is seen ejected to the left); (b) image of the fluorescence at the output plane of the channel for a 3.5 W input beam and a 27% duty cycle; (c) the ring in (b) focuses to a point after propagating ∼6 cm in free space.

Fig. 4 .
Fig. 4. (a) Classical evolution of the Lagrangian manifold in phase space for short times.(b)Corresponding density profiles-the red dashed line corresponds to the initial density, which is arbitrarily taken to be A 0 r 1 − r 2 ; the purple solid line is the semi-classical result Eq. (6); and the black dots are brute-force numerical solutions of Eq. (2); without loss of generality we take f r 1 cos πr∕2.

Fig
Fig. S3.A numerical simulation of beam propagation through and behind the channel, with the same parameters as in the simulation shown in Figure S1, except for substitution of the semi-circular static potential with a Gaussian form of similar width and with the same maximum value (at the center of the channel); (a) is a contour plot of the beam intensity as a function of the radial coordinate and the propagation distance inside the channel; (b) is the intensity profile at the output of the channel; (c) shows selfacceleration and focusing to a point of the output beam as it propagates in free space behind the channel.

Fig
Fig. S4.A sketch of some semi-classical trajectories obtained from Eq. (10) in the article, plotted as orange lines, and the (time dependent) location of the caustic, plotted as a red line; the curves show the evolution of the amplitude of the wave.