Dynamics of self-focused femtosecond laser

pulses in the near and far fields

Alex A. Zozulya
Department of Physics, Worcester Polytechnic Institute,
100 Institute Rd., Worcester, MA 01609-2280, USA

[email protected]

Scott A. Diddams
JILA, University of Colorado,
and National Institute of Standards and Technology
Boulder, CO 80309-0440, USA

[email protected]

Abstract: We investigate the propagation of femtosecond pulses in

a nonlinear, dispersive medium at powers several times greater than
the critical power for self focusing. The combined effects of diffraction,
normal dispersion and cubic nonlinearity lead to pulse splitting. We
show that detailed theoretical description of the linear propagation of
the pulse from the exit face of the nonlinear medium (near field) to the
measuring device (far field) is crucial for quantitative interpretation of
experimental data.
1998
c Optical Society of America
OCIS codes: (320.7110) Ultrafast nonlinear optics; (190.5530) Pulse propagation
and solitons

References
1. A. Braun, G. Korn, X. Liu, D. Du, J. Squier and G. Mourou, “Self-channeling of high-peak-power
femtosecond laser pulses in air,” Opt. Lett. 20, 73-75 (1995).
2. E. T. J. Nibbering, P. F. Curley, G. Grillon, B. S. Prade, M. A. Franco, F. Salin and A. Mysy-
rowicz, “Conical emission from self-guided femtosecond pulses in air,” Opt. Lett. 21, 62-64
(1996).
3. D. Strickland and P. B. Corkum, “Resistance of short pulses to self-focusing,” J. Opt. Soc. Am.
B 11, 492-497 (1994).
4. J. K. Ranka, R. W. Schirmer and A. L. Gaeta, ”Observation of pulse splitting in nonlinear
dispersive media,” Phys. Rev. Lett. 77, 3783-3786 (1996).
5. S. A. Diddams, H. K. Eaton, A. A Zozulya and T. S. Clement, “Amplitude and phase measure-
ments of femtosecond pulse splitting in nonlinear dispersive media,” Opt. Lett. 23, 379-381
(1998).
6. R. L. Fork and C. V. Shank and C. Hirlimann and R. Yen, “Femtosecond white-light continuum
pulses,” Opt. Lett. 8, 1-3 (1983).
7. P. B. Corkum, C. Rolland and T. Srinivasan-Rao, “Supercontinuum generation in gases,” Phys.
Rev. Lett. 57, 2268-2271 (1986).
8. N. A. Zharova, A. G. Litvak, T. A. Petrova, A. M. Sergeev and A. D. Yanukoviskii, “Multiple
fractionation of wave structures in a nonlinear medium,” JETP Lett. 44, 13-17 (1986).
9. P. Chernev and V. Petrov, “Self-focusing of light pulses in the presence of normal group-velocity
dispersion,” Opt. Lett. 17, 172-174 (1992).
10. J. Rothenberg,“Pulse splitting during self-focusing in normally dispersive media,” Opt. Lett.
17, 583-585 (1992).
11. G. G. Luther, J. V. Moloney, A. C. Newell and E. M. Wright, “Self-focusing threshold in normally
dispersive media’” Opt. Lett. 19, 862-864 (1994).
12. A. A. Zozulya, S. A. Diddams, A. G. Van Engen and T. S. Clement, “Propagation dynamics of
intense femtosecond pulses: Multiple splittings, coalescence, and continuum generation,” Phys.
Rev. Lett. 82, 1430-1433 (1999).

#9461 - $15.00 US Received March 22, 1999; Revised April 18, 1999
(C) 1999 OSA 26 April 1999 / Vol. 4, No. 9 / OPTICS EXPRESS 336
 13. R. Trebino, K. W. DeLong, D. N. Fittinghoff, J. N. Sweetser, M. A. Krumbügel and B. A.
Richman , “Measuring ultrashort laser pulses in the time-frequency domain using frequency-
resolved optical gating,” Rev. Sci. Inst. 68, 3277-3295 (1997).
14. A. A. Zozulya, S. A. Diddams and T. S. Clement, “Investigations of nonlinear femtosecond pulse
propagation with the inclusion of Raman, shock, and third-order phase effects,” Phys. Rev. A
58, 3303-3310 (1998). This work and Ref. [12] contain extensive references to related theoretical
work with modified nonlinear Schrödinger equations.
15. J. Rothenberg,“Space-time focusing: breakdown of the slowly varying envelope approximation
in the self-focusing of femtosecond pulses,” Opt. Lett. 17, 1340-1342 (1992).
16. J. K. Ranka and A. L. Gaeta, “Breakdown of the slowly varying envelope approximation in the
self-focusing of ultrashort pulses,” Opt. Lett. 23, 534-536 (1998).
17. F. DeMartini, C. H. Townes, T. K. Gustafson and P. L. Kelley, “Self-steepening of light pulses,”
Phys. Rev. 164, 312-323 (1967).
18. R. H. Stolen and W. J. Tomlinson, “Effect of the Raman part of the nonlinear refractive index
on propagation of ultrashort optical pulses in fibers, ” J. Opt. Soc. Am. B 9, 565-573 (1992)

