NLO

Download as pdf or txt
Download as pdf or txt
You are on page 1of 55

F RIEDRICH -S CHILLER-U NIVERSITÄT J ENA

P HYSIKALISCH -A STRONOMISCHE -FAKULTÄT

S OMMERSEMESTER 2021

Theory of Nonlinear Optics

U LF P ESCHEL

LATEX-Satz und Design von Martin Beyer


Contents
1 Phenomenological representation of non-resonant nonlinearities 3
1.1 The polarization as the origin of nonlinearity in optics . . . . . . . . . . . . . . . 3
1.2 The Taylor expansion of the polarization . . . . . . . . . . . . . . . . . . . . . . . 4
1.3 Polarization in frequency space . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5
1.4 Symmetries and conserved quantities . . . . . . . . . . . . . . . . . . . . . . . . . 10
1.5 Cubic nonlinearities in isotropic materials . . . . . . . . . . . . . . . . . . . . . . 11
1.5.1 Third harmonic generation (THG) . . . . . . . . . . . . . . . . . . . . . . . 12
1.5.2 Self phase modulation (SPM) or optical Kerr effect . . . . . . . . . . . . . 14

2 Resonant nonlinearities - a quantum mechanical representation 16


2.1 Equations of motion of a two-level-system . . . . . . . . . . . . . . . . . . . . . . 16
2.2 Rabi oscillations . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 21
2.3 The influence of damping . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 23

3 Nonlinear optics and field propagation - a plane wave approach 27


3.1 Perturbation theory for nonlinear solutions . . . . . . . . . . . . . . . . . . . . . 27
3.2 Propagation in z-direction only . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 28
3.3 Diffraction in linearly isotropic media . . . . . . . . . . . . . . . . . . . . . . . . . 30
3.4 Pulse evolution in linearly isotropic media . . . . . . . . . . . . . . . . . . . . . . 31

4 Field evolution for quadratic nonlinearities 34


4.1 Equations of motion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 34
4.2 Normalization . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 36
4.3 Conservation laws - the Hamiltonian . . . . . . . . . . . . . . . . . . . . . . . . . 37
4.4 General solution . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 38
4.5 Phase-matched SHG . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 39
4.6 Down conversion (phase matched case) . . . . . . . . . . . . . . . . . . . . . . . 40
4.7 Stability of the second harmonic field . . . . . . . . . . . . . . . . . . . . . . . . . 41
4.8 Visualization of the field dynamics using the Hamiltonian . . . . . . . . . . . . . 44
4.9 Low-depletion-limit . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 45

5 Diffraction, dispersion and nonlinearities 47


5.1 Beam evolution in isotropic Kerr materials . . . . . . . . . . . . . . . . . . . . . . 47
5.2 Pulse evolution in a two-level medium . . . . . . . . . . . . . . . . . . . . . . . . 51

2
1 Phenomenological representation of non-resonant
nonlinearities

1.1 The polarization as the origin of nonlinearity in optics

We start with M AXWELL’s equations in the optical domain without external charges, currents
and no magnetization

⃗ ∂B ⃗
∇×E =− ∇·D =0 (1.1)
∂t
⃗ ∂D ⃗
∇ × B = µ0 ∇ · B = 0, (1.2)
∂t
where E ist the electric field, B the magnetic induction and D the dielectric displacement
which is defined as

D = ε0 E + P , (1.3)

consisting of a vacuum contribution ε0 E and the polarization P of the material. We note that
M AXWELL’s equations are linear in all components.
The polarization of the material P is the source of nonlinearities and is driven by the elec-
tromagnetic fields. It cannot be deduced from electrodynamics, but we need quantum me-
chanics/solid state physics. Furthermore it has to obey symmetries which will be discussed
in this chapter.
We define the polarization P (r , t ) at a fixed point in space r and time t . In principle it de-
pends on electromagnetic fields acting on all other places r ′ and other times t ′ as P (r , t ) =
P [E (r ′ , t ′ ), B (r ′ , t ′ )]. However, we make several assumptions to simplify our calculations:
• There are only local interactions: r = r ′ . We explicitly exclude the propagation of a
polarization in the material. This means the polarization of one unit of the material
is influence by light acting on just that unit only. The reason for that is that nonlocal
effects on atomic scales take place at ∼ 10−10 m while the variations of the optical field
are ∼ 10−6 m.
• Causality: The polarization cannot depend on fields from the future t ′ ≤ t .
• No external action on the polarization within the memory time of the material. The
polarization depends on time differences only P (r , t ) = P [E (r , t − τ), B (r , t − τ)] with
τ > 0.
• No magnetic effects P (r , t ) = P [E (r , t − τ)]
The above assumptions are well justified for conventional materials. Unfortunately, a further
simplification is not straightforward.

3
1.2 The Taylor expansion of the polarization NLO

1.2 The Taylor expansion of the polarization

The general idea ist that for non-resonant nonlinearities the interaction between photons
via the material is weak. This implies that only a very limited number of photons hits one
atom of the material within the relaxation time/memory of the polarization. We can classify
the polarization with respect to the number of interacting photons

P (t ) = P (1) (t ) + P (2) (t ) + P (3) (t ). (1.4)

The linear polarization P (1) (t ) is induced by a single photon and is by far the dominant po-
larization. P (2) (t ) is the nonlinear polarization of quadratic/second order and is induced by
two-photon processes.

a.) Linear polarization


We start by only considering one vector component P i(1) (t ) which is driven by all other vector
components of the electric field E j having interacted with the atom at all times t − τ1 in the
past

3
ˆ∞
P i(1) (t ) = ε0 dτ1 R i(1) (τ1 )E j (t − τ1 ),
X
j
(1.5)
j =1
0

where R i(1)
j
(τ1 ) is the linear response function. It is a time dependent second rank tensor
which is decaying rapidly in time.

b.) Nonlinear polarization of second order


We can apply the previous procedure to higher order interactions in the same way. Now we
assume that two photons have interacted at times t − τ1 and t − τ2 in the past with the atom.
They both have to be taken into account leading to

3
ˆ∞ ˆ∞
P i(2) (t ) = ε0 dτ2 R i(2) (τ , τ )E j (t − τ1 )E k (t − τ2 ).
X
dτ1 jk 1 2
(1.6)
j ,k=1
0 0

The response function R i(2)


jk
is now a third rank tensor which depends on two times. In a
similar way we can represent all higher order interactions R (n) .
The response function R (n) has some special properties:
• It is real valued, because we have only considered real valued fields.
• It must obey causality corresponding R (n) (τk < 0) = 0.
• It must be symmetric with respect to an exchange of arguments
R i(n) (n)
... m ... k ... (...τm ...τk ...) = R i ... k ... m ... (...τk ...τm ...)
• Has the the symmetries of the respective material.

4
1.3 Polarization in frequency space NLO

Symmetry consideration for inversion symmetric materials


The latter property shall be demonstrated for a quadratic nonlinearity represented by R i(2) (τ , τ ).
jk 1 2
We assume the nonlinear material to have inversion symmetry as it is found in many crystals
as e. g. Si, NaCl or in all amorphous materials, liquids or gases.
We introduce a formal inversion operator Jˆ which does not act on material parameters of
inversion symmetric materials JˆR i(2) (τ , τ ) = R i(2)
jk 1 2
(τ , τ ) but performs a sign flip of vec-
jk 1 2
tor components JˆE i = −E i and JˆP i = −P i . We now apply the inversion operator to equa-
tion 1.6

ˆ∞ ˆ∞
 
3
Jˆ[P i(2) (t )] = Jˆε0 dτ2 R i(2) (τ , τ )E j (t − τ1 )E k (t − τ2 )
X
dτ1 jk 1 2
i , j =1
0 0
3
ˆ∞ ˆ∞
= ε0 dτ2 JˆR i(2) (τ , τ ) JˆE j (t − τ1 ) JˆE k (t − τ2 )
X
dτ1 jk 1 2
i , j =1
0 0
3
ˆ∞ ˆ∞
−P i(2) (t ) = ε0 dτ2 R i(2) (τ , τ )(−1)E j (t − τ1 )(−1)E k (t − τ2 )]
X
dτ1 jk 1 2
j ,k=1
0 0
−P i(2) (t ) = P i(2) (t ). (1.7)

This implies that P i(2) (t ) must be equal to zero for all times t which can only be fulfilled with
R i(2)
jk
= 0. In a similar way we can prove that all even order nonlinearities must vanish in cen-
trosymmetric and amorphous materials.

1.3 Polarization in frequency space

In optics we usually deal with monochromatic fields and are more interested in the spectral
amplitudes than in the fast temporal evolution. Often it is required to quantify the nonlinear
response on certain frequency components of the electric field.
a.) We start with the decomposition of the optical field into its frequency components via
Fourier transform which are defined in this lecture as follows
ˆ∞ ˆ∞
1
E (t ) = dω Ẽ (ω)e −iωt
and Ẽ (ω) = dt E (t )eiωt . (1.8)

−∞ −∞

The temporal field E (t ) is thus a sum of oscillating fields each having a spectral amplitude
Ẽ (ω). Note that E (t ) is real valued which implies that Ẽ (ω) is symmetric

ˆ∞ ˆ∞
 ∗
1 1
Ẽ (−ω) = dt E (t )e−iωt = dt E (t )eiωt  = Ẽ ∗ (ω). (1.9)
2π 2π
−∞ −∞

5
1.3 Polarization in frequency space NLO

This means that the spectral amplitudes at negative frequencies Ẽ (−|ω|) are nonzero, but
strictly correlated to their positive counterparts and thus need not be considered separately.
b.) Now we can find an expression for P̃ (1) (ω) by inserting (1.8) into (1.5)

3
ˆ∞
P i(1) (t ) = ε0 dτ1 R i(1) (τ1 )E j (t − τ1 )
X
j
j =1
0
3
ˆ∞ ˆ∞
= ε0 dτ1 R i(1) dω Ẽ j (ω)e−iω(t −τ1 )
X
j
(τ1 )
j =1
0 −∞
ˆ∞ 3
ˆ∞ ˆ∞
dω ε0 dτ1 R i(1) (τ1 )eiωτ1 Ẽ j (ω)e−iωt = dω P̃ i(1) (ω)e−iωt .
X
= j
(1.10)
j =1
−∞ 0 −∞

P̃ i(1) (ω) is the spectral amplitude of the linear polarization. It depends on fields oscillating
with the same frequency in a linear way
ˆ
 ∞ 
3 3
P̃ i(1) (ω) = ε0  dτ1 R (1) (τ1 )eiωτ1 Ẽ j (ω) = ε0 χ(1)
X X
ij ij
(−ω|ω)Ẽ j (ω). (1.11)
j =1 j =1
0

Here we introduce the complex valued linear susceptibility χ(1)


ij
(−ω|ω) which is the Fourier
transform of the real valued response function
ˆ∞
χ(1)
ij
(−ω|ω) = dτ1 R i(1)
j
(τ1 )eiωτ1 . (1.12)
0

The two arguments of χ(1) ij


(−ω|ω) are introduced for consistency with the following nota-
tion of complex coefficients. The first argument represents the negative of the frequency of
the induced polarization, where the second argument stands for the frequency of the elec-
tric field. Of course, in the linear case the frequency of the induced polarization and of the
electric field coincide and therefore the first argument is usually omitted in the linear case.
However, as we will see below such a simplification is no longer possible in the nonlinear
case.
c.) Similarly as before we can perform a spectral decomposition of the quadratic polarization
P (2) (t ) by inserting (1.8) into (1.6)

3
ˆ∞ ˆ∞
P i(2) (t ) = ε0 dτ2 R i(2) (τ , τ )E j (t − τ1 )E k (t − τ2 )
X
dτ1 jk 1 2
(1.13)
i , j =1
0 0
3
ˆ∞ ˆ∞ ˆ∞ ˆ∞
= ε0 dτ1 dτ2 R i(2) (τ , τ ) −iω1 (t −τ1 )
dω2 Ẽ k (ω2 )e−iω2 (t −τ2 )
X
jk 1 2
dω1 Ẽ j (ω1 )e
i , j =1
0 0 −∞ −∞
ˆ∞ ˆ∞ 3
ˆ ∞ ˆ∞

dω2 ε0 dτ2 R i(2) (τ , τ )ei(ω1 τ1 +ω2 τ2 ) Ẽ j (ω1 )Ẽ k (ω2 )e−i(ω1 +ω2 )t .
X
= dω1 dτ1 jk 1 2
i , j =1
−∞ −∞ 0 0

6
1.3 Polarization in frequency space NLO

Now we apply the Fourier transform to (1.13)

ˆ∞
1
P̃ i(2) (ωp ) = dt eiωp t P i(2) (t )

−∞
ˆ∞ ˆ∞ ˆ∞
1 h
= dt eiωp t dω1 dω2

−∞ −∞ −∞
3
ˆ∞ ˆ∞ i
ε0 dτ2 R i(2) (τ , τ )ei(ω1 τ1 +ω2 τ2 ) Ẽ j (ω1 )Ẽ k (ω2 )e−i(ω1 +ω2 )t
X
dτ1 jk 1 2
i , j =1
0 0
ˆ∞
1
= dt ei(ωp −ω1 −ω2 )t ...

−∞
ˆ∞ ˆ∞ ˆ∞ ˆ∞
 
3
dω2 ε0 dτ2 R i(2) (τ , τ )ei(ω1 τ1 +ω2 τ2 ) Ẽ j (ω1 )Ẽ k (ω2 ).
X
dω1 dτ1 jk 1 2
i , j =1
−∞ −∞ 0 0
(1.14)

We make use of the definition of the delta-distribution


ˆ∞ ˆ∞
1 i(ωp −ω1 −ω2 )t
dt e = δ(ωp − ω1 − ω2 ) with dx δ(x) f (x) = f (x = 0). (1.15)

−∞ −∞

Inserting the delta-distribution and integrating results in

ˆ∞ ˆ∞ 3
P̃ i(2) (ωp ) = χ(2) (−ωp |ω1 , ω2 )Ẽ j (ω1 )Ẽ k (ω2 )δ(ωp − ω1 − ω2 )
X
dω1 dω2 i jk
j ,k=1
−∞ −∞
ˆ∞ 3
dω1 ε0 χ(2) (−ωp |ω1 , ωp − ω1 )Ẽ j (ω1 )Ẽ k (ωp − ω1 )
X
= i jk
(1.16)
j ,k=1
−∞

with the nonlinear susceptibility of second order

ˆ∞ ˆ∞
χ(2)
i jk
(−ωp |ω1 , ω2 ) = dτ1 dτ2 R i(2) (τ , τ )ei(ω1 τ1 +ω2 τ2 ) .
jk 1 2
(1.17)
0 0

Several comments are in order:


• Nonlinearities generate new frequencies in the induce polarization, which are con-
verted to optical fields with new frequency components (new colors) by emission.
• The delta distribution δ(ωp −ω1 −ω2 ), which we have finally integrated, ensures energy
conservation.
• A photon of frequency ωp can only be created, if its frequency is equal to the sum of
the frequencies of the initial photons (including negative frequencies).

7
1.3 Polarization in frequency space NLO

