Density-based one-dimensional model potentials for strong-field simulations in $\text{He}$, $\text{H}_{2}^{+}$ and $\text{H}_{2}$

TL;DR

This paper develops and tests improved one-dimensional models for atomic and molecular systems under strong laser fields, achieving high accuracy in simulating dynamics and harmonic spectra compared to full 3D models.

Contribution

It introduces corrected 1D Coulomb potentials derived from 3D density corrections, enhancing the accuracy of strong-field simulations for He, H2+, and H2.

Findings

1D models closely match 3D simulation results.

High-order harmonic spectra agree well with 3D results.

Significant improvement in atomic quantity dynamics accuracy.

Abstract

We present results on the accurate one-dimensional (1D) modeling of simple atomic and molecular systems excited by strong laser fields. We use atomic model potentials that we derive from the corrections proposed earlier using the reduced ground state density of a three-dimensional (3D) single-active electron atom. The correction involves a change of the asymptotics of the 1D Coulomb model potentials while maintaining the correct ground state energy. We present three different applications of this method: we construct correct 1D models of the hydrogen molecular ion, the helium atom and the hydrogen molecule using improved parameters of existing soft-core Coulomb potential forms. We test these 1D models by comparing the corresponding numerical simulation results with their 3D counterparts in typical strong-field physics scenarios with near- and mid-infrared laser pulses, having peak intensities in the $10^{14}-10^{15}\,\mathrm{W/cm}^2$ range, and we find an impressively increased accuracy in the dynamics of the most important atomic quantities on the time scale of the excitation. We also present the high-order harmonic spectra of the He atom, computed using our 1D atomic model potentials. They show a very good match with the structure and phase obtained from the 3D simulations in an experimentally important range of excitation amplitudes.

Equations

$$i\frac{\partial {\Psi }^{3D}}{\partial t}=\left[T_z+T_{\rho }+V+z{\mathcal{E}}_z(t)+(N-1)V_{\mathrm{H}}\right]{\Psi }^{3D},$$

$$T_z=-\frac{1}{2}\frac{{\partial }^2}{\partial z^2},T_{\rho }=-\frac{1}{2}\left[\frac{{\partial }^2}{\partial {\rho }^2}+\frac{1}{\rho }\frac{\partial }{\partial \rho }\right].$$

$$V(z,\rho )=-\frac{1}{\sqrt{{\rho }^2+(z-\frac{d}{2}{)}^2}}-\frac{1}{\sqrt{{\rho }^2+(z+\frac{d}{2}{)}^2}}$$

$$V_{\mathrm{H}}(\mathbf{r},t;{\Psi }^{3D})=\int \frac{|{\Psi }^{3D}({\mathbf{r}}^{\prime },t){|}^2}{|\mathbf{r}-{\mathbf{r}}^{\prime }|}{\mathrm{d}}^3{\mathbf{r}}^{\prime }.$$

$$g(t)=1-{\left|⟨{\Psi }^{3D}(t=0)|{\Psi }^{3D}(t)⟩\right|}^2,$$

$${\varrho }^{3D}(z,\rho ,t)=N{\left|{\Psi }^{3D}(z,\rho ,t)\right|}^2,$$

$$H(t)=H_0+z{\mathcal{E}}_z(t)=T_z+V_0(z)+z{\mathcal{E}}_z(t).$$

$$i\frac{\partial }{\partial t}\Psi \left(z,t\right)=H(t)\Psi \left(z,t\right)$$

$${\varrho }_z(z)=2\pi \int {\varrho }^{3D}(z,\rho )\rho \mathrm{d}\rho ,$$

$$V_0(z)=E_0+\frac{1}{\sqrt{{\varrho }_z(z)}}{T}_z\sqrt{{\varrho }_z(z)},$$

$$V_{0,\mathrm{Sc}}(z)=-\frac{Z^{\ast }}{\sqrt{z^2+{\left[\alpha \left(Z^{\ast }\right)\right]}^2}}$$

$$Z^{\ast }=\frac{Z}{2},{\alpha }^2=\frac{1}{16(Z^{\ast }{)}^2}=\frac{1}{4Z^2}\text{with }{E}_0=-\frac{Z^2}{2}.$$

$$s(f)=\min \left(1+0.03{\left(100f-1\right)}^2,1+\left|100f-1\right|\right),$$

$$V_0^{(\mathrm{M})}(z_i;d)=E_0+\frac{1}{\sqrt{{\varrho }_z(z_i;d)}}{\tilde{T}}_z\sqrt{{\varrho }_z(z_i;d)},$$

$$V_{0,\mathrm{Sc}}^{(\mathrm{M})}(z;d)=V_{0,\mathrm{Sc}}\left(z-\frac{d}{2}\right)+V_{0,\mathrm{Sc}}\left(z+\frac{d}{2}\right)$$

$$H_{\mathrm{He}}^{1D}=\sum _{k=1}^2\left[T_{z_k}+V_{0,\mathrm{Sc}}(z_k)+z_k{\mathcal{E}}_z(t)\right]-\frac{1}{2}{V}_{0,\mathrm{Sc}}(z_1-z_2).$$

$${\varrho }_z(z,t)=2\int {\left|\Psi (z,z_2,t)\right|}^2\mathrm{d}{z}_2$$

$$g(t)=1-\frac{1}{2}{\left|\int \sqrt{{\varrho }_z(z,0)}\sqrt{{\varrho }_z(z,t)}\mathrm{d}z\right|}^2$$

$$H_{{\mathrm{H}}_2}^{1D}=\sum _{k=1}^2\left[T_{z_k}+V_{0,\mathrm{Sc}}^{(\mathrm{M})}(z_k;d)+z_k{\mathcal{E}}_z(t)\right]-\frac{1}{2}{V}_{0,\mathrm{Sc}}(z_1-z_2).$$

$${\mathcal{E}}_z(t)=F\cdot {\sin }^2\left(\frac{\pi t}{N_{\text{Cycle}}T}\right)\cos \left(\frac{2\pi t}{T}\right),$$

$$s(f)=\min \left(1+0.01{\left(100f-1\right)}^3,1+2.7\left|100f-1\right|\right),$$

