High order harmonics from relativistic electron spikes

High-order harmonic generation is one of the most fundamental processes of nonlinear optics. The coherency of high-frequency radiation associated with high-order harmonics is necessary for the tightest focusability and generation of ultrashort pulses required for applications. Numerous nonlinearities produce high-order harmonics in various laser–matter interaction regimes. At irradiances of the order of 10¹⁸ W cm⁻² and higher, corresponding to relativistic plasma dynamics [1, 2], several regimes of high order harmonics generation are known, including laser–solid target interaction [3–8], harmonics generation in a plasma channel [9], and the second harmonic from electro-optic shocks [10, 11]. Nonlinear Thomson scattering [12] by plasma electrons [13–17] and electron beams [18, 19] also produces harmonic spectra. At typical irradiances from ∼10¹³ to 10¹⁵ W cm⁻² harmonics are generated in atomic [20] or molecular [21] gases or in cold ablation plumes [22] via the nonlinearity of ionizing matter [23, 24]; a review of this process and related issues can be found in [25]. Recently [26] we discovered a new high-order harmonic generation regime in relativistic laser plasma, explained by a new mechanism based on catastrophe theory [27, 28].

The two most important parameters determining the relativistic laser–plasma interaction are the dimensionless laser pulse amplitude a₀ = eE₀/m_ecω₀ and ratio n_e/n_cr of the electron density n_e to the critical density n_cr = m_eω₀²/4πe² = π/r_eλ₀² ≈ 1.1 × 10²¹ cm⁻³ (μm/λ₀)². Here e and m_e are the electron charge and mass, c is the speed of light in vacuum, r_e = e²/m_ec² ≈ 2.8 × 10⁻¹³ cm is the classical electron radius, and λ₀, ω₀, and E₀ are the laser wavelength, angular frequency, and peak electric field, respectively. The laser irradiance in terms of the dimensionless amplitude is I_0L ≈ a₀² × 1.37 × 10¹⁸ W cm⁻² × (μm/λ₀)² for linear and I_0C = 2I_0L for circular polarizations.

In this paper we consider a relativistically strong laser pulse, a₀ ⩾ 1, and underdense plasma, n_e/n_cr ∼ 0.01 (see the details in the setup description below). Under these conditions plasma is essentially collisionless [29]. Propagating through the plasma, the laser pulse undergoes relativistic self-focusing [30–34] when its power P₀ exceeds the critical power [32, 34] P_SF:

$$\begin{eqnarray}&&{{P}_{0}}\gt {{P}_{{\rm SF}}}={{P}_{{\rm c}}}\frac{{{n}_{{\rm cr}}}}{{{n}_{{\rm e}}}}\approx 2{\rm TW}\left( \frac{{\rm 1}{{{\rm 0}}^{{\rm 19}}}{\rm c}{{{\rm m}}^{-3}}}{{{n}_{{\rm e}}}} \right){{\left( \frac{{\rm 1}\mu {\rm m}}{{{\lambda }_{0}}} \right)}^{2}},{{P}_{{\rm c}}}=2{{m}_{{\rm e}}}^{2}{{c}^{5}}/{{e}^{2}}\approx 17{\rm GW}.\end{eqnarray} \tag{ 1 }$$

Under our conditions, the typical values are n_e ∼ 10¹⁹ cm⁻³, and λ ∼ 1 μm, thus giving P_SF ∼ 2 TW. For mid-infrared lasers and high-density gas jets (n_e ∼ 10²⁰ cm⁻³) the critical power for self-focusing can be smaller than ∼ 100 GW. Due to self-focusing, the laser pulse acquires a much higher amplitude, which can be estimated using the stationary focusing condition [35]

$$\begin{eqnarray}&&{{a}_{{\rm SF}}}={{\left( 8\pi \frac{P}{{{P}_{{\rm c}}}}\frac{{{n}_{{\rm e}}}}{{{n}_{{\rm cr}}}} \right)}^{1/3}},\end{eqnarray} \tag{ 2 }$$

while the focal spot diameter is

$$\begin{eqnarray}&&{{d}_{{\rm SF}}}=2a_{{\rm SF}}^{1/2}\frac{c}{{{\omega }_{{\rm pe}}}}=\frac{\lambda }{\pi }{{\left( {{a}_{{\rm SF}}}\frac{{{n}_{{\rm cr}}}}{{{n}_{{\rm e}}}} \right)}^{1/2}},\end{eqnarray} \tag{ 3 }$$

where ω_pe = (4π n_e e²/m_e)^1/2 is the Langmuir frequency.

In the course of the relativistic laser–plasma interaction, nonlinear plasma waves are generated. The bow wave [36] is generated by the electrons expelled in the transverse direction. The bow wave formation condition, d_SF < 4a_SF^1/2 c/ω_pe, is satisfied during the stationary self-focusing, equation (3). The wake wave [37], emerging from the longitudinal electron oscillations, produces a moving potential, which is the basis of the laser wake-field accelerator (LWFA) concept [38, 39]. The electrons, spontaneously injected into and accelerated by the wake wave and oscillating in the plasma channel, can emit betatron x-ray radiation [40–42]. Nearly-breaking nonlinear wake waves act as partially reflective relativistic flying mirrors [43–47] for counter-propagating electromagnetic radiation. The flying mirror upshifts the frequency and shortens the duration of a reflected pulse due to the double Doppler effect, and also focuses the reflected pulse due to the close-to-parabolic mirror surface [48, 49].

In this paper we describe new features of the high-order harmonic generation regime in underdense relativistic laser plasma, first described in [26], including data obtained with a two-channel spectrograph, harmonic generation by circularly polarized laser pulses, spatial coherence and angular distribution. We compare the properties of the harmonics in the new generation regime with the previously known scenarios and show that several advantages make the new regime a possible candidate for a next-generation compact coherent x-ray source.

