# Calculation of the refractive index for plane waves propagating in ionized gas

Alain Bourdier, Jean-Claude Diels, Olivier Delage

► To cite this version:

Alain Bourdier, Jean-Claude Diels, Olivier Delage. Calculation of the refractive index for plane waves propagating in ionized gas. Results in Physics, 2020, 18, pp.103250. 10.1016/j.rinp.2020.103250. hal-03066613

HAL Id: hal-03066613

https://hal.univ-reunion.fr/hal-03066613

Submitted on 15 Dec 2020

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Calculation of the refractive index for plane waves propagating in ionized gas

<figure><img src="image_2.png" /></figure>

Alain Bourdier<sup>a,*</sup>, Jean-Claude Diels<sup>a</sup>, Olivier Delage<sup>b</sup>

## ARTICLE INFO

**Keywords:**

Refractive index

Laser pulse propagation in the atmosphere

Nonlinear optics

Relativistic effects

## ABSTRACT

The plasma refraction index can be calculated by different methods. In this paper, a physical approach is presented and applied to non-linear situations. We have particularly focused on the case where the mass of the electrons subjected to the wave becomes relativistic. In the case of a circularly polarized high intensity wave, ionization is considered.

## Introduction

We are interested in the propagation of laser radiation in the atmosphere in order to study it at very large distances. The long wavelengths ($\lambda$ = 10.6 µm) are of great interest although they require very high powers [1,2]. To date, research in the wavelength range of 0.8–1 µm has shown a fundamental limitation of guided energy to a few mJ in a single 100 µm channel. With a wavelength of 10.6 µm, we expect to propagate very large energies in filaments, we want to go from mJ (GW) to some Joules (TW = $10^{12}$ W) [2]. We know the essential role played by the index of refraction in the propagation of the wave. It is the result of a balance between the Kerr effect and the diffraction that can only be achieved when the power of the laser is above a certain threshold. According to the transverse direction of the gradient of the index, the light is either diffracted or focused in the propagation. Long wavelengths make the threshold power very high and require, for a CO₂ laser to be of the order of one terawatt. At these powers the ionized electrons can become relativistic. Even if for this wavelength we only meet weakly relativistic situations, taking these effects into account allows us to envisage more frankly relativistic situations that we could encounter with longer wavelengths generated for example by using of a free electron laser. This work is part of a wider field of research on wave propagation properties, an area in which the knowledge of the refractive index can be a major issue. In this area studied for many years, very innovative works have been published recently [3,4].

The refraction index concept seems to be a trivial subject introduced by Snell's law we learned in our elementary optics courses. Beyond this apparent simplicity, this work aims to better understand its physical

origin. Among the index calculations that can be presented, we have emphasized on a physical approach described in the courses of R. P. Feynman [5] that we have extended to complex non-linear situations. This method has been validated by comparing the result obtained with those obtained by other approaches. The originality of one of these methods lies in the fact that the impulse and the energy of an electron can under certain conditions become a 4-potential vector. This greatly simplifies the disturbance calculation which allows to find the dispersion relation of the wave and therefore the index of the medium. We thus extend the scope of the “Feynman method” and confirm a result already found in another way [6]. Thus, we draw the expression of the index when the electron dynamics becomes relativistic. We show that the plane wave approximation remains valid for laser pulses of a few wavelengths in length. The effect of ionization is also studied. Situations where the Keldysh factor [7] is less than unity have been considered, tunneling ionization has been preferred. In the case where the wave is circularly polarized the current associated with the ionization of the air is important and cannot be neglected. It is supposed that this ionization is instantaneous and that the free electron density is constant. We consider electrons have a constant drift velocity as they have zero velocity when they leave the tunnel [8–10]. This drift velocity is assumed to be in one direction only. Thus, we simplify physics in order to be able to achieve a simple analytical calculation using “Feynman’s method” conveniently to start roughing the problem. We introduce a Lorentz transformation that eliminates the electron drift due to ionization. Then, the field produced by the plasma electrons can be calculated simply. The index is calculated in the laboratory frame. A relativistic correction to the index of refraction is found in the direction of the ionization current which shows that ionization affects the index.

<footnote><sup>a</sup>Department of Physics and Astronomy, The University of New Mexico, Albuquerque, NM, USA</footnote>

<footnote><sup>b</sup>Laboratoire de l'atmosphère et des cyclones (LACy), CNRS, UMR 8105, Université de la Réunion, France</footnote>

<footnote>* Corresponding author.</footnote>

<footnote>E-mail address: alain.bourdier@gmail.com (A. Bourdier).</footnote>
Light propagates slightly more slowly in air than in vacuum. This effect is described by the index of refraction $n$. In a material, the phase velocity of a wave is expressed as $c/n$. The index of refraction when plasma effects are neglected can be expressed as [11–13]

$$n = n_0 + \bar{n}_2 I + \dots, \quad (1)$$

where $n_0$ is the usual weak field refractive index, $\bar{n}_2$ the second order index associated to bound electrons and $I$ the intensity of the optical field. A filament is created when spreading by diffraction does not occur, that is to say when $P$ the total beam power satisfies

$$P \ge P_{cr} = \frac{\lambda^2}{2\pi n_0 \bar{n}_2}, \quad (2)$$

where $\lambda$ is the light wavelength. In the case when $\lambda = 10.6\mu m$, we set $\bar{n}_2 \approx 1.08 \times 10^{-23} m^2/W$ for air [1], then $P_{cr} \approx 1.65TW$. The filament exists if $P \ge P_{cr}$. Considering a 1 cm diameter filament first [2], the intensity must be greater than $I_{cr} \approx 2 \times 10^{12} W/cm^2$. The relevance of considering the relativistic mass of one electron in the wave is expressed by solving the following simple non-relativistic equation

$$m \frac{d^2 x}{dt^2} = -e E_0 \cos \omega t, \quad (3)$$

where $x$ is the electron position, $m$ its mass, $\omega$ the frequency of the wave and $-e$ the charge of the electron. The average normalized square of the electron velocity is given by

$$\left\langle \frac{\nu^2}{c^2} \right\rangle = \frac{\varepsilon}{2} = \frac{\nu_0^2}{2} = e^2 E_0^2 / 2m^2 \omega^2 c^2, \quad (4)$$

where parameter $\varepsilon = \nu_0^2$ measures the importance of relativistic effects. It can be easily shown that $\varepsilon = 7.4 \times 10^{-19} I(W/cm^2)\lambda^2(\mu m)$. For $I = I_{cr}$, we have $\varepsilon \approx 1.7 \times 10^{-4}$ which means that relativistic effects can be neglected. If we consider a 1 mm diameter and that $P = 10 \times P_{cr}$, then $\nu_0 \approx 4 \times 10^{-1}$. Considering the same power and a 1 cm diameter filament we found $\nu_0 \approx 4 \times 10^{-2}$. Then, in a very large number of cases the relativistic mass of electron can play an important role. For very short pulses the motion of the electrons might be quite relativistic, as a consequence it is necessary to explore relativistic effects for wavelengths close to 10 µm.

In order to calculate the correction for the index of refraction due to free electrons, in every situation, a source called the external source(S) is assumed to be at a large distance from a thin plate of plasma [5] (Fig. 1).

The spherical wave is emitted by the source becomes a plane wave when reaching the plate. The total electric field in any physical circumstance can always be represented by the sum of the fields from all the charges in the universe

$$E = E_S + \sum_{\text{allothercharges}} E_{\text{eachcharge}}. \quad (5)$$

The effect of all the charges is assumed to be very small compared to the one of the source, the different charges are only driven by the

<figure><img src="image_3.png" /><figcaption>Fig. 1. Thin plasma layer.</figcaption></figure>

<figure><img src="image_4.png" /><figcaption>Fig. 2. The path of a moving charge.</figcaption></figure>

source. This means that.

After the plate, the field becomes the field from the source plus a field $\Delta E$ produced by the oscillating electrons of the plate. This field variation is also calculated by considering that the wave appears to propagate at the speed $c/n$ in the plate. Comparing these two approaches allows to determine the index of refraction $n$.

## Calculation of the field created by a moving electron far away from the plasma plate

The field produced at $P$ at large distance from a moving electron is given by [5,14]