The broad spectral bandwidths, high peak powers, and (3+1)-d nature of in-
tense femtosecond pulses result in a variety of complex and interesting propagation
effects. These include self-channeling[1,2], temporal break-up of the pulse[3,4,5] and
extreme spectral broadening–commonly called continuum generation[6,7]. To a large
extent, these phenomena rely on the basic process of self-focusing, due to an index
of refraction proportional to the intensity. However, other physical mechanisms turn
out to be of equal importance when a broadband ultrashort pulse is considered. Mate-
rial dispersion is typically one of the most important additional effects influencing the
dynamics of an ultrashort pulse. Indeed, at moderate powers normal group velocity dis-
persion acts to arrest catastrophic collapse with the result of temporal splitting of the
pulse into two[3,9,10,11]. However, recent observations of multiple splittings, reported in
Refs. [4,5,12], indicate that sharp temporal features, small beam diameters, and broad
bandwidths further require that space-time coupling, nonlinear shock, higher-order dis-
persion, and ionization effects be included in theoretical analysis aimed at quantitative
comparison with experiments. The importance of these higher order terms is elucidated
with “full-field” (amplitude and phase) measurements as provided by frequency-resolved
optical gating (FROG)[13]. It also becomes evident that quantitative comparison be-
tween theory and experiment requires careful theoretical description of the measurement
process, including propagation of the pulse from the nonlinear medium to a measuring
device. Due to the broad spectrum and small spatial dimensions, the vacuum propaga-
tion modifies the spatio-temporal characteristics of the pulse. In the present commu-
nication we will illustrate this point by presenting quantitative analysis of our results
regarding propagation of femtosecond pulses in a nonlinear, dispersive medium at pow-
ers several above the critical power for self focusing.
1. Model
We model the evolution of the complex envelope E(r, z, t) of the radially-symmetric
field
E(r, z, t) = E(r, z, t) exp(ikz − iω0 t) + c.c (1)
with the following modified nonlinear Schrödinger equation [12,14]:
 
∂ ∂ ∂2 ∂3
i E + 1 − iω ∇2 E − 2 E − i3 3 E
∂z ∂t ∂t ∂t
 
∂
+ 1 + iω g(|E|2 )E = 0 (2)
∂t
In Eq. 2, the time t is measured in the frame moving at the group velocity of the
pulse. The transverse Laplacian in cylindrical coordinates ∇2 = ∂ 2 /∂r2 + (1/r)∂/∂r

#9461 - $15.00 US Received March 22, 1999; Revised April 18, 1999
(C) 1999 OSA 26 April 1999 / Vol. 4, No. 9 / OPTICS EXPRESS 337
accounts for diffraction, while the second and third time derivatives describe group ve-
locity dispersion (GVD) and third order dispersion (TOD). The temporal, longitudinal,
and transverse coordinates are normalized to the characteristic pulse duration p τ , the
dispersion length lD = 2τ 2 /|k 00 |, and the characteristic transverse length l⊥ = lD /2k,
respectively. In addition, 3 = k 000 /(3k 00 τ ), and k = 2πn/λ, with n being the linear index
of refraction at the central wavelength λ. The dispersion coefficients k 00 and k 000 are the
second and third derivatives of k with respect to frequency, evaluated at the central
frequency ω0 . Space-time focusing [15] is given by the time derivative of the Laplacian,
while nonlinear shock [17] is given by the time derivative of the nonlinearity [g(|E|2 )].
As shown in Eq. (2), both of these terms are proportional to ω = 1/ω0 τ . These two
terms act together to shift energy towards the trailing edge of the pulse [14,15,16]. Due
to the short duration of the pulse, it is important to account for both instantaneous and
time-delayed Raman nonlinearities [18], such that
 Z t 
2 2πn2 lD 2 2
g(|E| ) = (1 − α)|E(t)| + α dτ f (t − τ )|E(τ )| ] , (3)
λ −∞