The paper is organized as follows. In section 2 we describe the experimental setup, followed by the experimental results in section 3. Section 4 contains a comparison of the observed harmonics with the previously suggested mechanisms, and a discussion about possible harmonics origin. Section 5 contains simulation results and their analysis, followed by a description of the new harmonics generation mechanism in section 6, a catastrophe theory based model in section 7, and simplified scaling in section 8. Section 9 contains our conclusions.

We have performed experiments with the J-KAREN laser in KPSI, JAEA, Japan [50] and the Astra Gemini laser in the CLF, RAL, UK [51]. The laser parameters used in these experiments are listed in table 1. The J-KAREN laser has the pulse energy E_L = 0.4 J, peak power P₀ ≈ 9 TW, and central wavelength λ₀ = 820 nm. The pulse duration measured with our transient grating frequency resolved optical gating (TG FROG) [52] is typically τ ≈ 27 fs (full width at half maximum, FWHM), with respect to power. The laser pulses are focused by an f/9 off-axis parabolic mirror down to the spot size of ∼10–12 μm FWHM and ∼20–30 μm at 1/e² intensity level measured by imaging the spot with magnification onto a charge-coupled device (CCD) using a laser pulse sufficiently attenuated with wedges. The derived peak irradiance in vacuum is I₀ ≈ 4 × 10¹⁸ W cm⁻². For the Astra Gemini laser, the pulse energy E_L is up to 10 J, P₀ ranges from 60–170 TW, central wavelength λ₀ = 804 nm, and typical duration τ ≈ 54 fs (FWHM) measured with spectral phase interferometry for direct electric field reconstruction (SPIDER). With an f/20 off-axis parabolic mirror, the typical focal spot size is 25 μm × 40 μm (FWHM) and ∼40 μm × 80 μm at 1/e²; these spot sizes are estimated from far field images of leakage through dielectric mirrors during the full-power shots. The derived peak irradiance in vacuum ranges from 2 × 10¹⁸–5 × 10¹⁸ W cm⁻².

The experimental setups are schematically shown in figure 1. The J-KAREN laser pulses have been focused into a supersonic high-purity helium gas jet streamed from a conical nozzle with the orifice diameter of 1 mm and calculated Mach number of M = 3.3 at a distance of 1 mm from the nozzle orifice. In the experiment with the higher-energy Astra Gemini laser, a conical nozzle with a smaller diameter 0.5 mm orifice with M = 2.0 has been used to decrease the energy of the accelerated electrons (by decreasing the acceleration length) and thus reduce the noise generated by hard x-ray photons. We noticed nozzle degradation after shots performed at the 1 mm distance between the focus and the nozzle orifice and, therefore, increased this distance to 1.5 mm. The gas density distribution has been characterized with interferometry using Ar gas, figure 2 (inset). The Ar gas clusterization does not significantly affect the density measurement because for Hagena's parameter [53] values of Γ* < 7 × 10⁴, the dryness is nearly at unity with an accuracy of a few percent [54] and, thus, the Ar density is close to the He density. We have not observed shocks at the gas jet–ambient gas interface. The electron density of the helium plasma is calculated to be twice the atomic density in the gas jet, because the laser irradiance significantly exceeds the threshold of the barrier suppression double ionization of helium [55]. Typical distributions of plasma density are shown in figure 2. The bright harmonic spectra are observed when the maximum electron density ranges from 1.7 × 10¹⁹–7 × 10¹⁹ cm⁻³, which corresponds to 0.01 ⪅ n_e/n_cr ⪅ 0.04 for λ₀ = 820 nm. Weaker spectra are detectable for densities down to 1.5 × 10¹⁹ cm⁻³ and up to ∼10²⁰ cm⁻³. The density gradient varies by a factor of five, from 2 × 10²⁰–11 × 10²⁰ cm⁻⁴. The parameters of the gas jet targets are listed in table 2.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

In our experiments the peak power in vacuum for both lasers significantly exceeds the critical power of the relativistic self-focusing, equation (1). Thus, during propagation in plasma the laser pulses self-focus, acquiring much higher irradiances than in vacuum. We note that the density gradient influencing the laser–target interaction [6, 8, 56] and the focusing optics f-number can be used for the optimisation of channel formation, laser pulse propagation and self-focusing [57–59]. As is typical for this parameter regime, we observe electron acceleration up to a few hundred MeV. To reduce the noise caused by bremsstrahlung x-rays, these electrons are deflected from the measurement direction by permanent magnets, and directed to the phosphor screens of the electron spectrometers.

The extreme ultraviolet (XUV) and soft x-ray harmonics have been observed in the forward direction with grazing-incidence flat-field spectrographs comprising a gold-coated collecting mirror, spherical varied-linespace grating, and back-illuminated CCD shielded by optical blocking filters. The parameters of the spectrographs are summarized in table 3. The observable spectral regions are 80–250 eV or 110–360 eV with the J-KAREN and 110 to 200 eV with the Astra Gemini lasers. The spectrograph used with the Astra Gemini laser consists of two channels with observation angles of α₀ ≈ 0° (forward direction) and α₀ + 0.53° (in the laser polarization plane for linear polarization); the spectrograph is operated in the 2nd and 3rd diffraction orders in order to reduce the scattered light arising from the zeroth order. The throughputs of the spectrographs, defined as products of the collecting mirror reflectivity, filter transmission, grating efficiency, and CCD quantum efficiency, are shown in figure 3(a). The mirror reflectivities and filter transmissions are calculated using the atomic scattering factors [60–62]. Measurements by the manufacturer at several wavelengths are in good agreement with the calculations, figure 3(b). The diffraction efficiencies at the incidence angles and diffraction orders used in these experiments have been measured using synchrotron radiation [63, 64]; comparison of the harmonic spectra obtained from the second and third diffraction orders (see later, in figure 6(c)), suggests that the diffraction efficiencies are valid with a relative accuracy of better than several percent. Quantum efficiency data and gain are provided by the CCD manufacturers.

Zoom In Zoom Out Reset image size