$$E = \frac{e}{4\pi\varepsilon_0 c^2} \frac{d^2 \hat{e}_{r'}}{dt^2}, \quad (6)$$

where $\hat{e}_{r'}$ is a unit vector in the apparent direction of the charge which stands at the retarded distance $r'(\tilde{x}', \tilde{y}', \tilde{z}')$ (Fig. 2). We have

$$\hat{e}_{r'} = \frac{PM'}{r'} = \frac{\tilde{x}'}{r'} i + \frac{\tilde{y}'}{r'} j + \frac{\tilde{z}'}{r'} z_p k, \quad (7)$$

where $z_p$ is the $z$-component of point $P$.

Further, calculating the field created by a small width plasma ring (Fig. 3), the charge will be supposed to move on a small distance compared to the distance $r'$ between the electron and $P$. The time delay, that is to say the time it takes to go, at speed $c$, from the electron to $P$ is $r'/c$ where $r'$ is the retarded distance. Assuming $r'$ is very large and neglecting terms in $1/r'^2$, we have

$$\frac{d}{dt} \left( \frac{\tilde{x}'}{r'} \right) = \frac{\tilde{x}'}{r'} - \frac{\tilde{r}'}{r'^2} \approx \frac{\tilde{x}'}{r'} \quad (8)$$

<figure><img src="image_5.png" /><figcaption>Fig. 3. Radiation produced by the plasma sheet.</figcaption></figure>
$$ \frac{d}{dt} \left( \frac{\tilde{z}'}{r'} - z_p \right) = \frac{\tilde{z}'}{r'} - \frac{\dot{r}' (\tilde{z}' - z_p)}{r'^2}, \qquad (9) $$

where $\dot{f} = df/dt$. As $r' \approx r_0 + \tilde{z}'$, we have

$$ \frac{d}{dt} \left( \frac{\tilde{z}'}{r'} - z_p \right) = \frac{\tilde{z}'}{r'} - \frac{\tilde{z}'}{r'} + \frac{\tilde{z} (r_0 + z_p)}{r'^2} \approx 0. \qquad (10) $$

We also have

$$ \frac{d^2}{dt^2} \left( \frac{\tilde{z}'}{r'} \right) \approx \frac{d}{dt} \left( \frac{\tilde{z}'}{r'} \right) = \frac{\tilde{z}'}{r'} - \frac{\dot{r}' \tilde{z}'}{r'^2} \approx \frac{\tilde{z}'}{r'}. \qquad (11) $$

Then, the electric field produced by the electron reads $E = (e/4\pi\epsilon_0c^2)[a(t - r'/c)/r']\hat{e}_{\infty}$, with $a(t - r'/c) = a[t(1 - v'.r/rc) - r/c]$ where $v'$ is the electron velocity. Far away from the plane, the relativistic correction $v'.r/rc$ can be neglected. We can replace $r'$ by $r$ which is a constant average distance in the denominator. In the numerator of the expression of the field relativistic effects can be ignored. Thus, the electron acceleration is just a function of the retarded time $t - r/c$. Consequently, the electric field created by an electron oscillating on a short distance is

$$ E = \frac{e}{4\pi\epsilon_0 c^2} \frac{a\left(t - \frac{r}{c}\right)}{r} \hat{e}_x, \qquad (12) $$

where $\hat{e}_x$ is a unit vector along the x-axis and $a$ the acceleration of the electron along this axis. This expression for the electric field produced by an electron when standing far away is derived in another way in Appendix A.

## Motion of an electron in the plasma plate when submitted to a linearly polarized external source

Let us derive the equation of one electron motion in the plasma sheet. When the external source is very intense, the electron relativistic mass will be considered by performing expansions and neglecting terms in $\gamma_0^3$. The plasma sheet is such as

$$ \Delta z \ll \lambda, \qquad (13) $$

where $\lambda$ is the wavelength of the external source.

The electromagnetic field is defined through the following 4-potential

$$ (V, A) = \left[ 0, -\frac{E_0}{\omega} \sin(\omega t - kz) \hat{e}_x \right]. \qquad (14) $$

It is assumed that the external source can be very intense, the relativistic electron mass can be significant. Thus, the following relativistic Hamiltonian for a free electron in the wave is taken

$$ H = \sqrt{\left[ P_x - \frac{eE_0}{\omega} \sin(\omega t - kz) \right]^2 + P_y^2 c^2 + P_z^2 c^2 + m^2 c^4}. \qquad (15) $$

Three constants of motion allow to integrate the system: $P_x$, $P_y$, and $C = H - (\omega/k)P_z$. It is assumed that $P_x = P_y = 0$. We choose $P_z = 0$ when $\omega t - kz = 0$, then $C = mc^2$.

The mechanical momentum along the x-axis is $p_x = m\gamma v_x = -(eE_0/\omega)\sin(\omega t - kz)$ where $\gamma$ is the Lorentz factor. Then, the equation of motion along this axis is

$$ \frac{dp_x}{dt} = m\dot{v}_x + m\gamma a = -eE_0\cos(\omega t - kz) + eE_0\frac{k}{\omega}\cos(\omega t - kz), \qquad (16) $$

where $a$ is the acceleration of the charged particle. The last term of this equation comes from the Lorentz force.

Let us make a number of reminders

$$ \dot{z} = \partial H / \partial P_z = \frac{P_z}{m\gamma} = \frac{P_z}{m\gamma}. \qquad (17) $$

As $C = mc^2$, we have

$$ p_z = mc^2 \frac{k}{\omega} (\gamma - 1). \qquad (18) $$

Keeping terms in $\epsilon$ only, Eqs. (17) and (18) imply that

$$ \dot{z} \approx \frac{1}{2} c \frac{v^2}{c^2} = c^2 \frac{k}{4\omega} \epsilon [1 - \cos 2(\omega t - kz)]. \qquad (19) $$

Eq. (19) shows that, in vacuum, the electron has a drift velocity along the z-axis. This drift velocity is not significant in plasma because of Coulomb forces.

As $\dot{\gamma} = (1 - v^2/c^2)^{-3/2}(v/c^2)a \approx (v/c^2)a$, the first term in the left-hand side of Eq. (17) is $m\dot{\gamma}v_x \approx m(v/c^2)v_x \approx m[v^2(1 - v^3/4c^3)/c^2]a \approx m(v^2/c^2)a$. Thus, the equation of motion along the x-axis is

$$ \frac{dp_x}{dt} = ma \left( 1 + \frac{3}{2} \frac{v^2}{c^2} \right) = -eE_0 \cos(\omega t - kz) + eE_0 \frac{k}{\omega} \dot{z} \cos(\omega t - kz). \qquad (20) $$

The last term in the right hand of Eq. (20), $(k/\omega)\dot{z} \approx \dot{z}/c \approx (1/2)v^2/c^2$ can be neglected in a non-relativistic approach, that is to say when terms in $v^2/c^2$ are neglected

Moreover, condition (13) allows to consider $E = E_0 \sin(\omega t) \hat{e}_x$. In the non-relativistic approach and, as in the relativistic one (when terms in $v^2/c^2$ are considered), the equation of motion of one electron in the layer will be given by

$$ \frac{dp_x}{dt} = -eE_0 \cos\omega t. \qquad (21) $$

## Calculation of the index of refraction for a plane wave propagating in plasma when the relativistic mass of electrons is neglected

### The wave is a pure plane wave

Only the “Feynman method” is applied in this preliminary part.

If the plate had no effect the field of the wave travelling to its right would be [5]

$$ E_S = E_0 \cos\omega \left( t - \frac{z}{c} \right), \qquad (22) $$

The wave takes an additional time when travelling in the plate. Without the plate, the wave would travel the distance $\Delta z$ with time of $\Delta z/c$. As the wave travels at the speed of $c/n$ in the plane, it takes the time $n \Delta z/c$. So, an additional time $\Delta t = (n - 1)\Delta z/c$ has to be considered. The extra delay in going through the plate is taken into account by replacing $t$ by $t - \Delta t$ in Eq. (22)

$$ \begin{aligned} E_{\text{after the plate}} &= E_0 \cos\omega \left[ t - (n - 1)\frac{\Delta z}{c} - \frac{z}{c} \right], \\ &= E_0 \cos\omega \left( t - \frac{z}{c} \right) + \frac{\omega(n-1)\Delta z}{c} E_0 \sin\omega \left( t - \frac{z}{c} \right), \\ &= E_0 \cos\omega \left( t - \frac{z}{c} \right) + \Delta E_{\omega}^{a}. \end{aligned} \qquad (23) $$

with

$$ \Delta E_{\omega}^{a} = \frac{\omega(n-1)\Delta z}{c} E_0 \sin\omega \left( t - \frac{z}{c} \right). \qquad (24) $$

The field $\Delta E_{\omega}^{a}$ which is added to the source field can be calculated in another way by calculating the field created by the electrons oscillating in the plate. There is no electron drift velocity acquired in ionization when considering linear polarization [8]. The Lorentz force can be neglected as we study a non-relativistic situation. Then, in the plasma sheet, at $z = 0$, the equation of motion of one electron is

$$ a = \frac{dv}{dt} = -\frac{eE_0}{m} \cos\omega t, \qquad (25) $$

as a consequence, the field generated by one electron is
<figure><img src="image_6.png" /><figcaption>Fig. 4. Fourier spectrum of the laser pulse. $\omega = 1.77 \times 10^{14} s^{-1}$. $\tau = 10 T$: solid line, $\tau = 3 T$: dashed line.</figcaption></figure>

$$E = -\frac{e^2 E_0}{4\pi\varepsilon_0 mc^2 r} \cos\omega \left(t - \frac{r}{c}\right) \hat{e}_x, \qquad (26)$$

The field created by a ring, that is to say the field produced by all the charges of the ring (Fig. 3), is obtained by adding the fields created by the infinitesimal bits of charge in surfaces $\delta^2S$

$$\delta^2 E_\omega = -\frac{Ne^2}{4\pi\varepsilon_0 c^2} \frac{E_0}{mr} \cos\omega(t - r/c), \qquad (27)$$

with $N = \eta\delta^2S$ where $\eta$ is the number of electrons per unit area. Then, the total field generated by the layer is obtained by integrating over all the rings

$$\Delta E_{\omega}^{b} = \int \delta E 2\pi \rho d\rho = -\frac{\eta e^2}{4\pi\varepsilon_0 c^2 m} \frac{E_0}{r} \int_0^{\infty} \cos\omega(t-r/c) \frac{2\pi\rho}{r} d\rho, \qquad (28)$$

as $r^2 = \rho^2 + z^2$ and as $z$ is independent of $\rho$ we have $2rdr = 2\rho d\rho$. Then, integrating over the $r$ - interval $[r = z, r = \infty]$ and neglecting the charge density when $r$ goes to infinity, we get [5]

$$\Delta E_{\omega}^{b} = -\frac{\eta e^2}{2\varepsilon_0 c} \frac{E_0}{m\omega} \sin\omega(t - z/c). \qquad (29)$$

Eq. (29) represents the part of the field which is not present in the original incoming one but produced through the mediation of the current in the plasma sheet.

Considering that $\eta = N\Delta z$ where N is the electron density and that plasma frequency is $\omega_p = \sqrt{e^2N/\varepsilon_0 m}$, comparing the two values of $\Delta E_\omega$ given by Eqs. (24), (29) and imposing their equality: $\Delta E_\omega^a = \Delta E_\omega^b$, we obtain

$$n = 1 - \frac{1}{2} \frac{\omega_p^2}{\omega^2}. \qquad (30a)$$

Laser pulses travel at the group velocity $v_G$ which is the speed of light is divided by $n_G$ the group index. We can deduce the group index $n_G$ from the relation $n_G n = 1$ [15]. Thus, in plasma, one has

$$n_G = 1 + \frac{1}{2} \frac{\omega_p^2}{\omega^2}. \qquad (30b)$$

These results (Eqs. (30a) and (30b)) are well known [15,16]. They can be derived in a different way by considering that the point of observation P is close to the electron and the mother ion [8] (Appendix B). Then, we focus on the static dipole field. These results could have also been obtained easily by introducing conductivity and using Maxwell equations, but we would have lost the insight into the origin of the fields.

### Calculation of the refractive index for a moderately short pulse

The source wave is chosen to be in the more realistic form

$$E_s = \hat{e}_x E_0 e^{-\frac{|t-z|}{\tau}} \cos\omega \left(t - \frac{z}{c}\right), \qquad (31)$$

where $\tau$ is the pulse duration.

As the thickness of the plasma layer is assumed to be very low, the wave phase is assumed to be function of $\omega t$ only inside. The Lorentz force is also neglected in this non-relativistic approximation. Thus, inside the very thin plasma sheet, the wave has the form

$$E_s = \hat{e}_x E_0 e^{-\frac{|t|}{\tau}} \cos\omega t. \qquad (32)$$

The length $\tau$ is assumed to be long enough so that no frequency of the Fourier spectrum is close to the plasma frequency. The frequencies of the Fourier spectrum are assumed to verify $\Omega \gg \omega_p$. Then, one can ignore nonlinear effects and the fact that the different frequency components move at different phase velocities. We assume that the pulse undergoes no distortion when propagating through the plasma sheet. The Fourier transform of this function is

$$\begin{aligned} E_{\omega} &= E_0 \int_{-\infty}^{+\infty} e^{-\frac{|t|}{\tau}} \cos\omega t e^{-i\Omega t} dt, \\ &= E_0 F(\Omega) = E_0 \tau \left[ \frac{1}{1 + (\omega - \Omega)^2 \tau^2} + \frac{1}{1 + (\omega + \Omega)^2 \tau^2} \right]. \end{aligned} \qquad (33)$$

Fig. 4 shows that the width at half maximum increases when $\tau$ decreases

The width at half maximum varies like $1/\tau: \Delta\Omega \sim 1/\tau$. We must have $\omega \pm 1/\tau > \omega_p$. For $\lambda = 10.6 \mu m$, one has $\omega = 1.77 \times 10^{14} s^{-1}$, if $\tau = 10T = 3.53 \times 10^{-13} s$ and the plasma density $10^{15} cm^{-3}$, then $\omega_p = 1.78 \times 10^{12} s^{-1}$. The condition to have no distortion mentioned just above is satisfied. When $\tau = T$, some significant frequencies are close to $\omega_p$ and distortion might start to take place (Fig. 5).

The acceleration of one electron, in the direction of the electric field, in the plasma layer is

$$a = -\frac{e}{m} E_0 e^{-\frac{|t|}{\tau}} \cos\omega t. \qquad (34)$$

The field produced by one electron of the plane at P is (Fig. 3)

$$E = \frac{e}{4\pi\varepsilon_0 c^2} \frac{a(t - r/c)}{r} = -\frac{e^2}{4\pi\varepsilon_0 mc^2} E_0 \frac{e^{-\frac{|t-r|}{\tau}}}{r} \cos\omega \left(t - \frac{r}{c}\right). \qquad (35)$$

The field created by a plasma ring is (Fig. 3)

$$\delta^2 E = -\frac{\eta(t) e^2 \delta^2 S}{4\pi\varepsilon_0 c^2} \frac{E_0}{mr} \left[ e^{-\frac{|t-r|}{\tau}} \cos\omega \left( t - \frac{r}{c} \right) \right]. \qquad (36)$$

The field created by all the planes is given by
<figure><img src="image_7.png" /><figcaption>Fig. 5. Fourier spectrum of the laser pulse. $\omega = 1.77 \times 10^{14} s^{-1}$, $\tau = T$.</figcaption></figure>

$$ \Delta E^b = \int \delta^2 E 2\pi\rho d\rho = -\frac{\eta e^2}{4\pi\varepsilon_0 c^2} \frac{E_0}{m} \int_0^\infty \left[ e^{-\frac{|t-r|}{\tau}} \cos\omega\left(t-\frac{r}{c}\right) \right] \frac{2\pi\rho}{r} d\rho. \quad (37) $$

Letting $r^2 = \rho^2 + z^2$ and $u = t - r/c$, $\Delta E$ becomes

$$ \Delta E^b = -\frac{\eta e^2}{2\varepsilon_0 c} \frac{E_0}{m} \int_{-\infty}^{\infty} e^{-\frac{|u|}{\tau}} \cos\omega u du. \quad (38) $$

Assuming $t - z/c \le 0$, one finds

$$ \Delta E^b = -\frac{\eta e^2}{2\varepsilon_0 c} \frac{E_0}{m} e^{\frac{1}{\tau}(t-\frac{z}{c})} \frac{\tau^2}{1+\omega^2\tau^2} \left[ \omega \sin\omega \left(t-\frac{z}{c}\right) + \frac{1}{\tau} \cos\omega \left(t-\frac{z}{c}\right) \right]. \quad (39) $$

As far as we focus on a part of the field which has the form of the source, we are interested in the following component of $\Delta E^b$

$$ \Delta E_{\cos\omega}^b = -\frac{\eta e^2}{2\varepsilon_0 c} \frac{E_0}{m} e^{\frac{1}{\tau}(t-\frac{z}{c})} \frac{\tau}{1+\omega^2\tau^2} \cos\omega \left(t-\frac{z}{c}\right). \quad (40) $$

Let us calculate $\Delta E$ by considering the fact that the pulse travels at the group velocity $v_G$ which is the speed of light is divided by the group index $n_G$. In the plate, the wave takes the additional time $\Delta t = (n_G - 1)\Delta z/c$. In this case, after the plasma sheet, $t$ has to be replaced by $t - \Delta t$ in expression of the wave

$$ E_{\text{after the plate}} = E_0 e^{\frac{t-(n_G-1)\Delta z/c-z/c}{\tau} + i\omega[t-(n_G-1)\Delta z/c-z/c]}. \quad (41) $$

As a consequence

$$ \Delta E_{\cos\omega}^a = -E_0 e^{\frac{t-\frac{z}{c}}{\tau} n_G - \frac{1}{\tau} \frac{\Delta z}{c} \cos\omega \left(t-\frac{z}{c}\right)}. \quad (42) $$

The fact that $\Delta E_{\cos\omega}^a = \Delta E_{\cos\omega}^b$ implies that

$$ \frac{n_G - 1}{\tau} \frac{\Delta z}{c} = \frac{\eta e^2}{2\varepsilon_0 mc} \frac{\tau}{1 + \omega^2 \tau^2}, \quad (43) $$

so

$$ n_G = 1 + \frac{1}{2} \frac{\omega_p^2}{\omega^2} K, \quad (44a) $$

and

$$ n = 1 - \frac{1}{2} \frac{\omega_p^2}{\omega^2} K. \quad (44b) $$

with

$$ K = \frac{\omega^2 \tau^2}{1 + \omega^2 \tau^2}. \quad (45) $$

Note that considering $t - z/c \ge 0$ leads to the same result.

In linear physics, this result can be found very simply by using Maxwell equation and introducing conductivity.

When $\tau$ becomes very large, the result of the non-relativistic plane wave is found again (Eq. (30b)). When the pulse duration, $\tau$, is 10 laser periods: $\tau = 10$ T, then the correction with respect to the case of the plane wave is $K \approx 0.999$, the correction is very weak and the plane wave approximation is valid.

## Calculation of the index of refraction for a strong linearly polarized relativistic plane wave

### One particle approach

We start by developing a one-particle electrostatic model that we will then insert into a plasma approach. The $v \times B$ force due to the wave applied to the electrons of the plasma sheet produces a density perturbation $\delta N$, having frequency $2\omega$, $\delta N$ will be calculated through a plasma approach. Then, the one-electron contribution to the field will be multiplied by $\delta N$ to obtain the contribution per unit plasma volume. The resulting current created, $\delta N e v_x$, is a source for $E$ and $B$.

When the external source is very intense, the electron relativistic mass will be considered by performing expansions and neglecting terms in $v_0^3$.

It is assumed that the electrons in the plasma layer remain at $z = 0$ and oscillate along the x-axis. The last term in Eq. (20), $eE_0(k/\omega)\dot{z}\cos(\omega t - kz)$, is ignored. In fact, it is considered when coupling this single particle motion to the density perturbation $\delta N$. These strong hypotheses will be justified later by comparing the result obtained through this approach to the one obtained with a rigorous plasma approach. Thus, the relativistic motion of one electron in the wave is described by Eq. (21)

$$ m\dot{\gamma}v_x + m\gamma a = -eE_0\cos\omega t, \quad (46) $$

while the plasma is coupled to the incident wave through density oscillations due to the longitudinal Lorentz force ($f = -ev \times B$).

Retaining second order terms in $v^2/c^2$ only, the electron acceleration along the x-axis is

$$ a = -\left(1 - \frac{3}{2} \frac{v^2}{c^2}\right) \frac{eE_0}{m} \cos\omega t, \quad (47) $$

Considering $v$ is the zero-order velocity $v_0 = -(eE_0/m\omega)\sin\omega t$, the acceleration formulation $a$ becomes

$$ a = \frac{-eE_0}{m} \left[ \left(1 - \frac{3}{8}\varepsilon\right) \cos\omega t + \frac{3}{8}\varepsilon\cos3\omega t \right]. \quad (48) $$

The field created by the slice of plasma due to all the electrons having an oscillating density due
to the incoming wave (Fig. 3) can be expressed as

$$ \Delta E^b = - \frac{\eta(t) e^2 E_0}{2\epsilon_0 c} \left[ \left( 1 - \frac{3}{8}\epsilon \right) \sin\omega(t - z/c) + \frac{1}{8}\epsilon \sin3\omega(t - z/c) \right], \quad (49) $$

where $\eta(t)$ oscillates due to the plasma oscillations.

### Calculation of the index by considering plasma density oscillations

Let us determine $\eta(t)$ through a plasma approach by solving Maxwell and Lorentz equation [17,18].

The electromagnetic field can be given by

$$ \begin{aligned} e\mathbf{E} &= -\frac{\partial\mathbf{p}}{\partial t} - \mathbf{grad}cp_0, \\ e\mathbf{B} &= \mathbf{rot}p, \end{aligned} \quad (50) $$

where $\mathbf{p}$ is the fluid mechanical momentum of the free electrons and where $p_0 = (p^2 + m^2c^2)^{1/2}$ is their energy. It means that, in some cases, $(cp_0/e, \mathbf{p}/e)$ is a four-vector potential. It is known that two four vectors potential $(\phi, \mathbf{A})$ and $(\phi', \mathbf{A}')$ give the same electromagnetic field when $\mathbf{A}' = \mathbf{A} + \mathbf{grad}f(\mathbf{r}, t)$, $\phi' = \phi - \partial f(\mathbf{r}, t)/\partial t$ where $f$ is some function depending on $\mathbf{r}$ and $t$ only [19]. We assume that the electric and magnetic fields of the wave are associated with the four-potential $(\phi, \mathbf{A})$, the four-vector potential $(cp_0/e, \mathbf{p}/e)$ can be allowed with to such a function. The initial position being $\mathbf{n}_0$, the Jacobi action $S(\mathbf{r}, \mathbf{n}_0, t)$ can be introduced by calculating the action of one electron on two adjacent very close actual trajectories [20]. We have $\mathbf{P} = \partial S(\mathbf{r}, \mathbf{n}_0, t)/\partial \mathbf{r}$, $\mathbf{P}_0 = -\partial S(\mathbf{r}, \mathbf{n}_0, t)/\partial \mathbf{r}_0$ where $\mathbf{P}$ and $\mathbf{P}_0$ are respectively the canonical momenta at some position $\mathbf{r}$ and at the initial position $\mathbf{n}_0$ of the electron. We have $\mathbf{p}/e = \mathbf{A} + \partial/\partial \mathbf{r}(S/e)$ and $cp_0/e = \phi - \partial/\partial t(S/e)$. Consequently $S/e$ is a function $f$ if $S$ does not depend on $\mathbf{n}_0$. Then, Eq. (57) can be used when $\mathbf{p}_0 - e\mathbf{A}_0 = -\partial S(\mathbf{r}, \mathbf{n}_0, t)/\partial \mathbf{n}_0 = 0$ which means that $\mathbf{p}$ is a vector field. For a constant magnetic field, if $\mathbf{A} = (1/2)\mathbf{B}_0 \times \mathbf{r}$, this condition is not satisfied and Eqs. (50) cannot operate.

Let us point out that Lorentz equation and Faraday's law are automatically satisfied by Eqs. (50) [17]. So, we only have to find a solution which verifies the two following equations in order to derive the wave equation

$$ \nabla \times \mathbf{B} = \frac{1}{c^2} \frac{\partial \mathbf{E}}{\partial t} - \mu_0 \mathbf{N} \mathbf{e}, \quad (51) $$

$$ \nabla \cdot \mathbf{E} = -\frac{e}{\epsilon_0}(N - N_0), \quad (52) $$

where $N$ is the free electron density and $N_0$ is the density of ionized atoms which are assumed to be a uniform background.

Considering $v/c = p/p_0 \sim p/mc$ and neglecting terms in $(p/mc)^3$, Eqs. (50)-(52) give the following wave equation

$$ \begin{aligned} & \frac{1}{c^2} \frac{\partial^2 p}{\partial t^2} + \nabla(\nabla \cdot \mathbf{p}) - \Delta p + \frac{\omega_p^2}{c^2} p \\ & = -\frac{1}{2mc^2} \nabla \left( \frac{\partial p}{\partial t} \right)^2 - \frac{\mathbf{p}}{mc^2} \frac{\partial}{\partial t} \frac{\partial}{\partial z} p_z + \frac{\omega_p^2}{2c^2 m^2 c^2} \mathbf{p}^2 - \frac{\mathbf{p}}{2m^2 c^2} \frac{\partial^2}{\partial z^2} \mathbf{p}^2, \end{aligned} \quad (53) $$

with $\omega_p^2 = e^2 N_0 / \epsilon_0 m$ where $N_0$ is the ion density which is assumed to be at rest. The order of magnitude of the terms are

$$ 1; \frac{k^2 c^2}{\omega^2}; \frac{k^2 c^2}{\omega^2}; \frac{\omega_p^2}{\omega^2}; \frac{2kc}{\omega} \frac{p}{mc}; \frac{kc}{\omega} \frac{p}{mc}; \frac{\omega_p^2}{2\omega^2} \left(\frac{p}{mc}\right)^2; 2 \frac{k^2 c^2}{\omega^2} \left(\frac{p}{mc}\right)^2 $$

The electron momentum is assumed to be in the form

$$ \begin{aligned} p_i &= \sum_n p_i^{(n)} = \sum_n \lambda_i^{(n)} f_i^{(n)}(\theta), \quad i = x, z, \quad \theta = \omega t - kz, \\ &= p_i^{(0)} + p_i^{(1)} + p_i^{(2)} + \dots \end{aligned} \quad (54) $$

Each successive term of this sum is assumed to be an order of

magnitude smaller than the preceding term.

The following initial conditions were chosen

$$ p_x = p_z = \frac{dp_z}{d\theta} = 0, \quad \frac{dp_x}{d\theta} = -\frac{eE_{0x}}{\omega}, \quad (55) $$

when $\theta = 0$

The Lindstedt-Poincaré method is applied, in addition to (55), it is assumed that [21]

$$ \omega^2 = \omega^{(0)2} + \alpha^{(2)} + \dots, \quad (56) $$

where $\alpha^{(2)}$ is relatively small (in $v^2/c^2$).

To zero order ($v/c = 0$)

$$ \begin{aligned} & \frac{1}{c^2} \frac{\partial^2 p_x^{(0)}}{\partial t^2} - \frac{\partial^2}{\partial z^2} p_x^{(0)} + \frac{\omega_p^2}{c^2} p_x^{(0)} = 0, \\ & \frac{1}{c^2} \frac{\partial^2 p_z^{(0)}}{\partial t^2} + \frac{\omega_p^2}{c^2} p_z^{(0)} = 0. \end{aligned} \quad (57) $$

We seek a solution in the form

$$ \begin{aligned} p_x^{(0)} &= \lambda_x^{(0)} \sin\theta, \\ p_z^{(0)} &= \lambda_z^{(0)} \sin\theta. \end{aligned} \quad (58) $$

The initial conditions imply $\lambda_z^{(0)} = 0$ and $\lambda_x^{(0)} = -eE_{0x}/\omega$. Then, Eq. (50) give

$$ \begin{aligned} E_x &= E_{0x} \cos\theta, \\ E_z &= 0. \end{aligned} \quad (59) $$

Ignoring the term $(\alpha^{(2)}/c^2)p_x^{(0)}$ from Eq. (57), the following dispersion relation is obtained

$$ \omega^2 = \omega^{(0)2} = \omega_p^2 + k^2c^2. \quad (60) $$

Let us go to first order in $v/c$ now to calculate the density oscillation of the plasma sheet. In this approximation, the wave equation is

$$ \begin{aligned} & \frac{1}{c^2} \frac{\partial^2 p_x}{\partial t^2} - \frac{\partial^2}{\partial z^2} p_x + \frac{\omega_p^2}{c^2} p_x = - \frac{p_x}{mc^2} \frac{\partial^2}{\partial z \partial t} p_z, \\ & \frac{1}{c^2} \frac{\partial^2 p_z}{\partial t^2} + \frac{\omega_p^2}{c^2} p_z = - \frac{p_z}{mc^2} \frac{\partial}{\partial t} \frac{\partial}{\partial z} p_z - \frac{1}{2mc^2} \frac{\partial^2}{\partial z \partial t} (p_x^{(0)} + p_z^{(0)}). \end{aligned} \quad (61) $$

As $\lambda_z^{(0)} = 0$, this set of equations leads to

$$ \begin{aligned} & \frac{1}{c^2} \frac{\partial^2 p_x^{(1)}}{\partial t^2} - \frac{\partial^2}{\partial z^2} p_x^{(1)} + \frac{\omega_p^2}{c^2} p_x^{(1)} = 0, \\ & \frac{1}{c^2} \frac{\partial^2 p_z^{(1)}}{\partial t^2} + \frac{\omega_p^2}{c^2} p_z^{(1)} = -\frac{1}{2mc^2} \frac{\partial^2}{\partial z \partial t} (p_x^{(0)} - p_z^{(0)})^2 = \frac{ek}{c} E_{0x} \nu_0 \cos 2\theta. \end{aligned} \quad (62) $$

The solution is

$$ \begin{aligned} p_x^{(1)} &= 0, \\ p_z^{(1)} &= \frac{eE_{0x}}{\omega^2} ck \nu_0 \frac{\omega^2}{\omega_p^2 - 4\omega^2} \left[ \cos 2\theta - \cos \frac{\omega_p}{\omega} \theta \right]. \end{aligned} \quad (63) $$

As a consequence (Eq. (50)), the electric field of the wave is

$$ \begin{aligned} E_x &= E_{0x} \cos\theta, \\ E_z &= \frac{kc}{\omega} \left( \frac{\Omega_L^2}{4 - \Omega_L^2} \right) \nu_0 E_{0x} \left[ \frac{\sin\Omega_L \theta}{\Omega_L} - \frac{\sin 2\theta}{2} \right], \end{aligned} \quad (64) $$

where $\Omega_L = \omega_p/\omega$. The dispersion law is

$$ \omega^2 = \omega^{(0)2} = \omega_p^2 + k^2c^2 \quad (65) $$

The plasma density is given by Poisson's equation

$$ N = N_0 - (\epsilon_0 / e) \frac{\partial E_z}{\partial z}, \quad (66) $$

that is to say

$$ N = N_0 + \frac{\epsilon_0 k^2 c}{e \omega} \left( \frac{\Omega_L^2}{4 - \Omega_L^2} \right) \nu_0 E_{0x} [\cos \Omega_L \theta - \cos 2\theta]. \quad (67) $$

In order to calculate the refractive index, we follow the form of the source field and focus on the $\omega$-component of $\Delta E$. As $\eta = N \Delta z$, the $\Delta E_\omega$
created by the plasma is

$$ \Delta E_{\omega}^{b} = -\Delta z \frac{E_{0x}}{2c\omega} \omega_{p}^{2} \left[ 1 - \frac{3}{8}\varepsilon + \varepsilon \frac{k^{2}c^{2}}{2\omega^{2}} \left( \frac{1}{4 - \Omega_{L}^{2}} \right) \right] \sin\omega \left( t - \frac{z}{c} \right), \quad (68) $$

as $E_0 = E_{0x}$.

Comparing this value of $\Delta E_{\omega}$ to the one given by (24), we obtain

$$ n = 1 - \frac{1}{2} \frac{\omega_p^2}{\omega^2} \left\{ 1 - \left[ \frac{3}{8} - \frac{\omega^2 - \omega_p^2}{2(4\omega^2 - \omega_p^2)} \right] \varepsilon \right\}. \quad (69) $$

### Fully plasma approach

The index has also been calculated through a purely plasma approach keeping terms in $\nu^2/c^2$. To second order in $\nu_0$, taking into consideration the term $(\alpha^{(2)}/c^2)p_x^{(0)}$ which was neglected in the non-relativistic approximation, the wave Eq. (53) becomes

$$ \begin{aligned} \frac{1}{c^2} \frac{\partial^2 p_x^{(2)}}{\partial t^2} - \frac{\partial^2}{\partial z^2} p_x^{(2)} + \frac{\omega_p^2}{c^2} p_x^{(2)} - \frac{\alpha^{(2)}}{c^2} p_{x0} &= - \frac{p_x^{(0)}}{mc^2} \frac{\partial^2}{\partial z \partial t} p_z^{(1)} + \frac{\omega_p^2}{c^2} \frac{p_x^{(0)3}}{2m^2c^2} \\ &\quad - \frac{\partial^2}{\partial z^2} (p_x^{(0)2} + p_z^{(0)2}) \frac{p_x^{(0)}}{2m^2c^2} \\ &= \frac{\lambda_x^{(0)3}}{m^2c^2} \left( \frac{2k^2}{\Omega_L^2 - 4} + \frac{3\omega_p^2}{8c^2} + \frac{k^2}{2} \right) \sin\theta \\ &\quad - \left( \frac{2k^2}{\Omega_L^2 - 4} + \frac{\omega_p^2}{8c^2} + \frac{k^2}{2} \right) \sin 3\theta \\ &\quad + \frac{1}{2} \frac{k^2\Omega_L^2}{\Omega_L^2 - 4} [\sin(1 + \Omega_L)\theta \\ &\quad + \sin(1 - \Omega_L)\theta)]. \end{aligned} \quad (70a) $$

$$ \frac{1}{c^2} \frac{\partial^2 p_z^{(2)}}{\partial t^2} + \frac{\omega_p^2}{c^2} p_z^{(2)} = 0. \quad (70b) $$

In order not to have a secular term, $\alpha^{(2)}$ is set to eliminate the term in $\sin\theta$ in the right-hand side of Eq. (70a). We have

$$ \alpha^{(2)} = -\varepsilon \left( \frac{2c^2k^2}{\Omega_L^2 - 4} + \frac{3\omega_p^2}{8} + \frac{c^2k^2}{2} \right). \quad (71) $$

As the dispersion law is $\omega^2 = \omega^{(0)2} + \alpha^{(2)} = \omega_p^2 + k^2c^2 + \alpha^{(2)}$;

$$ \omega^2 = \omega_p^2 \left(1 - \frac{3}{8}\varepsilon\right) + k^2c^2 \left\{ 1 - \varepsilon \left[ \frac{\omega_p^2}{2(\omega_p^2 - 4\omega^2)} \right] \right\}. \quad (72) $$

Then,

$$ n = \frac{ck}{\omega} = 1 - \frac{1}{2} \frac{\omega_p^2}{\omega^2} \left\{ 1 - \left[ \frac{3}{8} - \frac{\omega^2 - \omega_p^2}{2(4\omega^2 - \omega_p^2)} \right] \varepsilon \right\}, \quad (73) $$

which is the same result as the one given by equation (69), in the low plasma density limit. It is also the same result as the one previously derived by C.E. Max and J. Arons [6].

The momentum and the electric field have also been calculated (Appendix C).

### Calculation of the index of refraction in the case of a circularly polarized plane wave propagating in plasma

The wave propagates in a gas which has been pre-ionized

Ionization takes place during the propagation of the wave, it is ignored in this part. The optical field is given by

$$ \begin{aligned} E &= E_0 [\cos(\omega t - kz) \hat{e}_x + \sin(\omega t - kz) \hat{e}_y], \\ B &= \frac{kE_0}{\omega} [-\sin(\omega t - kz) \hat{e}_x + \cos(\omega t - kz) \hat{e}_y]. \end{aligned} \quad (74) $$

Considering the scalar potential is zero, the corresponding wave

vector is

$$ \mathbf{A} = -\frac{E_0}{\omega} [\sin(\omega t - kz) \hat{e}_x - \cos(\omega t - kz) \hat{e}_y]. \quad (75) $$

In this case $p_z = 0$. As $A_z = 0$, we have $P_z = 0$.

The Hamiltonian of one electron in this wave reads

$$ H = \sqrt{\left[P_x - \frac{eE_0}{\omega}\sin(\omega t - kz)\right]^2 + \left[P_y + \frac{eE_0}{\omega}\cos(\omega t - kz)\right]^2 + m^2c^4}, \quad (76) $$

where $P_x$ and $P_y$ are two constants of motion. It is assumed, in this paragraph, that $P_x = P_y = 0$, as $\gamma$ is a constant. The electron is at rest on average in the $(x, y)$ plane. The electron density fluctuations will not be considered as there is no longitudinal component of the electron motion when the wave propagates in plasma.

The relativistic case is approached first. The wave is still assumed to be passing through a thin layer of plasma. $\Delta E_{\omega}^{a}$, which is the field added to the source caused by the oscillating electrons, can be calculated considering the wave form (74) and the additional time necessary to travel the distance $\Delta z$ in the plasma. We find for the two components of $\Delta E_{\omega}^{a}$

$$ \begin{aligned} \Delta E_{\omega x}^{a} &= \frac{\omega(n-1)\Delta z}{c} E_0 \sin\omega \left(t - \frac{z}{c}\right), \\ \Delta E_{\omega y}^{a} &= -\frac{\omega(n-1)\Delta z}{c} E_0 \cos\omega \left(t - \frac{z}{c}\right). \end{aligned} \quad (77) $$

Let us calculate now $\Delta E_{\omega}^{b}$, with $\Delta E_{\omega}$ calculated by using “Feynman’s method”. The equations of motion of one electron in the plasma layer are

$$ \begin{aligned} \frac{dp_x}{dt} &= m\dot{\gamma}v_x + m\gamma a_x = -eE_0 \cos\omega t, \\ \frac{dp_y}{dt} &= m\dot{\gamma}v_y + m\gamma a_y = -eE_0 \sin\omega t, \\ \frac{dp_z}{dt} &= 0. \end{aligned} \quad (78) $$

One has $p = m\gamma\mathbf{v}$ and $\nu^2 = p^2/m^2\gamma^2$ as a consequence $\gamma = \sqrt{1 + p^2/m^2c^2} = \sqrt{1 + \nu_0^2}$ and $\dot{\gamma} = 0$. Then

$$ \begin{aligned} a_x &= -\frac{eE_0}{m\sqrt{1+\nu_0^2}} \cos\omega t, \\ a_y &= -\frac{eE_0}{m\sqrt{1+\nu_0^2}} \sin\omega t. \end{aligned} \quad (79) $$

Far enough away from the moving charge, the field is given by [5]

$$ \begin{aligned} \delta^2 E_x &= \frac{e}{4\pi\epsilon_0 c^2} \frac{a_x(t - \frac{r}{c})}{r} = -\frac{e}{4\pi\epsilon_0 c^2 m} \frac{1}{\sqrt{1 + \nu_0^2}} \frac{E_0}{r} \cos\omega(t - \frac{r}{c}), \\ \delta^2 E_y &= \frac{e}{4\pi\epsilon_0 c^2} \frac{a_y(t - \frac{r}{c})}{r} = -\frac{e}{4\pi\epsilon_0 c^2 m} \frac{1}{\sqrt{1 + \nu_0^2}} \frac{E_0}{r} \sin\omega(t - \frac{r}{c}). \end{aligned} \quad (80) $$

The field created by all the electrons of a ring are considered first (Fig. 3), then, an integration leads to the total field created by all the electron of the plasma layer

$$ \begin{aligned} \Delta E_{\omega x}^{b} &= \int \delta E 2\pi \rho d\rho = -\frac{\eta e^2}{4\pi\epsilon_0 c^2} \frac{E_0}{m} \frac{1}{\sqrt{1+\nu_0^2}} \int_0^\infty \cos\omega(t-r/c) \frac{2\pi\rho}{r} d\rho, \\ \Delta E_{\omega y}^{b} &= \int \delta E 2\pi \rho d\rho = \frac{\eta e^2}{4\pi\epsilon_0 c^2} \frac{E_0}{m} \frac{1}{\sqrt{1+\nu_0^2}} \int_0^\infty \sin\omega(t-r/c) \frac{2\pi\rho}{r} d\rho. \end{aligned} \quad (81) $$

The components of the resulting field are

$$ \begin{aligned} \Delta E_{\omega x}^{b} &= -\frac{\eta e^2 E_0}{2\epsilon_0 c m\omega} \frac{1}{\sqrt{1+\nu_0^2}} \sin\omega(t-z/c), \\ \Delta E_{\omega y}^{b} &= -\frac{\eta e^2 E_0}{2\epsilon_0 c m\omega} \frac{1}{\sqrt{1+\nu_0^2}} \cos\omega(t-z/c). \end{aligned} \quad (82) $$

We must have $\Delta E_{\omega}^{a} = \Delta E_{\omega}^{b}$, consequently

$$ n = 1 - \frac{1}{2\sqrt{1 + \varepsilon \omega^2}} \frac{\omega_p^2}{\omega^2}. \quad (83) $$
This result is in good agreement with the result previously found by A. I. Akhiezer and B.V. Polovin [22], by P. K. Kaw and J. Dawson [23] and by N.L. Tsintsadze and D.D. Tskhakaya [24] as part of a plasma approach. Thus, we have highlighted the simplicity of our method when we compare it to the different plasma approaches.

Let us point out that this expression of n is an exact solution that is accurate for high values of ε. The plasma effect is canceled at high intensity. In the non-relativistic limit, we find the same result as in the linear polarization case (Eq. (30a)).

Ionization which takes place during the propagation of the wave in the atmosphere is considered

Two channels of field ionization exist: multiphoton and tunneling. They are distinguished by the Keldysh parameter [7]

$$g = \sqrt{\frac{I_p}{2U_p}}, \qquad (84)$$

where $I_p$ is the ionization potential and $U_p$ is the ponderomotive energy. We chose $I_p = 12$ eV for air. As $U_p(eV) = 9.33 \times 10^{-14}I(W/cm^2)\lambda^2(\mu m)$, for $\lambda = 10.6 \mu m$ and $I \approx 5.2 \times 10^{11}W/cm^2$, we have: $g \approx 1$. If $I \approx 6 \times 10^{12}W/cm^2$, then $\approx 0.3 < 1$. When $I \approx 1 \times 10^{14}W/cm^2$, Eq. (84) yields $g \approx 7 \times 10^{-2} < 1$. As a consequence, tunneling ionization is prevailing.

The terms due to the relativistic mass of the electrons are considered. Feynman's method is used to estimate roughly the effect of ionization on the index of refraction. It is assumed that ionization is achieved instantaneously, and tunneled electrons have one average direction only. Ionization takes place in a very thin air layer before the optical wave is diffracted by the plasma. It starts at $t = 0$ and $z = 0$ and propagates with the plane wave. It is assumed here that, in the $\Delta z$ thick layer, the electron density is constant in the part traversed by the wave. These not very realistic assumptions are used to make a simple calculation and start exploring the effect of ionization on the index.

At the moment of ionization, the tunneled electron has zero velocity [8–10]. The optical field is assumed to be given by Eq. (74). The Hamiltonian of one electron in the wave is given by

$$H = \sqrt{\left[P_x - \frac{eE_0}{\omega}\sin(\omega t - kz)\right]^2 c^2 + \left[P_y + \frac{eE_0}{\omega}\cos(\omega t - kz)\right]^2 c^2 + p_z^2 c^2 + m^2 c^4}. \qquad (85)$$

Assuming that $P_x = 0$, the constant $P_y$ is set to $-eE_0/\omega$ so that the electron may have zero velocity at the moment of ionization. Thus, the momentum of this electron is given by

$$\begin{aligned} p_x &= -\frac{eE_0}{\omega}\sin(\omega t - kz), \\ p_y &= \frac{eE_0}{\omega}\cos(\omega t - kz) + P_y = \frac{eE_0}{\omega}\cos(\omega t - kz) - \frac{eE_0}{\omega}, \end{aligned} \qquad (86)$$

The equations of motion in the (x, y) plane are

$$\begin{aligned} \frac{dp_x}{dt} &= m\dot{\gamma}v_x + m\gamma a_x = -eE_0\cos(\omega t - kz) + eE_0 \frac{k}{\omega}v_z\cos(\omega t - kz), \\ \frac{dp_y}{dt} &= m\dot{\gamma}v_y + m\gamma a_y = -eE_0\sin(\omega t - kz) + eE_0 \frac{k}{\omega}v_z\sin(\omega t - kz). \end{aligned} \qquad (87)$$

An electron trajectory is shown in Fig. 6.

Fig. 7 shows a laser pulse whose envelope is symmetrical with respect to its direction of propagation. Electrons are ionized and they acquire an average velocity parallel to the y-axis.

We consider the average velocity along the y-axis, $V = \bar{v}_y$, then $V/c = -v_0$. We consider a Galilean frame ($L^*$) propagating in this direction with this velocity. In order to cancel this drift motion, a Lorentz transformation is achieved. Indeed, in the frame in translation, the drift velocity $V$ is transformed into $V' = 0$ according to the relativistic law transformation of velocity: $v' = (v - V)/(1 - vV/c^2)$ where $v'$ is the result of the Lorentz transformation on $v$ (everywhere in this paragraph where $g'$ is g after the Lorentz transformation). In ($L^*$) electrons can still

<figure><img src="image_8.png" /><figcaption>Fig. 6. Trajectory of one electron calculated ignoring the motion in the direction of propagation of the wave. $v_0 = 0.01$, $\omega = 1.77 \times 10^{14} s^{-1}$.</figcaption></figure>

<figure><img src="image_9.png" /><figcaption>Fig. 7. The average ionization current $j_i$ created along the y-axis.</figcaption></figure>

be considered to be locked in a ring. They move a small distance only and the delay time can be assumed to be almost constant. Far away from the ring one has $r' \approx r'_{ret}$ the subscript ret means that the quantity is evaluated at the retarded time. In this part, we avoid the premium symbol for these quantities so as not to confuse them with the quantities having undergone a transformation of Lorentz. Then, $r$ can be taken out the second derivative in the expression of the transverse electric field created by one electron [5]

$$\begin{aligned} \delta^2 E'_x &= \frac{e}{4\pi\epsilon_0 c^2} \frac{d^2}{dt'^2} \left(\frac{x'}{r'}\right)_{\text{ret}} = \frac{e}{4\pi\epsilon_0 c^2 r'} \frac{d^2 x'_{\text{ret}}}{dt'^2}, \\ \delta^2 E'_y &= \frac{e}{4\pi\epsilon_0 c^2} \frac{d^2}{dt'^2} \left(\frac{y'}{r'}\right)_{\text{ret}} = \frac{e}{4\pi\epsilon_0 c^2 r'} \frac{d^2 y'_{\text{ret}}}{dt'^2}. \end{aligned} \qquad (88)$$

As $\omega t - kz = \omega\Gamma[t' + (V/c^2)y'] - kz'$, the Lorentz transformations for the components of the optical field are

$$\begin{aligned} E'_{x} &= \Gamma(E_{x} + VB_{z}) = \frac{E_{0}}{\sqrt{1-v_{0}^{2}}} \cos(\omega't' - k'_{\perp}y' - kz'), \\ E'_{y} &= E_{y} = E_{0}\sin(\omega't' - k'_{\perp}y' - kz'), \\ E'_{z} &= \Gamma(E_{z} - VB_{x}) = -E_{0}v_{0}\sin(\omega't' - k'_{\perp}y' - kz'), \\ B'_{x} &= \Gamma(B_{x} - \frac{V}{c^2}E_{z}) = -\frac{E_{0}}{c} \frac{1}{\sqrt{1-v_{0}^{2}}} \sin(\omega't' - k'_{\perp}y' - kz'), \\ B'_{y} &= B_{y} = \frac{E_{0}}{c} \cos(\omega't' - k'_{\perp}y' - kz'), \\ B'_{z} &= \Gamma(B_{z} + \frac{V}{c^2}E_{x}) = -\frac{E_{0}}{c} v_{0} \cos(\omega't' - k'_{\perp}y' - kz'), \end{aligned} \qquad (89)$$

