Shocks in nonlocal media
Neda Ghofraniha,1 Claudio Conti,2,3 Giancarlo Ruocco,3,4 Stefano Trillo5∗
1 Research Center SMC INFM-CNR, Università di Roma “La Sapienza”, P. A. Moro 2, 00185, Roma, Italy
2Centro Studi e Ricerche “Enrico Fermi”, Via Panisperna 89/A, 00184 Rome, Italy
3Research center SOFT INFM-CNR Università di Roma “La Sapienza”, P. A. Moro 2, 00185, Roma, Italy
4Dipartimento di Fisica, Università di Roma “La Sapienza”, P. A. Moro 2, 00185, Roma, Italy
5 Dipartimento di Ingegneria, Università di Ferrara, Via Saragat 1, 44100 Ferrara, Italy
(Dated: October 1, 2018)
We investigate the formation of collisionless shocks along the spatial profile of a gaussian laser
beam propagating in nonlocal nonlinear media. For defocusing nonlinearity the shock survives
the smoothing effect of the nonlocal response, though its dynamics is qualitatively affected by the
latter, whereas for focusing nonlinearity it dominates over filamentation. The patterns observed in
a thermal defocusing medium are interpreted in the framework of our theory.
Shock waves are a general phenomenon thoroughly in-
vestigated in disparate area of physics (fluids and water
waves, plasma physics, gas dynamics, sound propagation,
physics of explosions, etc.), entailing the propagation of
discontinuous solutions typical of hyperbolic PDE mod-
els [1, 2]. They are also expected in (non-hyperbolic)
universal models for dispersive nonlinear media, such as
the Korteweg-De Vries (KdV) and nonlinear Schrödinger
(NLS, or analogous Gross-Pitaevskii) equations, since hy-
drodynamical approximations of such models hold true
in certain regimes (typically, in the weakly dispersive or
strongly nonlinear case) [3, 4, 5]. However, in the lat-
ter cases, no true discontinuous solutions are permitted.
The general scenario, first investigated by Gurevich and
Pitaevskii [3], is that dispersion regularizes the shock, de-
termining the onset of oscillations that appear near wave-
breaking points and expand afterwards. This so-called
collisionless shock has been observed for example in ion-
acoustic waves [6], or wave-breaking of optical pulses in
a normally dispersive fiber [7], and recently in a Bose-
Einstein condensate with positive scattering length [8].
In this Letter we investigate how nonlocality of the
nonlinear response affects the formation of a collisionless
shock in a system ruled by a NLS model. In fact nonlocal-
ity plays a key role in many physical systems due to trans-
port phenomena and finite range interactions (e.g. as in
Bose-Einstein condensation), and can be naively thought
to smooth and eventually wipe out steep fronts character-
istic of shocks. More specifically, we place this problem
in the context of nonlinear optics where nonlocality arises
quite naturally in different media [9, 10, 11, 12], study-
ing the spatial propagation of a fundamental (gaussian
TEM00) laser mode subject to diffraction and nonlocal
focusing/defocusing action (Kerr effect). In a defocus-
ing and ideal (local and lossless) medium, high intensity
portions of the beam diffract more rapidly than the tails
leading, at sufficiently high powers, to overtaking and os-
cillatory wave-breaking similar (in 1D) to what observed
in the temporal case [18]. We find that, while shock is
not hampered by the presence of (even strong) nonlocal-
ity, the mechanism of its formation as well as post-shock
patterns are qualitatively affected by the nonlocality. Ex-
perimental results obtained with a thermal defocusing
nonlinearity are consistent with our theory and shed new
light on the interpretation of the thermal lensing phe-
nomenon.
Importantly, our theory permits also to establish that
nonlocality allows the shock to form also in the focusing
regime where, contrary to the local case, it can prevails
over filamentation or modulational instability (MI).
Theory We start from the paraxial wave equation
obeyed by the envelope A of a monochromatic field
E = ( 2
)1/2A exp(ikZ − iωT ) (|A|2 is the intensity)
+ k0∆nA = −i
A. (1)
where k = k0n =
n is the wave-number, and α0
the intensity loss rate. A sufficiently general nonlocal
model can be obtained by coupling Eq. (1) to an equa-
tion that rules the refractive index change ∆n of nonlin-
ear origin. Introducing the scaled coordinates x, y, z =
X/w0, Y/w0, Z/L, and complex variables ψ = A/
and θ = k0Lnl∆n, where Lnl = (k0|n2|I0)−1 is the non-
linear length scale associated with peak intensity I0 and
a local Kerr coefficient n2 (∆n = n2|A|2), Ld = kw20 is
the characteristic diffraction length associated with the
input spot-size w0, and L ≡
LnlLd, such model can be
conveniently written as follows [12]
ψ + χθψ = −iα
εψ, (2)
−σ2∇2
θ + θ = |ψ|2, (3)
where α = α0L, ∇2⊥ = ∂2x + ∂2y , χ = n2/|n2| = ±1 is
the sign of the nonlinearity, and σ2 is a free parameter
that measures the degree of nonlocality. The peculiar
dimensionless form of Eqs. (2-3) where ε ≡ Lnl/L =
Lnl/Ld is a small quantity, highlights the fact that we
will deal with the weakly diffracting (or strongly non-
linear) regime, such that the local σ = 0 and lossless
α = 0 limit yields a semiclassical Schrödinger equation
with cubic potential (ε and z replace Planck constant
http://arxiv.org/abs/0704.0610v1
and time, respectively). We study Eqs. (2-3) subject
to the axi-symmetric gaussian input ψ0(r) = exp(−r2),
x2 + y2, describing a fundamental laser mode
at its waist. For ε ≪ 1, its evolution can be studied
in the framework of the WKB trasformation ψ(r, z) =
ρ(r, z) exp [iφ(r, z)/ ε] [4]. Substituting in Eqs. (2-3)
and retaining only leading orders in ε, we obtain
(D − 1)
ρu+ (ρu)r
= −αρ; uz + uur − χθr = 0,
θrr +
D − 1
+ θ = ρ. (4)
where u ≡ φr is the phase chirp, and D = 2 is the trans-
verse dimensionality. The 1D case described by Eqs. (4)
with D = 1 and r → x (∂y = 0) illustrates the ba-
sic physics with least complexity. In the defocusing case
(χ = −1) for an ideal medium (σ = α = 0, θ = ρ),
Eqs. (4) are a well known hyperbolic system of conser-
vation laws (Eulero and continuity equations) with real
celerities (or eigenspeeds, i.e. velocities of Riemann in-
variants) v± = u ±
−χρ, which rules gas dynamics (u
and ρ are velocity and mass density of a gas with pres-
sure ∝ ρ2). A gaussian input is known to develop two
symmetric shocks at finite z [4]. Importantly the diffrac-
tion, which is initially of order ε2, starts to play a major
role in the proximity of the overtaking point, and regu-
larize the wave-breaking through the appearance of fast
(wavelength ∼ ε) oscillations which connect the high and
low sides of the front and expand outwards (far from the
beam center) [3]. Such oscillations, characteristic of a col-
lisionless shock, appear simultaneously in intensity and
phase chirp (u) as clearly shown in Fig. 1(a,c).
In the nonlocal case, the index change θ(x) can be
wider than the gaussian mode (for large σ) and the shock
dynamics is essentially driven by the chirp u with ρ adia-
batically following. This can be seen by means of the fol-
lowing approximate solution of Eqs. (4): considering that
the equation for ρ is of lesser order [O(ǫ)], with respect
to those for θ and u [O(1)], we assume ρ = exp(−2x2)
unchanged in z and solve exactly the third of Eq. (4) for
θ(x) (though derived easily, its full expression is quite
cumbersome). Then, applying the theory of characteris-
tics [1], the second of Eqs. (4) is reduced to the following
ODEs, where dot stands for d/dz
ẋ = u ; u̇ = χθx, (5)
equivalent to the motion of a unit mass in the potential
V (x) = −χθ with conserved energy E = u(z)
+ V (x).
The solution of Eqs. (5) with initial condition x(0) =
s, u(0) = 0 yields x(s, z), u(s, z) in parametric form, from
which overtaking is found whenever u(x, z) (obtained
by eliminating s) becomes a multivalued function of x
at finite z = zs. The shock point corresponding to
|du/dx| → ∞ is found from the solution u(x, z) displayed
in Fig. 2(a) [ 2(b)], at positions x = ±xs 6= 0 (defocusing
FIG. 1: (Color online) 1D spatial profiles of phase chirp u(x)
(a-b) and amplitude |ψ(x)| =
ρ(x) (c-d), as obtained from
Eqs. (2-3) with ε = 10−3, α = 0, χ = −1 (defocusing), and
ψ0 = exp(−x
2), for different z as indicated: (a-c) local case,
σ2 = 0; (b-d) nonlocal case, σ2 = 5.
FIG. 2: (Color online) (a) u(x) for different z and χ = 1
(focusing), σ = 1; (b) as in (a) for χ = −1 (defocusing); (c)
shock distance zs (χ = −1 bold solid, χ = 1 thin solid) and
shock position xs in the defocusing case (dashed line).
case) or xs = 0 (focusing case). The shock distance zs
increases with σ in both cases, as shown in Fig. 2(c).
We have tested these predictions by integrating nu-
merically Eqs. (2-3). Simulations with χ = −1 [see Fig.
1(b,d)] show indeed steepening and post-shock oscilla-
tions in the spatial chirp u, which are accompanied by
a steep front in ρ moving outward. The shock location
in x and z is in good agreement with the results of our
approximate analysis summarized in Fig. 2.
Numerical simulations of Eqs. (2-3) validates also
the focusing scenario. The field evolution displayed in
Fig. 3(a) exhibits shock formation at the focus point
(xs = 0, zs ≃ 8, for σ = 5) driven the phase whose chirp
is shown in Fig. 3(b). This is remarkable because, in
the local limit σ = 0, the celerities become imaginary
(the equivalent gas would have pressure decreasing with
increasing density ρ), and no shock could be claimed to
exist. In this limit, the reduced problem (4) is elliptic and
the initial value problem is ill-posed [13], an ultimate con-
sequence of the onset of MI: modes with transverse (nor-
FIG. 3: (Color online) Level plot of the intensity in the
focusing case (χ = 1, ε = 0.01): (a) nonlocal case (σ2 =
25); (b) chirp profile for various z for (a); (c) quasi-local case
(σ2 = 10−5).
malized) wavenumber q < ∆q grow exponentially with z,
with both gain g and bandwidth ∆q scaling as 1/ε. How-
ever, the nonlocal response tends to frustrate MI (see also
Refs. [9, 12]), as shown by standard linear stability anal-
ysis which yields g =
d(2χ− d)/ε2 (we set d ≡ ε2q2/2
and χ ≡ χ/(1 + σ2q2)), in turn implying a strong reduc-
tion of both gain and bandwidth for large σ. In order
to emphasize the difference between the local and non-
local regime, we contrast Fig. 3(a) with the analogous
evolution [see Fig. 3(c)] obtained in the quasi-local limit
(σ2 = 10−5), which appears to be clearly dominated by
filamentation.
Thermal nonlinearity The physics of the defocusing
case can be experimentally tested by exploiting ther-
mal nonlinearities of strongly absorptive bulk samples,
that we show below to fit our model. In this case, the
system relaxes to a steady-state refractive index change
∆n = (dn/dT )∆T , where dn/dT is the thermal coef-
ficient, and ∆T the local temperature change due to
optical absorption. It is well known that this so-called
thermal lens distorts a laser beam propagating in the
medium [14, 15, 16]. However, only perturbative ap-
proaches to the problem have been proposed (ray optics
or Fresnel diffraction theory is applied after the lens pro-
file is worked out from gaussian ansatz [14]), while the
role of shock phenomena was completely overlooked.
FIG. 4: 2D evolution according to Eqs. (2-3) with σ2 = 1,
α = 1: (a) radial phase chirp at different z, as indicated,
showing steepening and shock formation for ε = 10−2; (b) cor-
responding intensity profile |ψ|2 (maximum scaled to unity)
at z = 4.9; (c) transverse intensity profile vs. x (at y = 0) at
z = 1/(4ε) and different values of ε (α0 = 62cm
−1, σ = 0.3).
We assume that the temperature field ∆T =
∆T⊥(X,Y ) obeys the following 2D heat equation
(∂2X + ∂
Y )∆T⊥ − C∆T⊥ = −γ|A|2 (6)
where the source term account for absorption pro-
portional to intensity through the coefficient γ =
α0/(ρ0cpDT ), where ρ0 the material density, cp the spe-
cific heat at constant pressure, and DT is the thermal
diffusivity (see e.g. [16]). Eq. (6) has been already em-
ployed to model a refractive index of thermal origin in
Ref. [10], and in Ref. [11] in the limit C = 0 which
is equivalent to consider the range of nonlocality (mea-
sured by 1/C, see below) to be infinite. Starting from
the 3D heat equation ∇2∆T = −γ|A|2, the latter regime
amounts to assume ∆T (X,Y, Z) = ∆T⊥(X,Y ), which is
justified when longitudinal changes in intensity |A|2 are
negligible as for solitary (invariant in Z) wave-packets
in the presence of low absorption [11]. Viceversa, in the
regime of strong absorption, we need to account for lon-
gitudinal temperature profiles that are known from solu-
tions of the 3D heat equations to be peaked at charac-
teristic distance Ẑ in the middle of sample and decay to
room temperature on the facets [14]. Since highly non-
linear phenomena occurs in the neighborhood of Ẑ where
the index change is maximum, we can use a (longitudi-
nal) parabolic approximation with characteristic width
Leff (∼ L) of the 3D temperature field ∆T (X,Y, Z) =
1− (Z−Ẑ)
∆T⊥(X,Y ) and consequently approximate
∇2∆T ≃ (∂2X + ∂2Y )∆T⊥ − L
eff∆T⊥, so that the 3D
heat equation reduces to Eq. (6) with C = 1/L2eff .
Following this approach, Eq. (6) coupled to Eq. (1)
can be casted in the form of Eqs. (2-3) by posing
θ = k0Lnl|dn/dT |∆T⊥ and σ2 = 1/(Cw20) = L2eff/w20.
The model reproduces the infinite range nonlocality for
negligible losses (Leff → ∞); while for thin samples
[|(∂2X + ∂2Y )∆T⊥| << |L
eff∆T⊥|], Leff can be related
to the Kerr coefficient n2 as
Leff =
γ|dn/dT |
DTρ0cp|n2|
α0|dn/dT |
which establishes a link between the degree of nonlocal-
ity and the strength of the nonlinear response (similarly
to other nonlocal materials [12]).
Having retrieved the model Eqs. (2-3), let us show
next that the scenario illustrated previously applies sub-
stantially unchanged in bulk (2D case) even on account
for the optical power loss (α 6= 0). An example of the
general dynamics is shown in Fig. 4, where we report a
simulation of the full model (2-3), with σ2 = 1 and rela-
tively large loss α = 1. In analogy to the 1D case, Fig.
4(a) clearly shows that the radial chirp u = φr steep-
ens and then develop characteristic oscillations after the
shock point (z ≃ 6, where |∂ru| → ∞). Correspond-
ingly the intensity exhibits also an external front which
is connected to a flat central region with a characteris-
tic overshoot [see Fig. 4(b)] corresponding to a brighter
ring [inset in Fig. 4(c)]. For larger distances this struc-
ture moves outward following the motion of the shock.
In the experiment such motion can be observed, at fixed
physical lenght, by increasing the power, which amounts
to decrease ε while scaling z and α accordingly (z ∝ 1/ε,
α ∝ ε), as displayed in Fig. 4(c) for σ = 0.3.
As a sample of a strongly absorbing medium we choose
a 1 mm long cell filled with an acqueous solution of
Rhodamine B (0.6 mM concentration). Our measure-
ments of the linear and nonlinear properties of the sam-
ple performed by means of the Z-scan technique gives
data consistent with the literature [17], and allows us
to extrapolate at the operating vacuum wavelength of
532 nm, a linear index n = 1.3, a defocusing nonlin-
ear index n2 = 7 × 10−7 cm2W−1, and α0 = 62 cm−1.
For our sample DT = 1.5 × 10−7 m2s−1, ρ0 = 103 kg
m−3, cp = 4 × 103Jkg−1K−1 and |dn/dT | = 10−4 K−1
(γ ∼= 104 K W−1), and exploiting Eq. (7) we estimate
Leff ∼= 10µm (Leff << L because of the strong ab-
sorption that causes strong heating of our sample near
the input facet), and correspondingly the degree of non-
locality σ ∼= 0.3. We operate with an input gaussian
beam with fixed intensity waist w0I = w0/
2 = 20 µm
(Ld ∼= 12 mm) focused onto the input face of the cell.
With these numbers, an input power P = πw20II0 = 200
mW yields a nonlinear length Lnl ∼= 8 µm (L ∼= 0.3 mm),
which allows us to work in the semiclassical regime with
ε ∼= 0.025. The radial intensity profiles together with the
2D patterns imaged by means of a 40×microscope objec-
tive and a recording CCD camera are reported in Fig. 5.
As shown the beam exhibits the formation of the bright
ring whose external front moves outward with increasing
power, consistently with the reported simulations. We
point out that, at higher powers, we observe (both ex-
perimentally and numerically) that the moving intensity
front leaves behind damped oscillations that correspond
to inner rings of lesser brightness, as reported in litera-
ture [15]. This, however, occurs well beyond the shock
point that we have characterized so far.
In summary, the evolution of a gaussian beam in the
strong nonlinear regime is characterized by occurence of
collisionless (i.e., regularized by diffraction) shocks that
survive the smoothing effect of (even strong) nonlocality.
While experimental results support the theoretical sce-
nario in the defocusing case, the remarkable result that
the nonlocality favours shock dynamics over filamenta-
tion requires future investigation.
∗ Electronic address: claudio.conti@phys.uniroma1.it
[1] G. B. Whitman, Linear and Nonlinear Waves (Wiley,
New York, 1974);
[2] L. D. Landau and E. M. Lifshitz, Fluid Mechanics (Perg-
amon, 1995); M. A. Liberman and A. L. Yelikovich,
Physics of Shock Waves in Gases and Plasmas (Springer,
FIG. 5: Radial profiles of intensity observed in the thermal
medium for different input powers. The insets show the cor-
responding 2D output patterns.
Heidelberg, 1986).
[3] A.V. Gurevich and L.P. Pitaevskii, Sov. Phys. JETP 38,
291 (1973); A.V. Gurevich and A. L. Krylov, Sov. Phys.
JETP 65, 944 (1987).
[4] J. C. Bronski and D. W. McLaughlin, in Singular Limits
of Dispersive Waves, NATO ASI Series, Ser. B 320, pp.
21-28 (1994); M. G. Forest and K. T. R. McLaughlin,
J. Nonlinear Science 7, 43 (1998); Y. Kodama, SIAM J.
Appl. Math. 59, 2162 (1999). M. G. Forest, J. N. Kutz,
and K. T. R. McLaughlin, J. Opt. Soc. Am. B 16, 1856
(1999).
[5] A. M. Kamchatnov, R. A. Kraenkel, and B. A. Umarov,
Phys. Rev. E 66, 036609 (2002).
[6] R. J. Taylor, D.R. Baker, and H. Ikezi, Phys. Rev. Lett.
24, 206 (1970).
[7] J. E. Rothenberg and D. Grischkowsky, Phys. Rev. Lett.
62, 531 (1989).
[8] M. A. Hoefer, M. J. Ablowitz, I. Coddington, E. A.
Cornell, P. Engels, and V. Schweikhard, Phys. Rev. A
74, 023623 (2006). V. M. Perez-Garcia, V.V. Konotop,
V.A. Brazhnyi, Phys. Rev. Lett. 92, 220403 (2004); B.
Damski, Phys. Rev. A 69, 043610 (2004);
[9] J. Wyller, W. Krolikowski, O. Bang, J. J. Rasmussen,
Phys. Rev. E 66, 066615 (2002).
[10] A. Yakimenko, Y. Zaliznyak and Y.S. Kivshar, Phys.
Rev. E 71, 065603(R) (2005).
[11] C. Rotschild, O. Cohen, O. Manela, M. Segev and T.
Carmon, Phys. Rev. Lett. 95, 213904 (2005).
[12] C. Conti, M. Peccianti and G. Assanto, Phys. Rev. Lett.
91 073901 (2003); Phys. Rev. Lett. 92 113902 (2004);
C. Conti, G. Ruocco and S. Trillo, Phys. Rev. Lett. 95
183902 (2005).
[13] P. D. Miller and S. Kamvissis, Phys. Lett. A 247, 75
(1998); J. C. Bronski, Physica D 152, 163 (2001).
[14] C. A. Carter and J. M. Harris, Appl. Opt. 23, 476 (1984);
S. Wu and N. J. Dovichi, J. Appl. Phys. 67, 1170 (1990);
F. Jürgensen and W. Schröer, Appl. Opt. 34 41 (1995).
[15] C. J. Wetterer, L. P. Schelonka, and M. A. Kramer, Opt.
Lett. 14, 874 (1989).
mailto:claudio.conti@phys.uniroma1.it
[16] P. Brochard, V. Grolier-Mazza and R. Cabanel, J. Opt.
Soc. Am. B 14, 405 (1997).
[17] S. Sinha, A. Ray, and K. Dasgupta, J. Appl. Phys. 87,
3222 (2000).
[18] paraxial diffraction in defocusing media is well known
to be isomorphus in 1D to propagation in a normally
dispersive focusing medium as considered in Ref. [7]

  We investigate the formation of collisionless shocks along the spatial