The wavelength calibration of both spectrographs is carefully performed using the spectra of Ar and Ne plasmas excited in-place by the same laser. Due to the finite resolving power of the spectrographs and the few hundred μm size of the emitting plasma region, the spectral lines of Ar and Ne ions are two to four pixels wide on the CCDs. Using many spectral lines, and determining line centres by multi-peak Gaussian fitting, we achieved sub-pixel calibration accuracy. To illustrate the method, the calibration of the spectrograph used with the Astra Gemini laser is shown in figure 4. In experiments with both setups (for J-KAREN and Astra Gemini lasers), the relative calibration error estimated from spectral line wavelength errors and difference between shots is δλ/λ ∼ a few × 10⁻⁴. Additional confidence in the spectral calibration accuracy follows from the correct wavelengths of silicon and carbon absorption edges visible in the data (e.g. figure 5 at 100 eV and figure 13 at 284 eV) and identical harmonic structures obtained from the second and third diffraction orders, figure 6(c).

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

High-frequency radiation with photon energy from ћω = 75–360 eV (corresponding to the full range of the spectrographs used) has been observed in the interaction of laser pulses (either linearly or circularly polarized) with power from P₀ = 9–170 TW, and helium gas jets with an estimated maximum electron density from n_e = 1.7 × 10¹⁹–7 × 10¹⁹ cm⁻³. The representatives of the variety of observed single-shot spectra are shown in figures 5– 7, 9, 10, 12 and 13. Simultaneously the electrons with broad spectra and energies up to a few hundred MeV have been measured. These electrons exhibit no apparent correlation with the high-frequency radiation properties. At lower plasma densities, higher-energy quasi-monoenergetic electron acceleration can be achieved, however the high-frequency radiation disappears. We also note that only usual thermal plasma emission has been observed in the cases of Xe, Ar, Ne, and CO₂ gas jets.

Zoom In Zoom Out Reset image size

The comb-like spectra comprising odd and even harmonics of similar intensity and shape are generated in the specified broad range of laser and plasma parameters, demonstrating the effect's robustness. This represents the majority of cases observed with the high-resolution spectrograph [65]. The spectra obtained with linearly polarized laser pulses are shown in figure 5 for the laser power of 9 TW, and in figure 6 for 120 TW. The spectrum obtained with the circularly polarized laser pulse is shown in figure 7 for 170 TW. In the case of circular polarization, the harmonic generation typically requires higher laser power or greater plasma density; however, additional experimental data are required to quantify the differences in optimum conditions for linear and circular polarizations.

Thanks to the high quality of the data and high calibration accuracy, several important features of the harmonic spectra are revealed. The base frequency of the spectra, ω_f, is downshifted from the laser frequency, ω₀. We attribute this to the well-known gradual frequency downshift of a laser pulse propagating in tenuous plasma [1, 39, 66, 67]. When the laser pulse transfers its energy to plasma waves in the adiabatic limit, the photon number is approximately conserved. This leads to the corresponding photon energy decrease, i.e. the frequency red-shift. This carrier frequency downshift of the laser pulse appears in the transmitted laser spectra, figure 8(a), in accordance with the observed base harmonic frequency.

Zoom In Zoom Out Reset image size

Another feature is a large number of resolved harmonics. For example, harmonics up to n_H* ∼ 126th order are clearly distinguishable in figure 5 and up to n_H* ∼ 370th in figure 7. In both cases, these numbers are close to the limits determined by the resolving power of the spectrographs. Figure 8(b) shows the resolving power of the spectrograph, which recorded the spectrum shown in figure 5, as estimated from the widths of spectral lines in neon plasma. The spectrograph enables resolving of at most ∼140 harmonics for the downshifted base frequency of ω_f = 0.885ω₀. As for figure 7, the separation between peaks in the spectrum approaches a two CCD pixels span, i.e. also the resolution limit. Despite this, the spectrum periodicity is clearly visible in the raw data figure 7(a), lineouts (b, d), and Fourier transformations (c, e).

The large number of distinctly resolved harmonic orders n_H* together with the gradual frequency downshift allow estimation of the harmonic's emission length δx. During the emission, the relative frequency drift, δω/ω, cannot be larger than (1/2n_H*), otherwise the harmonics n_H* and n_H* + 1 overlap. The relative frequency change is approximately equal to the relative energy change [66, 67], δω/ω ≈ δx/L_dep, where L_dep is the depletion length. Therefore, δx ≈ L_dep/2n_H*. We estimate L_dep ∼ 2.7 mm from the 70% energy transmission through the 0.9 mm long plasma, figure 8(a). For the harmonics spectrum shown in figure 5 this gives δx ≦̸ 12 μm, which is the upper bound, since the number of resolved harmonic orders may be limited by the spectrograph resolving power. Here we note that the longitudinal size of the harmonic's source can be much shorter than the emission length.

In the case of linearly polarized laser pulses, in some shots the spectra exhibit deep equidistant modulations, figure 9. There are even cases with the modulation depth approaching 100%. Such modulations are also visible in the case of unresolved or nearly unresolved harmonics, figure 10. The spectra with unresolved harmonics typically contain a larger photon number. These features can be explained by a longer harmonics emission time, during which the greater continuous downshift of the laser frequency leads to the observed harmonic structure blurring. In a total of 90 shots with circularly polarized 50–200 TW pulses, the large-scale spectral modulations have not been observed.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Two spectrograph channels observing the harmonics on-axis and 0.53° off-axis, figures 6 and 7, show spectra with similar intensities, figure 11, taking into account the ≈ 2.5-fold difference in acceptance angles determined during calibration using a nearly isotropic neon plasma emission. However, in some shots one of the channels shows a several times greater signal than another. These observations indicate that the angular distribution of the harmonic emission has a characteristic size of several degrees with possible direction fluctuations of the order of a few degrees. The high-resolution 2D particle-in-cell (PIC) simulations performed with the code REMP [68] (see below), suggest that the harmonics are emitted slightly off-axis and for harmonic orders of about a hundred their angular span is ∼3°.