• Because the spectral amplitudes of the electric field can also have arguments with neg-
ative frequencies, equation (1.16) incorporates difference as well as sum frequency
generation.
d.) Higher order nonlinear polarizations are introduced in the same way as to derive equa-
tion (1.16)
ˆ∞ ˆ∞ 3 3 h
dωn ε0 χ(n)
αp α1 ... αn (−ωp |ω1 , ..., ωn )Ẽ α1 (ω1 )
X X
P̃ αp (n)(ωp ) = dω1 ... ...
α1 =1 αn =1
−∞ −∞
n i
... Ẽ αn (ωn )δ ωp − ωk )
¡ X
(1.18)
k=1

with the nonlinear susceptibility


ˆ∞ ˆ∞
χ(n)
αp α1 ... αn (−ωp |ω1 , ..., ωn ) = dτ1 ... dτn R α(2)p α1 ... αn (τ1 , ..., τn )ei(ω1 τ1 ) ...ei(ωn τn ) , (1.19)
0 0
Pn
and ωp = i =1 ωi ensuring energy conservation.
Of course one of the frequencies found in the integral in equation (1.18) can be removed by
performing the integration and evaluating the delta-distribution. We also want to list some
properties of χ(n)
αp α1 ... αn (−ωp |ω1 , ..., ωn ):

• Only defined/different from zero for ωp = ω1 + ... + ωn


• Arguments originating from the driving fields (frequencies and indices of polarization)
can be exchanged together:
χ(2) (2)
αp α1 ... αk ... αl ... αn (−ωp |ω1 ...ωk ...ωl ...ωn ) = χαp α1 ... αl ... αk ... αn (−ωp |ω1 ...ωl ...ωk ...ωn )

• The susceptibility carries the symmetries of the respective material. The higher the
symmetry of the material the lower the number of independent coefficients.
• A K RAMERS -K RONIG transformation can be performed with respect to one of the fre-
quencies. It links real and imaginary parts of the coefficients as
ˆ∞ χ(n)
1 αp α1 ... αn (−ωp (ϖ)|ω1 , ...ϖ...ωn )
χ(n)
αp α1 ... αn (−ωp (ω0 )|ω1 , ...ω0 ...ωn ) = P dϖ
iπ ϖ − ω0
−∞

with P being Cauchy’s principal value integral.


e.) Furthermore we want to take a look at the nonlinear interaction of several continuous
waves (cw) which are discrete frequency components instead of a continuous spectrum
1 X¡
E λ e−iωλ t + c.c. .
¢
E (t ) = (1.20)
2 λ

The quantity |E λ | is the peak amplitude of the electric field component with frequency ωλ .
Fourier spectrum of sum of cw-waves is equal to a sum of delta functions
1 X£
E λ δ(ω − ωλ ) + E λ∗ δ(ω + ωλ ) .
¤
Ẽ (ω) = (1.21)
2 λ

8
1.3 Polarization in frequency space NLO

The nonlinear polarization of n-th order induced by cw-waves is a sum of harmonic oscilla-
tions
1 X ³ (n) −iωp t ´
P (n) (t ) = P ωp e + c.c. . (1.22)
2 ωp

We now want to determine the amplitudes of harmonic oscillations of the nonlinear polar-
ization of n-th order induced by cw-waves P ω(n)
p
by following these steps:
1. Formal integration of equation (1.18) using the Fourier spectrum from (1.21).
2. Inserting the nonlinear coefficient χ(n)
αp α1 ... αn (−ωp |ω1 , ..., ωn ) defined in equation (1.19).

3. Evaluating the Fourier spectrum of nonlinear polarization of n-th order induced by


cw-waves
1 X ³ (n) ´
P̃ (n) (ω) = P ωp δ(ω − ωp ) + P ω(n)∗
p
δ(ω + ω p ) . (1.23)
2 ωp

4. Evaluating permutation symmetry (incorporated into the K -coefficient).


The spectral amplitude of the nonlinear polarization can now be written as

3
(∗)
P ω(n) = ε0 K (−ωp | ± ω1 , ... ± ωn ) χ(n) (∗)
αp α1 ... αn (−ωp | ± ω1 , ... ± ωn )E 1,α1 E n,αn
X X
p αp
(1.24)
{±ωλ } α1 ... αn =1

• ±ωλ : one of the frequencies of the interacting cw-waves, either positive or negative
(∗)
• E λ,α : vector component αλ of the field amplitude E λ of the continuous wave oscillat-
λ
ing with ωλ . In case of the negative frequency the complex conjugate field profile E λ,α

λ
is chosen.
• {±ωλ }: all combinations of interacting ±ωλ including negative frequencies, which gen-
erate ωp , but no permutations.
The factor K (−ωp |±ω1 , ...±ωn ) incorporates the effect of permutations. It includes all factors
1
2 used to express the amplitudes of the electric fields and the polarization. We can generalize
it in the following way:

K (−ωp | ± ω1 , ... ± ωn ) = 2l +m−n p. (1.25)

• n: order of the nonlinear interaction


• p: number of permutations of the ±ωλ (p is equal to n! if all frequencies are different.)
1
• m: number of ωλ being zero (For those the factor 2
vanishes in the above expres-
sion (1.20) for the fields.)
• l : equal to 1 provided that ωp is not zero, otherwise l = 0

9
1.4 Symmetries and conserved quantities NLO

1.4 Symmetries and conserved quantities

a.) Complete permutation symmetry


We want to extend the permutation symmetry of the electric fields represented by ωλ to-
wards the induced polarization:

χαp α1 ... αk ... αn (−ωp |ω1 . . . ωk ωn ) = χαk α1 ... αp ... αn (ωk |ω1 . . . −ωp ωn ). (1.26)

This assumption is justified far away from resonance which means that nonlinear coeffi-
cients can be expressed by a Drude model

1 n 1
χ(n)
Y
αp α1 ... αn (−ωp |ω1 ...ωn ) ∼ , (1.27)
ω20 − ω2p λ=1 ω20 − ω2λ

where ω0 is the resonance frequency of the material. The frequency of the polarization does
not play an extraordinary role and can be mixed with the other components.

b.) Kleinmann symmetry (no frequency dependence)


The preconditions for the Kleinmann symmetry to hold are that we deal with frequencies in
the transparent domain of the material which are far away from resonances (ω20 ≫ ω2p and
ω20 ≫ ω2λ ). Then all frequency dependence can be neglected. Often Kleinmann symmetry is
assumed because of the lack of more detailed knowledge.
For quadratic nonlinearities we can introduce a matrix notation based on field amplitudes
of cw-fields defined in (1.20) and (1.21)
 
E x (ω1 )E x (ω2 )
 (2)   E y (ω1 )E y (ω2 ) 
Px

d 11 . . . d 16  
E (ω )E (ω )
 
 (2)  z 1 z 2
P y  = 2ε0 K (−ωp |ω1 , ω2 )d 21 . . . d 26  . (1.28)
 
(2)  E y (ω1 )E z (ω2 ) + E z (ω1 )E x (ω2 ) 
Pz d 31 . . . d 36  
 E y (ω1 )E z (ω2 ) + E z (ω1 )E x (ω2 ) 
E x (ω1 )E y (ω2 ) + E y (ω1 )E x (ω2 )

The d -coefficients are a very common representation of quadratic nonlinearities. They are
orientation dependent (x- and y-axis usually aligned parallel to the crystal axis). We can also
write down the K -factor for all cases possible:
• All frequencies different and nonzero: K = 1
• ω1 = ω2 or ω1 = −ω2 : K = 12
• ω1 = 0 or ω2 = 0: K = 2

10
1.5 Cubic nonlinearities in isotropic materials NLO

c.) Conservation laws


For arbitrary orders of nonlinearities we can find conservation laws which are valid in the
absence of dissipation or losses. We note Wλ as the energy flux of a field component of fre-
quency ωλ which points into a certain direction. Energy conservation means that λ Wλ =const.
P

For non-degenerate processes |ωα | ̸= |ωβ | and arbitrary α and β including the induced po-
larization we have photon flux conservation

Wα Wβ
− = const. (1.29)
ωα ωβ

To successfully create a photon at ωp a photon of each of the interacting waves must be


absorbed (photon number difference remains conserved).

1.5 Cubic nonlinearities in isotropic materials

Cubic nonlinearities can usually be neglected if quadratic nonlinearities are present. There-
fore they are only relevant if quadratic coefficients vanish (e. g. in inversion symmetric ma-
terials as glass) or if quadratically nonlinear processes are suppressed (e. g. due to the lack of
phase matching). The susceptibility tensor

χ(3)
i j kl
(−(ω1 + ω2 + ω3 )|ω1 , ω2 , ω3 ) (1.30)

possesses 34 = 81 elements.
From now on will will only consider isotropic materials with the highest possible degree of
symmetry, which means they have the minimum number of independent elements.

a.) All coefficients are equal j = k = l


These are coefficients such as χ(3) (3)
xxxx or χ y xxx . We assume that E points into the direction of
equal coefficients

E = (E j + E k + E l )êk .

Thus, the induced polarization must point into the same direction as E because the material
has no internal direction. Since E ||P we have

χ(3)
i kkk
= 0 for i ̸= k and χ(3) (3) (3)
xxxx = χ y y y y = χzzzz . (1.31)

The three non-vanishing elements must be equal as no direction is preferred. Thus we only
have one independent element.

11
1.5 Cubic nonlinearities in isotropic materials NLO

b.) One coefficient is different j = k ̸= l


Here we assume that E points into the two relevant directions

E = (E j + E k )êk + E l êl .

The vector components of P (3) at a fixed frequency are P i(3) ∼ χ(3) E E E . We assume that
i kkl j k l
the material is rotated around the singular axis (l ) by 180°. In the new coordinate system the
electric field is

E ′ = −(E j + E k )ê′k + E l ê′l .

As the material is isotropic the nonlinear coefficient has not changed by rotation χ(3)
i kkl
=

χ(3)
i kkl
. The vector components of P (3) in the rotated coordinate system are

P i(3) ∼ χ(3)
i kkl
(−E j )(−E k )E l = P i(3) . (1.32)

This is only possible if P (3) points along the axis of rotation. It follows that χ(3)
i kkl
= 0 for i ̸= l .
Thus we have 18 non-vanishing coefficients (all permutations are taken into account, for
different frequencies permutations give different coefficients)

χ(3) (3)
x y y x (−(ω1 + ω2 + ω3 )|ω1 , ω2 , ω3 ) ̸= χxx y y (−(ω1 + ω2 + ω3 )|ω1 , ω2 , ω3 ). (1.33)

However, we only have three independent coefficients e. g.

χ(3) (3) (3) (3) (3) (3) (3)


x y y x = χxzzx = χ y xx y = χ y zz y = χzxxz = χz y y z ̸= χxx y y . (1.34)

c.) All indices are different j ̸= k ̸= l ̸= j


With the same arguments as above we can perform a 180° rotation of the material and the
coordinate system around e. g. the j -direction and will find out that P (3) must point into that
direction. This we can do with every of direction as all are represented by a single field. As
P (3) can only point into a single direction it must vanish for that condition. All coefficients
are zero.
Concluding our results a general χ(3) tensor has 21 non-zero coefficients and only four in-
dependent coefficients. In the following we will introduce degeneracies in the interaction
waves leading to Third harmonic generation and self-phase modulation.

1.5.1 Third harmonic generation (THG)

For this process we consider the case j = k ̸= l = i with the nonlinear coefficient

χ(3)
i j kl
(−3ω, ω, ω, ω).

Since the arguments are equal the frequencies can be exchanged freely
(3)
χ(3) (3) (3)
xx y y = χx y x y = χx y y x = ... = χi i kk . (1.35)

12
1.5 Cubic nonlinearities in isotropic materials NLO

In total we have only two independent coefficients χ(3)


i i kk
and χ(3)
iiii
. In case of THG they are
not really independent which will be shown in the following.
We assume that we can rotate the material in such a way that the electric field E points into
the y- and z-direction
 
0
E  
E=p 1 .
2 1

The resulting components of the polarization are


E 2 ³ (3) ´
P y(3) = ε0 K (−3ω|ω, ω, ω) χ y zz y + χ(3)
yzyz + χ (3)
y y zz + χ (3)
yyyy Ey
2 (1.36)
(3) E 2 ³ (3) (3) (3) (3)
´
P z = ε0 K (−3ω|ω, ω, ω) χ + χz y z y + χzz y y + χzzzz E z .
2 zy yz
For THG there are no frequency permutations p = 1 which implies K = 1/4 (see equation (1.25)).
Then we have
E2
P (3) = ε0 K (−3ω|ω, ω, ω) (3χ(3) (3)
y zz y + χ y y y y )E . (1.37)
2
The polarization is parallel to E which is no surprise, as the material is isotropic.
Now we assume that the material (and the coordinate system) is rotated in such a way that
E points into the y-direction. This implies

P (3) = ε0 K (−3ω|ω, ω, ω)E 2 χ(3)


yyyyE. (1.38)

The relation between P (3) and E must be independent of the rotation. Therefore we can
equate (1.37) and (1.38)
3 (3) 1
χ y zz y + χ(3) = χ(3)
yyyy. (1.39)
2 2 yyyy
This means in an isotropic material THG is determined by a single independent coefficient
(3) (3)
χ(3) (3)
xx y y = χx y x y = ... = χTHG and χ(3) (3) (3)
xxxx = χ y y y y = χzzzz = 3χTHG . (1.40)

We can generalize the total third harmonic polarization by first looking at a single compo-
nent
1
P y(3) = ε0 (χ(3) 3 (3) 2 (3) 2
y y y y E y + 3χ y zz y E y E z + 3χ y xx y E y E x )
4
3 ³ ´
= ε0 χ(3) E 2
+ E 2
+ E 2
x Ey
4 THG y z

3
= ε0 χ(3) E 2E y . (1.41)
4 THG
Thus the total third harmonic polarization is

3
P (3) = ε0 χ(3) E 2E . (1.42)
4 THG
The third harmonic generated in an isotropic medium:

13
1.5 Cubic nonlinearities in isotropic materials NLO

• points into the direction of the exciting field


• vanishes for circularly polarized light.
2 2
For left circularly polarized light: E = pE (êx + iê y ) which means E 2 = E2 − E2 = 0.
2

We want to give an explanation for this phenomenon. In an isotropic medium the direction
of polarization is entirely determined by the direction of the exciting field. In case of circu-
larly polarized electric field the polarization can only rotate with the rotation speed of the
fundamental harmonic, hence one rotation per optical cycle of the fundamental harmonic.
Therefore, only polarization at fundamental harmonic frequency is generated.

1.5.2 Self phase modulation (SPM) or optical Kerr effect

Here we have the following nonlinear interaction: χ(3) i j kl


(−ω|ω, ω, −ω).
The negative frequency −ω corresponds to the complex conjugate of the field E ∗ . The non-
linear polarization and the driving field coincide in frequency. We observe a nonlinear phase
modulation Re χ(3)i j kl
and losses Im χ(3)
i j kl
of the initial wave only. We have more independent
coefficients than for THG, e. g.

χ(3) (3)
xx y y (−ω|ω, ω, −ω) ̸= χx y y x (−ω|ω, ω, −ω, ). (1.43)

The third order nonlinear polarization has the following form:


h ³ ´ i
P y(3) = ε0 |{z}
K 2 χ(3)y y xx |E x | 2
E y + χ (3)
y y zz |E z | 2
E y + χ (3)
yyyy |E y | 2
E y + χ (3)
E 2 ∗
E
y xx y x y + χ (3)
E 2 ∗
E
y zz y z y .
3/4
(1.44)

For SPM we have three frequency permutations (p = 3) and therefore K = 3/4 (see equa-
tion (1.25)). We use the same arguments as for THG

1 (3)
χ(3) (3) (3)
y y y y = (χ y zz y + 2χ y y zz + χ y y y y ) ⇒ χ(3) (3) (3)
y y y y = χ y zz y + 2χ y y zz . (1.45)
2

h i
P y(3) = ε0 K 2χ(3)
y y zz (|E x | 2
+ |E y | 2
+ |E z | 2
)E y + χ (3)
y zz y |E y | 2
E y + χ (3)
y zz y (E 2
x + E 2 ∗
z )E y
h i
= ε0 K 2χ(3) 2 2 2 (3) 2 2
y y zz (|E x | + |E y | + |E z | )E y + χ y zz y (E x + E y + E z )E y .
2 ∗
(1.46)

Neither x, y nor z direction are preferred:

3 ³ ´
P (3) (ω) = ε0 2χ(3) 2 (3) 2 ∗
y y zz |E | E + χ y zz y E E . (1.47)
4
We find two independent coefficients and want two discuss both nonlinear actions:
• First term: |E |2 E
depends on the intensity only and is a pure phase modulation or a nonlinear loss

14
1.5 Cubic nonlinearities in isotropic materials NLO

• Second term: E 2 E ∗
is sensitive to the state of polarization of the electric field and completely vanishes for
circularly polarized light as E 2 (see end of section 1.5.1). It bears half of the contribu-
tion as the first term for linearly polarized light.
We can further simply the result for K LEINMANN symmetry (no frequency dependence):
χ(3) (3)
y zz y = χ y y zz . Only a single coefficient is left and we obtain

3
P (3) = ε0 χ(3) 2 2 ∗
¡ ¢
y y zz 2|E | E + E E . (1.48)
4

Optical Kerr effect


Now we take linearly polarized light (e. g. E = E êx )

9 9 (3)
P (3) = ε0 χ(3) 2 2
y y zz |E | E = ε0 ∆χ with ∆χ = χ y y zz |E | . (1.49)
4 4
Thep
action of the nonlinear susceptibility can be expressed by a change of the refractive index
n = 1 + χ as
q ∆χ 9 (3)
∆n = n 02 + ∆χ − n 0 ≈ = χ |E |2 , (1.50)
2n 8n 0 y y zz

where n 0 is the refractive index in the absence of a nonlinear action. As there is no gen-
erally accepted definition of the amplitude of an oscillating field one prefers to refer to the
intensity
ε0
I= cn 0 |E |2 . (1.51)
2
We can replace the electric field in (1.50) and find
(3)
9 χ y y zz
∆n = n 2 I with n 2 = . (1.52)
4 ε0 cn 02

We want to give some examples of nonlinear refractive indices:

m2
silica n 2 = 2,3 · 10−20 (1.53)
W
2
m
Al0.18 Ga0.82 As(1,5 µm) n 2 = 2,1 · 10−17 . (1.54)
W

15
2 Resonant nonlinearities - a quantum mechanical
representation
In what follows we will shortly discuss the kind of nonlinear response, which shows up close
to a resonance of the material. If the energy of the incident photons just roughly matches
the energy of a quantum mechanical transition, the optical properties of the excited mate-
rial may change dramatically as a function of the incident intensity. A Taylor expansion of
the induced polarization with respect to powers of the acting electrical field as done in equa-
tion (1.4) is still applicable, however, the treatment applied in the last section is no longer
useful. The expansion can usually not be restricted to a few orders and resulting coefficients
turn out to be highly dispersive. Hence, a treatment taking into account the material re-
sponse more seriously has to be applied.

2.1 Equations of motion of a two-level-system

We first want to revise what a two-level-system actually Eb ¯ψ b


¯ ®
b
is: It is a quantum system with (at least) two energy
¯ ® lev-
els a and b corresponding to distinct states ¯ψa (ex-
cited state) and ¯ψb (ground state) with energies E a
¯ ®

and E b (E a > E b ). The energy of photons exciting the a Ea ¯ψa


¯ ®
ground state is ℏω0 ≈ E a − E b .
Real quantum states have many (usually infinite) possible states. A two-level-system is sin-
gled out by the excitation. The dynamics of the system are determined by the Hamiltonian
Ĥ for a single electron in an electromagnetic field

1
Ĥ = (p̂ + e A real )2 + V (r ). (2.1)
2m
The term V (r ) is the potential, which localizes the electrons. A real is the vector potential,
which represents the optical field by


E real = − A real . (2.2)
∂t
The subscript real is used to emphasize that observable fields are real valued and that they
contain components oscillating at +ω0 and −ω0 . We consider the vector potential as classical
which means that
• The strong field contains many photons
• quantum fluctuations can be neglected
• The field commutes with all operators
2
• The vector potential is weak compared to atomic fields ⇒ A real is neglected.

16
2.1 Equations of motion of a two-level-system NLO

We can now decompose the Hamiltonian into an unperturbed part Ĥ0 and a perturbation

1 2 e e 2 2

Ĥ = p̂ + V (r ) + p̂ A real + A . (2.3)
|2m {z } |m {z }  2m real
Ĥ0 Ŵ

The unperturbed part Ĥ0 is the Hamiltonian of the unperturbed two-level-system in absence
of any electromagnetic ¯ field.
® It defines ® energy levels E a and follows the S­CHRÖDINGER
the
equation Ĥ0 ψa = E a ψa , where ψa/b is a orthonormal basis of eigenstates ψa |ψb ¯ψa |ψb =
¯ ® ¯ ¯ ®
¯ ¯ ¯
δab . The evolution is restricted to the two eigenstates, which means that the system is either
in state a or b. The wave function of the whole system can be expressed as

¯ψ(t ) = a(t ) ¯ψa + b(t ) ¯ψb with |a(t )|2 + |b(t )|2 = 1,
¯ ® ¯ ® ¯ ®
(2.4)

where the second conditions follows from ψ(t )|ψ(t )¯ψ(t )|ψ(t ) = 1.
­ ¯ ®

The interaction operator Ŵ drives the nontrivial dynamics and induces changes. Using the
time dependent S CHRÖDINGER equation leads to

d ¯¯
ψ(t ) = Ĥ ¯ψ(t )
® ¯ ®
iℏ
dt

a(t ) ¯ψa + b(t ) ¯ψb = Ĥ (a(t ) ¯ψa + b(t ) ¯ψb )
¯ ® ¯ ®¢ ¯ ® ¯ ®
iℏ
dt
¯ ® d ¯ ® d
i ℏ ¯ψ a a(t ) + iℏ ¯ψb b(t ) = a(t ) Ĥ ¯ψa + b(t ) Ĥ ¯ψb .
¯ ® ¯ ®
(2.5)
dt dt

Multiplying from the left with ψa/b ¯ and using ψa |ψb ¯ψa |ψb = δab leads to
­ ¯ ­ ¯ ®

d
a = E a a + ψa |Ŵ |ψa a + ψa |Ŵ |ψb b
­ ® ­ ®
iℏ
dt
d (2.6)
iℏ b = E b b + ψb |Ŵ |ψb b + ψb |Ŵ |ψa a .
­ ® ­ ®
dt | {z }
interaction elements

Matrix elements of the interaction operator


We first calculate the symmetric elements:
® D e E e
ψa |Ŵ |ψa = ψa | p̂ A real |ψa = A real ψa |p̂|ψa
­ ­ ®
m ˆ m
e ℏ
= A real dV ψ∗a (r )⃗
∇ ψa (r )
m i
ˆ
e ℏ
dV ⃗
¢2
∇ ψa (r ) = 0.
¡
= A real (2.7)
m 2i

The last integral just gives us the values at the boundary. For a decaying wave functions the
boundaries are zero, leading to vanishing interaction terms. In the last line we assume a real
valued ψa (r ).

17
2.1 Equations of motion of a two-level-system NLO

®∗
The second matrix element ψa |Ŵ |ψb = ψb |Ŵ |ψa is responsible for transitions between
­ ® ­

the two states. It is usually nonzero , provided that an external field is present and both states
have different parity. For the evaluation we need a helpful commutator relation
· ¸ · ¸
1 2 1 2
[r̂ , Ĥ0 ] = r̂ p̂ + V (r ) − p̂ + V (r ) r̂
2m 2m
1 1 iℏ
= (r̂ p̂ 2 − p̂ 2 r̂ ) = ([r̂ , p̂]p̂ + p̂[r̂ , p̂]) = p̂. (2.8)
2m 2m m
Now we can calculate the matrix element as
e ® e D m E
ψa |Ŵ |ψb = A real ψa |p̂|ψb = A real ψa | (r̂ Ĥ0 − Ĥ0 r̂ )|ψb
­ ® ­
m m iℏ
ie
= A real (E a − E b ) ψa |r̂ |ψb .
­ ®
(2.9)

The optical field couples to the states via the dipole moment of the transition
ˆ
ψa |r̂ |ψb = dV ψ∗a (r )r ψb (r ).
­ ®
(2.10)

In the following we want to transform a and b into macroscopic observables:

a.) Inversion I = |a|2 − |b|2


The inversion describes the state of excitation of the two-level-system. It can vary between
−1 (no excitation) and 1 (maximum excitation). For its evolution we use equation (2.6)

dI d¡ 2 ¢ (2.4) d ¡ d
|a| − |b|2 = 2|a|2 − 1 = 2 |a|2
¢
=
dt dt dt dt
d 2
= 2a ∗ a + c.c. = ψa |Ŵ |ψb a ∗ b + c.c.
­ ®
dt iℏ
2e
= 2 A real (E a − E b ) ψa |r̂ |ψb a ∗ b + c.c.
­ ®
(2.11)

Now we want to do a frequency analysis of ab ∗ for a weak perturbation Ŵ and a resonant
excitation (E a − E b )/ℏ ≈ ω0 . If the energy of the state is much larger than the interaction
energy the solutions of a(t ) and b(t ) using (2.6) are
µ ¶ µ ¶
iE a iE b
a(t ) ∼ exp − , b(t ) ∼ exp − . (2.12)
ℏt ℏt
This leads to the following expression for ab ∗

E a − Eb
µ ¶

ab ∼ exp −i t ≈ exp(−iω0 t ), (2.13)

which means that ab ∗ is related to the optical polarization.

18
2.1 Equations of motion of a two-level-system NLO

b.) Polarization
The classical change of inversion is a result of absorbed optical power
dI absorbed power per system 2 1 ∂
=2 = E real P real . (2.14)
dt energy required to excite a system E a − E b N ∂t
This can be understood by noting that power is Voltage times current, where E real corre-

sponds to the voltage and ∂t P real to current. The latter one can be understood that the polar-
ization describes a deviation of the electrons from their rest position in the lattice. A change
of polarization describes the electron oscillations around their rest position and thus a cur-

rent. Then E real ∂t P real is the absorbed power per volume. It contains fast parts (oscillations
with ≈ 2ω0 ), however, they average out or can be neglected. In the following we only consider
the slow/almost stationary parts.
We decompose the fields into fast parts oscillating with the carrier frequency of the optical
field ω0 and slowly varying amplitudes

P (t ) exp(−iω0 t ) + P ∗ (t ) exp(iω0 t )
¤
P real (t ) =
2

E real (t ) = E (t ) exp(−iω0 t ) + E ∗ (t ) exp(iω0 t )
¤
(2.15)
2

A real (t ) = A(t ) exp(−iω0 t ) + A ∗ (t ) exp(iω0 t ) .
¤
2
The amplitudes are assumed to be slow that their time derivatives are negligible compared
with the derivatives of the fast-varying exponents
¯∂ ¯ (2.2) 1
|ω0 A(t )| ≫ ¯ A(t )¯ ⇒ A(t ) ≈ E (t ) (2.16)
∂t iω0
¯∂ ∂ iω0 £
P (t ) exp(−iω0 t ) − P ∗ (t ) exp(iω0 t ) .
¯ ¤
|ω0 P (t )| ≫ ¯ P (t )¯ ⇒ P real (t ) ≈ − (2.17)
∂t ∂t 2
Then we can formulate the classical evolution of inversion
dI (2.14) 2 1 ∂
= E real P real
dt E a − Eb N ∂t

2 1 E (t ) −iω0
= P (t ) + c.c. + fast oscillations
E a − Eb N 2 2
ω
= E ∗ (t )P (t ) + c.c. + fast oscillations. (2.18)
2i(E a − E b )N
Similarly we can formulate the quantum evolution of inversion
dI (2.11) 2e
A real (E a − E b ) ψa |r̂ |ψb a ∗ b + c.c.
­ ®
=
dt ℏ2
(2.13) e i
(E a − E b ) ψb |r̂ |ψa E ∗ (t ) exp(iω0 t )ab ∗ + c.c. + fast oscillations.
­ ®
= (2.19)
(2.15) ℏ ω0
2

If we now compare the slow parts containing E ∗ (t ) in the classical and quantum formula-
tions
e i ω0
(E a − E b ) ψb |r̂ |ψa E ∗ (t ) exp(iω0 t )ab ∗ ≈ E ∗ (t )P (t ).
­ ®
(2.20)
ℏ ω0
2 2i(E a − E b )N

19
2.1 Equations of motion of a two-level-system NLO

The slowly varying envelope of the optical polarization is then

E a − Eb 2 ­
µ ¶
ψb |r̂ |ψa ab ∗ exp(iω0 t ).
®
P = −2eN (2.21)
ℏω0
| {z }
≈1

Here we assume a resonant excitation which means that the frequency detuning ∆ω is small
E a − Eb
∆ω = − ω ≪ ω0 . (2.22)

The total or real valued polarization then is

P (t ) exp(−iω0 t ) + c.c. = −eN ψb |r̂ |ψa ab ∗ + c.c.
¤ ­ ®
P real = (2.23)
2
Now we also want to discuss the dynamics of the nonlinear polarization by using (2.6)

db ∗
µ ¶ µ ¶
d ¡ ∗¢ da ∗
iℏ ab = iℏ b − a iℏ
dt dt dt
= (E a − E b )ab − ψa |Ŵ |ψb |a|2 − |b|2

­ ®¡ ¢

ie
= (E a − E b )ab ∗ − A real (E a − E b ) ψa |r̂ |ψb I
­ ®

1 £ ¤ ie
≈ (E a − E b )ab ∗ − E (t )e−iω0 t − E ∗ (t )eiω0 t (E a − E b ) ψa |r̂ |ψb I
­ ®
2iω0 ℏ
e(E a − E b )
= (E a − E b )ab ∗ − E (t )e−iω0 t − E ∗ (t )eiω0 t ψa |r̂ |ψb I .
£ ¤­ ®
(2.24)
2ℏω0

The driving term of ab ∗ with the wrong frequency (∼ e+iω0 t ) is neglected. The evolution of
the polarization can be derived by using (2.24), (2.21) and (2.22)

dP Ne 2 ­
= ∆ωP + I ψb |r̂ |ψa ( ψa |r̂ |ψb E ).
® ­ ®
i (2.25)
dt ℏ
For an isotropic medium with linear polarization of the electric field, the polarization P is
parallel to E which leads to