Full text

I Introduction

The interpretation of typical experiments in attosecond and strong-field physics, including the pioneering results in Refs. Hentschel et al. (2001); Kienberger et al. (2002); Drescher et al. (2002); Baltuška et al. (2003); Uiberacker et al. (2007); Hommelhoff et al. (2009); Schultze et al. (2010); Haessler et al. (2010); Pfeiffer et al. (2012); Shafir et al. (2012); Ranitovic et al. (2014); Ciappina et al. (2017); Azoury et al. (2017), often relies on the quantum description of the involved atomic system driven by a strong laser pulse Keldysh (1965); Varró and Ehlotzky (1993); Lewenstein et al. (1994); Protopapas et al. (1997); Ivanov et al. (2005); Gordon et al. (2005); Krausz and Ivanov (2009); Frolov et al. (2012); Peng et al. (2015). Despite recent progress in analytic and numerical solution techniques Clarke et al. (2018); Jiang and Tian (2017); Kidd et al. (2017); van Dijk et al. (2017); Lackner et al. (2017); Kasza et al. (2018); Efimov et al. (2018); Patchkovskii and Muller (2016), the exact solution of the corresponding true 3D Schrödinger equation is beyond reach in this non-perturbative range (except for the simplest cases), which justifies the importance of good approximations.

If the strong driving laser pulse is linearly polarized then the most important features of the resulting quantum dynamics can usually be captured by a one-dimensional (1D) approximation Javanainen et al. (1988); Su and Eberly (1991); Bauer (1997); Chirilă et al. (2010); Silaev et al. (2010); Sveshnikov and Khomovskii (2012); Gräfe et al. (2012); Czirják et al. (2000); Benedict et al. (2012); Czirják et al. (2013); Geltman (2011); Szabó et al. (2013); Baumann et al. (2015); Teeny et al. (2016); Nishi et al. (2019). Such 1D models have increasing importance for longer laser wavelengths where typical 3D simulations become inefficient Liu and Thumm (2014); Miyagi and Madsen (2016). These typically use various 1D model potentials to account for the motion of the atomic system along the direction of the laser polarization. However, the particular model potential can strongly influence some of the 1D results and their quantitative comparison with the true three-dimensional (3D) results is usually non-trivial Bandrauk et al. (2009); Gräfe et al. (2012); Majorosi et al. (2017); Ivanov et al. (2019).

We addressed this problem for a single active electron atom in Ref. Majorosi et al. (2018): we introduced the density-based 1D model potential and, based on its features, we also found improved parameters for other well-known 1D model potentials. The promising strong-field simulation results inspired us to extend our modeling approach to simple atomic systems like the hydrogen molecular ion, the helium atom and the hydrogen molecule, which is the subject of the present paper. Our key idea is to require the Coulomb asymptotes in the 1D model potentials to be equal to those obtained from the corresponding reduced 3D ground state single-electron density along the direction of the laser polarization. Then we present the results of careful numerical simulations of strong-field ionization scenarios using these 1D model potentials, considering nearly single-cycle laser pulses with carrier wavelengths of 725 nm and 3045 nm. Comparing them with the corresponding 3D simulation results, we make a conclusion about the recommended use of these 1D model systems. We use atomic units in this paper.

II 3D reference systems

In this section, we specify in more detail the strong-field modeling of the selected three-dimensional systems: the helium atom, the hydrogen molecular ion and the hydrogen molecule driven by a linearly polarized laser pulse. We also outline the underlying numerical simulations, the results of which we use later as reference when we compare the corresponding one-dimensional results.

Although a suitable laser pulse may create also vibrations and rotations in a diatomic molecule, we assume the nuclear motion to be frozen throughout this paper, and we set the molecular axis parallel to the polarization of the laser pulse.

For the H${}_{2}^{+}$, we solve the three-dimensional Schrödinger equation, both to compute the ground state, and to obtain the time-evolution when driven by the laser pulse. For the two-electron systems, He and H2, we chose the time-dependent Hartree-Fock (TDHF) approach as the reference model, using a single atomic orbital in real three-dimensional space. According to Ref. Sato et al. (2016), this method provides a good approximation of the more elaborate and numerically demanding two-electron results for He, driven by a laser pulse with parameters very close to those in the present work, since the effect of electron correlation is relatively small for He. We assume that the ground states of He and H2 are spin singlets and that the laser pulse does not interact with the spin degrees of freedom, thus the orbital part of the two-electron wave function remains symmetric during the time evolution.

The governing equation of the electrons’ motion can be cast for all of the above cases into the following form, using cylindrical coordinates $\rho=\sqrt{x^{2}+y^{2}}$ and $z$:

[TABLE]

where the kinetic energy operator is split as

[TABLE]

The one-electron potential

[TABLE]

contains the Coulomb interaction with the nuclei. The parameter $d$ is the internuclear distance for $\mathrm{H}_{2}^{+}$ and $\mathrm{H}_{2}$, while for $d=0$ we get the helium atom (with its nucleus in the origin). The $z\mathcal{E}_{z}(t)$ term of (1) corresponds to the interaction with the laser field, polarized along the $z$-axis, using dipole approximation and length gauge Griffiths (2005); Bransden and Joachain (2003); Földi (2017). The $\mathcal{E}_{z}(t)$ denotes the electric field of the laser pulse evaluated in the origin, and we assume that it is present only after $t>0$, i.e. $\mathcal{E}_{z}(t\leq 0)=0$.

In (1), we distinguish the one- and the two-electron cases by the parameter$\ N$. For $N=2$ the electron-electron interaction is described in (1) by the time-dependent Hartree-potential which is given by

[TABLE]

The presence of this potential makes equation (1) nonlinear in $\Psi^{\mathrm{3D}}$. In actual computations, we obtain this potential by solving the corresponding discretized Poisson equation $\nabla^{2}V_{\mathrm{H}}=-4\pi|\Psi^{\mathrm{3D}}(\mathbf{r},t)|^{2}$  in cylindrical coordinates, to avoid the high dimensional integration.

This time-dependent Hamiltonian in (1) has axial symmetry around the polarization axis of the electric field of the laser pulse which makes the use of cylindrical coordinates practical, and provides efficient calculation of the reduced dynamics along the $z$-axis.