profile of a gaussian laser beam propagating in nonlocal nonlinear media. For
defocusing nonlinearity the shock survives the smoothing effect of the nonlocal
response, though its dynamics is qualitatively affected by the latter, whereas
for focusing nonlinearity it dominates over filamentation. The patterns
observed in a thermal defocusing medium are interpreted in the framework of our
theory.

Shocks in nonlocal media
Neda Ghofraniha,1 Claudio Conti,2,3 Giancarlo Ruocco,3,4 Stefano Trillo5∗
1 Research Center SMC INFM-CNR, Università di Roma “La Sapienza”, P. A. Moro 2, 00185, Roma, Italy
2Centro Studi e Ricerche “Enrico Fermi”, Via Panisperna 89/A, 00184 Rome, Italy
3Research center SOFT INFM-CNR Università di Roma “La Sapienza”, P. A. Moro 2, 00185, Roma, Italy
4Dipartimento di Fisica, Università di Roma “La Sapienza”, P. A. Moro 2, 00185, Roma, Italy
5 Dipartimento di Ingegneria, Università di Ferrara, Via Saragat 1, 44100 Ferrara, Italy
(Dated: October 1, 2018)
We investigate the formation of collisionless shocks along the spatial profile of a gaussian laser
beam propagating in nonlocal nonlinear media. For defocusing nonlinearity the shock survives
the smoothing effect of the nonlocal response, though its dynamics is qualitatively affected by the
latter, whereas for focusing nonlinearity it dominates over filamentation. The patterns observed in
a thermal defocusing medium are interpreted in the framework of our theory.
Shock waves are a general phenomenon thoroughly in-
vestigated in disparate area of physics (fluids and water
waves, plasma physics, gas dynamics, sound propagation,
physics of explosions, etc.), entailing the propagation of
discontinuous solutions typical of hyperbolic PDE mod-
els [1, 2]. They are also expected in (non-hyperbolic)
universal models for dispersive nonlinear media, such as
the Korteweg-De Vries (KdV) and nonlinear Schrödinger
(NLS, or analogous Gross-Pitaevskii) equations, since hy-
drodynamical approximations of such models hold true
in certain regimes (typically, in the weakly dispersive or
strongly nonlinear case) [3, 4, 5]. However, in the lat-
ter cases, no true discontinuous solutions are permitted.
The general scenario, first investigated by Gurevich and
Pitaevskii [3], is that dispersion regularizes the shock, de-
termining the onset of oscillations that appear near wave-
breaking points and expand afterwards. This so-called
collisionless shock has been observed for example in ion-
acoustic waves [6], or wave-breaking of optical pulses in
a normally dispersive fiber [7], and recently in a Bose-
Einstein condensate with positive scattering length [8].
In this Letter we investigate how nonlocality of the
nonlinear response affects the formation of a collisionless
shock in a system ruled by a NLS model. In fact nonlocal-
ity plays a key role in many physical systems due to trans-
port phenomena and finite range interactions (e.g. as in
Bose-Einstein condensation), and can be naively thought
to smooth and eventually wipe out steep fronts character-
istic of shocks. More specifically, we place this problem
in the context of nonlinear optics where nonlocality arises
quite naturally in different media [9, 10, 11, 12], study-
ing the spatial propagation of a fundamental (gaussian
TEM00) laser mode subject to diffraction and nonlocal
focusing/defocusing action (Kerr effect). In a defocus-
ing and ideal (local and lossless) medium, high intensity
portions of the beam diffract more rapidly than the tails
leading, at sufficiently high powers, to overtaking and os-
cillatory wave-breaking similar (in 1D) to what observed
in the temporal case [18]. We find that, while shock is
not hampered by the presence of (even strong) nonlocal-
ity, the mechanism of its formation as well as post-shock
patterns are qualitatively affected by the nonlocality. Ex-
perimental results obtained with a thermal defocusing
nonlinearity are consistent with our theory and shed new
light on the interpretation of the thermal lensing phe-
nomenon.
Importantly, our theory permits also to establish that
nonlocality allows the shock to form also in the focusing
regime where, contrary to the local case, it can prevails
over filamentation or modulational instability (MI).
Theory We start from the paraxial wave equation
obeyed by the envelope A of a monochromatic field
E = ( 2
)1/2A exp(ikZ − iωT ) (|A|2 is the intensity)
+ k0∆nA = −i
A. (1)
where k = k0n =
n is the wave-number, and α0
the intensity loss rate. A sufficiently general nonlocal
model can be obtained by coupling Eq. (1) to an equa-
tion that rules the refractive index change ∆n of nonlin-
ear origin. Introducing the scaled coordinates x, y, z =
X/w0, Y/w0, Z/L, and complex variables ψ = A/
and θ = k0Lnl∆n, where Lnl = (k0|n2|I0)−1 is the non-
linear length scale associated with peak intensity I0 and
a local Kerr coefficient n2 (∆n = n2|A|2), Ld = kw20 is
the characteristic diffraction length associated with the
input spot-size w0, and L ≡
LnlLd, such model can be
conveniently written as follows [12]
ψ + χθψ = −iα
εψ, (2)
−σ2∇2
θ + θ = |ψ|2, (3)
where α = α0L, ∇2⊥ = ∂2x + ∂2y , χ = n2/|n2| = ±1 is
the sign of the nonlinearity, and σ2 is a free parameter
that measures the degree of nonlocality. The peculiar
dimensionless form of Eqs. (2-3) where ε ≡ Lnl/L =
Lnl/Ld is a small quantity, highlights the fact that we
will deal with the weakly diffracting (or strongly non-
linear) regime, such that the local σ = 0 and lossless
α = 0 limit yields a semiclassical Schrödinger equation
with cubic potential (ε and z replace Planck constant
http://arxiv.org/abs/0704.0610v1
and time, respectively). We study Eqs. (2-3) subject
to the axi-symmetric gaussian input ψ0(r) = exp(−r2),
x2 + y2, describing a fundamental laser mode
at its waist. For ε ≪ 1, its evolution can be studied
in the framework of the WKB trasformation ψ(r, z) =
ρ(r, z) exp [iφ(r, z)/ ε] [4]. Substituting in Eqs. (2-3)
and retaining only leading orders in ε, we obtain
(D − 1)
ρu+ (ρu)r
= −αρ; uz + uur − χθr = 0,
θrr +
D − 1
+ θ = ρ. (4)
where u ≡ φr is the phase chirp, and D = 2 is the trans-
verse dimensionality. The 1D case described by Eqs. (4)
with D = 1 and r → x (∂y = 0) illustrates the ba-
sic physics with least complexity. In the defocusing case
(χ = −1) for an ideal medium (σ = α = 0, θ = ρ),
Eqs. (4) are a well known hyperbolic system of conser-
vation laws (Eulero and continuity equations) with real
celerities (or eigenspeeds, i.e. velocities of Riemann in-
variants) v± = u ±
−χρ, which rules gas dynamics (u
and ρ are velocity and mass density of a gas with pres-
sure ∝ ρ2). A gaussian input is known to develop two
symmetric shocks at finite z [4]. Importantly the diffrac-
tion, which is initially of order ε2, starts to play a major
role in the proximity of the overtaking point, and regu-
larize the wave-breaking through the appearance of fast
(wavelength ∼ ε) oscillations which connect the high and
low sides of the front and expand outwards (far from the
beam center) [3]. Such oscillations, characteristic of a col-
lisionless shock, appear simultaneously in intensity and
phase chirp (u) as clearly shown in Fig. 1(a,c).
In the nonlocal case, the index change θ(x) can be
wider than the gaussian mode (for large σ) and the shock
dynamics is essentially driven by the chirp u with ρ adia-
batically following. This can be seen by means of the fol-
lowing approximate solution of Eqs. (4): considering that
the equation for ρ is of lesser order [O(ǫ)], with respect
to those for θ and u [O(1)], we assume ρ = exp(−2x2)
unchanged in z and solve exactly the third of Eq. (4) for
θ(x) (though derived easily, its full expression is quite
cumbersome). Then, applying the theory of characteris-
tics [1], the second of Eqs. (4) is reduced to the following
ODEs, where dot stands for d/dz
ẋ = u ; u̇ = χθx, (5)
equivalent to the motion of a unit mass in the potential
V (x) = −χθ with conserved energy E = u(z)
+ V (x).
The solution of Eqs. (5) with initial condition x(0) =
s, u(0) = 0 yields x(s, z), u(s, z) in parametric form, from
which overtaking is found whenever u(x, z) (obtained
by eliminating s) becomes a multivalued function of x
at finite z = zs. The shock point corresponding to
|du/dx| → ∞ is found from the solution u(x, z) displayed
in Fig. 2(a) [ 2(b)], at positions x = ±xs 6= 0 (defocusing
FIG. 1: (Color online) 1D spatial profiles of phase chirp u(x)
(a-b) and amplitude |ψ(x)| =
ρ(x) (c-d), as obtained from
Eqs. (2-3) with ε = 10−3, α = 0, χ = −1 (defocusing), and
ψ0 = exp(−x
2), for different z as indicated: (a-c) local case,
σ2 = 0; (b-d) nonlocal case, σ2 = 5.
FIG. 2: (Color online) (a) u(x) for different z and χ = 1
(focusing), σ = 1; (b) as in (a) for χ = −1 (defocusing); (c)
shock distance zs (χ = −1 bold solid, χ = 1 thin solid) and
shock position xs in the defocusing case (dashed line).
case) or xs = 0 (focusing case). The shock distance zs
increases with σ in both cases, as shown in Fig. 2(c).
We have tested these predictions by integrating nu-
merically Eqs. (2-3). Simulations with χ = −1 [see Fig.
1(b,d)] show indeed steepening and post-shock oscilla-
tions in the spatial chirp u, which are accompanied by
a steep front in ρ moving outward. The shock location
in x and z is in good agreement with the results of our
approximate analysis summarized in Fig. 2.
Numerical simulations of Eqs. (2-3) validates also
the focusing scenario. The field evolution displayed in
Fig. 3(a) exhibits shock formation at the focus point
(xs = 0, zs ≃ 8, for σ = 5) driven the phase whose chirp
is shown in Fig. 3(b). This is remarkable because, in
the local limit σ = 0, the celerities become imaginary
(the equivalent gas would have pressure decreasing with
increasing density ρ), and no shock could be claimed to
exist. In this limit, the reduced problem (4) is elliptic and
the initial value problem is ill-posed [13], an ultimate con-
sequence of the onset of MI: modes with transverse (nor-
FIG. 3: (Color online) Level plot of the intensity in the
focusing case (χ = 1, ε = 0.01): (a) nonlocal case (σ2 =
25); (b) chirp profile for various z for (a); (c) quasi-local case
(σ2 = 10−5).
malized) wavenumber q < ∆q grow exponentially with z,
with both gain g and bandwidth ∆q scaling as 1/ε. How-
ever, the nonlocal response tends to frustrate MI (see also
Refs. [9, 12]), as shown by standard linear stability anal-
ysis which yields g =
d(2χ− d)/ε2 (we set d ≡ ε2q2/2
and χ ≡ χ/(1 + σ2q2)), in turn implying a strong reduc-
tion of both gain and bandwidth for large σ. In order
to emphasize the difference between the local and non-
local regime, we contrast Fig. 3(a) with the analogous
evolution [see Fig. 3(c)] obtained in the quasi-local limit
(σ2 = 10−5), which appears to be clearly dominated by
filamentation.
Thermal nonlinearity The physics of the defocusing
case can be experimentally tested by exploiting ther-
mal nonlinearities of strongly absorptive bulk samples,
that we show below to fit our model. In this case, the
system relaxes to a steady-state refractive index change
∆n = (dn/dT )∆T , where dn/dT is the thermal coef-
ficient, and ∆T the local temperature change due to
optical absorption. It is well known that this so-called
thermal lens distorts a laser beam propagating in the
medium [14, 15, 16]. However, only perturbative ap-
proaches to the problem have been proposed (ray optics
or Fresnel diffraction theory is applied after the lens pro-
file is worked out from gaussian ansatz [14]), while the
role of shock phenomena was completely overlooked.
FIG. 4: 2D evolution according to Eqs. (2-3) with σ2 = 1,
α = 1: (a) radial phase chirp at different z, as indicated,
showing steepening and shock formation for ε = 10−2; (b) cor-
responding intensity profile |ψ|2 (maximum scaled to unity)
at z = 4.9; (c) transverse intensity profile vs. x (at y = 0) at
z = 1/(4ε) and different values of ε (α0 = 62cm
−1, σ = 0.3).
We assume that the temperature field ∆T =
∆T⊥(X,Y ) obeys the following 2D heat equation
(∂2X + ∂
Y )∆T⊥ − C∆T⊥ = −γ|A|2 (6)
where the source term account for absorption pro-
portional to intensity through the coefficient γ =
α0/(ρ0cpDT ), where ρ0 the material density, cp the spe-
cific heat at constant pressure, and DT is the thermal
diffusivity (see e.g. [16]). Eq. (6) has been already em-
ployed to model a refractive index of thermal origin in
Ref. [10], and in Ref. [11] in the limit C = 0 which
is equivalent to consider the range of nonlocality (mea-
sured by 1/C, see below) to be infinite. Starting from
the 3D heat equation ∇2∆T = −γ|A|2, the latter regime
amounts to assume ∆T (X,Y, Z) = ∆T⊥(X,Y ), which is
justified when longitudinal changes in intensity |A|2 are
negligible as for solitary (invariant in Z) wave-packets
in the presence of low absorption [11]. Viceversa, in the
regime of strong absorption, we need to account for lon-
gitudinal temperature profiles that are known from solu-
tions of the 3D heat equations to be peaked at charac-
teristic distance Ẑ in the middle of sample and decay to
room temperature on the facets [14]. Since highly non-
linear phenomena occurs in the neighborhood of Ẑ where
the index change is maximum, we can use a (longitudi-
nal) parabolic approximation with characteristic width
Leff (∼ L) of the 3D temperature field ∆T (X,Y, Z) =
1− (Z−Ẑ)
∆T⊥(X,Y ) and consequently approximate
∇2∆T ≃ (∂2X + ∂2Y )∆T⊥ − L
eff∆T⊥, so that the 3D
heat equation reduces to Eq. (6) with C = 1/L2eff .
Following this approach, Eq. (6) coupled to Eq. (1)
can be casted in the form of Eqs. (2-3) by posing
θ = k0Lnl|dn/dT |∆T⊥ and σ2 = 1/(Cw20) = L2eff/w20.
The model reproduces the infinite range nonlocality for
negligible losses (Leff → ∞); while for thin samples
[|(∂2X + ∂2Y )∆T⊥| << |L
eff∆T⊥|], Leff can be related
to the Kerr coefficient n2 as
Leff =
γ|dn/dT |
DTρ0cp|n2|
α0|dn/dT |
which establishes a link between the degree of nonlocal-
ity and the strength of the nonlinear response (similarly
to other nonlocal materials [12]).
Having retrieved the model Eqs. (2-3), let us show
next that the scenario illustrated previously applies sub-
stantially unchanged in bulk (2D case) even on account
for the optical power loss (α 6= 0). An example of the
general dynamics is shown in Fig. 4, where we report a
simulation of the full model (2-3), with σ2 = 1 and rela-
tively large loss α = 1. In analogy to the 1D case, Fig.
4(a) clearly shows that the radial chirp u = φr steep-
ens and then develop characteristic oscillations after the
shock point (z ≃ 6, where |∂ru| → ∞). Correspond-
ingly the intensity exhibits also an external front which
is connected to a flat central region with a characteris-
tic overshoot [see Fig. 4(b)] corresponding to a brighter
ring [inset in Fig. 4(c)]. For larger distances this struc-
ture moves outward following the motion of the shock.
In the experiment such motion can be observed, at fixed
physical lenght, by increasing the power, which amounts
to decrease ε while scaling z and α accordingly (z ∝ 1/ε,
α ∝ ε), as displayed in Fig. 4(c) for σ = 0.3.
As a sample of a strongly absorbing medium we choose
a 1 mm long cell filled with an acqueous solution of
Rhodamine B (0.6 mM concentration). Our measure-
ments of the linear and nonlinear properties of the sam-
ple performed by means of the Z-scan technique gives
data consistent with the literature [17], and allows us
to extrapolate at the operating vacuum wavelength of
532 nm, a linear index n = 1.3, a defocusing nonlin-
ear index n2 = 7 × 10−7 cm2W−1, and α0 = 62 cm−1.
For our sample DT = 1.5 × 10−7 m2s−1, ρ0 = 103 kg
m−3, cp = 4 × 103Jkg−1K−1 and |dn/dT | = 10−4 K−1
(γ ∼= 104 K W−1), and exploiting Eq. (7) we estimate
Leff ∼= 10µm (Leff << L because of the strong ab-
sorption that causes strong heating of our sample near
the input facet), and correspondingly the degree of non-
locality σ ∼= 0.3. We operate with an input gaussian
beam with fixed intensity waist w0I = w0/
2 = 20 µm
(Ld ∼= 12 mm) focused onto the input face of the cell.
With these numbers, an input power P = πw20II0 = 200
mW yields a nonlinear length Lnl ∼= 8 µm (L ∼= 0.3 mm),
which allows us to work in the semiclassical regime with
ε ∼= 0.025. The radial intensity profiles together with the
2D patterns imaged by means of a 40×microscope objec-
tive and a recording CCD camera are reported in Fig. 5.
As shown the beam exhibits the formation of the bright
ring whose external front moves outward with increasing
power, consistently with the reported simulations. We
point out that, at higher powers, we observe (both ex-
perimentally and numerically) that the moving intensity
front leaves behind damped oscillations that correspond
to inner rings of lesser brightness, as reported in litera-
ture [15]. This, however, occurs well beyond the shock
point that we have characterized so far.
In summary, the evolution of a gaussian beam in the
strong nonlinear regime is characterized by occurence of
collisionless (i.e., regularized by diffraction) shocks that
survive the smoothing effect of (even strong) nonlocality.
While experimental results support the theoretical sce-
nario in the defocusing case, the remarkable result that
the nonlocality favours shock dynamics over filamenta-
tion requires future investigation.
∗ Electronic address: claudio.conti@phys.uniroma1.it
[1] G. B. Whitman, Linear and Nonlinear Waves (Wiley,
New York, 1974);
[2] L. D. Landau and E. M. Lifshitz, Fluid Mechanics (Perg-
amon, 1995); M. A. Liberman and A. L. Yelikovich,
Physics of Shock Waves in Gases and Plasmas (Springer,
FIG. 5: Radial profiles of intensity observed in the thermal
medium for different input powers. The insets show the cor-
responding 2D output patterns.
Heidelberg, 1986).
[3] A.V. Gurevich and L.P. Pitaevskii, Sov. Phys. JETP 38,
291 (1973); A.V. Gurevich and A. L. Krylov, Sov. Phys.
JETP 65, 944 (1987).
[4] J. C. Bronski and D. W. McLaughlin, in Singular Limits
of Dispersive Waves, NATO ASI Series, Ser. B 320, pp.
21-28 (1994); M. G. Forest and K. T. R. McLaughlin,
J. Nonlinear Science 7, 43 (1998); Y. Kodama, SIAM J.
Appl. Math. 59, 2162 (1999). M. G. Forest, J. N. Kutz,
and K. T. R. McLaughlin, J. Opt. Soc. Am. B 16, 1856
(1999).
[5] A. M. Kamchatnov, R. A. Kraenkel, and B. A. Umarov,
Phys. Rev. E 66, 036609 (2002).
[6] R. J. Taylor, D.R. Baker, and H. Ikezi, Phys. Rev. Lett.
24, 206 (1970).
[7] J. E. Rothenberg and D. Grischkowsky, Phys. Rev. Lett.
62, 531 (1989).
[8] M. A. Hoefer, M. J. Ablowitz, I. Coddington, E. A.
Cornell, P. Engels, and V. Schweikhard, Phys. Rev. A
74, 023623 (2006). V. M. Perez-Garcia, V.V. Konotop,
V.A. Brazhnyi, Phys. Rev. Lett. 92, 220403 (2004); B.
Damski, Phys. Rev. A 69, 043610 (2004);
[9] J. Wyller, W. Krolikowski, O. Bang, J. J. Rasmussen,
Phys. Rev. E 66, 066615 (2002).
[10] A. Yakimenko, Y. Zaliznyak and Y.S. Kivshar, Phys.
Rev. E 71, 065603(R) (2005).
[11] C. Rotschild, O. Cohen, O. Manela, M. Segev and T.
Carmon, Phys. Rev. Lett. 95, 213904 (2005).
[12] C. Conti, M. Peccianti and G. Assanto, Phys. Rev. Lett.
91 073901 (2003); Phys. Rev. Lett. 92 113902 (2004);
C. Conti, G. Ruocco and S. Trillo, Phys. Rev. Lett. 95
183902 (2005).
[13] P. D. Miller and S. Kamvissis, Phys. Lett. A 247, 75
(1998); J. C. Bronski, Physica D 152, 163 (2001).
[14] C. A. Carter and J. M. Harris, Appl. Opt. 23, 476 (1984);
S. Wu and N. J. Dovichi, J. Appl. Phys. 67, 1170 (1990);
F. Jürgensen and W. Schröer, Appl. Opt. 34 41 (1995).
[15] C. J. Wetterer, L. P. Schelonka, and M. A. Kramer, Opt.
Lett. 14, 874 (1989).
mailto:claudio.conti@phys.uniroma1.it
[16] P. Brochard, V. Grolier-Mazza and R. Cabanel, J. Opt.
Soc. Am. B 14, 405 (1997).
[17] S. Sinha, A. Ray, and K. Dasgupta, J. Appl. Phys. 87,
3222 (2000).
[18] paraxial diffraction in defocusing media is well known
to be isomorphus in 1D to propagation in a normally
dispersive focusing medium as considered in Ref. [7]