ψb |r̂ |ψa ( ψa |r̂ |ψb E ) = | ψb |r̂ |ψa |2 E .


­ ® ­ ® ­ ®
(2.26)

Then the final set of equations can be summarized to

Equations of motion

d Ne 2 ­
P = ∆ωP + | ψb |r̂ |ψa |2 E I
®
i (2.27)
dt ℏ
dI 1
= E ∗ P + c.c. (2.28)
dt 2iℏN

Note that up to now the equations of motion contain no damping, because we treated the
optical field classically thus omitting spontaneous emission. Equations (2.27) and (2.28)
are strictly valid only in time domains which are small compared with classical relaxation
times.

20
2.2 Rabi oscillations NLO

Conservation law of Bloch equations and the Bloch sphere

As stated before we can describe the state of the two- a


level-system as

¯ψ(t ) = a(t ) ¯ψa + b(t ) ¯ψb .


¯ ® ¯ ® ¯ ®
(2.29) ¯ψ
¯ ®

ϑ
Furthermore we have the conservation of the norm
|a(t )|2 + |b(t )|2 = 1 and no interest in absolute phases Im b
which means we can set a to a real value. Now we may ϕ
plot the space of the two-level state as a function of
Re b
a, Re(b) and Im(b) on the surface of the Bloch-sphere
(sphere with radius 1).
We can now express the norm in terms of slowly varying
envelope of the polarization and inversion

(|a|2 + |b|2 )2 = |a|4 + |b|4 − 2|a|2 |b|2 + 4|a|2 |b|2 = (|a|2 − |b|2 )2 + 4|ab ∗ |2
1
= I2 + 2 2 ­ ® |P |2 = 1 (2.30)
e N | ψb |r̂ |ψa |2

and find this conserved quantity connecting inversion and polarization. Here we used

P = −2eN ψb |r̂ |ψa ab ∗ exp(iω0 t ).


­ ®

Inversion and the scalar component of the normalized slowly varying envelope of the polar-
ization are on the surface of the Bloch sphere.

2.2 Rabi oscillations

In what follows we solve the system of equations (2.27) and (2.28) for a linearly polarized
optical field oscillating with the frequency of the transition ∆ω = 0 with constant amplitude,
but a step-like envelope (rectangular pulse)

E (t ) = E 0 ê Θ(t ) Θ(t 1 − t ), (2.31)

where Θ(t ) is the H EAVISIDE step function. We assume that for the initial conditions (t = −∞)
the system is purely in the ground state I = −1 and P = 0.

a.) At first we look at t < 0 with E = 0. Here we find

d d
I = 0 ⇒ I = −1, P = 0 ⇒ P = 0. (2.32)
dt dt
The system remains in its ground state.

21
2.2 Rabi oscillations NLO

b.) Now we look at 0 < t < t 1 with |E | = E 0 = const. (real valued). We can decompose the
polarization into its real and imaginary parts

P = [P R (t ) + iP i (t )]ê. (2.33)

Using this we find

d
P R = 0 ⇒ P R (t ) = 0
dt
d Ne 2 ­
| ψb |r̂ |ψa |2 E 0 I
®
Pi = − (2.34)
dt ℏ
d E0
I= Pi . (2.35)
dt ℏN
Differentiating the third equation we find an equation for a harmonic motion

d2 e
I + Ω2 I = 0 Ω = E 0 | ψb |r̂ |ψa |
­ ®
Rabi frequeny. (2.36)
dt 2 ℏ

The whole system performs oscillations with the Rabi frequency Ω. The solution of the sys-
tem is then
ℏN dI ℏN
I (t ) = − cos(Ωt ), Pi = = Ω sin(Ωt ). (2.37)
E 0 dt E0

c.) The last part is t > t 1 with E = 0. Here we find again

d ℏN
P i = 0 ⇒ P i (t ) = Ω sin(Ωt 1 ) stationary envelope
dt E0
d
I = − cos(Ωt 1 ). (2.38)
dt

The real valued polarization is then

ê £
P i (t ) exp(−iω0 t ) − P i∗ (t ) exp(iω0 t ) .
¤
P real (t ) = i (2.39)
2

After the end of the excitation the polarization keeps on oscillating forever. But, any real
system is damped, an effect we have omitted by describing the optical field classically. Hence
above solution is only valid in a time range short compared with the relaxation time, which
can range from 100 fs in solids up to ms in dilute gases.
We conclude the results from the solution:
• Polarization and electric field need not to be proportional.
• The maximum strength of the polarization need not depend on the electric field strength
(see equation (2.37)).
• The polarization can be present, although the optical field is absent.

22
2.3 The influence of damping NLO

I /P
1 I

t
t1

Fig. 1: Visualization of the the inversion (blue) and polarization (orange) for an electric field present
at 0 ≤ t ≤ t 1 .

• Any kind of susceptibility becomes meaningless in the time range of Rabi-oscillations.


• A two-level-system can be brought to complete inversion by a single very short laser
pulse.
• For a correctly tailored optical pulse the system can finally return to the ground state
which means that no energy is absorbed (transparent medium). This behaviour can
be reproduced for more complex pulse shapes leading to self-induced-transparency-
solitons.

2.3 The influence of damping

Damping is always present in real world systems. It is caused by the coupling to a continuum
of modes, which could be modes of the electromagnetic radiation, phonons or plasmons. Its
exact calculation requires a quantum mechanical description of the respective field. Here
we just introduce some phenomenological constants into respective evolution equations.

Damped equations of motion

Ne 2 ­
µ ¶
d i
i P = ∆ω − | ψb |r̂ |ψa |2 E I
®
P+ (2.40)
dt T2 ℏ
dI 1 I +1
= (E ∗ P − E P ∗ ) − , (2.41)
dt 2iℏN T0

where T2 describes the relaxation time of the polarization and T0 the life time of the upper
state. Note that the relaxation term of the inversion was chosen such that I relaxes towards
the equilibrium state −1.
We can describes different scenarios using these relaxations constants:
• Free atom/atom in a trap:
Here, relaxation is due to spontaneous emission with rather long relaxation times (ms).
The relaxation is strictly correlated T2 = 2T0 .

23
2.3 The influence of damping NLO

• Dense systems:
These are gases under normal pressure, liquids or solids which will be considered in
the following. They show a strong mutual interaction with many decay channels. The
polarization is dominantly affected by so-called “dephasing”. Therefore T2 is orders of
magnitude shorter than T0 , e. g. in a direct semiconductor: T2 ∼ 100 fs, T0 ∼ ns.
We can now distinguish different cases according to the actual pulse duration T p :

a.) T p < T2 < T0 (coherent regime)


During the optical pulse we observe Rabi-oscillations. After the optical pulse the polarization
still oscillates, but decays. The system finally relaxes to the ground state.

b.) T2 < T p < T0


Here the dephasing of the polarization is faster than the actual pulse propagates through the
medium. This means the Rabi oscillations have already disappeared and the evolution of
P (t ) follows that of E (t ). We can estimate the absolute values of the terms in equation (2.40)
as
¯ ¯ ¯ ¯ ¯ ¯
¯ d ¯ ¯P ¯ ¯P ¯
¯ P ¯ ≈ ¯ ¯ ≪ ¯ ¯. (2.42)
¯ dt ¯ ¯ T ¯ ¯ T ¯
p 2

We can now neglect the time derivative of P (t ) which leads to an algebraic equation for the
polarization

Ne 2 ­
µ ¶
i
0 = ∆ω − | ψb |r̂ |ψa |2 E I
®
P+
T2 ℏ
Ne 2
| ψb |r̂ |ψa |
­ ®2
⇒ P (t ) = − ℏ I E (t ). (2.43)
∆ω − Ti2

The polarization is determined by the actual optical field. Now we can also find the effective
inversion dependent susceptibility

Ne 2 | ψb |r̂ |ψa |2
­ ®
P (t ) = ε0 χ(t )E (t ) with χ(I ) = − I. (2.44)
ε0 ℏ ∆ω − Ti2
=: χ0
| {z }

The polarization has lost its independence. It has become enslaved. We can now also de-
scribe the evolution of the inversion
dI 1 I +1 ε0 I +1
= (E ∗ P − E P ∗ ) − = Im(χ(I ))|E |2 − (2.45)
dt 2iℏN T0 ℏN T0
2 | ψ |r̂ |ψ |2
­ ®
Ne b a 1
with Im(χ(I )) = − I ∼ I.
ε0 ℏ ∆ω + 2 T2
2 1
T2

24
2.3 The influence of damping NLO

We can rewrite the differential equation as a rate equation of Im[χ(t )] by multiplying (2.45)
with the prefactor of Im[χ(t )]
d ε0 Im(χ) − Im(χ0 )
Im(χ) = − Im(χ0 ) Im(χ)|E |2 − , (2.46)
dt | ℏN {z }| T
{z0 }
(1) (2)

with χ0 being the susceptibility in the ground state. The first term (1) describes absorption
(Im χ0 > 0) induced bleaching of the transition by reducing Im χ. It switches off, when Im χ →
0. The second term (2) describes the relaxation towards the ground state.
The complete susceptibility induced by the two-level-system is then

∆ω2 + T12
χ= 2
i
T2 Im(χ). (2.47)
∆ω − T2

Finally lets look at the situation when the optical field is constant after the Rabi oscillations
have decayed (t → ∞):
• We are in a steady state with Im χ > 0.
• Almost zero inversion (I → ∞) can be reached where the system is almost transparent.
A complete transparency is impossible, as absorption is switched off close to trans-
parency. Carrier relaxation represented by T0 prevents the system from reaching zero
inversion.
• A state of positive inversion cannot be attained.

c.) T2 < T0 ≪ T p (inversion enslaved)


In this situation the time derivative the the rate equations can be neglected dI
dt
= 0 ⇒ dtd
Im(χ) =
0. The imaginary part of the susceptibility is a direct function of the acting optical field. Us-
ing equation (2.46) we find
Im(χ0 )
Im(χ) = replace with intensity
1 + Tℏ0Nε0 Im(χ0 )|E |2
Im(χ0 ) ε0
= 2T0
I0 = cn|E |2
1+ Im(χ0 )I 0 2
c ℏN n
Im(χ0 ) Nℏc n
= I0
saturation intensity I S = . (2.48)
1 + IS 2T0 Im(χ0 )

Real and imaginary part of the susceptibility are proportional to each other (c. f. equation (2.47))
χ0 (ω) Nℏc n
χ(ω) = IS = . (2.49)
1 + IIS0 2T0 Im(χ0 )

This is the simplest formulation describing the response of an optical transition. It is only
valid for time scales much longer than all relaxation times of the material. It describes a
single transition, to which the optical field is resonant. We only require two quantities for
the whole description:

25
2.3 The influence of damping NLO

• The ground state susceptibility χ0 (ω) of the two-level-system (usually only a small frac-
tion of the complete susceptibility)
• The saturation intensity I S .
Equation (2.48) is not equivalent to the expressions derived in chapter (1.5.2) (see e. g. equa-
tion (1.52)). It is not only highly dispersive, but also a Taylor expansion with respect to the
optical intensity yielding such quantities as e. g. n 2 will only converge in a limited range of
intensities. If the expansion is done at I 0 = 0

I0
χ ≈ χ0 − χ 0 (2.50)
IS

it will diverge for intensities higher than I s and will give wrong results much earlier.

26
3 Nonlinear optics and field propagation - a plane
wave approach
In the last chapter we have learned about the nonlinear response of the material. Together
with M AXWELL’s equations we can now in principle determine any kind of field propagation
in the presence of a nonlinear response. Unfortunately, the complete set of equations in-
cluding all nonlinear interaction is much too difficult to be solved in a self-consistent way.
Hence, approximations have to be made which are usually based on the assumption that
the nonlinear response of the material is small compared with the linear one and that some
peculiar features of the geometry allow for further simplifications. A few examples of such
treatments will be presented in that chapter.

3.1 Perturbation theory for nonlinear solutions

First we start with some assumptions:


• Monochromatic fields: The fields consist of a single frequency ω. If required, e. g. for
frequency conversion processes, several discrete frequencies are taken into account.
E real , B real , P real take the same form namely

E (r )e−iωt + E ∗ (r )eiωt .
¢
E real (r , t ) = (3.1)
2

• The material polarization at a particular frequency is split into a linear and a nonlinear
part representing the linear part with the dielectric tensor ε̂

P (r ) = P L (r ) + P NL (r ) = ε0 (ε̂(r ) − 1)E (r ) + P NL (r ). (3.2)

• We further assume a linearly homogeneous, but not yet isotropic system

ε̂(ω, r ) = ε̂(ω). (3.3)

Remark: An index profile as it is typical for waveguides can be treated in a similar man-
ner by taking waveguide modes instead of plane waves into account.
• Small perturbation due to nonlinear action: |P NL (r )| ≪ |ε0 ε̂(r )E (r )|.
We want to briefly state M AXWELL’s equations for monochromatic fields

∇ × E = iωB , ⃗ ∇ × B = −iωµ0 (ε0 ε̂E + P NL ), ⃗
∇ · ε0 ε̂E + P NL = 0, ⃗
¡ ¢
∇ · B = 0.

Taking the curl on FARADAY’s law leads to the wave equation

ω2
µ ¶
⃗ 1 1
∇⃗
³ ´
∇ · E − ∆E = 2 ε̂E + P NL with ε0 µ0 = 2
c ε0 c
ω 2 ¶
ω 1
2
µ
∆ −⃗
∇div + 2 ε̂ E = − 2 P NL . (3.4)
c c ε0
| {z } | {z }
linear operator small perturbation

27
3.2 Propagation in z-direction only NLO

First we want to look at solutions E L (r ) of the linear system. Then we can write E L (r ) as an
eigenvector of the eigenvalue problem

ω2
ω2
ε̂−1 (∆ − ⃗
∇ div)E L (r ) = − 2 E L (r ) with eigenvalue − . (3.5)
c c2

For a homogeneous infinite system we find the solution E L (r ) = E k eik·r . This is the starting
point to determine the nonlinear solution via perturbation theory.
Remark: One might be tempted to drop the term ⃗ ∇ · E L , but this is only possible in homo-
geneous isotropic media. In anisotropic media as they are typical for a quadratic nonlinear
response, the term ⃗
∇ · E L is not necessarily zero and longitudinal electric fields may occur

εx
 
0 0
∂E x ∂E y ∂E z
ε̂ = 0
 εy 0  ⇒⃗
∇ · (ε̂E L ) = εx + εy + εz =0
∂x ∂y ∂z
0 0 εz

⃗ ∂E x ∂E y ∂E z
=⇒
̸ ∇ · EL = + + = 0. (3.6)
∂x ∂y ∂z

We want to state requirements on nonlinearly induced processes (e. g. frequency conver-


sion):
• The spatial scales are much longer compared with linear process (e. g. conversion
lengths ≫ wavelength)
• The complex amplitudes of linear solutions may change but not their structure (e. g.
polarization).
The general structure of a nonlinear solution is

E (r ) = E k u(r ) exp(ik · r ). (3.7)

In what follows we discuss several cases for which further assumptions can be made.

3.2 Propagation in z -direction only