For actual simulations, we use the efficient numerical method described in Majorosi and Czirják (2016), which incorporates the singularity of the Hamiltonian directly, using the required discretized Neumann and Robin boundary conditions. We compute the ground states via imaginary time propagation with high-order split-operator approximations Chin et al. (2009), then we compute the time evolution up to a specified time $T_{\max}$.

To characterize the effects of the external field, we shall  use the ground state population loss in a single-electron wavefunction or in an electron-orbital, defined as

[TABLE]

which is to be compared with the corresponding function in the 1D models  we consider below.

In our considerations, the electron density given by

[TABLE]

plays an important role as we shall construct our one dimensional model potentials by an appropriate reduction of this quantity.

III 1D model systems

In order to model the above described 3D strong-field process in 1D, it is customary to use the following form of the time-dependent 1D Hamiltonian

[TABLE]

where the effects of a strong few-cycle laser pulse are to be modeled by the very same electric field $\mathcal{E}_{z}(t)$ as in 3D.

We want to verify the physical correctness of the above models by numerically solving the time-dependent Schrödinger equation

[TABLE]

and compare the time-dependent physical response of this system with that of the original 3D TDSE.

III.1 Overview of 1D density-based atomic model potentials

The main question in (7) is the form of $V_{0}(z)$ Javanainen et al. (1988); Su and Eberly (1991); Bauer (1997); Chirilă et al. (2010); Silaev et al. (2010); Sveshnikov and Khomovskii (2012); Gräfe et al. (2012); Czirják et al. (2000, 2013); Geltman (2011); Baumann et al. (2015); Teeny et al. (2016). One of the possibilities proposed in our previous work Majorosi et al. (2018) for one electron atomic systems was that we have introduced the reduced one dimensional ground state density depending only on the $z$ coordinate as

[TABLE]

where $\varrho^{\mathrm{3D}}(z,\rho)$ is the 3D ground state density. This made it possible to calculate the density based model potential

[TABLE]

where $\sqrt{\varrho_{z}(z)}$ stands as the reduced ground state wavefunction, while $E_{0}$ is the corresponding exact 3D ground state energy. This construction ensures that the reduced problem yields exactly the same properties for the ground state as does the original 3D atomic calculation. It is also an important feature that this form preserves the ground state energy in our case, since the original 3D problem has a long-range Coulomb asymptotics.

In addition to these important physical properties, in the cases of single-electron atoms an analytic expression could be calculated for $V_{0}(z)$ Majorosi et al. (2018), which is of the form of a short-range correction plus a 1D regularized Coulomb potential. This long range Coulomb-part had an asymptotic form of $-\frac{1}{2}Z/z$, where $Z$ is the nuclear charge in 3D. This inspired us to develop Majorosi et al. (2018) an alternative, improved soft-core Coulomb potential form, which is smooth and easy-to use:

[TABLE]

Here the parameter $Z^{\ast}$ in the numerator has been determined by requiring to obtain the same asymptotic behavior as that of (10), while $\alpha$ is a fitting parameter depending on $Z^{\ast}$, set by demanding the 1D atomic system to have the same ground state energy as the corresponding 3D system. In case of single-electron atoms, the parameters turned out to be Majorosi et al. (2018)

[TABLE]

By using these potentials in the solutions of the corresponding one dimensional time-dependent Schrödinger equation (8) the results reproduced the 3D system’s strong-field response quantitatively correctly for several physical quantities. A detailed comparison can be also found in Majorosi et al. (2018). In those tests the soft-core Coulomb form of the atomic model potential gave the best physical response, despite yielding a slightly less accurate ground state density.

It has been also found that even the power spectrum $p^{\mathrm{3D}}(f)$ of coherent high-order-harmonic generation obtained from 3D simulations could be recovered from the corresponding 1D spectrum $p(f)$ by scaling the latter as $p(f)/s(f),$ where $s(f)$ is the simple scaling function

[TABLE]

which turned out to be essentially independent from the strength and the form of the exciting pulse $\mathcal{E}_{z}(t).$

The one-dimensional models we introduced in Majorosi et al. (2018) proved to be more than simple toy models, by their capability of providing quantitatively comparable results to the 3D ones. In this article we test the physical relevance of these one-dimensional models by extending them to simple composite atomic and molecular systems like the one-electron diatomic molecule H${}_{2}^{+}$, as well as to two two-electron systems namely the He atom and the $\mathrm{H}_{2}$ molecule.

III.2 1D hydrogen molecular ion model

In our 1D model of the H${}_{2}^{+}$ we assume that the molecular axis is set along the $z$-direction and the internuclear distance $d$ is fixed (i.e. we do not consider nuclear motion). We are interested in the strong-field dynamics of the electron according to the 1D Schrödinger equation with the Hamiltonian (7). For the 1D model potential $V_{0}(z)$ in (7), we suggest and test two candidates, $V_{0}^{(\mathrm{M})}$ and $V_{0,\mathrm{Sc}}^{(\mathrm{M})}$ as follows.

Based on our earlier results, summarized in the previous section, we define the density-based hydrogen molecular ion potential on an equidistant grid $z_{i}$ as

[TABLE]

where $\varrho_{z}(z_{i};d)$ is the reduced density of the $1\sigma_{g}$ ground state of a 3D hydrogen molecular ion (with a fixed internuclear distance $d$). This potential is calculated numerically with the finite-difference version of the kinetic energy operator, denoted by $\widetilde{T}_{z}$ in (14), and using the numerically exact ground state energy which equals that of the 3D reference system Majorosi et al. (2018).

The other 1D model potential we propose to use for the H${}_{2}^{+}$ is the soft-core molecular model potential

[TABLE]

where we assume an implicit dependence of the parameters $Z^{\ast}$ and $\alpha^{2}$ in $V_{0,\mathrm{Sc}}$ on the parameter $d$. The value of $Z^{\ast}$ determining the Coulomb asymptotics of (15) is calculated from the potential given by (14). Then, $\alpha^{2}$ can be determined by setting the correct single-electron energy from the reference 3D hydrogen molecular ion calculation. We plot the shape of the corresponding 1D potentials with $d=2$, along with $\sqrt{\varrho_{z}(z_{i};d=2.0})$ in Fig. 1 as well as the values of the parameters $Z^{\ast}$ and $\alpha^{2}$ obtained in this way.

