Comparison of strong-field ionization models in the wavelength-scaling of high harmonic generation

: We report the use of wavelength-tuneable laser pulses from an optical parametric ampliﬁer to generate high-order harmonics in a range of noble gases. The variation of the harmonic cut-oﬀ wavelength and phasematching pressure with gas species and fundamental wavelength were recorded. The experimental results are compared to a phenomenological model of the harmonic generation process, incorporating two separate models of photo-ionization. While the calculated phasematching pressure is generally insensitive to the ionization model, for the harmonic cut-oﬀ we obtain superior agreement between experiment and theory when the Yudin-Ivanov (YI) ionization model is used, compared to the commonly utilised Ammosov-Delone-Krainov (ADK) model.


Introduction
The interaction between intense, femtosecond duration laser pulses and matter is a widely studied phenomena, with the process underpinning many cutting-edge techniques in optical and atomic physics, including laser-induced electron diffraction [1] and high harmonic generation (HHG) [2]. As ultrafast laser technology has developed to encompass wavelengths spanning the vacuum ultraviolet [3] to mid-infrared [4] spectral regions, the challenge of accurately modelling light-matter interaction at high intensities has concurrently increased.
In the case of high harmonic generation, ionization by ultrafast laser pulses is the first step of the well-known, semi-classical, three-step model [5]. In this model, a linearly polarized laser field first ionizes an atom, with the subsequent laser-driven electron dynamics giving rise to the possibility of photo-electron -ion recombination, resulting in the emission of high energy photons at odd harmonics of the fundamental driving frequency. While the ionization rate can be calculated through numerically solving the time-dependent Schrödinger equation, several approaches have been developed to simplify the calculation and produce analytic expressions for the instantaneous ionization fraction [6][7][8]. Typically, HHG experiments have been modelled on the basis of ionization via tunnelling through a quasi-static barrier, and the corresponding rate of ionization is commonly calculated using the Ammosov-Delone-Krainov (ADK) model [8]. However, this approximation is not valid for all combinations of fundamental wavelength, peak intensity and gas species. Hence, alternate ionization mechanisms must be considered, such as multiphoton ionization, where multiple photons of the driving field are absorbed, promoting an electron into the continuum. In this case the instantaneous ionization rate can be calculated using the Yudin-Ivanov (YI) model [6], which incorporates both quasi-static tunnelling and multiphoton ionization, without resorting to averaging over the laser cycle, as commonly found in older approaches [7].
Typically, the Keldysh parameter is used to distinguish whether tunneling or multiphoton ionization is the dominant ionization mechanism. It is defined as: (1) where I p is the ionization potential, and U p is the ponderomotive potential: U p = e 2 λ 2 0 I 0 8π 2 c 3 0 m e (2) where I 0 is the laser intensity, λ 0 is the laser wavelength, e and m e are the electron charge and mass, respectively, c is the speed of light in vacuum and 0 is the permittivity of free space. The Keldysh parameter can be understood intuitively as the ratio of the electron tunnelling time with the laser period. Accordingly, tunnel ionization will dominate when γ 1, whereas multiphoton ionization will play a primary role when γ 1. Although the Keldysh parameter offers a convenient classification based on experimental conditions, it has been noted by Reiss that the simple dichotomy between multiphoton and tunneling can be misleading [9]. He showed that there exist scenarios where the Keldysh parameter may correspond to multiphoton ionization (γ > 1), yet physically, ionization can only occur through tunneling. Indeed, Reiss's results demonstrate that a general theory of ionization must account for both ionization mechanisms.
Recent experiments investigating the wavelength-scaling of HHG have shown that the choice of ionization model can be crucial for successful interpretation of experimental results. For example, Gkortsas et al., [10] reported that for HHG driven by 400 nm wavelength laser pulses, the YI model gave superior agreement with experimental data, compared to the ADK model. Further, Shiner et al. [11] used the YI model to relate the measured ion yield to the focused laser intensity, allowing comparison of the measured harmonic cut-off wavelength with calculated values in the case of a 1800 nm driving laser wavelength. In both examples, ionization is studied via detection of secondary emission, be it photons from HHG or the cation yield, allowing for a simplified experimental arrangement compared to that required for direct detection of photo-electrons.
In this paper we investigate the variation of the high harmonic phasematching pressure and shortest detectable harmonic wavelength (cut-off) as a function of fundamental wavelength and gas species. We compare the ADK and YI ionization models via a 1-D phenomenological model of HHG. We find that for 522 nm and 1300 nm wavelength driving pulses, the choice of ionization model has little impact on the predicted phasematching pressure. However, significant differences can arise when comparing the calculated cut-off wavelength with its experimental counterpart. In nearly all cases considered, using the YI model, instead of ADK, yields a more accurate estimate of the harmonic cut-off.