where $\Gamma = (\sqrt{1 - V^2/c^2})^{-1} = (\sqrt{1 - v_0^2})^{-1}$, $\omega' = \omega\Gamma = \omega/\sqrt{1 - v_0^2}$ and $k'_{\perp} = -\omega\Gamma V/c^2 = v_0\omega/c$.

In ($L^*$), assuming the particles remain very close to the plane $z = 0$,
<figure><img src="image_10.png" /><figcaption>Fig. 8. Radiation field of a sheet of oscillating charges in (L*).</figcaption></figure>

the transverse equations of motion of one electron are

$$
\begin{align}
\frac{dp_x'}{dt'} &= m\dot{\gamma}'v_x + m\gamma'a'_x = -\frac{eE_0}{\sqrt{1-v_0^2}}\cos(\omega't' - k'_\perp y'), \nonumber \\
\frac{dp_y'}{dt'} &= m\dot{\gamma}'v_y + m\gamma'a'_y = -eE_0\sin(\omega't' - k'_\perp y'), \nonumber \\
\frac{dp_z'}{dt'} &= m\dot{\gamma}'v_z + m\gamma'a'_z = -eE_0v_0\sin(\omega't' - k'_\perp y'), \tag{90}
\end{align}
$$

where $\dot{f} = df/dt'$. Here, $\gamma' = \sqrt{1+p'^2/m^2c^2} \approx \sqrt{1+v_0^2}$ as terms in $v_0^3$
are neglected. The transverse acceleration of the electron reads (Fig. 8)

$$
\begin{align}
a'_{x} &= -\frac{eE_0}{m} \cos[\omega't' - k'_{\perp} \rho' \cos\varphi'], \nonumber \\
a'_{y} &= -\frac{eE_0}{m} \frac{1}{\sqrt{1+v_0^2}} \sin(\omega't' - k'_{\perp} \rho' \cos\varphi'). \tag{91}
\end{align}
$$

