Tracking quantum clouds expansion in tunneling ionization

The notion of an electron trajectory proved itself extremely useful for the qualitative and, in many cases, quantitative description of various ionization phenomena. Even the simplest picture of classical electron motion in the field of an electromagnetic wave, the well-known simple man model (SMM) and its predecessors Brunel (1987, 1990); Corkum et al. (1989); Ehlotzky (1992); Krause et al. (1992); L’Huillier et al. (1993); Lewenstein et al. (1994); Corkum (1993); Lewenstein et al. (1995); Krausz and Ivanov (2009); Shvetsov-Shilovski et al. (2012); Arbo et al. (2015), provides a basis for understanding many ionization phenomena, such as above-threshold ionization (ATI) and high harmonic generation (HHG). The spectacular predictive power of the SMM gave rise to a variety of techniques in which the ionization is described quantum-mechanically and the subsequent electron motion is treated classically or semi-classically (classical trajectory Monte Carlo or CTMC method) Shvetsov-Shilovski (2021); Hofmann, Cornelia et al. (2021); Shvetsov-Shilovski et al. (2012); Pfeiffer et al. (2012); Arbo et al. (2015); Dimitrovski and Madsen (2015); Landsman and Keller (2015); Hofmann et al. (2013). In these approaches the role of quantum-mechanics (QM) consists in setting up initial conditions for the subsequent classical or semi-classical electron motion, and to assigning statistical weights to the trajectories originating at different times by employing the notion of the instantaneous ionization rate (IIR). Various analytical expressions for the IIR obtained in the framework of the quantum-mechanical strong field approximation (SFA) approach and its subsequent developments Keldysh (1965); Faisal (1973); Reiss (1980); Perelomov et al. (1966); Popov (2004); Ammosov et al. (1986); Ivanov et al. (2018) can be used for this purpose, such as the Ammosov-Delone-Krainov (ADK) Ammosov et al. (1986); Delone and Krainov (1991), or the Yudin-Ivanov Yudin and Ivanov (2001) formulas.

The initial values of the coordinates for the classical trajectories are determined by the tunnel exit position, which can be found using the field direction model (FDM) Landsman and Keller (2015) or a more refined approach based on the use of the parabolic coordinate system Hu et al. (1997) for the atomic systems governed purely by the Coulomb interaction. The initial velocities in the directions perpendicular to the orientation of the electric field vector at the time of ionization are typically assumed to be distributed according to the well-known SFA formula for the transverse velocity distribution Popov (2004). The initial velocity in the direction parallel to the electric field vector at the moment of ionization is typically assumed to be zero Hofmann et al. (2013).

The success and utility of this semi-classical picture of atomic and molecular ionization is ultimately due to the essentially semiclassical nature of many aspects of the ionization phenomena which can be explained quite satisfactorily using semiclassical trajectory simulations. The so-called low-energy structures in strong field ionization spectra W.Becker et al. (2014); Liu and Hatsagortsyan (2010), Coulomb focusing effect Brabec et al. (1996), nonadiabatic effects in strong field ionization Hofmann et al. (2014), frustrated tunneling ionization Ortmann et al. (2021), have been studied using the trajectory-based methods. The semiclassical approaches based on the classical trajectories can include the interference effects as well, which can be done by using the so-called Quantum Trajectory Monte Carlo (QTMC) approachLi et al. (2014); Fang et al. (2019), or the semiclassical two-step model for strong-field ionization Shvetsov-Shilovski et al. (2016), which allows to obtain angular photo-electron distributions in good agreement with fully quantum calculations based on the solution of the time-dependent Schrödinger equation (TDSE). Such semiclassical simulations usually require much less computational effort than the fully quantum calculations and consequently, for complex targets, when numerical solution of the TDSE becomes unfeasible, use of such methods may provide the only possibility to obtain quantitative description of the ionization process.