In many relevant case (e. g. frequency conversion) very broad collimated beams, which
travel into a single direction are used and any dependence on transverse coordinates (x, y)
can be neglected or taken into account parametrically only. We assume for our solution to
take the form

E (r ) = E β u(z) exp iβz .


¡ ¢
(3.8)

The related magnetic induction is obtained by the linear solution E β exp iβz via the curl
¡ ¢


∇ × (E β exp iβz ) = (iβêz ) × (E β exp iβz ) = iωB β exp iβz
¡ ¢ ¡ ¢ ¡ ¢

β
Bβ = êz × E β . (3.9)
ω

28
3.2 Propagation in z-direction only NLO

Now we want to insert our ansatz (3.8) into (3.4) where we use that u(z) is a slowly varying
envelope with
¯ ¯
¯d ¯
¯ u(z)¯ ≪ |βu(z)|. (3.10)
¯ dz ¯

This leads us to
ω2
µ ¶

= ∆ − ∇ div + 2 ε̂ [E k u(r ) exp(ik · r )]
c
d2 ω2 β
· ¸
β β
= (E − E z êz ) 2 + 2 ε̂E (u(z) exp iβz )
¡ ¢
dz c
2
ω2 β ¢ d2
· ¸ µ ¶
β β d β β d
= u(z) (E − E z êz ) 2 + 2 ε̂E exp iβz + (E − E z êz ) exp iβz
¡ ¢ ¡
+ 2iβ u(z)
dz c dz 2 dz
ω2
µ ¶ µ ¶
β ¢d d
⃗ β β
= u(z) ∆ − ∇ div + 2 ε̂ E exp iβz + (E − E z êz ) exp iβz
¡ ¢ ¡
+ 2iβ u(z)
c dz dz
| {z }
=0 linear solution
ω2 1
µ ¶
β β d d
+ 2iβ u(z) = − 2 P NL .
¡ ¢
= (E − E z êz ) exp iβz (3.11)
dz dz c ε0
Assuming the amplitude to vary slowly we can simplify
µ ¶
d
+ 2iβ u(z) ≈ 2iβu(z) (3.12)
dz
and multiplying a suitable vector from the left

β d ω2 1
(E β − E z êz ) exp iβz 2iβ u(z) = − 2 P NL
¡ ¢
dz c ε0
β
ω2 exp −iβz (E β − E z êz )∗ NL
¡ ¢
d
i u(z) = − 2 β
P . (3.13)
dz c 2βε0 |E β − E êz |2
z
q
W
As the amplitude u(z) has currently no unit, we change the unit to m2
that the squared
amplitude represents the intensity flowing in z-direction
rD E
β
a(z) = s z u(z) (3.14)

β
D E
with s z being the z-component of the time averaged P OYNTING vector of the linear solu-
tion which we want to calculate in the following:
D E
β 1
sz = Re(E × B ∗ )z
2µ0
1
Re[E β iβz × B β∗
¡ ¢ ¡  ¢
= exp exp
 −iβz ]z
2µ0
β 1 β h β2
· µ ¶¸
1 β β∗
i
= Re E × êz × E = Re |E | êz − E β (E β∗ êz )
2µ0 ω z 2µ0 ω z
1 β h
β
i 1 β¯ β ¯
β ¯2
¯
= |E β |2 − |E z |2 = ¯E − E z êz ¯ . (3.15)
2µ0 ω 2µ0 ω

29
3.3 Diffraction in linearly isotropic media NLO

We can summarize our findings of the nonlinearly driven evolution of a mode in z-direction

Evolution in z-direction
β
d ω2 1 (E β − E z êz )∗ NL ¡ ¢
i u(z) = − 2 P exp −iβz (3.16)
dz c 2βε0 |E β − E β êz |2
z
β
s
β
d ω3 (E − E z êz )∗ NL ¡ ¢
i a(z) = − P exp −iβz . (3.17)
dz 8βε0 c |E β − E β êz |
2
z

3.3 Diffraction in linearly isotropic media

To describe effects of a nonlinearly induced change of the beam profile one has to take into
account the effect of transverse coordinates to account for diffraction, self-focusing and spa-
tial soliton formation. To simplify the analysis, we restrict to isotropic media, which are most
relevant for cubic nonlinearities.
The starting point is again equation (3.4) with a scalar valued dielectric function ε̂ = ε

ω ´ ω2 1
∆ −⃗
³
∇ div + 2 ε E = − 2 P NL . (3.18)
c c ε0

Again we make the assumption of a weak nonlinear action |P NL | ≪ |ε0 εE | and homogeneity

∇ ε = 0 and use the first of M AXWELL’s equations

0 =⃗
∇·D =⃗
∇ · ε0 εE (r ) + P NL
£ ¤

≈⃗
∇ · [ε0 εE (r )] = ε0 ε⃗
∇ · E (r ) = 0. (3.19)

Then we can simplify (3.18) to

ω1 ω2 1
µ ¶
∆ + 2 E = − 2 P NL . (3.20)
c ε c ε0

Propagation mainly in z -direction


The field shape is determined by the respective linear wave E L (r ) = E β exp iβz with β2 =
¡ ¢
2
ε ωc 2 . The nonlinear wave takes the form

¯ ∂u ¯
¯ ¯
β
¡ ¢
E (r ) = E u(x, y, z) exp iβz with ¯¯ ¯¯ ≪ |β u(x, y, z)|, (3.21)
∂z

where u(x, y, z) is the slowly varying envelope. E β is the polarization of a linear wave propa-
gating in z-direction which is transverse (no z-component).

30
3.4 Pulse evolution in linearly isotropic media NLO

Derivation of the evolution of u(x, y, z)


We can find the evolution equation by inserting (3.21) into (3.20)

ω2
µ ¶
= ∆ + 2 ε [E β u(x, y, z) exp iβz ]
¡ ¢
c
à !
β 2 ∂ ∂2 ∂2 ∂2 ω2
= E exp iβz −β ε u(x, y, z)
¡ ¢
 + 2iβ ∂z + ∂z 2 + ∂x 2 + ∂y 2 + c

2
∂2 ∂2
· µ ¶ ¸
β
¡ ¢ d d
= E exp iβz  + 2iβ + + u(x, y, z)
dz dz ∂x 2 ∂y 2
∂2 ∂2 ω2 1
· ¸
β d
+ 2 + 2 u(x, y, z) = − 2 P NL .
¡ ¢
= E exp iβz 2iβ (3.22)
dz ∂x ∂y c ε0
We again assumed a slow varying amplitude. Similarly to section 3.2 we can multiply a suit-
able vector from the left and rearrange

1 ∂2 ∂2 ω2 1 ¢ E β∗ NL
· µ ¶¸
d ¡
i + + u(x, y, z) = − exp −iβz P . (3.23)
dz 2β ∂x 2 ∂y 2 c 2 2βε0 |E β |2
rD E
β
Furthermore we can apply the transformation a = s z u to arrive at intensities based on
β
D E
the z-component of the time averaged P OYNTING vector of the linear solution s z (c. f. (3.15))
D
β
E 1 β β2
sz = |E | . (3.24)
2µ0 ω
β
Note that E z is zero for isotropic media. We now obtain the evolution equation for nonlin-
earity driven diffraction

evolution for diffraction


s
1 ∂2 ∂2 ω3 E β∗ NL
· µ ¶¸
d ¡ ¢
i + + a(x, y, z) = − P (z) exp −iβz . (3.25)
dz 2β ∂x 2 ∂y 2 8βε0 c 2 |E β |

3.4 Pulse evolution in linearly isotropic media

To describe effects of a nonlinearly induced change of a pulse profile one has to take into
account a finite frequency spectrum and cannot restrict to monochromatic waves anymore
to account for dispersive spreading, self-compression and temporal soliton formation. To
simplify the analysis, we restrict to isotropic media and a propagation into a single direction
thus neglecting diffractive spreading.
We start with equation (3.20) for monochromatic waves propagating in an isotropic medium
with ε̂ = ε under the action of a weak nonlinear action at a fixed frequency
³ ω ´ ω2 1
∆ + 2 ε E (ω, z) = − 2 P NL (ω, z). (3.26)
c c ε0

31
3.4 Pulse evolution in linearly isotropic media NLO

As most of the quantities used in the further derivation are frequency dependent, we restrict
to a fixed frequency ω0 being in the middle of the investigated spectrum and representing
the carrier frequency of the pulse.
We assume that the field shape is similar to that of a linear wave propagating at ω0 for all
ω2
relevant frequency components E L (ω0 , r ) = E β0 exp iβ0 z with β20 = ε(ω0 ) c 20 . The nonlinear
¡ ¢

wave has the form


∂u
E (ω, r ) = E β0 u(ω, z) exp iβ0 z with | | ≪ |β u(ω, z)|
¡ ¢
(3.27)
∂z
with E β0 being the polarization of a linear (transverse) wave propagating in z-direction at
ω0 . u(ω, z) is again the slowly varying envelope. The spectrum is so narrow that the phase
evolution along z is still well represented by β0 .

Derivation of the evolution of u(ω, z)


We can find the evolution equation by inserting equation (3.27) into (3.26) and using the
slowly-varying-envelope-approach

ω2 ∂ ω2
µ ¶ · ¸
β0 iβ0 z 2
∆ + 2 ε(ω) E (ω, z) ≈ E e 2iβ0 + ε(ω) − β0 u(ω, z)
c ∂z c 2
ω2 1
= − 2 P NL (ω, z). (3.28)
c ε0
Now we can multiply a suitable vector from the left and rearrange
h ∂ β0 ∗
1 ω2 ω2 1 −iβ0 z E
µ ¶i
2
i + ε(ω) − β 0 u(ω, z) = − e β
P NL (ω, z). (3.29)
∂z 2β0 c 2 c 2β0 ε0
2
|E 0 |2
| {z }
(1)

2
We prepare for a F OURIER-back-transform under the assumption that β(ω) = ωc 2 ε(ω) is close
to β0
β(ω) + β0 1 GVD
(1) = (β(ω) − β0 ) ≈ β(ω) − β0 = (ω − ω0 ) + (ω − ω0 )2 + ... (3.30)
2β0 vg 2

where we use the expressions of group velocity v g and group velocity dispersion GVD
¯ ¸−1
∂ d2
· ¯
β(ω) ¯¯ β(ω)
¯ ¯
vg = , GVD = ¯ . (3.31)
∂ω ω0 dω2 ¯
ω0

Now we integrate the different spectral components/monochromatic waves of the pulse


ˆ∞ ˆ∞

dω E (z, ω)e−iωt + c.c. = Re dω E (z, ω)e−iωt
¤ £ ¤
E real (z, t ) =
2
0 0
ˆ
 ∞ 

= Re  dω E (z, ω)e−i(ω−ω0 )t e−iω0 t . (3.32)


0

32
3.4 Pulse evolution in linearly isotropic media NLO

Now we can insert our ansatz (3.27) and find

h ˆ ∞
i
β0 −i(β0 zω0 t )
E real (z, t ) = Re E dω u(ω, z)e−i(ω−ω0 )t e
| {z } | {z }
0 Fourier transform of SVE carrier wave
h i
= Re E β0 ũ(z, t )ei(β0 z−ω0 t ) , (3.33)

where we have the F OURIER transform of the slowly varying envelope as

ˆ∞
ũ(z, t ) = dω u(ω, z)e−i(ω−ω0 )t . (3.34)
0

Now we can apply the F OURIER transform (3.34) on equation (3.29) using approximation (3.30)

ˆ∞

· ¸
1 GVD
dω i + (ω − ω0 ) + (ω − ω0 ) + ... u(ω, z)e−i(ω−ω0 )t
2
∂z v g 2
0
∂ i ∂ GVD ∂2
· ¸
= i + − + ... ũ(t , z)
∂z v g ∂t 2 ∂t 2
ˆ∞
1 E β0 ∗
dω −ω2 P NL (ω, z)e−iωt e−i(β0 z−ω0 t )
¡ ¢
=
2β0 ε0 c 2 |E β0 |2
0
β0 ∗
E ∂2 NL
· ¸
1
= P̃ (z, t ) e−i(β0 z−ω0 t ) . (3.35)
2β0 ε0 c 2 |E β0 |2 ∂t 2

Remark: P̃ NL (z, t ) is not the real valued nonlinear polarization, but its complex counterpart,
the real part of which is the direct measurable polarization. It contains the positive frequency
part only, but not its complex conjugate counterpart on the negative frequency side.
rD E
β
Again we can apply the transformation a(z, t ) = s z 0 ũ(z, t ) to express the evolution in
β
D E
terms of intensities. As before s z 0 is the z-component of the time averaged P OYNTING
vector of the linear solution at the center frequency of the pulse ω0 . Again, we arrive at the
evolution equation for the slowly varying envelope of a pulse

evolution for a slowly varying pulse envelope

· β0 ∗ 2
∂ i ∂ GVD ∂2 E ∂
· ¸ ¸
1
i + − + ... a(t , z) = p β0 |2 ∂t 2
P̃ (z, t ) e−i(β0 z−ω0 t ) .
NL
∂z v g ∂t 2 ∂t 2 8β0 c ε0 ω0
2 |E
(3.36)

33
4 Field evolution for quadratic nonlinearities
The main subject of this chapter will be Second Harmonic generation (SHG). The basic prin-
ciple is that a fundamental harmonic (FH) wave at ω is frequency doubled to a second har-
monic (SH) wave oscillating at 2ω.
The geometry consists of plane waves propagating in z- ω 2ω
direction (no transverse structure). Note that the orien-
tation of the coordinates is defined by the propagation χ(2)
ω ω
direction and not by the nonlinear crystal. Because typ-
ical tensors of the nonlinear response are given with re-
spect to the crystal axis, an additional coordinate transformation might be required.

4.1 Equations of motion

The relevant types of nonlinear polarization are the following:

3
P i(2) (ω) = ε0 K (−ω|2ω, −ω) χ(2) (−ω|2ω, −ω)E 2ω j E ω∗k
X
i jk
with K (−ω|2ω, −ω) = 1
j ,k=1
3 1
P i(2) (2ω) = ε0 K (−2ω|ω, ω) χ(2) (−2ω|ω, ω)E ω j E ωk with K (−2ω|ω, ω) = .
X
i jk
(4.1)
j ,k=1 2

The interacting waves are one at FH frequency and one at SH frequency. We adapt (3.17) and
obtain
ω∗
s
d ω3 E⊥ NL
(z, ω) exp −iβ(ω)z
¡ ¢
FH i a ω (z) = − 2 ω∗ P
dz 8β(ω)ε0 c |E ⊥ |
s (4.2)
2ω∗
d ω3 E⊥ NL
¡ ¢
SH i a 2ω (z) = − P (z, 2ω) exp −iβ(2ω)z .
dz β(2ω)ε0 c 2 |E ⊥
2ω∗
|

ω/2ω
Here E ⊥ = E ω/2ω − E zω/2ω êz is the transverse field polarization at FH/SH frequency.
The total electric field driving the nonlinear polarization can be found by using (3.8), (3.14)
and (3.15)

a ω (z)
E ω (r ) = E ω u ω (z)eiβ(ω)z = E ω q­ ® eiβ(ω)z
s zω


s