An infinitesimal surface $\delta^2 S' = \rho' d\rho' d\rho'$) in (L*) which is part of a
plasma ring (Fig. 8)

The transverse field created by this surface $\delta^2 S$, far away from the
plasma sheet is

$$
\begin{align}
\delta^2 E'_{x} &= \frac{\eta \delta^2 Se}{4\pi\epsilon_0 c^2 r'} \frac{d^2 x'_{ret}}{dt^2} = - \frac{\eta e^2 E_0}{4\pi\epsilon_0 c^2 m} \cos \left[ \omega' \left(t' - \frac{r'}{c}\right) - k'_{\perp} \rho' \cos \varphi' \right] \frac{\varphi'}{r'} d\varphi' d\rho' \nonumber \\
&= - \frac{\eta e^2 E_0}{4\pi\epsilon_0 c^2 m} \left[ \sum_{n=-\infty}^{n=\infty} J_n(k'_{\perp} \rho') \cos n\left(\varphi' + \frac{\pi}{2}\right) \cos \omega' \left(t' - \frac{r'}{c}\right) \right. \nonumber \\
&\qquad \left. + \sum_{n=-\infty}^{n=\infty} J_n(\rho') \sin n\left(\varphi' + \frac{\pi}{2}\right) \sin \omega' \left(t' - \frac{r'}{c}\right) \right] \frac{\varphi'}{r'} d\varphi' d\rho', \nonumber \\
\delta^2 E'_{y} &= \frac{\eta \delta^2 Se}{4\pi\epsilon_0 c^2 r'} \frac{d^2 y'_{ret}}{dt^2} = - \frac{\eta e^2 E_0}{4\pi\epsilon_0 c^2 m} \frac{1}{\sqrt{1+v_0^2}} \sin \left[ \omega' \left(t' - \frac{r'}{c}\right) \right. \nonumber \\
&\qquad \left. - k'_{\perp} \rho' \cos \varphi' \right] \frac{\varphi'}{r'} d\varphi' d\rho' \nonumber \\
&= - \frac{\eta e^2 E_0}{4\pi\epsilon_0 c^2 m} \frac{1}{\sqrt{1+v_0^2}} \left[ \sum_{n=-\infty}^{n=\infty} J_n(k'_{\perp} \rho') \cos n\left(\varphi' + \frac{\pi}{2}\right) \sin \omega' \left(t' - \frac{r'}{c}\right) \right. \nonumber \\
&\qquad \left. - \sum_{n=-\infty}^{n=\infty} J_n(k'_{\perp} \rho') \sin n\left(\varphi' + \frac{\pi}{2}\right) \cos \omega' \left(t' - \frac{r'}{c}\right) \right] \frac{\varphi'}{r'} d\varphi' d\rho', \tag{92}
\end{align}
$$

where the $J_n$ are first kind Bessel functions. The field created by a ring is
obtained by integrating according to $\varphi'$

$$
\begin{align}
\delta E'_{x} &= -\frac{\eta e^2 E_0}{2\epsilon_0 c^2 m} J_0(k'_{\perp} \rho') \cos \left[ \omega' \left( t' - \frac{r'}{c} \right) \right] \frac{\varphi'}{r'} d\rho', \nonumber \\
\delta E'_{y} &= -\frac{\eta e^2 E_0}{2\epsilon_0 c^2 m} \frac{J_0(k'_{\perp} \rho')}{\sqrt{1 + \nu_0^2}} \sin \left[ \omega' \left( t' - \frac{r'}{c} \right) \right] \frac{\varphi'}{r'} d\rho'. \tag{93}
\end{align}
$$