If we are interested in a purely quantum mechanical description of the motion of ionized electrons and still want to be able to use some classical notions, one can apply the saddle-point method (SPM) to evaluate the SFA Lewenstein et al. (1994); Popruzhenko (2014a); Becker et al. (2002); Popov (2005); Popruzhenko (2014b); Lai et al. (2015, 2017) or the Feynman’s path-integral expressions for the ionization amplitude Salières et al. (2001). One obtains in this way a description of the ionization phenomena in terms of the so-called quantum trajectories (QT). QT are determined by the SPM equations and make the action in the integrals determining ionization amplitudes stationary. QT is a generally complex electron trajectory originating at the complex saddle point $t_{s}$ and propagating till the final moment of time $t_{f}$, when the electron arrives at the detector. This approach gives a very transparent and appealing view of the ionization process. It is, moreover, very flexible and allows to design a number of different generalizations and developments de Morisson Faria and Maxwell (2020); Popruzhenko (2014a). One may, for instance, start with the SFA ionization amplitude, ignoring effects of the ionic potential in the continuum Lohr et al. (1997); Kopold et al. (2000); Popruzhenko (2014a), evaluate it applying the SPM Kopold et al. (2000); Salières et al. (2001); Korneev et al. (2012); Popruzhenko (2014a) and obtain a description of the ionization process in terms of the complex QT which are solutions to Newton’s equations for a classical electron in the presence of the laser field. Such a description provides a link with the SMM. To include the effects of the ionic potential one can consider perturbatively the Coulomb effects on the QT used in the SFA (the so-called Coulomb corrected SFA or CCSFA method Popruzhenko and Bauer (2008)). Alternatively, one can consider the Coulomb and laser field effects on the trajectories on equal footing and find the QT as solutions to the Newton’s equations of motion in presence of the Coulomb and laser fields, still using the SFA equation defining the saddle-point (the Trajectory-based Coulomb SFA or (TCSFA) method Yan et al. (2010)). Yet more generally, one may apply the SPM to evaluate the Feynman’s path-integral representation of the ionization amplitude Salières et al. (2001), obtaining QT which are solutions to the Newton’s equations of motion in presence of the Coulomb and laser fields with a more complicated condition defining the starting time $t_{s}$ of the trajectory.

The path connecting $t_{s}$ and $t_{f}$ in the complex time-plane is often chosen to consist of two straight line segments: $(t_{s},Re\ t_{s})$ and $(Re\ t_{s},t_{f})$, with $Re\ t_{s}$ interpreted as the tunnel exit point. One should bear in mind, however, that QT are all but a convenient (albeit very useful and powerful) mathematical construct, arising as a result of the application of the SPM for evaluation of ionization amplitudes. In particular, the path connecting $t_{s}$ and $t_{f}$ described above is not unique and can be deformed to cross the real time-axis almost at any given point Popruzhenko (2014a). This remark applies equally, of course, to the notion of the tunnel exit used in the CTMC method Popruzhenko (2014a); de Morisson Faria and Maxwell (2020).

This path-dependence of the time and location of the tunnel exit does not affect the quantum amplitude since the integrals defining the amplitudes depend only on the end points $t_{s}$ and $t_{f}$ of the path (as long as deforming the path we do not cross singular points of the integrand in the complex time-plane). As it was mentioned in the review work Popruzhenko (2014a), no physical experiment can favor particular values for time or location of the tunnel exit event, which does not prevent these notions to be extremely useful for practical purposes. In practice, the path connecting $t_{s}$ and $t_{f}$, consisting of the two straight line segments: $(t_{s},Re\ t_{s})$ and $(Re\ t_{s},t_{f})$ is the most convenient choice. The exit time $Re\ t_{s}$ and the tunnel exit point are related to the corresponding sub-barrier part of the QT as the real part of the expression $\displaystyle{\bm{r}}=\int\limits_{t_{s}}^{Re\ t_{s}}{\bm{v}}(t)\ dt$, where ${\bm{v}}(t)$ is the complex-valued sub-barrier electron velocity. The sub-barrier part of the QT thus defined can be used to provide the necessary prerequisites for the CTMC simulations, such as position of the tunnel exit, the transverse and longitudinal momentum distributions at the tunnel exit and the instantaneous ionization rate Li et al. (2014, 2016). Taking the real part of the QT at $t=Re\ t_{s}$ as the tunnel exit position has the advantage that the imaginary part of the action which determines the ionization probability is accumulated in the sub-barrier region, while subsequent propagation of the QT in the real time only produces phase-shift due to the change of the real part of the action. This choice also allows to avoid complicated issues of branching points and branch cuts in the complex time-plane de Morisson Faria and Maxwell (2020).