= ω a ω (z)eiβ(ω)z
|E ⊥ | ε0 c 2 β(ω)
s
E 2ω 4ω
E 2ω (r ) = 2ω a 2ω (z)eiβ(2ω)z . (4.3)
|E ⊥ | ε0 c β(2ω)
2

34
4.1 Equations of motion NLO

We can now insert the fields (4.3) and the polarization (4.1) into the evolution equations (4.2)
s
d ω3 e−iβ(ω)z X ω∗ (2)
i a ω (z) = − ω E P i (ω)
dz 8β(ω)c 2 ε0 |E ⊥ | i =x,y i
s
ω3 e−iβ(ω)z X ω∗ X 3
=− ω E χ(2) (−ω|2ω, −ω)E 2ω j E ω∗k
8β(ω)c 2 |E ⊥ | i =x,y i j ,k=1 i j k
s
ω3 e−iβ(ω)z X ω∗ X 3
=− 2 ω Ei χ(2)
i jk
(−ω|2ω, −ω)
8β(ω)c |E ⊥ | i =x,y j ,k=1

E 2ω E ω∗
s s
j 4ω iβ(2ω)z k 2ω
· 2ω a 2ω (z)e ω a ω (z)∗ e−iβ(ω)z
|E ⊥ | ε0 c β(2ω)
2 |E ⊥ | ε0 c β(ω)
2

ω∗ ω∗ E 2ω∗
s
ω5 3 E ⊥i E j
i[β(2ω)−2β(ω)]z
χi j k ω kω
(2) ∗
X
= −e a 2ω a ω . (4.4)
β(ω) β(2ω)ε0 c i , j ,k=1
2 6 2ω
|E ⊥ | |E ⊥ | |E ⊥ |
| {z }
=:χ(2)
FH

Analogously we find the evolution equation for a 2ω (z)


s
d ω3 e−iβ(2ω)z X 2ω∗ (2)
i a ω (z) = − E P i (2ω)
dz β(2ω)c 2 ε0 |E ⊥ 2ω
| i =x,y i
s
ω3 e−iβ(2ω)z X 2ω∗ X 3
=− E χ(2) (−2ω|ω, ω)E ω j E ωk
4β(2ω)c 2 |E ⊥ 2ω
| i =x,y i j ,k=1 i j k
ω
E kω E j 2
s
2ω∗
−i[β(2ω)−2β(ω)]z ω5 3
(2)
E ⊥i
χ
X
= −e a (z).
2ω |E ω | |E ω | ω
(4.5)
β(ω)2 β(2ω)ε0 c 6 i , j ,k=1 i j k |E ⊥ | ⊥ ⊥
| {z }
=:χ(2)
SH

For further simplifications we introduce effective nonlinear coefficients of fundamental and


second harmonic as
ω∗ ω∗ E 2ω∗
s
(2) ω2 ω 3
(2)
E ⊥i Ek j
χFH = χ
X
(−ω|2ω, −ω) ω ω
β(ω)c 3 β(2ω)ε0 i , j ,k=1 i j k |E ⊥ | |E ⊥ | |E ⊥2ω
|
ω E ω (4.6)
2ω∗
s
ω2
ω 3 E E k j
χ(2) χ(2) (−2ω|ω, ω) ⊥i
X
= 2ω |E ω | |E ω |
.
SH β(ω)c 3 β(2ω)ε0 i , j ,k=1 i j k |E ⊥ | ⊥ ⊥

Now we assume that no further fields are involved and energy is conserved (Manley-Rowe
relations) which leads to the conservation of energy flux
µ ¶ µ ¶
d¡ 2 2
¢ ∗ d ∗ d
0= |a ω | + |a 2ω | = ia ω −i a ω + ia 2ω −i a 2ω + c.c.
dz dz dz
³ ´ ³ ´
∗ −i∆βz (2)
= ia ω e ∗
χFH a 2ω (z)a ω ∗
(z) + ia 2ω ei∆βz χ(2)
SH ω
a 2
(z) + c.c.
= ie−i∆βz χ(2) ∗2
a (z)a ω
FH 2ω
(z) − ie−i∆βz χ(2)∗ a a ∗ 2 (z) + c.c.
SH 2ω ω
h i
∗ 2 (2)
= −2 Im e−i∆βz a 2ω a ω (χFH − χ(2)∗
SH
) . (4.7)

35
4.2 Normalization NLO

Here we introduced the phase mismatch

∆β = 2β(ω) − β(2ω). (4.8)

Note that (4.7) must hold for every z. This can only be fulfilled by

χ(2)
FH
= χ(2)∗
SH
=: χ(2)
eff
. (4.9)

Now we may perform a phase transformation of the SH amplitude

a(z) = a ω (z) b(z) = a 2ω (z)e−i∆βz (4.10)

which gives us the final set of evolution equations

Evolution equations for SHG

d
FH: i a −∆βb + χ(2)
eff
a∗b = 0 (4.11)
dz
d
SH: i b − ∆βb + χ(2)∗
eff
a 2 = 0. (4.12)
dz

4.2 Normalization

We now aim towards a further simplification of the evolution equations for SHG which might
provide us with simple estimations of relevant scales. For that we introduce dimensionless
quantities:

z = Z0 Z , a = I 0 A and b = I 0 B eiϕ ,
p p
(4.13)

where Z0 is the characteristic length, I 0 the relevant intensity level and ϕ a constant phase
term. We insert these definitions into (4.11) and (4.12) and find

d
A −∆βZ0 B + Z0 χ(2)
p iϕ ∗
FH: i eff
I0e A B = 0 (4.14)
dZ
d
B − ∆βZ0 B + Z0 χ(2)
p −iϕ 2
SH: i eff
I 0 e A = 0. (4.15)
dZ
Now we can actually choose the normalization constants Z0 and I 0 in such a way, that all
constants become unity
p
1 I0
Z0 = , χ(2)
eff
eiϕ = 1. (4.16)
|∆β| |∆β|

Z0 now describes length scales for which the FH and SH waves run out of phase. For a per-
fectly phase-matched crystal this length scale goes to infinity. For the intensity we obtain

¯ ∆β ¯2
¯ ¯
1
I0 = = ¯ (2) ¯ . (4.17)
¯ ¯
(2)
|Z0 χ | ¯ χ ¯
eff eff

36
4.3 Conservation laws - the Hamiltonian NLO

The characteristic induces a noticeable nonlinear action (it is reduced by a growing non-
linear coefficient, but growing with the squared mismatch). This explains why second har-
monic generation is such a rare phenomenon. Particular care has to be taken to reduce the
mismatch. Otherwise required power levels grow to infinity.
The phase ϕ is chosen as
³ ´
ϕ = −arg χ(2)
eff
. (4.18)

It removes the phase of the nonlinear coefficient. This means the phase of the coefficient of
the quadratic nonlinearity has no physical relevance.
Remark: Normalization is not a unique procedure. For example, in case of vanishing mis-
match the sample length is much smaller than the inverse mismatch therefore becoming
the relevant scale.
The final set of normalized equations can be summarized to

Normalized evolution equations for SHG

d
FH: i A −σB + A ∗ B = 0 (4.19)
dZ
d
SH: i B − σB + A 2 = 0. (4.20)
dZ

σ = sgn(∆β) is called the scaled mismatch and can have values −1, 0, +1. The original fields
for a given A(Z ) and B (Z ) are


s
z iβ(ω)z−iωt
µ ¶
p h 2ω
E real (z, t ) = I 0 Re ω A Z = e
β(ω)ε0 c 2 |E ⊥ | Z0
s
E 2ω z 2iβ(ω)z−2iωt +iϕ i
µ ¶

+ B Z= e . (4.21)
β(2ω)ε0 c 2 |E ⊥2ω
| Z0

4.3 Conservation laws - the Hamiltonian

We have already enforced energy conservation onto our set of equations

Q = |A|2 + |B |2 = const. (4.22)

A second conserved quantity is the Hamiltonian of the system given by

H = σB ∗ B − A 2 B ∗ − A ∗2 B. (4.23)

We bring meaning to this Hamiltonian by considering that B ∗ corresponds to the creation


of a photon at the SH, whereas B describes the destruction of a photon. The same applies
for A. Then the first term describes the destruction of generated second harmonic light by
destructive interference in the mismatched case. The second term incorporates the SHG by
the annihilation of two photons from the FH wave creating a photon of the SH wave. The
last term describes down conversion of the SH wave into two photons of the FH wave.

37
4.4 General solution NLO

We now prove that the Hamiltonian is indeed conserved:

d d d d
H = σB ∗ B − A 2 B ∗ − 2AB ∗ A + c.c. (4.24)
dZ dz dz dz
= iσB ∗ (−σB + A 2 ) + iA 2 (−σB + A 2 )∗ − 2iAB ∗ A ∗ B + c.c.
= −iσ|B |2 + iσ|A|4 − 2i|A|2 |B |2 +  ∗2 ∗2
 
iσB A −
iσB A + c.c.
= 0. (4.25)

Since the other terms are purely imaginary they will cancel with they complex conjugated
terms. The quantum relation of the conserved quantities is the conserved norm of the wave
function which is equivalent to total intensity and the expectation value of the Energy (Hamil-
tonian operator) which is equivalent to the Hamiltonian.

4.4 General solution

Now we want to find a general solution for the set of normalized evolution equations for SHG.
We start with a transformation of complex variables into real ones by separating amplitudes
and phases

A(Z ) = |A(Z )|eiϕ A (Z ) and B (Z ) = |B (Z )|eiϕB (Z ) (4.26)

giving us four independent variables. We can insert this ansatz (4.26) into equations (4.19)
and (4.20) and divide by eiϕ A (Z ) and eiϕB (Z ) which yields

d ∂
FH i |A(Z )| − |A(Z )| ϕ A (Z ) + |A(Z )| · |B (Z )|eiϕB (Z )−2iϕ A (Z ) = 0
dZ ∂Z (4.27)
d ∂
SH i |B (Z )| − |B (Z )| ϕB (Z ) − σ|B (Z )| + |A(Z )|2 e2iϕ A (Z )−iϕB (Z ) = 0.
dZ ∂Z
Now we only consider the imaginary part of equation (4.27)

d
|A(Z )| + |A(Z )| · |B (Z )| sin ∆ϕ(Z ) = 0
¡ ¢
FH
dZ (4.28)
d
|B (Z )| − |A(Z )|2 sin ∆ϕ(Z ) = 0,
¡ ¢
SH
dZ
where ∆ϕ(Z ) = ϕB (Z ) − 2ϕ A (Z ) is the phase difference between SH and FH fields. We want
to investigate the evolution equation of the phase difference as well by using the real part of
equation (4.27)

d ∂ ∂
∆ϕ(Z ) = ϕB (Z ) − 2 ϕ A (Z )
dZ ∂Z ∂Z
|A(Z )|2
µ ¶
− 2|B (Z )| cos ∆ϕ(Z ) .
¡ ¢
= −σ + (4.29)
|B (Z )|

In all coherent systems without a fixed time reference the absolute phase does not matter
and only a phase difference plays a role. Thus we only have three independent real valued
variables ∆ϕ(Z ), |A(Z )|, |B (Z )|.

38
4.5 Phase-matched SHG NLO

We can reduce the number of free variables further by using the energy conservation (4.22)

d
|B (Z )| = Q − |B (Z )|2 sin ∆ϕ(Z )
¡ ¢ ¡ ¢
dZ
(4.30)
Q − |B (Z )|2
µ ¶
d
∆ϕ(Z ) = −σ + − 2|B (Z )| cos ∆ϕ(Z ) .
¡ ¢
dZ |B (Z )|

We can once more reduce the number of free variables by using the conserved Hamilto-
nian (4.23)

H = σ|B (Z )|2 − 2|A(Z )|2 |B (Z )| cos ∆ϕ(Z )


¡ ¢

= σ|B (Z )|2 − 2(Q − |B (Z )|2 )|B (Z )| cos ∆ϕ(Z ) .


¡ ¢
(4.31)

We can use this result to extract and eliminate sin ∆ϕ(Z ) = 1 − cos2 (∆ϕ(Z )) from equa-
¡ ¢ p

tion (4.30)
s
¶2
H − σ|B |2
µ
d 2 2
|B | = (Q − |B | ) sin ∆ϕ = ±(Q − |B | ) 1 −
¡ ¢
. (4.32)
dZ 2|B |(Q − |B |2 )

Now we transform to the SH power I = |B |2 by multiplying equation (4.32) with 2|B | because
dI
dZ
= 2|B | d|B
dZ
|

v
¶2
H − σI
u µ
d u
I = ±(Q − I ) 4I − 4I p
t
dZ 2 I (Q − I )
q
= ± 4I (Q − I )2 − (H − σI )2 . (4.33)

The formal solution can be obtained via integration. Using the initial conditions for given
SH power at Z = 0 : I (0) = I0

ˆ (L)
I ˆL
dI
p =± dZ = ±L. (4.34)
4I (Q − I )2 − (H − σI )2
I0 0

A practicable solution of equation (4.34) still requires solving the integral and inverting the
resulting function. This finally leads to so-called elliptic functions.
A nonlinear system with N real valued variables requires N − 1 conserved quantities or sym-
metries for an analytical solution. The existence of conserved quantities is essential for an
analytical solution.

4.5 Phase-matched SHG

Now we want to find solutions of (4.34) for the case of complete phase matching (σ = 0) and
no SH power at Z = 0

|B (Z = 0)|2 = I0 = 0. (4.35)

39
4.6 Down conversion (phase matched case) NLO

Looking at the conserved quantities we find that the total energy is Q = |A(Z = 0)|2 and the
Hamiltonian is simply zero H = 0 by comparing with (4.23). This simplifies (4.34) dramati-
cally and leads to
ˆ (L)
I
dI
p = ±L,
2(Q − I ) I
0
p !¯
Q + I ¯¯I (L)
Ãp
1
⇒ p ln p p = ±L. (4.36)
2 Q Q− I 0
¯

Inverting the results gives us the following expression

SH power in phase matched case


³p ´
I (L) = Q tanh2 QL . (4.37)

We conclude that total conversion is only reached for infinity propagation length. The con-
version efficiency η = I
p
Q depends only on the product QL. The required power scales in-
verse to the squared propagation length (see normalization in chapter 4.2).

1
SH power I

0.5

0
0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5
p (2) −1
propagation length in ( P in |χeff |)

Fig. 2: Up-conversion in the phase-matched case.

4.6 Down conversion (phase matched case)

A general property of Hamiltonian systems is the reversibility of their evolution. If SH gen-


eration is possible, as discussed in the previous chapter, also the opposite process called
parametric down conversion should be observed.
We assume the phase-matched case σ = 0 and all initial power to be in the SH wave

|B (Z = 0)|2 = I0 = Q and |A(Z = 0)|2 = 0. (4.38)

The Hamiltonian (4.23)

H = σ|B (Z )|2 − 2(Q − |B (Z )|2 )|B (Z )| cos ∆ϕ(Z ) = 0.


¡ ¢
(4.39)

40
4.7 Stability of the second harmonic field NLO

Only the initial value of I (I = Q instead of I = 0) has changed compared with SHG. The
derivation is analogous to section 4.5 but we have changed the initial conditions
p !¯
Q + I ¯¯I (L)
Ãp
1
p ln p p = ±L. (4.40)
2 Q Q− I Q
¯