In order to obtain the total field from all the charges of the plasma
sheet we integrate these values over all $r'$ using again the relation
$r'^2 = \rho'^2 + z'^2$. To do this, both integrations were performed in an ap-
proximate way by assuming that $J_0(X) \approx 1 - X^2/4$ for low values
of X (Appendix D). We focused on the term in $\sin[\omega'(t' - z'/c)]$ and in
$\cos[\omega'(t' - z'/c)]$ when integrating respectively $\delta E'_{x}$ and $\delta E'_{y}$

$$
\begin{align}
\Delta E'_{x} &= - \frac{\eta e^2 E_0}{2\epsilon_0 c^2 m} \int_{z}^{z_f} J_0(k'_{\perp} \rho') \cos \left[ \omega' \left( t' - \frac{r'}{c} \right) \right] dr', \nonumber \\
&= - \frac{\eta e^2 E_0}{2\epsilon_0 c^2 m} \int_{z}^{z_f} \left( 1 - \frac{k'^2}{4} (r'^2 - z'^2) \right) \cos \left[ \omega' \left( t' - \frac{r'}{c} \right) \right] dr', \tag{94}
\end{align}
$$

where $z_f$ is the value of $r$ for which the Bessel function is zero.

As the magnetic field is given by **B** = −**e**r' × **E**/c, the field compo-
nents are