In these computations, we have typically employed a grid spacing of $\Delta z=0.2$, which yields about 3-4 digit accuracy in the ground state density of the reference 3D method. For an intermolecular distance $d=2.0$ – that is near the equilibrium distance of the hydrogen molecular ion – the ground state energy is $E_{0}=-1.1026$ based on both our 3D reference calculation and the 1D model using (14). Schaad and Hicks Schaad and Hicks (1970) give a very accurate result for the equilibrium distance as $d=1.9972$, and an electron energy of $E_{0}=-1.10334$. (The binding energy they gave was $-0.602634619$ a.u. which also incorporated the $1/d$  Coulomb repulsion energy of the protons.)

III.3 1D helium atom model

We turn now to model a two-electron system, the helium atom in one spatial dimension denoting the electron coordinates by $z_{1}$ and $z_{2}$. The key element of our model is that we replace the true Coulomb potentials with 1D improved soft-core potentials $V_{0,\mathrm{Sc}}$ of the form (11), both for the electron-nucleus and for the electron-electron interaction. By setting $Z^{\ast}=1$ in $V_{0,\mathrm{Sc}}$, we have the following Hamiltonian:

[TABLE]

We solve the corresponding Schrödinger equation using this Hamiltonian to get the wave function $\Psi(z_{1},z_{2},t)$. Let us remind that the latter is a symmetric function in the two spatial variables. For the time-evolution, we set the initial state as the field-free ground state at $t=0$ (using $\mathcal{E}_{z}(t\leq 0)=0$). We handle the numerical time evolution and the imaginary time-propagation (for obtaining the ground state energy) with a split-step finite difference method Wang (2005); Majorosi and Czirják (2016); Majorosi et al. (2018).

In this two dimensional formalism we calculate the physical quantities that can be derived from the following form of the reduced density

[TABLE]

from which the mean values and the root-mean-square deviations of the spatial coordinate follow straightforwardly. We also make an approximate formula of

[TABLE]

that we call ground state population loss per electron orbital. This is directly comparable to the ground state population loss of the 3D Hartree-Fock formalism.

The model potential $V_{0,\mathrm{Sc}}(z)$ implicitly depends on the fitting parameter $\alpha^{2}$, which is determined by setting the initial two-electron energy. We use two types of parametrization. We call the first “ab initio” parametrization, in which we use the single-electron model parameters of (12) with $Z=2$, corresponding to $Z^{\ast}=1$, and $\alpha^{2}=0.0625.$ For the ground state energy of this 1D model system we have obtained $E_{0}=-3.02$ without additional parameter fitting. This is to be compared with the ground state energy of the real helium: $-2.903$ Drake (1999), indicating an error of about $3.4\%$ in the ground state energy of this “ab initio” 1D model.

In the second parametrization of $V_{0,\mathrm{Sc}}(z)$, the $\alpha^{2}$ was modified to reproduce the 1D ground state energy $E_{0}=-2.903$ (i.e. it matches that of the real helium), with $Z^{\ast}=1$, and $\alpha^{2}=0.0694$. Since the 3D Hartree-Fock method yields a ground state energy $E_{0}=-2.860$, the corresponding 3D reference simulations have also been modified to match the accurate ground state energy of $-2.903$. (To this end, the Hartree-potential was multiplied by a factor that is slightly different from 1.)

III.4 1D hydrogen molecule model

Finally, we create the 1D model of the hydrogen molecule in a similar way to the case of the He atom. We use again a wave function of two variables $\Psi(z_{1},z_{2},t)$ and we replace the true Coulomb potentials by our improved 1D soft-core potentials (11) and (15), which gives the following Hamiltonian for the 1D model hydrogen molecule:

[TABLE]

Note that this Hamiltonian reduces to the one used for the helium atom in Eq. (16) for the limiting value of the internuclear distance $d=0.$ The 1D potentials $V_{0,\mathrm{Sc}}^{(\mathrm{M})}$ and $V_{0,\mathrm{Sc}}$ depend implicitly on the parameters $Z^{\ast}$ and $\alpha^{2}$ but we use identical parametrization for these potentials. We are going to test the case when $Z^{\ast}=0.5$ which equals to that of the hydrogen atom, and alternatively the case when the asymptotics determined by $Z^{\ast}$ agree with the one derived from the 1D hydrogen molecular ion model.

Regarding the energy values, Doma Doma et al. (2016) gave $-1.173427$ a.u. for the dissociation energy of the real 3D hydrogen molecule with $d=1.4$, which means $-1.88729$ a.u. electronic ground state energy. Our reference Hartree-Fock computation gives $E_{0}=-1.848$ a.u. ground state energy which means a relative error of about $2.1\%$ in our 3D reference scheme. For consistency, the ground state energy of the 1D model system is adjusted to this Hartree-Fock energy, by setting $\alpha{}^{2}$ to the proper value, as to be explained in section IV.3.

IV Results and comparison of the 1D and 3D simulations

In this section, we present and compare the results of strong-field simulations of the 1D models of the previous sections to the results acquired from the 3D reference models.

In these simulations, we model the linearly polarized few-cycle laser pulse with a sine-squared envelope function. The corresponding time-dependent electric field has non-zero values only in the interval $0\leq t\leq N_{\text{Cycle}}T$ according to the formula:

[TABLE]

where $T$ is the period of the carrier wave, $F$ is the peak electric field strength and $N_{\text{Cycle}}$ is the number of cycles under the envelope function. In subsections IV.1-IV.3, we model a near-infrared laser pulse by setting $T=100$, corresponding to a ca. $725\,\mathrm{nm}$ carrier wavelength, and $N_{\text{Cycle}}=3$ which gives a pulse with its main peak at its center. In order to simulate the effect of a mid-infrared laser pulse in subsection IV.4, we set $T=420$, corresponding to a carrier wavelength of $3045\,\mathrm{nm}$, and $N_{\text{Cycle}}=1.5$ which gives a pulse with a zero at its center. We plot these pulse shapes in Fig. 2. From Fig. 2 on, the vertical dashed lines denote the zero crossings of the respective $\mathcal{E}_{z}(t)$ electric field.

