---
title: "Modulational instability"
firstLetter: "M"
publishDate: 2021-02-26
categories:
- Physics
- Fiber optics
- Optics
- Perturbation
- Nonlinear optics

date: 2021-02-26T20:36:22+01:00
draft: false
markup: pandoc
---

# Modulational instability

In fiber optics, **modulational instability** (MI)
is a nonlinear effect that leads to the exponential amplification
of background noise in certain frequency regions.
It only occurs in the [anomalous dispersion regime](/know/concept/dispersive-broadening/)
($\beta_2 < 0$), which we will prove shortly.

Consider the following simple solution to the nonlinear Schrödinger equation:
a time-invariant constant power $P_0$ at the carrier frequency $\omega_0$,
which is experiencing [self-phase modulation](/know/concept/self-phase-modulation/):

$$\begin{aligned}
    A(z,t) = \sqrt{P_0} \exp( i \gamma P_0 z)
\end{aligned}$$

We add a small perturbation $\varepsilon(z,t)$ to this signal,
representing background noise:

$$\begin{aligned}
    A(z,t) = \big(\sqrt{P_0} + \varepsilon(z,t)\big) \exp( i \gamma P_0 z)
\end{aligned}$$

We insert this into the nonlinear Schrödinger equation to get a perturbation equation,
which we linearize by assuming that $|\varepsilon|^2$ is negligible compared to $P_0$,
such that all higher-order terms of $\varepsilon$ can be dropped, yielding:

$$\begin{aligned}
    0
    &= - P_0 \sqrt{P_0} \gamma - P_0 \gamma \varepsilon + i \pdv{\varepsilon}{z}
    - \frac{\beta_2}{2} \pdv[2]{\varepsilon}{t}
    + \gamma \big(\sqrt{P_0} + \varepsilon\big)^2 \big(\sqrt{P_0} + \varepsilon\big)^*
    \\
    &= i \pdv{\varepsilon}{z}
    - \frac{\beta_2}{2} \pdv[2]{\varepsilon}{t}
    + \gamma \big( P_0 (\varepsilon + \varepsilon^*) + \sqrt{P_0} |\varepsilon|^2
    + \sqrt{P_0} \varepsilon (\varepsilon + \varepsilon^*) + \varepsilon |\varepsilon|^2 \big)
    \\
    &= i \pdv{\varepsilon}{z} - \frac{\beta_2}{2} \pdv[2]{\varepsilon}{t} + \gamma P_0 (\varepsilon + \varepsilon^*)
\end{aligned}$$

We split the perturbation into real and imaginary parts
$\varepsilon(z,t) = \varepsilon_r(z,t) + i \varepsilon_i(z,t)$,
which we fill in in this equation.
The point is that $\varepsilon_r$ and $\varepsilon_i$ are real functions:

$$\begin{aligned}
    0
    &= i \pdv{\varepsilon_r}{z} - \pdv{\varepsilon_i}{z}
    - \frac{\beta_2}{2} \pdv[2]{\varepsilon_r}{t} - i \frac{\beta_2}{2} \pdv[2]{\varepsilon_i}{t}
    + 2 \gamma P_0 \varepsilon_r
\end{aligned}$$

Splitting this into its real and imaginary parts gives two PDEs
relating $\varepsilon_r$ and $\varepsilon_i$:

$$\begin{aligned}
    \pdv{\varepsilon_r}{z} = \frac{\beta_2}{2} \pdv[2]{\varepsilon_i}{t}
    \qquad \quad
    \pdv{\varepsilon_i}{z} = - \frac{\beta_2}{2} \pdv[2]{\varepsilon_r}{t} + 2 \gamma P_0 \varepsilon_r
\end{aligned}$$

We [Fourier transform](/know/concept/fourier-transform/)
these in $t$ to turn them into ODEs relating
$\tilde{\varepsilon}_r(z,\omega)$ and $\tilde{\varepsilon}_i(z,\omega)$:

$$\begin{aligned}
    \pdv{\tilde{\varepsilon}_r}{z} = - \frac{\beta_2}{2} \omega^2 \tilde{\varepsilon}_i
    \qquad \quad
    \pdv{\tilde{\varepsilon}_i}{z} = \Big(\frac{\beta_2}{2} \omega^2 + 2 \gamma P_0 \Big) \tilde{\varepsilon}_r
\end{aligned}$$

We are interested in exponential growth, so let us make the following ansatz,
where $k$ may be a function of $\omega$, as long as it is $z$-invariant:

$$\begin{aligned}
    \tilde{\varepsilon}_r(z, \omega) = \tilde{\varepsilon}_r(0, \omega) \exp(k z)
    \qquad \quad
    \tilde{\varepsilon}_i(z, \omega) = \tilde{\varepsilon}_i(0, \omega) \exp(k z)