$$
\begin{align}
\Delta E'_{x} &= -\frac{\eta e^2 E_0}{2\epsilon_0 c m \omega} \sqrt{1 - \nu_0^2} \sin \left[ \omega' \left( t' - \frac{z'}{c} \right) \right], \nonumber \\
\Delta E'_{y} &= \frac{\eta e^2 E_0}{2\epsilon_0 c m \omega} \sqrt{\frac{1 - \nu_0^2}{1 + \nu_0^2}} \cos \left[ \omega' \left( t' - \frac{z'}{c} \right) \right], \quad \Delta E'_{z} = 0, \nonumber \\
\Delta B'_{x} &= -\frac{\eta e^2 E_0}{2\epsilon_0 c^2 m \omega} \sqrt{\frac{1 - \nu_0^2}{1 + \nu_0^2}} \cos \left[ \omega' \left( t' - \frac{z'}{c} \right) \right], \nonumber \\
\Delta B'_{y} &= -\frac{\eta e^2 E_0}{2\epsilon_0 c^2 m \omega} \sqrt{1 - \nu_0^2} \sin \left[ \omega' \left( t' - \frac{z'}{c} \right) \right], \quad \Delta B'_{z} = 0. \tag{95}
\end{align}
$$

Then, in the laboratory frame

$$
\begin{align}
\Delta E_x^b &= \Gamma(\Delta E'_x - V \Delta B'_z) = -\frac{\eta e^2 E_0}{2\epsilon_0 c m \omega} \sin\left[\omega\left(t - \frac{z}{c}\right)\right], \nonumber \\
\Delta E_y^b &= \Delta E'_y = \frac{\eta e^2 E_0}{2\epsilon_0 c m \omega} \left(1 - \nu_0^2\right) \cos\left[\omega\left(t - \frac{z}{c}\right)\right]. \tag{96}
\end{align}
$$

In the x-direction we have $\Delta E_x^a = \omega(n_x - 1)(\Delta z/c)E_0\sin\omega(t - z/c)$
and we must satisfy $\Delta E_x^a = \Delta E_x^b$. It implies

$$
n_x = 1 - \frac{1}{2} \frac{\omega_p^2}{\omega^2}. \tag{97}
$$

In the y-direction we have $\Delta E_y^a = -\omega(n_y - 1)(\Delta z/c)E_0\cos\omega(t - z/c)$, we must verify $\Delta E_y^a = \Delta E_y^b$. We obtain

$$
n_y = 1 - \frac{1}{2} \frac{\omega_p^2}{\omega^2} (1 - \nu_0^2). \qquad (98)
$$

Given the very strong hypothesis we made on the direction of the ionization current, only the refractive index, $n_y$, has a physical meaning. The relativistic correction due to ionization is in the y - direction, it is different from the one we have when the plasma is assumed to be pre-ionized (Eq. (83)). As the average velocity of the ionized electrons is in fact in all the directions of the (x, y) plane, we can conjecture that the index is isotropic and has a relativistic correction.

**Conclusions**

The key role played by the index of refraction in the propagation of
waves is well known. For instance, a good knowledge of the analytical
form of the index is necessary when the propagation of a wave in the
atmosphere is studied by using numerical codes. The fact that we are
interested in long wavelengths makes the threshold power very high. At
these powers, the ionized electrons of the atmospheres can become
relativistic. Special attention is given to this kind of non-linear situation
by applying a very intuitive method.

In this article, the refraction index is calculated mainly by enforcing
the physical approach described in the courses of R. P. Feynman. This
method gives insight to the characteristics origin fields which added to
the external source shifts its phase. The originality of this work lies in
applying this method to non-linear situations by comparing the result
obtained to the one obtained by plasma approaches before applying it
to a not yet studied problem. We insist on the simplicity of the method
although we show how to simplify the perturbation calculation in a
plasma model. The method was applied to the situation where the
electrons are relativistic due to the high intensity of the optical field.
First, when the electromagnetic field which propagates in the atmo-
sphere is linearly polarized, ionization is not considered. When the
relativistic electron mass is considered two approaches are im-
plemented. The first one consists in coupling the “Feynman’s method”
to a plasma approach, the second is a fully plasma approach. The good agreement between the two results confirms our interest of the first method which is much simpler to implement.

Then, the method is applied to the case of a circularly polarized optical wave while ionization and relativistic effects are considered. Situations corresponding to a Keldysh parameter lower than unity are considered. As a consequence, tunneling ionization dominates. In order to apply “Feynman’s method” conveniently, it was assumed that ionized electrons are ejected in one specific direction. Simplifying physics by

making a not very realistic assumption, a simple calculation allowed us to start exploring the effect of ionization on the index.

## Declaration of Competing Interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

## Appendix A