Of course, all the information we can hope to obtain about a physical system is encoded in its wave function. Any question about the motion of the ionized wave-packets for times within the laser pulse, should therefore be resolved from the solution of the TDSE, provided that this question has a physical meaning at all. The approach based solely on the information obtained by solving the TDSE, however, encounters the problem of identifying the contributions from different channels (ionization, excitation, and so on) for times within the laser pulse, when the wave-packets corresponding to these channels are not spatially separated yet. The splitting of the total wave-function of the system into the bound and ionized components seems to have been achieved in the SFA and the Perelomov-Popov-Terent’ev (PPT) approaches Keldysh (1965); Faisal (1973); Reiss (1980); Perelomov et al. (1966); Popov (2004); Delone and Krainov (1991). Such a splitting is not unique, however, and it is different in the two theories. Moreover, it is not gauge-invariant in the SFA or PPT approaches Popruzhenko (2014a). Similarly, the procedure based on projecting out contributions of the bound states of the field-free atomic Hamiltonian from the TDSE wave-function, which is sometimes used to define wave-packet describing ionized electron, is not gauge-invariant when applied for the times within the laser pulse duration.

In Tulsky and Bauer (2020) a method, allowing to identify the part of the wave-function describing ionized electron and relating the TDSE-based approach with the insight offered by the trajectory based approaches, has been proposed. In the framework of this method the photoionized part of the wave function is singled out by means of applying the time dependent surface flux (tSURFF) method Tao and Scrinzi (2012), which relies on the knowledge of the wave-function in the asymptotic region, when the photo-ionized part of the wave function is localized in space. By applying a short-time filter to the ionization amplitude, calculated using the tSURFF method, authors were able to identify the dominant pathways which form the photoelectron spectra.

Another group of methods allowing to connect the TDSE and the notion of trajectory is based on the Bohmian interpretation of the QM Bohm (1952). Bohmian view of the QM introduces a well-defined notion of the electron trajectory, exactly reproducing at the same time all the predictions of the conventional QM Bell (1987). This possibility of reintroducing trajectories in the QM framework has been exploited to describe ionization of atoms Botheron and Pons (2010); Jooya et al. (2015a) and molecules Sawada et al. (2014) driven by strong laser fields, and for the description of the HHG process Wu et al. (2013a, b); Jooya et al. (2015b). An approach to the problem of the tunneling time, based on the Bohmian QM, has been described in Zimmermann et al. (2016). In I.A.Ivanov et al. (2017) the notion of the coordinate distributions describing ionized electrons has been defined using the Bohmian trajectories, which allowed to look at the tunnel exit problem from the perspective offered by the Bohmian QM.