with
1 + (ωr τr )2
f (t) = exp(−t/τr ) sin(ωr t). (4)
ωr τr2
In Eq. (4), n2 is the nonlinear index of refraction, and α denotes the fractional amount
of the nonlinearity due to the Raman effect.
Unless noted, in the following numerical solutions of Eq. (2) we use parameters
typical of our recent experiments in fused silica[12]. The initial field is taken to be a
hyperbolic secant in time and a real Gaussian in space, having intensity full-width at
half maxiumum (FWHM) of 90 fs and 70 µm, respectively. Furthermore, the beam
waist is located at the entrance face of the sample. The linear index of refraction is
n = 1.45 at the center wavelength of λ = 0.8 µm. The nonlinear index of refraction
is n2 = 2.5 × 10−16 GW/cm2 , resulting in a self-focusing critical power of Pcrit =
(0.61λ)2 π/(8nn2 ) = 2.6 MW. The GVD and TOD coefficients are k 00 = 360 fs2 /cm,
and k 000 = 275 fs3 /cm. For the Raman response of fused silica we use α = 0.15,
τr = 50 fs, and ωr τr = 4.2[18]. Equation (2) is solved using a symmetric split-step
technique. During the linear part of each step, the time-dependent terms are solved
in the frequency domain, while the transverse Laplacian is computed using a finite
differences scheme. The time-dependent nonlinear part of the equation is also solved
using finite differences.
For accurate comparison with our experiments, it is important to account for
the propagation of the field from the output face of the nonlinear medium (fused silica)
to the measurement apparatus. In our experiments after propagation in the fused silica,
the entire field is allowed to diffract in air over 1.5 m. At this point an aperture selects
the on-axis portion of the field for characterization by the second-harmonic generation
(SHG) FROG apparatus[5,12]. Mathematically, we describe the linear propagation of
the pulse from the exit face of the nonlinear medium (near field) to the measuring device
(far field) by the equation:
 
1 ∂ 2 ∂
i E + ∇ 1 − iω E=0 (5)
n ∂z ∂t
where we are using longitudinal, transverse, and time scales introduced above for the
nonlinear medium. The appearance of the linear refractive index n in the first term of
Eq. (5) is due to the fact that the wave vector of the field in air (assumed to be vacuum)
is smaller than in the medium by a factor of n. The amplitude of the radially-symmetric

#9461 - $15.00 US Received March 22, 1999; Revised April 18, 1999
(C) 1999 OSA 26 April 1999 / Vol. 4, No. 9 / OPTICS EXPRESS 338
field at distance L from the nonlinear medium follows from the solution of Eq. (5) and
is given by the equation

Evac (r, L, t) =
Z    
−i dωdr0 r0 rr0 (r2 + r02 )
E(r0 , 0, ω)J0 exp i − iωt (6)
4πLn 1 − ω ω 2(1 − ω ω)Ln 4(1 − ω ω)Ln

where J0 is the zeroth-order Bessel function. Equation (6) is evaluated numerically

using the frequency domain representation E(r0 , 0, ω) of the field at the output of the
nonlinear medium. Finally we note that the on-axis in the far field (r = 0, L → ∞), the
amplitude of the pulse simplifies to the expression
Z ∞ Z Z ∞
E(r0 , 0, ω)
Evac (0, ∞, t) ∝ dr0 r0 dω exp(−iωt) ≈ dr0 r0 E(r0 , 0, t) (7)
0 1 − ω ω 0

This simplified expression tells us that the on-axis far field at time t can be approximated
by the integral over the cross-section of the near field at the same local time. This result
will be used below when discussing the physics of pulse splitting.
2. Results

Figure 1. (a) Animated surface plot of I(r, t) of a self-focusing femtosecond

pulse in fused silica for an input power of 4 MW (583kB QuickTime movie).
(b) Calculated intensity profile I(r, t) in the far field (z = 1.5 m) for the field
shown in (a).

The importance of detailed description of vacuum propagation is illustrated in Figure 1.

The animated surface plot of Figure 1(a) shows the basic dynamics of pulse splitting at
relatively small input powers as predicted by Eq. (2). The intensity I(r, t) is shown as a
function of the propagation distance z (z runs from 0 to 30 mm). The spatial dimension
is radially-symmetric, with position r = 0 at the center of the axis. The orientation of
the figure is such that early times are at the back of the figure. The peak intensity at
z = 0 cm is 72 GW/cm2 , corresponding to a peak power of 4 MW.
This input peak power is greater than the critical power, such that moderate
self-focusing occurs about t = 0; however, it is not sufficient to cause pulse splitting
inside the nonlinear medium. Nevertheless, splitting occurs as the pulse propagates in
vacuum to the far field [Figure 1(b)]. As discussed above, the free propagation of the field

