Electrons as probes of dynamics in molecules and clusters : a contribution from Time Dependent Density Functional Theory
EElectrons as probes of dynamics in molecules andclusters : a contribution from Time Dependent DensityFunctional Theory
P. Wopperer a,b , P. M. Dinh a,b ∗ , P.-G. Reinhard c , and E. Suraud a,b,d a CNRS, LPT (IRSAMC)118 route de Narbonne F-31062 Toulouse C´edex, France b Universit´e de Toulouse, UPS, Laboratoire de Physique Th´eorique (IRSAMC)118 route de Narbonne F-31062 Toulouse C´edex, France c Institut f¨ur Theoretische Physik, Universit¨at Erlangen,Staudtstrasse 7 D-91058 Erlangen, Germany d Physics Department, University at Buffalo, The State University New York, Buffalo, NY 14260,USA
Abstract
There are various ways to analyze the dynamical response of clusters and molecules to elec-tromagnetic perturbations. Particularly rich information can be obtained from measuring theproperties of electrons emitted in the course of the excitation dynamics. Such an analysis ofelectron signals covers observables such as total ionization, Photo-Electron Spectra (PES), Pho-toelectron Angular Distributions (PAD), and ideally combined PES/PAD. It has a long historyin molecular physics and was increasingly used in cluster physics as well. Recent progress inthe design of new light sources (high intensity, high frequency, ultra short pulses) opens newpossibilities for measurements and thus has renewed the interest on these observables, especiallyfor the analysis of various dynamical scenarios, well beyond a simple access to electronic den-sity of states. This, in turn, has motivated many theoretical investigations of the dynamics ofelectronic emission for molecules and clusters up to such a complex and interesting system asC . A theoretical tool of choice is here Time-Dependent Density Functional Theory (TDDFT)propagated in real time and on a spatial grid, and augmented by a Self-Interaction Correction(SIC). This provides a pertinent, robust, and efficient description of electronic emission includ-ing the detailed pattern of PES and PAD. A direct comparison between experiments and wellfounded elaborate microscopic theories is thus readily possible, at variance with more demandingobservables such as for example fragmentation or dissociation cross sections.The purpose of this paper is to describe the theoretical tools developed on the basis of real-time and real-space TDDFT and to address in a realistic manner the analysis of electronicemission following irradiation of clusters and molecules by various laser pulses. After a generalintroduction, we shall present in a second part the available experimental results motivatingsuch studies, starting from the simplest total ionization signals to the more elaborate PES and Preprint submitted to Elsevier 12 July 2018 a r X i v : . [ phy s i c s . a t m - c l u s ] J u l AD, possibly combining them and/or resolving them in time. This experimental discussionwill be complemented in a third part by a presentation of available theoretical tools focusing onTDDFT and detailing the methods used to address ionization observables. We shall also discussthe shortcomings of standard versions of TDDFT, especially what concerns the SIC problem,and show how to improve formally and practically the theory on that aspect. A long fourth partwill be devoted to representative results. We shall illustrate the use of total ionization in pumpand probe scenarios with fs lasers for tracking ionic dynamics in clusters. More challenging fromthe experimental point of view is pump and probe setups using attosecond pulses. The effortthere is more on the capability to define proper signals to be measured/computed at such a shorttime scale. TDDFT analysis provides here a valuable tool in the search for the most efficientobservables. PES and PAD will allow one to address more directly electronic dynamics itselfby means of fs or ns laser pulses. We shall in particular discuss the impact of the dynamicalregime in PES and PAD. We shall end this fourth part by addressing the role of temperaturein PES and PAD. When possible, the results will be directly compared to experiments. Thefifth part of the paper will be devoted to future directions of investigations. From the richchoice of developments, we shall in particular address two aspects. We shall start to discuss theinformation content of energy/angular spectra of emitted electrons in case of excitation by swiftand highly charged ions rather than lasers. The second issue concerns the account of dissipativeeffects in TDDFT to be able to consider longer laser pulses where the competition between directelectron emission and thermalization is known to play a role as, e.g., in experiments with C .Although such questions have been superficially addressed in the simple case of alkaline clustersby means of semi-classical methods, no satisfying quantum formulation, compulsory for mostrealistic systems, is yet available. First encouraging results will be presented on that occasion.We shall finally give a short conclusion. Key words:
Time-Dependent Density Functional Theory, Electronic observables, Ionization, Lasers,Charged projectiles, Photo-Electron Spectrum, Photoelectron Angular Distribution, Orientationaveraging, Self-interaction correction, Time-resolved observables, Temperature effects, Dissipationeffects
PACS:
Contents ∗ Corresponding author
Email-address : [email protected] .3 Dynamical aspects in photoelectron spectra 614.4 Photoelectron angular distributions (PAD): a sensitive tool 694.5 Impact of temperature in PES and PAD 775 Future directions 865.1 An excursion into irradiation by charged projectiles 865.2 Towards quantum dissipative electron dynamics 896 Conclusions 95References 99
1. General introduction and physical context
Irradiation of matter constitutes a key tool in physics, chemistry, and biology, foranalyzing structural and dynamical properties of atoms, molecules, clusters and bulkmaterial. Lasers offer here an especially flexible and powerful instrument which has beenwidely exploited, especially during the last decades with the enormous technologicalprogress reached in the manipulation of laser light [1, 2, 3]. We dispose now of a broadchoice of laser intensities, frequencies, pulse lengths, and pulse shapes. Collisions withcharged projectiles [4] are also used as sources of short electromagnetic pulses. However,they often require access to dedicated facilities.Radiation damage is the other side of irradiation studies and it is of high current inter-est, for example in connection with biological tissues (”human-controlled” as in a medicalcontext or ”natural” when referring to earth or space radiations) [5]. There are also otherinteresting domains of application. A typical example is the case of the irradiation of ma-terials (especially insulators) with applications to nuclear waste management. The fieldis rather unexplored from the microscopic dynamical point of view and any possibilityof treating irradiation scenarios on large systems would be here of invaluable help [6]. Inboth above examples, though, the lack of understanding of microscopic mechanisms callsfor dedicated studies on prototype, finite systems. Let us cite as an example the detailedstudies of irradiation of molecules of biological interest coated by a finite and well knownnumber of water molecules [7]. The study of the irradiation of finite molecular systemsand clusters is thus not only of interest for basic science but also for a wide range ofpractical applications.In all cases, the immediate electronic response of the irradiated system plays a keyrole as the doorway to all subsequent dynamical scenarios. A basic feature is here theoptical response corresponding to electronic oscillations [8, 9]. It delivers a first overviewof the coupling between irradiation and matter in a large variety of dynamical situations,from gentle to strong perturbations [10, 11, 12, 13]. Optical response related to photo-absorption is the leading signal in the case of gentle perturbations. It has been exploredin great detail for a large variety of electronic systems, from bulk down to atoms. For thecase of stronger perturbations, further response channels, especially ionization, becomehighly relevant [14, 12, 15]. Still, the optical response spectrum, which characterizes thestructural coupling of the system to light, provides a highly valuable information on anyensuing response mechanism, especially on ionization pattern. A typical example hereis the case of resonant ionization occurring when the laser frequency comes close to aneigenfrequency of the system [14].Equally important in energetic irradiation processes is electron transport, particularlyelectron emission. As typical examples, one can cite the many studies on irradiation of3lusters by short and intense laser pulses [15], providing invaluable information especiallythrough energy (Photo-Electron Spectra, PES [16]) and, more recently, angle-resolved [17]distributions of emitted electrons (Photoelectron Angular Distributions, PAD). Electronemission can also change the resonant ionization conditions in the course of time evolutionwhich, in turn, influences back again the optical response, making the whole scenarioextremely rich [15]. Secondary electrons in DNA damage [18] also provide a remarkableexample where a microscopic understanding of irradiation damage in biological systemswill only be achieved when including such complex non-linear electronic effects. A deeperunderstanding of the underlying mechanisms is highly desirable, as this example is ofgreat practical interest, especially in relation to oncology [5].The analysis and understanding of electronic emission from a finite system is thus akey issue in a wide range of physical, chemical and biological processes. Electrons areusually the first constituents to respond to an electromagnetic pulse. Strong excitationslead to immediate ionization of the system, often with dramatic long-time effects as,e.g., dissociation or Coulomb explosion [13]. It implies electronic transport and possibleindirect effects on neighboring species. A typical example of indirect effects is providedby Dissociative Electron Attachment (DEA) in biological systems [5] where electronsemitted somewhere else are attached to a target biological molecule which, in turn, leadsto the break up of the latter. Emitted electrons may also provide valuable insight intoreaction pathways when properly tracked. Typical examples are here PES and PAD.Moreover, Time-Resolved (TR) PES and PAD have been recorded in molecules [19] andmore recently in clusters, see e.g. [20]. Electrons are thus leading players at all stagesof an excitation of a system subject to an external electromagnetic perturbation (i.e.an irradiation). They are the first to respond at short time scales and distribute thenthe excitation more or less quickly to other degrees of freedom. They are finally usefulprobes along the whole dynamical process, especially when emitted from the system andproperly recorded.Analyzing the characteristics of emission properties of clusters and molecules is thusat the core of the understanding of irradiation processes. The numerous new experi-mental developments in analysis of electronic emission (PES, PAD) now allow an everimproving detailed access to electron dynamics in irradiated species. In turn, a theoreti-cal description of these highly involved dynamical scenarios calls for dedicated modeling.It is the aim of this paper to provide an overview of the theoretical description of observ-ables from electron emission on the basis of the well established theoretical frameworkof Time-Dependent Density Functional Theory (TDDFT) [21]. This will be done with aview on applications, as far as possible in direct relation to ongoing experiments. Beforegoing into the details, we will in this introductory section briefly remind the reader thetypical systems (and associated scales) that we aim at describing. It is also of relevanceto address here basic properties of laser pulses as presently accessible experimentally.1.1.
On the typical systems considered in this paper
In order to provide a basis for the forthcoming discussions, we shortly present here afew typical systems we shall consider in the following. This will be the occasion to remindtypical scales associated to these systems, especially in terms of times and energies.Fig. 1 provides four examples of systems computed with the tools described in Secs. 3,4.1.1 and 4.1.2. They cover several different binding types and properties. The figureshows single-particle (s.p.) energies, optical response and ionic structure. The four pre-sented systems are Na as an example of a simple metal cluster, C for its outstandingproperties and many applications, H O as a prototype of a covalent molecule, and C asa simple carbon chain, which displays interesting optical properties. Let us analyze eachsystem separately to extract typical properties.We start with the Na cluster (upper left block in Fig. 1) which is a medium sizealkaline cluster. It contains 40 valence electrons forming an electronic shell closure. Thisleads to a particularly abundant/stable species. The s.p. energies span an energy range oforder 2.6 eV and the Ionization Potential (IP) is of order 5.3 eV. Such values are typicalof alkaline clusters. The optical response displays a pronounced collective, essentiallysingle, peak around 2.6 eV. This is called the Mie surface plasmon and it is a typicalmode for simple metal clusters. In larger clusters, the density of s.p. states grows, whichleads to more Landau fragmentation and somewhat broadens the plasmon peak. The Mieplasmon frequency is related to a typical time scale of order 1.5 fs, again a characteristictime scale for simple metals. Ionic time scales (not shown in Fig. 1) are more sensitive tothe actual material due to the largely differing atomic weight. In Na clusters, vibrationalmodes typically lie in the 10 meV range and are associated to ionic motion in the 100 fsrange.The second example is C (upper right block) with 240 valence electrons (4 per Catom). The s.p. energies now ranges in a span of about 17.3 eV, much wider than inNa. The IP is 7.6 eV. The deeper binding and broader span of energies is typical ofcarbon, and more generally of organic systems with a covalent binding. Due to the highsymmetry close to sphericity, the optical response exhibits the same behavior in all spatialdirections. It has, however, a more complex structure than in Na . One can identifytwo prominent features, a strong resonance peak just below the IP and a much broadenedpeak centered around 20 eV. The latter part of the optical response lies well above theIP, whence its highly fragmented structure. It is considered to represent the Mie surfaceplasmon in C . The energies are higher and thus the associated time scales much smallerthan in Na, typically well sub-fs. Ionic vibration energies typically lie in the 40-200 meVrange with associated time scales of order 20-100 fs.The case of the small carbon chain C (lower left block) is complementing C in thesense that it has the same binding type, but a different geometry and thus differentoptical response. The s.p. energies span of the 20 valence electrons is of order 14 eV, andthe IP of order 9.9 eV. These values are of the same order of magnitude for larger chains.According to the linear geometry of the chain, the optical response shows a dominantresonance peak along the longitudinal direction at a frequency of 6.4 eV. The transversemodes are suppressed by at least one order of magnitude (mind that transverse strengthshave been multiplied by a factor of 20 to allow a better graphical comparison with thelongitudinal modes) and are significantly fragmented. There are three main peaks : one atthe same energy as the longitudinal plasmon peak, and two other ones at higher energiesnear the IP energy. The all dominant longitudinal mode lies well below the IP, a featurecommon to all carbon chains. Associated time scales are now typically ranging from sub-fs to fs. Ionic vibration energies lie again in the 0.15 eV range with associated time scalesof order 27.6 fs.We finally discuss the case of the prototypical water molecule H O (lower right block)which has 8 active valence electrons in our calculations (6 for O and 1 per each H).5 " i ( e V ) Na , 40 e -25-20-15-10-50 C , 240 e -30-25-20-15-10-50 " i ( e V ) C , 20 e -30-25-20-15-10-50 H O, 8 e = { ˜ D ( ! ) } Energy (eV) 0 5 10 15 20 25 30 35Energy (eV)0 5 10 15 20 25 30 35 = { ˜ D ( ! ) } Energy (eV) 5 10 15 20 25 30 35Energy (eV) xy dir. z dir. xyz dir. xy dir. ⇥ z dir. x dir. y dir. z dir. Fig. 1. Four typical examples of molecules and clusters explored theoretically in this review, namelythe metal cluster Na (upper left), the Buckminster fullerene C (upper right), the carbon chain C (lower left), and the covalent molecule H O (lower right). For each panel: Top row : ionic structures (allplotted at the same scale) and single particle energies of the valence electrons whose number is indicated;Bottom row : corresponding optical response. For C , the transverse spectrum is multiplied by 20 to easethe comparison with the longitudinal response. The vertical dashes indicate the position of the ionizationpotential in each case. , in spite of the muchsmaller number of electrons. The IP is of order 15.1 eV. Such large IP’s are typical ofcovalent systems of small to moderate size. The optical response, as well, is typical ofcovalent molecules with its highly fragmented structure above the IP, and some isolatedlow energy peaks below the IP. Associated time scales lie well below fs. Ionic vibrationsare more energetic than in other systems because of the especially light H species andthe strong covalent binding between H and O. The O-H ionic vibration energy is about0.5 eV with associated period of 8.3 fs.All in all, the four above examples point out the diversity and richness of the varioussystems nowadays accessible to both experimental and theoretical investigations. Thevarious cases also show that the range of energy and time scales to be investigatedis rather large from attosecond to several fs for electrons, and from fs to ps for ions. Inaddition, the optical spectra exhibit different pattern. Specific for metals is the especiallywell marked Mie surface plasmon with simple scaling properties with size [10]. The caseof covalent systems is more involved with basically no simple scaling properties, butnevertheless some generic trends. Optical spectra are generally much more fragmentedbelow and even more so above IP. Pure carbon systems contain besides covalent bindinga fraction of metallic binding which produces also plasmon structures amongst the highlyfragmented spectrum.Optical response is the key to understanding the coupling of the system to laser light,at least in the frequency-dominated regime (see Secs. 1.2 and 2.1). This will constitutea mostly used tool of investigation of dynamical scenarios in the following discussions.Before introducing the actual observables which can be attained that way, we will brieflydiscuss present days capabilities of lasers and the description of the electromagnetic fieldsthey deliver. This aspect is addressed in the following Sec. 1.2.1.2. On excitation mechanisms
Cluster dynamics requires excitation of the cluster formerly resting in its ground state.In this paper, we will exclusively addess excitation by electromagnetic fields, predom-inantly by laser pulses and in a few cases by short pulses from collisions with highlycharged ions. The corresponding excitation mechanisms are shortly explained in this sec-tion. Thereby, we focus on laser properties and finally address ion collisions in a shortparagrph.1.2.1.
Laser pulse characteristics
Laser science has experienced impressive progress during the last few decades [3]. Theversatility of laser pulses has increased remarkably, thus allowing one to shape a widerange of dynamical scenarios in the course of irradiation processes. We briefly remindhere key quantities of the laser pulses we are going to use in the following. Throughoutthis paper, we shall work in the dipole approximation which requires that the irradiatedsystem is much smaller than the laser wavelength λ = 2 πc/ω las . In practice, the dipoleapproximation is well justified in the optical domain ( λ ∼ µ m) for systems of nm size.It may become questionable for XUV photons and very large clusters in which fieldvariations inside the system itself should be accounted for. But we shall not considersuch cases here. In the non-relativistic regime, linearly polarized laser pulses acting on7toms, molecules or clusters can then be described as a homogeneous time-dependentelectric field of the form E ( t ) = e pol E f ( t ) sin( ω las t + ϕ ) . (1)In this expression, e pol denotes the (linear) polarization, E is the peak field strength, ω las is the photon frequency, and ϕ ( t ) is some phase shift, usually assumed to be zero.Finally f ( t ) is the pulse envelop. The peak laser intensity is I = c ε E / c beingvelocity of light in vacuum) usually expressed in W/cm . The net yield in a laser pulseis often characterized by the fluence F = (cid:82) d tI ( t ) ≈ I T FWHM , where the latter time T FWHM stands for the Full Width at Half Maximum of the pulse. This allows one tocompare the energy impact of laser pulses with different durations.For the sake of simplicity, we keep in the present discussion a fixed value of ω las butthe latter quantity can also be made time-dependent (”chirped”) which can induce inter-esting effects [3]. One may also render the phase ϕ time-dependent, which could produceinteresting phenomena. We shall not discuss these aspects here. The laser polarizationis usually taken linear but there also exists experiments/calculations using circularly po-larized light [3]. Again, we shall not discuss much such cases in the following and thusrecur to a linear polarization for the present discussion.The laser pulse envelop can be varied in a large range. Most flexible, and most widelyused, are optical lasers with pulse lengths from nano-seconds down to atto-seconds [22,23]. Free Electron Lasers (FEL) [24, 3] are yet on their way to comparable flexibility,with present pulse lengths down to 20 fs. It should also be noted that the actual shapeof f ( t ) is not exactly known experimentally. In many situations, the actual pulse ofinterest is built upon a (hopefully harmless) background of a long, low intensity, pulse.Moreover, the peak intensity has a spatial variation decreasing towards the edges of thepulse. This has to be kept in mind when assigning the observed signal to the laser pulsecharacteristics. Ignoring background, experimental short laser pulses have a pulse profileof Gaussian type. The theoretical situation is simpler as the pulse profile can be exactlyspecified. The Gaussian profile is theoretically not welcome since it never fully vanishesand requires unnecessarily long computation times to cover the pulse sufficiently well.Therefore, we mostly use for computations a sin pulse : f ( t ) = sin (cid:18) tπT pulse (cid:19) θ ( t ) θ ( T pulse − t ) . (2)where θ stands here for the Heaviside function. This pulse is limited to a finite timeinterval t ∈ [0 , T pulse ], but soft enough to deliver a clean frequency spectrum. It can besimply characterized by its FWHM which is in this case T FWHM = T pulse /
2. Note thatthe FWHM of the intensity I ( t ), which is proportional to the square of the field E ( t ),is rather T pulse /
3. The pulse maximum occurs at t = T pulse /
2. Note that the sin profileis written here for the laser field amplitude, which means that the time profile of theintensity time has a sin shape. Thus far, we have discussed simple one-peak pulses.More flexibility is conceivable. The next important tool are dual pulses as used in pump-and-probe experiments in which the laser irradiation is performed in two steps. We shallillustrate such cases at several places below.In practice, the effect of the laser field will be accounted for in our calculations as anexternal potential U ext ( r , t ) which delivers a time-dependent perturbation. In the longwavelength limit, the electric field is homogeneous and delivers the potential :8 ext ( r , t ) = − e E f ( t ) sin( ω las t ) e pol · r , (3)where f ( t ) is the time profile usually taken according to Eq. (2). This is, in fact, the laserfield in space gauge. Equivalently, one can use the velocity gauge for which the laser fieldis described by the interaction operator :ˆ U ( v )ext = − ec E F ( t ) e pol · ˆ p . (4)The rules of gauge transformation relate time profiles and wave functions by : F ( t ) = (cid:90) t −∞ d t (cid:48) f ( t (cid:48) ) sin( ω las t (cid:48) ) , (5) ϕ ( v ) i ( r , t ) = ϕ i ( r , t ) exp [i E F ( t ) e pol · r ] . (6)Both gauges are fully equivalent. Which one is to be preferred is a matter of the actualnumerical scheme. Most observables are not even sensitive to gauge. An exception is theevaluation of photoelectron spectra where the phase of the wave function plays a role. Inthis case, one has to consider gauges carefully. This will be addressed in more detail inSec. 3.3.1.2.2. Varying laser characteristics
As pointed out above, all laser parameters can be tuned in rather large ranges. Thepoint is illustrated in Fig. 2 which displays typical regions of interest in the intensity-frequency plane. One can notice the enormous large intensity range of optical lasers. Butthe range of available conditions is also dramatically extended by FEL, which exist forphotons in the IR, VUV and X-ray regime. Also indicated outside the axes are corre-sponding regions of relevance in atoms and molecules in terms of energy/frequency andfield/intensity. The lowest frequencies in the deep IR are associated with molecular vibra-tions, while the range around visible light belongs to the dynamics of valence electronsand core electrons move at much higher frequencies in X-ray regime. The gray box belowthe plot indicates typical atomic and molecular field strengths in terms of an equivalentlaser intensity.Laser characteristics have to be considered in relation to the electronic response. Thisis usually quantified via the ponderomotive potential U p and the associated Keldyshparameter γ . U p represents the electron kinetic energy (averaged over one photon cycle)of a freely oscillating electron (pure quiver motion, no drift velocity) in a laser field. Atpeak laser intensity, it reads : U p = e E m el ω = 9 . × − eV × I [W / cm ] ( λ las [ µ m]) , (7)where λ las is the photon wavelength. The other aspect concerns the electronic bindingin the system, which can be quantified by the ionization potential (IP) with associatedenergy E IP . What counts is the relation between E IP an U p , quantified by the Keldyshparameter [25] : γ = (cid:115) E IP U p = (cid:115) E IP ω I . (8)9 ig. 2. Schematic representation of various dynamical regimes as a function of laser intensity I and photonfrequency ω . The dashed diagonal line represents frequency-intensity combinations with constant Keldyshparameter γ = 1, see text for details. This line characterizes the transition from photon-dominated tofield-dominated regime for an assumed IP of a few eV. The blocks to the right side indicate typicalfrequency ranges as labeled. The block below the plot indicates typical atomic field strengths related togiven laser intensities. The value γ = 1 (see Fig. 2 in the case E IP = 1 eV) separates two regimes. For γ (cid:28) γ (cid:29) Not on lasers: collisions with fast ions
There is an alternative excitation mechanism by collisions with charged projectiles.We shall also marginally consider a few examples of collisions with fast ions and thuscomment briefly about this tool here. Experiments with charged, fast ions often requireaccess to large scale facilities. Thus there are much less experiments with irradiation bycharged projectiles than by the more easily accessible and versatile lasers. Although col-lisions with charged particles also provide a strong electromagnetic perturbation (oftenin form of a short pulse as soon as the projectile velocity is large enough), the charac-teristics of the perturbing field are significantly different from those delivered by a laserpulse. While lasers provide (up to details) an electromagnetic field with a well definedfrequency band (basically the laser frequency), collisions with charged projectiles delivera perturbation covering a very broad band of frequencies, the broader the shorter thepulse. This delivers useful, complementing information to that attained from lasers. It isimportant to note that collisions with charged projectiles also concern a wide range ofpotential applications of irradiation dynamics, especially in relation to radiation damage10nd applications thereof. The present review concentrates on laser excitations. Neverthe-less, we shall discuss a few cases with high energy projectiles. For them, the deliveredelectromagnetic perturbation can be modeled as an instantaneous boost ( ∝ δ ( t )) at theinitial time of the simulation. This is the way we shall treat this case in the following(see in particular Sec. 5.1).
2. From integrated to detailed observables
Electronic emission can be analyzed at various levels of sophistication, starting fromfully integrated quantities (total ionization) down to energy-resolved (Photo-ElectronSpectra, PES) and angle-resolved (Photo-Angular Distribution, PAD) quantities. Timeis also a key quantity as ionization signals can be followed in time, leading to Time-Resolved (TR) results. We briefly describe in this section the various types of observablesexperimentally accessible, starting from the simplest one, that is the total ionization, tothe most elaborate ones (TR-PES and PAD). In terms of cross sections, this means thatwe go from integrated ones to single-differential and even double-differential ones, allpossibly time-resolved. Before discussing these various observables, we briefly introducekey mechanisms of ionization, again focusing the discussion on laser induced ionization.2.1.
Ionization mechanisms
Basic ionization mechanisms are illustrated in Fig. 3. We start from the simplest case (b)
Fig. 3. Schematic view of ionization mechanisms in atoms and clusters. Occupied electron states areindicated by horizontal bars. (a) : multiphoton and optical field ionization in an atom, for which aredrawn the potentials of the unperturbed ion V ion (full line), of the laser V las (dots), and of the sum ofboth (dashes). (b) : inner and outer ionizations in a cluster, with the effective electron potential (withoutlaser) shown as a solid line. Adapted from [15]. of an atom (left panel in Fig. 3) to introduce two basic ionization mechanisms. Thefirst one corresponds to a vertical excitation of a bound electron by absorption of oneor several ( ν ) photons (Multi-Photon Ionization or MPI). This mechanisms may spreadover several laser cycles and prevails in weak and moderate fields, usually quoted pertur-bative regime. It is associated to large values of the Keldysh parameter ( γ (cid:29) γ (cid:46) inwhich the notion of collective plasmon is hard to disentangle from that of a moleculardipole transition. Na possesses two valence electrons. The mechanism actually remainsthe same and is thus illustrative of the role of the optical response. We shall considerplasmon effects on some examples later on (see in particular Sec. 4.2.1 and Fig. 30). Thedashed curve shows the optical response of the system with a well identified peak at2.12 eV. The full curves display the total ionization as a function of laser frequency for aset various laser intensities between 10 and 10 W/cm . One clearly observes that theionization signal directly follows the optical response : attaching a resonance peak leadsto enhanced ionization. The effect is especially visible at low intensity and vanishes withincreasing intensities. We gradually leave the photon-dominated regime (low intensity)to reach the field-dominated one. In terms of the Keldysh parameter γ , it decreases.In the present test case, γ takes values typically between 60 and 250 in the resonanceregion at low intensity, and reaches values between 2 and 10 in the high intensity case.The role of resonance peaks is thus crucial here and it should be noted that it does notreduce to the linear regime of excitation. The total ionization may reach rather largevalues (more than half of the available valence electrons) with increasing laser intensity,and still, the resonance enhancement remains very clear. This indicates that it will have12 .10.40.711.3 Na T pulse = 290 fs opt. resp. − − − N e s c ω las (eV) W/cm W/cm W/cm W/cm Fig. 4. Total ionization N esc of Na after irradiation by laser pulses polarized along the cluster axis,with duration of 290 fs, as a function of laser frequency ω las , for four different intensities as indicated.Dashed curve : optical response (power spectrum) of Na . Mind that for N esc < .
1, the vertical axis isin logarithmic scale, while a linear scale is used above 0.1. to be considered whatever the dynamical regime in the following, especially in the caseof metals. Although the basic enhancement mechanisms remain similar in non-metallicsystems (see Sec. 4.2.1), resonances are usually less collective and more narrow so thattheir impact is somewhat different. Still, in many systems such as for example C , oneobserves a wide bunch of resonances above continuum threshold which very clearly playa key role in the dynamics.2.2. Total ionization
Total ionization is the simplest ionization signal one can measure. Still, it alreadybrings interesting information, although not highly detailed, on irradiation mechanisms,as we just discussed in the previous section. We here illustrate the point on two exam-ples taken from rather original scenarios. The first case results from an irradiation withextremely large frequencies obtained from a FEL, and the second one directly addressesthe dynamical evolution of the system in a time-resolved experiment.Fig. 5 shows time-of-flight (TOF) spectra of Xe clusters irradiated by a FEL of fre-quency around 12 eV [32]. The TOF gives access to the various charge states attainedafter irradiation by a laser of intensity 2 × W/cm and pulse length of 100 fs. Thestriking point of the figure is the differences observed between the various cluster sizesin terms of attained charge states. While the atomic gas, under the present laser con-ditions, only allows to access singly charged cations, increasing cluster size allows toprogressively reach larger and larger charge states, clearly up to 8+ in the largest system13 ig. 5. Time-of-flight mass spectra of ionization products of Xe atoms (bottom curve) and clusters ofvarious sizes N as indicated. Irradiation was performed by a free-electron laser of wavelength of 98 nmand an average intensity of 2 × W/cm . The line splitting of the atomic spectrum (bottom curve) isdue to different isotopes. Inset : ion kinetic energies as a function of ion charge, in the case of 1500-atomclusters. From [32]. of about 3,000 atoms. The case very nicely illustrates the well known difference betweenenergy absorption by single atoms and clusters, as discussed on many occasions in thepast (see for example [15] and references therein). The mass peaks are rather broad. Theyare furthermore displaced with respect to the calculated flight times indicated by thinvertical lines (corresponding to the different charge states) in the top of the figure. Thisis an indication that ions have high kinetic energies. Not surprisingly, one can also notethat the higher the charge state, the higher the ejection energy (see inset in Fig. 5) andthe larger the above mentioned peak displacements.Analysis of total electron emission (or alternatively of charge state of ionized clusters)also gives information on the dynamics of the charging process. One of the most strikingearly example of such an analysis can be found in a series of experiments led by theRostock group on large size Pb clusters [33, 34, 35]. These experiments have shownstrong enhancement of cluster ionization for optimal pulse durations. More specifically,one observed Pb ions with very large state states, much larger than those attained inan atomic gas. Moreover, the attained charge state q strongly depends on the pulseduration. The shortest and most intense pulses of duration 150 fs yield ions up to chargestate q = 20. When increasing pulse duration, both the maximum charge state and thesignal intensity do grow towards a maximum attained for an optimal pulse width of 800fs. Charge states up to q = 28 can then be identified. For longer pulses, both maximum q and signal decrease again. Although other mechanisms can be envisioned, the efficientcharging for a certain pulse duration was in most cases attributed to resonant heating(plasmon-enhanced ionization) [36, 37, 38, 13].The above case of ionization enhancement was attained with a single laser pulse and14nly provides a rather indirect indication on the ionization mechanism. More detailedinvestigations were led with dual (or pump-and-probe) pulses, especially in the case ofAg clusters of about 20 000 atoms. An example of such experiments is shown in Fig. 6where one focuses on the yield of Ag and the maximum of the emitted electrons as afunction of pulse separation. One observes a strong variation of ionic signal as a function Fig. 6. Ionic charge state Ag yield (diamonds, left axis) in relation to the maximum kinetic energy ofthe emitted electrons (dots, right axis) following a laser excitation of Ag clusters of about 20 000 atomswith dual 100 fs laser pulses of intensity of 8 × W/cm and wavelength of 800 nm. From [39]. of pulse separation (delay between pump and probe) with a clear maximum around 5ps [39]. Such a behaviour indicates that cluster activation and enhanced ionization canbe clearly disentangled, which is also found in numerical simulations [40, 39, 41]. Thisagain provides an interesting insight the dynamical evolution of the system. There areeven clear indications that a sequence of two pulses might constitute an optimal pulseprofile for the production of very high charge ions [42], provided a proper tuning of pulseparameters. And pushing again the argument, one can even envision a route for targetedcontrol of the cluster dynamics [38]. Finally a word on the maximal electron energy shownin Fig. 6 is in place. The coincidence of high ionization yield and maximal electron energyagain points out the leading role of collective excitations, and this in both channels. Thisis compatible with other observations [43, 44].2.3. Energy- and angular-resolved ionization
The next step in the analysis of electronic emission consists in characterizing the prop-erties of the emitted electrons in terms of kinetic energy and angular distribution. Thisleads to Photo-Electron Spectra (PES) for the energy analysis and Photo-Angular Distri-bution (PAD) for angular signals. The terms come from laser irradiation but the signalsthemselves can as well be recorded in any ionization scenario, for example from collisionswith highly charged ions (see Sec. 5.1). PES and PAD signals turn out to provide ex-tremely rich information, both from a structural and from a dynamical viewpoint. Webriefly discuss their properties in this section and illustrate them on a few examplescovering several dynamical situations. 15inetic energies of emitted electrons can be measured in several ways. TOF devicesprovide here a versatile tool, but this time applied to the electrons themselves (whilethey are traditionally used for ions). Because of the well defined mass to charge ratiofor an electron, the arrival delay directly maps the electron kinetic energy, provideda carefully guiding of the electron flow e.g., by a magnetic mirror [11]. Another veryinteresting technique is provided by photo-imaging spectroscopy, also known as VelocityMap Imaging (VMI). This technique is more and more routinely used and providesa remarkable tool of investigation. It is based on a static electrical field which allowsone to map the distribution of electron velocities onto definite positions on a detectionscreen [45]. This is a polar representation of a velocity-resolved (or a momentum one) andangular-resolved photoelectron spectrum. Two experimental examples are presented inFig. 7, one in the metal cluster Na − in the monophoton regime [46], and another one inC in the multiphoton regime [47]. In these two examples, the vertical direction stands C ! θ Fig. 7. Left : Raw velocity map image of Na − irradiated by a laser pulse of frequency of 4.02 eV andpolarized vertically [46]. The length of the arrow stands for the velocity of the photoelectron, and θ itsangle to the laser polarization axis. Right : Inverted momentum map image of C irradiated by a laserpulse of FWHM of 150 ± . × W/cm , frequency of 3.68 eV, and polarizationalong the vertical axis [47]. for the laser polarization axis. If one draws an arrow from the origin of the circle, its lengthrepresents the norm of the velocity (or the momentum), while the angle θ to the verticaldirection is the angle of the photoelectron with respect to the laser polarization axis, andthe lighter the extremum of the arrow, the higher the yield at this point. The observedcircles correspond to peaks in the PES (not shown here). left panel of Fig. 7 is a rawimage, while the right one is obtained after some inversion analysis. An approach such asVMI allows a simultaneous determination of PES together with PAD, which is extremelyinteresting. From the thus combined PES/PAD distribution (double differential, energy-and angle-resolved, cross section) it is then easy to recover PES or PAD separately byproper energy or angular integration. Still, mostly because of signal intensity, the doubledifferential cross section can rarely be used as a whole. Therefore, energy or angularintegration usually allow a simpler access to the data. It is also simpler and usuallymore quantitative to compare theory to experiments in simpler representations, whererather than the double differential cross section d σ/ d Ed Ω PES/PAD, one considers singly16ifferential PES d σ/ d E or PAD d σ/ dΩ cross sections. We shall thus explore now in moredetail integrated PES and PAD.2.3.1. Photoelectron spectroscopy
A PES typically results from a multiphoton ionization (MPI) mechanism (see Sec. 2.1).Electrons absorb a certain number of photons to reach the continuum and be emitted.They can absorb more than the number of photons required to reach the continuumthreshold, which leads to copies of the signal (although much reduced in intensity). Thekinetic energies of the emitted electrons are then directly related to the single electronenergies ε i of the initially occupied electron states i inside the cluster through the simplerelation : ε kin = ε i + ν (cid:126) ω las , (9)where ν is the number of photons involved in the process.Fig. 8 illustrates the principle of PES in terms of a scheme (left part), both in the caseof mono- and multi-photons. A hypothetical system with two accessible valence states is Fig. 8. Left : schematic view of photoelectron spectroscopy (PES), including multiphoton ionizationscenario for the most bound state. The sample system has two single electron states, (cid:15) < (cid:15) <
0. Theemission threshold is taken as the reference of zero energy, here the ionization potential IP= − (cid:15) . Themeasured kinetic energies of emitted electrons are then recorded from threshold on upwards involvinga varying number ν of photons which produces successive copies of the single electron energies, sepa-rated by the laser frequency. Right : experimental example of PES measurement for the Na − cluster(monophoton regime), obtained with photons of energy 4.02 eV [46]. considered (levels 1 and 2) whose electrons can reach the continuum via 1 or 2 photonabsorption. In the multiphoton case, the resulting PES displays copies of the original PES,separated by the laser frequency. The PES is furthermore illustrated on an experimentalexample (right part) from Na − , in the monophoton case. The case of anionic clusters isemblematic of one-photon PES. Indeed, in such clusters, valence electron states are littlebound so that they can easily be turned to continuum electrons according to Eq. (9) withone photon in the visible. These measurements basically provide a structural informationon the system. In the present case, the PES exhibits well resolved peaks associated tothe single electron states, as indicated in standard spectroscopic notation. One can note17hat the degeneracy of the 1 g state is split into a series of sub-peaks because of symmetrybreaking of the ionic configuration. Such one-photon measurements on anionic clusterswere thus already performed in the early 1990’s [48]. More recent measurements nowadaysallow one to access PES for neutral or even cationic clusters, as those from neutralfullerenes [16, 47] and from positively charged metal clusters [49].MPI, as already indicated in Eq. (9) with ν >
1, is also possible thanks to the highcoherence of laser pulses. The impact on PES will be discussed at length in Secs. 4.3.3and 4.3.4. Let us however give a few words here. For moderate laser intensities, theMPI maps in the PES further copies of the occupied electron spectrum with increasingkinetic energy, each copy separated by (cid:126) ω las . For larger intensities, the regular pattern ofcopies of the single electron spectrum is blurred because of large ionization affecting thespectrum itself. At even higher intensities the signal mostly becomes exponential withbasically no structure left [15, 50].Fig. 9 shows a typical example of a PES measurement, in this case performed oncationic species. The chosen material is sodium in which it is well known that electronic d e n s i t y o f s t a t e s -6 -5 -4 -3 -2 P E S [ a r b . u . ] E kin [eV] P E S [ a r b . u . ] Na + Na + Na Na
1f 2p 1g -6 -5 -4 -3 -2 2 8 18 20 34 40 -6 -5 -4 -3 -2 E kin (eV) PES ( a r b . u . ) PES ( a r b . u . ) D en s i t y o f s t a t e s spherical cluster exp. ADSIC Fig. 9. Lower and middle panels: Experimental photoelectron spectra for Na (bottom) and Na (middle) obtained by irradiation from an ArF − excimer laser of frequency (cid:126) ω las = 6.42 eV [49]. Verticallines in the middle panel : static Kohn-Sham single particle energies of Na calculated in ADSIC.Upper panel: Kohn-Sham density of states for a spherical neutral Na calculated in LDA [51]. shell closure leads to especially stable configurations [10]. In turn, the PES is expected todisplay the corresponding shell structure. The figure focuses on the region of 40 electrons(which corresponds to a shell closure). For comparison, the expected shell sequence,as computed in the Clemenger-Nilsson approach [11], is indicated in the upper panel.Note that in that case only the two least bound shells altogether containing 20 electronswere measured. The figure exhibits several interesting features. First the comparisonbetween Na and Na (which contains 41 electrons) very clearly points out theshell closure at 40 electrons with the appearance of one single electron in the 1g levelaround 5.2 eV. This also complies with the expected level sequence displayed in theupper panel. Finally, for the sake of completeness, we have also indicated the results18f a DFT calculation performed with the ADSIC correction (see section 3.2 for details).The agreement obtained without any adjustment is remarkable. It should nevertheless benoted that in that case the PES mostly provides a structural information on the systemby giving access to the sequence of energies of occupied single electron levels. As we shallsee below, PES, especially in the MPI regime, can also provide valuable information onthe dynamics of the electron cloud.As a first example of study of electron dynamics by PES, an example on which weshall come back later (see section 4.5.2), we consider the case of C irradiated by laserpulses of various fluences, but fixed pulse duration of 150 fs [47], see Fig. 10. At variance C , ω las =1.59 eV, Δ t=150 fs P ho t oe l e c t r on s pe c t r a [ a r b . un i t s ] Kinetic energy [eV] ns laser
Fig. 10. Photoelectron spectra of C irradiated by a laser of frequency of 1.59 eV, pulse duration of150 fs, and various fluences of (from top to bottom) 2.19, 1.84, 1.70, 1.56, 1.42, 1.27, and 1.13 J/cm .Adapted from [47]. with the spectroscopic character of the PES in Fig. 9, the PES presented here display analmost monotonous exponential shape with little structures on top. The latter structuresare interpreted as signals from single-photon ionization of Rydberg states [16]. The ex-ponential slope is explained as reflecting thermal electron emission [52]. In this picture,the energy deposited by the laser is concerted into thermal electron energy. Concludingon the nature of the energy conversion on the single basis of the PES is nevertheless a bitquestionable as exponential PES are also naturally obtained by considering higher andhigher MPI processes [50]. On the other hand, the experiments of [47] also measured thePAD of emitted electrons and clearly identified a strong isotropic component which mightindeed be associated to thermal emission. The interpretation of [47] is thus certainly tobe considered very seriously. We shall come back on that point in Sec. 4.5.2 when dis-cussing effects of dissipation on electronic observables in more detail. At present stage,it is sufficient to conclude that PES clearly opens the door to the analysis of electrondynamics. And that PAD offers for sure an invaluable complement to such studies (seenext section 2.3.2). 19inally, and before discussing PAD we would like to discuss another possible applica-tion of the PES now involving rather long time scales. Fig. 11 shows an example of atime-resolved photoelectron spectrum (TRPES) measured in (H O) − [53]. The irradia- Kinetic energy (eV)
Fig. 11. Time-resolved photoelectron spectra of (H O) − irradiated by a pump of 1.00 eV and a probeof 1.57 eV. The spectrum below 1 eV has been multiplied by 15. See text for the explanation of thefeatures A, B, C and D. Inset: integrated intensity of features B and D as a function of pump-probedelay. From [53]. tion process is performed with a pump-and-probe setup of laser frequencies of 1 and 1.57eV respectively, and similar intensities (50-100 µ J per pulse). The PES exhibit a cleardependence on the delay between the pump and the probe. The four major structures areindicated by capital letters. The low energy structure (A) is associated to excited-stateautodetachment, while direct probe detachment from the ground state (B) is observedaround 0.25 eV. Structure around 0.6 eV (C) is attributed to resonant two-photon de-tachment from the pump , and finally transient excited-state signal (pump-probe, D)appears in the 1.00–1.50 eV kinetic energy range. Integrating the intensities of thesestructures provides the associated population dynamics, which is indicated in the insetof Fig. 11 for structures (B) and (D). Both exhibit a similar decay time. To summarize,the above result clearly shows that a TRPES provides an extremely rich tool of inves-tigation of details of electron dynamics. Such measurements, possibly complemented bytheoretical investigations, should thus help to reveal crucial information on irradiationscenarios. Even more so PAD bring an invaluable complement to PES, as we shall see inthe next section.2.3.2.
Photo-Angular Distributions (PAD) and PES/PAD
Photo-Angular Distributions bring an invaluable complement to Photo-Electon Spec-tra. An example of PAD is shown in Fig. 12, adapted from [17]. The PAD are plotted as a20 lectron kinetic energy (arb. units) Electron kinetic energy (eV)
Fig. 12. Top left : schematic view of various types of photoemission distributions (from left to right :oblate, isotropic, and prolate). Bottom left : ideal photoangular distributions corresponding to each case,with the respective value of the anisotropy parameter β . Right : experimental PES (top) and PES/PAD(bottom) of Na − cluster, irradiated by linearly polarized pulses from a dye laser (pulse width of about10 ns, peak intensity below 10 W/cm and photon energy 4.02 eV). Adapted from [17]. function of electron kinetic energies, so that they in fact represent a combined PES/PAD.The mere PAD can then be obtained by integrating over kinetic energies. It should beimmediately noted that the notion of PAD requires a proper definition of a referenceframe. The reference direction is given by the laser polarisation axis and angular distri-butions are thus measured with respect to this axis. But it should also be noted that,in the gas phase, the actual orientation of clusters or molecules with respect to this po-larisation axis is unknown so that one has at best access to an orientation averaged (ofthe molecule with respect to the laser polarisation) signal. This in particular reduces theangular distribution to a dependence on the angle between the laser polarization axisand the detection angle, because of angular averaging around the polarization axis. Forthen, in the case of single photon absorption, the cross section takes the simple form :d σ dΩ ∝ β P (cos θ ) , (10)where θ is the direction of the emitted electron measured with respect to the laser po-larisation, P is the second order Legendre parameter and β is known as the anisotropyparameter. In the simple case of one photon processes, the angular distribution is thusfully characterised by the anisotropy parameter β which takes values between -1 and2. Three values of β are thus special, as illustrated in the left part of Fig. 12 : β =2corresponds to a prolate-like form of the (orientation averaged) emission cloud along thelaser polarisation, so that signal will gather around 0 ◦ and 180 ◦ ; β =-1 corresponds to apurely transverse emission, oblate-like shape, with signal gathering around 90 ◦ and 270 ◦ ;finally β =0 corresponds to a fully isotropic emission.A realistic measurement is shown in the right part of Fig. 12 to complement theschematic part. The measurement has again been performed in Na − , thus complement-21ng the PES example of Fig. 8. We nevertheless indicate the latter PES for completenessand to ease the explanation of the features of the combined PES/PAD. As already notedin the discussion of Fig. 8, the case demonstrates a clear dependence of the photoemissionon the nature of the electronic wave functions (indicated with spectroscopic notationsin the figure). Comparing the PES (upper right panel) to the PES/PAD (lower rightpanel), one can see that the 2 p and 2 d electrons are emitted parallel to the laser polar-ization. On the contrary, the emission from all the 1 g states occurs preferentially alignedin the transverse direction. This demonstrates that PAD certainly adds further usefulinformation on the spatial structure of the emitting states.Another example of PES/PAD, this time in the standard VMI presentation, is shownin Fig. 13 in the case of C − [54]. In this representation, the VMI provides a polar image (a) (b) (c) (d) t <50 ns 90< t <150 ns C - Fig. 13. Velocity map images (left panels) of C − , and corresponding photoelectron spectra (rightpanels), after irradiation by laser pulses with frequency of 4.025 eV and duration of a few ns. Top row :yields accumulated for t <
50 ns. Bottow row : yields accumulated for 90 < t <
150 ns. Adaptedfrom [54]. of the directions (angle) and kinetic energies (radius) of the emitted electrons, again witha well defined reference axis provided by the laser polarization. The example of Fig. 13 isfurthermore time-resolved, or at least allows one to separate well separated time scales ofemission. Irradiation is performed with photons of frequency 4.025 eV and pulse durationsin the ns range. Panel (a) provides an image of photoelectrons emitted during the first 50ns after excitation, and the corresponding PES is plotted in panel (c). The PES exhibittwo maxima, also visible as rings in panel (a), one at high energy and one at low energy.22he high energy signal is associated to direct emission from the photo-excitation itself.The low energy component is attributed to thermionic emission in which the originallaser energy has been partly equipartitioned between vibronic degrees of freedom of thecluster, prior to electron emission. The time scale and the typical energies associated tothermoionic emission are thus much larger than the ones associated to direct emission.In the present experiment the typical time scales of thermoionic emission lie in the tensto hundreds of ns, which in that case can be identified experimentally. The scenario isconfirmed in panels (b) and (d) which present the VMI and the PES recorded in a latetime window, that is between 90 and 150 ns. The PES is now fully concentrated atlow energies, with no sign of a high energy, direct emission, component. This confirmsthe thermoionic nature of this late, low energy, emission. Therefore, even at a coarsetime level, such an analysis exemplifies the capabilities of PES and PES/PAD to analyseelectron dynamics in detail. We shall come back on those aspects later, see in particularSecs. 4.3, 4.4 and 4.5.
3. Theoretical approaches
Many-particle systems such as molecules or clusters, are highly correlated, and exactcalculations of their properties are extremely involved, mostly beyond feasibility for fi-nite systems without major symmetries. The main issue concerns here the treatment ofelectrons. Except for some specific cases, ions can be treated as classical particles. Thiswill always be the case in the following. To deal with the electronic problem, a varietyof approaches has been developed, each one being a compromise between precision andexpense. In this section, we present the most widely used schemes, paying a particularattention to density-functional theory (DFT) which is one of the most efficient tools incluster dynamics.Before going into details of the theoretical treatment, we schematically summarizein Fig. 14 the most widely used theoretical approaches and sketch the regimes of theirapplicability in a plane of excitation energy and particle number. The boundaries ofthe regimes are to be understood as very soft with large zones of overlap between themodels because the choice of a method also depends on several aspects (e.g., demandon precision, material, time span of simulation). The most elaborate models are the“ab-initio” methods which deal in a systematic manner with a Hamiltonian as exact aspossible. The simplest example is the Hartree-Fock approximation which, however, missesthe crucial electronic correlations. A typical example of the more elaborate approachesis the Configuration Interaction (CI) method which relies on an expansion of the exactmany-body wave function into a superposition of Slater states [55, 56]. The limitations forCI (and other ab-initio methods) are purely a matter of practicability. The limitation isnevertheless even more severe for dynamical applications of such theories, which are thuspresently restricted to rather small system sizes and small excitation energies. The rangeof applicability will slowly grow with the steadily increasing computer power. Densityfunctional theory (DFT) describes a system effectively in terms of a set of single-electronstates (see Sec. 3.1.2). It is limited in system size for practical reasons and in excitationenergy for physical reasons, because of the missing dynamical correlations from electron-electron collisions. Nevertheless, DFT and even more so its time dependent extensionTDDFT (especially when realized in full real time) nevertheless provide a most robust and23 ig. 14. Schematic view of applicability of different approaches (see text for details) in a landscape ofsystem size versus excitation energy per atom. The excitation energy can be loosely related to typicallaser intensities in the optical range. This is indicated by the intensity scales on top, which are, however,also strongly dependent on the response of the particular system, i.e. resonant or non-resonant. versatile tool in the field. A semi-classical mean-field description is provided by the Vlasovequation originally designed for plasma physics [57]. This approach ignores quantumeffects such as shell structure or tunneling and thus becomes questionable at low energies.It is furthermore reasonably tuned to metal electrons because of their ressemblance toa Fermi gas, but more difficult to apply in other materials, especially in covalent boundsystems. On the other hand, the semi-classical treatment allows one to include dynamicalcorrelations due to electron-electron collisions, leading to the Vlasov-Uehling-Uhlenbeck(VUU) approach [58, 59, 60], which extends the applicability to larger energies thanthose allowed by TDLDA. Even higher excitations and system sizes are the realm ofelectronic Molecular Dynamics approaches and rate equations which, however, are evenmore limited than VUU for low energies and small systems [61]. The upper limit inenergy is given by the onset of the relativistic regime, where retardation effects withinthe coupling begin to severely influence the dynamics.In the following, we shall use real-time TDDFT as the basic theory to describe ioniza-tion dynamics. We shall occasionally use VUU in order to discuss electronic temperatureeffects, as observed in some experiments. We thus briefly describe in this section ba-sics of TDDFT and practical implementations thereof. We discuss in some detail the24elf-interaction correction strategy to be developed to properly account for ionizationin a dynamical way within standard approximations of DFT. We next present in detailthe tools developed to access PES and PAD in TDDFT. We in particular discuss thedemanding inclusion of orientation effects of the irradiated clusters and molecules withrespect to laser polarization. We finally remind basics of VUU for completeness.3.1.
Basic formalism
Handling of the ionic background
The interaction between the ions in a cluster and the electrons is usually describedby pseudopotentials. This allows one to eliminate the inert, deep lying electron statesaround each ion and to concentrate on the relevant valence electrons. For a detaileddiscussion of pseudopotentials see, e.g., [62]. We go here a pragmatic way and take pub-lished pseudopotentials. For simple metals, we consider the soft, local pseudopotentialsof [63]. In more general cases, we employ mostly the local and non-local pseudopotentialsin separable form as introduced in [64]. More precisely, for each type of atom with Z valence electrons, we use the pseudopotential V PsP of the following form : V PsP ( r ) ϕ j ( r ) = V loc ( r ) ϕ j ( r ) + (cid:90) d r (cid:48) V nloc ( r , r (cid:48) ) ϕ j ( r (cid:48) ) , (11a) V loc ( r ) = − Zr erf (cid:16) x/ √ (cid:17) + exp (cid:0) − x / (cid:1) (cid:2) C + C x (cid:3) , x = rr loc , (11b) V nloc ( r , r (cid:48) ) = p ( r ) h p ( r (cid:48) ) , (11c) p ( r ) = √ r nloc3 / (cid:112) Γ (3 /
2) exp (cid:18) − r r nloc2 (cid:19) . (11d)Here ϕ j denotes the wave function of state j , erf the error function, Γ the Gammafunction, and x = r/r loc . The refitted parameters are C , C , r loc = r loc , and h . Thestandard parameters are given in [64]. However, for the results presented in Section 4,we use often refitted parameters which employ larger radii r loc and r nloc , thus softerpseudopotentials for more robust numerical handling, see Section 3.1.3.There also exists a commonly used alternative to pseudopotentials for metallic systems.In particular, simple metals have valence electrons with long mean free path throughout.The fine details of the ionic background can thus be seen by the electrons only in anaverage manner. This motivates the jellium approximation in which the ionic backgroundis smeared out to a constant positive background charge. This is a standard approachin bulk metals [65]. It has been generalized to finite clusters. In its simplest form, onecarves from the bulk a finite, homogeneously and positively charged sphere of radius R = r N / , whose total ion charge reproduces the given ionic charge eN ion . A moreflexible approach is achieved when allowing for a finite surface width, yielding the softjellium model ρ jel ( r ) = 34 πr s (cid:20) (cid:18) | r | − R jel σ jel (cid:19)(cid:21) − , (12)with R jel being defined by normalization to the total particle number (cid:82) d r ρ jel = N ion . The central density reproduces the bulk density ρ = 3 / (4 πr s ). The parameter σ jel ac-counts for a smooth surface transition from ρ to zero. The surface width (transition25rom 90% to 10% bulk density) is about 4 σ jel . The model can be extended to also de-scribe deformations which can have a considerable influence in metal cluster spectroscopydepending on the system [66, 67].3.1.2. Density Functional Theory and its Time-Dependent version
The goal of DFT is to develop self-consistent equations which employ effective po-tentials for the contributions from exchange and correlation. These potentials are to beexpressed in terms of the total local electron density ρ ( r ) of the system. The success ofDFT depends on a diligent choice of these effective potentials. For the brief review ofDFT, we take here a practitioners approach and discuss the Kohn-Sham (KS) schemefrom a given energy functional. We do not address the theoretical foundations of DFTin terms of the much celebrated Hohenberg-Kohn theorem [68] and Kohn-Sham formal-ism [69]. The many aspects of foundation and derivation can be found, e.g., in [70, 71, 72].3.1.2.1. The energy functional
The starting point is an energy functional for the totalelectronic energy E total , el . In the Kohn-Sham (KS) scheme, one represents the N (va-lence) electrons, by N non-interacting Kohn-Sham (KS) orbitals (or s.p. states) ϕ i ( r ), i ∈ { , . . . , N } . The total energy is then separated into kinetic energy (which then takesa simple form) and interaction energy (associated to the above mentioned effective pseu-dopotentials). The total electronic density is expressed from the KS orbitals as : ρ ( r ) = N (cid:88) i =1 (cid:12)(cid:12) ϕ i ( r ) (cid:12)(cid:12) . (13)Note that DFT schemes allow one to treat spin-up and spin-down density separately.For simplicity of presentation, we drop the spin dependence in the following. The totalelectronic energy is then composed as E total , el [ ρ ] = E kin ( { ϕ i } ) + E H [ ρ ] + E xc [ ρ ] + E coupl + E ext , (14a) E kin ( { ϕ i } ) = − (cid:126) m (cid:90) d r N (cid:88) i =1 ϕ ∗ i ( r ) ∇ ϕ i ( r ) , (14b) E H [ ρ ] = e (cid:90) (cid:90) d r d r (cid:48) ρ ( r ) ρ ( r (cid:48) ) | r − r (cid:48) | = 12 (cid:90) d r ρ ( r ) U H [ ρ ] , (14c) E coupl = (cid:90) d r N (cid:88) i =1 ϕ ∗ i ( r ) ˆ V coupl ϕ i ( r ) , (14d) E ext = (cid:90) d r ρ ( r ) U ext ( r ) . (14e)The kinetic energy is a functional of the s.p. orbitals ϕ i which serves to maintain thequantum shell structure in the KS calculations. The non-trivial correlation part of theexact kinetic energy is summarized in the interaction energy. The interacting term ismapped to the density functionals E H [ ρ ] + E xc [ ρ ]. The first term E H is the standard(direct) Coulomb Hartree energy, which naturally is a functional of ρ . We have introducedhere the notation U H for the corresponding Hartree potential. Conceptually simple are E coupl from the coupling to the ions ( ˆ V coupl is the potential operator built from the26seudopotentials) and the energy E ext modeling an external electromagnetic field U ext ( r ).Both these contributions couple to single electrons and are naturally well represented byan independent particle picture in terms of ϕ i .Finally, there is the exchange-correlation energy E xc which accumulates all pieces ofthe exact energy not yet accounted for. This is the most problematic part in the scheme,since its functional expression is not exactly known. Many approximations thereof doexist, among which the simplest and most robust one is the Local Density Approxima-tion (LDA). The construction of LDA is simple. One computes the ground state of thehomogeneous electron gas as exactly as possible and obtains the exchange-correlationenergy per volume E xc /V = ρ (cid:15) xc ( ρ ). Here this energy is still a function of the (ho-mogeneous) density ρ . The crucial point is to allow now for an inhomogeneous, andtime-dependent if needed, density ρ ( r , t ) in that expression. It amounts to consideringthe energy as composed piecewise from an infinite electron gas of densities ρ ( r , t ) whichis a bold approximation. Nonetheless, LDA provides a robust description for a widevariety of systems. There is an enormous body of literature pondering successes andfailures, for a more detailed discussion, see e.g. [71]. Note that a functional dependingon ρ ( r , t ) employs the instantaneous density and thus excludes any memory effect. Thistime-dependent generalization is often called adiabatic LDA (ALDA). Again, we use inthe following the generic notation LDA.The validity of LDA depends very much on the system under consideration. One ofthe major problems is the self-interaction error : the single particle state ϕ i is includedin the density ρ , and thus contributes to the mean-field Hamiltonian ˆ h KS (see Eq. (15b)below) which acts on ϕ i . This yields a wrong asymptotics for the Coulomb mean field.For example, for a neutral system described in LDA, it decays exponentially at largedistances instead of ∝ e /r as it should. An attempt to reduce the self-interaction error isthe generalized-gradient approximation (GGA) which augments LDA with an additionaldependence on ∇ ρ [73, 74]. GGA yields a significant improvement in the computation ofatomic and molecular binding. For example, it lifts the description of dissociation energiesto a quantitative level. However, GGA does not fully remove the self-interaction error.Thus there are various attempts for further improvement as, e.g., adding kinetic termsto DFT [75]. Another line of development is to explicitly implement a Self-InteractionCorrection (SIC). This helps to deliver correct ionization properties, which is crucial indescribing PES and PAD dynamically. We will therefore discuss this approach in moredetail in Sec. 3.2.As indicated above, time-dependent DFT, effectively using ALDA, makes also an adi-abatic approximation. In order to account for dynamical effects, a Current DFT (CDFT)has been developed, which is based on LDA augmented by a dependence on electronic cur-rents evaluated in the linear response [76, 77, 78]. The response kernels in the extendedfunctional include memory effects and allow one to describe relaxation [79]. CDFT israther involved and thus there exist so far only applications to symmetry restricted sys-tems as, e.g., in solids [80]. The underlying linear response modeling makes CDFT anextension for low excitation energies and/or amplitudes. Real electron-electron collisionsbecome important for more energetic processes. These are often treated by a quantumgeneralization of the Boltzmann collision term. We will address this extension of DFT inSecs. 3.5 and 5.2. 27.1.2.2. The Kohn-Sham equations
The stationary KS equations are derived by vari-ation of the total energy with respect to the s.p. wave functions ϕ ∗ i , yielding :ˆ h KS [ ρ ] ϕ i ( r ) = ε i ϕ i ( r ) , (15a)ˆ h KS [ ρ ] = − (cid:126) ∇ m + U KS [ ρ ] + ˆ V coupl + U ext , (15b) U KS [ ρ ] = U H [ ρ ] + U xc [ ρ ] . (15c)The local and density-dependent Kohn-Sham potential U KS consists in the direct Coulombterm U H and the exchange-correlation potential, which is a standard functional derivative U xc = δE xc /δρ . Coupling potentials to ions and to the external field are trivially given.The time-dependent KS equations analogously read :i (cid:126) ∂ t ϕ i ( r , t ) = ˆ h KS [ ρ ] ϕ i ( r , t ) , (16)where ˆ h KS is composed in the same manner as above, provided that one replaces ρ ( r ) by ρ ( r , t ). This assumes an instantaneous adjustment of the total electronic density, althoughmemory effects can play in some cases an important role, especially in E xc [72].The stationary KS equations (15a) pose an eigenvalue problem. They provide theelectronic ground state of a system. This is a highly non-linear problem due to theself-consistent feedback of the local density in the KS hamiltonian. It is usually solvedby iterative techniques [14]. The time-dependent KS equations imply an initial valueproblem. The natural starting point is the ground state obtained from the stationaryKS equations. The time-dependent KS system can then be solved by standard methodsof first order differential equations [14]. We finally remind that we wrote spinless KSequations. One can easily include the electron spin in Eqs. (15). We refer the readerto [70, 71, 72] for more details.3.1.3. A few words on numerical implementation
A representation of the s.p. wave functions ϕ α ( r , t ) and the fields ρ ( r , t ) and U KS ( r , t )on a coordinate-space grid is strongly recommended if one aims at computing electronicemission properties. Conceptually straightforward is a Cartesian 3D grid with equallyspaced mesh points. This leaves two choices for the description of the kinetic energy,that is finite difference schemes [81, 82, 83, 84] or the Fourier definition exploiting theextremely efficient fast Fourier transformation (FFT) [85, 86, 14]. The Coulomb problemis solved either by iterative methods (e.g., successive over-relaxed iterations) in con-nection with finite-difference schemes or by Fourier techniques in case of FFT. In thelatter case, one can produce an exact solution on the grid by using a double grid (ineach direction) [87] or, somewhat faster, by eliminating the long-range terms by a sep-arate analytical treatment [88]. A fast scheme for the static solution is provided by thedamped gradient iteration [89]. The scheme for time evolution depends on the representa-tion of the kinetic energy. For finite-difference schemes, one typically uses a second orderpredictor-corrector with a Taylor expansion of the KS time-evolution operator while forFFT schemes, the time-splitting (also called ˆ T - ˆ V splitting) technique is preferable in con-nection with the FFT representation [90]. An extensive comparative study of the variousgriding and iteration techniques can be found in [91].Let us end with a few words in the case of symmetric systems. For instance, the jelliummodel allows a description in higher symmetry, that is a representation on an axial 2D28rid [92]. In the case of explicit ions in simple metal clusters, they can be described bysoft, local pseudopotentials. This allows an averaging over axial angle, leading to thecylindrically averaged pseudopotential scheme (CAPS) [93] which is extremely efficientand thus has been used in many explorative studies. The CAPS allows one to treat explicitionic structure in full 3D with pseudopotentials while keeping electrons with cylindricalsymmetry [14]. This turns out to be a very good approximation for metals and it is evenexact for linear molecules such as carbon chains [94]. It has even allowed to step forth torather complex systems as, e.g., embedded clusters [95, 96]. The Fourier representationof kinetic energy cannot be applied in this geometry. Finite differences are the methodof choice in axial grids. The static solutions use the same iterative schemes as in 3D. Aparticularly suitable time-stepping scheme for axial 2D is the Peaceman-Rachford step,which is a separable version of the well known Crank-Nicholson step [97, 98].3.2. The self-interaction problem in DFT and TDDFT
As outlined above, LDA is plagued by the self-interaction error which is particularlyharmful for ionization properties. The safest way to deal with that is to introduce anexplicit Self-Interaction Correction (SIC). A conceptually simple and robust SIC wasintroduced by J. Perdew and A. Zunger in which all single-particle self-interactions aresubtracted from the DFT energy [99] : E SIC [ ρ ] = E LDA [ ρ ] − N (cid:88) α =1 E LDA [ ρ α ] , ρ α ( r ) = | ψ α ( r ) | , (17)where E LDA = E H + E xc is the LDA functional for the Coulomb-Hartree term as well asexchange and correlations. Note that we have changed here our standard notation for thesingle electron KS orbitals from ϕ i to ψ α . This is done on purpose as will become clearbelow in Eq. (20). The self-interaction corrected KS equations are then derived againby variation. A problem is that the emerging SIC-KS hamiltonian ˆ h SIC then becomesstate-dependent because δE LDA [ ρ α ] /δψ ∗ β ∝ δ αβ is state selective. We formulate this interms of a projector and obtain :ˆ h SIC = ˆ h LDA − (cid:88) α U α | ψ α (cid:105)(cid:104) ψ α | , U α = δE LDA [ ρ α ] δρ α . (18)It becomes apparent that the state dependence leads to a non-Hermitian SIC hamiltonian.This leads to a violation of orthonormality of the ψ α . To restore it, we have to add aconstraint (cid:80) αβ λ αβ (cid:104) ψ β | ψ α (cid:105) to the SIC energy with the (hermitian) Lagrange multiplier λ αβ . The SIC mean-field equations thus become :ˆ h SIC | ψ α (cid:105) = (cid:88) β λ βα | ψ β (cid:105) for the static case and (cid:16) ˆ h SIC − i (cid:126) ∂ t (cid:17) | ψ α (cid:105) = (cid:88) β λ βα | ψ β (cid:105) for the dynamic case. In both cases, these equations have to be complemented by the”symmetry condition” (cid:104) ψ β | U β − U α | ψ α (cid:105) = 0 (19)29hich is the crucial ingredient in the scheme stemming from the orthonormality con-straint [100, 101]. These SIC equations are hard to solve directly in the static case andnear to impossible in dynamics. The key to success is to introduce a second s.p. basisset { ϕ l , l = 1 , . . . , N } which is connected to the set { ψ α , α = 1 , . . . , N } by a unitarytransformation [102] | ϕ j (cid:105) = N (cid:88) α =1 u jα | ψ α (cid:105) (20)tuned to diagonalize the matrix of λ βα . This simplifies the static and dynamic SIC-KSequations to ˆ h SIC | ϕ j (cid:105) = ε j | ϕ j (cid:105) , (cid:16) ˆ h SIC − i (cid:126) ∂ t (cid:17) | ϕ j (cid:105) = 0 . (21)Eqs. (18–21) formulate the SIC problem in the “2setSIC” scheme. They are solved byinterlaced iterations. One performs a step of static or dynamics mean-field problem (21)and then adjusts the unitary transformation (20) to accomodate the symmetry condition(19). A detailed representation of the scheme can be found in [102, 100].Although we dispose with 2setSIC of a powerful technique to solve the static anddynamical equations with SIC, it remains a tedious task. There are several interestingsimplifications around. With the help of the optimized effective potential method (OEP)[103], one has developed implementations of SIC in terms of state-independent localpotentials, which lead to the Krieger-Li-Iafrate (KLI) method [104] and, one step simpler,to the Slater approximation to SIC [105] (for a detailed discussion in connection withclusters, see [106]). Metal clusters are special in the sense that their valence electronshave all very similar spatial extension and stay close in energy (see for instance the s.p.energies of Na in Fig. 1). This allows one to replace the detailed N s.p. densities ρ α in SIC by a single averaged representative ¯ ρ e = ρ/N which then defines the energyfunctional for Average-Density SIC (ADSIC) as : E ADSIC ( ρ ) = E LDA ( ρ ) − N E
LDA ( ρ/N ) . (22)ADSIC is simple, robust, and reliable. It provides the correct asymptotics of the KS fieldwhile it is formally as simple to handle as LDA. The correct asymptotics and simplicityrenders ADSIC very useful in calculating electron emission and its observables. Mostexamples in this article were computed with ADSIC. Although motivated by metal elec-trons, ADSIC also performs surprisingly well for covalent molecules, see [107, 108] andFig. 15. Within ADSIC, the concept of s.p. densities is not needed anymore such thatADSIC is also applicable to semi-classical schemes [106]. In fact, it was first proposed byFermi [109] in a semi-classical context.Fig. 15 demonstrates the effect of ADSIC for a couple of basic organic molecules. Ion-ization potentials (IP) presented in this figure have been computed using the energy ofthe Highest Occupied Molecular Orbital (HOMO) obtained from the ground state con-figuration of each system. The wrong asymptotic Coulomb Kohn-Sham potential of pureLDA leads to less binding and thus to much reduced IP. The deviation is uncomfortablylarge. Correcting the self-interaction error even with the simple ADSIC suffices to obtaina very satisfying reproduction of the experimental IP.Fig. 16 demonstrates the effect of SIC in more detail. As test case, we consider thecluster anion K − where the failure of LDA is particularly apparent. We use here aspherical soft jellium ionic background, see Eq.(12), with a Wigner-Seitz radius r s = 530 ig. 15. Ionization potentials calculated from the energy of the HOMO, for a selection of conjugatedmolecules. Compared are results from LDA and ADSIC with experimental data. Adapted from [107]. -3-2-101 -20 -10 0 10 20 U K S , ε i ( e V ) radial coordinate r (a ) LDA ADSIC s s p p K − Fig. 16. Kohn-Sham potential and single particles energies for K − described with spherical jelliumbackground. Compared are results from LDA (left) with those from ADSIC (right). a and surface parameter σ = 1 . . The cluster as a whole has a negative charge.Consequently, the total Coulomb potential as it is used in LDA has an asymptotics ∝ + e /r and produces a Coulomb barrier between inside and outside. ADSIC, on the otherhand, sees asymptotically the Coulomb potential of all electrons minus the one whichis departing. This is the potential of a neutral system which converges exponentially tozero from below. The different asymptotic potentials mostly lead to a global shift of the31.p. energies, while the energy differences between the occupied states are less affected.This global shift is particurlarly disastrous for this anion. The system is hardly boundwith LDA, while ADSIC produces comfortable and realistic binding, although weak.Finally, we consider an example of time-dependent SIC (TDSIC) solved in the 2setSICframework and directly analyze the ionization dynamics of a molecule. To that end, weuse a simple 1D model for a H dimer molecule with a smoothed Coulomb potential [110].It is well adapted for a consistent test of SIC [102]. The two electrons have aligned spins(triplet state) to make a non-trivial test of SIC. We work at the level of “exchange only”,so that the benchmark becomes time-dependent Hartree-Fock (computed with exactexchange). To test LDA consistently, a density functional for exchange has been developedfor this 1D model within LDA. This density functional is also used as a basis for SIC [102].As SIC has a large impact on the IP, we take the time evolution of ionization as a criticaltest. A result for an instantaneous boost is shown in Fig. 17. The failure of TDLDA Fig. 17. Time evolution of total ionization after an instantaneous boost for a 1D model of H in tripletstate, calculated in LDA (blue dots), Hartree-Fock (red dashes) and 2setSIC (black full curve). Adaptedfrom [102]. (blue dots) for this observable becomes obvious. The IP is grossly underestimated andconsequently, the ionization is too high. The 2setSIC (full line) cures the problem almostperfectly, as is visible in the excellent agreement with the Hartree-Fock calculation (reddashes).3.3. Total ionization, PES and PAD in TDDFT
In this section, we discuss detailed observables from direct electron emission. By directemission, we mean those processes which are caused without delay by the electronicexcitation process. They dominate at moderate excitations and short laser pulses withduration of some tens of fs where the competing process, that is thermalization andsubsequent thermal emission, is less important. At this short time scale, we can alsooften neglect explicit ionic motion. 32.3.1.
An example of PES and PAD as preview
As a preview, we show in Fig. 18 PES and PAD for a simple example, Na with sphericaljellium background. This system has only two occupied levels 1 s (twice degenerate) and1 p (six-fold degenerate) which simplifies the analysis. The left panel shows the PES, that -5 0 5 10 15 20 10 − − − − − − PE S ( a r b . un i t s ) E kin (eV)Na spherical jellium ε s = − . e V ω las = 12 .
24 eV ε p = − . e V ϑ (deg.) E kin = 6 .
74 eV β (1 s )2 = 2 .
00 0 30 60 90 120 150 18000.10.20.3 P A D ( a r b . un i t s ) Angle ϑ (deg.) E kin = 8 .
17 eV β (1 p )2 = 1 . s p β (1 s )2 P (cos ϑ ) 1 + β (1 p )2 P (cos ϑ )6.76.86.9 8.18.28.3 E k i n ( e V ) Fig. 18. Photoelectron spectrum (left), photoangular distributions (bottom middle and right) and com-bined PES/PAD (top middle and right) from Na with spherical jellium background, see Eq. (12), using r s = 3 .
65 a and σ = 1 a after excitation with a linearly polarized laser pulse, see Eq. (1), of intensity I = 6 . × W/cm , frequency ω las = 12 .
24 eV, and pulse duration T pulse = 90 fs. The two occupieds.p. states lie at ε s = − . ε p = − .
07 eV. The total ionization is N esc = 0 .
003 electron. Theangle ϑ is measured with respect to the laser polarization. is the distribution of kinetic energies of emitted electrons. Left of the vertical axis, thetwo originally occupied s.p. states are indicated. One photon adds 0.9 Ry energy andso places a peak at − . .
24 eV, or − .
07 + 12 .
24 eV respectively. The energy shiftby the photon is indicated by horizontal arrows. The peaks in the PES directly mapthe occupied states. Further 12.24 eV higher, one sees another two peaks. These are dueto two-photon processes moving the electrons 2 (cid:126) ω las up in energy. The two upper rightpanels of Fig. 18 show combined PES/PAD with an energy window around the first(middle) and second (right) peaks of the PES. The lower panels show the correspondingPAD for emission from the 1 s states (middle) and 1 p states (right). The PAD have avery simple structure. This test case with closed electron shells and spherical jelliumbackground is spherical throughout. The laser defines a preferred direction thus leavingaxial symmetry for the process. Therefore the PAD depend only on the angle ϑ relativeto the laser polarization axis. Moreover, they consist out of a constant contribution plus acos profile. We will see later on that this is the only possible structure for PAD from one-photon processes. The figure also indicates the anisotropy β , as defined later in Eq. (37),for each case. The PAD from a perfectly spherical 1 s state has the maximal possible value β = 2, which corresponds to strong alignment with the laser polarization and vanishingemission perpendicular to it. The less symmetrical 1 p states yield a somewhat loweranisotropy. This simple example already demonstrates the richness of PES and PAD. Italso shows that one needs get acqaueinted with technical details to better understandthe content and behavior of both PES and PAD.3.3.2. Absorbing boundary conditions
A grid representation naturally leads to reflecting or periodic boundary conditions.Reflection emerges for finite difference schemes. A representation of the kinetic energy33y complex Fourier transformation is associated with periodic boundary conditions whereflow leaving the box at one side is re-fed at the opposite side. Both can lead to artifactsif a sizable fraction of electronic flow hits the boundaries. There are several ways to solvethe problem. The conceptually simplest approach is to enhance the size of the numericalbox. However, this is not a realistic option as the expense grows cubically with thebox length in 3D and quadratically in 2D. Very recently, a multi-grid method has beenproposed [111] which renders the use of enlarged boxes feasible (although still at the edgeof present days computer capabilities). Perfect removal of escaping particles is achievedby or exact boundary conditions [112, 113] which, again, are not yet practicable in 3Dcalculations. Robust and efficient are absorbing boundary conditions by an especiallytailored imaginary potential [114] or by applying a mask function during time evolution[115]. The latter technique is particularly easy to implement and has been widely used inthe past. Its robustness and efficiency allow one to develop advanced analyzing techniqueson the grid as, e.g., the computation of PES and PAD [116]. A detailed description anddiscussion of this approach and its proper choice of numerical parameters is found in [117].In the following, we will only address the mask technique for absorbing bounds.Fig. 19 sketches the implementation of absorbing boundary conditions with computa-tion of PES and PAD on a coordinate space grid. Proper handling of electron emission
Fig. 19. Schematic view of a coordinate-space grid with absorbing bounds (ring zone), a sampling direc-tion for accumulating PAD, and measuring points r M for the PES. requires absorbing boundary conditions. These are indicated by the ring area in the fig-ure covering here 3 grid points in each direction (actual calculations typically use 6–8points.) The absorption is performed in each time step as :34 ( r , t ) −→ ˜ ϕ ( r , t + δt ) = ˆ U KS ( t + δt, t ) ϕ ( r , t ) , (23a) ϕ ( r , t + δt ) = M ( r ) ˜ ϕ ( r , t + δt ) , (23b) M ( r ) = | r | < R in , cos (cid:18) | r | − R in R out − R in π (cid:19) γ M for R in < | r | < R out , R out < | r | . (23c)First comes the standard step (23a) in terms of the TDLDA (or TDSIC) propagator ˆ U KS ,which yields the intermediate wave function ˜ ϕ ( r , t + δt ). This is followed by the actionin Eq. (23b) of the mask function M defined in Eq.(23c), which removes gradually anyamplitude towards the bounds. We use here a spherically symmetric mask. The sphericalprofile is helpful to minimize griding artifacts when computing angular distributions [118].The absorbing bounds steadily reduce the norm of the wave functions from the inner maskradius R in to the outer one R out . The mask technique is however not perfect. One willalways encounter a small amount of reflected flow, particularly for electrons with lowkinetic energy. One can minimize the back-flow by proper choice of the exponent γ M entering the mask profile, see Eq. (23c). This depends, however, on the actual numerics(number of absorbing points, size of time step), for a detailed discussion see [117]. Typicalvalues of γ M are of order 1 / Ionization
The first observable which can be computed using working absorbing boundaries is thetotal ionization, i.e. the number of escaping electrons N esc . This can be computed simplyfrom the, now decreasing, single-particle norms as : N esc ( t ) = N (cid:88) i =1 N esc ,i ( t ) , N esc ,i ( t ) = 1 − (cid:104) ϕ i ( t ) | ϕ i ( t ) (cid:105) . (24)This shows that we have access to even more than the mere net ionization. Indeed each N esc ,i yields the depletion of s.p. state i separately. Both, total ionization and detailedlevel depletion are very instructive observables, see e.g. Fig. 36.3.3.4. Photoemission angular distributions (PAD)
The angular distributions d σ/ dΩ( ϑ, φ ) are evaluated in angular segments labeled bythe azimuthal angle ϑ and the polar angle φ . The reference frame for these two angles isusually one axis, called the z axis, of the Cartesian 3D grid designed to be identical withthe laser polarization axis, for details see Sec. 3.4. We collect all probability which wasremoved by the absorption step (23b) and accumulate it. A straightforward collection ofgrid points in a segment tends to produce noisy results because the number of grid pointsper segment fluctuates. We therefore associate with each grid point a smoothing function W ( r ) which distributes the strength over a vicinity of order of grid spacing. This sufficesto produce acceptable smooth distributions. The PAD is thus computed as :35 ( ϑ, φ ) = N (cid:88) i =1 A ( i ) ( ϑ, φ ) , (25a) A ( i ) ( ϑ, φ ) = (cid:88) n ∈ abs.b.c. (cid:90) d r r W ( r e r − r n ) n esc ,i ( r n ) , (25b) W ( r ) = max(∆ x − | x | , x max(∆ y − | y | , y max(∆ z − | z | , z , (25c) n esc ,i ( r n ) = (cid:90) d t | ˜ ϕ i ( r n , t ) | [1 − M ( r n , t )] , (25d)where e r = (sin ϑ cos φ, sin ϑ sin φ, cos ϑ ) is the unit vector in the direction of the wantedangles. The smoothing is done by simple tent functions which comply with the integrationrule used in the normalization. The angular segments in Fig. 19 try to symbolize thissmoothing which collects (weighted) information in the vicinity of a ray. The above recipeapplies to state specific PAD A i as well as the total PAD A . Alternatively, one uses across-section-like notation d σ/ dΩ for PAD. The present choice is more flexible for thepresentation of orientation averaging, see Sec. 3.4.3.3.5. Photoemission spectra (PES)
The PES can be deduced from the temporal phase oscillations of the wave functions atmeasuring points r M close to the absorbing bounds. This technique had been introducedin [119, 120]. A detailed discussion of the method and extension to strong laser fields isfound in [121]. We summarize it here briefly. To explain the computation of PES, we firstconfine the considerations to 1D in order to keep things simple and extend it finally tothe general 3D case. The measuring point is thus denoted for a while z M . Starting pointis the solution of the Schr¨odinger equation for one electron in a laser field in velocity-gauge, see Eq. (4) in section 1.2.1. In this gauge, the electronic wave function ϕ ( v ) at thesampling point z M reads ϕ ( v ) ( z M , t ) = (cid:90) d k √ π e i kz M (cid:101) ϕ ( v )0 ( k ) e − i ω k t +i kδq − i δ Ω , (26a) ω k = k ←→ k = √ ω k , (26b) δq ( t ) = E (cid:90) t d t (cid:48) F ( t (cid:48) ) , (26c) δ Ω( t ) = E (cid:90) t d t (cid:48) F ( t (cid:48) ) , (26d)where F ( t ) is the time integrated laser pulse introduced in Eq. (5). It should also benoted that Eq. (26b) exploits the fact that z M is near the absorbing bounds, such thatonly outgoing waves with k > (cid:12)(cid:12)(cid:12) (cid:101) ϕ ( v )0 ( k ) (cid:12)(cid:12)(cid:12) from the sampled ϕ ( v ) ( z M , t ). This is straightforwardfor weak laser fields where we can neglect δq and δ Ω. For then, a time to frequencyFourier transformation of ϕ ( v ) , namely (cid:82) dt e i ωt ϕ ( v ) ( z M , t ), produces in the right-handside of Eq. (26a) a δ ( ω − ω k ) which, in turn, reduces the k integration to the point k = √ ω thus delivering the wanted (cid:101) ϕ ( v )0 ( √ ω ).36rying to directly apply this time-frequency transformation when the field strengthis not small runs into trouble due to the non-trivial time dependencies induced by thefactors δq ( t ) and δ Ω( t ) in Eqs. (26). A simple solution is to counter-weight the disturbingphase factors by a phase-correction factor e iΦ before the transformation. We thus consider (cid:90) d t √ π e i ωt − i √ ωδq +i δ Ω ϕ ( v ) ( z M , t ) = (cid:90) d k π (cid:90) d t e i kz M (cid:101) ϕ ( v )0 ( k ) e i( ω − ω k ) t − i( √ ω − k ) δq ≈ (cid:90) d k δ ( ω − ω k ) e i kz M (cid:101) ϕ ( v )0 ( k )= e i √ ωz M (cid:101) ϕ ( v )0 ( √ ω ) . We finally evaluate the PES yield Y M at measuring point z M by : Y M ( E kin ) ∝ (cid:12)(cid:12)(cid:12) (cid:101) ϕ ( v )0 ( √ ω ) (cid:12)(cid:12)(cid:12) = (cid:12)(cid:12)(cid:12)(cid:12)(cid:90) dt √ π e i ωt − i √ ωδq +i δ Ω ϕ ( v ) ( z M , t ) (cid:12)(cid:12)(cid:12)(cid:12) , ω ≡ E kin . The approximation here consists in assuming that the time integration, although com-plicated by the time profile in δq ( t ), will still deliver ω ≈ k / e − i( √ ω − k ) δq ≈
1. This approximation is valid anyway for weak fields. It extends the applicability ofthe time-frequency transformation to stronger fields. Fourier transformation of a phase-augmented wave function thus allows one to deduce the wanted momentum amplitude (cid:101) ϕ ( v )0 for a wide range of field strengths. The recipe may fail, however, for very strongfields where the temporal variation of ( √ ω − k ) δq dominates over ( ω − ω k ) t .For the extension to 3D, we have to take into account that there is a whole vector of k rather than just two directions k = ±√ ω . We exploit the fact that the 3D analyzingpoint r M is close to the absorbing bounds (no reflection) and sufficiently far from theemitting zone. Thus the prevailing outgoing momentum k at this point has the directionof r M , i.e. e k = k k = r M r M = e M . (27)The frequency-momentum relation Eq. (26b) is then generalized to k = + e M √ ω andthe PES yield at r M is identified as Y M ( E kin ) ∝ (cid:12)(cid:12)(cid:12)(cid:12)(cid:90) d t √ π e i ωt − i √ ωδq +i δ Ω ϕ ( v ) ( r M , t ) (cid:12)(cid:12)(cid:12)(cid:12) . (28)This is the straightforward 3D generalization of the 1D formula above. The phase Φ( t ) = − i √ ωδq ( t ) + i δ Ω( t ) is negligible in case of weak fields E (cid:28) (cid:112) ω/
2. It extends theapplicability of the method to stronger fields.A word is in order about the choice of gauge. The above evaluation of PES is formulatedin velocity gauge because the exact solution of the Schr¨odinger equation in the laser field ismuch simpler in this gauge, see Eq. (26). If one prefers to perform numerical calculationsin space gauge, one just has to apply the transformations Eq. (5) and Eq. (6) beforeusing Eq. (28).The effect of the phase correction in the recipe Eq. (28) is demonstrated in Fig. 20.The field strength for I = 3 × W/cm is obviously sufficiently small such that there ispractically no effect from the phase correction. The stronger field with I = 3 × W/cm clearly needs the phase correction and including it yields still a reliable PES. Results forsignificantly larger field strengths cannot be trusted.37 − − − − − − − PE S E kin (eV) I = 3 × W/cm Na , jellium 0 5 10 15 20 E kin (eV) I = 3 × W/cm rawphase-augmented Fig. 20. Ionization properties for Na with jellium background under the influence of a laser pulse havingfrequency ω las = 1 . T pulse = 12 fs, and intensity as indicated in each panel. Sphericalabsorbing bounds were used covering at least 16 grid points. The “phase-augmented” results (brownlines) are obtained from Eq. (28), while the “raw” results (light green lines) only use the time-frequencyFourier transform of ϕ ( v ) ( r M , t ), see Eq. (26a) in the 1D case. Alternative routes to compute PES/PAD
An alternative theoretical approach to evaluate PES/PAD has been proposed veryrecently in [111]. The principle is schematically explained in the left panel of Fig. 21.This method is very similar to ours in the sense that in region A , the TDDFT equations Fig. 21. Example of calculated combined PES/PAD. Left : 1D view of the mask function applied in region A where time-dependent Kohn-Sham equations are solved numerically in real-space, while wave functionsare propagated analytically in region B . Region C serves a a buffer region where wave functions, evaluatedin A , overlap with those evaluated in B [111]. Right : Experimental (top) combined PES/PAD [122] andtheoretical ones (bottom) [111] in logarithmic scale for N molecules in a 6 cycle infrared laser pulse( λ las = 750 nm, I = 4 . × W/cm ). The theoretical PES/PAD is averaged over four orientationsof the N molecule. M ( r )in region C . The new feature comes with region B where Volkov states are analyticallypropagated in momentum space to account for the interaction with the laser field. Moreprecisely, state i is described by the following wave functions ϕ i : ϕ A,i ( r , t ) = η A,i ( r , t ) + η B,i ( r , t ) , (29a)˜ ϕ B,i ( p , t ) = ˜ ξ A,i ( p , t ) + ˜ ξ B,i ( p , t ) , (29b)where η A,i ( r , t ) = M ( r ) U ( t, t (cid:48) ) ϕ A,i ( r , t (cid:48) ) , (30a) η B,i ( r , t ) = M ( r ) (cid:90) d p (2 π ) d/ e i p · r U V ( t, t (cid:48) ) ˜ ϕ B,i ( p , t ) , (30b)˜ η A,i ( p , t ) = (cid:90) d r (2 π ) d/ , e − i p · r [1 − M ( r )] U ( t, t (cid:48) ) ϕ A,i ( r, t (cid:48) ) , (30c)˜ ξ B,i ( p , t ) = U V ( t, t (cid:48) ) ˜ ϕ B,i ( p , t (cid:48) ) (cid:90) d r (2 π ) d/ e − i p · r η B,i ( r , t ) . (30d)In Eqs. (30a) and (30c), U ( t, t (cid:48) ) stands for the propagator from time t (cid:48) to t with the fullHamiltonian including external fields, while in Eqs. (30c) and (30d), only the laser fieldenters the Volkov time propagator U V ( t, t (cid:48) ). Finally the momentum distribution of thephotoelectrons is approximated to :d σ dΩ p (cid:39) (cid:88) i ˜ ϕ B,i ( p , t → ∞ ) . (31)This procedure has been applied to N irradiated by a 6 cycle laser pulse of wavelength750 nm and intensity of 4 . × W/cm , with an averaging over 4 orientations (0 ◦ ,30 ◦ , 60 ◦ , 90 ◦ , see right bottom panel of Fig. 21. The grid spacing in region A of radius R A = 35 a is 0.38 a , and the buffer region has a radius R C = 25 a . These results arecompared with experimental measurements [122] (see top right panel). The latter showthat photoelectrons are preferentially emitted parallel to the laser polarization axis (0and 180 ◦ ). Note that the signal at 0 ◦ is slightly different from that at 180 ◦ . This can beexplained by the fact that the laser pulse is so short that the symmetry along the laserpolarization axis is broken [111, 123]. The comparison with the theoretical PES/PAD isfairly good, especially at high kinetic energies. More discrepancies are observed at lowkinetic energies, probably because of the limited size of the numerical box.3.4. Orientation Averaging PAD (OAPAD)
The PAD obtained by Eqs. (25) determine the distribution for a fixed orientation ofthe molecule/cluster relative to the laser polarization axis. However measurements of freeclusters are usually done in the gas phase covering an isotropic distribution of clusterorientations. There are techniques for aligning molecules by strong laser pulses, for areview see [124] and for proposals using chains of pulses see [125, 126, 127]. Nonetheless,these have not yet been used in connection with measuring PAD on clusters. Thus wehave to perform orientation averaging of the TDDFT results to establish contact with39xisting measurements. Efficient techniques for orientation averaging of clusters havebeen developed in [128, 129]. We summarize them here briefly.3.4.1.
Direct averaging scheme
The cross-section detailed in Eqs. (25) stems from a TDDFT calculation with one fixedconfiguration in which the cluster orientation relative to laser polarization is known. Whatwe are looking for is the average of the cross-section over all possible cluster orientationswith equal weight. For its evaluation, we distinguish the laboratory frame and the clusterframe (in which all quantities are primed). The laboratory frame is defined by the laserpolarization axis, such that the polarization vector points along the z -axis of the 3DCartesian coordinate system, i.e. e pol = e z . The observed emission angles ( ϑ, φ ) aredefined with respect to this laboratory frame, where ϑ is the angle with respect to the z -axis and φ the angle in the x - y -plane. A cluster has three principle axes. The clusterorientation is defined by the Euler angles ( α, β, γ ) of these three axes with respect to thelaboratory frame [130]. Thus we deal, in fact, with an ensemble of PAD A ( i ) ( ϑ, φ ; α, β, γ )of the same cluster with different orientations. The orientation averaged one-photon PADfor emission from the s.p. state ϕ i then becomes A ( i ) ( ϑ, φ ) = (cid:90) d α dcos β d γ π A ( i ) ( ϑ, φ ; α, β, γ ) ≈ M (cid:88) m =1 λ m A ( i ) ( ϑ, φ ; α m , β m , γ m ) . (32)where the sum from 1 to M runs over a discrete set of points on which the integral willbe discretized (see below). Since spherical absorbing boundaries are used, the rotationby α about the laser axis can be done a posteriori and does not require any additionalTDDFT calculation. This leaves averaging over β and γ which is approximated by afinite-element representation of the integral, see rightmost part of Eq. (32). The chosenvalues for β m and γ m can be illustrated on a unit sphere. Figure 22 shows a samplingover 34 orientation points. The weight factors λ m are determined by cutting the surface
34 45 ◦ ◦ ◦ γβx yz Fig. 22. Example of a set of orientation points (black dots) used in direct orientation averaging. Thecluster is first rotated by γ about the z -axis (laser polarization) and then by β about the y -axis of thelaboratory frame. of the unit sphere into segments S m around unit direction e m of element m . The S m is40efined as the collection of points on the unit sphere which are closer to e m than to anyother e m (cid:48) . The summation weights are then simply the areas of S m divided by the areaof the whole sphere, i.e. λ m = Area( S m ) / (4 π ).3.4.2. OAPAD from one-photon processes in the perturbative regime
The direct averaging scheme is conceptually simple and can be applied in any situation.However, it requires a considerable amount of reference orientations (segments) to achievesufficiently reliable results. A great simplification can be worked out for one-photonprocesses in the perturbative regime. In that case, one can deduce, using perturbationtheory formally, that the rotated PAD can be represented as : A ( i ) ( ϑφ, αβγ ) = (cid:88) µµ (cid:48) ,lmm (cid:48) D (1) ∗ µ (cid:48) ( αβγ ) D (1) µ ( αβγ ) D ( l ) m (cid:48) m ( αβγ ) a ( i ) µµ (cid:48) ,lm (cid:48) Y lm ( ϑφ ) . (33)in terms of a couple of expansion coefficients a ( i ) µµ (cid:48) ,lm (cid:48) , rotation functions D ( l ) µµ (cid:48) [130],and spherical harmonics. In this form, the integration over Euler angles ( α, β, γ ) canbe worked out analytically (for details, see [128, 129]), yielding finally the OrientationAveraged PAD asd σ i dΩ ≡ A ( i ) ( ϑφ ) = C ( i )0 Y ( ϑφ ) + C ( i )2 Y ( ϑφ ) , (34) C ( i )0 = 13 (cid:88) µ a ( i ) µµ, , (35) C ( i )2 = (cid:88) µ a ( i ) µµ (cid:48) , µ − µ (cid:48) ( − µ × − µ µ (cid:48) µ − µ (cid:48) . (36)The corresponding total PAD is obtained by summing over the s.p. PAD, i.e. A = (cid:80) i A ( i ) . It remains to determine the expansion coefficients a ( i ) µµ (cid:48) ,lm (cid:48) . Fortunately, thereare only very few, that is three for l = 0 and six for l = 2. They can be determinedfrom the PAD in only six properly chosen orientations because the different l can beproduced from one PAD, see Eqs. (33). The six reference orientations should be chosensuch that Eqs. (33) can be solved in a stable manner as a linear equation for the a ( i ) µµ (cid:48) ,lm (cid:48) .A recommended set is e (1) = e x , e (2) = e y , e (3) = e z , e (4) = ( e x + e y ) / √ e (5) =( e x + e z ) / √
2, and e (6) = ( e y + e z ) / √ with detailed ionic structure which is, unlike the jellium case, not rotational invariant.The laser frequency is sufficiently high such that one-photon emission dominates. The leftpanel shows the PAD for a fixed cluster orientation where the cluster symmetry axis isaligned with the laser polarization axis. One sees a pronounced pattern, in particular thefour-fold structure of the ionic rings and the rotation of the upper ring by 45 ◦ relative tothe lower ring. The emission maxima are located at ϑ ≈ ◦ and ϑ ≈ ◦ , i.e. sidewardsto the laser polarization. The orientation averaged PAD becomes independent of φ , as itshould be, and the emission is forward/backward dominated with maxima at ϑ = 0 and ϑ = 180 ◦ . Altogether, the effect of orientation averaging is dramatic. Calculations for asingle orientation have thus little predictive value.41 .2x10 − ⇒ averaging φ ( ◦ )04590135180 ϑ ( ◦ ) φ ( ◦ ) Fig. 23. Total PAD A ( ϑ, φ ) for Na with explicit ionic background for one fixed orientation (left panel)and orientation-averaged (right panel), after laser excitation with frequency ω las = 7 . I = 10 W/cm , and pulse length with T pulse = 60 fs. The emission angles ϑ and φ are measured withrespect to the laser polarization. As already pointed out, the OAPAD only depend on the angle ϑ , and not on φ anymore.Moreover, they are symmetric with respect to the transformation ϑ ↔ − ϑ . Thus theycan be expanded in a standard manner in terms of even Legendre polynomials P l (cos ϑ )as d σ dΩ ∝ β P (cos ϑ ) + β P (cos ϑ ) + β P (cos ϑ ) + ... (37)Most important is β , called the anisotropy parameter. It is the only remaining parameterin the perturbative regime where all β > = 0. Note that we had already recognized thesimple 1 + β P (cos ϑ ) structure in Fig. 18.3.5. Thermal effects
The previous discussions were led assuming that electrons remain at zero temperature.Introducing an ionic temperature in the formalism raises no difficulty as the ions aretreated classically. An ionic temperature can thus either naturally build up in the courseof time, or be initially introduced by giving the proper velocities to ions. The system,of course, remains globally microcanonical, and the temperature is thus subject to thecorresponding fluctuations. The electronic temperature, in turn, is explicitly set to zeroby construction. TDLDA does not allow occupation numbers of KS orbitals to vary intime, which prevents an account of thermal effects at the side of the electrons. This iscertainly a major formal limitation of the formalism. We have seen that (possibly sizable)electronic temperatures have been observed experimentally (see for example Figure 10).These effects may play a significant role in the understanding of the dynamics of thesystem, and on the analysis of experimental results. A theoretical account of such thermaleffects is thus to be developed.A standard path to accommodate thermal effects in dynamical systems is to recurto kinetic theory. This theory was originally formulated in the framework of classicalmechanics with the Boltzmann equation as a prototypical example [131]. Quantum sys-tems add two more complications. The Pauli principle prevents collisions into occupiedstates. The uncertainty principle, in turn, raises difficulties when trying to treat colli-sions locally (as done in the Boltzmann equation). Fully quantal kinetic equations are42evertheless conceivable but much more involved than the classical ones [132, 133] andvery hard to apply in finite quantum systems with their discrete spectra. In practice,there is thus no practical quantum theory of such collisional correlations available yet.In turn, semi-classical approaches were developed over the years and allow to cover somesituations. The basics is then a version of the Boltzman equation, adapted to accountfor Pauli principle and known as the Boltzmann-Uehling-Ulhenbeck or Vlasov-Uehling-Ulhenbeck (VUU) equation. It was formally introduced in the early 1930’s [134], thenwidely used in nuclear dynamics [58, 135], and explored more recently in simple metalclusters [136, 137, 60]. By construction, the VUU equation is semi-classical, which meansthat details of quantum shell effects are washed out. It is thus applicable, at best, onlyin highly dissipative situations where the latter shell effects do indeed disappear, or inhomogeneous fermionic systems such as electrons in solids [138] or nucleons in a neutronstar [139]. In finite fermion systems like clusters and molecules, they are thus bound tohigh excitation energies, which certainly strongly limits their range of application.Another point is worth being mentioned here. In realistic calculations, and in spite ofthe semi-classical treatment, one wants to recover an acceptable description of groundstate properties of the studied systems. Experience shows that in electronic systems,this is practically viable only in metal clusters such as sodium clusters [137, 60]. Thisagain strongly restricts the range of applicability of such methods. In particular, it doesnot allow to attack such widely studied systems as C in a fully realistic manner. Still,results obtained in simple metals are interesting and demonstrate the importance of thisinclusion of collisional correlations. We will thus briefly present VUU and discuss someof these results here.There are various ways of introducing VUU but in the case of finite fermionic systems,the simplest and probably best founded presentation is to recur to a semi-classical ap-proximation on top of TDLDA. To perform such a semi-classical limit, we first recast thetime-dependent KS (TDKS) equations in a matrix form, by introducing the one-bodydensity matrix ˆ ρ KS associated to KS states, which reads in real space representation :ˆ ρ KS ( r , r (cid:48) ) = N (cid:88) i =1 ϕ ∗ i ( r (cid:48) ) ϕ i ( r ) . (38)This one-body density matrix fulfills the TDKS equations which, in a matrix form, read :i ∂ ˆ ρ KS ∂t = [ˆ h KS , ˆ ρ KS ] , (39)with ˆ h KS given by Eq. (15b). The semi-classical limit can then be attained by performinga Wigner transform (or, even better, an Husimi transform [140]) of Eqs. (38) and (39).This leads to the introduction of the phase space distribution f ( r , p , t ) which then fulfills,at lowest order in (cid:126) , the Vlasov equation : ∂f∂t = { h KS ( r , p , t ) , f ( r , p , t ) } , (40)where we have introduced Poisson brackets. At the LDA level, for which the xc functionalis local in density, and for simple metals, for which the pseudopotential is local as well,the semi-classical KS hamiltonian takes the following simple form : h KS ( r , p , t ) = p m + U KS [ ρ ] + V coupl ( r ) + U ext ( r , t ) , (41)43here the local density is now obtained from momentum space integration of f ( r , p , t ): ρ ( r ) = (cid:90) d p (2 π (cid:126) ) f ( r , p , t ) . (42)The Vlasov(-LDA) equation can now be complemented by the effect of ”two-body”collisions as usually done in kinetic theory. This leads to the Vlasov-Uehling-Ulhenbeck(VUU) equation : ∂f∂t = { h KS ( r , p , t ) , f ( r , p , t ) } + I UU [ r , p ] , (43)obtained as the Vlasov-LDA equation complemented by a UU collision term : I UU [ r , p ] = (cid:90) dΩ d p | p − p | m d σ dΩ (cid:104) f p (cid:48) f p (cid:48) (1 − ˜ f p )(1 − ˜ f p ) − f p f p (1 − ˜ f p (cid:48) )(1 − ˜ f p (cid:48) ) (cid:105) . (44)which exhibits a local gain-loss balance for elastic electron scattering ( p , p ) ↔ ( p (cid:48) , p (cid:48) ).The associated cross-section for electron-electron ”two-body” collisions is d σ/ dΩ anddepends on the relative momentum | p − p | and possibly on the scattering angle. Wehave furthermore used the shorthand notations ˜ f p = 2 π (cid:126) f p / f p = f ( r , p , t ),the collision term being local in space and time. To avoid a double counting with themean field h KS the differential cross-section can be evaluated with a screened Coulombinteraction following standard scattering theory [136, 137, 60, 141].The Vlasov and VUU equations are best solved using test particle methods ratherthan grid methods. Practically, this amounts to represent the phase space distribution f ( r , p , t ) by a swarm of numerical classical particles { ( r i , p i ) , i = 1 , ... } , each one with agiven weight ω and following classical equations of motion with forces derived from themean field hamiltonian h KS ( r , p , t ) [142]. Finally, the coupling to classical ionic motion,similarly as in quantal TDLDA, is performed, leading altogether to a Vlasov/VUU-LDA-MD approach. The theory was introduced in cluster physics in the late 1990’s [136] andfurther used since then [60], but, as already mentioned, only applied to simple metalclusters.The analysis of Vlasov or VUU dynamics proceeds in a rather simple manner as theelementary degrees of freedom are classical test particles. The total ionization is directlygiven by the number of test particles ( × ω ) outside a large given box in which the Coulombfield is computed. It is equivalently obtained as the integral of the semi-classical densityEq. (42) over this computational box. The latter density also enters the expression of thedipole moment (see Eq. (45) in the next section) to provide an analysis of optical response.The distribution of kinetic energies and the angular distribution are directly extractedfrom the velocities of the ”emitted” test particles. Note that the semi-classical nature ofthe approach makes the PES always exponentially decreasing, whatever the excitation, atvariance with the quantum case which exhibits the electronic single particles energies. Thequantum/classical comparison is then meaningful only at sufficiently large excitations forwhich the quantum PES is also exponentially decreasing with no more shell structure.The PAD, in turn, can be compared in a meaningful way whatever the excitation energy.For a detailed discussion, see Sec. 4.5.2. 44 . Illustrative results The coordinate-space and real-time technique to solve TDDFT as presented in Sec. 3offers a powerful tool to describe a broad range of scenarios for the dynamics of clustersand molecules. We present in this section a few illustrative examples of such theoreticalstudies thereby concentrating on electron ionization and related observables. As a starter,Fig. 24 shows a calculated combined PES/PAD (left column) for C , compared withan experimental one (right column). The experimental spectrum was recorded at the Fig. 24. Top : Combined PES/PAD of C obtained at ω las = 20 eV given by (a) TDLDA-ADSICcalculations ( I = 7 . × W/cm and duration of 60 fs, C radius of 6.44 a ) and (b) experimentalmeasurements using synchrotron radiation with duration of about 1 ps. Bottom : the correspondingvelocity map image of (a) is presented in panel (c), and that of (b) in panel (d). Maxlab Synchrotron facility, using an oven to produce the C molecular beam [143].Experiments and calculations are both performed at a photon energy of 20 eV. The toppanels present the angle- and energy-resolved distributions of photoelectrons for threeselected energies close to the HOMO level in a way similar to Fig. 18, while the bottompanels show the distribution in a polar representation called velocity map images (VMI),see introductory discussion in Sec. 2.3.2 and Fig. 7. The upper panel shows nice agreementbetween experimental and theoretical results. The VMI in the lower panels also nicelyagree for larger energies near the IP (outer part of the circle). But they differ in thelow-energy region because it is filled with electrons from thermal emission not accountedfor in TDDFT, see also Fig. 36 later on. However, although the physical content of such45MI is very rich, one cannot easily read off a quantitative comparison from such a figure.Therefore, we will in the following sections preferably discuss integrated VMI : PESobtained by integrating over the angle, PAD obtained by integrating over the kineticenergy, and total ionization obtained by integrating over both, energy and angle.As already mentioned in our first example in Sec. 3.3 and Fig. 18, the peaks observedin a PES are fingerprints of the single particle (s.p.) energies of the electrons before theyhad been ejected by the photon. However, s.p. energies are usually not well described byLDA. It has to be complemented by a self-interaction correction (SIC) to attain realisticionization potentials (IP) and thus ionization properties (see Sec. 3.2). Before presentingdetailed results on PAD and PES, we first illustrate in Sec. 4.1.1 the capabilities of SIC toproperly describe ionization dynamics. We remind again that we describe the processesfully dynamically (see the methodology of Sec. 3.3.5). This means, e.g., that discussingPES does not amount to a mere comparison of a computed static s.p. spectrum with ameasured PES. Since our access to PES, PAD, and ionization is fully dynamical, it isthus applicable to any dynamical regime and free of any adjustable parameter. This iswhy having a proper IP is a crucial step in our dynamical description, especially at lowenergy where ionization occurs mostly close to threshold, whence the importance of aSIC. Even when dealing with appropriate IP, a word of caution is to be added concerningthe interpretation of PES as map of s.p. spectra. There are cases where one finds slightdeviations for deep lying s.p. states. These can be understood as correction from finalstate interaction [144]. They should be considered for a high precision analysis of data.We will ignore the effect in the more principle considerations of this section.This section is thus organized as follows. In Sec. 4.1, we will first demonstrate theimpact of SIC on structural properties as single-electron energies and optical response.In Sec. 4.2, we will discuss ionization as a signal from laser-induced dynamics. Then,in Sec. 4.3, we will address PES, and in Sec. 4.4 PAD, both in the one-photon and inthe multi-photon regime. We finally discuss the impact of temperature, either ionic orelectronic, on PES and PAD in Sec. 4.5.4.1. Impact of the self-interaction correction on electronic emission
As just mentioned, ionization properties are very sensitive to the s.p. energies, whencethe importance of SIC. We will here demonstrate the impact of SIC on electronic prop-erties, such as s.p. spectrum or optical response, and finally discuss a typical example ofPES.4.1.1.
Ionization potentials and single electron spectra
Before starting a discussion of the ionization potential (IP), we have to specify that inmore detail. For a system with N electrons, the IP is defined by the difference I adia = E ( N − − E ( N ). This is, in fact, called the “adiabatic IP” if the final values of E ( N − I ∆ = E fix ( N − − E ( N ) which is obtained from the energy E fix of the ionized system while still maintaining the ions in their original configuration.This vertical IP is a purely electronic observable which renders it very instructive. Inpractice, induced ionization processes (laser pulse, ion collision) are so fast that the ionic46onfiguration is almost inert during the emission process. This motivates the use of thevertical IP, henceforth called simply “the IP”.As already discussed in Sec. 3.2, correct s.p. energies are strongly related to the correctasymptotic behavior of the KS potential, see the example of K − in Fig. 16. And as justexplained above, the IP is defined as difference of two stationary energies. However,this definition is not well suited to dynamical calculations where one aims at trackingionization ”on the fly”. For then, it becomes crucial to fulfill Koopmans’ theorem [145]which states that the vertical IP should be identical with the s.p. energy of the lastbound electron (HOMO level), i.e. I ε = − ε HOMO . Koopmans’ theorem is violated inLDA and usually recovered when invoking SIC. This rules out LDA for a fully dynamicalapproach to ionization and calls for some SIC. It is thus crucial to test the performanceof TDDFT in this respect. As discussed in section 3.2, we have basically a full SIC(practically implemented via the 2setSIC scheme and simply denoted in the following bySIC) and ADSIC at our disposal. ADSIC is orders of magnitude simpler than SIC andthus certainly worth being considered very seriously. As a first step, we shall thus checkthe capabilities of both approaches with respect to reproducing IP’s along the lines of alarge systematic study in [108].Fig. 25 depicts the difference between the calculated IP and the experimental one,for a selection of molecules. Both ways are used for the calculations : from the HOMOlevel as − ε HOMO (upper panel, open symbols) or from the difference of binding energies I ∆ (middle panel and filled symbols). The upper panel shows that, as expected, LDAperforms badly when one considers the IP I (cid:15) from the energy of the HOMO. SIC and evenmore so ADSIC come much closer to the experimental IP. The middle panel comparesthe IP I ∆ from energy differences. Here we see better agreement for all three methods,demonstrating that the I ∆ is the more robust definition. To emphasize the discrepancybetween both estimates of IP’s, we plot their difference in the lower panel of Fig. 25.Vanishing difference signifies fulfillment of Koopmans’ theorem. LDA produces largeerrors while SIC and ADSIC do well. It is a bit of a surprise that the simpler ADSICoften performs better than the more elaborate SIC.By construction, ADSIC is a priori well suited to systems with metallic binding [106].It was nevertheless soon realized that it also performs well in covalent systems [107]. Thevery systematic study of [108] led to the unexpected result that for an enormous rangeof atoms, molecules, carbon chains, and fullerenes ADSIC leads very often to smallernon-Koopmans errors than SIC. It was also shown that, when comparing theoretical IP’sto experimental ones, again ADSIC was providing the best results for a huge range ofmolecules. That however does not mean that ADSIC is the ultimate solution to the self-interaction problem. We have already mentioned its intrinsic limitations (essentially dueto the fact that the functional explicitly depends on the number of electrons) in Sec. 3.2and we should add here that the scaling properties of ADSIC with increasing system sizeare also raising some problems [108]. Furthermore, there are a few specific cases whichraise difficulties. These are molecules where very different types of bonding coexist, suchas covalent and metallic bonding. The point is illustrated in Fig. 26 which presents thes.p. energies of the two complexes NaH O and Na(H O) , calculated in LDA, ADSIC andSIC. The right part of the figure complements the picture by showing the correspondingnon-Koopmans errors. Not surprisingly the water molecule with one prevailing bondingtype is well described in ADSIC, while SIC yields a small but somewhat larger non-Koopmans’ error. The situation becomes different in the mixed complexes Na(H O) n ,47 − − H L i H C H NH O H H O H F L i L i F B e C H C H H C N C O N N O O F P C l I ∆ − I ε ( e V ) − ∆ I ∆ ( e V ) − − − ∆ I ε ( e V ) ADSIC SIC LDA
Fig. 25. Calculated ionization potentials (IP) of a selection of molecules with different bonding types forvarious level of SIC (see text for details). Top : difference ∆ I ε = I ε − I exp of the IP I ε deduced fromthe energy of the HOMO with the experimental (vertical) IP I exp . Middle : difference ∆ I ∆ = I ∆ − I exp for the IP I ∆ calculated as the difference of binding energies. Lower : non-Koopmans’ error I ∆ − I ε .Adapted from [108]. -35-30-25-20-15-10-50 SIC ADSIC LDA s . p . e n e r g i e s ( e V ) NaH O SIC ADSIC LDA
Na(H O) H O NaH O Na(H O) -20246 non - K oop m a n s e rr o r( e V ) LDAADSICSIC
Fig. 26. Single particle energies of NaH O (left) and Na(H O) (middle) in LDA, ADSIC and SIC, forspins up (full lines) and down (dashes). Right : non-Koopmans’ error calculated in various levels of DFTas indicated. From [146]. where both ADSIC and LDA exhibit similarly large errors (although with opposite signs),while the Koopmans’ theorem is perfectly fulfilled in SIC. Thus ADSIC should not beused for Na(H O) n . This may not be a total surprise if one reminds that ADSIC appliesthe same correction to every orbital while, in such a mix of bonding types, one encountersboth highly delocalized metallic orbitals and much more localized covalent orbitals. Asvisible in the left panel of Fig. 26, the HOMO in NaH O at 5.14 eV arises from the Na48tom whereas the deeper orbitals mostly come from the H O molecule, which has an IPof 12.6 eV (see Fig. 1). The electronic density of the HOMO is thus clearly different fromthat of the other orbitals, so that self-interaction for different electronic states becomesvery different.The example of a metal-covalent complex shows that ADSIC is not a safe fire solutionto the self-interaction problem. It is likely to fail whenever the s.p. states cover grosslydifferent regions of space. This also occurs in the fragmentation of a molecule and inprocesses with high ionization. Nonetheless, ADSIC provides very often a remarkablyaccurate and simple approximation to SIC. It is always worthchecking whether a givenproblem allows one to employ ADSIC. In the following examples, we will often recur toADSIC.4.1.2.
Optical response
One of the most prominent observables of electronic dynamics is the optical responsemeasured in terms of the photo-absorption strength. It gives insight into the spectrumof dipole transitions and provides useful information on the collective modes and theparticle-hole excitations of the system. Note that the word ”optical” is generic in thesense that it covers also the spectrum outside the optical range of frequencies. Consider,for instance, covalent systems where the dominant peaks rather lie in the UV domain.Here, we want to explore the impact of SIC on the optical response.The optical response can be calculated in various ways. It is often evaluated by comput-ing the response function directly in linearized TDLDA. Having a fully fledged TDLDAcode at hand, it is technically and conceptually simpler to employ spectral analysis forthat purpose [86, 147, 148]. To that end, we initialize the electronic dynamics by apply-ing an instantaneous dipole boost to the electronic wave functions. We then record thetime-dependent dipole momentum : D ( t ) = (cid:90) d r ( r − R cm , ion ) ρ ( r , t ) , (45)where R cm , ion denotes the center of mass of the ions. The dipole strength S D ( ω ) isobtained by Fourier transforming D ( t ) −→ ˜ D ( ω ) yielding finally S D ( ω ) ∝ (cid:61){ ˜ D ( ω ) } . Onecan alternatively look at the power spectrum | ˜ D ( ω ) | which basically contains the sameinformation as the dipole strength.We already presented in the bottom right part of Fig. 1 the optical response calculatedin ADSIC of H O in the three spatial directions. Fig. 27 now displays the power strengthaveraged over the three spatial directions, and compares the calculations done in LDAwith those in ADSIC. Each photo-absorption spectrum is complemented by the sequenceof one-particle-one-hole (1 ph ) states with dipole character calculated from the static s.p.energies of occupied and empty states. First, we notice that the distribution of 1 ph statesin LDA is very different to that in ADSIC. This reflects the different s.p. spectra of LDAand ADSIC. ADSIC tends to localize the wave functions more than LDA which resultsin a more compact electron cloud (associated to a larger IP) and a wider span of dipoletransitions. It is also interesting to note that most LDA dipole transitions lie in thecontinuum (IP (cid:39) . ph spectra, the dipolespectra look very similar. The reason is that the recoupling of the pure 1 ph states tothe true excited states is dominated by the Coulomb Hartree term [149, 150] which isthe same in LDA and SIC. This defines the overall position of dominant dipole strength49
10 15 20 25 30 P o w e r s p ec t r u m ( a r b . un it s ) Frequency (eV)
LDA
IP5 10 15 20 25 30 P o w e r s p ec t r u m ( a r b . un it s ) Frequency (eV)
LDA IP P o w e r s p ec t r u m ( a r b . un it s ) H O ADSIC IP P o w e r s p ec t r u m ( a r b . un it s ) H O ADSIC IP Fig. 27. Optical response of H O calculated in LDA (bottom) and ADSIC (top). The dashes indicate theionization potential (IP) calculated as the opposite of the HOMO energy. The vertical full lines standfor the possible static dipole transitions. which is more or less robust. The underlying 1 ph structure has an impact on the detailedfragmentation pattern which can depend more sensitively on the level of SIC treatment.The payoff between Coulomb interaction and 1 ph structure depends, of course, on thesystem. Metal clusters have the pronounced Mie plasmon mode which is dominated bythe Coulomb interaction, thus very robust against SIC. On the other hand, systems withfuzzy dipole spectra are more critical. The example H O is somehow in between.4.1.3.
SIC-revisited photoelectron spectra
Since the early days of DFT, the interpretation of Kohn-Sham (KS) orbitals has beena matter of debate. A direct comparison of experimental PES spectra to s.p. spectra oreven better so to our dynamically computed PES (see Sec. 3.3.5) is an a posteriori proofof the meaning to be given to KS s.p. energies, following the basic multiphoton ionization(MPI) relation (9). The point is rather easy to accept at the LDA or ADSIC level to theextent that the KS hamiltonian ˆ h KS (see Eq.(15b)) is well defined and common to all KSorbitals. The situation is more involved in the case of SIC where the KS hamiltonian (18)becomes state-dependent. The ”2setSIC” solution scheme (18–21) allows, nonetheless, todefine unambiguously s.p. energies. It is thus interesting to see how SIC performs forcomputing PES, to see to which extent the SIC s.p. energies have a similar meaning asin LDA or ADSIC. Rather than making a comparison to experiment, which will not testthe internal capabilities of the theory, it is here more interesting to check the evolution50f the PES peaks within varying the laser frequency and see whether they follow theMPI rule (9). This would give an indication on their possible interpretation. The pointis illustrated for the case of the planar metal cluster Na . Two different laser pulses havebeen used with frequencies 8.16 and 10.9 eV. The laser intensities have been adjustedin each case to obtain about the same low total ionization around 0.006, thus well inthe perturbative regime where PES signals are not yet blurred by Coulomb shift (seeSecs. 4.3.4 and 4.3.5.2). The laser polarization is taken normal to the cluster’s plane.Fig. 28 displays both PES. For a given laser frequency ω las , one clearly observes copies P E S ( a r b . un it s ) E kin (eV) ν = 1 ν = 2 ν = 3 I = 2 × W/cm , ω las = 8 . eV I = 2 . × W/cm , ω las = 10 . eV Fig. 28. Photoelectron spectra for Na (ionic configuration in the inset) irradiated by two different laserpulses with pulse duration of 60 fs, intensity I and frequency ω las as indicated, and polarized along thedirection normal to the Na plane. The static Kohn-Sham orbital energies, shifted by nω las where ν isthe number of involved photons, see Eq. (9), are also indicated as vertical dashed lines in both cases.Adapted from [151]. of the same pattern which are separated by ω las . Each pattern exhibits peaks whichare positioned at values of the kinetic energy E kin following the standard MPI relationEq.(9), i.e. ε kin ,j = ε j + ν (cid:126) ω las . The ε j entering this equation are the eigenvalues of thestationary equation h SIC | ϕ j (cid:105) = ε j | ϕ j (cid:105) , see Eq. (21), while ν corresponds to the numberof photons involved. The remarkable fact that the peaks of the PES, obtained from thecalculation of a time propagation of the 2 sets { ϕ j } and { ψ α } , coincide with the static ε j validates the interpretation (and the definition) of these energies as sound s.p. energies.This also supports the identification of the ϕ j ’s as the physical wave functions of theassociated s.p. states whose characteristics are measurable via the PES.4.2. Using ionization as an observable
In this section, we discuss basic mechanisms in the laser irradiation of an electronicsystem leading to significant electronic emission and their analysis in terms of totalionization as an observable. For moderate laser intensities, a major issue is the relationof the laser frequency with the optical response peaks, especially collective ones. Wewill present here two generic scenarios for this resonance effect. We will first explorelaser irradiation of a water molecule demonstrating off- and on-resonant ionization. In51ec. 4.2.2, we will take advantage of resonant enhancement of ionization to explore theionic dynamics via the use of a pump-and-probe (P&P) setup. In Sec. 4.2.3, we will finallyconsider again a P&P setup, but this time within using a train of attosecond pulses.4.2.1.
Off- and on-resonant ionization
As already discussed in Sec. 1.2, the great versatility of lasers through the choice offrequency, intensity, pulse duration and shape, offers experimentalists and theoreticians aworld of dynamical scenarios. To gather orientation in this huge landscape of options, wefirst explore the impact of laser frequency. To that end, we consider the dynamics of laserexcitation of a H O molecule, with techniques similar to those used in [152]. The resultsare shown in Fig. 29. The duration of the laser pulses is in all cases 20 fs. Four frequencies . . . .
002 0 10 20 30 40 50 60 N e s c Time (fs) − . − . . . D y ( a ) H O
13 15 17 19 ω = 12 . ω = 13 . ω = 15 . ω = 16 . | f D y | Fig. 29. Resonant and off-resonant irradiation of H O by laser pulses polarized along the symmetry axisof H O, denoted by y , and of various intensities ( I = 10 , I = I = 10 , I = 10 W/cm ) andvarious frequencies ( ω = 12 . ω = 13 . ω = 15 . ω = 16 . y direction. Bottom : time evolution of the total ionization N esc . Inset : optical strength of H O in y direction with horizontal energy axis in eV. The vertical full lines indicate the chosen laser frequencies,and the dotted one corresponds to the IP, here at 15.1 eV. have been explored, namely ω = 12 . ω = 13 . ω = 15 . ω = 16 . I = 10 , I = I = 10 and I = 10 W/cm to keep the52aximal dipole amplitude similar. The laser polarization is along the symmetry axis ofthe water molecule, denoted here by y .Let us start with ω and ω . Both frequencies are below the IP of H O (15.1 eV). Theoptical response of H O in y direction is shown in the inset in the bottom panel : ω lieson a double peak and, as we will see, corresponds to a resonant frequency, whereas ω does not match any dipole transition and is thus off-resonant. The top panel of Fig. 29shows the time evolution of the electronic dipoles. As expected [14], the red curve forthe off-resonant ω las = ω nicely follows the laser pulse profile and dies out with thelaser signal at 20 fs. The total ionization N esc shown in the bottom panel (red line)stays very close to zero for this low frequency case. The resonant ω las = ω proceedsdifferently : During the first 15 fs, the dipole signal (light green curve in top panel) stillfollows the laser profile, but then continues to oscillate with large amplitude long afterthe laser pulse is switched off. Such resonant oscillations come along with a larger depositof energy in the molecule and thus stronger ionization. This is clearly demonstrated inthe bottom panel of Fig. 29, where N esc (green line) steadily increases with visible stepsperfectly correlated to the maxima in the dipole oscillations. The total ionization N esc in the resonant case is orders of magnitude larger than in the off-resonant one, althoughthe laser intensity in the resonant case is 10 times smaller. This can be understood interms of the Keldysh parameter (8) which is, in both cases, much larger than 1. Thisindicates that we are in the frequency-dominated regime where such differences matter,see Sec. 1.2.2.Note also the initial profile of the ionization is determined by the laser pulse suchthat the maximum slope coincides with the laser peak amplitude. But this profile isdelayed by the time it takes for the escaping electrons to reach the box bounds. The finalnon-vanishing slope is related to the non-vanishing dipole oscillations.The next two cases consider an off-resonant frequency ( ω ) and a resonant one ( ω ),now both above the IP. For the sake of clarity, the corresponding dipole signals are notdisplayed. They show again the same typical pattern of resonant (long standing after-oscillations) and off-resonant (signal dies out with the laser) response. We only showthe associated total ionization N esc in the bottom panel of Fig. 29 (black thick line for ω and blue thin line for ω ). The off-resonant case increases with a slope following theamplitude of the dipole oscillations and levels off to a plateau after the laser pulse isover. The case ω yields much higher ionization than the case ω , even if its intensity I is an order of magnitude smaller than I . This happens because ω stays above the IPand can ionize directly with one-photon processes. Finally we compare the two resonantcases ω < IP with ω > IP. Although I /I = 100, both N esc are very similar. Two effectscooperate here : i) ω > IP and ii) the strength of the mode excited at ω is at least 3times stronger than that at ω . This demonstrates, once again, the importance of thelaser frequency in relation to the optical spectrum in order to drive large ionization.The point is again illustrated, this time in a more systematic manner, in Fig. 30 whichdisplays the dependence on laser frequency of the total ionization N esc of irradiatedC (left) and Na (right). The figure is a continuation to the pedagogical Fig. 4discussed in Sec. 2.1, this time for two more complex systems though. For C , with thechosen laser parameters ( I = 7 . × W/cm and ω las =14–28 eV), we are again inthe frequency-dominated regime. The ionization N esc ( ω las ) (light green curve) exhibitsstrong oscillations with ω las which match remarkably well the optical response of C (black dashed line). This shows that the signal of photoemission N esc ( ω las ) is close to the53 N e s c ω las (eV) C T pulse = 60fs I = 7 . × W/cm opt. resp. − − − − N e s c ω las (eV) Na T pulse = 300fs10 W/cm W/cm W/cm W/cm Fig. 30. Total ionization N esc as a function of laser frequency ω las . Left : case of C , calculated in full3D with a radius of 6.44 a , irradiated by a laser pulse with T pulse = 60 fs and I = 7 . × W/cm .Right : case of Na irradiated by laser pulses with T pulse = 300 fs and intensities I as indicated (usingCAPS). The black dashed line represents the optical response in both cases. signal of photo-absorption, at least above the emission threshold.The situation is similar in Na at the lowest laser intensity (light green curve in rightpanel) where N esc exhibits strong oscillations with ω las , once again fitting fairly well thoseof the optical spectrum. However, if the laser intensity is increased, we progressively leavethe frequency-dominated regime to enter the field- (intensity-)dominated domain (seeSec. 1.2.2). And indeed, the fragmented structure of N esc steadily broadens to be finallywashed out at the highest intensity (see top black curve). At the same time, the valuesof N esc for a given ω las also increase, since there are more and more photons pulling onthe valence electrons of the cluster.4.2.2. Pump and probe (P&P) analysis of ionic dynamics
The emergence of fs lasers allowed the development of time-resolved studies of molecu-lar reactions through pump-and-probe (P&P) experiments. The typical strategy of sucha fs spectroscopy is simple. An initial short laser pulse (pump) excites the electronic sys-tem which leads to subsequent ionic motion. This motion in turn changes the electronicresponse according to the actual ionic configuration. This change is explored by the re-sponse (e.g., ionization) of the system to a second laser pulse, the probe, sent after acertain time delay. Scanning the reaction strength as a function of delay time allows oneto map the time evolution of the molecular system. There are, of course, many variants ofthis generic strategy according to the variety of molecules and flexibility of laser pulses.Altogether, fs spectroscopy has become an extremely powerful analyzing tool in physicsand chemistry, for early reviews, see [153, 154].Of course, P&P analysis is also an extremely interesting tool in cluster physics. Verysmall clusters allow scenarios very similar to those in simple molecules, see e.g. experi-ments on trimers [155, 156] and associated theory [157]. Larger clusters are too complexfor the very subtle and detailed pathways followed in small molecules. One better looksfor global properties of the ionic background as, e.g., radius or deformation, and one needsprominent signal in the dense electronic spectrum. Metal clusters are distinguished by54he dominant Mie surface plasmon resonance [9], whose peak frequency is predominantlydetermined by cluster radius and deformation. Thus there are many P&P studies onmetal clusters, either free, deposited on a surface, or embedded in a substrate, see e.g.Sec. 5.3.4 of [12].As P&P experiments on clusters are rather demanding, early studies achieved a com-parable, although coarser, effect by varying the temporal width of a single laser pulse.For an early example on Pt clusters, see [33]. We exemplify this type of analysis here forthe case of Ag clusters embedded in a He droplet for better handling [158]. The clusters
Fig. 31. Yields for selected Ag q + ions after irradiation of a Ag cluster with a laser pulse of wavelength800 nm, drawn as a function of the width τ of the the laser pulse. For the shortest pulse of 130 fs, thepeak intensity is I = 1 . × W/cm . For other τ , the fluence ∝ τ I has been kept constant. Thedata are normalized for better comparison, and the fit curves serve as a guide to the eye. From [158]. are irradiated with laser pulses of fluence τ I = 156 W fs/cm (where τ is the pulsewidth and I the peak intensity). This strong pulse leads to a disintegration of clustersproducing all sorts of fragments and highly ionized Ag atoms. The charge state q ofthe emerging Ag q + ions is an indicator for the violence of the reaction and thus for thestrength of the laser-cluster coupling. Fig. 31 shows the ion yield as a function of pulsewidth. All ionization stages q show a strong dependence on τ with a distinct maximumfor a certain τ . This optimum pulse width, which was already observed in [33], resultsfrom an interplay of (ongoing) laser pulse, ionic expansion, and plasmon frequency. TheIR laser pulse first triggers ionization. The Coulomb pressure thus generated leads to aslow expansion of the cluster. And the plasmon frequency (originally in the visible range)decreases with increasing radius, until the laser comes into resonance with the plasmonwith subsequently strong energy absorption and violent reaction. If the laser pulse is tooshort, it is over before resonant conditions are reached. If it is too long, it becomes tooweak (remind the constant fluence thus implying decreasing intensity with increasingpulse duration) to trigger sufficient expansion. Such a maximum is seen in Fig. 31 foreach charge state, however at different delay times τ . The interpretation given in [158],furthermore, addresses a subtle point in laser experiments. The laser intensity is not con-stant over the spatial width of the beam. It decreases when going away from the focus.55hus clusters outside the focus receive a weaker signal than those right in the focus ofthe beam. It is assumed that these lower charge states are related to lower intensities outof focus. Therefore, the experiment so to say produces at once results for different laserfluences.The rather involved P&P analysis is easier for clusters in/on a substrate because thisallows a higher density of reactive centers. Therefore, most P&P studies on clustersare performed in/on substrate. The typical setup is that of a chromophore in an inertsubstrate. The latter thus serves mainly as a support for the cluster. The principles andthe richness of P&P analysis remain unaffected by the inert substrate. There is a coupleof measurements of an electronic property, the electronic relaxation time, for clusterson surfaces in a variety of material combinations [159, 160, 161]. More typical for P&Panalysis is the study of ionic oscillations which has been performed, e.g., for Ag clustersembedded in glass matrix [162, 163]. A much more gentle support is provided by liquid Heclusters, which were already used as useful laboratory for studying molecular propertiesunder well controlled conditions [164]. The He environment couples such softly to anyother material that one can consider the embedded system as being practically free. Thereare then several instructive P&P experiments of Ag clusters in He droplets, e.g. [38, 165](called “dual pulse” experiments in these publications). A detailed description of the largescale dynamics of Ag clusters is very expensive. Theoretical investigations are thus oftenperformed for Na clusters as practicable model systems for metal clusters [166, 167, 165].The dominance of the Mie plasmon peak in metal clusters allows a particular P&Pstrategy which does not rely on directly hitting the resonance but uses just the distanceof the Mie plasmon frequency to the laser frequency to map the underlying ionic dynam-ics. This strategy has been studied in detail for the global breathing (radius oscillations)of Na clusters in [166] and for the dynamics of cluster deformation in [167]. We illustratethe scheme here for the case of breathing. Fig. 32 shows the result of a theoretical explo-ration for the cluster Na using TDLDA for electronic dynamics coupled to moleculardynamics for the ionic motion [168, 14, 12]. The cluster Na is nearly spherical. Thepump pulse ionizes it quickly to charge state Na . This produces a Coulomb pressurewhich triggers slow breathing oscillations of the whole cluster, while deformation is neg-ligible along the whole dynamics. The radius oscillations after the mere pump pulse areshown in panel (a) of the figure. The Mie plasmon resonance depends on charge state andcluster radius. An estimate is shown in panel (b). One sees the fast initial blue-shift dueto the fast initial ionization to q = 4 + . After that, one finds oscillations which perfectlyfollow the radius oscillations according to ω ( t ) ∝ R ( t ) − / [150]. The laser frequency forthe probe pulse is also indicated. It was chosen safely below the Mie resonance such thatthe actual Mie frequencies never cross. The electronic response to the probe pulse is smallif the Mie frequency is far from the laser, and large if it comes close. This can be seenin panel (c) from the maximum amplitude of the dipole signal during the probe pulse.The strong dipole response leads to further ionization shown as additional number ofescaped electrons ∆ N esc in panel (d). It is, of course strongly correlated with the dipoleamplitude. Tracing the chain of correlations back to panel (d), we can conclude that theextra ionization directly maps the global ionic radius using the scheme with the remotelaser frequency as an “observer”. Net ionization thus provides a direct (time resolved)analysis of ionic motion, in that case dominated, on the rather short times considered,by radial oscillations. The actual long term evolution of the system is, in fact, Coulombexplosion. The interesting aspect is that the long path to explosion is accompanied by56 ig. 32. Pump and probe spectroscopy of ionic breathing vibrations of Na using the Mie plasmonresonance as indicator. Panel (a): time evolution of the ionic r.m.s. radius, (cid:112)(cid:80) I R I , after the pumppulse. Panel (b): time evolution of the Mie plasmon frequency after the pump pulse. The laser frequency ω las = 2 . N esc induced by the probe pulses as a function of time delay. Pump and probe pulses have the same properties:photon frequency ω las = 2 . I = 1 . × W / cm , and a sin shape with FWHM= 24 fs.The pump laser produces very quickly an initial emission of N esc = 3 electrons, thus delivering a totalcharge state q = 4 + . After [166]. monopole (radial) oscillations which are directly visible in the ionization signal (actuallyvia the plasmon peak).Almost all P&P studies on clusters have used the ionization yield as an observable.One expects that one could learn more from more detailed observables, particularlyfrom time-resolved PES and PAD. Such experiments are, of course, much more complexand still more demanding than the, already intricate, traditional P&P measurements.Nonetheless, first studies in that direction have been published [169, 170] so far in theregime of hefty excitations. The field is widely open for further studies in more moderate57xcitation regimes.All the P&P studies, including the above example, show that the ionization yieldcan depend sensitively on the choice of the laser pulse characteristics. As the temporalprofile and all other laser parameters can be tuned in an extremely flexible manner, thequestion naturally arises whether one could tune laser pulses for maximum yield (or otherdesired reaction properties). This is the idea of ”optimal control” which is of particularimportance in chemistry and molecular physics, see e.g. [171, 172]. Again, the applicationto large clusters is more involved and allows more strategies to be tracked. An interestingstudy using optimal control, e.g., to trigger the yield of highly charged Ag atoms fromAg clusters can, nevertheless, be found in [173].While addressing promising future developments of P&P studies, we ought to mentionthe upcoming availability of attosecond pulses. These allow P&P studies which resolvefeatures of electronic dynamics. We will discuss that in the next section 4.2.3.4.2.3. Towards P&P experiments with attosecond pulses
Electron dynamics can also be analyzed at its own pace if one is able to handle pulsesmuch shorter than typical electronic time scales in the fs regime. This is nowadays ex-perimentally accessible down to some hundreds of attoseconds, at least in the form of atrain of attosecond pulses. This yields access to details of electronic dynamics subject toelectromagnetic perturbations. Early convincing tests were performed in simple atomssuch as He and Ar [174] on the basis of a P&P setup involving a UV atto-train on top ofan IR field. It was shown that the total ionization may exhibit marked oscillations as afunction of the delay between UV train and IR signals, once the repetition rate of sevenattopulses per train is chosen to be half the IR period. The analysis of these experimentswas supported by simple simulations using the Time-Dependent Schr¨odinger Equation(TDSE) with a single active electron. The TDSE also served as a basis for further theo-retical investigations, either directly [175, 176, 177] or in perturbation theory [178]. Theseapproaches gave convincing clues on the origin of the modulation of the ionization butso far, no robust many-electron theory is available to explain the observations.More recently, experiments were generalized to simple molecules with qualitativelysimilar results as in the atomic case, although with a slightly different combination of IRand UV pulses [179, 180, 181]. The first fully microscopic calculations were performed onthis occasion and led to results remarkably compatible with experiments [182]. The actualinterpretation of the underlying mechanism is nevertheless, again, to be understood inmore detail. Still, the remarkable agreement between theory and experiments is worthbeing presented and discussed in one test case.We consider here the N molecule as a test case [182]. The laser pulse consists in anIR component : V IR ( t ) = E IR f ∆ T ( t ) sin( ω IR t ) , (46)with a frequency ω IR = 1 .
58 eV and with a sin pulse profile f ∆ T of FWHM= ∆ T / (cid:39)
29 fs (see Eq.(2)). At the same time is superimposed a train of n attopulses, each of themlabeled by i , which reads : 58 atto ( t ) = n (cid:88) i =1 E atto (˜ t ) f δt (˜ t − t i ) sin (cid:2) ω UV (˜ t − t i ) (cid:3) e − t / ∆ T , (47a) t i = ∆ t + ( i − δt + δt (cid:48) ] , (47b)˜ t = t − ∆ t , (47c)where δt (cid:39) .
29 fs is the actual duration of each individual attopulse, and δt (cid:48) the timeseparation between two successive attopulses. The attopulse train (APT) is delayed by avariable delay ∆ t with respect to the IR pulse (starting at t=0). Both the IR pulse andthe APT are linearly polarized. The individual attopulse shape is again a sin profile (seeEq. (2)) of FWHM= δt/
2. A key point of the setup is to fix the time interval betweentwo successive attosecond signals. We choose here δt + δt (cid:48) = T IR / (cid:39) . f s , which isexactly half the IR period. The amplitude of the APT is further modulated by a Gaussianenvelop of width such that the total APT duration is about half (29 fs) the total durationof the IR pulse. This fixes the number n of attopulses which is, in this case, 22. The peakintensity of the IR pulse is chosen to be I IR = 10 W/cm , while that of the attopulsesis, at maximum of the Gaussian envelope, I atto , = 10 W/cm . Finally, the frequencyof the APT lies in the UV domain, ω UV = 20 . is around 16.3 eV < ω UV . The remarkableresult of the experiments is that combining the IR and the attopulses leads to a significantenhancement of ionization (while it only adds 1.58 eV on top of the UV photons alreadyabove the continuum). Moreover, this enhanced ionization is strongly modulated by thedelay ∆ t . Ionization actually exhibits marked oscillations as a function of delay with aperiod equal to half the IR period. Maxima of oscillations are attained for delays suchthat the attopulses are in phase with maxima or minima of the IR pulse, which explainsthe doubled frequency of ionization maxima as compared to the IR frequency.The case is illustrated in Fig. 33 where we have plotted the total ionization, the averagekinetic energy of emitted electrons and the anisotropy β characterizing the PAD (seeEq. (37) and Sec. 4.4). As one UV photon suffices for ionization, β provides a completecharacterization of the PAD. The average kinetic energy is defined as (cid:104) E kin (cid:105) = (cid:126) m (cid:90) d r j ( r ) ρ ( r ) (cid:2) − M ( r ) (cid:3) , (48)where m is the electron mass, j ( r ) is the local current, ρ ( r ) the local density and M ( r )the mask function used to evaluate emitted electrons (see Eq. (23c)). It provides a simplemeasure of the PES in terms of one number. All three signals in Fig. 33 display remarkableoscillations as a function of delay time ∆ t with a period equal to half the IR period.Both (cid:104) E kin (cid:105) and β oscillate in phase and in opposite phase with the total ionizationrespectively. Indeed, if the deposited energy content is about the same, the higher theionization, the lower their average kinetic energy. And not surprisingly, the more energeticthe emitted electrons, the more aligned the emission along the laser polarization axis andthe larger the β . These oscillations of the total ionization perfectly match those observedexperimentally. A comparison with PES and PAD has to wait until these quantities areexperimentally available. 59 β ω IR = 1 .
58 eV, I IR = 10 W/cm ω UV = 20 . I UV = 10 W/cm N h E k i n i ( e V ) T IR / N e s c ( × ) Delay (fs)
Fig. 33. N ionization properties irradiated by an IR pump and an attopulse train (see text for details)as a function of delay time between attopulse train and IR pulse. Bottom : total ionization. Middle :average kinetic energy per emitted electron, see Eq. (48). Top : anisotropy parameter β , see Eq. (37).The faint dashed lines indicate the sequence of maxima and minima regularly separated by half the IRperiod. Adapted from [182]. We finally end this discussion by mentioning recent experimental P&P dynamics whichuse XUV pulses for both the pump and the probe. This is at variance with the examplediscussed above, where the pump is an IR pulse and the XUV probe is constructedfrom some of its high order harmonics using the so-called RABITT (Reconstruction ofAttosecond Beating by Interference of Two-photon Transitions) technique [183, 184].To distinguish the two kinds of P&P, one sometimes quotes them as IR-pump-XUV-probe and XUV-pump-XUV-probe experiments respectively. One type of setup uses acoherent splitting of a XUV light produced by a FEL. For instance, this technique hasbeen successfully applied to small molecules as N and O [185], or C H [186]. In theselatter examples, the XUV frequency is 38 eV, the photon intensity between 10 and10 W/cm , and the XUV pulse duration of 30 fs. The delay time resolution is 1 fs,and the whole delay time can vary over the ±
350 fs range. The advantage of using XUVlight for the pump and the probe is to ionize the species under study by absorption ofa few photons, at variance with an IR pulse. Therefore, this P&P setup can follow theinduced Coulomb explosion at a time scale of a few fs. Very recently, some experimentswent even further by taking advantage of high harmonic generation from an IR pulse60rradiating an atomic gas jet : the produced XUV attopulses are then separated from theIR pulse, filtered to keep only one pulse which is at the end split into two coherent XUV aspulses. This brand new technology has been applied in the irradiation of Xe atoms [187]and the H dimer [188]. The XUV intensities are about 10 − W/cm , their durationabout 600 as, their frequency is 14 eV, and the delay time is well below 1 fs. Such anexperimental apparatus thus enables to track the Coulomb explosion dynamics over thewhole reaction path and on a time scale never attained before. In both XUV-XUV P&Pexperiments, fragment or ion yields are measured as a function of delay time. To the bestof our knowledge, no electronic observable has been measured so far. There are also veryfew real-time calculations of such a dynamics [188, 189]. They employ a time-dependentSchr¨odinger equation of the full electronic and nuclear wave function (note however thatonly vibronic modes – no rotational mode – are considered in these calculations). Forlight atoms as those in H , a quantal description of the nuclei is probably compulsory.The necessity of such a fully quantal treatment is more questionable for heavier atoms,as in N . Anyway, the computational cost of such calculations becoming too prohibitivefor larger covalent systems, this probably calls for a classical treatment of the ions, evenif H + nuclei come into play.4.3. Dynamical aspects in photoelectron spectra
We have discussed above how (static) s.p. spectra can be extracted from the peaksobserved in PES using Eq. (9), see for instance Fig. 28. This identification can be workedout by perturbation theory and requires a laser with moderate intensity and high fre-quency resolution. True dynamical processes exploit more of the versatility of lasers. Theaim of this section is thus to discuss the impact of the laser parameters frequency, inten-sity, and pulse duration on PES, and find out how dynamical aspects can be analyzedthrough PES.4.3.1.
Impact of pulse duration
Any laser pulse of finite duration delivers a distribution of frequencies about its meanfrequency ω las . The longer the pulse, the sharper this distribution. One can evaluate thewidth of this distribution by calculating σ las = (cid:82) ∞ ( ω − ω las ) | ˜ I ( ω ) | d ω where ˜ I ( ω ) isthe Fourier transform of the time-dependent laser intensity I ( t ). As an illustration, wegive in the following table some widths related to pulse duration at ω las = 20 eV [190].One consequence of the finite frequency width of the laser pulse is that a PES (in the T pulse (fs) 10 30 60 75 200 1000 σ las (eV) 0.44 0.16 0.083 0.068 0.027 0.0058Table 1Width of frequency distribution for different laser pulse durations, around a mean value of ω las = 20 eV. perturbative regime) does not display a sharp spike exactly at E kin = ε i + nω las , but amore or less soft peak around this E kin . Finite pulse duration thus produces a broadeningof the PES peaks, the larger the shorter the pulse. To exemplify this effect, Fig. 34displays PES of Na irradiated by laser pulses of intensity 10 W/cm , frequency ω las =6 . PE S ( a r b . un i t s ) E kin (eV) Na ω las = 6 . I = 10 W/cm
50 fs100 fs200 fs400 fs
Fig. 34. Photoelectron spectra of Na after irradiation by laser pulses of indicated characteristics. SomePES have been up-shifted for the sake of clarity. as the intensity of the laser pulse remains sufficiently low. Higher intensities can blur thispicture because ongoing ionization induces a drift of the peaks due to a Coulomb shiftof the levels. This will be addressed in Secs. 4.3.4 and 4.3.5.2.4.3.2. Impact of laser polarization
In section 2.3.2, we have seen that a combined PES/PAD, that is an energy- and angle-resolved photoelectron spectrum, can deliver a lot of information on the dynamics of thephotoemission. In the perturbative regime, it reveals the angular distribution for emissionfrom specific s.p. states. A simplified view can be obtained by restricting the analysis totwo specific directions : one parallel to the laser polarization axis ( θ = 0 , ◦ ) and oneperpendicular to it ( θ = 90 ◦ , ◦ ). The point is illustrated in Fig. 35 where calculatedparallel and perpendicular PES of Na − are plotted. The sodium cluster has a fixedorientation here (orientation averaging will be discussed in Sec. 4.4). The calculated s.p.energies are ε = − .
82 eV, ε = − .
72 eV, and a pair of almost degenerate ε , = − .
43 eV. As in the case of Na (see Fig. 28), the peaks of the PES perfectly fulfill therelation E kin = ε i + νω las with ν = 1 ,
2. The first group of peaks between 1 and 3 eVcorresponds to ν = 1. States 1 and 2 predominantly emit along the laser polarization axis,while states 3 and 4 show a clear preference of emission in perpendicular direction. The2-photon process (between 5 and 7 eV) suppresses even more strongly the perpendiculardirection, and thus the parallel photoemission dominates. This indicates a general featureof multiphoton emission : the higher the photon number ν , the larger the anisotropy β corresponding to an increasing dominance of emission parallel to the laser polarizationaxis. 62 − − − − PE S ( a r b . un i t s ) E kin (eV) Na − ω las = 4 . eV I = 10 W/cm T pulse = 60 fs k⊥ Fig. 35. Photoelectron spectra of Na − for a given orientation, irradiated by a laser with parametersas indicated. Black full curve: PES along the laser polarization axis. Red dashes: PES in the directionperpendicular to it. Calculations were done in cylindrical approximation for the electronic potentials. Impact of laser frequency
The estimate (9), i.e. E kin = ε i + νω las , of PES peaks establishes correctly the relationbetween peaks at E kin and corresponding s.p. energies ε i . However, it does not tell any-thing about the strength with which the peaks appear. And here, we can have dramaticdifferences between one-photon processes and multiphoton ones. As an illustration, weshow in Fig. 36 ionization pattern of C in the one-photon (left panels) and in themulti-photon (right panels) regimes [191]. Orientation averaging presented in Sec. 3.4has been applied in the theoretical calculations, to allow a comparison with experimen-tal measurements. The upper panels compare theoretical and experimental PES, andthe lower panels show as complementing information the depletion (blue full lines) andthe occupancy (green dotted lines) of the corresponding s.p. levels. The one-photon casewas already presented in Fig. 24 which compares theoretical and experimental combinedPES/PAD and VMI. The laser pulse in this case experimentally stems from synchrotronradiation at 20 eV with a pulse duration of about 1 ps. Theoretical calculations were donefor the same frequency, intensity I = 7 . × W/cm , but shorter pulse T pulse = 60 fsfor practical reasons. In the one-photon regime ( ω las (cid:29) IP), we do not expect that thepulse duration is essential. Calculations yield a total ionization of about 0.006. In themulti-photon case, the theoretical parameters are chosen as in the experiment, that is T pulse = 60 fs, I = 7 . × W/cm , and ω las = 1 .
55 eV. Here, the calculated total ion-ization is about 0.07. Comparing both cases (one- and multi-photon), PES and depletionpattern are completely different. The one-photon case also shows a marked differencebetween experiment and theory. This has to be discussed in detail.The bottom panels of Fig. 36 display the s.p. levels (with occupation weight) and theirdepletion. For ω las = 20 eV (lower left panel), all states can emit by a one-photon process.And indeed, we observe that most states contribute to the total ionization. As expected,the theoretical PES resembles the s.p. depletion pattern very much [192].We now turn to the multi-photon case (right panels of Fig. 36). With the laser frequencychosen here, at least 6 photons are needed to bring electrons from the HOMO intothe continuum. Since the probability of ejection decreases with the number of photonsrequired, only the least bound states can be significantly depleted. And this is what is63
123 0 3 6 9 12 PE S ( a r b . un i t) Kinetic energy (eV) C ω las = 20 eV I theo = 7 . × W/cm s . p . d e p l e t i o n ( × ) s.p. energies (eV) Exp. ( ’ × occupation 10 − − − PE S ( a r b . un i t s ) Kinetic energy (eV) ω las = 1 . eV I = 1 . × W/cm T pulse = 60 fs00.10.20.30.40.5 -20 -17 -14 -11 -8 00.10.20.30.4 s . p . d e p l e t i o n s.p. energy (eV) ω las = 1 . eV I = 1 . × W/cm T pulse = 60 fs ω las ω las Fig. 36. Top : Theoretical and experimental photoelectron spectra of C with radius of 6.763 a andorientation averaging. Bottom : Calculated single particle depletion (blue full lines) compared with staticoccupation numbers (green dots). Left : one-photon regime ( ω las = 20 eV), laser pulse length of 60 fsand intensity of 7 . × W/cm (theory) or a synchrotron irradiation of duration of about 1 ps. Right :multi-photon regime with ω las = 1 .
55 eV, laser pulse length of 60 fs and intensity of 1 . × W/cm ,both in theory and experiment. In the bottom left panel, the single particle depletions are multiplied by100 for a better comparison. Adapted from [191]. observed in the right bottom panel : the states which emit the most are the HOMO,HOMO − −
3, precisely separated by about ω las = 1 .
55 eV. Therefore, wecannot expect that the PES maps the whole s.p. energy spectrum. Indeed, the theoreticalPES exhibits oscillations consisting in a broad peak repeated several times, separated by ω las from one copy to the other. These oscillations constitute the typical multiphotonionization (MPI) pattern (see Sec. 2.1). Each one of these MPI peaks is rather wellbundled due to the fact that the few emitting states line up rather well with the photonfrequency.In both cases (one- and multi-photon), the experimental PES differ from the theo-retical ones. The difference looks particularly large for the one-photon case (upper leftpanel). Here, the PES are still fairly comparable at higher energies (7–13 eV). Whiletheoretical calculations were done at the ground state configuration (zero temperature),the experimental peaks are broadened due to ionic vibrations of C which are ratherlarge at the experimental temperature of 900 K. A huge discrepancy between theoryand experiment is observed at low electron kinetic energy. The experimental PES growsalmost an order of magnitude above the theoretical one. We think that these low-energy64lectrons stem from electron-electron collisions which hinder part of the electrons frombeing directly emitted, but rather lead to auto-ionization mechanisms or thermal electronemission. Such dynamical electron-electron correlations are not included in TDLDA. Sowe are missing here most of the low-energy electrons. This could be cured with theoriesbeyond TDLDA which will be discussed in Sec. 4.5.2. For the time being, the comparisonbetween experiment and theory is relevant only for the highest photoelectron energiesdominated by direct electron emission. And there, the agreement is very satisfying.A difference between experimental and theoretical PES is also seen in the multi-photonregime (upper right panel). We first note a shift between the position of the theoreticalpeaks and that of the experimental ones. This might be due to a slight uncertainty in thedetermination of the experimental intensity and/or pulse duration, which can then pro-duce a different ionization stage. And we will see that this can cause a sizeable (Coulomb)redshift, see Secs. 4.3.4 and 4.3.5.2 below. Moreover, the amplitude of oscillations of PESdecreases with increasing kinetic energy in the measurements, while it remains moreor less constant in the theoretical calculations. Once again, the dynamical correlations,which are missing in the theory, are most probably the mechanism responsible for thedamping of the oscillations in the experimental data.4.3.4. Impact of laser intensity
A couple of experimental PES had already been shown in Fig. 10 for the case of C irradiated by laser pulses of various intensities. The results seemed to be all in the sameregime marked by smooth, exponentially decreasing PES throughout. Another examplewas just given above in Fig. 36 where the top right panel shows typical MPI patternof repeated peaks which have at least an exponentially decreasing envelope. For a moresystematic survey, we now discuss computed PES for two different clusters. As a firstexample, we discuss a set of PES of Na (see left panel of Fig. 37) obtained with thesame laser frequency (3.8 eV) and pulse duration (300 fs) but with varying intensities(from bottom to top) starting from I = 10 W/cm and up to 300 I . The chosen − − − − − PE S ( a r b . un i t s ) E kin (eV) Na , T pulse = 300 fs, ω las = 3 . eV I = 10 W/cm I I I − − − PE S ( a r b . un i t s ) E kin (eV)C , T pulse = 75 fs, ω las = 1 . I = 1 . × W/cm . I . I . I Fig. 37. Calculated photoelectron spectra for various intensities. Left : case of Na irradiated by laserpulses of duration of 300 fs, frequency of 3.8 eV, and I = 10 W/cm (calculations done with CAPS).Right : case of C (calculated in full 3D with radius of 6.44 a and with orientation averaging) irradiatedby laser pulses of duration of 75 fs, frequency of 1.5 eV, and I = 1 . × W/cm . The inset zooms inthe 10–14 eV range with the 13-photon ionization of the HOMO and the 14-photon one of the HOMO − frequency of 3.8 eV is basically off-resonant, as the dominant plasmon response of Na is 5.3 eV, thefirst peaks at low kinetic energies in Fig. 37 stem from two-photon processes. The pulseduration is 300 fs to allow a high spectral resolution of the PES peaks, as is observedin the lowest PES (black curve). The laser intensity is also small enough to stay in theperturbative regime, since the total ionization N esc is 0.004 in this case. When the laserintensity is increased by a factor 30 (light green curve), one can spot a slight broadening ofthe peaks, although the PES still shows clear signatures of the underlying s.p. spectrum.The broadening develops to the side of lower kinetic energies, in total yielding a weak red-shift of the peaks. It is to be noted that the total ionization amounts to N esc = 0 . ∝ I and so the absorption of morephotons becomes more probable.At the next stage, I = 100 I , the PES is already smeared to broad steps. Here, N esc =0 .
54 and the Coulomb shift significantly blurs the PES. But still, we can distinguish blocksof one-, two-, and three-photon processes. Finally, the highest intensity of I = 300 I produces an ionization of N esc = 4 . and the laser has a pulse duration of 75 fsand a frequency of 1.5 eV. Five laser intensities I have been considered : I = 1 . × W/cm , 1 . I , 1 . I , 1 . I , and 2 . I . The total ionization of course increaseswith I : N esc = 0 .
03, 0.06, 0.12, 0.31, and 0.60. With the chosen frequency, we need atleast 6 photons to extract electrons (the IP is here 8 eV). At the lowest intensities, weclearly observe the typical MPI patterns of repeated copies of the peak. The peaks aregradually red-shifted when we increase the laser intensity which is, again, the Coulombshift. Moreover, the red-shift increases with laser intensity due to increasing ionizationand finally washes out all structures at the highest I . The disappearance of the MPIpeaks is illustrated in the inset zooming into the 13-photon ionization of the HOMO andthe 14-photon ionization of the HOMO −
1. The net conclusions from that case are thesame as those from Na . But here, it becomes even more obvious that the envelopeof MPI follows in any regime an exponential decrease. Taking this together with thesmoothing due to high ionization delivers then the purely exponential profile resemblingthermal emission.4.3.5.
More on the role of plasmon
Competition between laser and plasmon frequencies
Thus far, we have dis-cussed PES emerging from the interplay between s.p. energies and the pulse frequency.This simple view has to be modified in the vicinity of strong excitation modes of the66ystem. In particular, the dominant Mie surface plasmon in metal clusters can also havea large impact on the PES. To illustrate this point, we discuss the irradiation of Na by laser pulses with duration of 48 fs, intensity of 10 W/cm , and six different laserfrequencies ω las in the vicinity of the Mie plasmon frequency of Na , ω pl = 2 . Fig. 38. Photoelectron spectra in the energy range of the 4-photon process (red vertical dots at ε s +4 ω las )of the 1 s state of Na , after irradiation by laser pulses with duration of 48 fs, intensity of 10 W/cm ,and six different frequencies ω las as indicated. The blue solid vertical line indicates the the 4-plasmonprocess located at ε s + 4 ω pl , with ω pl = 2 . at ε s + 4 ω las , is indicated for each laser frequency by vertical dots. The position of thispeak moves to the blue with increasing ω las . Additionally, one observes a peak whoseposition does not depend on the laser frequency, indicated by the solid vertical line. Itsposition matches the energy of a 4-plasmon process, i.e. ε s + 4 ω pl . When ω las is suffi-ciently separated from ω pl (see the two lowest and the two uppermost curves), one caneasily disentangle the four-photon process from the four-plasmon one. The four-plasmonpeak actually dominates the PES in most of the cases, a further indication of the alreadydiscussed resonant ionization mechanism (see Fig. 29 and 30 in Sec. 4.2.1). In most cases,one can even conceive a coexistence of plasmon and photon excitations. For instance, theuppermost curve shows three prominent peaks, that is the four-plasmon peak, the four-photon peak, and in between a two-plasmon–two-photon peak. For reasons not yet wellunderstood, we cannot find significant signals of a mix of one-plasmon–three-photon orthree-plasmon–one-photon processes. 67.3.5.2. Towards time-resolved PES
The above example dealt with the excitationof a resonance mode through a close, although not exactly matching, laser frequency.Resonant modes may also be excited by laser pulse whose frequency ω las is far away fromthe resonance, but where a multiple of ω las coincides with the mode. Such a situationshould also leave traces in the PES. We exemplify that for the case of Na irradiated bya laser polarized along the symmetry axis of the cluster, denoted by z , a pulse durationof 120 fs, a frequency of 1.1 eV and an intensity of 3 . × W/cm . The Na clusterconsists in two squares parallel to a plane (denoted by x and y ) and which are twistedby 45 ◦ around the z axis. It possesses three states of energies ε s = − .
75 eV, ε p xy = − . ε p z = − . ω sat = 3 . ω las = 1 . ω sat . It persists PE S ( a r b . un i t s ) E kin (eV) ε s + ω s a t ε p + ω s a t ε s + ω s a t ε p + ω s a t ε s + ω l a s ε p + ω l a s ε p + ω l a s -0.4-0.3-0.2-0.100.10.20.30.40.5 0 50 100 150 200 D z ( a ) Time (fs)Na , I = 3 . × W/cm ω las = 1 . eV ω sat = 3 . eV N esc ’ . full t > fs t < fs Fig. 39. Electronic dynamics of Na after irradiation by a laser with polarization along z , duration of120 fs, frequency of ω las = 1 . . × W/cm (calculations with pseudopotentialsand in full 3D). Left : Time evolution of electronic dipole along the symmetry axis of Na , denotedby z , and of ionization N esc . The horizontal bars emphasize time spans dominated by the indicatedfrequencies, that is by ω las below 90 fs and by a satellite frequency ω sat = 3 . s , the degenerate 1 p xy and the1 p z energies shifted by multiples of ω las (lower lines) or of ω sat (upper lines). All single electron energieshave been down-shifted by 0.18 eV to account for the Coulomb shift due to the total ionization of 0.06at the end of the simulation time. Adapted from [194]. even after 120 fs when the laser is switched off, since sizable oscillations remain at thishigher frequency. One should moreover notice that ω sat (cid:39) ω las , which indicates thatabsorption of three photons from the laser pulse triggers this excitation. It is a typical68xample for higher harmonic generation. The persistence of dipole oscillations provokesa continuous electronic emission. Thus the total ionization N esc , also plotted in the leftpanel of Fig. 39 as a red curve, does not level off after 120 fs but rather increases with aconstant slope, reaching the value of 0.06 at the end of the simulation time.It is interesting to observe how the PES builds up in time in such a case. To this end,we plot in the right panel of the PES calculated for the full time span, and compareit to that calculated for the first 120 fs (green lower curve) and that after 120 fs (blueupper curve). One observes a slight down-shift of the peaks from the early to the latetime windows. This provides a time-resolved illustration of the Coulomb shift (discussedin Sec. 4.3.4). Due to the final N esc = 0 .
06, an average Coulomb shift of the s.p. energiesof δ(cid:15) = − .
18 eV emerges. Hence the latter energies, indicated by vertical lines, havebeen shifted by δ(cid:15) = − .
18 eV to achieve a better matching with the PES peaks. Thefull PES (black middle curve) shows patterns repeated with equal spacing which are, atfirst glance, MPI peaks, as seen before for C (see right panel of Fig. 36). The first peaknear zero kinetic energy is related to a four-photon process emitting out of the 1 p z state.Most of the peaks in this full PES can be identified with E kin = ε i + νω las , as indicatedby the bottom vertical lines. There are, however, further peaks not explained in termsof photon frequency. To disentangle the peaks, we have also evaluated the PES in twotime windows, an early one during the laser pulse, i.e. 0-120 fs, and a late one after thepulse is over. The PES for the early window (lower green curve) is fully explained byMPI with the laser frequency. The PES from the late window (upper blue curve) showssharp peaks which can be identified as multi-resonance peaks located at ε i + µω sat , with µ = 2 , δ(cid:15) ). No MPI peak from the lasershows up in the late window. Note also that the multi-resonance peaks already slightlydevelop during the laser pulse (see vertical lines from above).In this example, we have thus demonstrated that a dynamical competition betweenvarious frequencies, here that from the laser pulse and that from a higher resonancematching the third harmonics of the laser, can provide mixed mappings of the PES andcan thus give rise to a complex structure of the PES. A time-resolved PES analysis canbe a way to disentangle the different contributions. To that end, even a coarse timeresolution as performed here may be sufficient.4.4. Photoelectron angular distributions (PAD): a sensitive tool
This section is devoted to PAD. We are using free clusters as examples. Thus weconsider orientation averaged PAD throughout, see Sec. 3.4. Remind that these can beexpanded in terms of Legendre polynomials P k (cos θ ) according to Eq. (37). The ex-pansion parameters β k carry all information about the orientation averaged PAD. Thelargest non-vanishing β ν is related to the number ν of photons involved in the pro-cess. The most important parameter is the anisotropy β which is also the only relevantparameter for one-photon processes. Therefore, β will play a key role in the followingpresentations.4.4.1. One-photon regime
Stationary state picture and Bethe-Cooper-Zare formula
The study of PADhas a long standing history, especially in atoms. Early works by Bethe [195] and Cooper69nd Zare [196, 197] are still routinely used in today’s cluster literature [17, 198]. TheBethe-Cooper-Zare formula delivers a compact expression for the evaluation of β inspherical potentials. At variance with our standard way of evaluating PAD (see Sec. 3.3.4),it provides a stationary state picture which thus requires evaluation of both bound and continuum electronic states to describe initial and final electronic states. By construction,it does not include possible electronic rearrangement effects following electronic emission,as was demonstrated in [194]. This limits the applicability to cases where electronic rear-rangement (through Coulomb residual interaction) can be neglected. Furthermore, beingdeveloped for atomic physics, the Bethe-Cooper-Zare is strictly limited to spherical ex-ternal potentials. Nonetheless, it may be useful as zeroth order estimate and reference.The Bethe-Cooper-Zare formula was first derived in first-order perturbation theory for one-electron atoms in [195]. But it can also be applied to many-electron systemsin an independent-state picture where the many-body wave function is a simple anti-symmetrized product of s.p. orbitals [196, 197]. For a given electronic level i , the totalcross-section σ ( i ) for emission and its anisotropy β ( i )2 are given by [199, 200]: σ ( i ) = (4 π ) N · L R − + ( L + 1) R (2 L + 1) , (49) β ( i )2 = L ( L − R − + ( L +1)( L +2) R − L ( L +1) R − R + cos ∆(2 L + 1)[ L R − + ( L + 1) R ] , (50)with R ± = (cid:90) ∞ d r r R ( f ) L ± ( r ) R ( i ) L ( r ) and ∆ = ∆ L +1 − ∆ L − , (51) N = 4 π e ω las (cid:126) c . where L is the angular momentum of the initial state. Once given the (spherical) po-tential, the radial wave functions of bound state R ( i ) L and continuum state R ( f ) L ± can becalculated by solving the associated radial Schr¨odinger equation. The phases ∆ L ± en-tering the continuum states can be obtained in a standard manner from the asymptoticbehavior of the outgoing wave R ( f ) l [201].Note that in the particular case of a spherical wave function ( L = 0), the Bethe-Cooper-Zare formula exactly delivers β ( s )2 = 2 (maximum possible value of β in theone-photon domain). In this case, the angular distribution of s states is not influencedby the radial form of the outgoing wave, so that the potential does not affect the angulardistribution; it only impacts the cross-section (49).In spite (or maybe because) of its simple compact form, the Bethe-Cooper-Zare for-mula has thus to be taken with a grain of caution in realistic cases, because of its strongdependence on the shape of bound and unbound electronic wave functions (see Fig. 40 be-low). Furthermore, it remains a stationary state picture thus well inside the perturbativeregime, when applicable (spherical potential). Its range of application is thus limited.4.4.1.2. On the sensitivity of β to model assumptions There are several approxima-tions around in the description of clusters and molecules. The above mentioned Bethe-Cooper-Zare formula for instance forces spherical symmetry and neglects electronic re-arrangement. Fully dynamical calculations often employ reduced symmetries as, e.g., in70APS or by using a jellium model for the ionic background. Many of these approxima-tions are validated for describing spectra and global emission properties. However, PADis very sensitive to this kind of theoretical details, and one has to check carefully theimpact of approximations.We test the sensitivity for the simplest case of one-photon processes which are fullycharacterized by the anisotropy β . To have a systematic test, we study variations of β as a function of laser frequency. We take as a test case Na and consider β (1 p )2 ( ω las ),that is, the anisotropy parameter for emission out of the 1 p state [202]. To explore theimpact of dynamical rearrangement effects in the PAD, the left panel of Fig. 40 comparesresults of a fully dynamical TDLDA-ADSIC calculation (full blue line) with the resultof the Bethe-Cooper-Zare formula or, alternatively, with a TDLDA calculation in whichthe electrons are propagated in the frozen ground-state Kohn-Sham potential (dashedred line). Note that the result of the Bethe-Cooper-Zare formula is basically identical to -0.500.511.52 5 10 15 20 25 β ( p ) ω las (eV) Na Jellium β (1 s )2 β ( i ) ω las (eV) Na Pseudopotential full TDLDA pot.Bethe-Cooper-Zare jellium, 1 p p z p x , y s Fig. 40. Anisotropy parameters β ( i )2 of single particle states i in Na , as a function of laser frequency ω las . Left : β of the 1 p shell in Na described by a spherical jellium ionic background (Wigner-Seitzradius r s = 3 .
65 a , surface thickness σ = 1 a ), calculated in TDLDA-ADSIC (full blue line), and in aBethe-Cooper-Zare approach, see Eq. (50) (pink dotted curve). For the 1 s level, β (1 s )2 = 2 at any ω las .Right : β of the 1 s , 1 p x,y , and 1 p z shells in Na described by explicit ions and pseudopotentials, and av-eraged over six orientations. For a better comparison, β (1 p )2 from the jellium calculation is superimposed.Adapted from [202]. a dynamical calculation with the KS potential kept fixed at its static form and driven atvery low laser intensity to stay safely in the one-photon regime. The laser pulse durationis of 60 fs and its intensity is scaled with the frequency ( I = 10 W/cm × ω las ) to keepthe total ionization in the range between 10 − and 0.1, and thus to stay in a perturbativeregime. The laser frequency is varied between 4.1 and 29 eV, so that only one photon isneeded to promote the electron from the 1 p state into the continuum ( ε p = − .
08 eV).The ionic background of Na is treated by a spherical jellium (Wigner-Seitz radius r s =3 .
65 a , surface thickness σ = 1 a ) which provides exactly the atomic situation for whichthe Bethe-Cooper-Zare was developed. Remind that β (1 p )2 can vary between − β (1 p )2 is close to2 for most frequencies. That means that the photoelectrons are mostly emitted alongthe polarization axis of the laser. There are however frequencies ω las at which we findpronounced dips down to negative values. Qualitatively, the pattern of the two casesare similar. However, the deep dips in β (1 p )2 ( ω las ) occur at very different places. Thisshows that dynamical effects, as the interaction of the photoelectrons with the residual71luster, strongly impact the PAD. A purely perturbative formula is therefore dangerousfor molecules which develop remarked rearrangement effects, as e.g. metal clusters.Moreover, the anisotropy parameter is very sensitive to the ionic background itself.This is demonstrated in the right panel of Fig. 40 where the jellium results are comparedto calculations with explicit ions and pseudopotentials (see Sec. 3.1.1). The jellium modelwas tuned such that both models for the ionic background provide about the same IP(4.08 eV for the jellium and 4.28 eV for the pseudopotentials). The spherical jelliumdelivers two occupied states, a 1 s states with two occupancies and a degenerate 1 p stateholding six electrons. The non-spherical ionic structure breaks the degeneracy into a1 p z state and two still degenerated 1 p xy states, and delivers a 1 s which is not perfectlyspherical anymore. To extract a sound β from the PAD, we apply orientation averagingwith the six reference orientations appropriate for one-photon processes (see Sec. 3.4).The β of these three states exhibit only faint dependence on ω las , and stay close to 1.7.The jellium model, on the contrary, systematically delivers higher values of β , with theexception of a few marked dips. It seems that the marked structures from the highlysymmetric jellium model are averaged out to a nearly constant anisotropy. This can beexplained by the rescattering of the photoelectrons on the ionic structure before leavingthe cluster.4.4.1.3. More on the dependence of β on laser frequency The results discussed inSec. 4.4.1.2 suggest that the ionic structure washes out strong variations in the frequencydependence of the anisotropy. We will address here two exceptions from this generalobservation.As a first example, we consider Na − for which experimental PAD exist [203]. Weshould mention that this case is numerically extremely demanding, because of the nega-tive charge and subsequently low IP (1.43 eV). We had to use a huge numerical box (160 mesh points and an overall box length of 280 a ). The left panel of Fig. 41 comparesthe calculated anisotropy for emission out of the group of 1 p states, β (1 p )2 ( ω las ), with theexperimental data. All laser frequencies shown in the figure correspond to emission fromthe 1 p states in the one-photon regime. The experimental data show three curves. Theseare associated with the three sub-peaks of the 1 p states found in the experimental PES,see lower right panel. The theoretical calculations did not disentangle these sub-peaksand show the β (1 p )2 from the PAD averaged over the whole 1 p group. The theoretical and(averaged) experimental curves nicely agree with each other. For higher frequencies, wesee again the smooth trend as was already observed in the example of Na in Fig. 40.The data stay systematically a bit below the theoretical β (1 p )2 . This is probably due toelectronic collisions not accounted for in TDLDA. The great surprise is the deep dip atlow frequencies which is not an artifact because it is also clearly seen in the experimentalresults. It is due to a very special situation for this loosely bound anion. Indeed, weare near threshold, and the electron cloud is thus emitted with near zero momentum.The KS potential seen by the escaping electron is extremely shallow and the outgoingelectronic wave function has an extremely long wavelength throughout. Thus it cannotresolve the ionic structure and the rescattering mechanism which wipes out the dips (seediscussion of Fig. 40) becomes obsolete. This is demonstrated in the upper right panelof Fig. 41 where we compare β (1 p )2 ( ω las ) for jellium and explicit ionic background over alarger frequency span. As in the case of Na , the jellium model produces values near 272 β ( p ) ω las (eV) Na − ω las (eV) β (1 p )2 -3 -2 -1 PE S binding energy (eV)1 s p γ p β p α exp. 1 p α exp. 1 p β exp. 1 p γ theo. ionic jelliumionic ’na7m-beta2.gnu’ i 2 u 1:2 Fig. 41. Top right : comparison of calculated β (1 p )2 of Na − with a jellium background and an explicitionic one, after irradiation by a laser with I = 10 W/cm and T pulse = 60 fs. Bottom right : experimentalPES of Na − irradiated by a laser of intensity < W/cm and duration of about 10 ns [203]. Left :Experimental (open symbols, [203]) and theoretical (full curve, [128]) anisotropy parameter of the 1 p states of Na − , as a function of laser frequency. 1 p α , 1 p β and 1 p γ correspond to state assignments ofpeaks observed in the experimental PES (bottom right). and a pronounced dip around 11 eV, while ionic structure delivers a generally smoothercurve with a maximum around 1.5. But both, jellium and ionic background, deliver thesame deep dip towards threshold. This confirms nicely that the extremely long wave-length of the outgoing electron state reduces the spatial resolution such that the softjellium and detailed ions cannot be distinguished anymore. This is also the reason whyan older calculation of β for Na − at low frequencies using a jellium model and linearresponse could provide realistic results [204].In contrast to the PES which exhibits a strong dependence on laser intensity, anorientation-averaged PAD seems to be not very sensitive to it. Indeed, when irradiatedby laser pulses of duration of 30 fs and frequency of 34 eV but with two different intensities(10 and 10 W/cm ), the obtained PAD both delivers the same β = 0 .
38, while thePES at the highest intensity is blurred and red shifted [190]. Note that the anisotropyparameter is here much smaller than for small Na clusters, see Fig. 43. This is again dueto the influence of the ionic structure : for sodium clusters, β decreases by about 25 %when going from jellium to ionic background. This effect should be even stronger in C ,since the number of ions is much higher than for the considered Na N with N = 3 − β . An orientationaveraging procedure is applied here (see Sec. 3.4.1). We focus on the photo-emission fromHOMO and HOMO − β are plotted in the right panel. The total β (blackcurve) does not depend on ω las very much. The anisotropy parameter of the HOMO73 β ω las (eV) C , I = 7 . × W/cm , T pulse = 60 fs ω las = 14 eV P h o t oa n g u l a r d i s t r i bu t i o n s ( a r b . un i t s ) Angle ( ◦ ) ω las = 26 eV Total HOMO HOMO − Fig. 42. Left : Calculated orientation-averaged PAD from C with radius of 6.763 a irradiated by alaser pulse of intensity I = 7 . × W/cm , duration of 60 fs, and frequency ω las of 14 eV (top) or26 eV (bottom). Right : anisotropy parameter β as a function of ω las . Red curves: total PAD and β from all single particle states ; light green curves : the same but from the HOMO only ; blue curves :the same but from the HOMO − and HOMO − β of the HOMO − β of the HOMO (green curve) steadilydecreases with ω las and changes of sign a bit before 24 eV. Below 22 eV, the β of theHOMO is also higher than the one of the HOMO −
1. The case once more demonstratesthe extreme sensitivity of the β as an observable characterizing a dynamical scenario.Mind that, in this monophoton domain, the anisotropy parameter is bound between − β iscertainly a very rich quantity to be measured and computed in a highly refined manner.4.4.1.4. Dependence of β on cluster deformation To explore the impact of clusterdeformation on the PAD, we now turn to a series of small neutral and cationic metalclusters which cover planar (Na ), prolate (Na , Na ), oblate (Na and Na ), andtriaxial (Na , Na ) systems. We consider detailed ionic background as well as a de-formed jellium approach to it. The jellium deformation is tuned in each case to reproducethe global deformation of the ionic configuration. The shape can be quantified by thequadrupole deformation α defined by α = (cid:113)(cid:80) m = − α m with α m = 4 πr Y m / N R , R rms the ionic root mean square radius, N the total number of ions, and Y m the spher-ical harmonics for l = 2. In Fig. 43, we compare α with the total anisotropy parameter β . As in the case of Na previously discussed, β is extracted from orientation aver-aged PAD, each PAD obtained after irradiation by a laser in the mono-photon regime,that is ω las = 7 . I = 10 W/cm . The total ionization always remains between 10 − and 10 − . We first note that β shows only small variations, particularly for the jelliummodel. We cannot spot any correlation of the anisotropy β with the deformation α . Thedifference between neutral clusters and cations is also very small. However, as observedpreviously, there can be a large sensitivity to the structure of the ionic background. The74 .41.51.61.71.81.92 3 6 9 12 15 18 A n i s o t r o p y β number of ions N D e f o r m a t i o n α Jellium: Na N + PsP: Na N + Jellium: Na N PsP: Na N Fig. 43. Comparison of quadrupole deformation α (top) and total anisotropy parameter β (bottom)for small neutral (circles) and cationic (squares) clusters, obtained in a jellium description of the ionicbackground (open symbols) or with an explicit ionic structure (closed symbols). Adapted from [128]. β is systematically smaller (more isotropic) when detailed ions are considered. It is par-ticularly interesting to note that the discrepancy grows systematically with increasingcluster size N . This trend is corroborated by results from larger clusters. For example,the total anisotropy for C comes close to zero, see Fig. 44 for low intensities. Thetrend complies with the interpretation that rescattering with ions enhances the isotropicbackground: the more scatterers, the closer to isotropy.4.4.2. PAD in the multiphoton regime
The orientation averaged PAD in the multiphoton regime develops more detailed an-gular dependence as indicated by Eq. (37). For clarity, we recall it here :d σ dΩ ∝ β P (cos ϑ ) + β P (cos ϑ ) + . . . In a strictly perturbative regime (low laser intensity), the series terminates at P ν (cos ϑ )where ν is the order of the multiphoton process which is determined by the relation of IPto photon frequency. This changes with increasing intensity where always all amounts ofphotons could be possible. Thus the series is, in principle, unterminated and we expectthat the contributions of higher β n increase with increasing intensity I . We thus brieflyanalyze the impact of intensity on the PAD in two emblematic cases, C and Na . Aswe are going beyond the perturbative regime, we will consider higher β n beyond theanisotropy β . It is to be noted that orientation averaging in the multiphoton regime hasto be done explicitly by integration over orientations. Due to the high symmetry of thetwo test cases, only 18 integration points need to be computed.We start with the case of C irradiated by a laser of frequency of 1.55 eV and pulseduration of 75 fs, with the same set of increasing intensities as in use for the systematics ofPES in the right panel of Fig. 37. The IP of C being 8 eV, it requires at least 6 photons75o extract electrons from occupied states. The PAD obtained in C with different laserintensities are shown in the left panel of Fig. 44, while the first anisotropy parameters β n extracted from these PAD are plotted in the right panel. The total ionization N esc (red P A D ( a r b . un i t s ) Angle (deg) C , ω las = 1 . eV, T pulse = 75 fs I = 1 . × W/cm . I . I . I . I W/cm ) β β β N esc Fig. 44. Anisotropy parameters of high orders (right) extracted from PAD (left) of orientation-averagedC irradiated by laser pulses of frequency of 1.55 eV, duration of 75 fs, for different laser intensities.For the sake of completeness, the total ionization N esc is plotted in the right panel. Adapted from [190]. boxes and dashed line in the right panel) increases rapidly with I as expected. The PADshown in the left panel become more and more aligned along the laser polarization withincreasing I . Accordingly, all β l , shown in the right panel, increase with I . Note that β > − ≤ β ≤ cluster. Differently as in previous sections, we now runit for laser frequencies below ionization threshold, namely ω las = 3 . I . Fig. 45 shows the intensity dependence of β (leftpanel) and β (right panel). The total ionization ranges from 0.001 to 0.1, which means β laser intensity (W/cm ) Na -0.09-0.06-0.0300.03 β laser intensity (W/cm ) Na Fig. 45. Anisotropy parameters β (left) and β (right) of Na extracted from PAD with an averagingprocedure over 18 orientations, after irradiation of laser pulses of duration of 60 fs and frequencies of 3.7(red full curves) or 3.9 eV (blue dashes). Adapted from [128]. that all cases constitute rather moderate excitation dynamics. We see again generally76n increase of total anisotropy β with intensity. Still, the growth of β is slower thanin the previous example C because we deal here only with two-photon processes. Thenext coefficient β also shows a marked trend to larger negative values with increasingintensity. It is interesting to note that both coefficients, β and β , are very sensitive tothe laser frequency. This is also a feature of MPI while frequency dependences are moremoderate for one-photon processes, see e.g. Sec. 4.4.1.2.4.5. Impact of temperature in PES and PAD
Effect of ionic motion on PES and PAD
So far, the presented calculations of PES and PAD were performed at ionic ground-state configuration, i.e. at a temperature of 0 K. This ideal situation is hardly everfeasible in an experiment. Depending on the production conditions, cluster beams havetemperatures in the range of several 100 K. This means that we encounter usually anensemble of ionic configurations fluctuating around the ground-state configuration. Inthe following, we will discuss the impact of thermal shape fluctuations on PES/PAD.Experimentalists are well aware of the temperature problem and have developed severaltechniques for dedicated cooling of cluster beams. Ion traps are particularly powerfuldevices for a clean handling of cluster beams [206]. We show here results from recentexperiments which used a trap and cooling with a He buffer gas to produce beams of Naanionic clusters with well defined temperatures between 6 and 265 K [203]. The uppertemperature is above the melting point of ≈
250 K for small Na clusters [207]. The testcase Na − is a cluster anion which has naturally a low IP. Thus one can easily realizeone-photon processes with standard laser pulses. A selection of combined PES/PAD forNa − is shown in the left panel of Fig. 46. The upper right panel shows the PES atdifferent temperatures T and the lower right panel the anisotropy β as a function of T .The PES/PAD in the left panels and the PES in the upper right panel (both analyzed atlow temperature) allow one to identify five emitting states, labeled A, B, C, D, and E. Thestructures produced by these states in the PES (and PES/PAD) are gradually blurredwith increasing temperature. This is easily understood from the fact that single-electronenergies can be very sensitive to changes in the cluster shape. The thermal ensemblethus represents a more or less broadened distribution of s.p. energies which, in turn, ismapped into PES and PES/PAD. On the other hand, the anisotropies β exhibit onlya weak dependence on T . This probably reflects the fact that the angular momentumcharacteristics (stemming from their wave functions) of the s.p. states are more robustthan their energies.The numerical simulation of a thermal ensemble is conceptually straightforward, al-though somewhat cumbersome. One starts an ionic dynamics from the ground-stateconfiguration by initializing ionic velocities stochastically according to a Maxwellian dis-tribution for the given temperature T . This state is then propagated by TDLDA-MDfor a few ps. About each 100 fs, a snapshot of the actual configuration is taken. Theset of all snapshots constitutes the thermal ensemble of cluster configurations. Now, in alaser-induced dynamics propagated for each sample, the wanted observables (e.g., PESand PAD) are evaluated and incoherently superimposed. This altogether yields the ob-servable for the ensemble. We have performed such a study for Na irradiated by laserpulses of intensity of 10 W/cm , and FWHM of 232 fs. This pulse duration allows a77 inetic energy (eV)
6K 100K 265K Na - A B C D E 0.0 0.2 0.4 0.6 0.8 1.0 A B C D E 0 o o o o o a ng l e P E S ( a r b . un it s )
0 50 100 150 200 250 temperature (K) 0 -1 1 2 a n i s o t r opy p a r a m e t e r Fig. 46. Left : Combined PES/PAD of Na − obtained by fit of five Gaussian peaks to the raw distri-butions at three different ionic temperatures as indicated. Right top : extracted PES for temperaturesranging from 6 K to 265 K. Right bottom : extracted anisotropy parameter β as a function of ionictemperature. Adapted from [203]. high resolution of the PES peaks (see the effect of pulse duration in Fig. 34) such thatline broadening comes predominantly from thermal effects. The pulse duration lies withinthe time scale of ionic motion. Thus the dynamical propagation is done at the level ofTDLDA-MD to include properly ionic motion. For reasons of simplicity, we are not per-forming orientation averaging such that we see exclusively the thermal effects. The laserpolarization is chosen along the symmetry axis of the T = 0 configuration. For a firsttest, we use a laser frequency at 6.8 eV just below the IP. The 1 s state ( ε s = − . p states ( ε p = − . T = 0 and at three increasingtemperatures T = 158, 315 and 473 K. As expected, the higher T , the broader the PESpeaks. There is also a faint red shift of the peaks with increasing temperature. Nonethe-less, it is surprising how well the structures survive in the PES even high above meltingtemperature (about 250 K). What the PAD is concerned, remind that it is computedwithout orientation averaging. It thus shows more structure and is asymmetric (reflect-ing the asymmetry of Na ). The thermal effects on the PAD are a bit larger than forthe PES, but remain still small, in accordance with the experimental results for Na − in Fig. 46. Fig. 47 also shows through the error bars the uncertainty associated withthermal fluctuations. These are computed in standard manner as the variance of PADand of the logarithm of the PES yield from the statistical ensemble. The uncertainties78 − − − − − − − − − − − − P h o t o e l ec t r o n s p ec t r a ( a r b . un i t s ) − − − − − − E kin (eV) 00.20.40.60.8 Na I = 10 W/cm FWHM=232 fs ω las = 6 . eV 00.20.40.60.8 P h o t oa n g u l a r d i s t r i bu t i o n s ( a r b . un i t s ) ◦ ) T = 158 K T = 0 K T = 315 K T = 0 K T = 473 K T = 0 K Fig. 47. Left column : Photoelectron spectra of Na of fixed orientation, irradiated by laser pulsespolarized along the symmetry axis of the cluster, with FWHM of 232 fs, intensity of 10 W/cm , andfrequency of 6.8 eV, at an ionic temperature T of 158 K (top), 315 K (middle) or 473 K (bottom). Theblue curve shows the PES at T = 0 K. Right column : corresponding PAD. grow with temperature. They stay rather small for the PES, showing once more thatthese structures are rather robust. The error bars are larger for the PAD in forward (0 ◦ )and backward (180 ◦ ) direction. Fortunately, these forward and backwards cones have asmall integration weight, such that global measures as, e.g., the anisotropy β are againrobust.We now concentrate on the case of T = 315 K and complement the ω las = 6 . W/cm for the two lower ω las and 10 W/cm for ω las = 13 . T = 0 case). For the highest frequency (right panels), the impact oftemperature is quite weak: we still observe a slight broadening of the peaks in the PES.However, the effect on the PAD is negligible. This result is intuitive because for such ahigh laser frequency, the photoelectrons are extracted by absorption of a single photonand basically follow the laser field, as is visible from the fact that the PAD is peaked79 − − − − − − PE S ( a r b . un i t s ) E kin (eV) 0 5 10 15 Na , FWHM=232 fs, T = 315 K E kin (eV) 0 5 10 15 20 E kin (eV)-0.400.40.81.21.622.42.83.2 0 30 60 90 120 150 P A D / Y Angles ( ◦ ) Y = 0 . I = 10 W/cm ◦ ) Y = 0 . I = 10 W/cm ◦ ) Y = 0 . I = 10 W/cm ω las = 3 . T = 0 K ω las = 6 . T = 0 K ω las = 13 . T = 0 K Fig. 48. Top row : Photoelectron spectra of Na with initial ionic temperature of 315 K, irradiatedby laser pulses of FWHM of 232 fs, intensity of 10 W/cm , and frequency of 3.4 eV (left), 6.8 eV(middle), and 13.6 eV (right). The blue curve shows the PES at T = 0 K. Bottom row : correspondingPAD. along ϑ = 0 ◦ and 180 ◦ . On the contrary, the lowest frequency (left panels) lies deeplyin the multi-photon regime. And the uncertainties produced by the ionic motion at thistemperature are extremely large. Apart from that, MPI peaks are still visible in the PES,although much broadened. At the side of the PAD, the uncertainties are larger than thesignal. This means that in the multi-photon regime, the PAD is very sensitive to ionictemperature. The intermediate frequency (middle panels) lies in between in all respects.Temperature effects are already well visible, but not as disastrous as for low frequency.This example shows that the ionic temperature should better be well controlled and keptat a sufficiently low value to allow a quantitative analysis of PAD.4.5.2. Impact of electronic dissipation on PES and PAD
In Sec. 4.3.4, we found that a smooth exponential PES develops for sufficiently highintensities, see Fig. 37. Figure 10 did also show a series of exponential PES from ameasurement with high intensities (fluences) and a long pulse. The question is to whatextent this could be a signature of thermal emission after full thermalization of thecluster.We start the discussion with looking back at Fig. 37. From the inverse slope of theexponential, we would read off an apparent temperature of T app = 1 . E ∗ (cid:39) . E ∗ / (cid:39) .
01 eV which is twoorders of magnitude smaller than T app . This, together with the fact that TDLDA doesnot account for electronic thermal effects, clearly rules out the thermal origin of theobserved exponential slope.The experimental results from [47] in Fig. 10 (see Sec. 2.3.1) did show also a remarkableseries of smooth exponential PES. Fig. 49 gathers the apparent temperatures (inverse Fig. 49. Experimental apparent electronic temperature as a function of laser fluence as extracted fromPES for C (squares) and C (circles). The dotted line is a fit to the temperatures vs laser fluenceobtained in the case of C in [52]. Note that the laser fluence has an overall systematic uncertainty of10%, which is not included in the present error bars. From [47]. slopes) extracted from these PES. Again we find rather large values in the range from 1up to 1.6 eV, near the 1.8 eV of the previous example. The remaining excitation energy isnot available for these experiments. Nevertheless, it is questionable that these exponentialslopes should correspond to real temperatures of the system in thermal equilibrium. Thisis why one wisely has coined the notion “apparent temperature” for the slope of thePES, see Fig. 10 from [47]. The example of Fig. 37 has shown that a smooth exponentialpattern can also be explained by TDLDA calculations. In fact, the exponential profilecan be nicely fully explained in terms of multiphoton perturbation theory [50]. Any MPIyields an exponentially decreasing slope. The Coulomb shift increasing with ionizationstage increasingly washes out the MPI peaks to yield eventually a purely exponentialPES.In order to check how far pure TDLDA (free of any thermalization) can describeexponential PES, we consider Na for which calculations can be compared with ex-perimental data. In this experiment both total ionization and PES have been measured.The exponential shape of the PES was attributed to thermal effects [208]. The resultingapparent temperature (inverse slope) and total ionization are plotted as a function oflaser intensity in Fig. 50. The experimental results are compared to standard TDLDAcalculations performed under the same laser conditions. The theoretical results are sur-81 s l o p e ( e V − ) log ( I in W / cm )110 N e s c Na ω las = 3 . Fig. 50. Ionization characteristics of Na irradiated by lasers of increasing intensities (in W/cm ) butfixed photon frequency at 3.1 eV and pulse FWHM at 200 fs. Calculations (red circles) are comparedto experimental data (green squares) from [208]. Ionic dynamics has been included with a Maxwelliandistribution of velocities according to a temperature of 100 K. Upper panel: total ionization in logarithmicscale. Lower panel: slope of the PES. Adapted from [50]. prisingly close to the data, which indicates that a purely thermal interpretation is notcompelling. Consider the apparent temperatures (inverse of the slopes). They are in range0.7–2 eV for TDLDA and 0.6–0.8 eV for the experiment. This would amount to about300 eV intrinsic energy in case of full thermal equilibrium ( ≈
150 eV if only electronswere thermalized) : this is too much as compared to the typical cluster binding in Na .Thus we are surely far from full thermalization in this case. On the other hand, there aresmall, but systematic, differences which may be a trace of thermal effects. The apparenttemperature (lower panel) and the ionization (upper panel) is somewhat lower in thedata than in TDLDA. This indicates that the data represent, in fact, a mixed situation,not fully thermalized yet, but somewhere on the way.There still remains the task to distinguish direct from thermal electron emission. Wehave argued above that exponential PES are only a necessary condition, not a suffi-cient proof. Additional information for a better discrimination is delivered by the PAD.Isotropic PAD are another necessary, usually much more conclusive, condition for thermalemission. For example, the experimental PES/PAD in Fig. 24 shows a larger inner spot ofisotropic emission at low energies which can be associated with thermal electrons. At thetheoretical side, a relevant description of thermalization requires dynamical correlationsbeyond TDLDA. This has been achieved at the semi-classical level, see Sec. 4.5.3. It is82till a great challenge for a fully quantum mechanical modeling, see Secs. 5.2.2 and 5.2.3.In the following, we will briefly present a simple estimate for the contribution of thermalelectrons which provides at least a first impression of the impact of thermal electrons onPAD.The case of C , illustrated in Fig. 51, addresses the complementing PAD signal, whichis expected to contain a significant isotropic component if strong thermal effects arepresent. We consider here the MPI regime, irradiated by a laser of frequency 1.55 eV, P A D [ a r b . un i t s ] Angle ϑ ( ◦ ) Exp.Theo., direct emissionTheo., direct+thermal
Fig. 51. Comparison of experimental photoangular distribution from C HOMO and HOMO − pulse FWHM of 20 fs and intensity of 1 . × W/cm . The experiments deliver aPAD energy-integrated over the HOMO, HOMO −
1, and HOMO − is 7.8 eV. The 1.8 eV excitation energy thus suffices to emit about0.22 electrons (when neglecting the possible C dissociation channel which requires largerenergy). We then add up the contribution of the extra 0.22 emitted electrons as a thermal,isotropic background to the PAD. This leads to the curve labeled “direct+thermal” inFig. 51, which now agrees fairly well with the experimental curve. Of course, the reasoningbasically provides an argument and does not constitute a theory by itself, but it onceagain confirms that bare TDLDA underestimates electronic collisions, while they areobviously non negligible in the MPI regime.83.5.3. The semi-classical route
As already discussed in Sec. 3.5, thermal effects in finite systems can presently onlybe attained via a semi-classical approximation, leading to semi-classical kinetic equationsuch as VUU. The underlying Vlasov equation is the semi-classical limit of TDLDA.VUU employs additionally the Uehling-Uhlenbeck collision term which accounts for thedynamical electron-electron correlations. Ionization, as a basic dynamical mechanism,has been discussed within VUU in several papers [209, 210, 60, 169, 141]. The interestingpoint for the present discussion is to analyze the impact of VUU as compared to Vlasov.We shall discuss the point on the example of Na . We fix the pulse duration at aFWHM of 25 fs. We consider two laser intensities, 10 W/cm and 6 × W/cm , andthree frequencies, ω las = 2 .
7, 3.0, and 3.3 eV, around the plasmon frequency of the system(3 eV). The time evolution of the ionization is shown in Fig. 52. It strongly depends on Na , FWHM=25 fs I = 10 W/cm ω las = 2 . eV02468 T o t a li o n i z a t i o n ω las = 3 . eV00.51.0 0 30 60 90 120 150 180Time (fs) ω las = 3 . eV 0246810 I = 6 × W/cm T o t a li o n i z a t i o n Fig. 52. Comparison of the time evolution of the total ionization in Vlasov (brown lines) and VUU (lightgreen lines) approaches. The system is Na irradiated by lasers of cos profiles with FWHM of 25 fs.Intensities are 10 (left panels) and 6 × W/cm (left panels). Frequencies are varied between 2.7and 3.3 eV, which covers the dominant optical response peak located around 3 eV for this system. ω las , as expected at the passage of a resonance [211]. The most interesting feature is theshape of the emission profile. In the Vlasov calculations, ionization grows very quicklyin the early stages and levels off once the laser pulse is switched off. Note that the hugeionization blue-shifts the plasmon resonance such that the case ω las = 3 eV becomesoff-resonant and subsequently emission is terminated after the pulse (see Sec. 4.2.1).The VUU results look much different. Emission is suppressed in early stages. This isovercompensated by a steadily continuing emission later on. The pattern are very similarfor all three laser frequencies, while the amplitude of the effect (and the relative valuesof Vlasov and VUU ionization at early times) depends sensitively on frequency. This isa long known effect that field amplification by the plasmon resonance naturally leadsto enhanced emission [212] in TDLDA and correspondingly in Vlasov [211, 210]. VUU84hows the same resonant behavior, but modulates the time profile of emission. It reducesionization in early stages because electron-electron collisions remove energy from thedirect emission channel, and it enhances emission in later stages by releasing gently thestored energy.The analysis of PES is of limited interest in Vlasov and VUU as at low energy neithercan identify electron single particle energies because of their semi-classical nature. Theobtained PES are thus always more or less exponentially decreasing, whatever the laserconditions. More interesting is the PAD which can be easily evaluated in Vlasov and VUUand which does not suffer so much from the semi-classical approximation. Furthermore,it is an ideal observable to identify thermal effects in terms of isotropy of the PAD, aswe have seen in Sec. 4.5.2. The point is illustrated in Fig. 53 for the same test case asin Fig. 52, and for a few laser frequencies, again around plasmon frequency of Na .Orientation averaging has not been performed, which is an acceptable approximation Na , FWHM=25 fs I = 10 W/cm P h o t oa n g u l a r d i s tr i bu t i o n ◦ ) 00.511.52 I = 6 × W/cm ω las = 2 . eV 01234 P h o t aoa n g u l a r d i s tr i bu t i o n ω las = 3 . eV0 30 60 90 120 150 18000.511.52Angles ( ◦ ) ω las = 3 . eVVlasovVUU Fig. 53. Comparison of angular distributions in Vlasov (brown lines) and VUU (light green lines) ap-proaches. Same system and laser conditions as in Fig. 52. here because Na is large and close to sphericity. For all laser frequencies, the PADfrom Vlasov exhibit strongly oriented emission along the laser polarization. The effectbecomes larger when ω las comes closer and closer to the plasmon frequency. This is againa consequence of field amplification near the resonance. A similar trend with frequencyis observed in the VUU calculations, although less pronounced. More striking is the factthat the PAD from VUU are much less peaked than the Vlasov ones. Very clearly, theVUU results have a strong isotropic component. The effect is especially clear far from theplasmon resonance but remains very visible close to it. We observe then a competitionbetween field amplification (from resonance) and thermalization [211]. In any case, wesee a significant enhancement of ionization perpendicular to the laser polarization, ascompared to Vlasov, a clear signature of isotropy.85 . Future directions We have seen in the previous sections the richness and variety of observables todayaccessible in experiments and theory. Of course, there remain many open questions forfuture research. We want to illustrate in this section two lines of development whichwe consider as especially promising and which require dedicated efforts both from theexperimental and theoretical side.The first aspect again focuses on the analysis of dynamics in terms of PES and PADfrom which we have seen that it is a powerful tool. We will now consider the analysisin terms of PES and PAD in connection with the short pulses delivered by a bypassingionic projectile. As we shall see, PES and PAD can again provide useful insights into theunderlying dynamics.The second topic concerns electronic thermalization which was already addressed tosome extent in Sec. 4.5.2. The excitation energy deposited originally by the laser pulseis released in the first stages of the dynamics by direct electron emission. However, partof the deposited energy is progressively converted into incoherent electronic excitation of”thermal” nature. This takes place on a moderate time scale of some tens of fs. Analysisof such effects is difficult and requires detailed experiments (see Figs. 10 and 49). Fromthe theoretical point of view, the situation is even worse as it requires the development ofdeep extensions of available current theories such as real-time TDDFT. In the following,we address both these directions in more detail.5.1.
An excursion into irradiation by charged projectiles
We briefly discuss in this section PES and PAD of ”photo”-electrons emitted aftercollision with a fast charged projectile. We put the word ”photo” in quotation marksbecause the electromagnetic pulse has not a well defined frequency here. Although rare,there exist data measured on atoms and mono-atomic dimers with a special focus atvery high kinetic energies ( E kin >
40 eV) of the projectiles. This was first motivatedby the observation of non-monotonous patterns in the PES after irradiation of O andN with photons in the 30–60 eV range [213], which were explained theoretically oneyear later [214] by Young-type interferences between electronic wave functions of elec-trons coherently emitted from identical atomic centers. Various experiments have beenperformed on H bombarded by He + and He of 20 and 40 keV [215] or by 8 keVelectrons [216, 217], on N colliding with 1–5 MeV H + [218], and on O colliding with3.5 MeV/ u C ions [219] or 30 MeV O and O ions [220]. Very recently, collisionsof 3.5 MeV/ u C ions on uracile [221] and of 4.5 MeV/ u O ions on H O [222] havebeen reported. The kinetic energy of the ejected electrons ranges from a few eV up to600 eV, and the emission is measured between 20 ◦ and 150 ◦ .The pulses from fast projectiles are extremely short and cover a very broad band offrequencies. At first glance, this looks like a disadvantage as there is thus no specificfrequency information in the pulse. However, it has the advantage that it enables toextract unambiguously effects from the system’s modes. To illustrate this point, we startwith two pedagogical examples which have a dominant dipole mode, namely the case ofNa and Cs . The Na is described by explicit ions and pseudopotentials, while we use adeformed jellium background for the description of Cs . The interaction of the irradiated86ystem with a very fast charged projectile can be modeled by an instantaneous boost ofthe electronic wave functions at t = 0 [148]. Note that this procedure is the same as theone we use for the calculation of an optical response, see Sec. 4.1.2. More precisely, weapply a boost p to each occupied s.p. wave function ϕ j, gs of the ground state, and takethe obtained wave functions as the initial states : ϕ j ( r , t = 0) = exp(i p · r ) ϕ j, gs ( r ) (52)This mimics the effect of the Coulomb field caused by a fast by-passing charged projectile.For simplicity, let us consider a boost in the z direction only. If the projectile is fastenough, we can assume that it travels on a straight line with constant velocity v proj .Therefore, one can evaluate the net force integrated over the collision, and the latterexhibits a component only in direction to the point (here the z axis) of closest impact[223]. The induced boost for a projectile of charge Z and an impact parameter b thenreads p = 4 Ze bv proj . (53)One can see that, for a given value of the boost, larger b and/or v proj can be compensatedby increasing the charge Z of the projectile. To derive Eq. (53), the passage time of theprojectile, which reads b/v proj , should be much smaller than a typical electron reactiontime ω − . If one uses the value ω el = 1 Ry and an impact parameter b = 10 a fora rough estimate, we have the constraint v proj (cid:29) bω el = 200 a /fs. This lower valuecorresponds to a kinetic energy of 700 keV for a colliding proton. Inserting this value of v proj in Eq. (53) with b = 10 a and Z = 1 yields a maximal value max = 0 . , whichis in the range of the boosts used in the following.We now come back to the first two test cases that we have studied, namely Na andCs . Both systems exhibit a very clean plasmon peak at 2.1 eV and 1.4 eV respectively.Fig. 54 shows the obtained PES. On top of the exponential decrease emerge some peaks.One can identify the dominant ones as multi-plasmon excitations, similar to those alreadydiscussed in Sec. 4.3.5.1 : the energy given by the boost to the system is mainly stored inthe dominant dipole mode, let us call it for simplicity the “plasmon”. Since the plasmonfrequency is below the ionization threshold for each case, two or more plasmons are neededin order to ionize the system. And indeed, the double and triple plasmon processes areclearly visible in the PES. It becomes more difficult to disentangle higher orders fromthe background.We now turn to a more involved case, that is the C chain. It is described by an ex-plicit ionic structure and the 20 valence electrons are shared among 8 different electroniclevels. Note that the HOMO − − in this direction is dominatedby one single, strong and sharp resonance at ω pl = 6 . − − − − − − − − − PE S ( a r b . un i t s ) E kin (eV)Boost= 0 .
01 a − Na pseudopotential ε a + ω ( a ) p l ε a + ω ( a ) p l ε a + ω ( a ) p l E kin (eV) Cs jellium ε b + ω ( b ) p l ε b + ω ( b ) p l ε b + ω ( b ) p l ε b + ω ( b ) p l Fig. 54. Photoelectron spectra of Na described by pseudopotentials (left) and Cs described by a jelliumbackground, excited by an instantaneous boost of 0.01 a − applied to the wave functions at t = 0. Thevertical lines are positioned at energies corresponding to the single particle energy of Na (of Cs ), ε a = − .
23 eV ( ε b = − .
18 eV), blue-shifted by multiples of the plasmon frequency ω ( a )pl = 2 . ω ( b )pl = 1 . − − − − − PE S ( a r b . un i t s ) E kin (eV) ε i + 2 ω pl ε i + 3 ω pl Totalstate 8state 7state 6 +2 ω pl +2 ω pl +2 ω pl E kin (eV) 2.7 3 04590135180 A n g l e ( d e g ) Fig. 55. Electron emission from C chain after excitation by an instantaneous boost of 0.07 a − along thechain direction of the electronic wave functions at t = 0. Left panel : PES (in the longitudinal direction)of all states (black), of the HOMO (state 8, blue), HOMO − − . The strikingfeature is that, although the boost is performed along the chain, the electrons are notexclusively emitted in this direction. For instance, states 6 and 8 exhibit in addition asizable emission at 60 ◦ and 120 ◦ . On the contrary, state 7 preferentially emits at 45 ◦ , 90 ◦ and 135 ◦ . Clearly, these combined PES/PAD allow one to relate the emission behaviorto the symmetry of the depleted wave functions as in cases with lasers, see e.g. Fig. 46.5.2. Towards quantum dissipative electron dynamics
From VUU to quantum world
Although VUU (see Sec. 4.5.3) provides a way to describe dissipation in dynamicalscenarios, it is limited to large excitation energies and to simple materials as, e.g., alkalineclusters. Both limitations are direct consequences of the semi-classical nature of VUU.Simple metals can be described because their electron cloud comes close to a Fermi gas.Reducing the excitation energy or considering other systems (as e.g. C ) requires toaccount for quantum effects which renders VUU inapplicable.This calls for a quantum kinetic theory. This can be seen two ways, either as quantumcounterparts of the VUU equation (43) or as TDLDA complemented by a quantumgeneralization of the UU collision term (44). Anyway, such a theory deals with impurequantum states which are described at the level of the one-body density matrix . Thecorresponding dynamical equation for ˆ ρ readsi ∂ ˆ ρ∂t = [ˆ h, ˆ ρ ] + I coll [ˆ ρ ] . (54)The commutator with the mean-field Hamiltonian ˆ h describes the mean-field evolutionaccording to TDLDA. It is complemented by a collision term I coll [ˆ ρ ] which, however,becomes awfully involved in the quantum case [224]. There does not yet exist any routinesolution to the problem in finite systems, in spite of the many investigations, particularlyin nuclear physics [58, 135, 225]. So far, most practical solutions rely on a (partial orfull) semi-classical treatment [58]. The many detailed experiments on cluster dynamicsdiscussed in the previous sections revive the call for a manageable quantum kinetic theory.We discuss in this section two promising directions of research along that line.89.2.2. A relaxation time ansatz
VUU in the semi-classical domain and stochastic TDHF/TDLDA (see Sec. 5.2.3) de-scribe dissipation in a very detailed, thus expensive, manner. In cases of moderate fluc-tuations, the system as such remains intact and the outcome is rather obvious : thedissipative dynamics drives steadily towards thermal equilibrium of a still compact sys-tem. This suggests a simplification in terms of the relaxation-time approximation whichhad been used since long in the homogeneous electron gas [65]. An implementation forfinite fermion systems had been proposed in the nuclear context in [226]. But the com-putational limitations at that time did not allow realistic applications. Just recently, wehave taken up this old idea of a relaxation time approximation and started to implementfor cluster dynamics. We give here a brief preview of this ongoing work.The relaxation-time approximation starts from Eq. (54). The collision term I [ˆ (cid:37) ] isapproximated by i ∂ t ˆ (cid:37) − (cid:104) ˆ h, ˆ (cid:37) (cid:105) = 1 τ relax (ˆ (cid:37) − ˆ (cid:37) equil [ ρ ( r ) , j ( r )]) . (55a)The right-hand-side is the effective collision term which forces the system to convergetowards the equilibrium. Note that this employs the local equilibrium ˆ (cid:37) equil [ ρ ( r ) , j ( r ) , E ]which depends on the instantaneous local density, current and energy E from the givenˆ (cid:37) ( t ). It readsˆ (cid:37) equil = (cid:88) α | ϕ α (cid:105) n (equil) α (cid:104) ϕ α | , n (equil) α = 11 + exp(( ε α − (cid:15) F ) /T ) (55b)and can be computed with density- and current-constrained TDLDA [227, 228]. Thetemperature T is tuned iteratively such that the total energy matches the wanted value E . The key parameter is the local relaxation time τ relax for which we need a reliable choice.To that end, we recur to a semi-classical estimate of relaxation time [229, 230, 231, 232].It is based on the Fermi gas model in which the relaxation time becomes simply (cid:126) τ relax = 1615 m (cid:126) σ ee T , (55c)where k B = 1 and σ ee is the effective in-medium cross section for electron-electroncollisions. For metal clusters, we find σ ee ≈ πr s [209]. Eqs. (55) together constitute thedissipative TDLDA in relaxation-time approach.Fig. 56 shows a first result from the newly developed dissipative TDLDA scheme. Weconsider as a test case Na after instantaneous boost of its electron cloud at variousboost energies. The s.p. entropy S = (cid:80) α ( n α log n α + (1 − n α ) log(1 − n α )) is shown onthe right panel. It demonstrates most clearly the evolution towards thermal equilibrium.The global relaxation time shrinks visibly with increasing excitation. This is a generalfeature already well known from VUU. It is related to the fact that the phase spacefor transition opens up with increasing energy. The relaxation times deduced from thisfigure range from about 40 fs for low excitation down to few fs for very energetic cases.This is in range of measured values [159, 233]. The time evolutions of the dipole envelope(top left) from pure TDLDA (dashed lines) are only very slowly decaying and look atthe scale of this figure nearly constant. Activating dissipation leads to a clear decay ofthe signal : the higher the excitation, the stronger the decay (in accordance with theright panel). But note that this decay starts only after some delay while the evolution inthe early stages is very similar to TDLDA. Finally, the left lower panel shows the time90 N e s c / E ∗ ( e V ) time (fs)0.0010.010.11 d i p o l e a m p li t ud e ( a ) Na , jellium e n t r o p y S / S a s y time (fs) dissip.TDLDA, E ∗ = 0 .
068 eVdissip.TDLDA, E ∗ = 2 E ∗ dissip.TDLDA, E ∗ = 4 E ∗ dissip.TDLDA, E ∗ = 10 E ∗ TDLDA, E ∗ = E ∗ TDLDA, E ∗ = 2 E ∗ TDLDA, E ∗ = 10 E ∗ Fig. 56. Time evolution of basic observables in Na with soft spherical jellium background (12) using r s = 3 .
65 a and σ jel = 1 a , after instantaneous boost with corresponding excitation energy E ∗ asindicated : Ionization N esc normalized to E ∗ (bottom left), envelope of the dipole signal (top left), andsingle particle entropy S normalized to the asymptotic entropy S asy (bottom right). Compared are resultsfrom pure TDLDA (dashes) with those from TDLDA with dissipation in relaxation-time approximation(full curves). evolution of ionization. The ongoing dipole oscillations in case of TDLDA leads to ongoingelectron emission. The case with dissipation shows a leveling off for the ionization. Thedipole signals has been damped away and the excitation energy is converted to intrinsic,thermal energy. This energy later on leads to a thermal electronix emission at a muchslower time scale, thus not visible here. After all, we see that dissipative TDLDA canprovide a pertinent of the thermalization of an excited electron cloud.5.2.3. Stochastic Time-Dependent Hartree-Fock
An alternative route to kinetic theory is provided by stochastic methods describingthe system as an ensemble of (pure) mean-field states. This leads to Stochastic Time-Dependent Hartree-Fock (STDHF), or Stochastic TDLDA (STDLDA) when combinedwith density functionals. It was originally formulated in the context of nuclear colli-sions [234], whence the acronym TDHF, but it can formulated for whatever system inwhich a quantum mean field provides a sound description of the ground state and to lowenergy properties. In the case of clusters and molecules TDLDA provides the obviouseffective mean field theory as a starting point, thus coming to STDLDA. For simplicity,91e will use the notion STDHF further on. STDHF contains all the ingredients of a stan-dard kinetic equation, complemented by proper statistical fluctuations. It accounts forcollisional correlations (from electron-electron collisions). They are treated in incoherentmanner and should not be mixed with coherent correlations as they typically dominatein low-energy processes.The original formulation of STDHF started from the quantum Liouville equation fordensity matrices. The early studies could show that the ensemble description of STDHFcan be reduced to a quantum Boltzmann equation complemented by a to the quantumBoltzmann Langevin equation [234]. The latter was introduced in [235, 236] and hasbeen particularly studied in the nuclear context [237, 238, 225, 239]. STDHF is thus awell founded theory containing all the ingredients necessary for a description of dissipa-tive electronic features. It has unfortunately never been explored at full quantum levelbecause of its complexity. It is only recently that the first calculations were performedin model systems [240] with a proper reformulation of the theory. The first results arequite promising and we will thus discuss here briefly the formalism and typical resultsobtained in a simple system.The STDHF describes the system by an ensemble of N Slater states {| Φ α (cid:105) , α =1 , . . . , N } . Each state | Φ α (cid:105) is associated with a set of single-particle (s.p.) states { ϕ αi , i =1 , . . . , Ω } . The labels i = 1 , . . . , N with N <
Ω stand for the occupied (hole) states.We also include in the description a sufficient amount of unoccupied (particle) states i = N +1 , , . . . , Ω which will serve as a ”reservoir” of levels for transitions to come. Withthis ensemble of s.p. states we can now unfold hierarchy of n -particle- n -hole ( nph ) ex-citations. The first ones to be appear are 2 ph excitations because 1 ph excitations arealready accounted for in the mean-field propagation (TDHF or TDLDA). The correlatedwave function (starting from an uncorrelated situation) can then be expanded as | Ψ α ( t ) (cid:105) = | Φ α ( t ) (cid:105) + (cid:88) pp (cid:48) hh (cid:48) c αpp (cid:48) hh (cid:48) ( t ) | Φ αpp (cid:48) hh (cid:48) ( t ) (cid:105) , (56)with | Φ αpp (cid:48) hh (cid:48) (cid:105) = ˆ a † p ˆ a † p (cid:48) ˆ a h (cid:48) ˆ a h | Φ α (cid:105) . (57)Note that the 2 ph states | Φ αpp (cid:48) hh (cid:48) (cid:105) are also Slater states. Starting from an uncorrelatedsituation, one then propagates a correlated state | Ψ α ( t ) (cid:105) up to a certain time τ at whichit is sampled in terms of an ensemble {| Φ ακ (cid:105) , w ακ } , where κ ∈ { , pp (cid:48) hh (cid:48) } . The weight w ακ = | c αpp (cid:48) hh (cid:48) | is the probability with which | Φ ακ (cid:105) appear. It is evaluated by meansof time-dependent many-body perturbation theory which finally leads to a transitionprobability following Fermi’s golden rule as [234] w αpp (cid:48) hh (cid:48) = τ (cid:12)(cid:12)(cid:12) (cid:104) Φ ακ | ˆ W | Φ α (cid:105) (cid:12)(cid:12)(cid:12) δ ( ε αp + ε αp (cid:48) − ε αh − ε αh (cid:48) ) , (58)where we have introduced the residual interaction ˆ W complementing the mean-fieldhamiltonian ˆ h . The original state | Φ α (cid:105) itself has a weight w α = 1 − (cid:80) w αpp (cid:48) hh (cid:48) (attributing w α to the ”no transition” case) which is the complement of all the other transition prob-abilities. The Dirac δ -function has to be taken with a word of caution. The full expressioninvolves an operator δ function of the mean-field Liouvillean [234]. The approximation(58) involves s.p. energies taken as expectation values over the s.p. states. These, however,are ambiguous to the extent that one has always the freedom of a unitary transforma-tion amongst the occupied states. We define the ε αi uniquely by diagonalizing the actual92ean-field hamiltonian ˆ h α separately amongst occupied states and unoccupied states. Inpractice, the Dirac δ -function has to be augmented by a finite width to account for thediscrete nature of spectra in finite systems [240]. The choice of τ and ˆ W also requiressome caution as the scheme, being based on time-dependent perturbation theory, has toremain in the weak coupling limit as typical of standard kinetic theory. In particular thesampling interval τ should be long enough to allow a sufficient number of ”collisions”to take place, which then justifies stochastic reductions and loss of coherence, but shortenough to remain perturbative, namely with w α (cid:28) | Φ (cid:105) and each member of the ensemble is initially set to | Φ α (0) (cid:105) = | Φ (cid:105) . We then propagate each | Φ α ( t ) (cid:105) individually, first from t = 0 to τ byTDHF/TDLDA. At time τ , all 2 ph states about | Φ α ( τ ) (cid:105) are evaluated as well as the asso-ciated jump probabilities w ακ following Eq.(58). One state | Φ ακ (cid:105) is then randomly selectedaccording to its weight w ακ . | Φ ακ (cid:105) is then again propagated according to TDHF/TDLDAfrom τ to 2 τ up to 2 τ at which a similar sampling takes place, and so on. The aboveprocedure is then restarted from initial time for each member of the ensemble separately.This altogether provides the STDHF ensemble {| Φ α ( t ) (cid:105) , α = 1 , . . . , N } : | Φ α (cid:105) TDHF −→ {| Φ ακ (cid:105) , w ακ } (cid:124) (cid:123)(cid:122) (cid:125) Sampling | Φ ακ (cid:105) TDHF −→ {| Φ ακ (cid:48) (cid:105) , w ακ (cid:48) } (cid:124) (cid:123)(cid:122) (cid:125) ...t = 0 τ τ ... α = ,..., N The ensemble {| Φ α (cid:105) , α = 1 , . . . , N } allows one to compute any observable by standardstatistical averages. In particular, one-body or two-body operators (both correlated) canbe directly constructed from the ensemble. The one-body density matrix readsˆ ρ = 1 N N (cid:88) α =1 ˆ ρ α = N (cid:88) α =1 N (cid:88) i =1 | ϕ αi (cid:105)(cid:104) ϕ αi | ≡ Ω (cid:88) ν =1 | ϕ ( nat ) ν (cid:105) n ν (cid:104) ϕ ( nat ) i | (59)where the second representation employs the natural s.p. orbitals | ϕ ( nat ) ν (cid:105) diagonalizing ˆ ρ and immediately delivers the associated (fractional) occupation numbers n ν . The latterquantities provide a natural tool for analyzing thermal effects.For a first test of STDHF, we use a simple 1D model simulating a dimer molecule. Themean-field hamiltonian reads (in x representation and taking (cid:126) =1) :ˆ h α = − ∆2 m + V ext ( x ) + λ ( (cid:37) α ( x )) . (60)It contains a self consistent term λ ( (cid:37) α ( x )) (with λ = 27 . ) involving the localone-body density (cid:37) α ( x ) (cid:37) ( x ) = (cid:80) N | ϕ αi ( x ) | . This terms stands for the effect of a simpledensity functional. The external potential V ext ( x ) has a Woods-Saxon shape : V ext ( x ) = V / (1 + exp(( x − x ) /a )) with V = −
68 eV, x = 15 a , a = 2 a . It is complemented,outside the well, by a confining harmonic oscillator ensuring soft reflecting boundaryconditions. These boundaries allow one to avoid direct emission and to focus the analysis93n the building up of thermal effects. Altogether, the model mimics a typical situation inclusters and molecules with a fixed external potential delivered by the ions and an energyscale typical of organic systems. The residual interaction is, in the present study, chosenschematically as a simple zero-range force W ( x, x (cid:48) ) = W δ ( x − x (cid:48) ) with W = 40 . τ = 1 fs (0.5 and 1.5 fsgive similar results) and the dynamics is followed over 100 fs, which is much larger thanthe optical period (1.15 fs) and long enough to study thermal relaxation. We propagatean ensemble of N =100 events.Fig. 57 shows the time evolution of s.p. energies for a typical STDHF event (bottom)and compares it to the corresponding pure TDHF evolution (top). Each member of the -60-50-40-30-20-100 s . p . e n e r g y ( e V ) ph excitation – E ∗ = 25 . s . p . e n e r g y ( e V ) Time (fs)STDHF
Fig. 57. Time evolution of single particle energies in a TDHF calculation (upper panel) and in oneSTDHF event (lower panel) for an initial excitation energy of 25.8 eV. Red lines correspond to occupiedstates, green ones to unoccupied ones (see text for details). In the STDHF case, five 2 ph transitionsactually occurred (were actually sampled) for this event, transitions which are indicated by faint verticaldashed lines. ensemble will actually deliver a different sequence of transitions which will finally lead tothe mixed state representing the correlated system. One can identify five 2 ph transitions(indicated by faint dashed lines). It is also interesting to note that, up to minor mean fieldrearrangements, the TDHF evolution does preserve the arrangement of particle and hole94tates in the course of time evolution. Pure TDHF evolution more or less preserves initialoccupations in time, hindering possible relaxation to a thermal state. STDHF overcomesovercomes this limitation as can be seen from the rearrangements in the lower panel.Fig. 58 displays snapshots of occupation numbers, extracted from the one body den-sity matrix, see Eq.(59), at several times along a STDHF propagation. The occupation O cc up a t i o nnu m b e r s Energy (eV)1 ph excitation – E ∗ = 25 . t = 5 fs t = 25 fs t = 50 fs t = 100 fs Fig. 58. Snapshots of occupation numbers as a function of time in an example of STDHF calculation,with initial excitation E ∗ of 25.8 eV and for a 100 event ensemble. numbers keep for some time a trace of the original excitation, in particular the initialhole around −
37 eV (see red curve). This hole is gradually filled and the occupationnumbers are soon washed out, leading asymptotically to an energy profile typical of ther-mal equilibrium. Note some unavoidable statistical fluctuations, still visible at very lowenergy. This figure therefore demonstrates the capability of STDHF to account for relax-ation effects at the side of electrons in a purely quantum mechanical manner. One shouldalso stress that STDHF enables to estimate fluctuations around average values as com-puted from one or two-body density matrix, which again represents a remarkable stepforward. The major obstacle of STDHF lies in the high cost of handling large ensembleswhich becomes particularly demanding for low excitation energies where smaller transi-tion probabilities require better statistics. This still more development work is needed toput STDHF fully into action.
6. Conclusions
In this paper we have reviewed the analysis of electron emission following irradiationof clusters and molecules by light pulses. Observables from electron emission give de-tailed insight into the dynamical response of the irradiated species. Understanding theirradiation and emission process is also essential in view of the many applications in ma-terials science, biology, and medicine. High-resolution studies of electron emission havemade tremendous progress over the past few years, both experimentally and theoreti-cally. In experiments, new developments in light sources now provide a broad choice ofelectromagnetic pulses with widely variable frequency, intensity, and time profile down to95ttosecond resolution in the range of electronic time scales. There is also great progress atthe side of the measurement giving access to increasingly detailed properties of emittedelectrons, high-resolution photo-electron spectra (PES) or angular distributions (PAD),often combined to velocity map imaging. As latest achievements, time-resolved PES/PADare waiting in the wings. This remarkable experimental progress calls for elaborate theo-retical treatments at the most microscopic level of description. In this respect, TDDFT,especially when solved in real time, constitutes an invaluable tool to simulate the variousdynamical scenarios of irradiation of clusters and molecules. Therefore, many groups allover the world are heavily working on such approaches. In this review, we have presenteda series of developments and results, mostly from the last decade, on irradiation of clus-ters and molecules by light pulses and subsequent detailed analysis of electron emission.Concerning the observables from emission, we consider total ionization, PES and PAD,also in connection with time-resolved measurements. At the theory side, we focus on amicroscopic description in terms of time-dependent density functional theory (TDDFT).This is practically handled at the level of the time-dependent local density approxima-tion (TDLDA) augmented by a self-interaction correction (SIC). We have illustrated thecapabilities of our approach on several systems, ranging from simple molecules like N tofashionable nano-clusters such as C , and also studying archetypal metal systems suchas Na clusters. We list below a few results that we consider as being emblematic of thesestudies.The theoretical modeling is based on TDDFT which is known to provide a robustmicroscopic description of the system dynamics. It enables to include ionic motion ata classical level. The fully coupled dynamics is needed in cases of long laser pulses andfor thermal ensembles. In most cases, we consider short pulses and keep ions frozen.Starting level for TDDFT is the time-dependent local-density approximation (TDLDA).However, in order to describe ionization dynamics properly, one needs a theory fulfillingKoopmans’ theorem which states that the ionization potential (IP) has to be identicalwith the single particle (s.p.) energy of the least bound state. This is violated by theself-interaction error in LDA. It can be cured by augmenting LDA with a SIC. Thelatter has been turned manageable by a handling in terms of two sets of occupied s.p.states, the 2setSIC scheme. A less expensive alternative is offered by averaged densitySIC (ADSIC) which performs surprisingly well as long as the dynamics stays off theregime of fragmentation and/or huge ionization. ADSIC has thus been used here in mostcases and it allowed us to obtain remarkably accurate results in good agreement withexperiments.The numerical handling of ionization dynamics is most efficient in a coordinate-spacerepresentation. To describe ionization, we augment the coordinate-space grid by absorb-ing boundary conditions. They allow one to trivially compute the observables of totalionization and ionization out of each s.p. state separately (level depletion). PAD arecomputed by collecting the electron loss in angular segments on the grid. A compli-cation arises when comparing PAD with measurements, since clusters or molecules ingas phase have an undefined orientation. They represent, in fact, an isotropic ensembleof orientations. Theoretical calculations need thus to be complemented by orientationaveraging. For the one-photon domain, we have worked out a compact formula whichcan live with only six reference orientations to be computed. Multi-photon processes re-quire direct integration where we find that one can obtain reliable results with typically18–36 integration points, depending on the symmetry of the cluster. PES are computed96y recording the phase oscillations of outgoing wave functions close to the onset of theabsorbing bounds, and by finally Fourier transforming the temporal oscillations to theenergy domain. Special care has to be taken to cope with strong laser pulses. They maymodify the phase of the wave functions at the sampling point. Fortunately, this effect canbe evaluated analytically and allows us to derive a phase correction, thus rendering thescheme for computing PES reliable up to rather large laser intensities (typically 10 –10 W/cm . Altogether, we have thus at hand powerful and versatile tools to simulateionization dynamics and to evaluate the observables deduced thereof. In the following,we will briefly summarize the results for each observable separately.The simple signal of total ionization is already useful when combined with system-atics. For example, the frequency dependence of ionization maps the underlying dipoleresponse. Ionization becomes the key signal in pump and probe (P&P) scenarios whichconstitute a well established tool for a time-resolved measurement of ionic motion. Clus-ters are rather complex systems where the motion of single ions is hard to track. P&Pmeasurements at least enable to identify global properties of the ionic configuration, asradius and quadrupole deformation, which is already very useful information in studies ofCoulomb explosion of clusters. On a first example, we could show how the now upcomingattosecond pulses allow time-resolved analysis at electronic pace.The main body of the paper dealt with results on PES and PAD. Combined PES/PADas obtained from velocity map imaging (VMI) contain a very rich amount of information,but are usually hard to appreciate as such. The energy- or angular-integrated versionsthereof, delivering PAD and PES, are better suited for detailed analysis and comparisonsbetween experiments and theory. A PES in a strictly one-photon domain delivers a mapof the clusters s.p. spectrum. In experiments, this was limited previously to negativelycharge cluster anions due to a low ionization potential. The availability of coherent high-frequency light sources now allows one to employ the one-photon analysis for neutralclusters and even cations. Multi-photon processes make PES more involved and richer.In the low-intensity regime, one can identify multiple copies of the s.p. spectra separatedby the photon frequency. But PES goes beyond just mapping s.p. spectra. It deliversa picture of the whole dynamical processes. We have illustrated that by working outthe impact of plasmon resonances which can directly drop its signatures in the PESthemselves. Increasing intensity produces more ionization which Coulomb shifts the s.p.states gradually downwards, thus broadening the peaks in the PES. This eventually leadsto totally smoothed PES with straightforward exponential decrease. This kind of patternsuggests at first glance an interpretation as purely thermal electron emission. A closerinspection from energetic considerations and TDLDA simulations reveals that it cannotbe fully a thermal process. We probably encounter a mixed situation. Direct emissionstill prevails and electronic re-collisions add first thermal effects to the picture. We havealso seen that PES alone cannot distinguish unambiguously between direct and thermalemission.The unavoidable orientation averaging of PAD wipes out many details which werecontained in before averaging. Fortunately, there remains a lot of useful information.Orientation averaged PAD in the one-photon regime can be characterized by one singleparameter, the anisotropy β . A systematic survey of β and its frequency dependencerevealed that PAD are extremely sensitive to every detail of the modeling. This holdseven more so for state-resolved anisotropies β ( i )2 (where i stands for a s.p. state or de-97enerated group thereof). No compromises are allowed. One must invoke the full ma-chinery of TDLDA+SIC and a careful description of the detailed ionic structure to havea chance for a relevant description. For example, the β ( ω las ) computed with smoothjellium background shows marked fluctuations which disappear when ionic backgroundis used. However, there remains one remarkable exception in the low-frequency tail of β ( ω las ) for the loosely bound anion Na − . This shows a deep dip towards one-photonemission threshold and it does so independent of the model for the ionic background. Inthe multi-photon domain, a PAD provides crucial information which helps to distinguishdirect from thermal emission. This was nicely visible in the combined PES/PAD of C where one could associate uniquely the region of low kinetic energy with thermal elec-trons while higher kinetic energies show clear sign of direct emission with the PAD beingforward/backward dominated.The collection of results presented in this review has several open ends which call forfurther development and investigation. We briefly quote a few of them which we considerto be important next steps. It was already mentioned above that it becomes increasinglypossible to extend time-resolved analysis to the attosecond domain. Theoretical stud-ies are required to explore the huge space of new possibilities and to find out the mostpromising experimental conditions. A further interesting perspective emerges if combin-ing time-resolved analysis with PES and PAD. This enables, e.g., to track the branchingbetween direct and thermal emission in the course of time. First studies in this directionare very promising. The discussion of PES and PAD left as a yet unsolved problem thedistinction between direct and thermal processes in electron emission. This calls for aproper theoretical modeling of electron-electron collisions (dynamical electron correla-tions) and subsequent dissipation effects. Two promising development lines for such anextended TDLDA have been presented, that is the relaxation time approach which mod-els dissipation phenomenologically and a stochastic mean field model which describes thesystem as an ensemble of mean-field states incorporating the electron-electron collisionsas stochastic jumps between these states.The results summarized above prove that observables from electron emission are anextremely powerful tool to analyze irradiation processes on clusters and molecules. Theyallowed one in the past to reveal many interesting aspects of structure and dynamics ofthese systems. The future directions sketched in the previous paragraph show that we arenot nearly at an end of the investigations. The field remains lively and highly interestingand challenging for the future. Acknowledgments
The authors gratefully thank K. Andrae, C. Bordas, F. Calvayrac, J.-M. Escartin, B.and M. Farizon, B. Faber, T. Fennel, C. Z. Gao, M. Ivanov, P. Kl¨upfel, F. L´epine, K.-H.Meiwes-Broer, J. Messud, A. Pohl, P. Romaniello, J.-M. Rost, N. Slama, O. Smirnova,J. Tiggesba¨umker, M. Vincendon, B. von Issendorff, Z. P. Wang, and F.-S. Zhang, forfruitful discussions during the realization of this work. The authors would like to ac-knowledge financial support from ANR-10-BLAN-0428, ANR-10-BLAN-0411, ANR-11-IS04-0003, ITN-CORINF, and the Institut Universitaire de France. The theoretical workwas granted access to the HPC resources of IDRIS under the allocation 2013–095115made by GENCI (Grand Equipement National de Calcul Intensif), of CalMiP (Calcul en98idi-Pyr´en´ees) under the allocation P1238, and of RRZE (Regionales RechenzentrumErlangen).
References [1] U. Keller, Nature 424 (2003) 831.[2] C. Rulli`ere (Ed.), Femtosecond Laser Pulses: Principles and Experiments, 2nd ed., Advanced Textsin Physics, Springer, New-York, 2005.[3] R. Paschotta, Encyclopedia of Laser Physics and Technology, volumes 1 and 2, Wiley-VCH, Berlin,2008.[4] N. Ishikawa, Y. Saitoh, T. Yamaki, F. Hori, M. Sasase, Nucl. Inst. & Meth. B 314 (2013) 1.[5] G. Garcia, M. C. Fuchs, Series in Biological and Medical Physics, Biomedical Engineering, Springer-Verlag, Berlin, 2012.[6] Y. Zhang, Y. Wang, W. J. Weber, Nucl. Inst. & Meth. B, vol 286, 2012.[7] B. Liu, S. B. Nielsen, P. Hvelplund, H. Zettergren, H. Cederquist, B. Manil, B. A. Huber, Phys.Rev. Lett. 97 (2006) 133401.[8] G. Mie, Ann. Phys. (Leipzig) 25 (1908) 377.[9] U. Kreibig, M. Vollmer, Optical Properties of Metal Clusters, Vol. 25, Springer Series in MaterialsScience, 1993.[10] M. Brack, Rev. Mod. Phys. 65 (1993) 677.[11] W. A. de Heer, Rev. Mod. Phys. 65 (1993) 611.[12] P.-G. Reinhard, E. Suraud, Introduction to Cluster Dynamics, Wiley, New York, 2003.[13] U. Saalmann, C. Siedschlag, J. M. Rost, J. Phys. B 39.[14] F. Calvayrac, P.-G. Reinhard, E. Suraud, C. A. Ullrich, Phys. Rep. 337 (2000) 493.[15] T. Fennel, K.-H. Meiwes-Broer, J. Tiggesb¨aumker, P. M. Dinh, P.-G. Reinhard, E. Suraud, Rev.Mod. Phys. 82 (2010) 1793.[16] E. E. B. Campbell, K. Hansen, K. Hoffmann, G. Korn, M. Tchaplyguine, M. Wittmann, I. V.Hertel, Phys. Rev. Lett. 84 (2000) 2128.[17] C. Bartels, C. Hok, J. Huwer, R. Kuhnen, J. Schw¨obel, B. von Issendorff, Science 323 (2009) 132333.[18] B. Bouda¨ıffa, P. Cloutier, D. Hunting, M.-A. Huels, L. Sanche, Science 287 (2000) 1658.[19] I. Seideman, Ann. Rev. Phys. Chem. 53 (2002) 41.[20] M. Rohmer, M. Bauer, T. Leissner, C. Schneider, A. Fischer, G. Niedner-Schatteburg, B. vonIssendorff, M. Aeschlimann, Phys. Status Solidi B 247 (2010) 1132.[21] M. A. L. Marques, N. T. Maitra, F. M. S. Nogueira, E. K. U. Gross, A. Rubio, Lect. Notes in Phys.vol 837, Springer-Verlag, Berlin, 2012.[22] T. Brabec, F. Krausz, Rev. Mod. Phys. 72 (2000) 545.[23] F. Krausz, M. Ivanov, Rev. Mod. Phys. 81 (2009) 163.[24] J. Feldhaus, J. Arthur, J. B. Hastings, J. Phys. B 38 (2005) 799.[25] L. V. Keldysh, Sov. Phys. JETP 20 (1964) 1307.[26] S. Augst, D. Strickland, D. D. Meyerhofer, S. L. Chin, J. H. Eberly, Phys. Rev. Lett. 63 (1989)2212.[27] T. Seidemann, M. Y. Ivanov, P. B. Corkum, Phys. Rev. Lett. 75 (1995) 2819.[28] T. Zuo, A. D. Bandrauk, Phys. Rev. A 52 (1995) R2511.[29] V. V´eniard, R. Ta¨ıeb, A. Maquet, Phys. Rev. A 65 (1) (2001) 013202.[30] C. Siedschlag, J. M. Rost, Phys. Rev. Lett. 89 (2002) 173401.[31] I. Last, J. Jortner, Phys. Rev. A 60 (1999) 2215.[32] H. Wabnitz, L. Bittner, A. R. B. de Castro, R. Dohrmann, P. Gurtler, T. Laarmann, W. Laasch,J. Schulz, A. Swiderski, K. von Haeften, T. Moller, B. Faatz, A. Fateev, J. Feldhaus, C. Gerth,U. Hahn, E. Saldin, E. Schneidmiller, K. Sytchev, K. Tiedtke, R. Treusch, M. Yurkov, Nature 420(2002) 482.[33] L. K¨oller, M. Schumacher, J. K¨ohn, S. Teuber, J. Tiggesb¨aumker, K.-H. Meiwes-Broer, Phys. Rev.Lett. 82 (1999) 3783.[34] L. Schweikhard, A. Herlert, S. Kr¨uckeberg, M. Vogel, C. Walther, Phil. Mag. B 79 (1999) 1343.[35] T. D¨oppner, S. Teuber, M. Schumacher, J. Tiggesb¨aumker, K. Meiwes-Broer, Appl. Phys. B 71(2000) 357.
36] E. Suraud, P.-G. Reinhard, Phys. Rev. Lett. 85 (2000) 2296.[37] U. Saalmann, J. M. Rost, Phys. Rev. Lett. 91 (2003) 223401.[38] T. D¨oppner, T. Fennel, T. Diederich, J. Tiggesb¨aumker, K. Meiwes-Broer, Phys. Rev. Lett. 94(2005) 013401.[39] T. D¨oppner, T. Fennel, P. Radcliffe, J. Tiggesb¨aumker, K.-H. Meiwes-Broer, Phys. Rev. A 73 (2006)031202.[40] C. Siedschlag, J. M. Rost, Phys. Rev. A 71 (2005) 031401(R).[41] T. Bornath, P. Hilse, M. Schlanges, Laser Phys. 17 (2007) 591.[42] S. Zamith, T. Martchenko, Y. Ni, S. A. Aseyev, H. G. Muller, M. J. J. Vrakking, Phys. Rev. A 70(2004) 011201.[43] Y. L. Shao, T. Ditmire, J. W. G. Tisch, E. Springate, J. P. Marangos, M. H. R. Hutchinson, Phys.Rev. Lett. 77 (1996) 3343.[44] E. Springate, S. A. Aseyev, S. Zamith, M. J. J. Vrakking, Phys. Rev. A 68 (2003) 053201.[45] C. Bordas, F. Paulig, H. Heln, D. L. Huestis, Rev. Sci. Instr. 67 (1996) 2257.[46] O. Kostko, C. Bartels, J. Schwobel, C. Hock, B. v Issendorff, J. Phys. : Conf. Ser. 88 (2007) 012034.[47] M. Kjellberg, O. Johansson, F. Jonsson, A. V. Bulgakov, C. Bordas, E. E. B. Campbell, K. Hansen,Phys. Rev. A 81 (2010) 023202.[48] O. Cheshnovsky, K. J. Taylor, J. Conceicao, R. E. Smalley, Phys. Rev. Lett. 64 (1990) 1785.[49] G. Wrigge, M. A. Hoffmann, B. von Issendorff, Phys. Rev. A 65 (2002) 063201.[50] A. Pohl, P.-G. Reinhard, E. Suraud, J. Phys. B 37 (2004) 3301.[51] A. Rytk¨onen, H. H¨akkinen, M. Manninen, Phys. Rev. Lett. 80 (1998) 3940.[52] K. Hansen, K. Hoffmann, E. E. B. Campbell, J. Chem. Phys. 119 (2003) 2513.[53] R. M. Young, M. A. Yandell, S. B. King, D. M. Neumark, J. Chem. Phys 136 (2012) 094304.[54] J. Wills, F. Pagliarulo, B. Baguenard, F. L´epine, C. Bordas, Chem. Phys. Lett. 390 (2004) 145.[55] Y. Yamaguchi, Y. Osamure, J. D. Goddard, H. Schaefer, A New Dimension to Quantum Chemistry:Analytic Derivative Methods in Ab Initio Molecular Electronic Structure Theory, Oxford UniversityPress, Oxford, 1994.[56] P. Krause, T. Klamroth, P. Saalfrank, J. Chem. Phys. 127.[57] A. A. Vlasov, Many Particle Theory and Its Applications to Plasma, Gordon & Breach, New York,1950.[58] G. F. Bertsch, S. Das Gupta, Phys. Rep. 160 (1988) 190.[59] A. Domps, E. Giglio, P.-G. Reinhard, E. Suraud, J. Phys. B 33 (2000) L333.[60] T. Fennel, G. F. Bertsch, K.-H. Meiwes-Broer, Eur. Phys. J. D 29 (2004) 367.[61] C. Sieber, W. Harbich, K.-H. Meiwes-Broer, C. F´elix, Chem. Phys. Lett. 433 (2006) 32.[62] L. Szasz, Pseudopotential Theory of Atoms and Molecules, Wiley, New York, 1985.[63] S. K¨ummel, M. Brack, P.-G. Reinhard, Eur. Phys. J. D 9 (1999) 149.[64] S. Goedecker, M. Teter, J. Hutter, Phys. Rev. B 54 (1996) 1703.[65] N. W. Ashcroft, N. D. Mermin, Solid State Physics, Saunders College, Philadelphia, 1976.[66] W. Ekardt, Z. Penzar, Phys. Rev. B 43 (1991) 1331.[67] B. Montag, T. Hirschmann, J. Meyer, P.-G. Reinhard, M. Brack, Phys. Rev. B 52 (1995) 4775.[68] P. Hohenberg, W. Kohn, Phys. Rev. 136 (1964) 864.[69] W. Kohn, L. J. Sham, Phys. Rev. 140 (1965) 1133.[70] R. G. Parr, W. Yang, Density-Functional Theory of Atoms and Molecules, Oxford University Press,Oxford, 1989.[71] R. M. Dreizler, E. K. U. Gross, Density Functional Theory: An Approach to the Quantum Many-Body Problem, Springer-Verlag, Berlin, 1990.[72] E. K. U. Gross, J. F. Dobson, M. Petersilka, Top. Curr. Chem. 181 (1996) 81.[73] A. D. Becke, Phys. Rev. A 38 (1988) 3098.[74] J. P. Perdew, K. Burke, M. Ernzerhof, Phys. Rev. Lett. 77 (1996) 3865.[75] J. P. Perdew, S. Kurth, A. Zupan, P. Blaha, Phys. Rev. Lett. 82 (1999) 2544.[76] A. K. Dhara, S. K. Ghosh, Phys. Rev. A 35 (1987) 442(R).[77] G. Vignale, W. Kohn, Phys. Rev. Lett. 77 (1996) 2037.[78] C. Ullrich, Oxford University Press, Oxford, 2012.[79] R. D’Agosta, G. Vignale, Phys. Rev. Lett. 96 (2006) 016405.[80] N. T. Maitra, I. Souza, K. Burke, Phys. Rev. B 68 (2003) 045109.[81] P. Bonche, B. Grammaticos, S. E. Koonin, Phys. Rev. C17 (1978) 1700.
82] J. Terasaki, P.-H. Heenen, H. Flocard, P. Bonche, Nucl. Phys. A 600 (1996) 371–386.[83] A. Castro, H. Appel, M. Oliveira, C. Rozzi, X. Andrade, F. Lorenzen, M. Marques, E. Gross,A. Rubio, Phys. Stat. Sol. B 243 (2006) 2465.[84] X. Andrade, J. Alberdi-Rodriguez, D. A. Strubbe, M. J. T. Oliveira, F. Nogueira, A. Castro,J. Muguerza, A. Arruabarrena, S. G. Louie, A. Aspuru-Guzik, A. Rubio, M. A. L. Marques, J.Phys.: Cond. Matt. 24 (2012) 233202.[85] R. Y. Cusson, R. K. Smith, J. A. Maruhn, Phys. Rev. Lett. 36 (1976) 1166–1169.[86] F. Calvayrac, P.-G. Reinhard, E. Suraud, Phys. Rev. B 52 (1995) R17056.[87] J. W. Eastwood, D. R. K. Brownrigg, J. Comp. Phys. 32 (1979) 24.[88] G. Lauritsch, P.-G. Reinhard, Int. J. Mod. Phys. C 5 (1994) 65.[89] P.-G. Reinhard, R. Y. Cusson, Nucl. Phys. A 378 (1982) 418.[90] M. D. Feit, J. A. Fleck, A. Steiger, J. Comp. Phys. 47 (1982) 412.[91] V. Blum, G. Lauritsch, J. A. Maruhn, P.-G. Reinhard, J. Comp. Phys 100 (1992) 364.[92] K. T. R. Davies, S. E. Koonin, Phys. Rev. C23 (1981) 2042–2061.[93] B. Montag, P.-G. Reinhard, Phys. Lett. A 193 (1994) 380.[94] T. Berkus, P.-G. Reinhard, E. Suraud, Int. J. Mol. Sci. 3 (2002) 69.[95] P. M. Dinh, F. Fehrer, P.-G. Reinhard, E. Suraud, Eur. Phys. J. D 45 (2007) 415.[96] P. M. Dinh, P.-G. Reinhard, E. Suraud, Phys. Rep. 485 (2010) 43.[97] R. S. Varga, Iterative Matrix Analysis, Prentice-Hall, Englewood Cliffs, 1962.[98] W. H. Press, S. A. Teukolsky, W. T. Vetterling, B. P. Flannery, Numerical Recipes, CambridgeUniversity Press, Cambridge, 1992.[99] J. P. Perdew, A. Zunger, Phys. Rev. B 23 (1981) 5048.[100] J. Messud, P. M. Dinh, P.-G. Reinhard, E. Suraud, Ann. Phys. 324 (2009) 955.[101] M. R. Pederson, R. A. Heaton, C. C. Lin, J. Chem. Phys. 80 (5) (1984) 1972.[102] J. Messud, P. M. Dinh, P.-G. Reinhard, E. Suraud, Phys. Rev. Lett. 101 (2008) 096404.[103] S. K¨ummel, L. Kronik, Rev. Mod. Phys. 80 (2008) 3.[104] J. B. Krieger, Y. Li, G. J. Iafrate, Phys. Rev. A 45 (1992) 101.[105] R. T. Sharp, G. K. Horton, Phys. Rev. 30 (1953) 317.[106] C. Legrand, E. Suraud, P.-G. Reinhard, J. Phys. B 35 (2002) 1115.[107] I. Ciofini, C. Adamo, H. Chermette, Chem. Phys. 309 (2005) 67–76.[108] P. Kl¨upfel, P. M. Dinh, P.-G. Reinhard, E. Suraud, Phys. Rev. A 88 (2013) 052501.[109] E. Fermi, E. Amaldi, Accad. Ital. Rome 6 (1934) 117.[110] Y. H. Yu, T. Zuo, A. D. Bandrauk, J. Phys. B 31 (1998) 1533.[111] U. De Giovannini, D. Varsano, M. A. L. Marques, H. Appel, E. K. U. Gross, A. Rubio, Phys. Rev.A 85 (2012) 062515.[112] K. Boucke, H. Schmitz, H.-J. Kull, Phys. Rev. A 56 (1997) 763.[113] M. Mangin-Brinet, J. Carbonell, C. Gignoux, Phys. Rev. A 57 (1998) 3245.[114] T. Nakatsukasa, K. Yabana, Phys. Rev. C 71 (2005) 024301.[115] J. L. Krause, K. J. Schafer, K. C. Kulander, Phys. Rev. A 45 (1992) 4998.[116] P.-G. Reinhard, E. Suraud, Cluster dynamics in strong laser fields, in: M. A. L. Marques, C. A.Ullrich, F. Nogueira (Eds.), Time-dependent density functional theory, Vol. 706 of Lecture Notesin Physics, Springer, Berlin, 2006, p. 391.[117] P.-G. Reinhard, P. D. Stevenson, D. Almehed, J. A. Maruhn, M. R. Strayer, Phys. Rev. E 73 (2006)036709.[118] A. Pohl, P.-G. Reinhard, E. Suraud, Phys. Rev. A 70 (2004) 023202.[119] A. Pohl, P.-G. Reinhard, E. Suraud, Phys. Rev. Lett. 84 (2000) 5090.[120] A. Pohl, Ph.D. thesis, Friedrich-Alexander-Universit¨at, Erlangen/N¨urnberg (2003).[121] P. M. Dinh, P. Romaniello, P.-G. Reinhard, E. Suraud, Phys. Rev. A 87 (2013) 032514.[122] A. Gazibegovic-Busuladzic, E. Hasovic, M. Busuladzic, D. B. Milosevic, F. Kelkensberg, W. K.Siu, M. J. J. Vrakking, F. L´epine, G. Sansone, M. Nisoli, I. Znakovskaya, M. Kling, Phys. Rev. A84 (2011) 043426.[123] D. B. Miloˇsevi´c, G. G. Paulus, D. Bauer, W. Becker, J. Phys. B 39 (14) (2006) R203.[124] H. Stapelfeldt, T. Seideman, Rev. Mod. Phys. 75 (2003) 543.[125] S. Pabst, R. Santra, Phys. Rev. A 81 (2010) 065401.[126] S. Pabst, R. Santra, Phys. Rev. A 82 (2010) 044901(E).[127] P.-G. Reinhard, E. Suraud, Eur. Phys. J. D 34 (2005) 145.206] L. Schweikhard, S. Kr¨uckeberg, K. L¨utzenkirchen, C. Walther, Eur. Phys. J. D 9 (1999) 15.[207] H. Haberland, T. Hippler, J. Donges, O. Kostko, M. Schmidt, B. von Issendorff, Phys. Rev. Lett.94 (2005) 035701.[208] R. Schlipper, R. Kusche, B. von Issendorff, H. Haberland, Appl. Phys. A 72 (2001) 255.[209] E. Giglio, P.-G. Reinhard, E. Suraud, Ann. Phys. (Leipzig) 11 (2002) 291.[210] E. Giglio, P.-G. Reinhard, E. Suraud, Phys. Rev. A 67 (2003) 43202.[211] E. Giglio, P.-G. Reinhard, E. Suraud, J. Phys. B 34 (2001) 1253.[212] P.-G. Reinhard, E. Suraud, Eur. Phys. J. D 3 (1998) 175.[213] J. A. R. Samson, R. B. Cairns, J. Opt. Soc. Am. 55 (1965) 1035.[214] H. D. Cohen, U. Fano, Phys. Rev. 150 (1966) 30–33.[215] F. Fr´emont, A. Hajaji, A. Naja, C. Leclercq, J. Soret, J. A. Tanis, B. Sulik, J.-Y. Chesnel, Phys.Rev. A 72 (2005) 050704.[216] S. Chatterjee, D. Misra, A. H. Kelkar, L. C. Tribedi, C. R. Stia, O. A. Foj´on, R. D. Rivarola, Phys.Rev. A 78 (2008) 052701.[217] S. Chatterjee, A. N. Agnihotri, C. R. Stia, O. A. Foj´on, R. D. Rivarola, L. C. Tribedi, Phys. Rev.A 82 (2010) 052709.[218] J. L. Baran, S. Das, F. J´arai-Szab´o, K. P´ora, L. Nagy, J. A. Tanis, Phys. Rev. A 78 (2008) 012710.[219] S. Nandi, A. N. Agnihotri, S. Kasthurirangan, A. Kumar, C. A. Tachino, R. D. Rivarola, F. Mart´ı n,L. C. Tribedi, Phys. Rev. A 85 (2012) 062705.[220] M. Winkworth, P. Fainstein, M. Galassi, J. Baran, B. Dassanayake, S. Das, A. Kayani, J. Tanis,Nucl. Instr. Meth. Phys. Res. B 267 (2) (2009) 373, proceedings of the Fourth InternationalConference on Elementary Processes in Atomic Systems.[221] A. N. Agnihotri, S. Nandi, S. Kasthurirangan, A. Kumar, M. E. Galassi, R. D. Rivarola,C. Champion, L. C. Tribedi, Phys. Rev. A 87 (2013) 032716.[222] S. Nandi, S. Biswas, A. Khan, J. M. Monti, C. A. Tachino, R. D. Rivarola, D. Misra, L. C. Tribedi,Phys. Rev. A 87 (2013) 052710.[223] M. B¨ar and B. Jakob and P.–G. Reinhard and C. Toepffer, Phys. Rev. A 73 (2006) 022719.[224] K. Goeke, P.-G. Reinhard, H. Reinhardt, Ann. Phys. 166 (1986) 257.[225] Y. Abe, S. Ayik, P.-G. Reinhard, E. Suraud, Phys. Rep. 275 (1996) 49.[226] C. Wong, K. Davies, Phys. Rev. C 28 (1983) 2321.[227] R. Cusson, P.-G. Reinhard, J. Maruhn, W. Greiner, M. Strayer, Z. Phys. A 320 (1985) 475.[228] A. Umar, V. Oberacker, Phys. Rev. C 74 (2006) 021601R.[229] G. Bertsch, Z. Phys. A 289 (1978) 103.[230] P. Danielewicz, Phys. Lett. B 146 (1984) 168.[231] P. Danielewicz, Ann. Phys. (N.Y.) 152 (1984) 239.[232] P. Danielewicz, Ann. Phys. (N.Y.) 152 (1984) 305.[233] J. Lehmann, M. Merschdorf, W. Pfeiffer, A. Thon, S. Voll, G. Gerber, Phys. Rev. Lett. 85 (2000)2921.[234] P.-G. Reinhard, E. Suraud, Ann. Phys. (NY) 216 (1992) 98.[235] M. Bixon, R. Zwanzig, Phys. Rev. 187 (1969) 267.[236] R. Zwanzig, J. Stat. Phys. 9 (1973) 215.[237] S. Ayik, C. Gregoire, Phys. Lett. B 212 (1988) 269.[238] J. Randrup, G. F. Burgio, P. Chomaz, Nucl. Phys. A 538 (1992) 393.[239] P. Napolitani, M. Colonna, Phys. Lett. B 726 (2013) 382.[240] N. Slama, P. M. Dinh, P. G. Reinhard, E. Suraud, J. Phys. Conf. Ser., in press.[241] L. van Hove, Physica 21 (1955) 517.