In our calculations, we set the typical step sizes as $\Delta z=0.2$ and $\Delta t=0.01$ since these are sufficient for the numerical errors to be within line thickness. We use box boundary conditions and we set the size of the box to be sufficiently large so that the reflections are kept below $10^{-8}$ atomic units in the wave function.

The results belonging to the correct reference 3D simulation of a given system are plotted in blue and labeled as “3D-reference”. We also plot dashed blue overlays for these reference curves for clarity.

IV.1 Hydrogen molecular ion driven by a near-infrared laser pulse

In this section we present our results for the hydrogen molecular ion with two selected external field amplitudes $F$ and three selected molecular distance parameters $d$. We have investigated the following model potentials, the results of which are plotted on each of the figures belonging to this section (from Fig. 3 to Fig. 4): the density-based molecular 1D model potential (14) (in purple), and the 1D soft-core model (15) using the correct asymptotes (in gold). In these figures we plot the time-dependence of the following physical quantities: the mean value $\left\langle z\right\rangle(t)$, the standard root-mean-square deviation $\sigma_{z}(t)$, and the ground state population loss $g(t)$, as determined by the time-dependent wave function.

Let us begin the discussion with a hydrogen molecular ion frozen in its equilibrium distance ($d=2.0$). From time-independent calculations of the reduced ground-state density, we inferred that the correct 1D density-based model does not behave like the composition of two 1D density-based atomic binding potentials having Coulomb asymptotes with $Z^{*}=0.5$. Instead, in this case, these long range 1D Coulomb asymptotes have the value of around $Z^{*}=0.625$, that we obtained by numerical calculations based on (14). This means that for this separation distance, the Coulomb cores exhibit a screening-like effect in the 1D model. Accordingly, we parametrized the 1D improved soft-core model using $d=2.0$ with $Z^{*}=0.625$ and $\alpha^{2}=0.277$.

In Fig. 3 we can see the time-evolution results for this model system under the influence of the external fields with $F=0.1$ (left panels) and $F=0.15$ (right panels), corresponding to a peak intensity of $3.51\times 10^{14}\,\mathrm{W/cm}^{2}$ and $7.89\times 10^{14}\,\mathrm{W/cm}^{2}$, respectively. These results have comparable accuracy to that of our previous calculations in Ref. Majorosi et al. (2018) using single-electron atom models (although we use larger $F$ values here). In the case of $F=0.1$, hardly any ground state population loss occurs after the laser pulse, even though temporarily it raises to relatively larger values, in sync with the laser field. This behavior is quantitatively correctly predicted by our 1D models. The curves of $\sigma_{z}(t)$ show that the 1D improved soft-core model has somewhat lower accuracy for the 1D ground state ($\sigma_{z}(t)$ has a somewhat larger initial value in the soft-core case compared to 3D case), but its laser-field induced dynamics follows better the 3D results in overall. This latter is even more pronounced for the stronger field value $F=0.15$. The curves of the mean value $\left\langle z\right\rangle(t)$ also show good quantitative agreement.

Although these 1D models qualitatively reproduce even the finer details of the 3D ionization dynamics of $\mathrm{H}_{2}^{+}$, we should note that the 1D models yield high-order harmonics with larger amplitudes than the 3D results. This is mainly caused by the fact that the interference effects in the 1D density are stronger than in the reduced 3D density (since these latter are suppressed by the spatial integration along the $\rho$ coordinate). This feature of the HHG spectra was successfully accounted for by the scaling function (13) for a hydrogen atom Majorosi et al. (2018), which turned out to be suitable also for the present case of $\mathrm{H}_{2}^{+}$.

Next, let us discuss hydrogen molecular ion configurations where the parameter $d$ is different from its equilibrium value $d=2.0$. We chose for the case of closer nuclei $d=1.4$, which may be regarded as an instantaneously ionized static hydrogen molecule. Based on the calculation of the 1D density-based potential (14), we computed the parameter of the effective Coulomb asymptotes as $Z^{*}=0.597$, having a ground state energy of $E_{0}=-1.284$. Accordingly, the fitting parameter has to be $\alpha^{2}=0.2023$. For the case of nuclei with a larger distance, we chose the value $d=2.6$, which is closer to limit of the molecular dissociation. From the calculation of the 1D density based potential (14) we get its ground state energy as $E_{0}=-0.975$, and the parameters of the Coulomb asymptotes as $Z^{*}=0.647$ and $\alpha^{2}=0.351$. Note that the dependence of these parameters on $d$ is not negligible, and the parameter $Z^{*}$ appears to be increasing as $d$ is increasing. However, for nuclei very far apart, the asymptotic value of $Z^{*}=0.5$ should hold to describe a 1D H-atom and a proton. For the other limiting case, $d=0$, we get the 1D $\mathrm{He}^{+}$ with the value of $Z^{*}=0.5$ again.

In Fig. 4 we plot the results for $d=1.4$ (left panels) and $d=2.6$ (right panels) using the external electric field with $F=0.15$. We can see that the induced dynamics with $d=1.4$ is similar in behavior to the case of $d=2.0$ using $F=0.1$ (see Fig. 3). In the case of $d=2.6$, the ground state population loss is much larger than in the previous cases with the same $F$ value. Most importantly, both of these 1D models reproduce even the finer details of the curves of $g(t)$, $\sigma_{z}(t)$, $\left\langle z\right\rangle(t)$ of the true 3D dynamics.

We can also observe that the accuracy of the results improves as the value of $d$ increases, especially for the 1D improved soft-core potential, which is mainly due to the fact that the probability concentrated between the nuclei decreases with increasing $d$. This increasing accuracy suggests that these 1D potentials are suitable also to model strong field processes leading to molecular dissociation.

IV.2 Helium atom driven by a near-infrared laser pulse, high-order-harmonic spectra