#9461 - $15.00 US Received March 22, 1999; Revised April 18, 1999
(C) 1999 OSA 26 April 1999 / Vol. 4, No. 9 / OPTICS EXPRESS 339
from the output of the fused silica to the measuring apparatus involves the evaluation
of Eq. (6). Figure 1(b) shows the result of this calculation of the intensity distribution
I(r, t) at z = 1500 mm beyond the output of the 30 mm fused silica media. This far-field
spatio-temporal shape of the pulse clearly demonstrates splitting. Numerical simulations
show that for the above input intensity the pulse splitting becomes noticeable several
centimeters after the exit face of the nonlinear medium.
We can understand the far-field splitting with reference to Eq. (7), which shows
that the axial far-field amplitude can be approximated by the integral over the near-field
spatial cross section at the same local time. As such, a maximum in the on-axis field
can come either from a cross section with the largest values of the field, or from a cross
section with possibly smaller values of the field, but nearly constant phase. This second
possibility is quite generally created during the process of self-focusing, where the spatial
phase written onto the field varies rapidly about the peak intensity (t = 0), but is more
constant in the wings of the field. As a result, the differently-phased portions of the
field about t = 0 in the near field may destructively interfere in the far field leading to
a local minimum and the observed temporal pulse splitting of Figure 1(b).

Figure 2. (a) Animated surface plot of I(r, t) for input power of 4.7 MW
(583kB QuickTime movie). (b) Calculated intensity profile I(r, t) in the far
field for the field shown in (a). Early times are at the back of the figure.

The animated surface plot of Figure 2(a) illustrates the basic dynamics of pulse
splitting at moderate input powers. The peak intensity at z = 0 cm is 85 GW/cm2 ,
corresponding to a peak power of 4.7 MW. In Figure 2(a), we observe that an initially
uniform Gaussian input self-focuses in both space and time before splitting into two
separate pulses. The pulse splitting results when self-focusing moves off-axis energy
towards the peak of the pulse, while positive dispersion acts to pull this energy away
from t = 0. As this process continues, the peak intensity drops, stopping the collapse
at t = 0. However, off-axis energy continues to focus at t 6= 0 such that two pulses are
resolved.
The temporal asymmetries seen in Figure 2(a) are predominantly the result of
the space-time focusing and nonlinear shock terms of Eq. (2)[16,14]. Both of these terms
act to move energy towards the back edge of the pulse, creating a trailing shock front
that is seen around z = 19 mm in Figure 2(a). The resulting shock front causes extensive
broadening on the blue side of the spectrum via self phase modulation. Although this

#9461 - $15.00 US Received March 22, 1999; Revised April 18, 1999
(C) 1999 OSA 26 April 1999 / Vol. 4, No. 9 / OPTICS EXPRESS 340
shock formation initially results in a trailing pulse with higher intensity and shorter
duration, this pulse spreads much faster than the longer and lower intensity leading
pulse as seen for z > 20 mm in Figure 2(a). We also note that the Raman nonlinearity,
which is included Eq. (2), counteracts the nonlinear shock and space-time focusing terms
by transferring energy to the leading (red-shifted) pulse. However, this is a smaller effect.

10
1.0 (a) 1 (c)
0
0.1

Phase (rad)

Spectrum
Intensity

0.5

-10 0.01

0.0 1E-3
-20
5
1.0
(b) 1 (d)
0
Phase (rad)
0.1

Spectrum
Intensity

0.5

-5 0.01

0.0 -10 1E-3

-200 0 200 700 750 800 850 900

Tim e (fs) W avelength (nm)

Figure 3. (a) Axial intensity (black) and phase (red) measured with
SHG-FROG. (b) Calculated axial intensity and phase (c) Measured axial
spectrum of field shown in (a). The blue line is the square modulus of the
Fourier transform the data of (a), while the red points are measured with a
spectrometer. (d) Calculated axial spectrum. Frame (a) of this figure is linked
to an animated graphic with sound (330 kB QuickTime movie with sound)
that allows the reader to see and hear the frequency variations associated
field.