Zoom In Zoom Out Reset image size

A remarkable feature of the harmonic spectra is their high brightness, so each of the spectra shown in this paper is the result of a single-shot acquisition. The harmonic pulse energy and photon numbers are conservatively estimated using spectrograph throughputs, figure 3(a). Examples of the spectra in the units of μJ/(eV·sr) are shown in figures 6(c), 12 and 13. These data are obtained within the respective acceptance angles of the spectrographs; thus, without the necessity of assumptions on the angular distribution and spatial beam profile. The photon yield scales favourably with the laser power, as follows from the comparison of the harmonics spectra obtained with the 9 TW and 120 TW laser powers, figure 12.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

At 120 TW laser power, the harmonic pulse energy and photon number in a unit solid angle reach 40 μJ sr⁻¹ and 2 × 10¹² photons sr⁻¹ in one harmonic at 120 eV. Assuming the 3° angular width derived from the simulations, these amount to 90 nJ pulse energy and 4 × 10⁹ photons.

The spectrograph used with the 9 TW laser has a throughput cut-off at ∼360 eV with the Pd filters, figure 3(a). The harmonic emission extends up to this photon energy, figure 13, thus reaching the 'water window' spectral range enclosed between the K absorption edges of carbon at 284 eV and oxygen at 543 eV. The carbon absorption edge, originating from the hydrocarbon optics contamination, is visible in figure 13 as a dip in the spectrum. The spectrum in the 'water window' range contains $0.8_{-0.2}^{+0.3}$ μJ sr⁻¹ pulse energy or $1.5_{-0.4}^{+0.5}$ × 10¹⁰ photons sr⁻¹. For the assumed 3° angular width, these correspond to 1.7 nJ and 3 × 10⁷ photons. The uncertainties are due to the CCD noise and 10% filter thickness tolerance.

The spectrograph throughput estimates do not take into account the absorption of the hydrocarbon contamination of optics, nor the effect of the 200 μm slit present in one of the spectrographs; measurement after the experiment shows that the geometrical throughput of the slit is ∼0.25. Thus the actual throughputs are lower, and photon numbers are larger than the conservative estimates presented here.

The temporal coherence of the radiation is deduced from the spectral fringes at the harmonic and large-scale modulation frequencies. To study the spatial coherence of the harmonics in the experiment with the J-KAREN laser with the peak power increased up to 20 TW, we use the diffractive imaging method [69]. This method is based on the analysis of the intensity distribution formed by diffraction of the investigated radiation off test objects. We employ a submicron spatial resolution LiF film detector, which is sensitive to photons with energies greater than 14 eV and has a large dynamic range in the EUV, XUV and soft x-ray spectral regions [70, 71]. The LiF detector is placed at the observation angle of 8° from the laser beam propagation, the minimum allowed by the setup. As a test object we use a square Ni mesh with a density of 70 wires/inch and a wire size of 41 μm. The distance between the source and mesh is 550 mm and between the mesh and detector is 27.2 mm. The LiF detector consists of 0.5 μm thick LiF film evaporated on glass substrates with a diameter of 20 mm and thickness of 2 mm. As LiF is not sensitive to radiation with a photon energy of less than 14 eV (λ > 88 nm), an optical blocking filter is not used. At the same time the maximum penetration depth in LiF in the spectral range under investigation, i.e. from 14–360 eV, is shorter than 500 nm, and radiation is fully absorbed by the LiF film. After the accumulation of 11 shots, the irradiated LiF film is observed using a fluorescence microscope. A portion of the photoluminescent mesh image is presented in figure 14(a). Despite the fact that many laser shots are accumulated in the LiF film, the diffraction pattern is clearly seen in the image and lineout of intensity distribution, figure 14(b).

Zoom In Zoom Out Reset image size

To model the intensity distribution in the diffractive image (for details of the algorithm see [72]) we used a harmonic spectrum at the observation angle of 8°, figure 14(c), obtained in the 2D PIC simulation, see later in section 5. Only harmonic orders above ten are taken into account, corresponding to the LiF sensitivity threshold. Each harmonic has a Gaussian spectrum with an 8% bandwidth corresponding to the initial laser pulse spectrum. The harmonics originate from a point source with a divergence of 0.034 rad corresponding to the full detector size. The comparison of modelling with experimental intensity distribution, figure 14(b), suggests that the harmonics are indeed spatially coherent. However, additional experiments are necessary to fully characterize the spatial coherence properties. We note that the modelling of the diffractive image produced with helium plasma lines (not shown) gives fringe periodicity that is not consistent with the experimentally observed intensity distribution.

The properties of the observed harmonics cannot be explained by the previously suggested generation mechanisms. Atomic harmonics cannot produce odd and even orders with linearly or circularly polarized driving pulses with the observed low sensitivity to the plasma density. Betatron emission is determined by the plasma frequency and electron energy and does not produce harmonics of optical frequencies.

The nonlinear Thomson scattering cannot produce the obtained photon number even under the most favourable assumptions. Indeed, for the shot shown in figure 5, equation (2) gives a_SF ≈ 6. Assuming an over-estimated self-focusing up to a₀ = 7, τ = 30 fs, and observation at the angle of 15°, where the nonlinear Thomson scattering exhibits a maximum at high frequencies, numerical calculations [17, 73] for a single electron at the photon energy of 100 eV give at most 2 × 10⁻¹⁰ nJ eV⁻¹ sr. The high-frequency part of the nonlinear Thomson scattering near the laser propagation direction where the experimental spectra have been measured is much weaker. Assuming the maximum electron density and neglecting nonuniform irradiance distribution, the number of emitting electrons can be estimated as N_e ≈ π d_SF² δx n_e/4 ≈ 6 × 10⁹, where δx = 12 μm is the emission length estimated above, and d_SF ≈ 5 μm is the self-focusing channel diameter, equation (3). These electrons radiate 1 nJ eV⁻¹ sr at most, i.e. 200 times smaller than the observed value. An additional difficulty arises if we take into account the inhomogeneity of the laser field across and along the self-focusing channel [74]. The emission region, where the relative emission frequency variation is smaller than 1/2n_H* = 0.004 (required to produce resolved harmonic orders up to the order n_H* = 126), is much smaller than the whole laser spot, even after self-focusing. Finally we note that other spectra shown here are even stronger than in figure 5, so the nonlinear Thomson scattering corresponds to an even smaller fraction of those spectra.