Experiment
We investigate high-order harmonics generated by laser pulses produced by a custom-built, three stage, synchronously pumped, optical parametric amplifier (OPA) [12]. A schematic of the OPA is shown in Fig. 1. The OPA is pumped with pulses from a Ti:sapphire regenerative amplifier (wavelength of 800 nm, pulse duration of 40 fs, pulse energy of 3 mJ, and repetition rate of 1 kHz). The final stage of the OPA can be configured in one of two ways: from the amplified signal beam, pulses with a centre wavelength tunable in the short-wave infrared region (λ 0 = 1200 − 1550 nm) can be produced; or, alternatively, the final stage can be set to sum frequency generation of the signal and pump, to produce visible pulses with centre wavelength tunable across λ 0 = 485 − 530 nm. In all cases the resultant pulses have a peak power in excess of 1 GW direct from the OPA, without additional pulse compression. The duration of the pulses produced by the OPA was measured using a home-built Frequency Resolved Optical Gating (FROG) device, comprised of a wavefront division interferometer, silver-coated off-axis parabolic mirror, and fibre coupled spectrometer (Ocean Optics USB 4000). For the short-wave infrared output, the SHG-FROG technique was used [13], with a BBO crystal placed at the focus of the off-axis parabola in the FROG apparatus. For the visible output pulses, the SD-FROG technique was used, with the BBO crystal replaced by a thin, glass plate [14].
For the experiments described in this article, high order harmonics were generated by focusing laser pulses from the OPA into a gas cell backed continuously with different noble gases. The gas cell was housed within a vacuum chamber, with a background pressure below 0.07 mbar for the highest gas cell pressures utilised. The gas cell was made by pressing a hollow, thin-walled, nickel tube to a thickness of < 1 mm. The gas cell thickness was less than the Rayleigh range of the focused laser for all but one case considered below (λ 0 = 522nm and neon gas). The focused laser beam drilled entrance and exit holes into the cell. After the gas cell, the residual fundamental beam was filtered using metallic foils (either Al or Zr, depending on spectral region of interest) before the harmonic spectrum was recorded on a home-built flat-field spectrometer, comprised of a variable line-spaced grating and x-ray sensitive CCD (Andor DO440-BN).
In this study we investigate the driver wavelength and gas species dependence of two experimental quantities: the phasematching pressure (P m ) and effective harmonic cut-off wavelength (λ min ). We define the phasematching pressure as the gas cell backing pressure for which the spatially and spectrally integrated harmonic intensity (hereafter "the harmonic signal") is first maximised. It was found that comparing spectral integration over a single harmonic order or the entire detected harmonic bandwidth does not appreciably alter the recorded phasematching pressure, since P m is only weakly dependent on harmonic order. The gas pressure was measured near to the gas cell, to ensure that the measured value was close to the pressure in the cell itself. In Fig. 2(a) the harmonic signal is plotted as a function of backing pressure for a fundamental wavelength of 522 nm, for the case of argon (orange squares) and krypton (green circles). The recorded harmonic signal is clearly maximised at different backing pressures (P) for the two gas species. To extract the phasematching pressure, a function I = a × sinc k(P − P m ) 2 is fitted to the data, where a, k and P m as fitting parameters and I is the measured harmonic intensity. The phasematching pressure corresponds to P = P m . Owing to the low absorption of the short, cut-off wavelengths we consider, reabsorption is neglected in the determination of P m . For a given set of experimental conditions, we define λ min as the shortest harmonic wavelength Variation of the measured harmonic signal with gas cell pressure, for λ 0 = 522 nm and krypton (green circles) and argon (orange squares) as the generating gas. b) High harmonic spectra recorded for λ 0 = 1300 nm with either xenon (blue line), krypton (green line) or argon (yellow line) as the generating gas and Al foils for filtering the fundamental. The orange line is the HHG spectra recorded in argon with Zr foils used in place of the Al foils. c) High harmonic spectra recorded for λ 0 = 522 nm, and with argon (blue line) and neon (green line) as the generating gas. Al foils were used for filtering in both cases. d) Calculated ionization rate in argon, for λ 0 = 522 nm, using either the ADK (solid blue line) or YI (dashed red line). In both cases the peak intensity was 2.8 × 10 14 W/cm 2 and the pulse duration was 150 fs.
recorded. For consistency, λ min is evaluated at a backing pressure equal to P m . We have measured both P m and λ min for two different driving wavelengths and a variety of noble gas species.
With the output wavelength of the OPA tuned to 1300 nm, the pulse duration was measured using a SHG-FROG device and found to be 108 fs in duration. Pulses with an energy of 240 µJ were focused to a spot size of ≈ 25.1 µm inside the vacuum chamber using a silver coated, spherical mirror with focal length f = 150 mm operated at near normal incidence. The gas cell was placed close to the focal plane, with the fundamental beam drilling entrance and exit holes in the gas cell whilst under vacuum. Before data was recorded, the cell was repeatedly translated longitudinally to ensure the holes were sufficiently large to avoid clipping of the focused fundamental. Recorded high harmonic spectra are shown in Fig. 2(b), where the gas cell was backed, separately, by xenon, krypton and argon. The same laser parameters were used for the three different gas species. To observe the harmonic spectrum beyond the aluminum L-edge (≈ 17.1 nm), Zr foils were used in place of Al. With Zr foils and with argon as the generating gas, the harmonic cut-off extended to 15.3 nm (the 85th harmonic order).
When the output of the OPA was tuned to a wavelength of 522 nm, the pulse duration was measured to be ≈ 150 fs using the SD-FROG technique. The pulse energy was 135 µJ, measured immediately after the focusing optic. For Kr and Ar the beam was focused to a size ≈ 15.2 µm using an un-coated, fused silica lens with focal length f = 75 mm, while for Ne, tighter focusing using an f = 50 mm focal length lens resulted in a spot size of ≈ 10.1 µm. High harmonic spectra produced by argon (blue line) and neon (green line) are shown in Fig. 2(c). In the case of neon, the cut-off extended to 21.2 nm (the 25th harmonic order). In the specific case of krypton, the high harmonic spectrum was not recorded for this fundamental wavelength since the harmonic cut-off wavelength was too long to be detected by the spectrometer. The phasematching pressure was determined in this specific case by recording the backing pressure-dependence of the harmonic signal detected from the zero-order grating reflection in the spectrometer.

