% File: nlopt/lect1/lect1.tex [pure TeX code]
% Last change: January 7, 2003
%
% Lecture No 1 in the course ``Nonlinear optics'', held January-March,
% 2003, at the Royal Institute of Technology, Stockholm, Sweden.
%
% Copyright (C) 2002-2003, Fredrik Jonsson
%
\input epsf
\font\ninerm=cmr9
\font\twelvesc=cmcsc10
\def\lecture #1 {\hsize=150mm\hoffset=4.6mm\vsize=230mm\voffset=7mm
\topskip=0pt\baselineskip=12pt\parskip=0pt\leftskip=0pt\parindent=15pt
\headline={\ifnum\pageno>1\ifodd\pageno\rightheadline\else\leftheadline\fi
\else\hfill\fi}
\def\rightheadline{\tenrm{\it Lecture notes #1}
\hfil{\it Nonlinear Optics 5A5513 (2003)}}
\def\leftheadline{\tenrm{\it Nonlinear Optics 5A5513 (2003)}
\hfil{\it Lecture notes #1}}
\noindent\epsfxsize 100pt\epsfbox{../info/kthtext.eps}
\vskip-26pt\hfill\vbox{\hbox{{\it Nonlinear Optics 5A5513 (2003)}}
\hbox{{\it Lecture notes}}}\vskip 36pt\centerline{\twelvesc Lecture #1}
\vskip 24pt\noindent}
\def\section #1 {\medskip\goodbreak\noindent{\bf #1}
\par\nobreak\smallskip\noindent}
\def\subsection #1 {\smallskip\goodbreak\noindent{\it #1}
\par\nobreak\smallskip\noindent}
\lecture{1}
Nonlinear optics is the discipline in physics in which the electric
polarization density of the medium is studied as a nonlinear function
of the electromagnetic field of the light. Being a wide field of
research in electromagnetic wave propagation, nonlinear interaction
between light and matter leads to a wide spectrum of phenomena, such
as optical frequency conversion, optical solitons, phase conjugation,
and Raman scattering. In addition, many of the analytical tools applied
in nonlinear optics are of general character, such as the perturbative
techniques and symmetry considerations, and can equally well be applied
in other disciplines in nonlinear dynamics.
\section{The contents of this course}
This course is intended as an introduction to the wide field of
phenomena encountered in nonlinear optics.
The course covers:
\smallskip
\item{$\bullet$}{The theoretical foundation of nonlinear interaction
between light and matter.}
\item{$\bullet$}{Perturbation analysis of nonlinear interaction
between light and matter.}
\item{$\bullet$}{The Bloch equation and its interpretation.}
\item{$\bullet$}{Basics of soliton theory and the inverse
scattering transform.}
\smallskip
It should be emphasized that the course does not cover state-of-the-art
material constants of nonlinear optical materials, etc.~but rather
focus on the theoretical foundations and ideas of nonlinear optical
interactions between light and matter.
A central analytical technique in this course is the perurbation
analysis, with its foundation in the analytical mechanics.
This technique will in the course mainly be applied to the
quantum-mechanical description of interaction between light and
matter, but is central in a wide field of cross-disciplinary
physics as well.
In order to give an introduction to the analytical theory of
nonlinear systems, we will therefore start with the analysis
of the nonlinear equations of motion for the mechanical pendulum.
\section{Examples of applications of nonlinear optics}
Some important applications in nonlinear optics:
\medskip
\item{$\bullet$}{Optical parametric amplification (OPA) and
oscillation (OPO), $\hbar\omega_{\rm p}\to\hbar\omega_{\rm s}
+\hbar\omega_{\rm i}$.}
\item{$\bullet$}{Second harmonic generation (SHG),
$\hbar\omega+\hbar\omega\to\hbar(2\omega)$.}
\item{$\bullet$}{Third harmonic generation (THG),
$\hbar\omega+\hbar\omega+\hbar\omega\to\hbar(3\omega)$.}
\item{$\bullet$}{Pockels effect, or the linear electro-optical
effect (applications for optical switching).}
\item{$\bullet$}{Optical bistability (optical logics).}
\item{$\bullet$}{Optical solitons (ultra long-haul communication).}
\vfill\eject
\section{A brief history of nonlinear optics}
Some important advances in nonlinear optics:
\medskip
\item{$\bullet$}{Townes et al. (1960), invention of the
laser.\footnote{${}^1$}{Charles H.~Townes was in 1964 awarded
with the Nobel Prize for the invention of the ammonia laser.}
}
\item{$\bullet$}{Franken et al. (1961), First observation ever of
nonlinear optical effects, second harmonic generation
(SHG).\footnote{${}^2$}{Franken et al. detected
ulvtraviolet light ($\lambda=347.1$ nm) at twice the frequency of a
ruby laser beam ($\lambda=694.2$ nm) when this beam traversed a
quartz crystal; P.~A.~Franken, A.~E.~Hill, C.~W. Peters, G.~Weinreich,
Phys.~Rev.~Lett. {\bf 7}, 118 (1961).
Second harmonic generation is also the first nonlinear effect ever
observed where a coherent input generates a coherent output.}
}
\item{$\bullet$}{Terhune et al. (1962), First observation of
third harmonic generation (THG).\footnote{${}^3$}{In their experiment,
Terhune et al. detected only about a thousand THG photons per pulse,
at $\lambda=231.3$ nm, corresponding to a conversion of one photon
out of about $10^{15}$ photons at the fundamental wavelength at
$\lambda=693.9$ nm; R.~W. Terhune, P.~D. Maker, and
C.~M. Savage, Phys. Rev. Lett. {\bf 8}, 404 (1962).}
}
\item{$\bullet$}{E.~J.~Woodbury and W.~K.~Ng (1962),
first demonstration of stimulated Raman
scattering.\footnote{${}^4$}{E.~J.~Woodbury and W.~K.~Ng,
Proc.~IRE {\bf 50}, 2347 (1962).}
}
\item{$\bullet$}{Armstrong et al. (1962), formulation of general
permutation symmetry relations in nonlinear
optics.\footnote{${}^5$}{The general permutation symmetry
relations of higher-order susceptibilities were published
by J.~A.~Armstrong, N.~Bloembergen, J.~Ducuing, and
P.~S. Pershan, Phys.~Rev.~{\bf 127}, 1918 (1962).}
}
\item{$\bullet$}{A.~Hasegawa and F.~Tappert (1973),
first theoretical prediction of soliton generation in optical
fibers.\footnote{${}^6$}{A.~Hasegawa and F.~Tappert,
``Transmission of stationary nonliner optica pulses in dispersive
optical fibers: I, Anomalous dispersion; II Normal dispersion'',
Appl. Phys. Lett. {\bf 23}, 142--144 and 171--172
(August 1 and 15, 1973).}
}
\item{$\bullet$}{H.~M.~Gibbs et al. (1976), first demonstration
and explaination of optical
bistability.\footnote{${}^7$}{H.~M. Gibbs, S.~M. McCall,
and T.~N.~C. Venkatesan, Phys. Rev. Lett. {\bf 36}, 1135 (1976).}
}
\item{$\bullet$}{L.~F.~Mollenauer et al.~(1980),
first confirmation of soliton generation in optical
fibers.\footnote{${}^8$}{L.~F.~Mollenauer, R.~H. Stolen,
and J.~P. Gordon, ``Experimental observation of picosecond
pulse narrowing and solitons in optical fibers'',
Phys. Rev. Lett. {\bf 45}, 1095--1098 (September 29, 1980);
the first reported observation of solitons was though made in 1834 by
John Scott Russell, a Scottish scientist and later famous
Victorian engineer and shipbuilder, while studying water waves
in the Glasgow-Edinburgh channel.}
}
\smallskip
Recently, many advances in nonlinear optics has been made,
with a lot of efforts with fields of, for example, Bose-Einstein
condensation and laser cooling; these fields are, however,
a bit out of focus from the subjects of this course, which
can be said to be an introduction to the 1960s and 1970s
advances in nonlinear optics.
It should also be emphasized that many of the effects observed in
nonlinear optics, such as the Raman scattering, were observed
much earlier in the microwave range.
\section{Outline for calculations of polarization densities}
\subsection{Metals and plasmas}
From an all-classical point-of-view, the calculation of the
electric polarization density of metals and plasmas, containing
a free electron gas, can be performed using the model of free
charges acting under the Lorenz force of an electromagnetic field,
$$m_{\rm e}{{d^2{\bf r}_{\rm e}}\over{dt^2}}
=-e{\bf E}(t)-e{{d{\bf r}_{\rm e}}\over{dt}}\times{\bf B}(t),$$
where ${\bf E}$ and ${\bf B}$ are all-classical electric and magnetic
fields of the electromagnetic field of the light.
In forming the equation for the motion of the electron, the origin was
chosen to coincide with the center of the nucleus.
\subsection{Dielectrics}
A very useful model used by Drude and Lorentz\footnote{${}^9$}{R.~Becker,
{\it Elektronen Theorie}, (Teubner, Leipzig, 1933).} to calculate
the linear electric polarization of the medium describes the
electrons as harmonically bound particles.
For dielectrics in the nonlinear optical regime, as being the focus
of our attention in this course, the calculation of the electric
polarization density is instead performed using a nonlinear spring
model of the bound charges, here quoted for one-dimensional motion as
$$m_{\rm e}{{d^2{x}_{\rm e}}\over{dt^2}}
+\Gamma_{\rm e}{{d{x}_{\rm e}}\over{dt}}+\alpha^{(1)} x_{\rm e}
+\alpha^{(2)} x^2_{\rm e}+\alpha^{(3)} x^3_{\rm e}+\ldots
=-eE_x(t).$$
As in the previous case of metals and plasmas, in forming the equation
for the motion of the electron, the origin was also here chosen to
coincide with the center of the nucleus.
This classical mechanical model will later in this lecture be applied
to the derivation of the second-order nonlinear polarzation density
of the medium.
\section{\bf Introduction to nonlinear dynamical systems}
In this section we will, as a preamble to later analysis of
quantum-mechanical systems, apply perturbation analysis to a
simple mechanical system.
Among the simplest nonlinear dynamical systems is the pendulum,
for which the total mechanical energy of the system, considering the
point of suspension as defining the level of zero potential energy,
is given as the sum of the kinetic and potential energy as
$$E=T+V={1\over2m}|{\bf p}|^2 - mgl(\cos\vartheta-1),\eqno{(1)}$$
where $m$ is the mass, $g$ the gravitation constant, $l$ the length,
and $\vartheta$ the angle of deflection of the pendulum, and where ${\bf p}$
is the momentum of the point mass.
\medskip
\centerline{\epsfxsize=50mm\epsfbox{../images/pendulum/pendulum.1}}
\centerline{Figure 1. The mechanical pendulum.}
\medskip
\noindent
From the total mechanical energy~(1) of the system, the equations of
motion for the point mass is hence given by Lagranges
equations,\footnote{${}^{10}$}{Herbert Goldstein, {\it Classical Mechanics},
2nd ed.~(Addison-Wesley, Massachusetts, 1980).}
$$
{{d}\over{dt}}\bigg({{\partial L}\over{\partial p_j}}\bigg)
-{{\partial L}\over{\partial q_j}}=0,\qquad j=x,y,z,\eqno{(2)}
$$
where $q_j$ are the generalized coordinates, $p_j=\dot{q}_j$ are the
components of the generalized momentum, and $L=T-V$ is the Lagrangian of
the mechanical system. In spherical coordinates $(\rho,\varphi,\vartheta)$,
the momentum for the mass is given as its mass times the velocity,
$$
{\bf p}=ml(\dot{\vartheta}{\bf e}_{\vartheta}
+\dot{\varphi}\sin\vartheta{\bf e}_{\varphi})
$$
and the Lagrangian for the pendulum is hence given as
$$
\eqalign{
L&={{1}\over{2m}}(p^2_{\varphi}+p^2_{\varphi})+mgl(\cos\vartheta-1)\cr
&={{ml^2}\over{2}}(\dot{\vartheta}^2
+\dot{\varphi}^2\sin^2\vartheta)
+mgl(\cos\vartheta-1).\cr
}
$$
As the Lagrangian for the pendulum is inserted into Eq.~(2),
the resulting equations of motion are for $q_j=\varphi$ and $\vartheta$
obtained as
$$
{{d^2\varphi}\over{dt^2}}
+\bigg({{d\varphi}\over{dt}}\bigg)^2\sin\varphi\cos\varphi=0,
\eqno{(3{\rm a})}
$$
and
$$
{{d^2\vartheta}\over{dt^2}}+(g/l)\sin\vartheta=0,
\eqno{(3{\rm b})}
$$
respectively, where we may notice that the motion in ${\bf e}_{\varphi}$ and
${\bf e}_{\vartheta}$ directions are decoupled.
We may also notice that the equation of motion for $\varphi$ is
independent of any of the physical parameters involved in te Lagrangian,
and the evolution of $\varphi(t)$ in time is entirely determined by
the initial conditions at some time $t=t_0$.
In the following discussion, the focus will be on the properties of the
motion of $\vartheta(t)$.
The equation of motion for the $\vartheta$ coordinate is here
described by the so-called Sine-Gordon equation.\footnote{${}^{11}$}{The
term ``Sine-Gordon equation'' has its origin as an allegory over the
similarity between the time-dependent (Sine-Gordon) equation
$${{\partial^2\varphi}\over{\partial z^2}}-{{1}\over{c^2}}
{{\partial^2\varphi}\over{\partial t^2}}=\mu^2_0\sin\varphi,$$
appearing in, for example, relativistic field theories, as compared
to the time-dependent Klein-Gordon equation, which takes the form
$${{\partial^2\varphi}\over{\partial z^2}}-{{1}\over{c^2}}
{{\partial^2\varphi}\over{\partial t^2}}=\mu^2_0\varphi.$$
The Sine-Gordon equation is sometimes also called ``pendulum equation''
in the terminology of classical mechanics.}
This nonlinear differential equation is hard\footnote{$\dagger$}{Impossible?}
to solve analytically, but if the nonlinear term is expanded as a
Taylor series around $\vartheta=0$,
$$
{{d^2\vartheta}\over{dt^2}}+(g/l)
(\vartheta-{{\vartheta^3}\over{3!}}+{{\vartheta^5}\over{5!}}+\ldots)=0,
$$
Before proceeding further with the properties of the solutions to
the approximative Sine-Gordon, including various orders of nonlinearities,
the general properties will now be illustrated.
In order to illustrate the behaviour of the Sine-Gordon equation,
we may normalize it by using the normalized time $\tau=(g/l)^{1/2}t$,
giving the Sine-Gordon equation in the normalized form
$${{d^2\vartheta}\over{d\tau^2}}+\sin\vartheta=0.\eqno{(4)}$$
The numerical solutions to the normalized Sine-Gordon
equation are in Fig.~2 shown for initial conditions (a) $y(0)=0.1$,
(b) $y(0)=2.1$, and (c) $y(0)=3.1$, all cases with $y'(0)=0$.
\medskip
\centerline{\epsfxsize=45mm\epsfbox{sg010.eps}
\hskip 7.5mm\epsfxsize=45mm\epsfbox{sg210.eps}
\hskip 7.5mm\epsfxsize=45mm\epsfbox{sg310.eps}}
\centerline{Figure 2. Numerical solutions to the normalized Sine-Gordon
equation.}
\medskip
\noindent
First of all, we may consider the linear case, for which the approximation
$\sin\vartheta\approx\vartheta$ holds. For this case, the Sine-Gordon
equation~(4) hence reduces to the one-dimensional linear wave-equation,
with solutions $\vartheta=A\sin((g/l)^{1/2}(t-t_0))$.
As seen in the frequency domain, this solution gives a delta peak
at $\omega=(g/l)^{1/2}$ in the power spectrum $|\tilde{\vartheta}(\omega)|^2$,
with no other frequency components present.
However, if we include the nonlinearities, the previous sine-wave
solution will tend to flatten at the peaks, as well as increase in
period, and this changes the power spectrum to be broadened as well
as flattened out.
In other words, the solution to the Sine-Gordon give rise to a wide
spectrum of frequencies, as compared to the delta peaks of the
solutions to the linearized, approximative Sine-Gordon equation.
From the numerical solutions, we may draw the conclusion that whenever
higher order nonlinear restoring forces come into play, even such a
simple mechanical system as the pendulum will carry frequency components
at a set of frequencies differing from the single frequency given
by the linearized model of motion.
More generally, hiding the fact that for this particular case the
restoring force is a simple sine function, the equation of motion
for the pendulum can be written as
$${{d^2\vartheta}\over{dt^2}}+a^{(0)}+a^{(1)}\vartheta+a^{(2)}\vartheta^2
+a^{(3)}\vartheta^3+\ldots=0.\eqno{(5)}$$
This equation of motion may be compared with the nonlinear wave
equation for the electromagnetic field of a travelling optical
wave of angular frequency $\omega$, of the form
$$
{{\partial^2 E}\over{\partial z^2}}
+{{\omega^2}\over{c^2}}E
+{{\omega^2}\over{c^2}}(\chi^{(1)}E+\chi^{(2)}E^2+\chi^{(3)}E^3+\ldots)=0,
$$
which clearly shows the similarity between the nonlinear wave propagation
and the motion of the nonlinear pendulum.
Having solved the particular problem of the nonlinear pendulum,
we may ask ourselves if the equations of motion may be altered
in some way in order to give insight in other areas of nonlinear
physics as well.
For example, the series~(5) that define the feedback that tend to
restore the mechanical pendulum to its rest position clearly defines
equations of motion that conserve the total energy of the
mechanical system.
This, however, in generally not true for an arbitrary series
of terms of various power for the restoring force.
As we will later on see, in nonlinear optics we generally have a
complex, though in many cases most predictable, transfer of energy
between modes of different frequencies and directions of propagation.
\section{The anharmonic oscillator}
Among the simplest models of interaction between light and matter
is the all-classical one-electron oscillator, consisting of a
negatively charged particle (electron) with mass $m_{\rm e}$, mutually
interacting with a positively charged particle (proton)
with mass $m_{\rm p}$, through attractive Coulomb forces.
\medskip
\centerline{\epsfxsize=90mm\epsfbox{../images/spring/udspring.1}}
\smallskip
\centerline{Figure 5. Setup of the one-dimensional undamped spring model.}
\medskip
In the one-electron oscillator model, several levels of approximations
may be applied to the problem, with increasing algebraic complexity.
At the first level of approximation, the proton is assumed to
be fixed in space, with the electron free to oscillate around
the proton.
Quite generally, at least within the scope of linear optics,
the restoring spring force which confines the electron
can be assumed to be linear with the displacement distance
of the electron from the central position.
Providing the very basic models of the concept of refractive index and
optical dispersion, this model has been applied by numerous authors,
such as Feynman~[R.~P.~Feynman, {\sl Lectures on Physics} (Addison-Wesley,
Massachusetts, 1963)], and Born and Wolf~[M.~Born and E.~Wolf,
{\sl Principles of Optics} (Cambridge University Press, Cambridge, 1980)].
Moving on to the next level of approximation, the bound proton-electron
pair may be considered as constituting a two-body central force
problem of classical mechanics, in which one may assume a fixed center
of mass of the system, around which the proton as well as the electron
are free to oscillate.
In this level of approximation, by introducing the concept of reduced
mass for the two moving particles, the equations of motion for the
two particles can be reduced to one equation of motion, for the
evolution of the electric dipole moment of the system.
The third level of approximation which may be identified is when
the center of mass is allowed to oscillate as well, in which case
an equation of motion for the center of mass appears in addition to
the one for the evolution of the electric dipole moment.
In each of the models, nonlinearities of the restoring central force
field may be introduced as to include nonlinear interactions as well.
It should be emphasized that the spring model, as now will be introduced,
gives an identical form of the set of nonzero elements of the
susceptibility tensors, as compared with those obtained using a
quantum mechanical analysis.
Throughout this analysis, the wavelength of the electromagnetic
field will be assumed to be sufficiently large in order to neglect
any spatial variations of the fields over the spatial extent of the
oscillator system.
In this model, the central force field is modelled by a mechanical
spring force with spring constant $k_{\rm e}$, as shown schematically
in Fig.~5, and the all-classical Newton's equations
of motion for the electron and nucleus are
$$
\eqalign{
m_{\rm e}{{\partial^2 x_{\rm e}}\over{\partial t^2}}
&=\underbrace{-eE(t)}_{\rm optical}
-\underbrace{k_0(x_{\rm e}-x_{\rm n})
+k_1(x_{\rm e}-x_{\rm n})^2}_{\rm spring},\cr
m_{\rm n}{{\partial^2 x_{\rm n}}\over{\partial t^2}}
&=\underbrace{+eE(t)}_{\rm optical}
+\underbrace{k_0(x_{\rm e}-x_{\rm n})
-k_1(x_{\rm e}-x_{\rm n})^2}_{\rm spring},\cr
}
$$
corresponding to a system of two particles connected by a spring with
spring ``constant''
$$
k=-{{\partial F^{({\rm spring})}_{\rm e}}
\over{\partial (x_{\rm e}-x_{\rm n})}}
={{\partial F^{({\rm spring})}_{\rm n}}
\over{\partial (x_{\rm e}-x_{\rm n})}}
=k_0-2k_1(x_{\rm e}-x_{\rm n}).
$$
By introducing the reduced mass\footnote{${}^{12}$}{Herbert Goldstein,
{\it Classical Mechanics}, 2nd ed.~(Addison-Wesley, Massachusetts, 1980).}
$m_{\rm r} = m_{\rm e} m_{\rm n} / (m_{\rm e} + m_{\rm n})$ of the system,
the equation of motion for the electric dipole moment
$p=-e(x_{\rm e}-x_{\rm n})$ is then obtained as
$$
{{\partial^2 p}\over{\partial t^2}}+{{k_0}\over{m_{\rm r}}}p
+{{k_1}\over{e m_{\rm r}}}p^2
={{e^2}\over{m_{\rm r}}}E(t).\eqno(6)
$$
This inhomogeneous nonlinear ordinary differential equation for the
electric dipole moment is the primary interest in the discussion
that now is to follow.
The electric dipole moment of the anharmonic oscillator is now
expressed in terms of a perturbation series as
$$
p(t)=p^{(0)}(t)+\underbrace{p^{(1)}(t)}_{\propto E(t)}
+\underbrace{p^{(2)}(t)}_{\propto E^2(t)}
+\underbrace{p^{(3)}(t)}_{\propto E^3(t)}
+\ldots,
$$
where each term in the series is proportional to the applied electrical
field strength to the power as indicated in the superscript of repective
term, and formulate the system of $n+1$ equations for $p^{(k)}$,
$k=0,1,2,\ldots,n$, that define the time evolution of the electric
dipole. By inserting the perturbation series into Eq.~(6), we hence
have the equation
$$
\eqalign{
{{\partial^2 p^{(0)}}\over{\partial t^2}}&
+{{\partial^2 p^{(1)}}\over{\partial t^2}}
+{{\partial^2 p^{(2)}}\over{\partial t^2}}
+{{\partial^2 p^{(3)}}\over{\partial t^2}}
+\ldots\cr
&+{{k_0}\over{m_{\rm r}}}p^{(0)}
+{{k_0}\over{m_{\rm r}}}p^{(1)}
+{{k_0}\over{m_{\rm r}}}p^{(2)}
+{{k_0}\over{m_{\rm r}}}p^{(3)}
+\ldots\cr
&+{{k_1}\over{e m_{\rm r}}}
(p^{(0)}+p^{(1)}+p^{(2)}+\ldots)(p^{(0)}+p^{(1)}+p^{(2)}+\ldots)
={{e^2}\over{m_{\rm r}}}E(t).\cr
}
$$
Since this equation is to hold for an arbitrary electric field $E(t)$,
that is to say, at least within the limits of the validity of the
perturbation analysis, each set of terms with equal power dependence
of the electric field must individually satisfy the relation.
By sorting out the various powers and identifying terms in the left
and right hand sides of the equation, we arrive at the system
of equations
$$
\eqalign{
&{{\partial^2 p^{(0)}}\over{\partial t^2}}
+{{k_0}\over{m_{\rm r}}}p^{(0)}
+{{k_1}\over{e m_{\rm r}}}p^{(0)}{}^2=0,\cr
&{{\partial^2 p^{(1)}}\over{\partial t^2}}
+{{k_0}\over{m_{\rm r}}}p^{(1)}
+{{k_1}\over{e m_{\rm r}}}2p^{(0)}p^{(1)}
={{e^2}\over{m_{\rm r}}}E(t),\cr
&{{\partial^2 p^{(2)}}\over{\partial t^2}}
+{{k_0}\over{m_{\rm r}}}p^{(2)}
+{{k_1}\over{e m_{\rm r}}}(2p^{(0)}p^{(2)}+p^{(1)}{}^2)=0,\cr
&{{\partial^2 p^{(3)}}\over{\partial t^2}}
+{{k_0}\over{m_{\rm r}}}p^{(3)}
+{{k_1}\over{e m_{\rm r}}}(2p^{(0)}p^{(3)}+2p^{(1)}p^{(2)})=0,\cr
}
$$
where we kept terms with powers of the electric field up to and including
order three.
At a first glance, this system seem to suggest that only the first
order of the perturbation series depends on the applied electric
field of the light; however, taking a closer look at the system,
one can easily verify that all orders of the dipole moment
is coupled directly to the lower order terms.
The system of equations for $p^{(k)}$ can now be solved for
$k=0,1,2,\ldots$, in that order, to successively provide the
basis of solutions for higher and higher order terms, until reaching
some $k=n$ after which we may safely neglect the reamaining terms,
hence providing an approximate solution.\footnote{${}^{13}$}{It should
though be emphasized that in the limit $n\to\infty$, the described
theory still is an exact description of the motion of the electric
dipole moment within this model of interaction between light and matter.}
The zeroth order term in the perturbation series is decribed
by a nonlinear ordinary differential equation of order two,
a so-calles {\sl Riccati equation}, which analytically can be
solved exactly, either by directly applying the theory of Jacobian
elliptic integrals of by applying the Riccati
transormation.\footnote{${}^{14}$}{For examples of the application of the
Riccati transformation, see Zwillinger, {\it Handbook of Differential
Equations}, 2nd ed.~(Academic Press, Boston, 1992).}
However, by considering a system starting from rest, at a state of
equilibrium, we can immediately draw the conclusion that $p^{(0)}(t)$
must be identically zero for all times $t$.
This, of course, only holds for this particular model; in many
molecular systems, such in water, a permanent static dipole moment
is present, something that is left out in this particular spring model
of ours. (Not to be confused with the static polarization induced
by the electric field, which by definition of the terms in the
perturbation series is included in higher order terms, depending
on the power of the electric field.)
The first order term in the perturbation series is the first and only
one with an explicit dependence of the electric field of the light.
Since the zeroth order perturbation term is zero, the differential
equation for the first order term is linear, which simplifies the
calculus.
However, since it is an inhomogeneous differential equation, we must
generally look for a total solution to the equation as a sum of
a homogeneous solution (with zero right hand side) and a particular
solution (with the electric field in the right hand side present).
The homogeneous solution, which will contain two constants of integration
(since we are considering second-order ordinary differential equation)
will though only give the part of the solution which depend on initial
conditions, that is to say, in this case a harmonic natural oscillation
of the spring system which in the presence of damping terms rapidly
would decrease to zero.
This implies that in order to find steady-state solutions, in which
the oscillation of the dipole moment directly follows the oscillation
of the electric field of the light, we may directly start looking for
the particular solution.
For a time harmonic electric field, here taken as
$$E(t)=E_{\omega}\sin(\omega t),$$
the particular solution for the first order term is after some
straightforward algebra given as\footnote{${}^{15}$}{For the sake
of self consistency, the general solution for the first order
term is given as
$$p^{(1)}=A\cos((k_0/m_{\rm r})^{1/2}t)+B\sin((k_0/m_{\rm r})^{1/2}t)
+{{e^2/m_{\rm r}}\over{k_0/m_{\rm r}-\omega^2}}E_{\omega}\sin(\omega t),$$
where $A$ and $B$ are constants of integration, determined by initial
conditions.}
$$
p^{(1)}={{(e^2/m_{\rm r})}\over{(k_0/m_{\rm r}-\omega^2)}}E_{\omega}
\sin(\omega t),\qquad\omega^2\ne(k_0/m_{\rm r}).
$$
For a material consisting of $N$ dipoles per unit volume, and by
following the conventions for the linear electric susceptibility
in SI units, this corresponds to a first order electric polarization
density of the form
$$
\eqalign{
P^{(1)}(t)
&=P^{(1)}_{\omega}\sin(\omega t)\cr
&=\varepsilon_0\chi^{(1)}(\omega)E_{\omega}\sin(\omega t),\cr
}
$$
with the {\sl first order (linear) electric susceptibility} given as
$$
\chi^{(1)}(\omega)
=\chi^{(1)}(-\omega;\omega)
={{N}\over{\varepsilon_0}}{{(e^2/m_{\rm r})}\over{(\Omega^2-\omega^2)}},
$$
where the resonance frequency $\Omega^2=k_0/m_{\rm r}$ was introduced.
The Lorenzian shape of the frequency dependence is shown in Fig.~6.
\vfill\eject
\medskip
\centerline{\epsfxsize=90mm\epsfbox{chi1.eps}}
\smallskip
\centerline{Figure 6. Lorenzian shape of the linear susceptibility
$\chi^{(1)}(-\omega;\omega)$.}
\medskip
Continuing with the second order perturbation term, some straightforward
algebra gives that the particular solution for the second order term
of the electric dipole moment becomes
$$
\eqalign{
p^{(2)}(t)&=-{{k_1e^3}\over{2k_0 m^2_{\rm r}}}
{{1}\over{(\Omega^2-\omega^2)}}E^2_{\omega}
+{{k_1e^3}\over{2m^3_{\rm r}}}
{{1}\over{(\Omega^2-\omega^2)(\Omega^2-4\omega^2)}}
E^2_{\omega}\cr
&\qquad\qquad+{{k_1e^3}\over{m^3_{\rm r}}}
{{1}\over{(\Omega^2-\omega^2)(\Omega^2-4\omega^2)}}
E^2_{\omega}\sin^2(\omega t).\cr
}
$$
In terms of the polarization density of the medium, still with $N$
dipoles per unit volume and following the conventions in regular
SI units, this can be written as
$$
\eqalign{
P^{(2)}(t)
&=P^{(0)}_{0}+P^{(0)}_{2\omega}\sin(2\omega t)\cr
&=\underbrace{\varepsilon_0\chi^{(2)}(0;\omega,-\omega)
E_{\omega}E_{\omega}}_{\rm DC\ polarization}
+\underbrace{\varepsilon_0\chi^{(2)}(-2\omega;\omega,\omega)
E_{\omega}E_{\omega}\sin(2\omega)}_{\rm second
\ harmonic\ polarization}\cr
}
$$
with the {\sl second order (quadratic) electric susceptibility} given as
$$
\eqalign{
\chi^{(2)}(0;\omega,-\omega)&=
{{N}\over{\varepsilon_0}}{{k_1e^3}\over{2m^3_{\rm r}}}
\bigg[{{1}\over{(\Omega^2-\omega^2)(\Omega^2-4\omega^2)}}
-{{1}\over{\Omega^2(\Omega^2-\omega^2)}}\bigg],\cr
\chi^{(2)}(2\omega;\omega,\omega)&=
{{N}\over{\varepsilon_0}}{{k_1e^3}\over{m^3_{\rm r}}}
{{1}\over{(\Omega^2-\omega^2)(\Omega^2-4\omega^2)}}.\cr
}
$$
From this we may notice that for one-photon resonances, the nonlinearities
are enhanced whenever $\omega\approx\Omega$ or $2\omega\approx\Omega$,
for the induced DC as well as the second harmonic polarization density.
The explicit frequency dependencies of the susceptibilities
$\chi^{(2)}(-2\omega;\omega,\omega)$ and $\chi^{(2)}(0;\omega,-\omega)$
are shown in Figs.~7 and 8.
\vfill\eject
\medskip
\centerline{\epsfxsize=90mm\epsfbox{chi2a.eps}}
\smallskip
\centerline{Figure 7. Lorenzian shape of the linear susceptibility
$\chi^{(2)}(-2\omega;\omega,\omega)$ (SHG).}
\bigskip
\centerline{\epsfxsize=90mm\epsfbox{chi2b.eps}}
\smallskip
\centerline{Figure 8. Lorenzian shape of the linear susceptibility
$\chi^{(2)}(0;\omega,-\omega)$ (DC).}
\medskip
A well known fact in electromagnetic theory is that an electric dipole
that oscillates at a certain angular frequency, say at $2\omega$,
also emits electromagnetic radiation at this frequency.
In particular, this implies that the term described by the
susceptibility $\chi^{(2)}(-2\omega;\omega,\omega)$ will generate
light at twice the angular frequency of the light, hence generating
a second harmonic light wave.
\bye