In the case of Figure 2(a), the transformation from the near to the far field can
result in multiple splitting of the field as shown in Figure 2(b). While there were two
pulses at the output of the nonlinear medium [z = 30 mm frame of Figure 2(a)], there
are now three pulses near r = 0 in the far field. This is shown in more detail by the
data of Figure 3. In this figure, we present a comparison of the measured and calculated
axial field for the same conditions shown in Figure 2. Figure 3(a) is the measured
axial (r = 0) intensity and phase of the complex envelope E(r = 0, z = 1500 mm,t).
This data was acquired using the SHG-FROG technique under experimental conditions
corresponding to those of the calculation of Figure 2. For comparison, Figure 3(b)
is the calculated axial intensity and phase, demonstrating good agreement with the
measurement. The measured and calculated frequency domain representation of the
field is shown in Figure 3(c) and (d), respectively. In Figure 3(c) the blue line is the
square modulus of the Fourier transform of the measured field of Figure 3(a), while
the red points are an independently measured spectrum (using a 0.27 m spectrometer).
The discrepancy seen on the short wavelength side of the spectrum is due to bandwidth
limitations in the SHG-FROG measurement. Nonetheless, this data verifies the fidelity

#9461 - $15.00 US Received March 22, 1999; Revised April 18, 1999
(C) 1999 OSA 26 April 1999 / Vol. 4, No. 9 / OPTICS EXPRESS 341
of the FROG measurement over several orders of magnitude.

1.5
(a) (b)
1.0

Axial Intensity
0.5

t r 0.0
-200 0 200
Time (fs)

1.0
(c) (d)
Axial Intensity

0.5

t r
0.0
-200 0 200
Time (fs)
Figure 4. (a) Intensity distribution I(r, t) after 30 mm of propagation in
fused silica. The input power in this case is 5.5 MW, corresponding to an
intensity of 100GW/cm2 . (b) Axial intensity corresponding to (a). (c) I(r, t)
in the far field for z = 1.5 m beyond the exit face of the fused silica. (d) Axial
intensity corresponding to (c)

The broadening of the pulse spectrum (to ∼ 150 nm) as seen in Figure 3(c) and
(d) is the result of intensity induced modulation of the phase of the field. For reference,
the near-transform-limited input field has a spectral width of ∼ 45 nm at the 10−3 level.
As means of further illustrating the self-phase modulation that results in the spectral
broadening, we present animated video with sound that is linked to Figure 3(a). In this
multimedia clip the intensity (blue dots) and phase (black line) of Figure 3(a) are again
shown; however, the pulse phase is additionally presented as a “sonogram.” The pitch
of the sound one hears when playing the movie is proportional to the local frequency
(indicated by the moving red point). This local, or instantaneous, frequency is defined as
ωinst = ω0 −∂φ/∂t with φ(t) being the phase of the complex envelope. One will note that
sharp negative slopes of the phase correspond to higher pitch and therefore correspond
to blue shifts. On the other hand, positive slopes in the pulse phase correspond to red
shifts as indicated by the lower pitch.

#9461 - $15.00 US Received March 22, 1999; Revised April 18, 1999
(C) 1999 OSA 26 April 1999 / Vol. 4, No. 9 / OPTICS EXPRESS 342
 At yet higher input power the spatio-temporal evolution of the pulse is charac-
terized by further broadening of its spectrum, especially on the blue side. The resulting
temporal dynamics can be described as a far field coalescence of the pulse toward a
single broad pulse as illustrated by Figure 4. Figure 4(a) is the calculated intensity
profile I(r, t) at the output of 30 mm of fused silica for an input power if 5.5 MW. The
axial intensity in this case is shown in Figure 4(b). Of significance here is the absence of
repeated multiple splittings of the form seen in Figure 2. This implies that the original
two pulses of Figure 2(a) do not simply create two additional daughter pulses of their
own. In the far field, shown in Figure 4(c), we see the field beginning to coalesce to-
wards a single dominant peak. Our measurements at this and higher powers show good
agreement with these calculations[12].
The fact that the axial intensity is larger in the far field at times where minimal
intensity is seen in the near field can be qualitatively explained by Eq. (7). It is from
the smaller intensity region about +75 fs in the near field that the large trailing peak
arises in the far field. As already described above, smaller nonlinear phase variations are
written onto the less intense regions of the near field, resulting in reduced destructive
interference in the far field at the same local time.
3. Acknowledgement
This work was supported in part by the National Science Foundation and the National
Institute of Standards and Technology. S. A. Diddams thanks T. S. Clement for valuable
input, financial support and the use of laboratory equipment.

#9461 - $15.00 US Received March 22, 1999; Revised April 18, 1999
(C) 1999 OSA 26 April 1999 / Vol. 4, No. 9 / OPTICS EXPRESS 343