In our earlier works Ivanov and Kim (2022a, b); Ivanov et al. (2023) we described an alternative method that allowed us to extract from the TDSE information about the time-development of the ionization process for times within the laser pulse duration. The method is based on the analysis of two-time correlation functions, computed using the time-dependent wave-function describing evolution of the system, which was obtained by solving the TDSE numerically. In essence, this procedure allows us to use the notion of the conditional probability, where the condition is imposed at an instant when the laser pulse is gone. In other words, we formulate the questions about different observables characterizing the electron motion in the following way: What would be the probability of observing a given value of a certain observable during the laser pulse, provided that the electron is found in a given state at the end of the laser pulse? We applied this technique to study the evolution of the electron velocity distribution in strong field ionization Ivanov and Kim (2022b) and to study trajectories of the electron wave-packets for the process of the frustrated tunneling ionization (FTI) Ivanov et al. (2023). Here we report a study of evolution of the ionized wave-packets for the process of strong field ionization, based on the analysis of the information obtained from the numerical solution of the TDSE.

Atomic units with $\hslash=1$, $e=1$, $m=1$ where $e$ and $m$ being the charge and the mass of the electron are used throughout the paper.

We consider an atom interacting with a linearly polarized laser pulse which we define in terms of its vector potential: $\displaystyle\bm{E}(t)=-{\partial{\bm{A}}(t)\over\partial t}$, where:

Here $T_{1}=N_{c}T$ is the total pulse duration and $T=2\pi/\omega$ is the optical cycle (o.c.) corresponding to the central frequency $\omega=0.057$  a.u. (the wavelength of 800 nm). In the calculations below we use pulses with $N_{c}=4$. The target system is described using the single-active electron (SAE) approximation and a spherical potential $V(r)$. As targets, we will consider the hydrogen atom with $V(r)=-1/r$, a model atom with a short range (SR) Yukawa-type potential $V(r)=-1.903e^{-r}/r$ and the Ar atom described by means of an effective potential Sarsa et al. (2003). The target atom is initially in the ground state $|\phi_{0}\rangle$, which is an $s$ state for the hydrogen and Yukawa atoms (both with the ionization potential of 0.5 a.u.) and a $p$ state with the energy -0.59 a.u. in the case of the Ar atom.

We have shown in earlier works Ivanov and Kim (2022a, b); Ivanov et al. (2023) that tunneling ionization dynamics can be studied in detail by analyzing suitably chosen two-time correlation functions:

where the operators $\hat{A}^{H}(t)$ and $\hat{B}^{H}(t)$ in (2) are taken in the Heisenberg representation, and $|\phi_{0}\rangle$ is the initial state of the system. The particular choice of the operators $\hat{A}$ and $\hat{B}$ in Eq. (2) is dictated by the nature of the problem under consideration. We have shown in Ivanov and Kim (2022b); Ivanov et al. (2023) that by choosing for $\hat{B}^{H}(t)$ the Heisenberg form $\hat{Q}^{H}(T_{1})$ of a suitable Schrödinger projection operator $\hat{Q}$ one can study the dynamical development of various ionization processes. The reason why we may expect the correlation function (2) to provide a useful dynamical information with such a choice of $\hat{B}$ can be easily understood if in Eq. (2) we transform the operators to the more familiar Schrödinger picture:

where $\hat{U}(t,0)$ is the operator describing quantum evolution of the system, so that the wave function of the system at time $t$ is $\Psi(t)=\hat{U}(t,0)\phi_{0}$. Applying the transformation (3) we rewrite Eq. (2) as:

Let us assume, for instance, that $\hat{B}=\hat{P}$, where $\hat{P}$ is the projection operator on the continuous spectrum of the field-free atomic Hamiltonian, and $t_{2}=T_{1}$, where $T_{1}$ is the moment of time when the laser pulse is gone. Then, according to the well-known projection postulate of QM Lambropoulos and Petrosyan (2007), the ket-vector $\hat{P}\Psi(t_{2})\rangle$ represents, apart from an unimportant normalization factor, the wave-function of the system immediately after the measurement that has found the electron in an ionized state. Eq. (4) can therefore be interpreted as a quantum-mechanical amplitude of finding an electron in the state $|\hat{A}\Psi(t_{1})\rangle$ at the moment $t=t_{1}$ provided that the electron has been found in an ionized state after the end of the pulse. With a suitable choice of the operator $\hat{A}$ (we will discuss this choice in more detail below) we can now have a glimpse of the dynamical characteristics of the ionized electrons for the moments of time $t<T_{1}$ within the laser pulse duration.