In the following, we present our results regarding the 1D model helium atom based on Sec. III.3 under the influence of an external laser pulse, for selected values of $F$. In Figs. 5-7, we plot results of the following 1D model soft-core potentials: $(i)$ the 1D improved soft-core potentials with “ab initio” parametrization ($Z^{*}=1$, $\alpha^{2}=0.0625$) from (12) in red, $(ii)$ the 1D improved soft-core potentials with energy fitting ($Z^{*}=1$, $\alpha^{2}=0.0694$) in gold, and $(iii)$ the 1D usual soft-core potentials Bauer (1997) using the normal Coulomb asymptotes and with the same form of energy fitting as the previous case ($Z^{*}=2$, $\alpha^{2}=0.5474$) in green. Except for $(i)$, the time evolution of the 1D two-electron model systems starts from a ground state that has the same energy as the reference time-dependent Hartree-Fock simulation.

We show in Fig. 5 the time-dependence of the mean values $\left\langle z\right\rangle(t)$, the standard root-mean-square deviations $\sigma_{z}(t)$ and the ground state population losses per electron orbital $g(t)$, with $F=0.15$ (left panels) and $F=0.2$ (right panels), corresponding to a peak intensity of $7.89\times 10^{14}\,\mathrm{W/cm}^{2}$ and $1.40\times 10^{15}\,\mathrm{W/cm}^{2}$, respectively. For $F=0.15$, we can see that the dynamics induced by the laser field are weak, which is true both for the 3D reference and the 1D improved soft-core models, while the results of the usual 1D model are significantly off. Our “ab initio” parametrization $(i)$ of the 1D improved soft-core Coulomb forms is still quantitatively acceptable, and corresponds to a stronger bound than the real helium atom. This latter seems to be in agreement with the fact that it has a lower bound state energy by $0.1$ a.u. If we now look at the results corresponding to $F=0.2$, we can see that the results of our improved models are similar, albeit slightly off with stronger ionization. Based on the similarity to the results of $\left\langle z\right\rangle(t)$ and $\sigma_{z}(t)$, we also note that the quantity $g(t)$ defined in (18) indeed bears the same meaning as in (5) calculated from the Hartree-Fock orbital.

We also plot in Fig. 6 the $\left\langle z\right\rangle(t)$ and $g(t)$ results for $F=0.25$. We can see now that the 1D results are relatively up shifted compared to the previous figures in such a way that now the curves of 1D “ab initio” parametrization much better match the 3D reference. Since such an effect does not occur in the 1D models for the hydrogen molecular ions (or hydrogen atoms) under similar conditions, this may indicate the inaccuracy of our approximations at such high intensities, which presumably affects also the Hartree-Fock calculations.

Let us next investigate the high-order harmonic response of the particular 1D model helium. The accurate computation of the high-order harmonic spectrum is especially important in strong-field physics: its well-known characteristic features McPherson et al. (1987); Ferray et al. (1988); Harris et al. (1993); Krausz and Ivanov (2009); Gombkötő et al. (2016) represents the highly nonlinear atomic response to the strong-field excitation, and its suitable phase relations enable the generation of attosecond pulses of XUV radiation Farkas and Tóth (1992); Paul et al. (2001); Hentschel et al. (2001); Drescher et al. (2002); Kienberger et al. (2002); Carrera et al. (2006); Sansone et al. (2006). In accordance with Ref. Majorosi et al. (2018) we expect that the structure of the spectra based on the 1D and the 3D simulations is similar, but the amplitudes are larger in the 1D results. The reason for the latter is that rescattering on the ion-core is a much stronger effect in a 1D dynamics than in a 3D dynamics (amplifying the high-frequency oscillations in the 1D results), and on the other hand, the integration over the transverse directions decreases the effect of the oscillations of the 3D wave function on the mean values, like $\left\langle z\right\rangle(t)$. (These also cause the small oscillations are visible only on the 1D curves, but not on the corresponding 3D curves, of some of the plots in Figs. 5 and 6). This feature of the 1D models can be handled again by introducing a frequency dependent scaling function. In Fig. 7 we can see the scaled $p(f)/s(f)$ power spectrum of the second derivative of $\left\langle z\right\rangle(t)$, where we applied the following scaling function:

[TABLE]

which we obtained by fitting its parameters to reproduce spectrum for the case $F=0.15$. In Fig. 7 (a) we can see that the structure of the scaled spectra is indeed a good match compared to the reference Hartree-Fock results. The spectrum for the “ab initio” case has a similar structure, and it is shown with the same scaling function. If we compare (21) to (13) used for a single electron atom, we see that the 1D helium atom needs a stronger scaling (by about $2.7$ times) especially for the higher harmonics. What makes this scaled spectrum interesting, that it works for different configuration of electric fields: in Fig. 7 (b) and 7 (c) we applied the same scaling. For the larger field value of $F=0.2$ the improved soft-core results still replicate the 3D spectra really well. (However, there seems to be some deviation near the 20th harmonic.) For $F=0.25$, we can see that the spectrum is scaled properly, but the positions of the harmonics are slightly shifted from the 3D reference. Here, the “ab initio” model better describes the structure of the spectra, which is connected to the closer matching to the 3D reference of the other quantities of Fig. 6.

Finally we note that the match of the spectral phase is also very good, especially in the higher frequency range, which is of fundamental importance for the generation of isolated attosecond pulses.

In overall we can say that, the 1D improved soft-core models replicate the strong-field response of a real helium atom comparing to the time-dependent Hartree-Fock approximation up to around $F=0.2$. This includes the low frequency response of the mean motion on the level of the reduced density. The low dimensionality however, gives a much larger high-order-harmonic amplitude, which can be converted to the corresponding 3D spectra using one scaling function. The improved models discussed here appear to be quantitatively correct 1D models of the helium atom.

The merits of these 1D model potentials may also pave the way to simulate properties of a dilute medium, like an atomic gas-jet, used in actual strong field or attosecond physics experiments with He Shafir et al. (2012); Chen et al. (2012); Stooß et al. (2018). Previously, this was done Tosa et al. (2016); Heyl et al. (2016) by calculating Lewenstein’s integral Lewenstein et al. (1994), but using the density-based model potentials and integrating the low dimensional TDSE has now become also an option Berman et al. (2019).