\end{aligned}$$

With this, we can write the system of ODEs for
$\tilde{\varepsilon}_r(z,\omega)$ and $\tilde{\varepsilon}_i(z,\omega)$
in matrix form:

$$\begin{aligned}
    \begin{bmatrix}
        k & \beta_2 \omega^2 / 2 \\
        \beta_2 \omega^2 / 2 \!+\! 2 \gamma P_0 & - k
    \end{bmatrix}
    \cdot
    \begin{bmatrix} \tilde{\varepsilon}_r(0, \omega) \\ \tilde{\varepsilon}_i(0, \omega) \end{bmatrix}
    =
    \begin{bmatrix} 0 \\ 0 \end{bmatrix}
\end{aligned}$$

This has non-zero solutions if the system matrix' determinant is zero,
which is true when:

$$\begin{aligned}
    k = \pm \sqrt{ - \frac{\beta_2}{2} \omega^2 \Big( \frac{\beta_2}{2} \omega^2 + 2 \gamma P_0 \Big) }
\end{aligned}$$

To get exponential growth, it is essential that $\mathrm{Re}\{k\} > 0$,
so we discard the negative sign,
and get the following condition for MI:

$$\begin{aligned}
    - \frac{\beta_2}{2} \omega^2 \Big( \frac{\beta_2}{2} \omega^2 + 2 \gamma P_0 \Big) > 0
    \quad \implies \quad
    \boxed{
        \omega^2 < -\frac{4 \gamma P_0}{\beta_2}
    }
\end{aligned}$$

Since $\omega^2$ is positive, $\beta_2$ must be negative,
so MI can only occur in the ADR.
It is worth noting that $\beta_2 = \beta_2(\omega_0)$,
meaning there can only be exponential
noise growth when the "parent pulse" is in the anomalous dispersion regime,
but that growth may appear in areas of normal dispersion,
as long as the above condition is satisfied by the parent.

This result has been derived using perturbation,
so only holds as long as $|\varepsilon|^2 \ll P_0$.
Over time, the noise gets amplified so greatly
that this approximation breaks down.

Next, we define the **gain** $g(\omega)$,
which expresses how quickly the
perturbation grows as a function of the frequency offset $\omega$:

$$\begin{aligned}
    \boxed{
        g(\omega)
        = \mathrm{Re}\{k\}
        = \mathrm{Re} \bigg\{ \sqrt{ - \frac{\beta_2}{2} \omega^2 \Big( \frac{\beta_2}{2} \omega^2 + 2 \gamma P_0 \Big) } \bigg\}
    }
\end{aligned}$$

The frequencies with maximum gain are then found as extrema of $g(\omega)$,
which satisfy:

$$\begin{aligned}
    g'(\omega_\mathrm{max}) = 0
    \qquad \implies \qquad
    \boxed{
        \omega_\mathrm{max} = \pm \sqrt{\frac{2 \gamma P_0}{-\beta_2}}
    }
\end{aligned}$$

A simulation of MI is illustrated below.
The pulse considered was a solition of the following form
with settings $T_0 = 10\:\mathrm{ps}$, $P_0 = 10\:\mathrm{kW}$,
$\beta = -10\:\mathrm{ps}^2/\mathrm{m}$ and $\gamma = 0.1/\mathrm{W}/\mathrm{m}$,
whose peak is approximately flat, so our derivation is valid there,
hence it "wrinkles" in the $t$-domain:

$$\begin{aligned}
    A(0, t)
    = \sqrt{P_0} \sech\!\Big(\frac{t}{T_0}\Big)
\end{aligned}$$

<a href="pheno-mi.jpg">
<img src="pheno-mi-small.jpg">
</a>

Where $L_\mathrm{NL} = 1/(\gamma P_0)$ is the characteristic length of nonlinear effects.
Note that no noise was added to the simulation;
what you are seeing are pure numerical errors getting amplified.

If one of the gain peaks accumulates a lot of energy quickly ($L_\mathrm{NL}$ is small),
and that peak is in the anomalous dispersion regime,
then it can in turn also cause MI in its own surroundings,
leading to a cascade of secondary and tertiary gain areas.
This is seen above for $z > 30 L_\mathrm{NL}$.

What we described is "pure" MI, but there also exists
a different type caused by Raman scattering.
In that case, amplification occurs at the strongest peak of the Raman gain $\tilde{g}_R(\omega)$,
even when the parent pulse is in the NDR.
This is an example of stimulated Raman scattering (SRS).



## References
1.  O. Bang,
    *Numerical methods in photonics: lecture notes*, 2019,
    unpublished.
2.  O. Bang,
    *Nonlinear mathematical physics: lecture notes*, 2020,
    unpublished.