Similarly, if we use $\hat{B}=\hat{I}-\hat{P}$ in Eq. (4) and again choose $t_{2}=T_{1}$, the expression for the correlation function can be interpreted as an amplitude of finding an electron in the state $|\hat{A}\Psi(t_{1})\rangle$ at the moment $t=t_{1}$ provided that the electron remains bound after the end of the pulse, which allows us to study dynamics of the FTI process for $t<T_{1}$.

We can concentrate on various aspects of the electron dynamics by choosing the operator $\hat{A}$ in Eq. (4) appropriately. We can choose, for instance, $\hat{A}$ to be a projection operator in momentum space. This choice together with $\hat{B}=\hat{P}$ allows us to study the development of the ionized electron velocity distribution Ivanov and Kim (2022b). The choice of $\hat{B}=\hat{I}-\hat{P}$ and a coordinate space projection operator for $\hat{A}$ allowed us to study evolution of the FTI electrons in coordinate space. We exploit below yet another possibility, using the following Schrödinger operators in the definition of the correlation function (4):

where the components of the ket vector $|\phi_{z_{0}}\rangle$ in the position representation are:

In Eq. (6) $N$ is the normalization factor. The ket vector $|\phi_{z_{0}}\rangle$ and its components in the coordinate basis given by (6) depend on the parameters $z_{0}$ and $a$, defining a point in space with the coordinates $(0,0,z_{0})$ and the resolution with which we look at the neighborhood of this point. In the calculations below we used $a=4\ln 2$ which gives us the spatial resolution of approximately one atomic unit. We use the position representation in all the calculations below. In this representation action of the projection operator $\hat{A}_{z_{0}}$ on a state vector $|\Phi\rangle$ with the components $\Phi({\bm{r}})=\langle{\bm{r}}|\Phi\rangle$ can be found as:

We choose $t_{2}=T_{1}$ in Eq. (4), where $T_{1}$ is duration of the laser pulse and we will be looking at various moments of time $t\leq T_{1}$. We will be thus studying the correlation function:

for $t\leq T_{1}$, with the operators $\hat{A}$ and $\hat{P}$ specified in Eq. (5). It is clear from the above discussion that with such a choice, Eq. (8) can be interpreted as giving us (apart from an unimportant normalization factor) a quantum-mechanical amplitude of finding the electron near the point with the coordinates $(0,0,z_{0})$ at the time $t$ provided that the electron will be found in an ionized state after the end of the pulse. In other words, this expression provides a means of studying trajectories of the ionized electrons during the laser pulse.

To calculate the correlation function (8) we use a procedure similar to the one we used previously in Ivanov and Kim (2022a, b); Ivanov et al. (2023), and we will only briefly describe the technical details. The calculation can be reduced to multiple solutions of the 3D time-dependent Schrödinger equation (TDSE):

where $\displaystyle H_{\rm atom}={\hat{\bm{p}}^{2}\over 2}+V(r)$ is the field free atomic Hamiltonian and $\hat{H}_{\rm int}(t)$ is the atom-field interaction Hamiltonian for which we use the length form:

We first propagate the TDSE forward in time on the interval $(0,T_{1})$, using ground atomic state as the initial state, obtaining position representation of the state vector $|\Psi(T_{1})\rangle$. Acting on $|\Psi(T_{1})\rangle$ with the projection operator $\hat{P}$ we obtain the wave-function $\Phi({\bm{r}})$ corresponding to the vector $|\Phi\rangle=\hat{P}|\Psi(T_{1})\rangle$. To find the vector $\hat{U}(t,T_{1})\hat{P}\Psi(T_{1})\rangle$, that we need to compute the matrix element in Eq. (8), we propagate the TDSE (9) backward in time using $\Phi({\bm{r}})$ as an initial (or rather final) wave-function, obtaining the vector $|\Phi(t)\rangle$ with the components $\Phi({\bm{r}},t)$ for the times $t$ within the laser pulse. Simultaneously, we propagate backward in time the vector $|\Psi(T_{1})\rangle$, obtaining the vector $|\Psi(t)\rangle$ and the wave-function $\Psi({\bm{r}},t)$- solution to the TDSE for the times $t$ within the laser pulse. Of course, $\Psi({\bm{r}},t)$ had already been computed during the first, forward run of the TDSE, but we cannot store it in memory for all the times $t$ we need as it would require too much memory space. We recompute it again, therefore, in the process of the back-propagation of the TDSE.

Calculating overlaps of $|\Phi(t)\rangle$ and $|\hat{A}_{z0}\Psi(t)\rangle$ for a given $z_{0}$ and for a given set of times $t$ (we use the grid of $t$ with twenty points for every optical cycle), we obtain the correlation function $C(z_{0},t)$ defined in Eq. (8). A single calculation using the backward propagation that we described above, gives us $C(z_{0},t)$ for the whole grid of $t$ and a single $z_{0}$. The procedure is repeated for different values of $z_{0}$. More specifically, we used a grid of a hundred $z_{0}$-values equally spaced on the interval $(-60$ a.u., $60$ a.u.).

The TDSE was solved numerically using the procedure tested and described in detail earlier Ivanov (2014); Ivanov and Kheifets (2013); I.A.Ivanov et al. (2016). The procedure relies on representing the coordinate wave-function as a series in spherical harmonics with the quantization axis along the laser polarization direction. Spherical harmonics with orders up to $L_{\rm max}=50$ were used. The radial variable was treated by discretizing the TDSE on a grid with a step-size $\delta r=0.05$ a.u. in a box of size $R_{\rm max}=400$ a.u. The initial ground state of the system was obtained by using a variational calculation employing the Slater basis set I.A.Ivanov and Y.K.Ho (1999) with subsequent propagation in imaginary time Bader et al. (2013) on the spatial grid we described above. The necessary convergence checks were performed. As we discussed above, to calculate the correlation function (4) we have to propagate the TDSE both forward and backward in time. That was achieved by using the matrix iteration method Nurhuda and Faisal (1999).

In summary, we devised a procedure based on the correlation function analysis and employed it to study ionization dynamics of the three atomic targets: the SR Yukawa, hydrogen and Ar atoms. The starting point of our analysis is the time-dependent wave function returned by a numerical solution of the TDSE. Our approach allows us to look closely at early stages of the photo-electron evolution and to separate various components of the wave function which, at later stages, will contribute to distinct outcomes of the laser-atom interaction. We achieve this result by reformulating the problem in terms of the conditional amplitudes, i.e., the amplitudes describing outcomes of measurements of different observables provided that electron is found in the ionized state after the end of the pulse. By choosing electron coordinate as such an observable we were able to track the motion of the ionized electron wave-packets basing on an ab initio TDSE calculation.

Our study demonstrates the somewhat limited character of the notion of a photo-electron trajectory for the description of the ionization process. The true photo-electron dynamics is more complex and resemble more a “quantum cloud” expansion. We demonstrate that the photo-electron wave-packets obtained in this way follow closely the CTMC trajectories at the initial stages of the evolution both for the SR and Coulomb systems. However, at the later stages of the evolution, the picture becomes more complicated due to the spread and interference of the wave-packets originated at different field maxima.

