1 Introduction
The transport of high-current relativistic electron beams driven by ultraintense laser interactions with plasmas is relevant to many applications of high energy density physics, particularly in areas of the fast ignition scheme for inertial confinement fusion[Reference Tabak, Hammer, Glinsky, Kruer, Wilks, Woodworth, Campbell, Perry and Mason1], laser-driven ion acceleration[Reference Macchi, Borghesi and Passoni2, Reference Yang, Ma, Shao, Xu, Yu, Gu, Yu, Yin, Tian and Kawata3] and production of ultrashort radiation sources[Reference Jin, Chen, Zhuo, Kon, Nakatsutsumi, Wang, Zhang, Gu, Wu, Zhu, Wang, Yu, Sheng and Kodama4–Reference Zhu, Yu, Sheng, Yin, Turcu and Pukhov7]. A target can be ionized by relativistic electrons both through field ionization and collisional ionization, inducing nonlinear and collective effects that can feed back to the transport of relativistic electrons. It is important to comprehensively investigate the transport process of relativistic electrons in such a target, especially in insulators that are without free electrons initially.
The transport of relativistic electrons can be inhibited by charge separation fields until a cold return current is generated or the electrons from ionization neutralize space charge fields. The transport of relativistic electrons through targets was found to be dependent on target materials, showing spatial disruption in insulators but more uniform propagation in metals[Reference Fuchs, Cowan, Audebert, Ruhl, Gremillet, Kemp, Allen, Blazevic, Gauthier, Geissel, Hegelich, Karsch, Parks, Roth, Sentoku, Stephens and Campbell8–Reference Yang, Xu, Ma, Ge, Zhuo and Shao10]. A collimated ionization channel[Reference Borghesi, Mackinnon, Bell, Malka, Vickers, Willi, Davies, Pukhov and Meyer-ter-Vehn11] or an ionization front behind the collimated jet[Reference Gremillet, Amiranoff, Baton, Gauthier, Koenig, Martinolli, Pisani, Bonnaud, Lebourg, Rousseaux, Toupin, Antonicci, Batani, Bernardinello, Hall, Scott, Norreys, Bandulet and Pépin12] was observed by optical shadowgraph during ultraintense laser interactions with silica targets with propagation velocities (
$0.4c$–
$0.5c$) much slower than the speed of light. One-dimensional (1D) particle-in-cell (PIC) simulations showed a transient longitudinal electric sheath field at the plasma–neutral interface in insulators, which was induced by relativistic electrons penetrating into the neutral matter and was more intense than that in metals, leading to a shorter penetration depth in the former[Reference Kemp, Sentoku, Cowan, Fuchs and Ruhl13]. In addition to the electric sheath field, a so-called ‘fountain’ field was observed by ultrafast interferometry during the propagation of a relativistic electron beam in a solid dielectric[Reference Sarkisov, Ivanov, Leblanc, Sentoku, Yates, Wiewior, Chalyy, Astanovitskiy, Bychenkov, Jobe and Spielman14], which was confirmed by numerical simulations. The authors attributed the fields to relativistic electrons that move faster than the ionization wave and turn back to the target. However, Bai et al. [Reference Bai, Tian, Zhou, Zeng, Sun, Wang and Liu15] showed that such fields can be described by the collective effect of free electrons and ions. Theoretical analysis of high-current relativistic electron beam interactions with the insulator showed that the ionization process can be separated into four regions, i.e., a charge accumulation region; a field ionization region; a current and charge neutralization region; and a collisional ionization region[Reference Debayle and Tikhonchuk16–Reference Tikhonchuk18]. The propagation velocity of the ionization front in insulators depends on the relativistic electron energy and current density. The previous numerical simulations on high-current relativistic electron beam inducing ionization in targets were mainly limited to 1D PIC simulations due to the highly computing resource requirement. However, to comprehensively understand the ionization process in the targets, detailed two-dimensional (2D) PIC simulations are required.
In this paper, ultraintense laser-driven relativistic electron transport in a dielectric silicon (Si) target is studied. An ionization wave with a velocity less than the speed of light is observed by including the field ionization and collisional ionization processes. Widely spread ‘fountain’ electric fields associated with the sheath field occur during the relativistic electron propagation. In addition, two-stream instability (TSI) behind the ionization front arises when the relativistic electron current is sufficiently high. The structure of this paper is as follows. The ionization models (Section 2) and the simulation model (Section 3) are introduced first. The results of relativistic electron transport in an initially unionized Si target are presented in Section 4.1. Section 4.2 shows the influence of laser intensity on the transport of the ionization wave. Finally, in Section 4.3, we present the results of relativistic electron transport in an initially ionized Si target for comparison. Conclusions are provided in Section 5.
2 The ionization models
The 2D3V PIC simulation code EPOCH[Reference Arber, Bennett, Brady, Lawrence-Douglas, Ramsay, Sircombe, Gillies, Evans, Schmitz, Bell and Ridgers19] is employed to model the relativistic electrons propagating in the insulator. Both field ionization induced by the intense laser electric field and self-generated electric field and collisional ionization induced by the relativistic electrons and cold return electrons are included. Following the Keldysh theory, we introduce the Keldysh parameter 
$\unicode[STIX]{x1D6E4}$ to separate the field ionization into multi-photon ionization and tunneling ionization[Reference Keldysh20],
$$\begin{eqnarray}\unicode[STIX]{x1D6E4}=\frac{\unicode[STIX]{x1D714}\sqrt{2m_{e}\unicode[STIX]{x1D716}_{i}}}{eE},\end{eqnarray}$$ where 
$\unicode[STIX]{x1D714}$ is the field frequency, 
$m_{e}$ is the electron mass, 
$\unicode[STIX]{x1D716}_{i}$ is the ionization energy, 
$e$ is the electron charge and 
$E$ is the magnitude of the electric field at the electron. (1) For 
$\unicode[STIX]{x1D6E4}\gg 1$, multi-photon ionization can occur though one photon’s energy is insufficient to cause ionization or excitation. This process is modeled by the semi-empirical WKB approximation[Reference Ammosov, Delone and Krainov21]. (2) When 
$\unicode[STIX]{x1D6E4}\ll 1$, the strong electric field would deform the atomic Coulomb potential and create a finite potential energy barrier that induces tunneling ionization. The ADK ionization rate equation averaged over all possible values of the magnetic quantum equation[Reference Delone and Krainov22] is applied to model the process. The barrier-suppression ionization as a special case of tunneling ionization, in which the electron may escape classically due to the potential energy barrier being lower than the electron binding energy, is also included by a correction to the ADK ionization rate by Posthumus et al. [Reference Posthumus, Thompson, Frasinski and Codling23].
Since we focus on the relativistic electrons propagating through a high-density solid target that is strongly collisional, collisional ionization is also considered. The MBELL equation[Reference Haque, Alfaz Uddin, Basak, Karim, Saha and Malik24] is used as the electron impact ionization cross-section for 
$q\leqslant 36$, while the relativistic modified binary encounter Bethe model[Reference Guerra, Parente, Indelicato and Santos25] is applied for 
$q>36$, where 
$q=A-N_{nl}$; here 
$A$ is the atomic number and 
$N_{nl}$ is the total number of electrons in all suborbitals up to the ionizing orbital. The ionization models have been verified in Ref. [Reference Arber, Bennett, Brady, Lawrence-Douglas, Ramsay, Sircombe, Gillies, Evans, Schmitz, Bell and Ridgers19], which presents a reasonable ionization evolution.
Note that the three-body recombination that is dependent on 
$n_{e}^{2}T_{e}^{-9/2}$ can be neglected for electron temperatures above 10 eV and on the time scale of hundreds of femtoseconds here[Reference Klimo, Tikhonchuk and Debayle26, Reference Hahn27], where 
$n_{e}$ and 
$T_{e}$ are the electron density and temperature, respectively.
3 The simulation model
To investigate the transport of relativistic electrons in insulators, the binary collision model of charged particles proposed by Sentoku et al. [Reference Sentoku and Kemp28] is also considered in addition to the aforementioned ionization models. The target is Si and is initially unionized with a density of 
$2.32~\text{g}/\text{cm}^{3}$, corresponding to an atomic density of 
$44n_{c}$, where 
$n_{c}=1.12\times 10^{21}~\text{cm}^{-3}$ is the critical density for the laser wavelength of 
$\unicode[STIX]{x1D706}_{L}=1.0~\unicode[STIX]{x03BC}\text{m}$. The target has an exponentially increasing density profile 
$n_{a}=0.357n_{c}\exp [(x-10~\unicode[STIX]{x03BC}\text{m})/1~\unicode[STIX]{x03BC}\text{m}]$ in the range 
$x=6{-}10~\unicode[STIX]{x03BC}\text{m}$ and then maintains the maximum atomic density of 
$44n_{c}$ in the range 
$x=10{-}35~\unicode[STIX]{x03BC}\text{m}$. The initial temperature of the target is 10 eV. The simulation box has dimensions 
$40~\unicode[STIX]{x03BC}\text{m}\times 50~\unicode[STIX]{x03BC}\text{m}$ with 
$2400\times 3000$ cells. Each cell contains 10 macro-particles for Si atoms in the simulations. It is worth mentioning that such a simulation is computationally expensive since a large number of electrons and ions will be ionized at later times, which can cause severe load imbalance. A p-polarized laser pulse with a peak intensity of 
$I_{L}=10^{20}~\text{W}/\text{cm}^{2}$ is incident normally from the left boundary, corresponding to 
$a_{L}=eE_{L}/m_{e}c\unicode[STIX]{x1D714}_{L}=8.54$ for the dimensionless maximum amplitude of the laser electric field, where 
$c$ is the speed of light and 
$E_{L}$ and 
$\unicode[STIX]{x1D714}_{L}$ are the electric field and frequency of the laser pulse, respectively. Both the spatial and temporal profiles of the laser are Gaussian, with a spot radius of 
$4~\unicode[STIX]{x03BC}\text{m}$ and a duration of 80 fs (full-width at half-maximum). In order to suppress numerical heating, a fifth-order interpolation scheme is employed to evaluate the currents and a fourth-order finite difference scheme is applied for solving Maxwell’s equations. Open boundary conditions are used for both the fields and particles.
Distributions of average ionization degree (
$\bar{Z}$) (a) and 
$\log _{10}$ of electron density (
$n_{e}$) (b) at 
$t=150~\text{fs}$. Profiles of 
$\bar{Z}$ and 
$n_{e}$ (c) and density of 
$\text{Si}^{4+}$ and 
$\text{Si}^{12+}$ (d) along the laser propagation axis at 
$t=100~\text{fs}$ and 150 fs, which are averaged over one wavelength around 
$y=0$. Both the electron and ion densities are in units of 
$n_{c}$ here and in other figures.
Distributions of the quasi-static magnetic field (
$B_{z}$) [(a) and (b)], the longitudinal electrostatic field (
$E_{x}$) [(c) and (d)] and the transverse electrostatic field (
$E_{y}$) [(e) and (f)] at 
$t=100~\text{fs}$ [(a), (c) and (e)] and 150 fs [(b), (d), and (f)]. The fields are averaged over two laser cycles and the fields in front of the solid target (
$z<10~\unicode[STIX]{x03BC}\text{m}$) are not shown for clarity. The magnetic field and electric field are in units of tesla and 
$\text{V}/\text{m}$, respectively.
4 Results and discussions
4.1 Relativistic electron transport in a solid Si target
Figure 1 shows distributions of the average ionization degree (
$\bar{Z}$) and the electron density (
$n_{e}$) of the target as the intense laser pulse irradiates on the Si target. The field ionization is the dominant ionization scheme in the preplasma[Reference Lang, Yang, Xu, Jin and Zhuo29], and the maximum ionization degree in the preplasma is 
${\sim}12$, which is consistent with the field ionization by the laser electric field. The final ionization of an atom by a given electric field can be estimated from[Reference Delone and Krainov22]
$$\begin{eqnarray}E_{L}=\frac{E_{h}}{4Z}\left(\frac{\unicode[STIX]{x1D716}_{i}}{2\unicode[STIX]{x1D716}_{h}}\right)^{2},\end{eqnarray}$$ where 
$Z$ is the ionization degree, 
$E_{h}=5.14\times 10^{11}~\text{V}/\text{m}$ is the atomic electric field, 
$\unicode[STIX]{x1D716}_{h}=13.6~\text{eV}$ is the ionization energy of hydrogen and 
$\unicode[STIX]{x1D716}_{i}$ is the ionization energy of the 
$i$th electron for a given atom. The electric field corresponding to the peak laser intensity here is 
$2.74\times 10^{13}~\text{V}/\text{m}$. Thus it can ionize the target easily to 
$\text{Si}^{12+}$. The collisional ionization occurs after the relativistic electrons are generated and propagate through the target. The target is ionized rapidly, and the Si atoms are ionized completely to 
$\text{Si}^{14+}$ in the laser–target interaction region, as shown in Figure 1(a). Filamented ionization distribution appears deep in the target (
$x>25~\unicode[STIX]{x03BC}\text{m}$) due to filamented propagation of the relativistic electrons (Figure 1(b)). For clarity, the profiles of 
$\bar{Z}$ and 
$n_{e}$ at two different times are presented in Figure 1(c). One can see two steep steps in the average ionization profiles, corresponding to 
$\text{Si}^{4+}$ and 
$\text{Si}^{12+}$ due to the sharp increase of the ionization energy for the states [from 45.14 eV (
$\text{Si}^{4+}$) to 166.77 eV (
$\text{Si}^{5+}$) and 523.42 eV (
$\text{Si}^{12+}$) to 2437.66 eV (
$\text{Si}^{13+}$)]. The electrons associated with the target ionization also have step distributions in the density profile. The electrons in front of the ionization front are very few (
${\sim}10^{-3}n_{c}$), but they are the source of the ‘fountain’ field (see the following paragraph). In Figure 1(d), the distribution of 
$\text{Si}^{4+}$, with a narrower width, is always ahead of that of 
$\text{Si}^{12+}$, and both widths increase with the penetration depth of the relativistic electrons. At the head of the relativistic electron beam, the field ionization is the main ionization mechanism due to the high energy and lower density of relativistic electrons, while the collisional ionization takes place after the return currents are turned back by the sheath field (see the discussion below) and the resistive electric field (
$\mathbf{E}=\unicode[STIX]{x1D702}\mathbf{j}_{c}$), which contributes most of the ionized electrons in the target[Reference Klimo, Tikhonchuk and Debayle26]. Here 
$\unicode[STIX]{x1D702}$ is the resistivity of the target and 
$\mathbf{j}_{c}$ is the return current density.
The self-generated magnetic field, longitudinal electrostatic field and transverse electrostatic field are shown in Figure 2. A significant resistive magnetic field (∼25 T) is observed at the periphery of the relativistic electron beam, generated by[Reference Yang, Borghesi and Robinson30]
$\unicode[STIX]{x2202}\mathbf{B}/\unicode[STIX]{x2202}t=\unicode[STIX]{x1D6FB}\times (\unicode[STIX]{x1D702}\mathbf{j}_{h})=\unicode[STIX]{x1D6FB}\unicode[STIX]{x1D702}\times \mathbf{j}_{h}+\unicode[STIX]{x1D702}\unicode[STIX]{x1D6FB}\times \mathbf{j}_{h}$, where 
$\mathbf{j}_{h}$ is the relativistic electron current density. Such a field has a pinch effect on the relativistic electron beam and can collimate the relativistic electrons, provided it has the appropriate amplitude and structure[Reference Xu, Yang, Sheng, McKenna, Ma, Zhuo, Yin, Ren and Zhang31]. Magnetic filaments are observed behind the ionization front due to the filamented current propagation of the relativistic electrons. The ionization front here is defined as the position where the electric field is equal to the threshold electric field of field ionization (
$E_{th}$) given by Equation (2), with 
$\unicode[STIX]{x1D716}_{i}$ being the ionization energy of the outermost electron. The ionization energy of the outermost electron for Si is 8.15168 eV; thus 
$E_{th}=1.154\times 10^{10}~\text{V}/\text{m}$.
Initially, the laser pulse can penetrate the target before the ionization takes place, as shown in the transverse electric field distributions averaged over two laser cycles (Figure 2(e)). After the target has been ionized, it can only penetrate into the region of relativistic critical density (
$a_{L}n_{c}$). A moving longitudinal electric field is observed in the target, with an amplitude approximately that of the threshold electric field of the ionization. (A weaker and widely spread field is also observed in front of the moving field.) The intense field is generated by charge separation in the target, similar to the sheath field near the interface of the target and vacuum in laser-driven ion acceleration schemes[Reference Yang, Ma, Shao, Xu, Yu, Gu, Yu, Yin, Tian and Kawata3]. The field will slow down the relativistic electrons and pull back the cold electrons to neutralize the charge separation, and collisional ionization is induced that contributes to most of the free electrons in the target. The sheath field and the ‘fountain’ field are clearly present in the transverse electric field distributions (Figures 2(e) and 2(f)); both are very weak near the laser propagation axis (
$y=0$) due to symmetry. Comparing these figures to Figures 1(a) and 1(b), we conclude that the fountain field is associated with the very rarefied relativistic electrons (of density 
$10^{-3}n_{c}$) that propagate faster than the ionization wave. The distribution of the electron longitudinal momentum (Figure 3(b)) shows that the electrons ahead of the ionization wave only have positive momentum, indicating that they are forward propagating in the target. This is consistent with the explanation that the fountain field is produced by the collective effect of the free electrons and ions[Reference Bai, Tian, Zhou, Zeng, Sun, Wang and Liu15]. The fountain field is purely radial and is too weak to induce breakdown of the target. The width of the sheath field remains constant during its propagation in the target, while the fountain field spreads with time due to the spreading of the relativistic electrons ahead of the ionization front.
The energy spectrum of the electrons (a) and the distribution of the electron longitudinal momentum (
$P_{x}$) along the laser axis at 
$t=150~\text{fs}$ (b).
The relativistic electrons are mainly accelerated by the 
$\mathbf{J}\times \mathbf{B}$ acceleration[Reference Wilks, Kruer, Tabak and Langdon32] during the laser and target interaction, as can be seen from the distribution of the electron longitudinal momentum (
$p_{x}$) along the laser axis (Figure 3(b)), which shows a 
$2\unicode[STIX]{x1D714}_{L}$ modulation. The maximum 
$p_{x}$ is 
${\sim}40m_{e}c$, which is close to the longitudinal momentum of a relativistic electron in a plane electromagnetic wave in vacuum (
$p_{x}=\frac{1}{2}a_{L}^{2}m_{e}c$), while the maximum transverse momentum (
$p_{y}$) reaches 
$24m_{e}c$, which is much higher than that in vacuum (
$p_{y}=a_{L}m_{e}c$)[Reference Yu, Bychenkov, Sentoku, Yu, Sheng and Mima33]. It also shows that in addition to the large population of forward propagating relativistic electrons, many electrons are accelerated backward by the reflected laser near the front target surface. In addition, the cold electrons are turned back by the sheath field and resistive electric field as the relativistic electrons propagate forward in the target (as indicated by the negative longitudinal momentum). The relativistic electron current is nearly neutralized by the cold return current, as shown in Figure 4, which presents the current densities for the relativistic electrons and cold electrons. This leads to the continuous forward propagation of the relativistic electrons in the target. Figure 3(a) shows the energy spectrum of the electrons, which has a three-temperature distribution. The temperature for the highest energy electrons (kinetic energy greater than ∼22 MeV) is not given in the figure since their population is too small to be reliably fitted with a Maxwellian distribution function. The temperature for the main relativistic electron population is approximately 
$T_{h}=2.85~\text{MeV}$, which is lower than the ponderomotive scaling for an ultraintense p-polarized laser pulse interacting with plasmas[Reference Wilks, Kruer, Tabak and Langdon32], i.e., 
$T_{h}=m_{e}c^{2}(\unicode[STIX]{x1D6FE}_{h0}-1)$, where 
$\unicode[STIX]{x1D6FE}_{h0}=\sqrt{1+a_{L}^{2}}$ is the Lorentz factor of the relativistic electrons. This can be attributed to the relatively short preplasma employed in the simulation, which is pushed into the solid target by the intense radiation pressure of the laser pulse at later times. Since the laser pulse can only penetrate to skin depth after the target has been ionized, the electrons cannot be accelerated by the full ponderomotive potential.
Current density distributions for the relativistic electrons (kinetic energy 
$K_{e}\geqslant 50~\text{keV}$) (a) and cold electrons (
$K_{e}\leqslant 10~\text{keV}$) (b) at 
$t=150~\text{fs}$. The current is in units of 
$en_{c}c$.
We now calculate the velocity of the ionization front using a reference frame co-moving with the ionization front. Assume that the velocity of the ionization front in the laboratory frame is 
$v_{f}$ and the corresponding Lorentz factor is 
$\unicode[STIX]{x1D6FE}_{f}$. Then the temperature of the relativistic electrons in the reference of the ionization front can be obtained by Lorentz transformation,
$$\begin{eqnarray}T_{h}^{\prime }=m_{e}c^{2}\left[\unicode[STIX]{x1D6FE}_{f}\left(\unicode[STIX]{x1D6FE}_{h0}-\unicode[STIX]{x1D6FD}_{f}\sqrt{\unicode[STIX]{x1D6FE}_{h0}^{2}-1}\right)-1\right],\end{eqnarray}$$ where 
$\unicode[STIX]{x1D6FD}_{f}=v_{f}/c$. The corresponding relativistic electron density is 
$n_{h}^{\prime }=n_{h}/\unicode[STIX]{x1D6FE}_{f}$, where 
$n_{h}$ is the density in the lab frame.
Fast Fourier transform of 
$E_{x}$ (a) and 
$E_{y}$ (b) behind the ionization front (i.e., 
$23~\unicode[STIX]{x03BC}\text{m}<x<26~\unicode[STIX]{x03BC}\text{m},-8~\unicode[STIX]{x03BC}\text{m}<y<8~\unicode[STIX]{x03BC}\text{m}$) at 
$t=150~\text{fs}$. 
$k_{0}$ is the wave number of the laser pulse.
The sheath field in the ionization front reference is estimated by 
$E_{sh}\approx (8\unicode[STIX]{x1D70B}n_{h}^{\prime }T_{h}^{\prime })^{1/2}$[Reference Bychenkov, Novikov, Batani, Tikhonchuk and Bochkarev34]. One can simply consider that the sheath field is equal to the threshold of the ionization field in order to obtain the velocity of the ionization front, i.e., 
$E_{sh}=E_{th}$[Reference Sarkisov, Ivanov, Leblanc, Sentoku, Yates, Wiewior, Chalyy, Astanovitskiy, Bychenkov, Jobe and Spielman14]; thus,
$$\begin{eqnarray}\frac{E_{th}^{2}}{8\unicode[STIX]{x1D70B}m_{e}c^{2}n_{h}}=\unicode[STIX]{x1D6FE}_{h0}-\unicode[STIX]{x1D6FD}_{f}\sqrt{\unicode[STIX]{x1D6FE}_{h0}^{2}-1}-\sqrt{1-\unicode[STIX]{x1D6FD}_{f}^{2}}.\end{eqnarray}$$ Given the density and temperature of the relativistic electrons, 
$\unicode[STIX]{x1D6FD}_{f}$ can be explicitly solved from Equation (4). Note that the velocity of the ionization front given by Equation (4) is very close to that given in Ref. [Reference Debayle and Tikhonchuk16], in which the current density and electron energy conservation, and the Poisson equation were applied together to derive the ionization front velocity. Since the threshold of the ionization electric field is typically low, 
$\unicode[STIX]{x1D6FD}_{f}$ depends weakly on 
$E_{th}$ and 
$n_{h}$, and it is mainly dependent on the relativistic electron energy 
$\unicode[STIX]{x1D6FE}_{h0}$. Also the ionization wave velocity is always smaller than the relativistic electron velocity. For the laser intensity employed here, Equation (4) gives 
$v_{f}\sim c$, which is greater than that obtained in the simulation (
${\sim}0.81c$). This can be attributed to two factors: (1) the average energy of the relativistic electrons given by the ponderomotive scaling is higher than that in our simulation; (2) the inhibition of the relativistic electrons by the resistive electric field and other linear or nonlinear instabilities, which could dissipate the electron energy, is not included in the theoretical analysis.
Both magnetic and electric fields in the target show instabilities in Figure 2. Here we mainly focus on the small-scale instabilities for the electric fields. To understand the mechanism of the instabilities, we plot in Figure 5 the fast Fourier transform of 
$E_{x}$ and 
$E_{y}$ behind the ionization front. Both spectra show a peak near 
$|k|=14k_{0}$, i.e., 
$|k_{x}|\approx 14k_{0}$ for the 
$|E_{xk}|$ distribution and 
$|k_{y}|\approx 14k_{0}$ for the 
$|E_{yk}|$ distribution, where 
$k_{0}$ is the wave number of the laser pulse. Such electric field instability corresponds to a wavelength of 71 nm, and this short scalelength instability can be explained by TSI. Combining the mass continuity equation for both species, momentum equation and Poisson equation, the dispersion equation of TSI is given as[Reference Klimo, Tikhonchuk and Debayle26]
$$\begin{eqnarray}1-\frac{\unicode[STIX]{x1D714}_{ph}^{2}}{\unicode[STIX]{x1D6FE}_{h0}^{3}(\unicode[STIX]{x1D714}-kv_{h})^{2}}-\frac{\unicode[STIX]{x1D714}_{pe}^{2}}{\unicode[STIX]{x1D714}(\unicode[STIX]{x1D714}+i\unicode[STIX]{x1D708}_{e})^{2}}=0,\end{eqnarray}$$ where 
$\unicode[STIX]{x1D714}_{ph,pe}=(4\unicode[STIX]{x1D70B}n_{h,e}e^{2}/m_{e})^{1/2}$ is the plasma frequency for the relativistic electrons and return electrons, 
$\unicode[STIX]{x1D708}_{e}=(8\unicode[STIX]{x1D70B}/3\sqrt{3})(n_{e}e^{4}\ln \unicode[STIX]{x1D6EC}/m_{e}^{1/2}T_{h}^{3/2})$ is the collision frequency of the relativistic electrons and return electrons and 
$\ln \unicode[STIX]{x1D6EC}$ is the Coulomb logarithm. For the relativistic electron density 
${\sim}0.4n_{c}$ and the return electron density 
${\sim}190n_{c}$ near 
$x=25~\unicode[STIX]{x03BC}\text{m}$, the collision frequency is 
$\unicode[STIX]{x1D708}_{e}\sim 2.45\times 10^{9}/\text{s}$, which is negligible compared to the growth rate of TSI discussed below. The maximum growth rate of TSI for collisionless electrons is given by[Reference Bret35]
$\unicode[STIX]{x1D6E4}=(\sqrt{3}/2^{4/3})(\unicode[STIX]{x1D6FC}^{1/3}/\unicode[STIX]{x1D6FE}_{h0})\unicode[STIX]{x1D714}_{pe}$, where 
$\unicode[STIX]{x1D6FC}=n_{h}/n_{e}$, leading to 
$\unicode[STIX]{x1D6E4}\sim 2.67\times 10^{14}/\text{s}$. It corresponds to the electron plasma wave with phase velocity equal to the velocity of the relativistic electrons 
$\unicode[STIX]{x1D714}_{pe}/k_{res}\approx v_{h}$, that is, the resonance mode near the return electron plasma frequency. The oscillating frequency of the electric field calculated from simulation is 
${\sim}14\unicode[STIX]{x1D714}_{0}$, close to the electron plasma frequency. Since the growth rate increases as 
$n_{h}^{1/3}n_{e}^{1/6}$, and recalling that 
$n_{e}\propto n_{h}$, such instability disappears for the case with laser intensity less than 
$5\times 10^{19}~\text{W}/\text{cm}^{2}$, which is consistent with the results of Klimo et al. [Reference Klimo, Tikhonchuk and Debayle26].
4.2 Influence of laser intensity on the transport of ionization waves
To study the influence of laser intensity on the transport of the ionization wave, we vary the laser intensity while keeping other parameters the same as those in the above case. The profiles of average ionization degree and electron density along the laser propagation axis for different laser intensities are shown in Figure 6. The ionization depth in the target decreases as the laser intensity decreases, due to the decrease of both the relativistic electron currents and kinetic energy for the lower laser intensities. The full ionization depth in the target is 
$x\sim 16.4~\unicode[STIX]{x03BC}\text{m}$ for a laser intensity of 
$10^{20}~\text{W}/\text{cm}^{2}$ at 
$t=150~\text{fs}$, while only the target near the solid interface (
$x\sim 10~\unicode[STIX]{x03BC}\text{m}$) is ionized completely for a laser intensity of 
$5\times 10^{18}~\text{W}/\text{cm}^{2}$. The steep steps in the average ionization profiles for 
$\text{Si}^{4+}$ and 
$\text{Si}^{12+}$ are present at all intensities. Since the ionization decreases for lower laser intensities, the electron density is decreased. The penetration depth for an electron density of 
$44n_{c}$ (corresponding to the ionization of 
$\text{Si}^{+}$) is 
$x\sim 26.9~\unicode[STIX]{x03BC}\text{m}$ for a laser intensity of 
$10^{20}~\text{W}/\text{cm}^{2}$ at 
$t=150~\text{fs}$ and 
$x\sim 17.6~\unicode[STIX]{x03BC}\text{m}$ for a laser intensity of 
$5\times 10^{18}~\text{W}/\text{cm}^{2}$.
Profiles of the average ionization degree 
$\bar{Z}$ (a) and electron density 
$n_{e}$ (b) along the laser propagation axis (averaged over one wavelength near 
$y=0$) at 
$t=150~\text{fs}$ for laser intensities of 
$5\times 10^{18}~\text{W}/\text{cm}^{2}$, 
$1\times 10^{19}~\text{W}/\text{cm}^{2}$, 
$5\times 10^{19}~\text{W}/~\text{cm}^{2}$ and 
$1\times 10^{20}~\text{W}/\text{cm}^{2}$.
The dependence of the ionization front velocity on laser intensity is presented in Figure 7, which also includes the results from the theoretical analysis and two experiments. The ionization front velocities are close to the speed of light when the laser intensity is greater than 
$10^{19}~\text{W}/\text{cm}^{2}$ from the theory, which are much higher than those obtained in the 2D PIC simulations. One of the reasons for this inconsistency could be the lower relativistic electron energy obtained in the numerical simulation compared to that given by the ponderomotive scaling (e.g., 2.85 MeV versus 3.88 MeV for a laser intensity of 
$10^{20}~\text{W}/\text{cm}^{2}$). However, we find that this is not the main reason. The collective and nonlinear effects should be considered appropriately to describe the transport of relativistic electrons in the insulator. Our results are higher than those from the experiments. The discrepancies could be due to the following factors. (1) A glass target, with oxygen (O), was used instead of pure silicon in Ref. [Reference Sarkisov, Ivanov, Leblanc, Sentoku, Yates, Wiewior, Chalyy, Astanovitskiy, Bychenkov, Jobe and Spielman14]. The ionization energy of the outermost electron of O is 13.618 eV, corresponding to 
$E_{th}=3.22\times 10^{10}~\text{V}/\text{m}$, which is three times that of Si. This indicates that, in order to reach the threshold, more relativistic electrons are required [since 
$E_{sh}\propto (n_{h}^{\prime }T_{h}^{\prime })^{1/2}$], leading to a slower velocity of the ionization wave in the experiment compared to that obtained in the pure Si target here. (2) For the experiment by Gremillet[Reference Gremillet, Amiranoff, Baton, Gauthier, Koenig, Martinolli, Pisani, Bonnaud, Lebourg, Rousseaux, Toupin, Antonicci, Batani, Bernardinello, Hall, Scott, Norreys, Bandulet and Pépin12], a silica target coated by an Al foil was used. An intense magnetic field could be generated due to the resistivity gradient near the interface of Al and the silica target. This can inhibit the relativistic electron penetration into the silica target, as shown in our previous work[Reference Yang, Xu, Zhuo, Ma, Shao, Yin and Borghesi36, Reference Yang, Zhuo, Ma, Xu, Yu, Zou, Ge, Xu, Zhu, Shao and Borghesi37]. (3) The velocity of the ionization wave decreases as it propagates in the target. We have applied a relatively short duration of diagnosis in the simulations, which also could induce the discrepancy between our simulations and the experiments. This indicates that 2D simulations are better than the 1D theoretical analysis to describe accurately the transport of a relativistic electron beam in an insulator since the 1D theoretical results might overestimate the ionization front velocity significantly.
The velocity of ionization wave as a function of laser intensity in the dielectric target, where 1D theory denotes the results from Equation (4); UNR and LULI denote the experimental results from Refs. [Reference Gremillet, Amiranoff, Baton, Gauthier, Koenig, Martinolli, Pisani, Bonnaud, Lebourg, Rousseaux, Toupin, Antonicci, Batani, Bernardinello, Hall, Scott, Norreys, Bandulet and Pépin12, Reference Sarkisov, Ivanov, Leblanc, Sentoku, Yates, Wiewior, Chalyy, Astanovitskiy, Bychenkov, Jobe and Spielman14], respectively.
4.3 Relativistic electron transport in a plasma target without ionization
In order to make a comparison with the transport of relativistic electrons in the dielectric target, relativistic electron transport in a preionized Si target is also investigated. The target is composed of neutralized electrons and 
$\text{Si}^{3+}$ initially. Except for the exclusion of the ionization process, particle collisions and other parameters are the same as those in the above cases. As expected, the sheath field and the fountain field are not observed in the plasma target (not shown for brevity) because free electrons are distributed throughout the target. Figure 8 shows the electron energy spectrum and the electron longitudinal momentum distribution along the laser axis. One can see that the electron spectrum has a similar profile to that of the ionization case, and the temperature of the relativistic electrons is nearly identical. The only difference is that the number of relativistic electrons generated in the ionization case is slightly lower than that of the case without ionization, which can be attributed to the energy loss in the ionization process. The phase distribution of 
$p_{x}$ in Figure 8(b) shows that the relativistic electrons propagate much faster and penetrate deeper than that in the ionization case because of the absence of the sheath field. However, the return currents pulled back by the charge separation field and the resistive electric field are still present. It is important to ensure the high-current relativistic electron beam forward propagation though the magnitude of 
$p_{x}$ and phase space density are much smaller than those in the ionization case. This suggests that the ionization process should be considered in PIC simulations to accurately describe the relativistic electron generation and propagation in dielectric targets.
The energy spectrum of electrons (a) and distribution of the electron longitudinal momentum (
$P_{x}$) along the laser axis at 
$t=150~\text{fs}$ (b) for the case with a 
$\text{Si}^{3+}$ target and without the ionization process. The electron energy spectrum for the ionization case is also presented for comparison.
5 Conclusion
In summary, ultraintense laser-driven relativistic electron transport in a dielectric Si target is studied using PIC simulations including the field and collisional ionization processes. An ionization wave propagating in the target is observed. In addition to the intense sheath fields (close to the threshold electric field of field ionization), a widely spread ‘fountain’ electric field occurs ahead due to the collective effect of free electrons and ions. The velocity of the ionization wave increases with laser intensity but is much less than the speed of light and also that from the 1D theoretical analysis, indicating that 2D3V numerical simulations are better to describe the relativistic electron transport in dielectric targets. The TSI behind the ionization front arises with increase in laser intensity. The results are useful for applications related to laser-driven relativistic electron transport in dielectric targets such as fast ignition and laser-driven ion acceleration.
Acknowledgements
This work was supported by the National Natural Science Foundation of China (Nos. 11775305, 11675264 and 11705282), Science Challenge Project (No. TZ2018001) and Open Fund of the State Key Laboratory of High Field Laser Physics (SIOM). X. H. Yang also acknowledges the support from the China Scholarship Council.
Open access