IV.3 Hydrogen molecule driven by a near-infrared laser pulse

We present now the results for the 1D model hydrogen molecule as introduced in Sec. III.4. We apply the same two-electron formalism as for the helium atom, and we compare again the results with the corresponding 3D Hartree-Fock simulation as the reference. We also make use of our results given in Sec. IV.1 regarding the dependence of the Coulomb asymptotes on the intermolecular distance $d$ in the case of the hydrogen molecular ion model. Now the presence of the extra electron poses the question, whether it affects the values of the Coulomb asymptotes. To get the answer, we have considered two different parametrizations of the improved soft-core model $V_{0,{\rm Sc}}^{({\rm M})}$ in Eq. (19): $(i)$ the model using the corresponding hydrogen molecular ion asymptotes ($d=1.4$, $Z^{\ast}=0.597$, $\alpha^{2}=0.235$) and $(ii)$ soft-core models using atomic asymptotes ($d=1.4$, $Z^{\ast}=0.5$, $\alpha^{2}=0.112$). The $\alpha^{2}$ parameters were determined by fitting the ground state energies, so that the time evolution of the 1D two-electron model systems started from a ground state that had the same energy value as the reference (3D) time-dependent Hartree-Fock simulation. It turned out that a rather good agreement could be obtained with the asymptotes from the H${}_{2}^{+}$ model, i.e. from the choice $(i)$, especially for lower excitation amplitudes. This is demonstrated in Fig. 8  where the results are shown for two specific values: $F=0.07$ (left panels) and $F=0.1$ (right panels), in yellow for model $(i).$ For comparison the functions obtained with model $(ii)$ are also shown in orange.

The two particular $F$ values for which we show the results for H2 were chosen because the induced ionization response is similar in magnitude to that of the He atom with the amplitudes $F=0.15$ and $F=0.2$ shown in Fig. 5. If we compare the left panels of of Fig. 8 and Fig. 5, we can see that the yellow curves of the physical quantities $\left\langle z\right\rangle(t)$ and $\sigma_{z}(t)$, behave similarly and they show a good match with the 3D reference results if the ionization is low. We also note that the quantity $g(t)$ behaves in the same way, compared to the reference.  This agreement is a consequence of the picture that considers the He atom as the $d=0$ limit of H${}_{2},$ but in the case of He, which has more strongly bound electrons, larger electric field amplitudes are required to achieve a similar effect. In the right panels of the respective figures, which show simulation results with higher intensities we can see that for model $(i)$ an up-shift occurs, the size of which in Fig. 8 is again consistent with the results of Fig. 5. This value of up-shift may still be quantitatively acceptable for a one dimensional molecular model, especially because the reference Hartree-Fock method may certainly become inaccurate at these larger field strengths. We also note, that for higher intensities this up-shift becomes larger, causing that already for $F=0.15$ the model $(ii)$ matches the 3D reference better (similar to Fig. 6, but not shown). From this we can conclude that using the asymptote value $Z^{\ast}=0.597$ from the respective static hydrogen molecular ion configuration behaves like the 1D improved soft-core model of the helium calculations, i.e. it gives quantitatively correct results for the physical quantities in a similar manner shown here.

Overall it is impressive that using improved asymptotes with one dimensional soft-core Coulomb potentials the physics of a 3D hydrogen molecule becomes quantitatively reproducible by the corresponding 1D model system. From the results we can say that indeed the correct model is based on the density based potential, and its asymptotes are surely between the respective hydrogen molecular ion asymptote with $Z^{*}=0.597$ and the hydrogen atom asymptote with $Z^{*}=0.5$ in the tested peak electric field strength range. For more accurate tests, more accurate 3D reference simulation methods are necessary.

IV.4 $\mathrm{H}_{2}^{+}$ and He driven by a mid-infrared laser pulse

Although the established laser technology for strong-field physics is mainly in the spectral range from ultraviolet to near-infrared currently, the use of longer carrier wavelengths has several important advantages Wolter et al. (2015), which makes the use of mid-infrared laser pulses a promising new research direction in strong-field physics Wolter et al. (2014); Kübel et al. (2017); Quan et al. (2017); Ortmann et al. (2017); Pullen et al. (2017); Shao et al. (2017); Maurer et al. (2018); Shaaran et al. (2019); Ghimire and Reis (2019). However, the continuum electron wavepackets may travel much larger distances with increasing wavelength in the case of gas targets, which rises the numerical demand of the corresponding 3D quantum simulations sharply Liu and Thumm (2014); Miyagi and Madsen (2016). Thus, accurate 1D simulations are very important, therefore we present below a few promising results obtained with our 1D models for the hydrogen molecular ion and the helium atom, driven by the mid-infrared laser pulse specified after Eq. (20). (Since the 1D He atom models show similar accuracy as the 1D $\mathrm{H}_{2}$ models do, we do not consider the latter here.)

In particular, we test our 1D models of $\mathrm{H}_{2}^{+}$ (with $d=2.0$) and He in the same manner as we did in Sec. IV.1 and IV.2, respectively, with some simulation parameters adjusted to the requirements due to the mid-infrared laser pulse. We perform the simulations from 0 to 630 atomic time units, and we ensure that the box size is large enough to keep the numerical errors within line thickness by setting $z_{\min}=-2000$, $z_{\max}=2000$, $\rho_{\min}=0$, $\rho_{\max}=1750$ for all simulations in this subsection. Note also that the mid-infrared pulse has a zero crossing at its center (at $t=315$), and the two main peaks in opposite directions at $t=229.6$ and $t=400.4$ have a magnitude of ca. $0.795F$. The anti-symmetry of this pulse with respect to its center helps to keep the electron’s motion confined which makes the 3D simulations less demanding and more accurate.