High Harmonic Generation Model
In order to calculate P m and λ min we use a 1-D, analytic, phenomenological model first described in [12]. Phasematched generation corresponds to the case of zero net dispersion: where k is the wavevector and ∆k G accounts for the geometric dispersion resulting from the Gouy phase in a free focus geometry. The possibility of transverse phase-matching [15] is not accounted for in our model.
In the loose focusing limit, where the longitudinal extent of the generation region is much less than the Rayleigh range of the fundamental, ∆k = 0 may be satisfied by balancing the positive (neutral gas) and negative (plasma) contributions to the dispersion. However, this can only be achieved up to a maximum ionization level (plasma density), beyond which the plasma contribution to the dispersion is larger than the neutral gas contribution. This maximum ionization level occurs at the critical ionization fraction, (η crit ), given by [16]: where N 0 is the number density at atmospheric pressure, r e is the classical electron radius, and ∆n is the difference between the refractive index evaluated at the fundamental and harmonic frequencies. Strictly, equation 4 is valid in the loose-focusing limit, where the Gouy phase contribution to the dispersion is negligible. Outside of this limit, equation 4 is an upper limit on the critical ionization fraction and is generally a good approximation of η crit for all but the tightest focusing geometries.
Since the ionization fraction η increases with time during the passage of the laser pulse, η crit occurs at a time within the pulse after which phasematching, and hence efficient harmonic generation, can no longer be achieved. This situation has been referred to as transient phasematching [17], since ∆k = 0 is only satisfied over a limited period of time prior to when η = η crit occurs.
From equation 3 it is possible to derive an expression for the gas pressure which satisfies the phasematching condition (∆k = 0), i.e. the phasematching pressure. For a driving laser with a Gaussian transverse profile the phasematching pressure can be written as [18,19]: where P 0 is the standard pressure and w 0 is the laser spot size. Equation 5 assumes that generation occurs at the focal plane of the fundamental and is valid within the paraxial approximation, for both loose and tight focusing geometries [19].
For a single atom, the harmonic cut-off wavelength λ min can be calculated using [5]: where h is Planck's constant, κ is a constant and I 0 is the peak laser intensity. It is known that the experimentally measured harmonic cut-off wavelength λ min is nearly always longer than the theoretical single atom cut-off, calculated using equation 6 (i.e. λ min ≥ λ min ). [20] One explanation for this difference is that the highest harmonic orders are generated by the highest laser intensities which occur at the temporal peak of the laser pulse. These conditions correspond to high ionization fractions, and hence phasematching is often not possible. This leads to a proportionally smaller and often undetectable signal from the non-phasematched harmonic orders. Consequently, an analogous expression for the effective cut-off can be written as: where I eff is an effective intensity, such that I eff ≤ I 0 . Within our model, λ min is found from equation 7, with the laser intensity evaluated at the moment in time that the ionization fraction is equal to the critical ionization fraction. In order to calculate P m and λ min from equations 5 and 7 we first need to define I eff . For a given harmonic wavelength λ q > λ min , we take I eff to be the laser intensity evaluated at the time before the temporal peak of the laser pulse, when the intensity is just high enough to generate wavelength λ q according to the expression for the single atom cut-off (equation 6), i.e. I eff (λ q ) = I 0 (t eff ), where t eff is the moment in time that λ q is first generated. An ionization fraction can be associated with λ q through: η q = η(t eff ). From this process, P m may be calculated using equation 5. Further, this approach allows us to define an effective Keldysh parameter: i.e. the Keldysh parameter evaluated at a time t eff .