Now we summarize the experimental observations, estimates, and theoretical considerations which can shed light on the possible source of the harmonic emission. The head of the laser pulse exhibits a frequency downshift, while the tail tends to be frequency upshifted [74]; the observation of a base harmonic frequency red-shift suggests that the electric field of the laser pulse head drives the harmonics. Therefore, the harmonic-emitting electric charge density is situated in the location of the laser pulse, moving together with it. The conservative estimate presented above shows that the number of electrons encountered by the laser pulse within the emission length is at least two orders of magnitude less than the numbers necessary for single electron scattering to produce the observed spectra, suggesting that the emission grows faster than linearly with the electron number, so the emitter is at least partially coherent. The observation of both even and odd harmonics indicates the absence of the electric field sign inversion invariance at the harmonic source location, which suggests that this source is off-axis. The large scale modulations seen in the spectra produced by linearly polarized laser pulses with the depth near 100% suggest that those modulations result from interference between two nearly identical coherent sources separated in time and/or space. This can be connected with the mirror symmetry of linear polarization, in contrast to the axial symmetry of circular polarization, taking into account that such modulations have not been observed with circularly polarized laser pulses.

In order to reveal the high-order harmonics generation mechanism, we performed 2D and 3D PIC simulations of the laser–plasma interaction and high-order harmonic generation using the code REMP [68]. The ions are assumed to be immobile due to their inertia. In a broad range of laser and plasma parameters, the interaction can be described by the following scenario. The laser pulse undergoes relativistic self-focusing [30–34] accompanied by amplitude increase and spot narrowing. The laser pulse pushes the plasma electrons out, thus forming a nearly electron-free cavity [39, 75]. Due to the finite transverse dimension of the self-focused laser pulse, the electrons are pushed transversely as well. When the laser pulse and cavity radius become sufficiently narrow, the bow wave detaches from the cavity [36]. Several simulations presented here demonstrate various aspects of these processes.

In order to simulate the laser pulse propagation in millimetre-scale plasma and its evolution we performed a 2D PIC simulation, figures 15(a)–(d), in the 122λ₀ × 100λ₀ moving window with absorbing and periodical boundaries in the longitudinal and transverse directions, respectively. The mesh size is (λ₀/256) × (λ₀/64), and the quasiparticle number is 5.8 × 10⁸. The p-polarized nearly Gaussian laser pulse with the original amplitude of a₀ = 5 and duration τ = 20λ₀/c is focused with a spot of 25λ₀ to the location where the electron density is one third of the maximum density of n_e,max = 3.8 × 10⁻³n_cr. The laser pulse transverse profile and the density longitudinal profile, figure 15(e), closely match the experiment with the Astra Gemini laser. The laser pulse expels electrons from its location, creating the cavity, as seen starting from frame (a). The wake wave excited by the laser pulse undergoes transverse wave breaking due to the relativistic nonlinearity, causing the wake wave length dependence on the laser pulse amplitude [76]. The laser pulse self-focuses, as seen in frames (c) and (d), causing a multi-stream plasma flow. Electrons expelled in the transverse direction form a bow wave, which detaches from the cavity walls at about t = 530, frame (c). The cavity and bow wave, and a strong electron density concentration at their joining boundaries are visible in frame (d). At some later time the electrons are injected into the first period of the wake wave (not shown). An animation showing the simulation results can be found in the supplemental movie S1 (available at stacks.iop.org/NJP/16/093003/mmedia).

Zoom In Zoom Out Reset image size

The typical electron density distribution at the stage of the detached bow wave is shown in figure 16. This distribution is obtained in a 3D PIC simulation of a laser pulse linearly polarized along the y-axis, with a₀ = 6.6, duration of 10λ₀/c, and FWHM waist of 10λ₀, interacting with a plasma of electron density n_e = 0.001 n_cr. The simulation is performed within the 125λ₀ × 124λ₀ × 124λ₀ box with periodical transverse and absorbing longitudinal boundaries, with a resolution of λ₀/32 along the x (propagation) axis and λ₀/8 along the y- and z-axes. The number of quasiparticles is 2.3 × 10¹⁰. The laser pulse modifies the originally uniform electron density, creating the density modulations at the head of the pulse, the electron-free cavity with dense and sharp wall, and the bow wave. The highest electron density is at the joining of the cavity wall and the bow wave front, denoted as electron spike. Figure 16(a) also shows the region of intense high order harmonics generation with ω > 4ω₀. This region coincides with the joining of the cavity wall and the bow wave, i.e. with the electron density spike. The generation of frequencies greater than ω₀ can also be seen at the cavity wall and wake wave breaking and injection point; however, in the simulations presented here, at the highest frequencies corresponding to the experimentally observed high-order harmonics, the radiation from the electron spike is much stronger.

Zoom In Zoom Out Reset image size

In order to simulate the emission of very high harmonic orders comparable to those obtained in the experiment, we performed a high-resolution 2D PIC simulation, figures 17– 20. The simulation was performed in a 87λ₀ × 72λ₀ window with a grid size of λ₀/1024 along the x-axis and λ₀/112 along the y-axis, with 6 × 10⁸ quasiparticles. The parameters of the laser pulse in the simulation are close to the parameters of the self-focused Astra Gemini laser pulse, i.e. a₀ = 10, τ = 16λ₀/c ≈ 43 fs, and a spot diameter of 10λ₀ = 8 μm. The density of the plasma is n_e = 0.01n_cr = 1.7 × 10¹⁹ cm⁻³. The electron density profile shown in figure 17(a) contains all of the features described above, i.e. the cavity, cavity walls, bow wave, and the electron density spike. The density lineout across the spike is shown in figure 17(b); this data is taken ten cycles earlier than in (a) in order to demonstrate the structure of the spike responsible for the emission of the high-frequency radiation within the spectral range from 60ω₀ to 100ω₀ shown in (a) by the red colour scale.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

