Table Of Contentsubmitted to Phys. Rev. A; 7 January 2011
Breakup of the aligned H molecule by xuv laser pulses: A time-dependent treatment
2
in prolate spheroidal coordinates
Xiaoxu Guan1, Klaus Bartschat1, and Barry I. Schneider2
1Department of Physics and Astronomy, Drake University, Des Moines, Iowa 50311, USA and
2Physics Division, National Science Foundation, Arlington, Virgina 22230, USA
(Dated: January 11, 2011)
We have carried out calculations of the triple-differential cross section for one-photon double
ionization of molecular hydrogen for a central photon energy of 75 eV, using a fully ab initio,
nonperturbative approach to solve the time-dependent Schr¨odinger equation in prolate spheroidal
1
coordinates. The spatial coordinates ξ and η are discretized in a finite-element discrete-variable
1
representation. The wave packet of the laser-driven two-electron system is propagated in time
0
through an effective short iterative Lanczos method to simulate the double ionization of the hy-
2
drogen molecule. For both symmetric and asymmetric energy sharing, the present results agree
n to a satisfactory level with most earlier predictions for the absolute magnitude and the shape of
a the angular distributions. A notable exception, however, concerns the predictions of the recent
J
time-independent calculations based on the exterior complex scaling method in prolate spheroidal
8 coordinates [Phys. Rev. A 82, 023423 (2010)]. Extensive tests of the numerical implementation
were performed, including the effect of truncating the Neumann expansion for the dielectronic in-
] teraction on thedescription of theinitial bound stateand thepredicted cross sections. Weobserve
h
that the dominant escape mode of the two photoelectrons dramatically depends upon the energy
p
sharing. Intheparallelgeometry,whentheejectedelectronsarecollected alongthedirectionofthe
-
m laserpolarizationaxis,back-to-backescapeisthedominantchannelforstronglyasymmetricenergy
sharing, while it is completely forbidden if thetwo electrons share the excess energy equally.
o
t
a PACSnumbers: 33.80.-b,33.80.Wz,31.15.A-
.
s
c
i I. INTRODUCTION center, multi-electron system. The single-center con-
s
y vergent close-coupling method was used to model the
h double ionization of H2 by Kheifets and Bray [8, 9].
Ameasurementofthe completebreakupofthe atomic
p Later McCurdy, Rescigno, Mart´ın and their collabora-
[ heliumtargetbyxuvradiationwasachievedover10years tors [10–12] implemented a formulation based on time-
ago [1]. Since then, rapid developments in strong xuv independent exterior complex scaling (ECS) in spheri-
1
light sources and momentum imaging techniques have
v calcoordinates,with the originof the coordinatesystem
made it possible to record all the reaction fragments,
7 placed at the center of the molecule, to treat the double
0 nuclei and electrons, in double photoionization of the photoionization at a photon energy of 75 eV. The radial
6 simplest two-electron hydrogen/deuterium molecule by coordinates of the two electrons are measured from the
1 one-photon absorption [2–5], and most recently also for
center, and the radial parts of the wave function were
. two-photon absorption [6]. For the double ionization
1 either expanded in B-splines or using a finite-element
0 of H2 by single-photon absorption, only randomly ori- discrete-variable representation (FE-DVR).
1 ented molecules were investigated in earlier experiments
The time-dependent close-coupling (TDCC)
1 (e.g. [7]). Using “fixed-in-space” molecules, more recent
method [13], again in spherical coordinates, was
: experimentaleffortsincludethemeasurementsofenergy-
v also extended to calculate the triple-differential cross
i andangle-resolveddifferentialcrosssectionsbyWeber et
X al. for either equal energy sharing [2, 3] or asymmetric section (TDCS) for double photoionization of the H2
molecule. While the agreement between the published
r energysharing[4],andbyGisselbrechtetal. [5]forequal
a TDCSs from the ECS [11, 12] and TDCC [13] calcu-
energy sharing at a photon energy of 76 eV. These ex-
lations is basically acceptable, noticeable discrepancies
perimental studies were at least partially stimulated by
remained for a few particular geometries. In the parallel
the goalofunderstanding the similarities anddifferences
geometry, for instance, where the molecular axis ζ is
between the hydrogen molecule and its atomic counter-
chosenalongthelaserpolarizationvectorǫ,thecoplanar
parthelium. However,allrecordedfullydifferentialcross
TDCSpredictionsfromtheECSandTDCCcalculations
sectionstodatesufferfromsomeexperimentaluncertain-
differ by up to 40 percent when one of the electrons (we
ties regarding the alignment angle of the molecule with
will refer to it as the “fixed electron” below) is observed
respect to the polarization vector of the laser and the
along the direction perpendicular to the ǫ-axis. In some
emission angles of the photoelectrons.
other cases, there exists a noticeable “wing” structure
From a theoretical point of view, the hydrogen in the published TDCC predictions for equal energy
molecule exhibits a significant complexity compared to sharing. Additional TDCC calculations [14] suggest
helium and, therefore, provides an enormous challenge that the agreement can be systematically improved,
to a fully ab initio description inherent in a multi- albeit the above-mentioned discrepancy still exists at a
2
somewhat reduced level. II. THE SCHRO¨DINGER EQUATION IN
Another independent approach [15] to this problem is PROLATE SPHEROIDAL COORDINATES
the very recent time-independent ECS treatment, for-
mulated – as in the current work – in prolate spheroidal The prolate spheroidal coordinates with the two foci
coordinates. Quitesurprisingly,the resultsofthatcalcu- separated by a distance R are defined by
lation differed from both the earlier ECS and also the
r +r r r
TDCC predictions, both obtained in spherical coordi- 1 2 1 2
ξ = , η = − , (1)
nates. Specifically, the ECS prolate spheroidal calcula- R R
tions showed differences from the earlier spherical coor-
and the azimuthal angle ϕ. Here r and r are the dis-
dinate calculation for the TDCSs, at a level of about 1 2
tancesmeasuredfromthetwonuclei,respectively. These
20% depending on the details of the energy sharing. As
coordinates are specified in the ranges ξ [1,+ ),
will be demonstrated below, we have gone to consider- ∈ ∞
η [ 1,+1], and ϕ [0,2π]. The volume element is
ablelengthsinanattempttoresolvethesediscrepancies. dV∈=−(R/2)3(ξ2 η2)∈dξdηdϕ. According to the asymp-
However, significant differences between the present re- −
totic behaviors as r ,r + , ξ and η approach 2r/R
sults and those of the ECS [15] still remain. 1 2 → ∞
and cosθ, respectively, where r and θ are the standard
Both the ECS [11, 12] and the TDCC [13] calcula-
spherical coordinates. Consequently, ξ is the “quasi-
tions made some attempt to deal with the experimental
radial” coordinate, while η is “quasi-angular”. The
uncertainties in the scattering angles. Given the experi-
Hamiltonian of a single electron, , (q = 1,2 for the
q
mental uncertainties and the differences in the previous H
two electrons in H below), is written as
2
calculations, however, it appeared worthwhile to inves-
tigate the computational effort required to obtain accu-
rate TDCSs before averaging over any experimental ac-
2 ∂ ∂ ∂ ∂
creepprteasnecnetsanagnleisn.deCpoenndseeqnuteinmtlpyl,emtheentparteisoenntofcatlhceultaitmioen- Hq =− R2(ξq2−ηq2)"∂ξq(ξq2−1)∂ξq + ∂ηq(1−ηq2)∂ηq
dependentFE-DVRapproachinprolatespheroidalcoor- (2)
dinates. Asintheotherapproachesmentionedabove,the 1 ∂2 1 ∂2 4ξ
q
internuclear separation (R) was held fixed at its equilib- + + .
riumdistanceof1.4bohr. Thetwo-centerprolatesystem, (ξq2−1)∂ϕ2q (1−ηq2)∂ϕ2q#− R(ξq2−ηq2)
with the foci located on the nuclei, provides a suitable
We solve the time-dependent Schr¨odinger equation
description for the two-center characteristics of the H
2
(TDSE) of the laser-driven H molecule (with two elec-
molecule. The formulation of the Schr¨odinger equation 2
trons) in the dipole length gauge:
in prolate spheroidal coordinates for diatomic molecules
is not new. The pioneering work of Bates, O¨pik, and
∂ 1
Poots [16] for the H+ ion, which is exactly separable in i Ψ(1,2,t)= + + +E(t) (r +r ) Ψ(1,2,t).
2 ∂t H1 H2 r · 1 2
prolate spheroidal coordinates, already revealed the ap- 12
h i (3)
pealing features of the prolate system. In particular,the
Here r is the coordinate of the q-th electron measured
electron-nuclearinteractionisrenderedbenigninthisco- q
relative to the center of the molecule and r = r r
ordinate system. A partial list of recent applications of 12 | 1− 2|
is the interelectronic distance. Without loss of general-
prolatespheroidalcoordinatestodiatomicmoleculescan
ity, we choose the molecular axis along the z axis, and
be found in [17–21].
the plane formed by the molecular axis and the polar-
As has been demonstrated in a number of recent pub-
ization vector as the xz plane. Generally, we can de-
lications, a grid-based approach provides a very ap-
composethepolarizationvectorintoitstwocomponents,
propriate description of laser-driven atomic and molec-
ǫ = cosθ e +sinθ e , where θ is the angle between
N z N x N
ular physics when combined with an efficient time-
theς andǫaxes,ande ande aretheunitvectorsalong
z x
propagation method such as the short iterative Lanczos
the z and x axes, respectively. The dipole interaction is
(SIL) method [22, 23]. In the present work, we employ
therefore given as
the FE-DVR/SIL approachin prolate spheroidal coordi-
nates to study the correlated response of a two-electron
E(t) (r +r )=E(t) (z +z )cosθ
molecule in the double ionization process. · 1 2 1 2 N
The remainder of this manuscript is organized as fol- h
+(x +x )sinθ , (4)
1 2 N
lows. After presenting the Hamiltonian of the hydrogen
i
molecule in Section II and providing some details about
where the rectangular coordinates x and z and the pro-
the discretization of the system in an FE-DVR basis in
late spheroidal coordinates are related through
Sec. III, the solution of the two-center Poisson equation
is presented in Sec. IV. This is followedby a description R R
of the procedure for extracting the cross sections of in- x= (ξ2 1)(1 η2)cosϕ, z = ξη. (5)
2 − − 2
terestin Sec.V. The results arepresentedanddiscussed
p
in Sec. VI, before we finish with a summary in Sec. VII. We expand the wave function for the H molecule in the
2
3
body-frame as spheroidal coordinates, we have chosen to work directly
in the FE-DVR basis. This differs from what is usually
Ψ(1,2,t)= Πm1m2(ξ1,η1,ξ2,η2,t) done for atoms in spherical coordinates, where spheri-
mX1m2 cal harmonics are used for the angular variables and an
Φ (ϕ ,ϕ ). (6) FE-DVR for the radial coordinates. Consequently, the
× m1m2 1 2
boundary conditions in prolate spheroidal coordinates
HereΦm1m2(ϕ1,ϕ2)=ei(m1ϕ1+m2ϕ2)/(2π)istheangular require some more discussion. Analyzing the asymp-
function,wherem1andm2denotethemagneticquantum totics reveals that in the region near the boundaries of
numbers of the two electrons along the molecular axis. ξ = 1 and η = 1, which correspond to the molecu-
Next, Πm1m2(ξ1,η1,ξ2,η2,t) is expanded in a prod- lar axis, the singl±e-electron wave function behaves like
uctofnormalized“radial” fi(ξ) and“angular” gk(η) (ξ2 1)|m|/2(1 η2)|m|/2. This indicates that the physi-
DVR bases: { } { } cal −wavefunctio−n is finite for m =0, whereas it goes to
| |
Π (ξ ,η ,ξ ,η ,t)= zerointhe regionclosetothe molecularaxisfor m =0.
m1m2 1 1 2 2 | |6
More importantly, the behavior of the wave function for
fi(ξ1)fj(ξ2)gk(η1)gℓ(η2)Cimjk1ℓm2(t). (7) odd m containsa square-rootfactor,givinga decidedly
| |
ijkℓ nonpolynomialbehaviorto the wavefunction that is im-
X
possible to capture in a straightforward fashion using a
Note that the basis is not symmetrized with respect to
DVR basis.
the coordinates of the two electrons. Since we always
The former problem is readily treated by using a
beginthecalculationwithaproperlysymmetrizedinitial
Gauss-Radauquadraturein the first DVR element for ξ,
state, however, that symmetry will be preserved in the
where only the right-most point is constrained to lie on
calculation.
theboundarybetweenthefirstandsecondfiniteelement.
To discretize this partial differential equation, we em-
The volume element ensures that the integrand is well
ploy the FE-DVR approach for both the ξ and η vari-
behaved near the end points and makes it unnecessary
ables [21, 24]. If we normalize the DVR bases according
to invoke a separate quadrature for different m values.
to
For all the other elements, a Gauss-Lobatto quadrature
dξf (ξ)f (ξ)=δ and dηg (η)g (η)=δ , is employed. This allows us to make the FE-DVR basis
i i′ ii′ k k′ kk′
continuous everywhere and to satisfy the m-dependent
Z Z (8) | |
boundary condition.
respectively, the overall (ξ,η) DVR basis is not normal-
To overcome the nonanalytic behavior of the basis for
ized with respect to the volume element in the prolate
odd m, Bachauand collaborators[17] explicitly factored
coordinatesystem. Thisiscorrectedbydefiningthetwo- out the (ξ2 1)|m|/2(1 η2)|m|/2 part before the wave
electron basis − −
function was expanded in terms of B-splines in the dis-
2 3 1 cretization approach. We have adopted a similar idea
bm1m2(1,2)= (9)
ijkℓ R (ξ2 η2)(ξ2 η2) in our FE-DVR treatment of the two-center problem to
(cid:16) (cid:17) i − k j − ℓ achieve much faster convergence, as was also done in
f (ξ )fq(ξ )g (η )g (η )Φ (ϕ ,ϕ ), Ref.[21]. Forthecaseofeven m,nochangesneedtobe
× i 1 j 2 k 1 ℓ 2 m1m2 1 2 | |
made to define the DVR basis, i.e., the normalized basis
which satisfies the desired normalization is written as
dV1dV2bmijk1ℓm2∗(1,2)bmi′j′1′mk′′2ℓ′(1,2)= f (ξ)= 1 ξ−ξk and g (η)= 1 η−ηk .
i i
ZZ ωi ξi ξk ωi ηi ηk
δii′δjj′δkk′δℓℓ′δm1m′1δm2m′2, (10) ξ kY6=i − η kY6=i −
q q (12)
to expand Π (ξ ,η ,ξ ,η ,t). Specifically, we have
m1m2 1 1 2 2 For odd m, however, we define the DVR basis as
| |
Ψ(1,2,t)= bm1m2(1,2)Xm1m2(t). (11)
ijkℓ ijkℓ 1 (ξ2 1)1/2 ξ ξ
k
Introducing the mnXo1rmm2aXijlkizℓed (ξ,η) DVR basis eliminates fi(ξ)= ωξi (ξi2−−1)1/2 kY6=iξi−−ξk (13)
thecomplexitiesofmatrixoperationsrelatedtotheover- q
and
lap matrix at each time, and hence makes the standard
SIL algorithm directly applicable to study the temporal 1 (1 η2)1/2 η η
k
response of the molecule to laser pulses. gi(η)= ωηi (1−−ηi2)1/2 kY6=iηi−−ηk. (14)
q
III. THE FE-DVR BASIS Hereωi andωi aretheweightfactorsrelatedtotheDVR
ξ η
bases f (ξ) and g (η), respectively. The goal of using a
i i
In our current implementation of the FE-DVR ap- uniquesetofmeshpoints,whichare m-independent,to
| |
proach for the time-dependent wave function in prolate discretize the (ξ,η) coordinates has now been achieved
4
in this scheme. The same technique was employed in After integrating over ϕ and ϕ , the matrix element of
1 2
recent calculations of one- and two-photon double ion- 1/r can be written as
12
ization of H [15, 24]. In principle, it is also possible
2
ttoheinDtVroRdubcaesetshetofaccirtcoursmv(ξe2nt−t1h)e|md|i/ffi2(c1u−ltiηes2)r|mel|a/t2edinttoo ijkℓm1m2|r1 |i′j′k′ℓ′m′1m′2 (17)
12
the nonanalytic behavior near the boundary. However, D ∞ E 2
this results in an m-dependence of the DVR bases and = 1 ( 1)|m|(2l+1) (l−|m|)! i′j′k′ℓ(l),
| | a − (l+ m)! Iijkℓ
quadrature points. This, in turn, leads to a number of l=X|m| (cid:18) | | (cid:19)
unnecessary complications in the practical implementa-
tionofthecomputationalmethodology. Onemightargue where the selection rule m = m1 −m′1 = m′2 −m2 has
that an m-dependent discretization procedure could be been used. Hence m is uniquely determined for a given
usefulfo|ra|systeminwhichthemagneticquantumnum- pair of angular partial waves. Above we introduced the
bermisconserved. AnexampleistheH+ ioninexternal reduced (ξ,η) integral
2
magnetic fields along the molecular axis [25]. However,
that is not the situation in the current calculation. i′j′k′ℓ′(l) (18)
Iijkℓ
= ijkℓ P|m|(ξ )Q|m|(ξ )P|m|(η )P|m|(η ) i′j′k′ℓ′ .
l < l > l 1 l 2
IV. THE ELECTRON-ELECTRON COULOMB
INTERACTION IN PROLATE SPHEROIDAL Since(cid:10) m (cid:12)(cid:12)is fixedinthe aboveequation,we om(cid:12)(cid:12)ittedit(cid:11)in
COORDINATES i′j′k′ℓ|′(l|) and will do so in the related quantities below
Iijkℓ
as well. It is now worthwhile to define the two electron
Similar tothe expansionofthe electron-electroninter- densities [29]:
action in terms of spherical coordinates, a counterpart
ρ (ξ,η)=f (ξ)g (η)f (ξ)g (η),
exists in prolate spheroidal coordinates [26] through the A i k i′ k′
(19)
Neumann expansion ρB(ξ,η)=fj(ξ)gℓ(η)fj′(ξ)gℓ′(η).
1 1 ∞ l (l m)! 2 After truncating the radial integral to the edge of the
= ( 1)|m|(2l+1) −| | (15) box, ξ , this yields
r a − (l+ m)! max
12 l=0m=−l (cid:18) | | (cid:19)
X X
×Pl|m|(ξ<)Q|lm|(ξ>)Pl|m|(η1)Pl|m|(η2)eim(ϕ1−ϕ2), Iiij′jk′ℓk′ℓ′(l)= dVξdVξ′ρB(ξ,η)Pl|m|(ξ<)Ql|m|(ξ>)
ZZ
awlhoenrgetahe≡zRa/x2is. aTnhdeξ>tw(<o)n=uclmeiaaxr(emilnoc)(aξt1e,dξ2a)t. ±BRo/th2 ×Pl|m|(η)Pl|m|(η′)ρA(ξ′,η′)
ittnihovenorslev[ge2ud7l]ai,nrwPthhl|eimche|(xξap)raeannsddieofiinrnraeesdgutihlnaert“Qhrea|lmdri|ea(gξli”)onLpae(rg1te,,n+wd∞rheil)fe,utanhrcee- =Z1ξmaxdVξPl|m|(η)ρB(ξ,η)Ul(ξ). (20)
“angular” part is only related to P|m|(η). Note that
l Here a convention for the volume element was made in
we chose to work in terms of an un-normalized “an-
such a way that, for any function F(ξ,η), we define
gular” basis rather than the usual spherical harmonics.
dV F(ξ,η) dξa3 +1(ξ2 η2)F(ξ,η)dη to simplify the
The matrix elements of 1/r12 in a traditional basis, for ξ ≡ −1 −
example, a B-spline or Slater-type basis, can be com- notation. MostimpRortantly,thefunctionUl(ξ)isdefined
putedthroughthewell-knownMehler-Ruedenbergtrans- by
formation [19, 28]. Due to the discontinuous derivative
ξ
along the line of ξ1 = ξ2 in the Neumann expansion, (ξ)=Q|m|(ξ) dV ρ (ξ′,η′)P|m|(ξ′)P|m|(η′) (21)
the straightforward computation of the matrix element Ul l ξ′ A l l
Z1
of 1/r12, using the value of this interaction potential at ξmax
the mesh points, is very slowly convergent. We seek a +Pl|m|(ξ) dVξ′ρA(ξ′,η′)Q|lm|(ξ′)Pl|m|(η′).
more robust representation of the 1/r12 matrix which Zξ
retains both the underlying Gauss quadrature and the
Instead of evaluating the above integrals directly, we
DVR property of all potentials being exactly diagonal
solve the differential equation satisfied by (ξ). As will
l
with respect to the highly localized DVR basis. U
becomeapparentlater,this equationcanbe shownto be
In the following, we use the simplified notation
the “radial” Poisson equation in the prolate spheroidal
ijkℓm m = f (ξ )f (ξ )g (η )g (η )Φ (ϕ ,ϕ )
| 1 2i | i 1 j 2 k 1 ℓ 2 m1m2 1 2 i coordinate system. The differential equations satisfied
to denote the basis. Essentially, we need the integral
by the Legendre functions P|m|(ξ) and Q|m|(ξ) suggests
l l
ijkℓm m P|m|(ξ )Q|m|(ξ )P|m|(η )P|m|(η ) (16) that we introduce the operator
1 2 l < l > l 1 l 2
D (cid:12) d d m2
(cid:12) eim(ϕ1−ϕ2) i′k′j′ℓ′m′m′ . 2 = (ξ2 1) l(l+1) (22)
×(cid:12) 1 2 ∇ξ dξ − dξ − − ξ2 1
(cid:12) E −
(cid:12)
(cid:12)
5
for given quantum numbers l and m. This is the one- at ξ =ξ , respectively.
max
dimensional Laplacian operator in the ξ coordinate. There are two important points to realize, namely:
After some algebra, we obtain First,onthe right-handboundary,the function (ξ) as-
l
U
sumesanonzerovalue,whichisgivenbyEq.(30),forall
d dQ|m|(ξ) ξ possible m values. Itsasymptoticbehaviorreliesonthe
dξUl(ξ)= ldξ dVtρA(t,τ)Pl|m|(t)Pl|m|(τ) (23) function|Q||m|(ξ) at large ξ, which behaves like 1/ξl+1.
Z1 l max
dP|m|(ξ) ξmax This indicates that it is nonzero generally, although it
+ l dV ρ (t,τ)Q|m|(t)P|m|(τ). could be small at the large ξ values used in practical
dξ t A l l max
Zξ calculations. Second, on the left-hand boundary, the sit-
uationdependsonthevalueof m. (ξ)takesanonzero
and l
| | U
value if m =0, while it becomes zero if m =0.
d2 d2Q|m|(ξ) ξ Follow|in|g the philosophy employed |to|h6 andle the
dξ2Ul(ξ)= dlξ2 Z1 dVtρA(t,τ)Pl|m|(t)Pl|m|(τ) Pspohisesroicnaleqcausaeti[o2n9s],(w26e)fi-(r3s0t)setehkatassaotliustfiieosn,thUel0(zξe)r,ot-voatluhee
+d2Pl|m|(ξ) ξdV ρ (t,τ)Q|m|(t)P|m|(τ)+ boundary condition at ξmax by using exactly the same
dξ2 t A l l ξ mesh points as those for the wave functions. In other
Z1
+W(Pl|m|(ξ),Q|lm|(ξ))a3Z−+11dτ(ξ2−η2)ρA(ξ,τ)Pl|m(|2(4τ)). wUUllo(00(r)ξd(msξ,a)xw=)e=haµ0vc.eµ∇fAµ2ξf(Utξel0)r(oξsf)utb=hset̺ist(ouξl)tuintwigoitnhthinUetl0oD(1tVh)Re=deUixfflpe(1rae)nnsatinioadnl
equation, we obtain a system of linear equations for the
Here the Wronskian of the Legendre functions P|m|(ξ) unknown cPoefficients cµ :
l { }
and Q|m|(ξ) is given by
l (l+ m)! 1
W(P|m|(ξ),Q|m|(ξ))= (−1)|m| (l+|m|)!. (25) Xµ′ cµ′Tµ|mµ′| =(−1)|m|(l−||m||)! ωξia3δµiδii′δkk′ (31)
l l (1 ξ2)(l m)!
− −| | (ξ2 η2)P|m|(η )q.
× i − k l k
Consequently, we obtain
The matrix T is defined by its elements
2 (ξ)=̺(ξ). (26)
∇ξUl ξmax d2 d
Thisisthesecond-orderinhomogeneousPoissonequation Tµ|mµ′| =− dξfµ(ξ) (ξ2−1)dξ2 +2ξdξ (32)
satisfied by l(ξ) with the “source” term given by Z1 (cid:20)
U m2
l(l+1) f (ξ).
̺(ξ)=( 1)|m|+1(l+|m|)!a3 (27) − − ξ2−1(cid:21) µ′
− (l m)!
−| | Therefore, the coefficient c can formally be written as
+1 µ
dη′(ξ2 η′2)ρ (ξ,η′)P|m|(η′).
After carryin×gZo−u1tthe inte−gralovAerη′ vialGaussquadra- cµ = [T|mω|]i−µi1(−1)|m|((ll+|mm|))!!a3δii′δkk′(ξi2−ηk2)Pl|m|(ηk),
ξ −| |
ture, we recast the source term as
q (33)
(l+ m)! where [T|m|]−1 denotes the inverse of the matrix T|m|.
̺(ξ)=δkk′( 1)|m|+1 | | a3fi(ξ)fi′(ξ) In this case 0(ξ) fulfills the left-hand boundary condi-
− (l−|m|)! (28) tion 0(1)=U0l, and so does (1). Recall, however, that
×(ξ2−ηk2)Pl|m|(ηk). its riUghlt-hand boundary conUdlition differs from those of
(ξ ). Thissuggeststhatthefinalanswertothefunc-
l max
As one might expect, the inhomogeneous equation re- tUion (ξ) can be constructed as (ξ) = 0(ξ)+F (ξ),
tdouctehsetPootishseonhoemquoagteinoenosu(s2o6)neanifdk(26=8)kc′a.nTbheeusnoilquutieolny i.e., wUel add the difference functioUnlFl(ξ) tUol Ul0(ξ). lThe
function F (ξ) is also a solution to the homogeneous
l
determined by enforcing the boundary conditions
Poisson equation, subject to the boundary condition
F (1)=0 and F(ξ )= (ξ ). After writing it as a
Ul(1)=Pl|m|(1) ξmaxdVξ′ρA(ξ′,η′)Q|lm|(ξ′)Pl|m|(η′) linlearcombinatilonmoafxPl|m|U(ξl)amnadx Q|lm|(ξ),andimposing
Z1 the boundary conditions, F (ξ) takes the form
(29) l
at ξ =1 and
F(ξ)=δ δ a3(ξ2 η2)P|m|(ξ )P|m|(η )
l ii′ kk′ i − k l i l k
ξmax
Ul(ξmax)=Q|lm|(ξmax) dVξ′ρA(ξ′,η)Pl|m|(ξ′)Pl|m|(η′) Q|lm|(ξmax)P|m|(ξ). (34)
Z1 (30) ×Pl|m|(ξmax) l
6
We finally arrive at m =0. Substituting Eq. (35) into Eq. (20) allowsus to
| |
obtain the kernel integral,
( 1)|m|(l+ m)!
Ul(ξ)= − ωξi (l−||m||)!a3δii′δkk′(ξi2−ηk2)Pl|m|(ηk) Iiij′jk′ℓk′ℓ′(l)=δii′δjj′δkk′δℓℓ′a6(ξi2−ηk2)(ξj2−ηℓ2)Pl|m|(ηℓ)
×Xµq[T|m|]−µi1fµ(ξ)+δii′δkk′(ξi2−ηk2)Pl|m|(ξi) ×"(−ω1)ξi|ωmξj|((ll−+||mm||))!![T|m|]−ji1Pl|m|(ηk)+
×Pl|m|(ηk)QPl||lmm||((ξξmmaaxx))Pl|m|(ξ). +Plq|m|(ξi)Pl|m|(ξj)Pl|m|(ηk)QP||lmm||((ξξmax))#.
(35) l max
(36)
At this point, the DVR version of the solution (ξ) is
l
U
ready for all possible values of m, either m = 0 or The matrix element of 1/r can finally be written as
12
| | | | 6
1 lmax (l m)!
ijkℓm1m2 r i′j′k′ℓ′m′1m′2 =δii′δjj′δkk′δℓℓ′a5(ξi2−ηk2)(ξj2−ηℓ2) (2l+1)(l+−|m|)!Pl|m|(ηk)Pl|m|(ηℓ)
12
(cid:10) (cid:12) (cid:12) (cid:11) l>X|m| | | (37)
(cid:12) (cid:12) 1 (l m)! Q|m|(ξ )
[T|m|]−1+( 1)|m| −| | P|m|(ξ )P|m|(ξ ) l max ,
×" ωξiωξj ji − (l+|m|)! l i l j Pl|m|(ξmax)#
q
where we truncated the l summation in the Neumann maximumvaluesof m and m ,determinesthelargest
1 2
| | | |
expansionto l . The above equationcanbe converted time step ∆t for the propagationinrealtime. Typically,
max
to the normalized (ξ,η) bases with the help of Eq. (9). E is about 6,000 atomic units (a.u.) in our calcu-
max
Thisresultsinadiagonalrepresentationofthematrixel- lations. Although the chances of electrons populating
ements of the electron-electronCoulomb interaction and states with such high energies are practically negligible
thus considerably simplifies the FE-DVR discretization forshorttime scalesofthelaser-moleculeinteraction,we
procedure. Theabovetreatmentofthe 1/r matrixwas generally require ∆t . 2π/E in order to resolve the
12 max
successfully applied to the two-photon double ionization mostrapidoscillationsinthetimeevolution. Thismeans
of H [24]. The implementation of this representation that at least a few time steps are needed during one pe-
2
will be illustrated below. riod of 2π/E . We refer readers to Refs. [30–32] for
max
further details and discussions behind the SIL method.
In order to ensure that the double-ionization wave
V. TIME EVOLUTION AND EXTRACTION OF
packet is sufficiently far away from the nuclei, and also
CROSS SECTIONS
that the two photoelectrons are well separated, we al-
low the system to evolve for a few more cycles in the
The time-dependent laser-driven electronic wave field-free Hamiltonian, i.e., after the laser pulse has died
packet in the hydrogen molecule is obtained by solving off. This is the wave packet we use to extract the phys-
the TDSE on the (ξ,η) grid. Launched from the pre- ical information. The ionization probabilities and the
viously determined ground state, the time evolution of corresponding cross sections are extracted by projecting
the system is achieved by using our recently developed the time-dependent wave packet onto uncorrelated two-
SIL method [30, 31]. The ground state is determined electroncontinuumstatessatisfyingthe standardincom-
by relaxing the system in imaginary time from an ini- ingboundaryconditions. ThelatterstatesofH arecon-
2
tial guess of the wave function on the grid points. At structedfromtheone-electroncontinuumstateoftheH+
2
eachtimestepweonlyneedtogeneratethevaluesofthe ion described in the following subsection.
discretized wave function on the selected grid points. If
desired, the information at arbitrary points within the
spatial box can be obtained from the interpolation pro-
A. Continuum states of H+
cedure in terms of the DVR bases. 2
Afewremarksseemappropriateregardingtheefficient
implementation of the SIL algorithm. The highest en- The field-free wave function Φ(ξ,η,ϕ) of the one-
ergy, E , which essentially depends on the smallest electron molecular ion is completely separable in pro-
max
separationbetweenthe(ξ,η)meshpointsandalsoonthe late spheroidal coordinates. For a given, and conserved,
7
magnetic quantum number m, the wave function takes satisfying the incoming boundary condition can be ex-
the form Φ(ξ,η,ϕ) = T (ξ)Ξ (η)Φ (ϕ), where the az- panded as
m m m
imuthal dependence, is given by Φ (ϕ) eimϕ/√2π.
The“radial”partTmq(ξ)andthe“anmgular≡”partΞmq(η) Φ(−)(r)=1 +∞ iℓe−i∆mq(k)
k
of the wave function satisfy the equations k
m=−∞ℓ>|m|
X X
∂ (ξ2 1) ∂ m2 +2Rξ+c2ξ2 A T (ξ)=0 ×Yℓ∗m(k)Yℓm(ηr,ϕr)Tm(kq)(ξ). (43)
∂ξ − ∂ξ −(ξ2 1) − mq mq
(cid:20) − (cid:21) This function is normalized in momentum space accord-
(38)
ing to Φ(−) Φ(−) =δ(k k′), provided the asymptotic
and h k | k′ i −
solution in Eq. (41) is satisfied.
∂ ∂ m2 Uncorrelated two-electron continuum states with to-
∂η(1−η2)∂η − (1 η2) −c2η2+Amq Ξmq(η)=0, tal spin angular momentum S (S = 0 in our case) can
(cid:20) − (cid:21) generally be constructed as follows:
(39)
respectively. Here c = kR/2 for the continuum state Φ(−) (r ,r )= (44)
whose momentum vector has the magnitude k. In addi- k1k2 1 2
1
tion, we need to introduce another quantum number q, Φ(−)(r )Φ(−)(r )+( 1)SΦ(−)(r )Φ(−)(r ) .
which denotes the number of nodes of Ξ (η) in the re- √2 k1 1 k2 2 − k2 1 k1 2
m
h i
gion η [ 1,+1], to label the states, and finally the
∈ − WiththehelpofEq.(43),itspartial-waverepresentation
separation constant A .
mq can be written as
When the angularfunction Ξ (η)Φ (ϕ) is discretized
m m
in terms ofthe relevantDVR bases,a few “spurious”so- Φ(−) (r ,r )= 2 3 1 1 iℓ1+ℓ2
lutionsmightbeencountered. Thisiscausedbytheresid- k1k2 1 2 R √2k1k2
ualerrorsassociatedwiththeGaussquadratures. Conse- (cid:16) (cid:17) ℓ1mX1ℓ2m2
quently, we expand the angular function, or “spheroidal bm1m2(1,2) (ξ2 η2)(ξ2 η2) (45)
× ijkℓ i − k j − ℓ
harmonics” function Yℓm(η,ϕ) ≡ Ξmq(η)Φmq(ϕ) with Xijkℓ q
ℓ= m+qinsteadintermsofsphericalharmonics. These
func|tio|ns are normalized according to ×"e−i ∆|m1|ℓ1(k1)+∆|m|2ℓ2(k2) Yℓ∗1m1(k1)Yℓ∗2m2(k2)
(cid:0) (cid:1)
+1 2π
dη dϕ ∗ (η,ϕ) (η,ϕ)=δ δ . (40)
Z−1 Z0 Yℓm Yℓ′m′ mm′ ℓℓ′ Ciℓj1kmℓ1ℓ2m2(k1,k2)+(−1)S(k1 ↔k2)#.
After obtaining the separation constant A by solv-
mq Here we introduced
ing Eq. (39), the “radial” function T (ξ) is again ex-
mq
panded in terms of the DVR bases. The last DVR point Cℓ1m1ℓ2m2(k ,k )= (46)
ijkℓ 1 2
at ξ = ξ needs to be kept for the continuum state.
max T˜(k1) (ξ )T˜(k2) (ξ )Ξ˜(k1) (η )Ξ˜(k2) (η ),
Asymptotically, the radial function behaves like ℓ1|m1| i ℓ2|m2| j ℓ1|m1| k ℓ2|m2| ℓ
byrepresentingtheradialandangularpartsonthe(ξ,η)
1 8 R ℓπ
T (ξ) sin cξ+ ln(2cξ) +∆ (k) grid points:
mq mq
→ ξR π c − 2
r (cid:20) (4(cid:21)1) T(k)(ξ)= f (ξ)T˜(k)(ξ ), (47)
ℓm i ℓm i
as ξ + . Here ∆ (k) is the two-center Coulomb
→ ∞ mq Xi
phase shift. The normalization factor on either the en-
ergyorthemomentumscaleandtheCoulombphaseshift Ξ(k)(η)= g (η)Ξ˜(k)(η ). (48)
can be determined by matching the numerical solution ℓm µ ℓm µ
µ
of T (ξ) according to its asymptotic behavior given in X
mq
Eq. (41). Here the exchange symmetry
The plane wave in prolate spheroidal coordinates can
Cℓ2m2ℓ1m1(k ,k )=Cℓ1m1ℓ2m2(k ,k ) (49)
be written as [33] jiℓk 2 1 ijkℓ 1 2
is satisfied.
eik·r =4π iℓ (η ,ϕ ) ∗ (η ,ϕ )R(k)(ξ), (42)
Yℓm r r Yℓm k k ℓm
ℓm
X
B. Extraction of double-ionization cross sections
where R(k)(ξ) 1/(cξ)sin cξ ℓπ/2 in the asymp-
ℓm → −
totic region. Note that η and η are related to the Ithasbeendemonstrated[24,30,34]thatusinguncor-
k r
(cid:2) (cid:3)
directions of k and r in spherical coordinates through related two-electroncontinuum states is a good approxi-
η = cosθ . The partial-wave expansion of the plane mation in a time-dependent propagation approach, pro-
k,r k,r
waveeik·r remindsusthatthetwo-centerCoulombwave vided the two ejected electrons are well separated from
8
each other. The probability amplitude of double ioniza- special case of the generalized N-photon effective inter-
tion is then given by action time T(N) [35] for a one-photon reaction. Gener-
eff
alized cross sections for two-photon double ionization of
hΦk(−1k)2|Ψ(t)i= (50) the hydrogen molecule were extracted in the same for-
1 malism [24].
( i)ℓ1+ℓ2ei ∆|m1|ℓ1(k1)+∆|m|2ℓ2(k2)
k k − In the TDCC treatment [13], different strategies were
1 2
m1Xℓ1m2ℓ2 (cid:0) (cid:1) employedto describe the linearone-photonandthe non-
×Yℓ1m1(k1)Yℓ1m2(k2)Fℓ1m1ℓ2m2(k1,k2), linear two-photon double ionization processes of atoms
and molecules. For the one-photon case, the cross sec-
where
tions were obtained through the time derivative of the
double-ionization probability, ∂P2+(t)/∂t. The laser
F (k ,k )= (51) ion
ℓ1m1ℓ2m2 1 2 field does not need to be turned off in this case. On the
√2 Cℓ1m1ℓ2m2∗(k ,k )Xm1m2(t). otherhand,atruelaserpulsewasusedforthetwo-photon
ijkℓ 1 2 ijkℓ
case and an effective time, defined as the time integral
ijkℓ
X
under a flat-top pulse with a smooth turn-on and turn-
Here the exchange symmetry off, was introduced. In the present work, we employed a
unified formulationthrough an effective interaction time
F (k ,k )=( 1)SF (k ,k ) (52)
ℓ2m2ℓ1m1 2 1 − ℓ1m1ℓ2m2 1 2 forbothone-photonandmulti-photonionizationinlaser
pulses.
is satisfied. We also see that the probability ampli-
For the one-photon double ionization initialized from
tude formulated in prolate spheroidal coordinates takes
thelowestX1Σ state,thetwoejectedelectronscanonly
a similar form as for the atomic case in spherical coordi- g
populate the final 1Σ and 1Π continuum states, with
nates. However,a subtle difference from the atomic case u u
the specifics depending on the relative orientationof the
is worth pointing out. Strictly speaking, the spheroidal
molecular axis and the laser polarization vector. Con-
harmonics involved in the probability amplitude gener-
sequently only partial waves with ungerade parity [i.e.,
ally depend on the magnitude of the momenta k and
1 ( 1)ℓ1+ℓ2 = 1] need to be included in Eq. (53).
k2, in addition to their directions. − −
The energy sharing of the two photoelectrons can
be specified by introducing the hyperangle α =
VI. RESULTS
tan−1(k /k ). This describes the double-ionization re-
2 1
action with kinetic energies E = E cos2α and E =
E sin2α for the two ionize1d elecetxrcons, respectiv2ely. A. Preparation of the initial electronic X1Σg state
exc
HereE istheavailableexcessenergyabovethedouble-
exc
ionization threshold. In the present work, specifically, For the nonsequential double-ionization process in-
E =23.6 eV for absorption of a 75-eV photon. duced by one- or two-photon absorption, electronic cor-
exc
For double ionization by one-photon absorption, the relationplaysadominantrole,asthetwophotoelectrons
triple differential cross section with respect to one of the mustsharethe availableexcessenergyEexc. Double ion-
kinetic energies and the two solid angles kˆ and kˆ can ization by a single photon would not occur at all if the
1 2
be extracted using the same formalism as in the corre- two-electronatom or molecule were approximated by an
sponding He case [30, 34], i.e., independent-electron model. Therefore, the quality of
thedescriptionofelectron-electroncorrelationinalaser-
d3σ 1 ω 1 kmax k1′ driven system is crucially important for accurate results
dE1dkˆ1dkˆ2 = k1k2cos2αI0Te(ff1) Z0 dk1′ Z0 dk2′ ttwoobeeleocbttraoinnsedh.asTthoebCeodueloscmrbibeindteinracaticoonnsbisettewnetenmtahne-
nerforboththe initialboundstateandthe time-evolved
×k1′δ(k2′ −k1′ tanα)(cid:12) (−i)ℓ1+ℓ2ei ∆|m1|ℓ1+∆|m|2ℓ2 wave packet. Before we go any further, it is worth dis-
(cid:12)m1Xℓ1m2ℓ2 (cid:0) (cid:1) cussing how we prepare the initial X1Σg state at the
(cid:12)
(cid:12) 2 equilibrium distance of R=1.4 bohr.
(k′,kˆ(cid:12)) (k′,kˆ )F (k′,k′) . As seen from Eq. (37), the magnetic quantum num-
×Yℓ1m1 1 1 Yℓ2m2 2 2 ℓ1m1ℓ2m2 1 2 (cid:12) ber in the Neumann expansion of the matrix element of
(cid:12)
((cid:12)53) 1/r isuniquelydeterminedbytheangularbases. How-
(cid:12) 12
(cid:12) ever, this is not the case for the index l, if we choose to
Here ω and I are the central photon energy and the
0 discretize the coordinate η, rather then expanding that
peak intensity of the laser pulse, respectively, while T(1) part of the wave function into spherical harmonics. In
eff
denotes the effective interaction time between the tem- practice, the summation over l must be truncated at a
poral electric laser field and the electrons in the one- finite value of l . In principle, the higher-order ex-
max
photon absorption process. For a laser pulse of time du- pansion terms always guarantee well-converged results.
rationτ with asine-squaredenvelopeforthe fieldampli- However, as mentioned earlier, we approximate the rel-
tude, T(1) = (3/8)τ. Note that T(1) corresponds to the evant η-integrals by using Gauss-Legendre quadrature.
eff eff
9
This is the price we have to pay for making the dielec-
tronicCoulombpotentialdiagonalintheDVRbases. As -1.886 (9,1)
)
a consequence, we need to determine how the approxi- u (11,1)
a
mationintroduced inthe η-integralsfor the two-electron ( (13,1)
y -1.887
integralsaffectstheresultsforthe crosssectionsofinter- g
r
est. ne
e
To answer this question, we first investigate the de- c -1.888
pendence of the energy obtained for the initial X1Σ ni
g o
state on the value of lmax used in the Neumann expan- ctr -1.889 (9,4)
sion. Figure 1 shows the variation of the initial-state Ele (11,4)
electronic energy of the hydrogen molecule with respect
E =−1.888761au (13,4)
to lmax, obtained with a ξ setup of ten elements in the -1.890 gs
region of 1 < ξ 6 15.82, five in the region 1 < ξ 6 5 0 10 20 30 40
and another five in 5 6 ξ 6 15.82. Each element, in
l
max
turn, contains five DVR points to further discretize the
configuration space. Furthermore, we employ 9th-order
DVR points for η. For a given number of η mesh points FIG. 1. (Color online) Energy of thelowest electronic X1Σg
state at R=1.4 bohr as a function of the lmax value used in
(n ) and m = m = m , we observe that
η | |max | 1|max | 2|max the Neumann expansion of 1/r12. The number of η points,
the resulting energy typically exhibits a plateau-like be-
nη, and the largest magnetic quantum number, |m|max, are
havior with increasing lmax. For given nη, when lmax labeled as (nη,|m|max). The open symbols correspond to
is relatively small, the η-integral can be computed very |m|max =1, while the filled symbols are for |m|max =4. The
accurately by using Gauss quadrature. However, when benchmark energy (Egs) from Ref. [36] is shown as well.
l is too large, the numerical errors introduced from
max
the Gauss quadrature cause the energy value to fluctu-
ate. Thisoccurswhenl approaches2n andisshown B. Convergence of the TDCS
max η
by the grey stripes in Fig. 1. In this region of l an
max
unphysically low energy can be produced. Beyond that Beforewe presentourresults forthe crosssections,let
point,thecalculatedenergyincreasestothenextplateau. us take a closer look at the survival probability
Ultimately, this is not too surprising, since any Gauss
P = Ψ Ψ(t) 2 (54)
quadrature is only reasonably accurate up to a limited surv |h gs| i|
polynomial order of the integrand. Consequently, if we of the aligned H molecule in its ground state Ψ . This
2 gs
want to keep more terms in the Neumann expansion, is shownin Fig.2. For homonuclearmolecules,the inde-
we have to increase nη correspondingly. This finding is pendent alignment angle θN between the molecular axis
further substantiated by the dependence of the energy andthe polarizationvectorcanbeconfinedtotheregion
found for nη = 11 and 13. The plateaus are indeed ex- from 0◦ to 90◦. In the xuv regime, we observe that the
tendedtothe correspondinglylargervaluesof2nη. Most hydrogen molecule shows a larger probability of being
importantly, the amplitude of the energy fluctuation is ionized or excited (i.e., a lower probability of staying in
systematically reduced with increasing nη. The error in the initial state) at the end of the pulse in an aligned
the energy is lowered from 2.13 10−3 to 1.45 10−3 geometry. This indicates that the perpendicular compo-
and finally 1.05 10−3 a.u., whe×n nη increases×from 9 nent of the temporal electric field exerts more influence
×
to 11 and then 13 for mmax = 4. We obtained the en- on the ionization process due to the larger dipole mo-
| |
ergy at R = 1.4 bohr as 1.8887324 a.u. for lmax = 10, mentum. Interestingly, at the earlier stages of the time
−
nη =9, and mmax =4, resulting in a double-ionization evolution(e.g.,t.9a.u.),whentheionizedwavepacket
| |
potentialof51.394eV.Keepingtheotherparametersun- is driven back by the change in direction of the electric
changed, we obtained an energy of 1.8887128 a.u. for field, the tilted molecule has a largerprobability of stay-
−
nη = 11. The benchmark energy in the literature is ing in its ground state. This happens near the various
1.888761428a.u. at the same R [36], after we take out minima in P . However,once the electric field has be-
− surv
the nucleus-nucleus interaction of 1/1.4 a.u. come sufficiently strong (t & 9 a.u.), the wave packet is
Tosummarize: Unlikeforotherexpansionparameters, driven out and spread into a larger space. This leads to
it is importantto be consistentin the size ofthe angular lower minima in P for the tilted molecule.
surv
quadrature and the largest l employed in the Neu- When the wave packet is driven back to the nuclear
max
mannexpansioninpracticalcalculations,ifwediscretize region and therefore has a chance to recombine with the
the coordinate η. However, this provides a way to ex- H+ ion, a maximum in P appears. Not surprisingly,
2 surv
amine a potential sensitivity of the physical observables the parallel geometry always has the largest probability
of interest (here the differential cross sections) to the for this to happen. Although the wave packet can also
ground-state wave functions generated by varying l be scattered for the untilted molecule in the plane per-
max
and other parameters. This will be further discussedbe- pendicular to the molecular axis, the probability is un-
low. doubtedlylargerifthelaserelectricfieldisperpendicular
10
to the molecular axis. A similar behavior of H+ in xuv
2 TABLEI.Thediscretizationandexpansionparametersofthe
pulses was observed in Ref. [37].
H2 wave function in prolate spheroidal coordinates. Here ξb
stands for the border between the inner and outer regions
1.000
in the ξ coordinate, while ξmax is the size of the ξ box. In
addition, n denotes the number of ξ mesh points in each
ξ
y element. The numbers of ξ elements in the inner and outer
t 0.998 (a)
bili region are ninn and nout, respectively. These ξ parameters
ba produce the total number of ξ mesh points Nξ. The ξ grid I
ro 0.996 and ξ grid II are used to examine the convergence of our
p results.
al θN = 0◦
viv 0.994 θN =30◦ ξb ninn nout ξmax nξ Nξ
ur θN =60◦ ξ grid I 5 5 67 150 5 288
S θN =90◦ ξ grid II 9 1 11 97 14 156
0.992
0 5 10 15
Time (au)
differential cross section, since it reveals the fine details
1.000 ofpossibleenergysharingsandpreferreddirectionsofthe
ejected electrons in the double-ionizationprocess. Given
y the discrepancies between results from various calcula-
t 0.998
bili tions found in the literature, we carried out comprehen-
a sive convergence tests for our predictions of the TDCSs.
b
ro 0.996 These tests are essentially divided into two groups. The
p
first group concerns the laser parameters, while the sec-
al
v ond one deals with the discretization and expansion pa-
vi 0.994
r rameters. Anexampleoftwodifferentparametersetsfor
u
S (b) theξ gridisgiveninTableIandwillbefurtherdiscussed
below.
0.992
In order to obtain a good handle on the sensitivity of
15 20 25
the results to the various parameters and the resulting
Time (au)
level of “convergence”, we try to only vary a single pa-
rameter while keeping all others fixed if possible. For
FIG. 2. (Color online) Survival probability of the hydrogen
the dependence on the laser parameters, we use the ξ
molecule subjected to a sine-squared laser pulse with a peak
grid I combined with (n , m ,l ) = (9,4,10). For
intensity of 1015 W/cm2. The laser pulse lasts for 10 optical η | |max max
the tests regarding the discretizations and expansions,
cyclesandthesystemisfollowedforaperiodofanother2cy-
the peak intensity of laser was fixed at 1015 W/cm2 and
cles of field-free propagation. The central photon energy of
thelaser pulse is 75 eV. a time scale of “10+2” optical cycles (o.c.) was used.
Here “10+2”refers to a10-cyclelaserpulse with a sine-
squared envelope for the field amplitude, followed by a
For most calculations performed in this study, we ex-
2-cycle field-free propagation.
pose the hydrogen molecule to a laser pulse with a peak
intensity of 1015 W/cm2. Looking at Fig. 2 we see that Figures 3 and 4 show the convergence pattern of our
the depletion of the initial ground state can be safely TDCS results for asymmetric energy sharing in the par-
neglected for our typical interaction times. Even for allel geometry (θN = 0◦). The energy sharing between
θN = 90◦, Psurv = 0.99677 remains very close to unity. electron1(observedatthe fixedangleθ1)andelectron2
The negligible depletion of the ground state suggests (observed at the variable angle θ2) is 20% : 80%. Only
that the concept of cross sections is valid and applica- the electron that takes away 20% of the excess energy is
ble. On the other hand, it also presents a numerical recordedat fixed positions either parallelor perpendicu-
challenge to predict the cross sections accurately from lar to the polarization axis.
a time-dependent treatment, due to the generally small Since the laser pulse is explicitly involvedin our time-
ionization probability. dependent treatment, we first have to be sure that the
At first glance, a peak intensity of 1015 W/cm2 might extractedcrosssectionsareessentiallyindependentofthe
seemveryintenseformostatomicandmoleculartargets. laserintensityandthetimescales. Onlythenarethecal-
Here, however, we consider an xuv rather than an IR culations of cross sections meaningful. This also allows
pulse. For an xuv pulse with central photon energy of us to compare the physical information extracted from
75eV,suchlaserfieldsdefinitelyfallintothe“weak-field” our time-dependent scenario to that obtained through
regime. The ponderomotive energy in the xuv regime is conventionaltime-independent treatments, which are ef-
much smaller than the photon energy of interest. fectively equivalent to the weak-field approximation and
In this work, we are mainly interested in the triple- “infinitely” longinteractiontimes. Rather thancomput-