In the present work we considered these effects using a quantum mechanical approach based on the correlation function analysis. Alternative description might be based on incorporating these effects into the trajectory based methods. Interference effects can be included in the consideration following the prescriptions of the QTMC method Li et al. (2014); Fang et al. (2019) or the semiclassical two-step model for strong-field ionization Shvetsov-Shilovski et al. (2016), which supply each trajectory with a phase accumulated along the trajectory and allow thus to describe the interference effects. The dispersion effect could be described analogously to the description we obtained in the SFA-based model we presented above, where this effect manifests itself as a spread of the wave-packet moving along the classical trajectory. In the simplest case given by the Eq. (24) the spread does not depend on the core potential and is described by a simple analytical formula. Such a trajectory based description of the interference and dispersion effects possesses the advantage of the trajectory based methods, since it can be applied for more complex targets such as molecules, for which use of the TDSE based technique becomes prohibitively computationally demanding.

Our approach also allowed us to obtain information about coordinate and velocity electron distributions at the tunnel exit.

This work was supported by the Institute for Basic Science grant (IBS-R012-D1) and the National Research Foundation of Korea (NRF), grant funded by the Korea government (MIST) (No. 2022R1A2C3006025). Computational works for this research were performed on the IBS Supercomputer Aleph in the IBS Research Solution Center. IAI wishes to thank the Australian National University for hospitality.

We use the well-known expression for the SFA ionization amplitude Popov (2004):

where ${\bm{A}}(t)$ is the vector potential (1) of the pulse.

The physical meaning of the amplitude (16) is that it gives us the momentum space wave-function of the ionized wave-packet. Fourier transform of $a_{{\bm{p}}}(t)$ will give us then coordinate wave-function of the ionized wave-packet:

Below, we will be interested in the absolute value $|\Psi_{ion}({\bm{r}},t)|$ of the coordinate wave-function (15). Choice of the length or velocity gauges to describe the atom-field interaction is, therefore, immaterial for our purposes and we use the velocity gauge in Eq. (14) which makes the formulas somewhat simpler.

We consider expression (14) for the ionization amplitude for times $t$ within the interval $(0,T_{1})$ of the laser pulse duration. To evaluate expression (14) we employ the SPM, supplemented with the rule used in Yudin and Ivanov (2001) that for $t<T_{1}$ we need to consider only the saddle points $t_{s}$ with $Re(t_{s})<t$ for the evaluation of the integral in Eq. (14). Following the standard prescriptions of the SPM we obtain:

where $t_{s}$ are saddle-points of the integrand in Eq. (14), satisfying the SPM equation $\displaystyle({\bm{p}}+{\bm{A}}(t_{s}))^{2}+2I_{p}=0$ and:

We can compute the Fourier transform defining the coordinate wave-function (15) using the SPM again. One can see that it is the region of small momenta ${\bm{p}}$ that dominate the integral in Eq. (15). It is sufficient, therefore, to expand the action in Eq. (17) in powers of ${\bm{p}}$ keeping only the constant, linear and quadratic terms:

A simple integration than will give us for the coordinate wave-function of the ionized wave-packet:

where we have absorbed all constant factors into the factors $C_{s}$.

To see the physical meaning of Eq. (19) we have to take a closer look at the coefficients of the expansion (18). The integration path in Eq. (17) can be chosen to consist of two segments: a vertical line $(t_{s},Re(t_{s}))$ descending on the real time axis and a horizontal segment $(Re(t_{s}),t)$. We can then represent the action in Eq. (17) as a sum $S(t_{s},t,{\bm{p}})=S(t_{s},Re(t_{s}),{\bm{p}})+S(Re(t_{s}),t,{\bm{p}})$, where

with

and:

with

It is customary to refer to the action (20) as describing under-the-barrier, and (22) as describing post-ionization motion. Of course, as we mentioned in the Introduction, this division is arbitrary to a degree, since it corresponds to a particular choice of the integration path in Eq. (17) which is not unique.