### Electric field created by an oscillating charge derived by using the lienard and Wiechert solution

The wave electric field emitted by a moving electron along the x-axis at an observation point P (Fig. A1) is calculated.

<figure><img src="image_11.png" /><figcaption>Fig. A1. Finding the electric field created by a moving electron at P.</figcaption></figure>

The 4-potential produced by a moving electron (Fig. A1) is given by the Liénard and Wiechert solution [14,25]

$$ 
\begin{aligned}
A &= -\frac{1}{4\pi\epsilon_0 c^2} \frac{ev_{ret}}{\left(r - \frac{v \cdot r}{c}\right)_{ret}} \approx -\frac{ev_{ret}}{4\pi\epsilon_0 c^2 r_{ret}} \left(1 + \frac{v_{ret} \cdot \tilde{\eta}_{ret}}{c}\right), \\
\phi &= -\frac{1}{4\pi\epsilon_0} \frac{e}{\left(r - \frac{v \cdot r}{c}\right)_{ret}} \approx -\frac{e}{4\pi\epsilon_0 \tilde{r}_{ret}} \left(1 + \frac{v_{ret} \cdot \tilde{\eta}_{ret}}{c}\right),
\end{aligned}
\quad (A1) $$

where the subscript “ret” means that the quantity in the brackets is evaluated at a retarded time and $\eta_{ret} = (r/r)_{ret}$ a unit vector along the $r_{ret}$ direction. As the electron travels a short distance, then $r_{ret} \approx r$. The term $(v \cdot \eta/c)_{ret}$ is assumed to be small as its trajectory remains far from P, the velocity becomes perpendicular to $r$. As $\partial[v \cdot (t - r/c)]/\partial r = -(v \cdot \eta)/c$, $\partial\phi/\partial r$ and $\partial A/\partial t$ are given by

$$ 
\begin{aligned}
\frac{\partial \phi}{\partial r} &= -\frac{ev}{4\pi\epsilon_0} \left[ -\frac{1}{r^2} \left(1 + \frac{v_{ret} \cdot \eta}{c}\right) - \frac{v_{ret} \cdot \eta}{c^2 r} \right], \\
\frac{\partial A}{\partial t} &= -\frac{e}{4\pi\epsilon_0 c^2 r} \left[ \dot{v}_{ret} \left(1 + \frac{v_{ret} \cdot \eta}{c}\right) + \frac{v_{ret}}{c} \dot{v}_{ret} \cdot \eta \right].
\end{aligned}
\quad (A2) $$

Thus the field is

$$ 
\begin{aligned}
E &= -\nabla\phi, -\frac{\partial A}{\partial t} \\
&= \frac{e}{4\pi\epsilon_0} \left\{ \left[ -\frac{1}{r^2} \left(1 + \frac{v_{ret} \cdot \eta}{c}\right) - \frac{v_{ret} \cdot \eta}{c^2 r} \right] \eta + \frac{1}{c^2 r} \left[ \dot{v}_{ret} \left(1 + \frac{v_{ret} \cdot \eta}{c}\right) + \frac{v_{ret}}{c} \dot{v}_{ret} \cdot \eta \right] \right\}.
\end{aligned}
\quad (A3) $$

Far away from the plate terms with $\dot{v}_{ret} \cdot \eta$ can also be neglected. As $\dot{v}_{ret} \approx a(t - r/c)$, the electric field is approximated by

$$ E = -\frac{e}{4\pi\epsilon_0} \left( \frac{\eta}{r^2} - \frac{1}{c^2 r} \dot{v}_{ret} \right) \approx \frac{e}{4\pi\epsilon_0 c^2} \frac{a\left(t - \frac{r}{c}\right)}{r} \hat{e}_x. \quad (A4) $$

The electrostatic term in (A4) can be neglected for large values of $r$. Thus, the electric field expression given by Eq. (12) is found again. This confirms that, if the observation position is far from the plate, the field expression given by Eq. (12) constitutes a good approximation

## Appendix B

### Electric field created by an oscillating dipole close to the charges

The ionized electron and its parent ion are considered here. The electron driven by the optical field oscillates around the ion which remains at rest. We consider the field created by an oscillating dipole [8,25] (Fig. B1).

The 4-potential produced by the moving charge $q$ (Fig. B1) is given by the Liénard and Wiechert solution [14,25]

$$ 
\begin{aligned}
\mathbf{A} &= \frac{1}{4\pi\epsilon_0 c^2} \frac{q v_{ret}}{\left(r - \frac{v \cdot r}{c}\right)_{ret}} \approx \frac{q v_{ret}}{4\pi\epsilon_0 c^2 \tilde{r}_{ret}} \left(1 + \frac{v_{ret} \cdot \tilde{\eta}_{ret}}{c}\right), \\
\phi &= \frac{1}{4\pi\epsilon_0} \frac{q}{\left(r - \frac{v \cdot r}{c}\right)_{ret}} \approx \frac{q}{4\pi\epsilon_0 \tilde{r}_{ret}} \left(1 + \frac{v_{ret} \cdot \tilde{\eta}_{ret}}{c}\right),
\end{aligned}
\quad (B1) $$

where the subscript “ret” is defined in Appendix A and $\tilde{\eta}_{ret} = \tilde{r}_{ret}/\tilde{r}_{ret}$. We have considered that the term $(v \cdot \tilde{\eta}/c)_{ret}$ is a small relativistic correction. Relativistic terms will be ignored further. The point of observation P is assumed to be very close to the dipole, then $r_{ret} \approx r$ as r/c can be neglected. Thus, the field in P is
<figure><img src="image_12.png" /><figcaption>Fig. B1. Finding the electric field created by an oscillating dipole.</figcaption></figure>

$$ 
\begin{aligned}
E_{dip} &= -\nabla(\phi + \phi^+) - \frac{\partial A}{\partial t}, \\
&= \frac{e}{4\pi\epsilon_0} \left\{ \left[ -\frac{1}{\bar{r}^2} \left( 1 + \frac{v_{\text{ret}} \cdot \tilde{\eta}}{c} \right) - \frac{v_{\text{ret}} \cdot \tilde{\eta}}{c^2 \bar{r}} \right] \tilde{\eta} + \frac{1}{c^2 \bar{r}} \left[ \dot{v}_{\text{ret}} \left( 1 + \frac{v_{\text{ret}} \cdot \tilde{\eta}}{c} \right) + \frac{v_{\text{ret}}}{c} \dot{v}_{\text{ret}} \cdot \tilde{\eta} \right] + \frac{1}{\bar{r}_+^3} \tilde{\eta}^+ \right\},
\end{aligned}
\quad (B2) $$

with $\tilde{\eta}^+ = \tilde{r}_{\text{ret}}^+ / \tilde{r}_{\text{ret}}^+$. The approximations which have just been described lead to the following electric field expression due to the dipole in P.

$$ E_{dip} \approx \frac{e}{4\pi\epsilon_0} \left( \frac{\tilde{\eta}^+}{\bar{r}^{+2}} - \frac{\tilde{\eta}}{\bar{r}^2} \right), \quad (B3) $$

Considering that the distance $d$ between the two charges is a small quantity and as $\tilde{r} \approx \tilde{r}^+ \approx r = \sqrt{x^2 + y^2 + z^2}$, we can write $\tilde{r}^2 = (x - d/2)^2 + y^2 + z^2 \approx r^2 - xd$ and $\tilde{r}^{+2} = (x + d/2)^2 + y^2 + z^2 \approx r^2 + xd$. When P is on the x-axis, the x-component of the electric field is given by

$$ E_{dip} \approx \frac{1}{2\pi\epsilon_0 r^3} \hat{e}_x, \quad (B4) $$

where $\hat{e}_x$ is a unit vector in the direction of the x-axis and $p$ is the dipole moment: $p = -ed$. This electric field adds to the incoming optical field. The dipole moment per volume unit of the plasma P is given by $P = Np = -Ned$ [25,26] and is proportional to $E_{dip}$

$$ 
\begin{aligned}
E_{dip} &= \alpha \frac{P}{\epsilon_0}, \\
\alpha \frac{P}{\epsilon_0} &= -\alpha \frac{Ned}{\epsilon_0} \hat{e}_x = -\frac{1}{2\pi\epsilon_0 r^3} \hat{e}_x,
\end{aligned}
\quad (B5) $$

where $\alpha$ is a constant coefficient which is set at $\alpha = 1/3$ [26,27]. Eq. (B5) is satisfied when $r^{-3} = 2\pi N/3$. The displacement of the charge is calculated from the equation of motion: $d = (eE_0/m\omega^2) \cos \omega t$. Then,

$$ 
\begin{aligned}
E_{dip} &= - \frac{e^2 N}{3\epsilon_0 m \omega^2} E_0 \cos \omega t, \\
&= - \frac{\omega_p^2}{3\omega^2} E_0 \cos \omega t.
\end{aligned}
\quad (B6) $$

We know that $P = (n^2 - 1)\epsilon_0 E_0 \cos \omega t$ where $n$ is the plasma index and $E_0$ the amplitude of the source. The following equality: $P = 3\epsilon_0 E_{dip}$ must be satisfied. Consequently, we have

$$ n^2 = 1 - \frac{\omega_p^2}{\omega^2}, \quad (B7) $$

or in the case of low-density plasma

$$ n = 1 - \frac{\omega_p^2}{2 \omega^2}, \quad (B8) $$

Finally, Eq. (30a) is found again by using a different approach.

## Appendix C

### Relativistic solution for the momentum and the field

In order to calculate the momentum, we let

$$ p_x^{(2)} = A_{\theta}^{(2)} \sin\theta + A_{3\theta}^{(2)} \sin3\theta + A_{(1+\Omega_L)\theta}^{(2)} \sin[(1 + \Omega_L)\theta] + A_{(1-\Omega_L)\theta}^{(2)} \sin[(1 - \Omega_L)\theta]. \quad (C1) $$

As $\alpha^{(2)}$ (Eq. (56)) is set to eliminate the $\sin\theta$ force term in the right-hand side of Eq. (70a). The amplitude $A_{3\theta}^{(2)}$ must satisfy

$$ \left( -\frac{\omega^{(0)2}}{c^2} \frac{\partial^2}{\partial t^2} + 9k^2 + \frac{\omega_p^2}{c^2} \right) A_{3\theta}^{(2)} \sin 3\theta = -\frac{\lambda_x^{(0)3}}{m^2 c^2} \left( \frac{2k^2}{\Omega_L^2 - 4} + \frac{\omega_p^2}{8c^2} + \frac{k^2}{2} \right) \sin 3\theta, \quad (C2) $$

this implies
$$A_{3\theta}^{(2)} = \frac{\lambda_x^{(0)3}}{m^2 c^2} \frac{1}{9 \left( \frac{\omega^2}{c^2} - k^2 \right) - \frac{\omega_p^2}{c^2}} \left( \frac{2k^2}{\Omega_L^2 - 4} + \frac{\omega_p^2}{8c^2} + \frac{k^2}{2} \right) \quad (C3)$$

Furthermore

$$A_{(1 \pm \Omega_L)\theta}^{(2)} = \frac{1}{2} \frac{\lambda_x^{(0)3}}{m^2 c^2} \frac{1}{(1 \pm \Omega_L)^2 \left(k^2 - \frac{\omega^2}{c^2}\right) + \frac{\omega_p^2}{c^2} \frac{k^2 \Omega_L^2}{\Omega_L^2 - 4}} \quad (C4)$$

with $A_{(a \pm b)}^{(2)} = A_{(a+b)}^{(2)} + A_{(a-b)}^{(2)}$.