We present the results for $\text{H}_{2}^{+}$ driven by a mid-infrared pulse of $F=0.15$ in the left panels of Fig. 9. Comparing these plots to those in Fig. 3, we see that the improved 1D soft-core potential gives again a better model of the 3D process than the density-based 1D potential does. However, the mid-infrared pulse creates somewhat different dynamics with respect to the near-infrared pulse: although $\left\langle z\right\rangle(t)$ in Fig. 3 (a) yet follows the near-infrared pulse shape, the $\left\langle z\right\rangle(t)$ in Fig. 9 (a) deviates from the mid-infrared pulse shape considerably more than expected based on Fig. 3 (b), due to the longer period. The longer period also enables the continuum wavepackets to spread for a longer time, thus the $\sigma_{z}(t)$ has considerably increased values in Fig. 9 (c). On the other hand, the peaks of the ground state population loss in Fig. 9 (e) are in accordance both in magnitude and in timing with those in Fig. 3 (e) and (f), but the final value of of $g(t)$ in Fig. 9 (e) is considerably less than that in Fig. 3 (f). This latter shows that strong-field ionization from the double-well potential of $\text{H}_{2}^{+}$, which has interesting internal and strong-field dynamics Chomet et al. (2019), is very sensitive to the peak value of the laser electric field which is ca. 20% lower for the mid-infrared than for the near-infrared pulse.

We present the results for He driven by a mid-infrared pulse with $F=0.2$ in the right panels of Fig. 9, which are to be compared to the right panels of Fig. 5. The 1D results corresponding to $V_{0,\mathrm{Sc}}(z)$ with $\alpha^{2}=0.0694$ give even better results than those for the near-infrared pulse, but the "ab initio" parameter for $V_{0,\mathrm{Sc}}(z)$ is clearly less accurate for the mid-infrared pulse. This suggests that the exact match of the 1D ground state energy with the 3D ground state energy becomes more important with increasing laser wavelength. The difference in the shapes of the $\left\langle z\right\rangle(t)$ for the near- and the mid-infrared pulse is clearly due to the different pulse shapes and the increased wavelength. The ground state population loss is very similar to Fig. 5 (f), both in behaviour and in magnitude, which is readily explained by taking into account that the increasing effect of the longer period of the mid-infrared pulse is largely compensated by the actually ca. 20% less peak electric field strength. However, the longer period of the mid-infrared pulse does considerably increase the values of the $\left\langle z\right\rangle(t)$ and $\sigma_{z}(t)$ curves, since the continuum part of the wave function has more time to travel and spread.

Summarizing this section, the carefully parametrized 1D soft-core model potentials work very well also in the case of a mid-infrared laser pulse. Note also, that the small oscillations visible on some of the 1D curves in Figs. 3-6 are absent now, due to the mid-infrared pulse shape which causes only one rescattering.

V Discussion and conclusions

In this paper we have shown that the density-based model potentials, developed in Majorosi et al. (2018) for strong field simulations in single-active-electron atoms can be extended to two-electron systems and simple molecules. Our results show that the modeling based on low dimensional Coulomb asymptotes selected by the reduced density works also for two-electron problems, not just for single-active electron systems. The molecular potentials we built in this way fulfilled the requirement that the 1D model should recover a quantitatively correct ionization response compared to the respective 3D molecular system under the influence of a linearly polarized external laser field. This improvement is mainly due to the accurate match of the 1D and 3D ground state energy and the correct 1D asymptotics which provides better 1D continuum states.

For H${}_{2}^{+}$ one of the molecular potentials has been obtained by using the reduced density according to Eq. (14). The other choice for this system was a sum of two soft-core Coulomb potentials as given by Eq. (11) with appropriate values of the fitting parameters $Z^{\ast}$ and $\alpha^{2}$. For the two-electron systems He and H2 the potentials were built by combining 1D improved soft-core Coulomb potentials from the corresponding single-electron density-based models. Here we have used two different sets of parameters to simultaneously reproduce the correct density based Coulomb asymptotics and ground state energies. We compared the results of numerical strong-field simulations for the complete 3D systems with our 1D molecular models excited with the same linearly polarized laser field. We have shown that our simpler 1D constructions provide impressive accuracy, for being 1D models, for the helium atom and the hydrogen molecular ion, driven by the near-infrared and even by the mid-infrared pulse. These two models performed exceptionally well, especially if the strong-field ionization response was relatively weak. For the case of the hydrogen molecule the correct asymptotic values of the potential turned out to be near to that of the hydrogen molecular ion with the same intermolecular distance. Experimental interest in strong field and attosecond processes of He Shafir et al. (2012); Chen et al. (2012); Stooß et al. (2018) inspired us to calculate high-order-harmonic generation spectra  by using our 1D model for He. It turned out that quantitatively correct spectra could be recovered with a simple scaling for different external electric fields.

These results overall provide more possibilities to explore. One can extend the description of the molecular systems and model them under strong-field circumstances with moving nuclei. In order to include the motion of the nuclei into the simulation using e.g. the Born-Oppenheimer or the Ehrenfest approximations Ullrich (2011); Bransden and Joachain (2003), the determination of $Z^{\ast}(d)$ and $\alpha^{2}(d)$ is required beforehand numerically. Such a $d$ dependent molecular model potential seems to be a promising way of modeling the true molecular dynamics in strong-field scenarios with electron wave functions in one spatial variable, which would make it especially effective. It is also a possibility to model a linear chain of atomic cores in 1D by reducing the proper 3D density in advance. Regarding the multiple electron systems, one can apply the time-dependent Hartree-Fock, or multiconfigurational Hartree-Fock approach Beck et al. (2000) to the 1D helium and hydrogen systems directly, which may enable massive performance gain, while providing quantitatively correct reduced dynamics and high-order-harmonic spectra. This would also enable to efficiently perform low dimensional strong-field calculations in gas-jets, or statistical mixtures.

Acknowledgements.

The authors thank Péter Földi, Katalin Varjú and Sándor Varró for stimulating discussions. This research was performed in the framework of the project Nr. GINOP-2.3.2-15-2016-00036 titled Development and application of multimodal optical nanoscopy methods in life and material sciences. The project has been also supported by the European Union, co-financed by the European Social Fund, Grant No. EFOP-3.6.2-16-2017-00005. Partial support by the ELI-ALPS project is also acknowledged. The ELI-ALPS project (GINOP-2.3.6-15-2015-00001) is supported by the European Union and co-financed by the European Regional Development Fund.