The tilted coefficients in Eq. (21) and Eq. (23) are themselves functions of ${\bm{p}}$. This dependence is due to the dependence of the saddle point position $t_{s}$ on the momentum. To obtain the untilted quantities in Eq. (18) we must renormalize these coefficients, by expanding $t_{s}$ in powers of momentum components and keeping only constant, linear and quadratic terms in the resulting expressions. Corresponding formulas become rather bulky and add little to understanding the physical picture, we will not present them here. In the actual calculation reported below we performed all the necessary re-expansions numerically.

Before presenting results of this calculation we will first illustrate the physical meaning of Eq. (19) by making a few simplifying assumptions leading to more transparent formulas. Let us suppose first that we can drop under-the-barrier part of the action and substitute the expressions for the tilted coefficients from Eq. (23) in Eq. (19). We obtain then:

with ${\bm{\beta}}_{s}(t)$ given by the second of equations (23). It is not difficult to see that for each $t_{s}$ exponential factor in Eq. (24) describes evolution of an electron prepared at the moment $t=Re(t_{s})$ in the state described by a delta function $\delta({\bm{r}})$, which evolves subsequently under the action of the laser field only. For each $t_{s}$ the corresponding evolving wave-packet is weighted by a factor $e^{-i\gamma_{s}}$. Since $\gamma_{s}$ is complex, the exponential function $e^{-i\gamma_{s}}$ is sharply peaked around the field maximum nearest to $t_{s}$ (it is this factor, in fact, that gives the characteristic exponential dependence of the ionization probability on the field strength in the ADK and similar formulas). We obtain thus a simple picture of very narrow electron wave-packets created at times near the local field maxima and propagating subsequently under the action of the laser field.

Including the under-the-barrier part of the action makes this picture more realistic. Assuming that $\alpha_{s}(t)=\alpha^{s}_{1}(t)+\alpha^{s}_{2}(t)$, where $\alpha^{s}_{1}(t)$ and $\alpha^{s}_{2}(t)$ are given by Eq. (21) and Eq. (23), we obtain from Eq. (19):

where $t_{s}^{0}$ is the zero order term in the expansion of $t_{s}$ in powers of momentum. The role that this correction plays in Eq. (25) is quite clear. It is easy to see that exponential factor in Eq. (25) describes now the spread and motion of the wave-packet prepared initially in a state describe by a Gaussian $e^{-cr^{2}}$ with $c={1\over 2Im(t_{s})}$, and evolving subsequently under the action of the laser field only. Role of the under-the-barrier part of the coefficient $\alpha$ consists, therefore, in giving a non-zero initial spread to the ionized wave-packet. Similarly, one can see, that under-the-barrier part of the coefficient $\beta$ leads to the non-zero initial value of the coordinate.

Results for the coordinate wave-function of the ionized electron provided by Eq. (19) by systematically obtaining terms of the expansion Eq. (18) from Eq. (21) and Eq. (23) are shown in Fig. 7. The re-expansions needed to compute coefficients in Eq. (18) from Eq. (21) and Eq. (23) were done numerically. The figure shows absolute value $|\Psi_{ion}({\bm{r}},t)|$ along the line ${\bm{r}}=(0,0,z)$ computed for the field parameters used in Fig. 2a. We note the slight offset of the $z$-coordinate of the maxima of the correlation function at the staring points of the classical trajectories. This offset is due to the fact that the classical trajectory calculations use the FDM initial coordinates values. The initial coordinates in the SFA calculation, on the other hand, are essentially determined by the second equation (21), which does not include ionic potential, and differs, therefore, from the FDM value. One can see, nevertheless, that the plot shown in Fig. 7 reproduces qualitatively the main features of the TDSE correlation pattern. The wave-packets follow initially the classical trajectories, broadening subsequently. We also see the structures appearing at the latest stages of evolution, for times near the end of the pulse, which are reminiscent of the structures seen in Fig. 7. These structures disappear if we retain in Eq. (19) only the contribution of the main field maximum, and are, therefore, a manifestation of the interference of the wave-packets born at different field maxima.