Home
Search
Collections
Journals
About
Contact us
My IOPscience
A classical model of H in an intense laser field
This content has been downloaded from IOPscience. Please scroll down to see the full text. 2013 J. Phys. B: At. Mol. Opt. Phys. 46 235601 (http://iopscience.iop.org/0953-4075/46/23/235601) View the table of contents for this issue, or go to the journal homepage for more
Download details: IP Address: 157.82.9.138 This content was downloaded on 08/11/2013 at 02:37
Please note that terms and conditions apply.
IOP PUBLISHING
JOURNAL OF PHYSICS B: ATOMIC, MOLECULAR AND OPTICAL PHYSICS
doi:10.1088/0953-4075/46/23/235601
J. Phys. B: At. Mol. Opt. Phys. 46 (2013) 235601 (9pp)
A classical model of H+ 3 in an intense laser field Erik L¨otstedt 1 , Tsuyoshi Kato and Kaoru Yamanouchi Department of Chemistry, School of Science, The University of Tokyo, 7-3-1 Hongo, Bunkyo-ku, Tokyo 113-0033, Japan E-mail:
[email protected]
Received 2 August 2013, in final form 1 October 2013 Published 7 November 2013 Online at stacks.iop.org/JPhysB/46/235601 Abstract Soft-core Coulomb potentials for the pairwise interactions among electrons and protons are + introduced so that H+ 3 , H2 , H2 and H can take stable geometrical configurations in classical mechanics, and responses of H+ 3 to an ultrashort intense laser pulse are investigated by integrating the classical equations of motion. The Monte Carlo simulation shows that the single and double ionization as well as the dissociation proceed, indicating that classical mechanics can give insight into dynamics of small molecular systems interacting with an intense laser field. (Some figures may appear in colour only in the online journal)
1. Introduction
the intensity dependence of the yield of nonsequential double ionization [17–20]. Classical calculations of systems in intense laser fields have been extended to small molecular systems such as H+ 2 [21, 22] and H2 [23–25]. Apart from purely classical models, there are also approaches in which quantum effects are included in an approximate way. One such class of methods is referred to as ‘fermionic molecular dynamics’ [26, 27], where the main idea is to add a momentum-dependent potential to the Coulomb interaction. The auxiliary potential, supposed to mimic the Heisenberg principle, is sufficient to render both many-electron atoms [28] and molecules [29] stable. This model was found to be helpful in the interpretation of experimental data when quantum mechanical calculations were not available, and has been applied to collision problems [29–31], and to laser-induced ionization of atoms [32–35] and molecules [36–39]. In this paper, we attempt to construct an entirely classical model of a triatomic H+ 3 molecule. In contrast to previous studies [37, 38], no quasiclassical potentials are employed. A fundamental difficulty that has to be overcome in modelling classical molecules is the stability of the ground state, having the equilibrium configuration of minimum energy. In the case of classical H+ 2 , where one electron and two nuclei interact with the Coulomb force, a few stable trajectories were found [4], but no stable particle trajectories are known for H2 or any other polyatomic molecules. A stable trajectory is defined in this
Classical models of atoms and molecules are interesting for a number of reasons. From a fundamental point of view, it should be worth investigating to what extent quantum systems can be understood in terms of classical mechanics. In particular, a classical helium atom, a fundamental three-body system, has been extensively studied in this perspective [1, 2]. Some efforts have also been made to understand classical, purely Coulombic H+ 2 [3, 4] and H2 [5–7]. Another motivation to study atomic and molecular systems in classical mechanics is rather practical. Since the classical equations of motion are straightforwardly integrated numerically, dynamics induced by a particle collision or by an external field can easily be studied. This kind of research was performed for the field ionization of a hydrogen atom [8, 9], and experimental results were found to be well reproduced. Also high-order harmonic spectra may be computed from classical trajectory ensembles [10]. In the relativistic regime of laser–matter interaction, even the numerical solution of the one-particle quantum wave equation (Dirac or Klein– Gordon) poses a considerable challenge [11–13], but classical simulations may provide practical alternatives [14–16]. A classical model of the helium atom convincingly explained 1
Present address: Laser Technology Laboratory, RIKEN, 2-1 Hirosawa, Wako, Saitama 351-0198, Japan.
0953-4075/13/235601+09$33.00
1
© 2013 IOP Publishing Ltd
Printed in the UK & the USA
J. Phys. B: At. Mol. Opt. Phys. 46 (2013) 235601
E L¨otstedt et al + having a stable H+ 3 , we also demand H2 and H2 to be stable, since we are interested in calculating dissociation probabilities, and we require H2+ 3 to be unstable, as in the real, quantum mechanical case. The H atom, with Hamiltonian p2 1 P2 HH = + − , (2) 2M 2 (R − r)2 + α 2
context as a trajectory that would not lead to autoionization or autodissociation of the molecule when all particles, both electrons and nuclei, are allowed to move. We show in this study that a suitable choice of softened Coulomb interaction potentials may be employed to create a stable H+ 3 molecule. This type of approach has been adopted before [23–25], but to the best of our knowledge, this study is the first showing that such a model may be made to work for polyatomic molecules. The choice of H+ 3 is partially inspired by recent experimental measurements [40–43]. It is expected that theoretical models like the one proposed in the current paper can provide useful feedback for interpretation of experimental data. Quantum mechanical studies of H+ 3 exposed to a pulsed laser field can be found in [44–46], where the protons were fixed, and in [47, 48] in which time-dependent density functional theory was adopted for describing electrons, and classical dynamics for nuclei. The paper is organized as follows. In section 2, we introduce our theoretical model of a classical polyatomic molecule, and examine its stability properties. After a brief description of the procedures of simulation in section 3, the results of a Monte Carlo simulation of the response of a model H+ 3 molecule to a few-cycle, intense laser pulse are presented in section 4. Conclusions are described in section 5, and a brief discussion of laser-induced dissociation of H2 and the effect of changing the model parameters can be found in the appendix.
ep
is stable for any value of αep , with the equilibrium configuration R0 −r0 = 0, P0 = p0 = 0. Note that by ‘stable’, we refer to the existence of a static equilibrium configuration, such that P j = pk = 0 for all particles, which is stable against small perturbations in all directions. For H+ 2 , which has the Hamiltonian 2 2 P2j 1 p2 + HH2 = + − 2M 2 2 (R j − r)2 + αep j=1 j=1 +
In this section, and throughout the paper, we work in atomic units (au). As mentioned in the introduction, the idea is to construct a stable H+ 3 molecule, consisting of two electrons and three protons, and follow the motion of all particles interacting with soft Coulomb potentials. The so-called softcore potential is a standard model potential of the electron– nucleus interaction in laser-driven atoms, both in quantum calculations in reduced dimensionality [49, 50] and in classical models [17–19]. The starting point of our investigation is the classical H+ 3 Hamiltonian, HH+3 =
j=1
+
j α pp > 22/3 − 2 + 2−1/3 αep , and R0 2 = 0 + H for 21/3 αep α pp . The case where R0 2 = 0 means that all three particles sit on top of each other. H+ 2 is unstable for √ α pp < 22/3 − 2 + 2−1/3 αep . The second case where analytical solutions are available is H2+ 3 ,
2. Theoretical model
3 P2j
1
,
H2+
stable with the equilibrium nuclear separation R0 3 = 0, so that all particles are located at the same point. Let us now turn to the H2 molecule. The H2 Hamiltonian reads 2 2 P2j p2k 1 + + HH2 = 2 2M k=1 2 (r1 − r2 )2 + αee j=1
(1) where lower case letters rk , pk , 1 k 2 are used to denote electron coordinates and momenta, respectively, and upper case letters R j , P j , 1 j 3, and M are used for the coordinates, momenta, and mass, respectively, of the protons. Throughout this paper, we take M = 1836.2 au as the numerical value of the proton mass. In equation (1), αee , αep and α pp are real, positive constants that specify the three different inter-particle interactions. Before fixing the value of these constants, we have to investigate the stability properties of the molecule as a function of αee , αep and α pp . In addition to
+
1 (R1 − R2 )2 + α 2pp
−
2 2 j=1 k=1
1 2 (R j − rk )2 + αep
.
(6)
Here, the minimum energy is restricted to configurations where the two electrons are either placed symmetrically on a line 2
J. Phys. B: At. Mol. Opt. Phys. 46 (2013) 235601 (a)
(b)
(c)
(d)
E L¨otstedt et al + demand the equilibrium configurations for H+ 2 , H2 and H3 to occur at some non-zero internuclear separation, so that two or three nuclei do not coalesce at the same point. Similarly to the quantum mechanical case, one point on the PES or PEC is obtained by minimizing the total energy of the + 2+ molecule (for H+ 2 , H2 , H3 and H3 ) while keeping the nuclear coordinates fixed. Repeating the minimization for several nuclear configurations, we obtain a PEC as a function of the coordinates specifying the positions of the nuclei. For comparative purposes, a few simulations were performed using two different sets of soft-core parameters, set 2:
αee = 1.4658,
αep = 1.1371,
α pp = 0.8444,
(8)
and set 3: αee = 2,
connecting the two nuclei, or placed symmetrically on a line perpendicular to the connecting line. When the electrons are on the perpendicular line, analytical solutions for the equilibrium distances may be obtained, but not when both electrons and protons are on the same line. Finally, for H+ 3 , no analytical results may be derived. In order to find the combinations of parameters (αee , αep , α pp ) that would make both H2 and H+ 3 stable, we have numerically searched for the particle positions that give the minimal values of the Hamiltonians (1) and (6), using the downhill simplex method [51], and the results are summarized in figure 1. Interestingly, there are areas where H+ 3 is stable, but not H2 . In this case the electron–proton attraction is slightly stronger than the proton–proton repulsion, and therefore the energy is lowered when a proton is added to the otherwise unstable H2 to form stable H+ 3. In the simulations, we use the following set of soft-core parameters, αep = 1.1560,
α pp = 0.8669,
α pp = αep = 3/2.
(9)
Set 2 is similar to set 1, and is employed to illustrate the effect of a slight shift of the model parameters. Set 3 results in a qualitatively different molecular structure, and consequently leads to rather different final spectra, which is discussed in more detail in the appendix. The PECs for the different parameter sets (7)–(9) are shown in figure 2. For comparison, also shown in figure 2(a) are the PECs obtained by quantum chemical methods. The data for H2 was taken from [52], for H+ 2 from [53], while GAMESS [54] was used to produce the one-dimensional cuts of the PESs for H+ 3 (full configuration interaction, 6-311++G(2p) basis set) and H2+ 3 (6-311++G(2p) basis set). The curves displayed in figure 2(a) resemble those of the molecules in quantum mechanics. We note that the energy difference between the H+ 3 curve and the H2 curve, as well as the dissociation energies + of H+ 3 , H2 and H2 are too small compared to the quantum mechanical case. In figure 2(b), the positions of the electrons in the minimum energy configuration are shown. At R = 0.4 au and R = 1.33 au (the equilibrium internuclear distance of H+ 3 ), the electrons in H+ 3 are at the same position, on top of each other. This is energetically favourable since αee > αep and therefore the electron–electron interaction is slightly weaker than the electron–proton interaction. We can see in figure 2(c) that while sets 1 and 2 give rise to similar PECs, set 3 leads to qualitative differences. In particular, the equilibrium internuclear distance for H+ 3 is H+ H+ 3 3 shorter (R0 ≈ 1.1 au for set 3 compared to R0 ≈ 1.3 au for 2 set 1), and the minimum of the PEC for H2 occurs at RH 0 = 0.
2+ + Figure 1. Stability of the H+ 3 , H3 , H2 and H2 molecules in the present model, as a function of the parameters αee , αep and α pp (expressed in au). Red dots: H+ 3 and H2 unstable. Blue diamonds: + H+ 3 unstable, H2 stable. Green squares: H3 stable, H2 unstable. + 2+ Black triangles: H3 and H2 stable. H3 is stable for values of α pp above the thick dashed line given by α pp = 31/3 αep , and unstable for values below. H+ 2 is stable for α pp above the thin, solid line given by √ α pp = 22/3 − 2 + 2−1/3 αep , unstable below. Acceptable values for the parameters are thus given by the black triangles in the area between the two lines.
αee = 1.4991,
αep = 3/2,
3. Simulation procedure We have calculated the response of the present classical H+ 3 molecule to a linearly polarized laser field E(t ) = zˆ E(t ), pointing in the z-direction. E(t ) is taken to be in the following form, π t sin(ω0t )[(t ) − (t − T )], (10) E(t ) = E0 sin2 T where (·) is the step function, E0 the peak field strength, ω0 the laser frequency, T = 2π N/ω0 the total pulse duration, with N being the number of cycles. Below, we take
(7)
referred to as set 1. This parameter set leads to stable H2 and 2+ H+ 3 , but unstable H3 , which is also the quantum mechanical case. This choice resulted from a limited search among the 2+ values leading to stable H2 and H+ 3 and unstable H3 , with the suitability being judged from manual inspections of the respective potential energy surface (PES), or potential energy curve (PEC) in the case of H2 and H+ 2 . In particular, we 3
J. Phys. B: At. Mol. Opt. Phys. 46 (2013) 235601
E L¨otstedt et al
the system vanishes. The system is then allowed to evolve in the absence of the laser field for 300 laser field periods, yielding a trajectory corresponding to a slightly excited, but non-dissociating H+ 3 molecule. The set of initial values rk (0), pk (0), R j (0) and P j (0) for the simulation with the laser field are then taken from points on this trajectory. All sets of initial values are therefore slightly different, but have the same total energy. A new set of initial values is sampled for each trajectory run, so that there are as many sets of initial values as simulated trajectories. In addition, we let the initial distribution be unaligned in space, by rotating the molecule randomly in space before starting each run, without adding an angular momentum. In this way, an ensemble of molecules randomly oriented in space with respect to the laser field polarization vector can be obtained. After running many trajectories with different initial conditions, total and differential probabilities are obtained by counting the trajectories belonging to a particular final channel or energy/angular bin.
4. Results As a result of the interaction with the laser pulse, H+ 3 may decompose via one of the following seven channels [40, 43]. + H+ (12a) 3 → H3 no ionization or dissociation, ⎫ + H+ ⎬ 3 → H + H2 + H+ → H + H dissociation only, (12b) 3 2 ⎭ + → H + H + H H+ 3 + + − H+ single ionization and 3 → H + H2 + e (12c) + + − H+ dissociation, 3 →H +H +H+e
Figure 2. (a) The minimum energy as a function of the nuclear + separation R, for H+ 3 (solid line), H2 (broken line), H2 (dash-dotted 2+ 2+ + line) and H3 (dotted line), in the case of H3 and H3 in the equilateral triangle configuration. The coloured, thick curves were obtained with soft-core parameter set 1 (see (7)), and the thin, grey curves show the PECs (one-dimensional cuts of the PESs in the case 2+ of H+ 3 and H3 ) obtained in quantum mechanics (shifted by −0.3 au in energy). The inset shows a magnification around the minima of the H2 and H+ 3 curves. (b) Positions of the protons (large red circles) and electrons (small blue circles) in H+ 3 , in the minimum energy configuration of the electrons, at fixed internuclear distance R as indicated, using soft-core parameter set 1. (c) Thick, coloured curves: PECs calculated with parameter set 1; thin, black curves: parameter set 2 (see (8)); thin, grey curves: parameter set 3 (see (9)). The curves obtained with parameter sets 1 and 2 are almost overlapping. The inset shows the difference between parameter sets 1 and 2 close to the minimum of the H+ 3 PEC.
+ + + − − H+ 3 →H +H +H +e +e
double ionization (12d) and dissociation.
In the computer code employed, the final channel is automatically recognized at the end of the run, based on the internuclear distances, the electron–proton distances, and the single-particle electron energies. The single-particle electron energies are used to determine whether or not an electron is ejected. We adopt the following definition for the singleparticle electron energy εk(e) , k = 1, 2,
ω0 = 5.8 × 10−2 au, which corresponds to the wavelength 790 nm, and N = 3 cycles. Simulations are performed for different values of the field amplitude E0 in the range 0.1 au E0 0.6 au, as described below. ini ini ini Given the initial positions Rini j , rk , and momenta P j , pk , the classical trajectory of the system is computed accordingly to the classical equations of motion, ∂HH+3 ∂HH+3 dP j dpk =− =− − E(t ), + E(t ), (11) dt ∂rk dt ∂R j
εk(e) =
p2k 1 1 + − , 2 2 2 2 (r1 − r2 )2 + αee (R j − rk )2 + αep j=1 3
(13) which implies that the electron–electron repulsion energy is split equally between the two electrons. An electron k is considered to be ejected if its energy εk(e) > 0. 4.1. Total probabilities We have calculated the total probability for H+ 3 to break up according to the seven possible pathways (12a)–(12d), as a function of the amplitude E0 of the incident laser pulse. The soft-core parameter set 1 (see (7)) was used. The results are shown in figure 3. For each field amplitude, 4×104 trajectories
with 1 j 3, 1 k 2. The initial values for the momenta and positions are drawn from a distribution constructed as follows: first, an additional amount = 0.03 au of the kinetic energy is given to the electrons in the lowest energy configuration, in such a way that the total momentum of 4
J. Phys. B: At. Mol. Opt. Phys. 46 (2013) 235601
E L¨otstedt et al
Figure 3. Total probability per pulse as a function of the field amplitude E0 , for each of the seven possible decomposition + + + pathways: H+ 3 (+), H2 + H (∗), H + 2H (×), H + H2 (• ), + + + H + H2 (), 2H + H () and 3H+ (). Note the logarithmic scale on the probability axis. Statistical error bars are shown when they exceed the size of the curve symbol.
were computed. In terms of the laser intensity I, E0 = 0.1 au corresponds to I = 3.5 × 1014 W cm−2 , and E0 = 0.6 au corresponds to I = 1.3 × 1016 W cm−2 by the conversion formula I/[W/cm2 ] = 3.5 × 1016 E20 /[au]. From figure 3, we see that for low field amplitudes, the neutral (non-ionizing) dissociation channels are dominant. At 0.2 au E0 0.3 au, the probability for the neutral dissociation channels becomes maximal, and drops to negligibly small values when E0 becomes larger. This is natural, since high field strengths inevitably lead to electron ejection. On the other hand, a laser intensity close to 1016 W cm−2 is required in this model for the double ionization to be the dominant pathway. Note that since the model is purely classical, ionization necessarily proceeds with the scheme called over-the-barrier ionization. Another thing worthy of attention in figure 3 is the negligibly small yield of the ionization channel H+ 3 → − H+ + H+ 2 + e (red-filled squares in figure 3). Final states in this channel were found only at two values of E0 at around E0 = 0.3 au (see figure 3). Simulations conducted at other choices of the values of E0 did not produce any events from this pathway. Even though the primary purpose of this paper is not to reproduce experimentally obtained data, we note that events in this channel were found to be rather abundant in the recent experimental results [40]. This illustrates the difficulty in constructing realistic, classical models of multi-particle dynamics, with many possible final products. It should be noted that the yield of the two neutral + dissociation channels H+ 2 + H and H + H2 have comparable magnitudes in the entire range of field strengths considered, which is consistent with the experimental data [40]. This can be compared with the results presented in [37, 38], where a quasiclassical model was used to study the interaction of H+ 3 with intense laser light. In this quasiclassical model, the model parameters were adjusted so that the ionization potentials for + H+ 3 , H2 and H2 agreed with the ionization potentials obtained by quantum chemical calculations. As a result, some features of the experimental spectra for the ionization channels (12c) and (12d) were reproduced. However, for the dissociation channels (12b), the simulation in [37, 38] produced much lower yields for the H+ + H2 channel as compared to the
Figure 4. Distributions of the KER Ekin for different values of the laser field amplitude. In panel (a), we have E0 = 0.32 au, in panel (b) E0 = 0.39 au and in panel (c) E0 = 0.60 au. The legend in panel (c) applies to panels (a) and (b) as well. All curves were normalized to unity at the largest value. Statistical error bars are shown when they exceed the size of the curve symbol. Note that the energy scale is not the same in the three panels.
H+ 2 +H channel, in disagreement with the experimental results. Finally, we comment that in this paper, we have not considered averaging of the results over the laser intensity distribution in the laser focus. Such averaging leads to total probability curves that are monotonically increasing with increasing laser intensity for all channels [38], due to the contribution from the regions of lower laser intensity. 4.2. Kinetic energy distributions A quantity often measured in experiments is the kinetic energy release (KER). In this study, the KER, which we label by Ekin , is defined as the sum of the kinetic energy of the protons at the end of the trajectory, when the fragments are sufficiently separated for the proton–proton repulsion term to be negligible,
fin 2 3 Pj . (14) Ekin = 2M j=1 The KER for three different values of E0 , calculated using soft-core parameter set 1, is displayed in figure 4. To obtain better statistics for the KER and angular distributions, 2 × 105 trajectories were run at each value of E0 . The KER for the two ionization channels H+ + H+ + H and H+ + H+ + H+ depends strongly on the field strength. For the intermediate value E0 = 0.39 au in figure 4(b), the two peaks are well separated, while at the higher value E0 = 0.6 au in figure 4(c), the peaks are overlapping. Also interesting is the peak at Ekin ≈ 0.6 au for the neutral dissociation channel H+ + H + H 5
J. Phys. B: At. Mol. Opt. Phys. 46 (2013) 235601
E L¨otstedt et al
in figure 4(b). This peak, which has the same KER as the ionization channel H+ + H+ + H, can be explained by a mechanism called ‘frustrated ionization’, first used to clarify the high-energy neutral H fragments measured in the laserinduced dissociation of H2 [55, 56]. In our calculation, the explanation is the following: one electron absorbs energy from the laser field, but lacks the necessary kinetic energy to escape the attractive potential formed by one of the protons. At the end of the pulse, it therefore recombines to form an excited H atom. Since a highly excited electron moves around the nucleus at a large radius, the proton is less screened from the other protons, and consequently the resulting KER of the H+ +H+H fragments may reach the same values as for the corresponding ionization channel H+ + H+ + H. We note that this process was also found in the laser-induced dissociation and ionization of H+ 3 , studied with a quasiclassical model [37, 38], and also in experiments with D+ 3 [57]. We therefore suggest that this kind of process could commonly exist in any unimolecular break-up processes induced by an intense laser field, as long as long-ranged Coulomb forces are present. 4.3. Angular distributions In the case of photodissociation of H+ 3 , it is natural to investigate three kinds of angular distributions [40]. The angle θ is defined as the angle between the polarization direction of the laser field and the normal vector of the plane spanned by the final momentum vectors of the protons. The total final momentum of the protons is in general non-zero, while the sum of the momenta the two electrons 3of thefinthree protons 2 and fin = P = − p = 0. This means vanishes, Pfin tot j=1 j k=1 k that the final momentum vectors of the protons do not lie exactly in the same plane. However, due to the large mass of fin the protons, we have |pfin k | |P j | for all j and k, and therefore
3 fin |Pfin tot |/ j=1 |P j | 1. It is therefore meaningful to speak about an approximate plane spanned by the final momentum vectors of the protons, we define the normal vector nˆ to this plane by the average b fin , b=− (−1)k− j Pfin (15) nˆ = j × Pk . |b| j