This leads to some problems. The function on the left side of equation (4.40) becomes in-
finity on the lower boundary. The intensity of the SH wave can only deviate from Q, if the
propagation length L tends to infinity. This is consistent with an inversion of the SH process
discussed above. The latter one also needs an infinite propagation length to obtain complete
conversion.
A direct solution of equation (4.36) yields

d p
I = ±2(Q − I ) I ⇒ I (L) = Q. (4.41)
dZ
The system can stay in a state with all the power in the SH field and formally requires an
infinite propagation length to obtain noticeable down-conversion.

1
SH power I

0.5

0
−5 −4 −3 −2 −1 0 1 2 3 4 5
p (2) −1
propagation length in ( P in |χeff |)

Fig. 3: Down-conversion for the phase-matched case. The down conversion process starts from
quantum fluctuations.

In real systems any small perturbation (e. g. quantum noise) starts the down-conversion pro-
cess after a finite length. In the phase-matched case complete down-conversion is obtained
(inverse behaviour to SHG) and SHG starts (see previous chapter).

4.7 Stability of the second harmonic field

We have already noticed that down-conversion is started by field noise. To investigate that
phenomenon in more detail we perform a stability analysis. We take the results from sec-
tion 4.2
d
FH: i A −σB + A ∗ B = 0
dZ (4.42)
d 2
SH: i B − σB + A = 0.
dZ

41
4.7 Stability of the second harmonic field NLO

p −iσZ
We want to look for stationary solutions B 0 = Qe and A 0 = 0 with small fluctuations
δA(Z ) and δB (Z )
σ
³p ´
A(Z ) = δA(Z )e−i 2 Z and B (Z ) = Q + δB (Z ) e−iσZ . (4.43)

Then the evolution equations for the small perturbations become

d σ ³p ´
FH: i δA + δA + δA ∗ Q + δB
 =0
dZ 2 (4.44)
d
SH: i δB δA2
+ = 0.

dZ
We assume δA and δB to be small, so that quadratic terms can be neglected. Then the SH
perturbation δB remains
¡ in first¢ order approximation a constant (as long as δA(Z ) remains
small). Now we apply −i dz + σ2 to the first part of (4.44)
d

d σ σ
µ ¶·µ ¶ ¸
d ∗
p
0 = −i + i + δA + δA Q
dz 2 dZ 2
¶ ¸∗
σ 2 d σ
µ 2 ³ ´ ¶ ·µ
d p
= + δA + Q i + δA
dz 2 2 dz 2
σ 2
µ 2 ³ ´ ¶
(4.44) d
= + −Q δA. (4.45)
dz 2 2

the general solution of this equation can be expressed in terms of exponential functions
Ãr ! Ã r !
³ σ ´2 ³ σ ´2
δA(Z ) = δA + exp Q − Z + δA − exp − Q − Z . (4.46)
2 2

The two coefficients δA ± are fixed and determined by the initial fluctuations δA(Z = 0) and
their first derivative
σ
¯
d
= i δA(Z = 0) + i QδA ∗ (Z = 0).
p
δA(Z ) ¯¯
¯
(4.47)
dZ Z =0 2

Now we can use this and the initial conditions

d
r ³ σ ´2
δA(Z = 0) = δA + + δA − and δA(Z = 0) = Q− [δA + − δA − ] (4.48)
dZ 2
to find an expression for δA ±
à !
σ/2
p
1 Q
δA ± = ±ip δA(Z = 0) ± i p δA ∗ (Z = 0). (4.49)
2 4Q − σ2 4Q − σ2

Now we can distinguish two cases:

42
4.7 Stability of the second harmonic field NLO

³ σ ´2
a.) SH power is below a threshold Q ≤
2
This occurs for a mismatched sample (σ = ±1) only: Q ≤ 41 . The arguments in the exponential
terms become imaginary
Ãr ! Ã r !
1 1
δA(Z ) = δA + exp i −Q Z + δA − exp −i −Q Z . (4.50)
4 4

This leads to an oscillatory behaviour. The initial perturbation is conserved (which is typ-
ical for Hamiltonic systems). We do not observe growth of the SH wave or further down-
conversion. The system is (marginally) stable.

³ σ ´2
b.) SH power is above a threshold Q >
2
Then the fluctuations will grow exponentially
Ãr ! Ã r !
³ σ ´2 ³ σ ´2
δA(Z ) = δA + exp Q − Z + δA − exp − Q − Z . (4.51)
2 2
| {z }
≈ 0 for large Z

Inserting the equation for the coefficients δA ± we find for large Z

"Ã ! # Ãr !
σ/2 ³ σ ´2
p
1 Q
+ip δA(Z = 0) + i p δA ∗ (Z = 0) exp Q − Z . (4.52)
2 4Q − σ2 4Q − σ2 2

We observe that this solution is unstable. Also the initial perturbation δA and its complex
conjugate δA ∗ are amplified (so-called parametric amplifier).
We also want to discuss the case of large SH power with 4Q ≫ σ2 = 1. Then we can approxi-
mate δA(Z ) as
1¡ ¢ p
δA(Z ) ≈ δA(Z = 0) + iδA ∗ (Z = 0) e Q Z . (4.53)
2
The phase of the amplified perturbation for large Z is
Im[δA(Z = 0) + iδA ∗ (Z = 0)]
µ ¶
ϕ = arctan
Re[δA(Z = 0) + iδA ∗ (Z = 0)]
π
µ ¶
Im[δA(Z = 0)] + Re[δA(Z = 0)]
= arctan = arctan(1) = . (4.54)
Re[δA(Z = 0)] + Im[δA(Z = 0)] 4
This means that after a certain propagation length the phase of the perturbation is fixed and
independent of δA(Z = 0). This means the phase fluctuations (even of quantum origin) are
suppressed which leads to squeezed (non-classical) light.
Lastly we want to mention the critical intensity for down conversion in real world units Q <
1
4
¯ ¯2
³ σ ´2 1 ¯¯ ∆β ¯¯
I thresh. = I 0 = ¯ (2) ¯ . (4.55)
2 4¯χ ¯
eff

43
4.8 Visualization of the field dynamics using the Hamiltonian NLO

4.8 Visualization of the field dynamics using the Hamiltonian

We aim for a demonstration of the field dynamics in an intuitive way. We start again with
the conserved quantities, namely the Hamiltonian H given in (4.31) and the total intensity Q
given by (4.22). Then we can find another conserved quantity

H σ |B (z)|2 |B (Z )|2 |B (z)|


µ ¶

H (σ, |B |, ∆ϕ) = p cos ∆ϕ(Z )
¡ ¢
=p −2 1 −
Q 3/2 Q Q Q Q
| {z } | {z ′ }
=σ′ =:|B |

H (σ , |B |, ∆ϕ) = σ |B | − 2(1 − |B ′ |2 )|B ′ | cos ∆ϕ ,


′ ′ ′ ′ ′ 2
¡ ¢
(4.56)

where |B ′ | is the scaled amplitude of the SH field with 0 ≤ |B |′ ≤ 1, ∆ϕ the phase difference
between SH and FH fields and σ′ the rescaled mismatch.
The system now evolves along the lines/contours of fixed rescaled Hamiltonian H ′ . The con-
served direction of motion on each contour line is determined by

Q
µ ¶
d
∆ϕ(Z ) = −σ + − 3|B (Z )| cos ∆ϕ(Z )
¡ ¢
dZ |B (Z )|
µ ¶
1 d ′ 1 ′
∆ϕ(Z ) = −σ + − 3|B | cos ∆ϕ(Z ) .
¡ ¢
p ′
(4.57)
Q dZ |B |

We can draw the contour lines for several values of the Hamiltonian H ′ and different mis-
matches in figure 4.

Fig. 4: Contour plots of the scaled Hamiltonian H ′ for zero mismatch (left) and σ′ = 0.5 (right). The
radial coordinate |B ′ | describes the scaled SH power. For our initial conditions we start at con-
tour lines going through the origin.

For SHG we find a contour line of H ′ = 0 which starts at the origin (at zero SH power). For a
zero mismatch 4 (left) we see that |B ′ | approaches 1 (complete conversion) after infinite time.
However, if we introduce a mismatch 4 (right) then the inversion is incomplete, because the

44
4.9 Low-depletion-limit NLO

1
0.5
H′

H′
0
0
−0.5 1 1
0.5 ¡ ¢ 0.5 ¡ ¢
−1 −0.5 0 0
sin ∆ϕ −1 −0.5 0 sin ∆ϕ
¡ 0 ¢ 0.5 −0.5
1 −1 ¡ ¢ 0.5 −0.5
1 −1
cos ∆ϕ cos ∆ϕ

Fig. 5: 3D visualization of the scaled Hamiltonian for zero mismatch (left) and σ′ = 0.5 (right).

system evolves on a contour line which never reaches the border of |B ′ | = 1. We observe a
periodic motion corresponding to up- and down-conversion.
If we start outside the origin (e. g. simultaneously inject FH and SH fields), the phase between
SH and FH fields matters. We see a periodic motion around the minimum or maximum of
the scaled Hamiltonian. If we start at one of the extrema of H ′ we will have a stationary state
consisting of SH and FH fields.
The scaled mismatch σ′ = pσ controls the dynamics of the system. Increasing the power
Q
corresponds to a reduction of mismatch.

4.9 Low-depletion-limit

We now want to consider a system with a large mismatch with the unscaled equations (4.12).
We will find a low conversion efficiency, which means that the FH intensity is nearly constant
(undepleted pump approximation). Then
p
a(z) ≈ I FH = const. (4.58)

Then we may remove a(z) from the equation of the SH field (4.12) yielding
d
SH: i b − ∆βb + χ(2)∗ I = 0.
eff FH
(4.59)
dz
Using the initial condition of no SH intensity b(z = 0) = 0 we find the solution

χ(2)∗
eff ¡ ¢
b(z) = − I FH (exp −i∆βz − 1). (4.60)
∆β
The intensity of the SH wave is then

2
|χ(2)∗
eff
|2 ¯ ∆βz ¯2 ¯ ∆βz
2 ¯ i 2 ¯ ¯ −i 2
∆βz ¯2
¯
I SH = |b(z)| = I FH − ei 2 ¯
∆β2
¯e ¯ ¯e
´ 2
∆β
³
sin 2 z
= 4|χ(2)∗ | 2 2 
I FH
 . (4.61)
eff ∆β

45
4.9 Low-depletion-limit NLO

I SH

distance z in 1/∆β
1 2 3 4 5 6 7 8 9 10 11 12
L

Fig. 6: SH power in the undepleted pump approximation with a(z) = const. L describes the oscillation
period of SH-generation.

We observe (figure 6) that the intensity decays quickly with growing mismatch and oscillates
periodically with a period of
¯ ¯
¯ ¯
λ λ
¯ ¯ ¯ ¯ ¯ ¯
2π ¯ 2π ¯ ¯ 1 ¯ ¯ FH SH ¯
L= = ¯¯ ¯=¯
n n
¯=¯ ¯
|∆β| β(2ω) − 2β(ω) ¯
¯ SH − 2 FH ¯
¯ ¯ ¯ λFH n SH − 2λSH n FH ¯
¯λ λFH ¯
SH
¯ λSH
¯ ¯
¯
= ¯¯ ¯, (4.62)
n −n ¯
FH SH

where n FH , n SH are the refractive indices at FH/SH frequency and λFH , λSH the vacuum wave-
length of the FH/SH wave.

46
5 Diffraction, dispersion and nonlinearities
Now we want to discuss effects of further (transverse) coordinates which leads us to partial
differential equations instead of ODEs. Usually we will find no analytical solutions. However,
we need the transverse coordinates in order to account for diffraction.

5.1 Beam evolution in isotropic Kerr materials

We make several assumptions concerning the field propagation:


• The propagation takes place mainly in z-direction.
• We use a linearly polarized monochromatic electric field.
• The Kerr medium is isotropic: β = n(ω) ωc .
We use the evolution equation for nonlinearity driven diffraction that was derived in chap-
ter 3.3 (equation (3.25))
s
∂ 1 ∂2 ∂2 ω3 E β∗ NL
· µ ¶¸
¡ ¢
i + + a(x, y, z) = − β
P (z) exp −iβz . (5.1)
∂z 2β ∂x 2 ∂y 2 2
8βε0 c |E |

The monochromatic field implies that no Third-Harmonic Generation takes place and we
only have to deal with self-phase modulation. Since the electric field is linearly polarized, the
nonlinearity induces a refractive index change ∆n = n 2 I (I -intensity). Then the nonlinear
induced polarization takes the form

P NL (z) ∼ n 2 |a(x, y, z)|2 a(x, y, z) exp iβz .


¡ ¢
(5.2)

The nonlinearly induced change of the propagation constant is then


ω ω
∆β = ∆n = n 2 |a(x, y, z)|2 . (5.3)
c c
Now we can include the effect of the nonlinear polarization as a change of the propagation
constant ∆β

∂ 1 ∂2 ∂2 ω
· µ ¶ ¸
2
i + + + n 2 |a(x, y, z)| a(x, y, z) = 0. (5.4)
∂z 2β ∂x 2 ∂y 2 c

Now we again perform the normalization process with z = Z0 Z , x = W0 X , y = W0 Y and


p
a(x, y, z) = I 0 A(X , Y , Z )
" #
∂ ∂ ∂2 ω
µ 2
1 Z0

2
i + + + n 2 Z0 I 0 |A(X , Y , Z )| A(X , Y , Z ) = 0. (5.5)
∂Z 2 βW02 ∂X 2 ∂Y 2 c

47
5.1 Beam evolution in isotropic Kerr materials NLO

We can now choose the normalization constants in such a way that the prefactors become
one

Z0 2πnW02
=1 ⇒ Z0 = βW02 =
βW02 λ0
ω λ0 1 λ0 λ0 λ0 2
µ ¶
1
n 2 Z0 I 0 = 1 ⇒ I0 = = = . (5.6)
c 2π |n 2 |Z0 2π 2πn|n 2 |W02 n|n 2 | 2πW0

We see that W0 is a free scaling parameter, so we choose it e. g. as the initial beam width. This
means that the principal effects we discuss will not depend on the beam width.
Since we chose W0 to be the beam width its advantageous to introduce the total beam power

ˆ∞ ˆ∞ ˆ∞ ˆ∞
2
P= dx dy |a(x, y, z)| = W02 I 0 dX dY |A(X , Y , Z )|2 , (5.7)
−∞ −∞ −∞ −∞

µ ¶2
1 λ0
with P 0 = W02 I 0 = as a new scaling quantity.
n|n 2 | 2π
Note that P 0 only depends on material properties and wavelengths, not the beam width.
Qualitative changes are now to be expected to depend on the beam power in relation to P 0 ,
but not on the beam width, as it can be scaled away.
The final equation we have to solve is a nonlinear 2+1 dimensional S CHRÖDINGER equation
(NLS)

∂ 1 ∂2 ∂2
· µ ¶ ¸
2
i + + + κ|A| A = 0 with κ = sgn(n 2 ) = ±1. (5.8)
∂z 2 ∂X 2 ∂Y 2

However, this equation has no general solution (not integrable). The nonlinearity creates an
effective potential κ|A|2 ∼ −V . For κ = +1(n 2 > 0) high intensities correspond to potential
dips.