The solution for $p$ becomes

$$P_x = \lambda_x^{(0)} \sin\theta + \frac{\lambda_x^{(0)3}}{m^2 c^2} \frac{1}{9 \left( \frac{\omega^2}{c^2} - k^2 \right) - \frac{\omega_p^2}{c^2}} \left( \frac{2k^2}{\Omega_L^2 - 4} + \frac{\omega_p^2}{8c^2} + \frac{k^2}{2} \right) \sin3\theta \\ + \frac{1}{2} \frac{\lambda_x^{(0)3}}{m^2 c^2} \frac{1}{(1 \pm \Omega_L)^2 \left(k^2 - \frac{\omega^2}{c^2}\right) + \frac{\omega_p^2}{c^2} \frac{k^2 \Omega_L^2}{\Omega_L^2 - 4}} \sin(1 \pm \Omega_L)\theta, \quad (C5a)$$

$$P_z = - \frac{\lambda_x^{(0)}}{\omega} c k \nu_0 \frac{\omega^2}{\omega_p^2 - 4\omega^2} \left[ \cos 2\theta - \cos \frac{\omega_p}{\omega} \theta \right] \quad (C5b)$$

The initial conditions lead to

$$\begin{aligned} p_x = & - \frac{e}{\omega} E_{0x} \sin\theta + \frac{e}{\omega} E_{0x} \nu_0^2 \left[ \frac{3}{9 \left( \frac{\omega^2}{c^2} - k^2 \right) - \frac{\omega_p^2}{c^2}} \left( \frac{2k^2}{\Omega_L^2 - 4} + \frac{\omega_p^2}{8c^2} + \frac{k^2}{2} \right) \right. \\ & \left. + \frac{1}{(1 \pm \Omega_L) \left( k^2 - \frac{\omega^2}{c^2} \right) + \frac{\omega_p^2}{c^2} \frac{k^2 \Omega_L^2}{\Omega_L^2 - 4}} \sin\theta \right. \\ & - \frac{e}{\omega} E_{0x} \nu_0^2 \frac{1}{9 \left( \frac{\omega^2}{c^2} - k^2 \right) - \frac{\omega_p^2}{c^2}} \left( \frac{2k^2}{\Omega_L^2 - 4} + \omega_p^2 + \frac{k^2}{2} \right) \sin3\theta \\ & - \frac{1}{2} \frac{e}{\omega} E_{0x} \nu_0^2 \frac{1}{(1 \pm \Omega_L)^2 \left( k^2 - \frac{\omega^2}{c^2} \right) + \frac{\omega_p^2}{c^2} \frac{k^2 \Omega_L^2}{\Omega_L^2 - 4}} \sin(1 \pm \Omega_L)\theta, \end{aligned} \quad (C6a)$$

$$P_z = \frac{e E_{0x}}{\omega} c k \nu_0 \frac{\omega^2}{\omega_p^2 - 4\omega^2} \left[ \cos 2\theta - \cos \frac{\omega_p}{\omega} \theta \right] \quad (C6b)$$

Eq. (50) give the expression of the wave electric field

$$\begin{aligned} E_x = & E_{0x} \cos\theta - \nu_0^2 E_{0x} \left[ \frac{3}{9 \left( \frac{\omega^2}{c^2} - k^2 \right) - \frac{\omega_p^2}{c^2}} \left( \frac{2k^2}{\Omega_L^2 - 4} + \frac{\omega_p^2}{8c^2} + \frac{k^2}{2} \right) \right. \\ & \left. + \frac{1}{(1 \pm \Omega_L) \left( k^2 - \frac{\omega^2}{c^2} \right) + \frac{\omega_p^2}{c^2} \frac{k^2 \Omega_L^2}{\Omega_L^2 - 4}} \sin(1 \pm \Omega_L) \right] \cos\theta \\ & + \nu_0^2 E_{0x} \frac{1}{3 \left( \frac{\omega^2}{c^2} - k^2 \right) - \frac{\omega_p^2}{c^2}} \left( \frac{2k^2}{\Omega_L^2 - 4} + \frac{\omega_p^2}{8c^2} + \frac{k^2}{2} \right) \cos3\theta \\ & + \frac{1}{2} \nu_0^2 E_{0x} \frac{1}{(1 \pm \Omega_L) \left( k^2 - \frac{\omega^2}{c^2} \right) + \frac{\omega_p^2}{c^2} \frac{k^2 \Omega_L^2}{\Omega_L^2 - 4}} \cos(1 \pm \Omega_L) \theta, \end{aligned} \quad (C7a)$$

$$E_z = \frac{kc}{\omega} \left( \frac{\Omega_L^2}{4 - \Omega_L^2} \right) \nu_0 E_{0x} \left[ \frac{\sin \Omega_L \theta}{\Omega_L} - \frac{\sin 2\theta}{2} \right] \quad (7b)$$

## Appendix D

Approximations used in integrating Eq. (93)

Let us justify the very strong approximations from $\delta E_x'$ in Eq. (93) in the calculation of

$$\Delta E_x = - \frac{\eta e^2 E_0}{2\epsilon_0 c^2 m} \int_z^\infty J_0(k_\perp \sqrt{r^2 - z^2}) \cos \left[ \omega \left( t - \frac{r}{c} \right) \right] dr, \quad (D1)$$

The primes were forgotten here for the sake of simplicity.

First, let us calculate numerically the quantity $f = J_0(k_\perp \sqrt{r^2 - z^2})$ versus $r$ and let us compare it to its approximate value $g \approx 1 - k_\perp^2 (r^2 - z^2)/4$ in
<figure><img src="image_13.png" /><figcaption>Fig. C1. f versus r (solid line), g versus r (long dashed line).</figcaption></figure>

the case when $\omega = 1.77 \times 10^{14}$ s$^{-1}$ ($\lambda = 10.6$ µm) and $z = 1$ m. We have $k = 5.9 \times 10^5$ m$^{-1}$, assuming $v_0 = 0.01$, we took $k_⊥ = 5.9 \times 10^3$ m$^{-1}$. Fig. C1 shows that g is a pretty good approximation for f.

In Eq. (D1), we integrate from z to $z_f$, where $z_f$ is the value of r for which the Bessel function is zero. The integration from $z_f$ to infinity is ignored as we focus on $\sin \omega(t - \frac{z}{c})$ or $\cos \omega(t - \frac{z}{c})$ terms.

It is assumed that: $J_0(X) \approx 1 - X^{2/4}$ when r is in the range $[z - z_f]$. Thus

$$ \Delta E_x = - \frac{\eta e^2 E_0}{2\epsilon_0 c^2 m} \int_z^{z_f} \left(1 - \frac{k_\perp^2}{4}(r^2 - z^2)\right) \cos \omega \left(t - \frac{r}{c}\right) dr. \qquad (D2) $$

As z and r remain very close to each other in the integration

$$ \begin{aligned} \Delta E_x &= - \frac{\eta e^2 E_0}{2\epsilon_0 c^2 m} \int_z^{z_f} \left(1 - \frac{k_\perp^2}{2}(r-z)\right) \cos \omega \left(t - \frac{r}{c}\right) dr, \\ &= - \frac{\eta e^2 E_0}{2\epsilon_0 c m \omega} \left( \left(1 + \frac{k_\perp^2}{2}z\right) \sin \left[ \omega \left(t - \frac{r}{c}\right) \right]_z^{z_f} - \frac{k_\perp^2}{2} \frac{c}{\omega^2} \left[ \cos \omega \left(t - \frac{r}{c}\right) - r\omega \sin \omega \left(t - \frac{r}{c}\right) \right]_z^{z_f} \right), \\ &= - \frac{\eta e^2 E_0}{2\epsilon_0 m \omega^2} \left[ \sin \omega \left(t - \frac{z}{c}\right) - \sin \omega \left(t - \frac{z_f}{c}\right) \right] + \frac{\eta e^2 E_0}{2\epsilon_0 c m \omega^2} \frac{k_\perp^2}{2} \left[ \cos \omega \left(t - \frac{z_f}{c}\right) - \cos \omega \left(t - \frac{z}{c}\right) \right]. \end{aligned} \qquad (D3) $$

Then, we consider the $\sin \omega(t - \frac{z}{c})$ term

$$ \Delta E_{x \sin \omega} = - \frac{\eta e^2 E_0}{2\epsilon_0 m \omega^2} \sin \omega \left(t - \frac{z}{c}\right). \qquad (D4) $$

In the same way, $\Delta E_y \cos \omega$ was calculated from the second equation giving $\delta E'_y$.

## References

[1] Geints YE, Zemlyanov AA. J Opt Soc Am B 2014;31:788.

[2] S. Tochitsky, E. Welch, M. Polyanskiy, I. Pogorelsky, P. Panagiotopoulos, J. Moloney, and C. Joshi, CLEO: QELS Fundamental Science, CLEO_QELS 2018 (Vol. Part F93-CLEO_QELS 2018).

[3] Wu L, Liu Y, Wu J, Hong S, Wu L, Chen Y, et al. Phys Plasmas 2019;26:092111.

[4] Huang J, Liang Z, Deng F, Yu W, Zhao R, Chen B, et al. J Opt Soc Am A 2015;32:2104.

[5] R.P. Feynman, Lectures on Physics 1, Addison-Wesley, 1977.

[6] Max CE, Arons J. Phys Rev Lett 1974;33:209.

[7] Keldysh LV. Soviet Phys JET 1965;20:1307.

[8] Arissian L, Diels JC. J Lasers Opt Photon 2014;1:107.

[9] Ivanov MY, Spanner M, Smirnova O. J Modern Opt 2005;52:165.

[10] Eckle P, Pfeiffer AN, Cirelli C, Staudte A, Dörner R, Muller HG, et al. Science 2008;322:1525.

[11] Chiao RY, Garmire E, Townes CH. Opt Lett 1964;13:479.

[12] Boyd RW. Nonlinear optics. Boston: Academic; 1992.

[13] Liu W, Petit S, Becker A, Aközbek N, Bowden CM, Chin SL. Opt Commun 2002;202:189.

[14] Jackson JD. Classical electrodynamics. New York: McGraw-Hill; 1975.

[15] Quémada D. Ondes dans les plasmas. Collections Méthodes: Hermann Paris; 1968.

[16] Chen FF. Introduction to Plasma Physics. New York: Plenum; 1974.

[17] Bourdier A, Fortin X. Phys Rev A 1979;20:2154.

[18] Bourdier A, Babonneau D, di Bona G, Fortin X. Phys Rev A 1978;18:1194.

[19] Landau LD, Lifshitz EM. The classical theory of fields. Oxford: Pegamon; 1962.

[20] Rax JM. Physique des Plasmas. Dunod; 2005.

[21] M; Tabor, Chaos and Integrability in Nonlinear Dynamics. John Wiley & Sons, New York, 1989.

[22] A. I. Akhiezer and B.V. Polovin, Zh. Uksp. Teor. Fiz. 30, (1956) 915 (Sov. Phys. JETP 3, (1956)1) 696.

[23] Kaw PK, Dawson J. Phys Fluids 1970;13:472.

[24] Tsintsadze NL, Tskhakaya DD. Sov J Plasma Phys 1975;1:150.

[25] Feynman RP. Lectures on Physics 2. Basic Books; 2013.

[26] Budden KG. The propagation of Radio Waves. Cambridge University Press; 1988.

[27] Ginzburg VL. The propagation of electromagnetic waves in Plasmas. Oxford: Pergamon Press; 1970.