The profile of the electron density spike in the log–log scale is shown in figure 17(c). It is remarkable that the density growths as ∝(y – y_C)^−2/3 when the transverse coordinate approaches the spike location, y_C. This power-law density dependence on the coordinate at the joining of two density walls, caused by a multi-stream plasma flow, is very characteristic, and gives an important clue to the mechanism of the electron density spike formation, see section 7.

The spectrum of harmonics in the region denoted by the dotted rectangle in figure 17(a) is shown in figure 18. The spectrum extends up to the ∼128th order, as is allowed by the simulation resolution. The emission is slightly off-axis, figure 19, with an emission cone size of several degrees. The spectrum of lower-order harmonics at the observation angle of 8°, obtained in the same simulation, is shown in figure 14(c). The off-axis angle decreases when the harmonic order increases. The harmonic energy yield obtained in the simulation is close to the experimental data obtained with the Astra Gemini laser, figure 12.

The electron density spikes at the joining of the cavity wall and bow wave move together with the laser pulse, and oscillate with the laser carrier frequency. With a linearly polarized pulse, the electron density spike is located at the crescent curves embracing the pulse head, figure 16, following the pulse electric field mirror symmetry. Circularly polarized pulses produce annular electron density spikes. At any time, the spike consists of different electrons, flowing in different directions through a curve embracing the laser pulse head, with momenta induced by the laser pulse field. As a strongly localized charge, oscillating with relativistic velocity, the spike emits electromagnetic radiation in accordance with classical electrodynamics [73, 77], in the form of harmonics with the base frequency equal to the local carrier frequency of the laser pulse. Owing to a strong spike localization, the short-wavelength radiation undergoes a constructive interference for electrons located at distances below a half-wavelength from each other. This radiation is coherent, and its intensity scales as a square of the number of emitting particles. We note that for constructive interference in the forward direction, it is sufficient that only the longitudinal size of the source is small enough in comparison with the emitted wavelength. In the simulations presented above, the number of electrons within the interval of width λ₀/100 near the spike is N_e ∼ 10⁶–10⁷. With the emission intensity growing as N_e², the spikes provide a signal level close to the experiment. The divergence of radiation is determined by the source size δ in the direction perpendicular to the observation direction, and by the source velocity with respect to the observer. The cusps move with relativistic velocity, which decreases the divergence angle by the Lorentz factor γ, i.e. Δθ ∼ λ/γδ.

In the region of the spike, the electron density is at a maximum in the direction of the electric field vector, since a greater amount of electrons are expelled in this direction. In the case of linear polarization, the maxima appear at a pair of opposite points, oscillating along the polarization axis. In the case of circular polarization, the maximum is located on a point travelling on a helical trajectory, following the electric field vector of the laser pulse. The coherent emission intensity strongly depends (as N_e²) on the number of electrons, therefore the higher the harmonic order, the smaller the emitting region near the electron density maximum.

The described mechanism is consistent with all of the features of the experimentally observed high-order harmonics. The strong localization of the electron density spikes removes the discrepancy between the estimated number of electrons which are encountered by the laser pulse, and the number of electrons emitted in a single-particle scattering regime required for the experimentally observed harmonic spectra. The electron density spikes are off-axis and therefore produce harmonics with both even and odd orders having the same shape and intensity.

The new harmonics generation mechanism is dominant when the laser pulse creates a bow wave while the accelerated electrons cannot produce coherent radiation. This is realized (1) in underdense plasma, n_e < n_cr, when (2) the laser dimensionless amplitude is much greater than 1, a₀ ≫⃒ 1; (3) the laser pulse length is comparable with the Langmuir wavelength (so that the wake wave undergoes a transverse break in the first period suppressing electron self-injection into the wake wave); (4) the laser is tightly focused, making a spot with the diameter of $d\leqslant \frac{2\lambda }{\pi }{{\left( {{a}_{{\rm 0}}}\frac{{{n}_{{\rm cr}}}}{{{n}_{{\rm e}}}} \right)}^{1/2}}$ [36]. Under our experimental conditions, the large amplitude and tight focusing are achieved owing to relativistic self-focusing. These conditions remove or supress previously known harmonics generation mechanisms: atomic harmonics do not occur in plasma; accelerated electrons have too small a density and too large an energy spread so betatron emission is not coherent; all the electrons encountered by the tightly focusing laser pulse are transversely expelled, forming the bow wave and flow through the cusps (density spikes), so that their coherent cusp emission is always stronger than the nonlinear Thomson scattering.

The formation of a strongly localized electron density spike and its peculiar properties, namely robustness to oscillations imposed by the laser field and the −2/3 power-law density dependence, can be explained with the help of catastrophe theory [27, 28].

The laser pulse expels electrons in both the longitudinal and transverse directions, while the electrostatic potential of the non-compensated background of slowly responding ions pulls the electrons back. This leads to a multi-stream motion of the electron fluid [76, 78], as is seen in figure 20, corresponding to the simulation shown in figure 17. As is shown in figure 20, the quasiparticles in the phase space (x, y, y – y₀) are well-ordered according to their initial position y₀ and are situated on a surface. Here $y-{{y}_{0}}=\int \left( {{p}_{\bot }}/{{m}_{{\rm e}}}{{\gamma }_{{\rm e}}} \right){\rm d}t=\left\langle {{p}_{\bot }}/{{m}_{{\rm e}}}{{\gamma }_{{\rm e}}} \right\rangle$ is the displacement along the y-axis from the initial coordinate y₀, and γ_e is the electron Lorentz factor. The electron density increases abruptly at the location of the cavity walls and bow wave boundaries, representing catastrophes of the plasma flow.