Conserved quantities
For our further discussions we want to investigate the conserved quantities. We start with

48
5.1 Beam evolution in isotropic Kerr materials NLO

the total power P


ˆ∞ ˆ∞
P = P0 dX dY |A(X , Y , Z )|2
−∞ −∞
ˆ∞ ˆ∞

· ¸
dP ∗
= P0 dX dY A A + c.c.
dz ∂Z
−∞ −∞
ˆ∞ ˆ∞
1 ∂ ∂2
µ · µ 2 ¶ ¸ ¶
∗ 2
= P0 dX dY A i + + |A| A + c.c.
2 ∂X 2 ∂Y 2
−∞ −∞
ˆ∞ ¶ ¯¯∞ ˆ∞
∂ 2
∂ ∂ ∗ ∂
µ
A∗ A = A∗ A ¯¯ − dX A A
∂X 2 ∂X −∞ ∂X ∂X
−∞ | {z } −∞
→0
ˆ∞ ˆ∞
i ∂ ∗ ∂ i ∂ ∗ ∂
µ ¶
4
= P0 dX dY − A A− A A + i|A| + c.c.
2 ∂X ∂X 2 ∂Y ∂Y
−∞ −∞
ˆ∞ ˆ∞
i ¯¯ ∂ ¯¯2 i ¯¯ ∂ ¯¯2
µ ¯ ¯ ¯ ¯ ¶
4
= P0 dX dY − ¯ A − ¯ A + i|A| + c.c. = 0. (5.9)
2 ∂X ¯ 2 ∂Y ¯
−∞ −∞
The last part cancels with its complex conjugate.
The second conserved quantity is the Hamiltonian which is given by
ˆ∞ ˆ∞ ³ ¯¯ ∂A ¯¯2 ¯¯ ∂A ¯¯2 ´
H= dX dY ¯¯ ¯ +¯ ¯ − κ|A|4 . (5.10)
∂X ¯ ¯ ∂Y ¯ | {z }
−∞ −∞ | {z } self-attraction
kinetic energy self-repulsion

The second term corresponds to a potential energy since |A|4 can be interpreted as a multi-
plication of particle density ∼ |A|2 times Potential ∼ |A|2 . We can proof that the Hamiltonian
is conserved in a similar way as for the total power.

Second moment
A third quantity we want to introduce is called the second moment defined by
ˆ∞ ˆ∞
S(Z ) = dX dY (X 2 + Y 2 )|A(X , Y , Z )|2 . (5.11)
−∞ −∞
Despite the other quantities
q the second moment is not conserved. It is a measure of the
beam width 〈W (Z )〉 = S(Z )
Q .

We can use the definition of S(Z ) to derive its evolution equation. However, we only want to
mention the result:
d2
S = 4H conserved. (5.12)
dZ 2
In the following we want to distinguish two cases:

49
5.1 Beam evolution in isotropic Kerr materials NLO

ˆ∞ ˆ∞
¯ ∂A ¯2 ¯ ∂A ¯2
µ¯ ¯ ¯ ¯ ¶
4
a.) H = dX dY ¯
¯ ¯ + ¯ ¯ − κ|A| > 0.
∂X ¯ ¯ ∂Y ¯
−∞ −∞

We fulfill this case always for κ = −1(n 2 < 0). Here n 2 corresponds to a defocussing nonlin-
earity. The case κ = +1 also works for prevailing kinetic energy (e. g. curved phase fronts).
Then we find

d2
S = 4H > 0. (5.13)
dZ 2
The solutions of this equations are parabolas which diverge for Z → ∞
s
S(Z )
〈W (Z )〉 = → ∞. (5.14)
Q

We can also draw the solutions as done in figure 7.

S S

Z
Z
ZC

Fig. 7: Solutions of the second moment S(Z ) for H > 0 (left) and H < 0 (right). Since S(Z ) > 0 per
definition, the evolution path is not possible after S(Z ) = 0.

ˆ∞ ˆ∞
¯ ∂A ¯2 ¯ ∂A ¯2
µ¯ ¯ ¯ ¯ ¶
4
b.) H = dX dY ¯
¯ ¯ + ¯ ¯ − κ|A| < 0.
∂X ¯ ¯ ∂Y ¯
−∞ −∞

This case only works for κ = +1(n 2 > 0). Using the evolution equation of the second moment
we find that the beam will collapse within a finite propagation distance

S(Z ) = 0 ⇒ W (Z ) = 0 ⇒ |A|2 → ∞. (5.15)

The collapse is unphysical because of zero width (paraxial approximation fails) and infinite
intensities (power expansion of P NL cannot be stopped after the third order). However, the
evolution is often very well described up to the damage threshold of the material.
Furthermore we want to calculate the critical power for self-focusing for κ = 1. We assume
that the initial field has a fixed beam profile f (X , Y ) but a varying amplitude: A(X , Y , Z =

50
5.2 Pulse evolution in a two-level medium NLO

0) = A 0 f (X , Y ). Then the Hamiltonian can be written as

ˆ∞ ˆ∞
¯ ∂A ¯2 ¯ ∂A ¯2
µ¯ ¯ ¯ ¯ ¶
H= dX dY ¯
¯ ¯ +¯
¯ ¯ − |A| = |A 0 |2 h kin − |A 0 |4 h int
4
(5.16)
∂X ¯ ∂Y ¯
−∞ −∞
ˆ∞ ˆ∞ ˆ∞ ˆ∞
¯ ∂ f ¯2 ¯ ∂ f
µ¯ ¯ ¯ ¯2 ¶
1 1
dY | f |4 .
¯
with h kin = dX dY ¯
¯ ¯ +¯ ¯ and h int = dX
2 ∂X ¯ ¯ ∂Y ¯ 2
−∞ −∞ −∞ −∞

Both h kin and h int are only shape dependent. We can easily find the critical power as

h kin
H = |A 0 |2 h kin − |A 0 |4 h int < 0 for |A 0 |2 > |A crit | = . (5.17)
h int

Every beam collapses above a certain (profile dependent) power level. The beam profile
with the minimum power for collapse is called the Townes soliton. It is a metastable (quasi-
stationary) state between spreading and self-focusing. It is cylindrically symmetric, has no
analytical solution and thus must be numerically determined. It looks similar to a Gaussian
profile (except of a simple exponential decay in the tails). Numerically we find Q Townes as

ˆ∞ ˆ∞
Q Townes = dX dY |A Townes (X , Y , Z )|2 = 5.85043. (5.18)
−∞ −∞

In real world units, the critical power for collapse would be


µ ¶2
1 λ0
P critical = Q Townes P 0 = 5.85043 . (5.19)
n|n 2 | 2π
2
For fused silica with n = 1.45, n 2 = 2,4 · 10−20 m
W
and a wavelength of λ0 = 800 nm we find

P critical = 2,8 MW. (5.20)

For a Gaussian profile the critical power is increased by 2 %.

5.2 Pulse evolution in a two-level medium

Lastly we want to consider the pulse evolution in a two-level-system. First we have to make
some assumptions:
• The propagation is entirely in z-direction (no diffraction).
• The medium is isotropic: β0 = n(ω0 ) ωc0 .
• The two-level-system is not damped (pulse duration ≪ damping time).
• The electric field is linearly polarized.
• Group velocity dispersion (GVD) of the background medium can be neglected.

51
5.2 Pulse evolution in a two-level medium NLO

The polarization of the two-level medium is


1h i
TLS
P real (z, t ) = P TLS (z, t )ei(β0 z−ω0 t ) + c.c. , (5.21)
2
where P TLS (z, t ) is the slowly varying envelope of the polarization of the two-level-system.
The total electric field is
1h i
E real (z, t ) = E (z, t )ei(β0 z−ω0 t ) + c.c. (5.22)
2 s
a(z, t ) a(z, t ) E β0 2
with E (z, t ) = E β0 ũ(z, t ) = E β0 rD E = E β0 q = β a(z, t ).
β ε0 c 2 β 0 β0 2 |E 0 | ε0 cn(ω0 )
sz 0 2 ω0 |E |

The evolution equation for the slowly varying envelope of a pulse is (see (3.36))
· β0 ∗ 2
∂ i ∂ E ∂
· ¸ ¸
1
i + a(t , z) = p β0 |2 ∂t 2
P̃ (z, t ) e−i(β0 z−ω0 t )
NL
∂z v g ∂t 8β0 c ε0 ω0
2 |E
· β0 ∗ 2 ³
1 E ∂ i(β0 z−ω0 t )
´¸
=p β0 |2 ∂t 2
TLS
P (z, t )e e−i(β0 z−ω0 t )
8β0 c ε0 ω0
2 |E
E β0 ∗
à !
1 2 ∂ ∂2 TLS
=p −ω0 − 2iω0  + 2 P (z, t )
8β0 c 2 ε0 ω0 |E β0 |2 ∂t ∂t
s
ω20 E β0 ∗ TLS
≈− P (z, t ). (5.23)
8ncε0 |E β0 |2
Now we use the evolution equations of the polarization for zero detuning ℏω0 = E a − E b
(equations (2.27) and (2.28))
d TLS Ne 2 ­
| ψb |r̂ |ψa |2 E I
®
i P =
dt ℏ (5.24)
dI 1
= E ∗ P TLS + c.c.
dt 2iℏN
q
Again we introduce a normalization z = Z0 Z , t = T0 T, a(Z , T ) = Iˆ0 A(Z , T ) and
E β0
P TLS (z, t ) = −i Pˆ0 P (Z , T ). (5.25)
|E β0 |
Then equation (5.23) becomes
s
∂ i Z0 ∂ Z0 ω20 E β0 ∗ E β0
· ¸ · ¸
i + A(Z , T ) = − q −i β P̂ 0 P (Z , T ) . (5.26)
∂z v g T0 ∂T ˆ 8ncε0 |E β0 | |E 0 |
I0
Inserting all the expressions for P TLS (5.25) and E (5.22) into the evolution equations we
find
s
E β0 ® 2 E β0
" #
Ne 2 ­
· ¸
d 2
q
i −i β Pˆ0 P (Z , T ) = T0 Iˆ0 | ψb |r̂ |ψa | A(Z , T ) I (5.27)
dT |E 0 | ℏ |E β0 | ε0 cn(ω0 )
q
#∗ ·
T0 Iˆ0 E β0
s
E β0
" ¸
dI 2
= A(Z , T ) ˆ
−i β P 0 P (Z , T ) + c.c.
dT 2iℏN |E β0 | ε0 cn(ω0 ) |E 0 |
(5.28)

52
5.2 Pulse evolution in a two-level medium NLO

We choose our scale parameters in such a way, that the constants vanish. This leads us to the
following conditions:
s
(5.26)Z0 Z0 ω20
⇒ = 1, q P̂ 0 = 1,
v g T0 8ncε0
Iˆ0
q
Iˆ0 Ne 2 ­
s
(5.27) T0 2
| ψb |r̂ |ψa |2
®
⇒ = 1, (5.29)
P̂ 0 ℏ ε0 cn(ω0 )
q
Iˆ0
s
(5.28) T0 2
⇒ P̂ 0 = 1.
ℏN ε0 cn(ω0 )

Some algebraic rearrangements result in the following solutions for the scale parameters
with µ := | ψb |r̂ |ψa |:
­ ®

s
1 2n c ℏ ε0 1
P̂ 0 = N eµ, T0 = , Iˆ0 = ℏω0 v g N . (5.30)
eµ v g ω0 N 4

The final set of normalized equation now looks like

Normalized equations in a two-level-system


P (Z , T ) = A I (5.31a)
∂T
∂ 1¡
I (Z , T ) = − A ∗ P + AP ∗
¢
(5.31b)
∂T µ2
∂ ∂

P (Z , T ) = + A(Z , T ). (5.31c)
∂Z ∂T

We use the initial conditions lim P (Z , T ) = 0, lim A(Z , T ) = 0 and lim Z (Z , T ) = −1


T →−∞ T →−∞ T →−∞
where the system is in the ground state.
We try to find solution using conserved quantities. One of these is

|P |2 + I 2 = 1 (5.32)
∂ ¡ 2 ∂P ∂P ∗
∂I
· ¸
1
|P | + I 2 = P ∗ = P ∗ AI + P A ∗ I + 2I − (A ∗ P + AP ∗ ) = 0.
¢
+P + 2I
∂T ∂T ∂T ∂T 2

In the following we restrict to real valued solutions for A, P, I which simplifies (5.31b) to


I (Z , T ) = −A(Z , T )P (Z , T ). (5.33)
∂T
Since (5.32) looks similar to the trigonometric P YTHAGOREAN, we introduce a substitution

∂Θ(Z , T )
P (Z , T ) = sin Θ(Z , T ), I (Z , T ) = − cos Θ(Z , T ), A(Z , T ) = − . (5.34)
∂T

53
5.2 Pulse evolution in a two-level medium NLO

This inherently fulfils the conservation law |P |2 +I 2 = 1, but also equations (5.31a) and (5.31b).
An evaluation of equation (5.31c) leads us to

∂ ∂ ∂
µ ¶
+ Θ(Z , T ) = − sin Θ(Z , T ). (5.35)
∂Z ∂T ∂T

For further simplification we transform into a moving spatial frame Z̃ = 2Z −T which changes
the derivatives accordingly

∂ ∂ ∂ ∂ ∂ ∂ ∂ ∂ ∂ ∂2
µ ¶ µ 2 ¶
⇒2 , ⇒ − , + ⇒ + . (5.36)
∂Z ∂ Z̃ ∂T ∂T ∂ Z̃ ∂Z ∂T ∂T ∂T 2 ∂ Z̃ 2

Using this, equation (5.35) transforms into the Sine-Gordon equation

∂2 ∂2
µ ¶
+ Θ( Z̃ , T ) = sin Θ( Z̃ , T ) (5.37)
∂ Z̃ 2 ∂T 2

which is actually solvable. We want to briefly state its physical origins:


• Pendulum under the action of gravity:
d2
Θ(T ) + sin Θ(T ) = 0
dT 2
• Set of coupled pendula under the action of gravity:
d2
Θ(Zi , T ) + sin Θ(Zi , T ) = [Θ(Zi +1 , T ) − Θ(Zi , T )] + [Θ(Zi −1 , T ) − Θ(Zi , T )].
dT 2 | {z }
∂2
∼ Θ(T, Z̃ )
∂ Z̃ 2
The sine-Gordon equation is integrable, however, it requires an infinite number of conserved
quantities. It bears conserved multi-soliton solutions, e. g. first order solitons like

Z̃ − vT
· µ ¶¸
Θ( Z̃ , T ) = 4 arctan exp ± p +δ , (5.38)
1 − v2
where 0 ≤ v < 1 is the velocity of the soliton. We observe that the pendula make a full turn
in angular space changing their phase from 0 to 2π (kink, the case of + sign) or from 2π to 0
(anti-kink, the case of − sign). For our description this means that the inversion goes from
−1 to +1 and returns to −1. The given solution is displayed in figure 8.

54
5.2 Pulse evolution in a two-level medium NLO

Fig. 8: Basic single soliton solution, simulated in Mathematica as systems of spring-coupled torsional
pendula, see V.G. Ivancevic and T.T. Ivancevig, JGSP 31(2013) 1–56

55

You might also like