Comparing Ionization Models
To illustrate the difference between the two ionization models, in Fig. 2(d) the calculated, instantaneous ionization rate in argon is plotted for both ADK (solid blue line) and YI (dot-dashed red line) models, assuming a peak intensity of 2.8 × 10 14 W/cm 2 , λ 0 = 522 nm and pulse duration of 150 fs. This corresponds to a Keldysh parameter of 1.05 at the temporal peak of the pulse. In the case of the ADK calculation, the ionization rate is localised to the half-cycle peaks of the laser pulse, separated by T 2 = 0.87 fs, where T is the laser period. The same peaks in the ionization rate are present in the YI calculation, however now there is also a non-zero ionization rate between half-cycle peaks, leading to a larger overall ionization fraction. For the YI model, the contribution between half-cycle peaks is attributed to multi-photon ionization. This mechanism is not accounted for in the ADK model, leading to an underestimation of the true ionization rate. Orange cross-hatched boxes show where both models give the same result. b) As a) but for λ 0 = 1300 nm and either Xe, Kr, or Ar backing the gas cell. c) Experimental harmonic cut-off (blue circles) for λ 0 = 522 nm and the gas cell backed with Kr, Ar or Ne. In the case of Kr only, the cut-off could not be measured and the data-point represents the lower bound for λ min according to the wavelength coverage of the XUV spectrometer. The calculated cut-off is shown for: the YI model (red box) or ADK model (yellow box). d) As c) but for λ 0 = 1300 nm and Xe, Kr or Ar backing the gas cell.
We compare the experimentally measured P m and λ min to values calculated using our phenomenological model. The laser peak intensity is a key parameter in the model. Therefore we calculate P m and λ min for a range of peak intensities, corresponding to the experimental uncertainty in this parameter. For λ 0 = 522 nm, the range is: I 0 = 2.1 − 2.9 × 10 14 W/cm 2 , except for when Ne is the gas species, where it is I 0 = 4.6 − 6.6 × 10 14 W/cm 2 instead. For λ 0 = 1300 nm, the range is: I 0 = 1.6 − 2.5 × 10 14 W/cm 2 . The experimental laser parameters themselves are stated in Table 1. Table 1. Experimental laser parameters and associated uncertainties. For the case of λ 0 = 1300 nm, the spot size was inferred from the collimated beam size and lens focal length, with the unknown M 2 beam parameter contributing to the uncertainty. The peak intensity is calculated from the experimental laser parameters.
In Figs 3(a) and (b) we plot the measured phasematching pressures for λ 0 = 522 nm and λ 0 = 1300 nm, respectively, for the different noble gases tested. In Figs 3(c) and (d) the gas species dependence of the observed high harmonic cut-off wavelength is shown for λ 0 = 522 nm and λ 0 = 1300 nm, respectively. For each driving wavelength, the generating gases are shown in order of increasing ionization potential. For both driving wavelengths, the phasematching pressure increases and cut-off wavelength decreases with increasing ionization potential. In Fig.  3 we overlay calculated values for P m and λ min using equations 5 and 7, respectively. For the calculated values of P m and λ min we consider two scenarios: the ionization fraction calculated using the YI model (red hatched boxes) and the ADK model (yellow hatched boxes).
For the phasematching pressure calculations, we choose λ q = 73 nm when λ 0 = 522 nm, and λ q = 31 nm when λ 0 = 1300 nm. In all cases λ q < λ min . Generally, the calculated values for P m do not strongly depend on the choice of λ q : for η η crit , the calculated P m varies by less then 1% for adjacent harmonic orders.
In spite of the simplicity of the calculation, we find that experimental and calculated values of P m agree well, shown in Figs 3(a) and (b), with little difference between YI and ADK calculations for the two highest ionization potential gases under consideration for each driver wavelength.
For the case of Kr, the cut-off wavelength was out of the range of the spectrometer used so the cut-off could not be determined experimentally. Hence the errorbar in Fig. 3(c), for the case of Kr only, represents the range of possible values of λ min , given the longest possible wavelength the spectrometer was able to resolve in first order (≈ 48 nm). Effective Keldysh parameters, calculated according to equation 8, for the laser parameters presented in table 1, are shown in table 2, for both YI and ADK models. The YI γ eff values are larger than those from the ADK calculation, owing to a smaller I eff , due to η crit being achieved earlier in time for the YI ionization fraction compared to ADK.