The formation of catastrophes is illustrated by the following one-dimensional model of electron motion along the transverse coordinate y, figure 21 and supplemental movie S2 (available at stacks.iop.org/NJP/16/093003/mmedia). Initially, the electrons with spatially varied momentum ${{p}_{\bot }}$ homogeneously fill the y-axis, figure 21(a). In the phase space $(y,{{p}_{\bot }})$, the phase 'volume' occupied by the electrons is represented in an implicit form by a curve $S\left( t,y,{{p}_{\bot }} \right)=0$. Due to the momentum inhomogeneity, the order of electrons eventually changes, thus a multi-stream flow arises and the curve S becomes folded. The electron density is obtained by a projection of the electron phase space onto the y-axis,

$$\begin{eqnarray}&&{{n}_{{\rm e}}}\left( t,y \right)=\int \limits_{-\infty }^{+\infty } f\left( t,y,{{p}_{\bot }} \right){\rm d}{{p}_{\bot }},\end{eqnarray} \tag{ 4 }$$

where $f\left( t,y,{{p}_{\bot }} \right)\propto \delta \left[ S\left( t,y,{{p}_{\bot }} \right) \right]$ is the particle density in the phase space $\left( t,y,{{p}_{\bot }} \right)$ and δ is the Dirac delta-function. At the position of a fold, the density rises sharply. The emergence of a singularity in a flow is called a catastrophe. Note that in the case of the fluid approximation, projection (4) gives an infinite density, while in the case of a finite particle number, as in the supplemental movie S2 (available at stacks.iop.org/NJP/16/093003/mmedia), the density, although exhibiting a sharp rise, remains finite.

Zoom In Zoom Out Reset image size

The singularity type is determined by the structure of the phase 'volume'. By analogy with the observation that any smooth function at a regular point can be approximated by a linear form, in the general case, without special conditions on the shape of the curve S, a fold can be approximated by a portion of a parabola,

$$\begin{eqnarray}&&{{y}_{{\rm f}}}-y\propto {{\left( {{p}_{\bot {\rm f}}}-{{p}_{\bot }} \right)}^{2}},\end{eqnarray} \tag{ 5 }$$

figure 21(c), while the portion near the curve S inflection can be approximated by a cubic,

$$\begin{eqnarray}&&{{y}_{{\rm c}}}-y\propto {{\left( {{p}_{\bot {\rm c}}}-{{p}_{\bot }} \right)}^{3}},\end{eqnarray} \tag{ 6 }$$

figure 21(b). Therefore, the projection (4) gives the density singularity of

$$\begin{eqnarray}&&{{n}_{{\rm eF}}}\propto \frac{\partial S/\partial y}{\partial S/\partial {{p}_{\bot }}}=\frac{{\rm d}{{p}_{\bot }}}{{\rm d}y}\propto {{\left( {{y}_{{\rm f}}}-y \right)}^{-1/2}}\end{eqnarray} \tag{ 7 }$$

at the fold, and of

$$\begin{eqnarray}&&{{n}_{{\rm eC}}}\propto {{\left( {{y}_{{\rm c}}}-y \right)}^{-2/3}}\end{eqnarray} \tag{ 8 }$$

at the inflection point. In a higher dimensionality space $(t,y,{{p}_{\bot }})$, the inflection point corresponds to a merger of two folds of a surface given by $S\left( t,y,{{p}_{\bot }} \right)=0$, figure 21(d). In the vicinity of the merger point, $({{t}_{{\rm c}}},{{y}_{{\rm c}}},{{p}_{\bot {\rm c}}})$, the surface S can be represented implicitly by

$$\begin{eqnarray}&&{{\left( {{p}_{\bot }}-{{p}_{\bot {\rm c}}} \right)}^{3}}+\alpha {\mkern 1mu} \left( t-{{t}_{{\rm c}}} \right){\mkern 1mu} \left( {{p}_{\bot }}-{{p}_{\bot {\rm c}}} \right)+\beta \left( y-{{y}_{{\rm c}}} \right)=0\end{eqnarray} \tag{ 9 }$$

