Open Access
Add to BibSonomy
Abstract
We present a computational model for the nonlinear response of molecular oxygen exposed to strong mid-wavelength and long-wavelength infrared optical fields. Based on a non-Hermitian approach utilizing metastable electronic states, the nonlinear polarization and strong-field ionization are described as intimately connected properties. Good agreement with the measured nonlinear index and ionization rates is shown, and parameterized response functions are provided to facilitate applications in large-scale simulations of infrared optical pulses interacting with gaseous media.
© 2020 Optical Society of America under the terms of the OSA Open Access Publishing Agreement
1. Introduction
Nonlinear optics of gaseous media in the mid-wavelength and long-wavelength infrared regions has recently attracted considerable interest. Thanks to the richness of the chromatic dispersion landscapes at longer wavelengths, various highly nonlinear regimes can occur. Similar to the optical filamentation [1–3] studied in near-infrared and visible, many situations could be characterized as envelope-driven nonlinear optics. On the other hand, when the nonlinearity dominates, carrier-wave driven interactions can reshape optical waves on sub-cycle time scales [4], and harmonic generation affects the propagation dynamics of the pulse even at the fundamental frequency [5]. Numerical modeling has played an important role in the search for the best way to describe the rich physics including nonlinear polarization [6–22], below threshold harmonic generation [23–27], supercontinuum [28–35], plasma dynamics [36], and the interplay between ionization and nonlinear polarization [37].
In comparison to the near-infrared regime, numerical simulations at high intensities at longer wavelengths must capture what is a comparatively delicate dynamic relationship between various effects that control light-matter interactions, and this is especially true for long-distance propagation of optical pulses. For example, while the typical plasma densities are lower, their effects may accumulate over longer propagation paths, resulting in a need for ever more accurate descriptions of a number of concomitant effects that govern the high-intensity pulse dynamics. Strictly quantum-level description [38–40] of light-matter interactions is still computationally demanding, although significant progress is being achieved in this field [41]. That is why for large-scale simulations, such as intense light delivery to distant targets, one requires an approach that honors the microscopic origin of the interaction, but allows for computationally efficient model implementation.
The purpose of this work is to present a model for molecular oxygen, as one of the main components of air important for many applications of laser pulse propagation. Specifically, we describe its nonlinear response which includes both the nonlinear polarization and ionization as two intimately connected aspects in a unified framework. The model is based on the Metastable Electronic State Approach (MESA) [42] in which quantum mechanical calculations are used for one-time characterization of the given species. While the initial simulation of an atom or a molecule may require a lot of computations, subsequent response evaluation is inexpensive and allows one to apply the model in large-scale simulations. The method we use in this work has been tested against exact solutions [43], further verified with the results of TDSE simulations [44,45], and quantitatively compared with space-and-time resolved measurements [11,46] of the nonlinear response in noble gases [47]. Since a MESA model for molecular nitrogen was developed utilizing previous multielectron calculations [48], molecular oxygen is the most important missing piece needed for applications involving atmospheric constituents. It is our aim here to fill in this gap. With the help of the nonlinear response of O$_2$ presented in this work, large-scale numerical simulations of long-distance laser pulse propagation in atmosphere will be able to eliminate their most uncertain aspect, which is the relative strength between the de-focusing mediated by free-electrons and the nonlinear self-focusing.
The present result is the first non-Hermitian description of the nonlinear field-induced polarization of a molecule that is based on the single active electron application. As such, it provides a much needed test bed for future application to other small molecules. Indeed, methods have been demonstrated [49] that allow the data needed to construct MESA models for molecules to be obtained.
We shall start with a brief description of a simulation framework for intense optical pulses, and identify the required descriptors for the light-matter interactions. We continue with our numerical approach to obtain meta-stable states of molecular oxygen that serve as a basis of our model. We then describe our main results in Section 4., and compare them to experiments. Finally, we want to provide parameterizations for the nonlinear dipole moment and the ionization rate that are meant to facilitate applications in numerical simulations.
2. Numerical modeling of optical filamentation
The model described in this work can be applied in various approaches designed to simulate optical pulse propagation [50] as long as the representation utilizes the electric field amplitude and not the envelope of the pulse [1,2]. As an example, Forward Maxwell Equation [51] represents one of the suitable methods. Here we summarize the application in the framework of the unidirectional pulse propagation equation (UPPE) [52] where our O$_2$ model can be implemented as a user-defined add-on for the the gUPPEcore simulator [53]. For simplicity, we restrict this illustration to the scalar approximation, and note that the vector generalization is straightforward. The evolution equation for the temporal and spatial spectral amplitude of the electric field reads
In the framework of the metastable electronic-state approach [42,47], a molecule exposed to a strong electric field (we shall use notation $F\equiv E(t)$ to indicate atomic units) is described by two functions or descriptors. We restrict this work to wavelengths roughly longer than two microns, where post-adiabatic corrections [45] of the MESA models can be omitted and both descriptors respond to the instantaneous value of the field $F(t)$.
The first function to represent the molecule is the nonlinear dipole moment, $P_{\textrm {nl}}(F(t))$, that gives rise to the nonlinear polarization in (1)
where $e$ is the elementary charge, $a_{\textrm {B}}$ stands for the Bohr’s radius, and $N$ is the number density of molecules in the gas. Note that only the nonlinear (with respect to the field strength) part of the dipole moment is used, This is because the linear polarizability is included in the pulse propagation simulation through the tabulated function of $\epsilon (\omega )$, and the removal of the linear-in-field dipole is needed to avoid double-counting.The second function, used to calculate the ionization rate in a strong field, is the Stark-shifted, complex-valued energy of the unstable (i.e. resonant) ground state, $E_g(F)$, tabulated as a function of the electric field strength $F$. It enters the evolution equation for the free electron density $\rho (t)$, from which the current density is obtained by integrating in time,
3. Resonant electronic states
3.1 Single active electron model
Nonlinear effects in laser driven atoms and molecules are often studied within the single-active electron (SAE) approximation. For example, the phenomenon of suppressed molecular ionization in a strong laser field has been first theoretically explained for the O$_2$ [54] molecule and later extended to complex molecules, such as C$_{60}$ [55], using the strong-field approximation with SAE wavefunctions. Similarly, the nonlinear dipole response in the form of high harmonic generation and the application to imaging of diatomic molecules has been successfully studied using the SAE approximation, starting with the work in Ref. [56]. Thus, SAE calculations often provide accurate quantitative estimates, which can be used as parameters in other complex calculations, e.g. in propagation and filamentation codes where a coupled numerical solver of the multielectron time-dependent Schrödinger equation along with the Maxwell equations is unfeasible. Furthermore they can be scaled to parameter regimes (e.g., long wavelength regime), that are very challenging to explore numerically with multielectron codes. In the present case we test the applicability of the SAE approximation against experimental data for the ionization yields.
In the SAE approximation it is assumed that the laser field only interacts with one electron in the atom or molecule (here, an electron in the highest occupied molecular orbital, (HOMO)), while all other electrons are treated as frozen. We utilize an approach to construct a SAE potential based on an effective Kohn-Sham potential from density functional theory (DFT), and fit an analytic form to this numeric potential [57]. This method has been shown to provide a key element for benchmark tests between calculations based on the time-dependent Schrödinger equation and multielectron simulations for atoms [57], which shows the high accuracy of the present method to construct SAE potentials. Here, we apply the same approach to the oxygen molecule. This overcomes a deficiency in our previous application of the MESA model, in which we used a simplistic model potential and did not achieve the same level of agreement with experimental data as in the present case (see Section 4).
The effective Kohn-Sham potential $V_{eff,\sigma }(\textbf {r}) = V_{ext}(\textbf {r}) + V_H(\textbf {r}) + V_{xc,\sigma }(\textbf {r})$ is calculated using the OCTOPUS code [58–60] on a cylindrical 2D ($\rho ,z$) grid with spacing $0.05$ a.u. and length of 10 a.u. in both dimensions. The potential comprises Coulomb terms from the nuclei,
For the analytic form, we use a pair of radially-symmetric functions centered at each nucleus plus a cylindrically-symmetric term along the molecular axis. This approximately captures the details of the potential from the core and valence electrons while retaining a simple form. Thus, for a molecule aligned with the $o$-axis, we fit an analytic potential of the form
The radial functions $f(r)$ are a truncated form of atomic SAE potentials, which we developed earlier [57]. The nuclear-electron Coulomb potential is given by $-Z/r$ (here, $Z=8$ for oxygen) while the electron-electron Hartree potential is modeled by a Yukawa shielding term. Note that the net +1 charge seen by the active electron is split evenly across the two nuclei. This ensures that the long- and short-range behavior of the potential is physically accurate. Additionally, an exponential decay term captures the character of the potential at intermediate distances. The length scales of the Yukawa and exponential terms as well as the exponential coefficient are left as fitting parameters. The cylindrical Gaussian term $g(x, \rho )$ with three fitting parameters coarsely captures the effect of molecular bonding orbitals.
The numerical effective Kohn-Sham potential from the OCTOPUS code is fit to the analytic form (8) using least-squares fitting. Fitting of the radial term is facilitated by (Fourier) transforming the 2D potential to a radial form. This improves the speed of the fit, but is not strictly necessary. The accuracy of the fit can be seen from the comparison of the analytical SAE fit with the numeric DFT potential in Fig. 1. Shown are cuts along (a) and perpendicular to (b) the internuclear axis as well as the relative error as a function of the coordinates (c). Since the electron-nucleus potential $V_{ext}$ is the same in the two forms of the potential, we have excluded it from the visualization to enhance the comparison. We further note that the ionization potential of the analytic fit ($I_p = 0.43$ a.u.) agrees well with that of the OCTOPUS calculations ($I_P=0.40$ a.u.).
Fig. 1. Top row: Comparison of the analytical SAE potential fit (dashed red line) and the numeric DFT potential (solid black line) (a) along the molecular axis ($\rho =0$) and (b) for a line bisecting the molecular axis at the origin ($o=0$). Bottom row: (c) Relative error between numeric DFT potential and analytic SAE fit (Eq. (8)). Note that the error, within a range of $\approx 1.5$ a.u. from the nuclei, is $\leq 5\%$.
3.2 Wavefunctions
Having selected the SAE potential, we search for the relevant wavefunctions, which are the Stark resonances born from the highest-energy populated states. The eigenvalue problem can be written as
In order to enforce the outgoing boundary conditions, we have used the technique described in Ref. [45,47]. It is based on external complex scaling realized by a complexified spatial axis. As the grid-axis departs from the real axis, the outgoing waves appear to decay at large distances, and the wavefunction can be terminated by (approximate) zero at the domain edge. However, in this case the problem is three-dimensional, and presents a challenge to solve. In order to keep the underlying matrices as small as possible we have improved the previous approach by implementing transparent boundary conditions that allow us to utilize a smaller computational domain.
Assuming that the domain boundary $x=a$ is far enough from the origin so that we can neglect the molecular potential, the solutions for $x>a$ can be expressed in terms of Airy functions, and the boundary condition reads
The obvious complication with this transparent boundary is that it depends on the unknown eigenvalue $W$, and as a result we have a nonlinear eigenvalue problem to solve. Fortunately, it turns out that a direct iteration procedure converges very fast when $W(F)$ calculated previously is utilized as a guess for $W(F+\delta F)$: as few as two iterations are enough to obtain a $W$-value converged to four or five digits. We have tested this method on the case of atomic Hydrogen for which there are many accurate results available in the literature (see e.g. [63]). The method is applicable in other contexts, and its detailed description will be presented elsewhere. For the present purpose, the simplified version of the boundary condition (12) with $k_\perp \to 0$, is used to obtain the discretized Laplacian operator for the grid points at the boundary of the domain. For example, a three-point scheme reads
The advantage of this approach is that we can make the domain smaller. Specifically, when external complex scaling is used to enforce the outgoing boundary conditions as in Ref. [44,64,65], the domain must be large enough to provide ample space for the wavefunction to decay to (approximate) zero. This is not needed here; we can “intercept” the oscillating tail of the outgoing resonance, and bring the boundary condition as close to the origin as the molecular potential allows (for a short-range potential the domain edge can even be placed under the barrier in the classically forbidden region).
As yet another modification aimed at keeping our grids as coarse as possible, we have implemented the Laplacian with a thirteen-point scheme giving us a higher, fourth-order discretization accuracy. While the resulting matrices are obviously denser and cause significantly higher memory demands from the eigen-solver, the higher-order discretization helps to reduce the numerical dispersion which on a coarse grid affects the oscillatory wavefunction tail representing an escaping particle. We have performed our calculation on different grids in order to check that our results are well behaved. Moreover, we have calculated the wavefunction in two different ways, rotating either the molecule or the direction of the field. We show illustrative comparison of data obtained on different grids in the discussion section.
To obtain the resonant energies, we have utilized the FEAST solver [66]. Two states, which are degenerate at zero field $F=0$, contribute comparably to the nonlinear polarization at stronger fields. They can be distinguished by their symmetry w.r.t. the plane given by the molecular axis and the field, and we utilize this symmetry to calculate them separately. We solved the eigenvalue problem for a range of field strengths $F\le 0.08$ and for a discrete set of inclination angles, $\alpha =k \pi /12$, between the molecular axis and the direction of the field. The results for a randomly oriented molecular ensemble are obtained by averaging over the inclination angle $\alpha$ with weight $\sin (\theta ) d\theta$, which reflects the assumption that all molecule orientations in the three dimensional space are equally likely. The outcome of the averaging, $E_g(F) = \sum w_k E_g(F, \alpha = k \pi /12)$, is expressed in terms of weights $w_k$ which are calculated analytically for an interpolation in the form $E_g(F, \alpha ) = \sum c_n(F) \cos (2 n \alpha )$, which reflects the symmetry of the tabulated function. Explicit form of $w_k$ can be readily obtained with Mathematica.
The nonlinear dipole can be obtained from the calculated meta-stable wavefunction in two ways. One is the numerical integration of the expectation value and the other utilizes the derivative of the energy as a function of the field strength [48]. It is the latter that we have used for this work, because while the two methods are equivalent on a sufficiently large grid, the former is not suited for grids with transparent boundaries which we utilize here. On the other hand, the disadvantage of the energy-based approach to obtain the nonlinear dipole is that a differentiable representation of $E_g(F)$ must be created, usually in terms of a fitting function. In order to minimize the number of fitting functions, we perform the angle averaging, as described above, before extracting the nonlinear dipole moment. Thus, our results describe an ensemble of randomly oriented molecules, but not a molecule with a specific inclination angle.
It is known that intense laser fields impart partial orientation order on molecules. However, the degree of alignment is small for the range of laser intensities and pulse durations typical in optical filamentation. While even small alignment among the molecules manifests in the well-known orientational contribution to the nonlinear index (Raman effect) caused by the anisotropy of the linear molecular polarizability, its effect on the nonlinear dipole moment is small, at least given the accuracy of experimental values of $n_2$ and/or of our model.
4. Nonlinear polarization and ionization rate
Let us turn to our central result. Figure 2 summarizes the two MESA descriptors for an ensemble of randomly oriented oxygen molecules. The nonlinear dipole per molecule is shown in atomic units together with the ionization rate, both as functions of the external-field strength also given in atomic units. These two functions constitute the core of our description as outlined in Section 2.
Fig. 2. Nonlinear dipole moment and the ionization rate of the molecular oxygen calculated from two highest-energy metastable states. Full and dashed lines represent results obtained on different discretization grids (see text for details).
To verify that our results are independent of the parameters of the numerical grid used for the eigenvalue problem, we have obtained solutions on a number of different discretizations. In Fig. 2 we compare two representative sets. The grid labeled A is our reference, and has the resolution of $\Delta x = 0.06$ in the central region. With 120 discretization points in the $x$-dimension, the grid is inhomogeneous with grid spacing increasing toward the boundaries. Grid B is denser ($\Delta x = 0.0375$), and used 200 points along the direction of the external field and has a comparable spatial extent, with the outer $x$-boundary around 30 a.u. While small deviations between the results obtained from different discretizations are evident for strong fields, in the region relevant for optical filamentation ($F < 0.05$) the results are sufficiently close. From the practical standpoint, both MESA descriptors, i.e. the nonlinear dipole and the ionization rate, appear to be sufficiently accurate for applications in simulation of optical pulses typical for extreme nonlinear optics in mid- and long-wavelength infrared.
Of course, the question of how these results compare to various experimental characterizations of molecular oxygen is most important. Let us start with nonlinear polarization. To facilitate an intuitive comparison, Fig. 3 depicts the nonlinear dipole moment transformed into what can be called the effective nonlinear index. It is formed as
and its value can be interpreted as the effective nonlinear index experienced by a harmonic wave. The shape of $n_{2\textrm {eff}}$ as a function of the field strength reflects the deviations of the full nonlinear dipole from the cubic behavior expected for the pure instantaneous optical Kerr effect. Figure 3(a) shows this quantity transformed into usual units of cm$^{2}/$W under the assumption that the field is a continuous wave with a fixed cycle-averaged intensity. For comparison in Fig. 3(b) are data from Ref. [67] showing the results of the measurements of the nonlinear index for a number of mid- and long-infrared wavelengths. The figure shows that the effective nonlinearity compares quite well with the experiment. Because the experimental data were obtained for lower intensities, these measurements are most relevant for the weak-field region of the theoretical results. Here we see that the model gives somewhat weaker nonlinearity.Fig. 3. Effective nonlinear index as function of intensity (a) calculated from the nonlinear-dipole function. The gray area indicates the uncertainty due to different fitting functions used to represent the nonlinear dipole. Symbols (b) show the experimental results for the nonlinear index measured for a number of mid- and long-infrared wavelengths (reproduced from Ref. [67]). The vertical axes in the two panels are aligned to facilitate comparison of the measured and calculated value of the nonlinear index
The nonlinear dipole is obtained from the meta-stable energy averaged over the inclination angle and fitted as a function of the external field $F$, and this is subsequently differentiated with respect to the field strength. Also shown in Fig. 3 is the gray area that approximately represents the uncertainty of the calculated effective nonlinear index due to different choices of the fitting functions. Naturally, the deviations are most pronounced in the low-field region, and they indicate that if we had to characterize our MESA model in terms of the nonlinear index $n_2$, we should ascribe to it an uncertainty of about twenty percent. However, it must be emphasized that it is never the data shown in Fig. 3 that is used in simulations. Rather, it is the curves in Fig. 2 that represent the nonlinear response, and there the accuracy is better.
As for the shape of the effective nonlinear index function versus field, what we have here is quite typical for the adiabatic regime in noble gases, and it also is in line with our previous results obtained for molecular nitrogen [48]. It suggests that the nonlinearity is slightly higher for stronger fields. For very strong fields, the effective nonlinearity goes down but this is concomitant with the increase in ionization which results in intensity clamping so that this regime is unlikely to be observed. To the best of our knowledge there are no measurements available at the time of this writing that could be used to verify the shape of this curve or the increasing trend.
Yet another way to verify that the new O$_2$ model gives reasonable values of the nonlinear index is a comparison with Argon. Absolute measurements over a range of wavelengths show that the two species have rather similar nonlinear indices [67]. That motivates Fig. 4, which compares the nonlinear dipole curve of Argon obtained in Ref. [45] and verified against measurements in Ref. [47] with its counterpart for O$_2$ calculated in this work. It shows that up to typical maximal fields that occur in spontaneously created optical filaments, the two curves are actually quite similar. This is in line with the fact that the nonlinear indices are comparable, and it further supports the new model. Deviations between the species do occur at stronger fields, which we think reflects the fact that the critical field is lower for oxygen. However, for these strong fields it is usually the freed electron response that largely dominates the dynamics, and it is difficult to measure the precise shape of these nonlinear dipole curves.
Fig. 4. Nonlinear dipole of O$_2$ (line) compared to that of Argon (symbols) [45]. The curves are close to each other for weak to moderated fields, reflecting similar nonlinear index values in the two species.
Let us turn to the ionization rate. Unfortunately, unlike for noble gases [47], comparison with experiments is not yet possible on absolute scale. In order to assess the functional dependence on the field strength, Fig. 5 depicts our calculated ionization yield compared with two experimental data sets, Refs. [68,69]. To facilitate this comparison, we have digitized the data from the relevant figures. Since the measurement results are in arbitrary units, we have re-scaled them vertically to match as well as possible our theoretical curve. From the theoretical ionization rate, we have calculated the yield for the pulse durations corresponding to the experiments, and integrated the result over the focal volume, under an assumption of an ideal Gaussian beam. One can see that the overall shape of the ionization-yield function fits reasonably well. We consider the set L (by Lin at al.) being the most relevant as this was measured at $\lambda =2000$ nm, while the other set is for shorter wavelengths, where post-adiabatic corrections [45] may be needed to improve the present adiabatic MESA-based theory. Such corrections could be also relevant for other data sets at 800 nm with a relatively long 220 fs duration [70].
Fig. 5. Ionization yield calculated for the MESA-model of O$_2$ (solid blue line) compared with measurements taken from two representative works, L: [69] (Lin et al.) and G: [68] (Guo et al.). Experimental data, originally in arbitrary units, were digitized from the publication graphs and scaled.
To get at least a rough idea about the accuracy of our ionization rate, we turn to the ratio between O$_2$ and its “companion” noble gas Xe. This ratio has been discussed in the literature in the context of observed ionization suppression in a molecular gas with respect to an atomic species with comparable ionization potential [54,68–71]. Here we utilize our previous calculations for Xenon to evaluate the ratio of the ionization yields. Figure 6 shows our result compared to the data obtained from Refs. [68,69]. We have digitized the measurements of the ionization yield signals (integrated over focal volume) for Xe and O$_2$ from each work, and calculated the ratio as a function of the intensity from the linear interpolation of the digitized data; the results are shown in the figure. Note that the noise in the data is partially due to our digitization and subsequent interpolation. However, the noise is not especially important for our current purpose. Given the significant spread between the experimental results, one can say that our calculated value might be similar to that seen in experiments in the range of low intensities. However, while the measurements do not precisely agree on the particular value of the ratio, they exhibit an increase with the intensity. The data shown in Fig. 6 is why we think that our model underestimates the actual ionization rate of molecular oxygen, perhaps by a factor between two to four. While this may seem like a lot, this accuracy is in fact quite sufficient for applications in optical filamentation; the highly nonlinear nature of the ionization-rate curve ensures that the resulting error is a small adjustment of the clamped intensity within a high-intensity filament.
Fig. 6. Ratio of the ionization yield between O$_2$ and Xe. The solid and dashed lines represent our calculation, while symbols are obtained from Refs. G: [68] and L: [69], as described in the text.
For the sake of practical implementation, we have designed fitting functions that can represent our results with an accuracy sufficient for the simulations typical in optical filamentation, namely in the range of the field strength of up to 0.08 atomic units. This corresponds to the range of light intensities of 225 TW/cm$^{2}$ which is adequate for the regime encountered in spontaneously formed optical filaments, both in the near infrared and even more safely in the mid infrared where typical peak intensities tend to be smaller. For the ionization rate we found
to represent the data adequately. Note that while the above functional form could be motivated by the expected functional behavior of the tunneling-regime ionization rate formulas, the coefficient values could change when terms linear and higher are included in the exponential, and the values are in general dependent on the interval fitted. For the nonlinear dipole, one can use this polynomial parameterization5. Conclusion
The main result of this work is the model for the nonlinear polarization and ionization rate of the molecular oxygen. It is based on the previously described metastable electronic state approach, and utilizes the single active electron approximation to characterize the response of the molecule. Together with our previous report on the MESA model for the molecular nitrogen [48], and our earlier results concerning noble gases [47], this completes the set of most important species for extreme nonlinear optics in the atmosphere and other gaseous media. The envisioned applications are mainly in the realm of mid- and long-wavelength infrared, which is the field that has attracted a lot of recent interest.
The metastable electronic state approach gives the functions for the nonlinear polarization and the ionization rate that are intimately connected; the two descriptors are derived from a single metastable state of a system, and this is important for practical applications. Namely, it eliminates the often-perceived freedom to separately adjust the ionization parameters and self-focusing parameters (nonlinear index in the standard model of optical filamentation) to a given experiment.
Of course it would be folly to expect that our parameterized model is accurate as is and needs no adjustment. Should absolute, accurate and simultaneous measurements of the ionization rate and nonlinear polarization as functions of the field strength become available, this model can be further calibrated for accuracy. Such a task could be achieved by adjusting two scaling parameters while preserving the functional shape(s) of the nonlinear response, following the scaling idea outlined in Ref. [72]. From this standpoint, the present work can serve as a basis for a future improved description.
Another possible improvement of our model would be its extension for shorter wavelengths. Inclusion of post-adiabatic corrections as described in [45] requires higher accuracy, for both the eigen-energies and wavefunctions, than that achieved in this work, but we trust it is not an unrealistic goal.
Finally, an interesting and important basic question that this work raises is related to our use of the single active electron approximation. As it is not a priori given or self-evident which quantities can be accurately described in this approximation, we must admit the possibility that the quantities we calculate could be biased due to the single active electron assumption. Thus, our model can also serve as a useful test bed that can help to shed additional light on this issue.
Last but not least, this work is intended to motivate future application to other small molecules. We would like to emphasize that methods for small molecules, e.g. carbon dioxide, do exist and have been in use that can produce the necessary data for construction of MESA-based models. All that is needed is to apply some version of absorbing boundaries or external complex scaling and collect the information regarding both the real and complex parts of the resonance energy (often only the imaginary part is reported). Alternatively, the dipole moment can be utilized as well. The challenge in such applications consist in the need to obtain a sufficiently dense tabulation as a function of the applied field, which can be averaged over molecular orientations. Needless to say, these are merely technical issues, and it is our hope that this paper will motivate further related work on the nonlinear response of molecules.
Funding
Air Force Office of Scientific Research (FA9550-16-1-0121, FA9550-18-1-0183).
Disclosures
The authors declare no conflicts of interest.
References
1. A. Couairon and A. Mysyrowicz, “Femtosecond filamentation in transparent media,” Phys. Rep. 441(2-4), 47–189 (2007). [CrossRef]
2. L. Bergé, S. Skupin, R. Nuter, J. Kasparian, and J.-P. Wolf, “Ultrashort filaments of light in weakly ionized optically transparent media,” Rep. Prog. Phys. 70(10), 1633–1713 (2007). [CrossRef]
3. S. L. Chin, Femtosecond Laser Filamentation (Springer, 2009).
4. P. Panagiotopoulos, P. Whalen, M. Kolesik, and J. V. Moloney, “Carrier field shock formation of long-wavelength femtosecond pulses in single-crystal diamond and air,” J. Opt. Soc. Am. B 32(8), 1718–1730 (2015). [CrossRef]
5. N. Aközbek, A. Iwasaki, A. Becker, M. Scalora, S. L. Chin, and C. M. Bowden, “Third-harmonic generation and self-channeling in air using high-power femtosecond laser pulses,” Phys. Rev. Lett. 89(14), 143901 (2002). [CrossRef]
6. V. Loriot, E. Hertz, O. Faucher, and B. Lavorel, “Measurement of high order Kerr refractive index of major air components,” Opt. Express 17(16), 13429–13434 (2009). [CrossRef]
7. V. Loriot, E. Hertz, O. Faucher, and B. Lavorel, “Measurement of high order Kerr refractive index of major air components: erratum,” Opt. Express 18(3), 3011–3012 (2010). [CrossRef]
8. P. Béjot, J. Kasparian, S. Henin, V. Loriot, T. Vieillard, E. Hertz, O. Faucher, B. Lavorel, and J.-P. Wolf, “Higher-order Kerr terms allow ionization-free filamentation in gases,” Phys. Rev. Lett. 104(10), 103903 (2010). [CrossRef]
9. P. Béjot, E. Hertz, B. Lavorel, J. Kasparian, J.-P. Wolf, and O. Faucher, “From higher-order Kerr nonlinearities to quantitative modeling of third and fifth harmonic generation in argon,” Opt. Lett. 36(6), 828–830 (2011). [CrossRef]
10. Z. Wang, C. Zhang, J. Liu, R. Li, and Z. Xu, “Femtosecond filamentation in argon and higher-order nonlinearities,” Opt. Lett. 36(12), 2336–2338 (2011). [CrossRef]
11. J. K. Wahlstrand, Y.-H. Cheng, and H. M. Milchberg, “High field optical nonlinearity and the Kramers-Kronig relations,” Phys. Rev. Lett. 109(11), 113904 (2012). [CrossRef]
12. A. Spott, A. Jaroň-Becker, and A. Becker, “Ab initio and perturbative calculations of the electric susceptibility of atomic hydrogen,” Phys. Rev. A 90(1), 013426 (2014). [CrossRef]
13. J. K. Wahlstrand, S. Zahedpour, A. Bahl, M. Kolesik, and H. M. Milchberg, “Bound-electron nonlinearity beyond the ionization threshold,” Phys. Rev. Lett. 120(18), 183901 (2018). [CrossRef]
14. M. Nesrallah, A. Hakami, G. Bart, C. R. McDonald, C. Varin, and T. Brabec, “Measuring the Kerr nonlinearity via seeded Kerr instability amplification: conceptual analysis,” Opt. Express 26(6), 7646–7654 (2018). [CrossRef]
15. M. Hofmann and C. Brée, “Adiabatic Floquet model for the optical response in femtosecond filaments,” J. Phys. B: At., Mol. Opt. Phys. 49(20), 205004 (2016). [CrossRef]
16. L. Wang and W. Lin, “The impact of the varied nonlinear refractive index of higher-order Kerr effect on the laser pulse’s propagation,” Optik (Munich, Ger.) 126(24), 5387–5391 (2015). [CrossRef]
17. L. Su-Yu, G. Fu-Ming, Y. Yu-Jun, and J. Ming-Xing, “Defocusing role in femtosecond filamentation: Higher-order Kerr effect or plasma effect?” Chin. Phys. B 24(11), 114207 (2015). [CrossRef]
18. M. Hofmann and C. Brée, “Femtosecond filamentation by intensity clamping at a Freeman resonance,” Phys. Rev. A 92(1), 013813 (2015). [CrossRef]
19. M. Tarazkar, D. A. Romanov, and R. J. Levis, “High-order nonlinear refractive indices for He, Ne, Kr, and Xe atoms,” Phys. Rev. A 90(6), 062514 (2014). [CrossRef]
20. M. Tarazkar, D. A. Romanov, and R. J. Levis, “Higher-order nonlinearity of refractive index: The case of Argon,” J. Chem. Phys. 140(21), 214316 (2014). [CrossRef]
21. D. Wang and Y. Leng, “Measuring high-order Kerr effects of noble gases based on spectral analysis,” Opt. Commun. 328, 41–48 (2014). [CrossRef]
22. T. C. Rensink, T. M. Antonsen, J. P. Palastro, and D. F. Gordon, “Model for atomic dielectric response in strong, time-dependent laser fields,” Phys. Rev. A 89(3), 033418 (2014). [CrossRef]
23. G. O. Ariunbold, P. Polynkin, and J. V. Moloney, “Third and fifth harmonic generation by tightly focused femtosecond pulses at 2.2 μm wavelength in air,” Opt. Express 20(2), 1662–1667 (2012). [CrossRef]
24. A. Nath, J. A. Dharmadhikari, A. K. Dharmadhikari, and D. Mathur, “Seventh-harmonic generation from tightly focused 2 um ultrashort pulses in air,” Opt. Lett. 38(14), 2560–2562 (2013). [CrossRef]
25. D. L. Weerawarne, X. Gao, A. L. Gaeta, and B. Shim, “Higher-order nonlinearities revisited and their effect on harmonic generation,” Phys. Rev. Lett. 114(9), 093901 (2015). [CrossRef]
26. A. Spott, A. Becker, and A. Jaroň-Becker, “Transition from perturbative to nonperturbative interaction in low-order-harmonic generation,” Phys. Rev. A 91(2), 023402 (2015). [CrossRef]
27. N. Garejev, V. Jukna, G. Tamošauskas, M. Veličkė, R. Šuminas, A. Couairon, and A. Dubietis, “Odd harmonics-enhanced supercontinuum in bulk solid-state dielectric medium,” Opt. Express 24(15), 17060–17068 (2016). [CrossRef]
28. M. A. Porras, A. Dubietis, A. Matijošius, R. Piskarskas, F. Bragheri, A. Averchi, and P. D. Trapani, “Characterization of conical emission of light filaments in media with anomalous dispersion,” J. Opt. Soc. Am. B 24(3), 581–584 (2007). [CrossRef]
29. C. Brèe, A. Demircan, and G. Steinmeyer, “Saturation of the all-optical Kerr effect,” Phys. Rev. Lett. 106(18), 183902 (2011). [CrossRef]
30. A. A. Lanin, A. A. Voronin, E. A. Stepanov, A. B. Fedotov, and A. M. Zheltikov, “Multioctave, sub-two-cycle supercontinua from self-compressing, self-focusing soliton transients in a solid,” Opt. Lett. 40(6), 974–977 (2015). [CrossRef]
31. H. Liang, P. Krogen, R. Grynko, O. Novak, C.-L. Chang, G. J. Stein, D. Weerawarne, B. Shim, F. X. Kärtner, and K.-H. Hong, “Three-octave-spanning supercontinuum generation and sub-two-cycle self-compression of mid-infrared filaments in dielectrics,” Opt. Lett. 40(6), 1069–1072 (2015). [CrossRef]
32. A. V. Mitrofanov, A. A. Voronin, S. I. Mitryukovskiy, D. A. Sidorov-Biryukov, A. Pugžlys, G. Andriukaitis, T. Flöry, E. A. Stepanov, A. B. Fedotov, A. Baltuška, and A. M. Zheltikov, “Mid-infrared-to-mid-ultraviolet supercontinuum enhanced by third-to-fifteenth odd harmonics,” Opt. Lett. 40(9), 2068–2071 (2015). [CrossRef]
33. A. V. Mitrofanov, A. A. Voronin, D. A. Sidorov-Biryukov, S. I. Mitryukovsky, M. V. Rozhko, A. Pugžlys, A. B. Fedotov, V. Y. Panchenko, A. Baltuška, and A. M. Zheltikov, “Angle-resolved multioctave supercontinua from mid-infrared laser filaments,” Opt. Lett. 41(15), 3479–3482 (2016). [CrossRef]
34. R. Šuminas, G. Tamošauskas, G. Valiulis, V. Jukna, A. Couairon, and A. Dubietis, “Multi-octave spanning nonlinear interactions induced by femtosecond filamentation in polycrystalline ZnSe,” Appl. Phys. Lett. 110(24), 241106 (2017). [CrossRef]
35. R. Šuminas, G. Tamošauskas, V. Jukna, A. Couairon, and A. Dubietis, “Second-order cascading-assisted filamentation and controllable supercontinuum generation in birefringent crystals,” Opt. Express 25(6), 6746–6756 (2017). [CrossRef]
36. E. M. Wright, S. W. Koch, M. Kolesik, and J. V. Moloney, “Memory effects in the long-wave infrared avalanche ionization of gases: a review of recent progress,” Rep. Prog. Phys. 82(6), 064401 (2019). [CrossRef]
37. M. Richter, S. Patchkovskii, F. Morales, O. Smirnova, and M. Ivanov, “The role of the Kramers-Henneberger atom in the higher-order Kerr effect,” New J. Phys. 15(8), 083012 (2013). [CrossRef]
38. E. Lorin, S. Chelkowski, and A. Bandrauk, “A numerical Maxwell-Schrödinger model for intense laser-mater interaction and propagation,” Comput. Phys. Commun. 177(12), 908–932 (2007). [CrossRef]
39. E. Lorin, S. Chelkowski, E. Zaoui, and A. Bandrauk, “Maxwell — Schrödinger — Plasma (MASP) model for laser — molecule interactions: Towards an understanding of filamentation with intense ultrashort pulses,” Phys. D (Amsterdam, Neth.) 241(12), 1059–1071 (2012). [CrossRef]
40. E. Lorin, M. Lytova, A. Memarian, and A. D. Bandrauk, “Development of nonperturbative nonlinear optics models including effects of high order nonlinearities and of free electron plasma: Maxwell — Schrödinger equations coupled with evolution equations for polarization effects, and the SFA-like nonlinear optics model,” J. Phys. A: Math. Theor. 48(10), 105201 (2015). [CrossRef]
41. N. Berti, P. Béjot, E. Cormier, J. Kasparian, O. Faucher, and J.-P. Wolf, “Ab initio calculations of laser-atom interactions revealing harmonics feedback during macroscopic propagation,” Phys. Rev. A 99(6), 061401 (2019). [CrossRef]
42. M. Kolesik, J. M. Brown, A. Teleki, P. Jakobsen, J. V. Moloney, and E. M. Wright, “Metastable electronic states and nonlinear response for high-intensity optical pulses,” Optica 1(5), 323–331 (2014). [CrossRef]
43. J. M. Brown, A. Lotti, A. Teleki, and M. Kolesik, “Exactly solvable model for non-linear light-matter interaction in an arbitrary time-dependent field,” Phys. Rev. A 84(6), 063424 (2011). [CrossRef]
44. A. Bahl, J. M. Brown, E. M. Wright, and M. Kolesik, “Assessment of the metastable electronic state approach as a microscopically self-consistent description for the nonlinear response of atoms,” Opt. Lett. 40(21), 4987–4990 (2015). [CrossRef]
45. A. Bahl, E. M. Wright, and M. Kolesik, “Nonlinear optical response of noble gases via the metastable electronic state approach,” Phys. Rev. A 94(2), 023850 (2016). [CrossRef]
46. S. Zahedpour, J. K. Wahlstrand, and H. M. Milchberg, “Measurement of the nonlinear refractive index of air constituents at mid-infrared wavelengths,” Opt. Lett. 40(24), 5794–5797 (2015). [CrossRef]
47. A. Bahl, J. K. Wahlstrand, S. Zahedpour, H. M. Milchberg, and M. Kolesik, “Nonlinear optical polarization response and plasma generation in noble gases: Comparison of metastable-electronic-state-approach models to experiments,” Phys. Rev. A 96(4), 043867 (2017). [CrossRef]
48. A. Bahl, V. P. Majety, A. Scrinzi, and M. Kolesik, “Nonlinear optical response in molecular nitrogen: from ab-initio calculations to optical pulse simulations,” Opt. Lett. 42(12), 2295–2298 (2017). [CrossRef]
49. A. H. Larsen, U. De Giovannini, D. L. Whitenack, A. Wasserman, and A. Rubio, “Stark Ionization of Atoms and Molecules within Density Functional Resonance Theory,” J. Phys. Chem. Lett. 4(16), 2734–2738 (2013). [CrossRef]
50. A. Couairon, E. Brambilla, T. Corti, D. Majus, O. d. J. Ramírez-Góngora, and M. Kolesik, “Practitioner‘s guide to laser pulse propagation models and simulation,” Eur. Phys. J.: Spec. Top. 199(1), 5–76 (2011). [CrossRef]
51. A. V. Husakou and J. Herrmann, “Supercontinuum generation of higher-order solitons by fission in photonic crystal fibers,” Phys. Rev. Lett. 87(20), 203901 (2001). [CrossRef]
52. M. Kolesik and J. V. Moloney, “Nonlinear optical pulse propagation simulation: From Maxwell’s to unidirectional equations,” Phys. Rev. E 70(3), 036604 (2004). [CrossRef]
53. M. Kolesik, http://acms.arizona.edu/FemtoTheory/MK_personal/guppelab/.
54. J. Muth-Böhm, A. Becker, and F. H. M. Faisal, “Suppressed molecular ionization for a class of diatomics in intense femtosecond laser fields,” Phys. Rev. Lett. 85(11), 2280–2283 (2000). [CrossRef]
55. A. Jaron-Becker, A. Becker, and F. H. M. Faisal, “Saturated ionization of fullerenes in intense laser fields,” Phys. Rev. Lett. 96(14), 143006 (2006). [CrossRef]
56. M. Lein, N. Hay, R. Velotta, J. P. Marangos, and P. L. Knight, “Interference effects in high-order harmonic generation with molecules,” Phys. Rev. A 66(2), 023805 (2002). [CrossRef]
57. R. Reiff, T. Joyce, A. Jaroń-Becker, and A. Becker, “Single-active electron calculations of high-order harmonic generation from valence shells in atoms for quantitative comparison with TDDFT calculations,” J. Phys. Commun. 4(6), 065011 (2020). [CrossRef]
58. M. A. Marques, A. Castro, G. F. Bertsch, and A. Rubio, “Octopus: a first-principles tool for excited electron–ion dynamics,” Comput. Phys. Commun. 151(1), 60–78 (2003). [CrossRef]
59. A. Castro, H. Appel, M. Oliveira, C. A. Rozzi, X. Andrade, F. Lorenzen, M. A. L. Marques, E. K. U. Gross, and A. Rubio, “Octopus: a tool for the application of time-dependent density functional theory,” Phys. Status Solidi B 243(11), 2465–2488 (2006). [CrossRef]
60. X. Andrade, J. Alberdi-Rodriguez, D. A. Strubbe, M. J. T. Oliveira, F. Nogueira, A. Castro, J. Muguerza, A. Arruabarrena, S. G. Louie, A. Aspuru-Guzik, A. Rubio, and M. A. L. Marques, “Time-dependent density-functional theory in massively parallel computer architectures: the octopus project,” J. Phys.: Condens. Matter 24(23), 233202 (2012). [CrossRef]
61. J. P. Perdew and A. Zunger, “Self-interaction correction to density-functional approximations for many-electron systems,” Phys. Rev. B 23(10), 5048–5079 (1981). [CrossRef]
62. J. B. Krieger, Y. Li, and G. J. Iafrate, “Exact relations in the optimized effective potential method employing an arbitrary Exc[}ψiσ{],” Phys. Lett. A 148(8-9), 470–474 (1990). [CrossRef]
63. F. M. Fernández and J. Garcia, “Highly accurate calculation of the resonances in the Stark effect in hydrogen,” Appl. Math. Comput. 317, 101–108 (2018). [CrossRef]
64. A. Bahl, E. M. Wright, and M. Kolesik, “Nonlinear optical response of noble gases via the metastable electronic state approach,” Phys. Rev. A 94(2), 023850 (2016). [CrossRef]
65. A. Bahl, A. Teleki, P. K. Jakobsen, E. M. Wright, and M. Kolesik, “Reflectionless beam propagation on a piecewise linear complex domain,” J. Lightwave Technol. 32(22), 4272–4278 (2014). [CrossRef]
66. E. Polizzi, “Density-matrix-based algorithm for solving eigenvalue problems,” Phys. Rev. B 79(11), 115112 (2009). [CrossRef]
67. S. Zahedpour, S. W. Hancock, and H. M. Milchberg, “Ultrashort infrared 2.5-11μm pulses: spatiotemporal profiles and absolute nonlinear response of air constituents,” Opt. Lett. 44(4), 843–846 (2019). [CrossRef]
68. C. Guo, M. Li, J. P. Nibarger, and G. N. Gibson, “Single and double ionization of diatomic molecules in strong laser fields,” Phys. Rev. A 58(6), R4271–R4274 (1998). [CrossRef]
69. Z. Lin, X. Jia, C. Wang, Z. Hu, H. Kang, W. Quan, X. Lai, X. Liu, J. Chen, B. Zeng, W. Chu, J. Yao, Y. Cheng, and Z. Xu, “Ionization suppression of diatomic molecules in an intense midinfrared laser field,” Phys. Rev. Lett. 108(22), 223001 (2012). [CrossRef]
70. A. Talebpour, C.-Y. Chien, and S. L. Chin, “The effects of dissociative recombination in multiphoton ionization of O2,” J. Phys. B: At., Mol. Opt. Phys. 29(18), L677–L680 (1996). [CrossRef]
71. X. M. Tong, Z. X. Zhao, and C. D. Lin, “Theory of molecular tunneling ionization,” Phys. Rev. A 66(3), 033402 (2002). [CrossRef]
72. M. Kolesik and E. M. Wright, “Universal long-wavelength nonlinear optical response of noble gases,” Opt. Express 27(18), 25445–25456 (2019). [CrossRef]