Discussion
Considering the laser parameters used in the calculations in section 4, the estimated peak laser intensity is 2.5 × 10 14 W/cm 2 for λ 0 = 522 nm (when Ar and Kr were used) and 2.2 × 10 14 W/cm 2 for λ 0 = 1300 nm (when Xe, Kr and Ar were used). In the case of argon, this corresponds to a single atom cut-off of λ min = 34.6 nm, for the 522 nm driver, and λ min = 9.9 nm, for the 1300 nm driver. As per Figs 3(c) and (d), the experimentally measured cut-off is much closer to the single atom value in the case of the visible driver (λ min = 35.0 nm) than is the case of the infrared driver (λ min = 15.3 nm). An explanation for this wavelength-dependent behaviour has been described previously [21,22]: visible pulses experience lower plasma dispersion, compared to infrared pulses because the fundamental frequency is further from the plasma frequency in the generating medium. Consequently, harmonic generation can be phasematched at a higher ionization fraction, occurring later in the pulse and hence at a higher intensities. Therefore the effective cut-off is closer to the single atom cut-off in the case of visible pulses. The wavelength-dependence of the phasematching process has been utilised recently to efficiently produce high harmonics from UV driving pulses in multiply ionized gas media [3]. For the case of λ 0 = 1300 nm, our results indicate that the ADK calculation underestimates the ionization fraction, compared to the YI model, even though γ eff <1 such that tunnel ionization is expected to dominate. This is evident in Fig. 3(d), where the ADK rate leads to a calculated harmonic cut off which is substantially smaller than both experiment and the calculation using the YI model. Considering the calculation for λ min , η crit occurs later in time but prior to the peak of the pulse. Hence a higher effective intensity is reached for ADK compared to YI. This leads to a calculated harmonic cut-off wavelength that is considerably shorter than measured in the case of ADK. Since the YI calculation better matches the actual ionization rate, a lower effective intensity is calculated, and hence a harmonic cut-off which more closely matches experiment. The difference between YI and ADK is less dramatic for the visible pulses because harmonic generation is phasematched at a comparatively higher ionization fraction (η crit ≈ 8.5%, compared to η crit ≈ 1.3% for the infrared pulse). We note in Fig. 3(d) that the experimentally measured cut-off for Kr and Ar, is at a shorter wavelength compared to our calculated values. This difference is not accounted by the uncertainties in the experimental or calculated values. Rather, we attribute this discrepancy to propagation effects in the ionizing gas medium, (e.g. pulse broadening/compression and/or self-focusing) which are neglected in our calculations. Such effects could have lead to a larger laser intensity than that used for the calculation and are more prevalent at longer wavelengths and moderately high gas pressures, as is the case in Fig.  3(d).
For the fundamental wavelengths under investigation, and the two highest ionization potential gases considered in each case, we find that the calculated range of phasematching pressures do not depend on the ionization model chosen [see Figs 3(a) and (b)]. This can be understood by rewriting equation 5 as: can be thought of as the phasematching pressure evaluated when η = 0 (i.e. when the Gouy and neutral gas contributions to the dispersion are balanced). In the limit where η η crit , equation 9 can be approximated with a binomial expansion as: In this limit, corresponding to the case where the plasma contribution to the dispersion is comparatively weak, changes to η are a small perturbation to the calculated phasematching pressure. Therefore, differences in η from the different ionization models do not lead to noticeable differences in P m . Alternatively, if η ≈ η crit , the plasma contribution to the dispersion cannot be treated as a perturbation and the calculation of P m is more sensitive to the ionization level and hence the ionization model. We see this reflected in our data for the lowest ionization potential gases (Xe for λ 0 = 1300 nm; Kr for λ 0 = 522 nm), where the choice of λ q was closer to λ min than for the other gas species, meaning η eff ≈ η crit in this case. This leads to a greater sensitivity of P m on η, with the result that the choice of ionization model has a noticeable impact on the calculated values of P m .

Conclusion
In summary, we have measured high harmonic spectra produced by two different fundamental wavelengths (522 and 1300 nm) and a range of noble gases (Xe, Kr, Ar and Ne). For each combination of fundamental wavelength and gas species investigated, the phasematching pressure and harmonic cut-off wavelength were measured. Using a 1-D phenomenological model, we calculated the high harmonic phasematching pressure and cut-off wavelength for two different ionization models: ADK and YI. The calculated phasematching pressure does not vary with choice of ionization model for harmonic wavelengths longer than the harmonic cut-off and reasonable agreement is observed between the experimental and calculated phase-matching pressures. In addition, we find better agreement between experimental and calculated harmonic cut-off wavelengths when the YI ionization model is used compared to the more commonly utilised ADK model. Our results show that even in situations where tunnel ionization is expected to dominate over multiphoton ionization (i.e. γ < 1), an ionization model which takes into account multiphoton contributions to the ionization rate is important for accurately predicting the experimental high harmonic cut-off wavelength. We anticipate that this result will further inform efforts to accurately model the wavelength-scaling behaviour of a range of strong-field phenomena.