with some constants α ≠ 0 and β ≠ 0 ([28], chapter 5 paragraph 2), in analogy with the approximations by a parabola and a cubic used above. Here the projection onto the (t, y) plane produces two singular lines joining at a point with a singularity of higher order. If we apply small smooth perturbations to the electron phase 'volume' represented by the curve S, the fold and the inflection still remain, and can be approximated locally by a parabola and a cubic, respectively. In other words, the corresponding singularities are robust with respect to the perturbations. The rigorous proof of this fact is given in catastrophe theory ([27], chapter 2; [28], chapter 6, paragraphs 4 and 5), which is a well-established branch of mathematics studying the formation of singularities in dynamical systems and in mappings. Here catastrophe theory establishes two structurally stable singularities, the fold catastrophe (A₂, according to V I Arnold's classification) and the cusp catastrophe (A₃).

Catastrophes of the same types as in the model in figure 21 appear in the laser plasma during cavity and bow wave formation [36], figure 22. Here the x coordinate of the laser pulse propagation corresponds to the reversed time in figure 21. The cusp catastrophe, appearing at the joining of the cavity walls and bow wave boundaries, figure 22(a), causes the formation of the electron density spike, figure 22(b), seen in the 3D PIC simulations, figure 16(b). The spike is located in a ring surrounding the cavity head, moving together with the laser pulse. Being pushed by the electromagnetic field of the laser pulse, the electrons form a flow modulated with the local laser frequency. Correspondingly, the electron spike undergoes periodic oscillations and emits electromagnetic radiation with high-order harmonics of the local laser frequency. Catastrophe theory ensures that the described singularities are structurally stable which means that they are insensitive to perturbations imposed by the laser pulse and other factors such as electron density fluctuations. Singularities stronger than the cusp can momentarily appear [78], however, they are not stable against the perturbations. Other singularities in collisionless plasma are discussed in [76, 78].

Zoom In Zoom Out Reset image size

Although the electron spike consists of different electrons constantly flowing through it, here we use a simplified model where the spike is considered as a point charge. The point charge approximation becomes better for higher order harmonics. According to classical electrodynamics [73], an oscillating charge emits harmonics up to the critical frequency ω_c, corresponding to the critical order n_Hc, proportional to the cube of the particle energy ε_e. Assuming mostly transverse motion of electrons, where ε_e ≈ a₀m_ec², we obtain

$$\begin{eqnarray}&&{{n}_{{\rm Hc}}}=\frac{{{\omega }_{{\rm c}}}}{{{\omega }_{{\rm f}}}}\tilde{\ }a_{0}^{3}.\end{eqnarray} \tag{ 10 }$$

At harmonic orders substantially higher than this, the spectrum exponentially vanishes, as observed in experiments. Assuming the stationary self-focusing condition, equation (2), we obtain the following scaling of the critical harmonic order, expressed via the laser peak power P₀ and electron density n_e:

$$\begin{eqnarray}&&{{n}_{{\rm Hc}}}\tilde{\ }\frac{{{P}_{0}}}{{{P}_{{\rm c}}}}\frac{{{n}_{{\rm e}}}}{{{n}_{{\rm cr}}}}\approx 500\left( \frac{{{P}_{0}}}{1{\rm PW}} \right)\left( \frac{{{n}_{{\rm e}}}}{{{10}^{19}}{\rm c}{{{\rm m}}^{-3}}} \right){{\left( \frac{{{\lambda }_{{\rm 0}}}}{1\mu {\rm m}} \right)}^{2}},\end{eqnarray} \tag{ 11 }$$

where P_c is given by equation (1). We see that a 100 TW laser pulse with a wavelength of 1 μm, propagating in plasma with the density of 5 × 10¹⁹ cm⁻³ can produce at least about 250 harmonic orders.

The total energy ε_S emitted by the spike can be obtained using formulae from section 14.2 of [73]. It is proportional to the forth power of the charge energy, (a₀m_ec²)⁴, and the charge squared, i.e. e²N_e², and the emission time τ_H:

$$\begin{eqnarray}&&{{\varepsilon }_{{\rm S}}}\approx \frac{{{e}^{2}}N_{{\rm e}}^{{\rm 2}}a_{0}^{4}\gamma \omega _{{\rm f}}^{2}{{\tau }_{{\rm H}}}}{8c}\propto N_{{\rm e}}^{{\rm 2}}P_{{\rm 0}}^{{\rm 4}/3}n_{{\rm e}}^{{\rm 5}/6}\omega _{{\rm f}}^{{\rm 1}/3}{{\tau }_{{\rm H}}}.\end{eqnarray} \tag{ 12 }$$

Here the harmonic base frequency ω_f ≈ ω₀. We take into account the relativistic motion of the charge, together with the laser pulse, and introducing factor γ, the Lorentz factor associated with the spike velocity, which is close to the laser pulse group velocity in plasma. The latter can be estimated in the first approximation as γ ≈ (n_cr/n_e)^1/2. The stationary self-focusing condition, equation (2), and the relation between the emitting charge Lorentz factor and plasma density allows one to reveal the harmonic emission energy scaling. Equation (12) for P₀ = 120 TW, n_e = 1.9 × 10¹⁹ cm⁻³, N_e = 10⁷ gives the total emitted energy of ∼1 mJ. Here we estimate a₀ = a_SF = 12, γ = 10, τ_H = 0.1 fs, using equation (2) and PIC simulation results. This energy corresponds to ∼60 μJ sr⁻¹ in one harmonic (assuming that the divergence angle is ∼1/γ and the energy per harmonic is about total energy/a₀³), which compares well to the experimental value of 40 μJ sr⁻¹ at 120 eV and the simulation result of 30 μJ sr⁻¹ at harmonic order 100. For P₀ = 9 TW, n_e = 2.7 × 10¹⁹ cm⁻³, N_e = 10⁶ the total emitted energy is ∼3 μJ. Here a₀ = a_SF = 6, γ = 8, τ_H = 0.7 fs. This energy corresponds to ∼1 μJ sr⁻¹, in agreement with the experimental observation (∼1 μJ sr⁻¹ at 80 eV).

We have experimentally discovered a new regime of high-order harmonic generation using multi-terawatt femtosecond lasers irradiating gas jet targets. The harmonics are generated by linearly as well as circularly polarized laser pulses in wide ranges of laser power and plasma density. The comb-like harmonic spectra containing even and odd orders extend to the 'water window' spectral region. The harmonics are emitted slightly off-axis with the emission angle reducing at higher orders. The harmonic yield grows faster than linearly with the laser power. We have proposed the harmonic generation mechanism based on self-focusing, following the nonlinear wake and bow waves excitation which produce electron density spikes coherently emitting high-order harmonics. Catastrophe theory explains the electron spike sharpness and stability. The maximum harmonic order scaling with the laser intensity will allow the attainment of the keV range. Our results open the way to a compact coherent x-ray or XUV source built on a university laboratory scale repetitive laser and an accessible, replenishable, and debris-free gas jet target. This will impact many areas requiring a bright compact x-ray or XUV source for pumping, probing, imaging, or attosecond science.

We acknowledge the financial support from MEXT, Japan (Kakenhi 20244065, 21604008, 21740302, 23740413, 25287103, 25390135 and 26707031), the JAEA Research and Development Adjustment Funding for Exploratory Research, the STFC facility access fund (UK), RFBR-JSPS collaboration program, Russia and Japan (RFBR Grant 14-02-92107) and Program for Basic Research No. 2 of the Presidium of the RAS, Russia.