Semi-classical Electron Dynamics in Metal Clusters beyond Mean-Field

Annals of Physics 280, 211235 (2000) doi:10.1006aphy.1999.5990, available online at http:www.idealibrary.com on

Semi-classical Electron Dynamics in Metal Clusters beyond Mean-Field A. Domps Laboratoire de Physique Quantique, Universite Paul Sabatier, 118 Route de Narbonne, F-31062 Toulouse Cedex, France

P.-G. Reinhard Institut fu r Theoretische Physik, Universitat Erlangen, Staudtstrasse 7, D-91058 Erlangen, Germany

and E. Suraud 1 Laboratoire de Physique Quantique, Universite Paul Sabatier, 118 Route de Narbonne, F-31062 Toulouse Cedex, France E-mail: suraudirsamc2.ups-tlse.fr Received March 1, 1999; revised June 18, 1999

We discuss the importance of dynamical correlations in the electronic response of metal clusters to violent excitations. The theoretical description is based on the Vlasov equation, the semi-classical variant of the time dependent local density approximation (TDLDA) to density functional theory (DFT). Dynamical correlations are added by a quantum-kinetic collision term leading to the VlasovUehlingUhlenbeck (VUU) equation. As test cases we consider excitations of sodium clusters either by a rapid collision with a highly charged projectile or by intense femtosecond laser pulses. We find that the collective response and direct ionization (following immediately the excitation process) are insensitive to dynamical correlations for the case of short excitation processes (taking at most a few fs) where two-body collisions have no time to become active. For longer excitations, however, the VUU collision term plays an important role. Whatever the excitation time, dynamical correlations, as described in VUU, are always playing a crucial role in the energy deposition mechanism and subsequent heating of the electron cloud. Such a phenomenon is thus likely to have strong impact on ionic motion inside the cluster on even longer time scales.  2000 Academic Press

1. INTRODUCTION Dynamical properties of metal clusters have been much studied over the past two decades. In the first stage, spectroscopic experiments and slow excitations leading 1

Membre de l'Institut Universitaire de France.

211 0003-491600 35.00 Copyright  2000 by Academic Press All rights of reproduction in any form reserved.

212

DOMPS, REINHARD, AND SURAUD

eventually to Coulomb explosion or fragmentation had been in the foreground [13]. Increased understanding and new experimental devices have raised interest in a non-linear regime where the valence electrons are violently excited within very short time scales, thus placing the cluster in strongly off-equilibrium situations. Processes of that kind are, e.g., the bombardment of clusters with highly charged projectiles [4] or irradiation with intense femtosecond laser pulses [5, 6]. In any regime, the cluster response is dominated in the first stage by excitations of sizable collective oscillations of electrons with respect to ions (Mie plasmon). The subsequent relaxation of the plasmon energy towards other degrees of freedom is a complicated many-body problem. According to the generally accepted picture, the plasmon decays by transferring first its energy to single electrons. These electrons may then be emitted or exchange energy with the remaining electron cloud, therefore leading to a gradual thermalization of the electron distribution, usually at times below 100 fs. At larger times, energy starts being transferred to the ions, possibly leading much later to the evaporation of monomers or dimers or fragmentation of the cluster. The details of electronic and ionic thermalization in clusters are not yet fully understood, neither experimentally nor theoretically. It is the aim of this paper to investigate from a theoretical point of view the electronic relaxation processes in connection with electronelectron collisions. Thanks to pump-probe femtosecond laser experiments, such electronic relaxation processes have already been studied in bulk metals, leading to a better knowledge of the time scales involved [7, 8] for the relaxation of hot electrons [9]. Similar questions have also been much discussed during the 1980s in heavy-ion physics: after being excited by an heavy-ion collision, the energy of collective nuclear modes is redistributed among the nucleons, leading to the formation of considerably hot nuclei [10, 11]. In both cases, it is well known that a proper theoretical description of the heating process needs to go beyond the mean-field level because two-body dissipation, from electronelectron or nucleonnucleon scattering, respectively, plays a decisive role. Cluster physics is just at the beginning of these investigations [12]. The experimental difficulties are substantially larger than in the case of bulk material. Nonetheless, with the steady development of femtosecond laser devices, more information on time resolved electron dynamics at short time scales will come up soon. Clusters as probes add a further interesting aspect in that the subsequent ionic motion dramatically depends on the initial electron dynamics, especially where the ionization rate of the cluster is concerned. In some circumstances, the ionic and electronic time scales might even interfere [13]. The description of this fast and large-amplitude electron dynamics is a challenge for a theoretical description. The enormous complexity of the many-body problem is much simplified in the framework of density functional theory (DFT) [14, 15] at the level of the time dependent local density approximation (TDLDA). The (TD)LDA has long been used in metal clusters with success for calculations of structure or optical response [1, 16]. The use of these methods in truly (i.e., nonlinear) dynamical problems [13, 1719], possibly with non-adiabatic coupling to ionic motion [13, 19] or within a semi-classical approach [20, 21], is more recent

SEMI-CLASSICAL ELECTRON DYNAMICS

213

and goes somewhat beyond the safely checked grounds of DFT. But the complexity of the problem and the simplicity of TDLDA leaves little other choice. Experience shows that TDLDA as a mean-field theory gives probably a correct view of the gross structure of the dynamical evolution [22]. Correlation effects, however, become important for a more detailed description. One can consider them as a perturbation to the effective mean-field theory and thus derive the correction in terms of a two-body collision integral [23] which can still carry all memory effects [24]. Furthermore, the time for the two-body collision shrinks dramatically with the excitation of the system [25] which soon allows us to neglect the memory effects and to simplify the correlations to a Boltzmann collision integral. Of course, one deals with a generalized version of the Boltzmann kernel which takes properly into account the blocking of occupied Fermion states and which is called the Uhling Uhlenbeck collision integral [23]. The limitation of the present approach to very energetic excitations allows the semi-classical approximation at the side of the TDLDA, the Vlasov-TDLDA. Together with the collision integral, we thus end up at the VlasovUhlingUhlenbeck (VUU) approach for which well developed numerical schemes exist [10]. This is the level of treatment which we will investigate here. We mention in passing that there are recent developments to describe dynamical correlations including memory effects and at a fully quantal level from a systematic expansion within DFT [26]. These are, however, still limited to the linear domain of small excitations. The paper is organized as follows. In Section 2, we introduce the theoretical framework of TDLDA, Vlasov-LDA, and VUU. We then turn to practical examples in Section 3, emphasizing the effect of electronelectron collisions by systematic comparisons to the mean-field description. Conclusions are drawn in Section 4.

2. THEORETICAL FRAMEWORK 2.1. Ionic Background Metal clusters are considered as bound systems of ions and valence electrons. A detailed description of structural properties requires the use of (possibly sophisticated) pseudo-potentials [27] accounting for the electron-ion interactions. Long time evolution will involve ionic motion and thus requires an explicit treatment of ionic dynamics. On the other hand, the huge difference between ionic and electronic masses allows us to neglect ionic motion on short time scales (typically up to t100 fs). In the situations considered below we can thus safely keep ions frozen [13]. Furthermore, in the case of violent excitations, as considered here, the details of ionic structure play only a minor role and become less and less important beyond the linear regime [28]. We thus further approximate the ionic background by a jellium droplet [1]. We use the ``soft'' version of the jellium model, which

214

DOMPS, REINHARD, AND SURAUD

provides the correct position of the Mie plasmon in TDLDA calculations [29]. The jellium background for a spherical cluster is then represented by the density * jel(r)=

* bulk 1+exp((r&R)(r c - 3))

,

(1)

where the bulk density is given by * bulk =3(4?r 3s ), and R is close to the sharp jellium radius r s N 13 (R is adjusted so that * jel(r) provides the correct total ionic charge and r s is the Wigner Seitz radius of the material). The finite surface width r c - 3 is defined via the Ashcroft core radius r c [30] (r c t1.73 a 0 for sodium). 2.2. Mean-Field Description of Electrons 2.2.1. From quantal to semi-classical TDLDA. We base our description of valence electrons on the time dependent version of density functional theory (TDDFT) [15], which has been used in cluster physics, both in the linear [1, 2] and non-linear domain [28, 31]. Furthermore, it has been shown that the simplest version of TDDFT, namely the time dependent local density approximation (TDLDA), turns out to provide a reliable and robust tool for exploratory investigations of highly non-linear dynamics [22]. The basic dynamical constituents of TDLDA are the time-dependent single-electronic wave functions , i (r) (i=1, ..., n) which obey the KohnSham equations [32], representing effective one-electron Schrodinger equations. A most compact notation is reached when introducing the one-body density matrix \^ = j |, j )( , j |, which reads, e.g., in coordinate space representation \^(r, r$)= nj=1 , j*(r$) , j (r). The KohnSham equation then reads i\*^ =[h(\^ ), \^ ]

(2)

in operator form, where h(\^ ) stands for the effective one-body Hamiltonian. TDDFT shows that under certain conditions h( \^ ) can be written as a functional of the diagonal part of the density matrix operator *(r)=\^(r, r), which makes Eq. (2) particularly handy. The operator form (2) of the KohnSham equations is better suited to formal manipulations, e.g., for going ``beyond'' mean-field in the realm of kinetic equations, or for deducing semi-classical approximations [33]. Both aspects will be explored here. 2.2.2. The Vlasov equation. We start with the semi-classical limit of TDLDA, called the Vlasov-LDA equation. It can be obtained formally from the quantal equation Eq. (2) through replacing the density operator \^ by a one-body phase space distribution f (r, p, t) and the commutator by Poisson brackets [33], \^(r, r$) Ä f (r, p, t) [ . , . ] Ä [ . , . ].

(3)

SEMI-CLASSICAL ELECTRON DYNAMICS

215

This leads to the Vlasov equation f p f V f + & =0, t m r r p

(4)

where V stands for the self-consistent (KohnSham) potential

|

V(r, t)=V ext(r, t)+ d 3r$

*(r$, t)&* jel(r$) +V xc[*](r, t). |r&r$|

(5)

The electronic density is now computed from the phase space distribution as

|

*(r, t)= d 3pf (r, p, t).

(6)

The KohnSham field Eq. (5) is composed of the external field V ext describing the excitation process (see below), the Coulomb field with both the direct Hartree contribution (electronelectron) and electron-ion part (Eq. (1) for * jel ), and the exchange correlation term (xc). The xc potential V xc[*](r, t) is a functional of the density * which is approximated in LDA (actually using here the functional of [34]). The Vlasov-LDA method is numerically less demanding than the full solution of the KohnSham equations for energetic cases and yields satisfying results for the most prominent aspects of the electron dynamics [20, 21]. A warning has nevertheless to be added here. The semi-classical approach wipes out any information on the Fermionic properties of the electrons. The Vlasov equation in finite numerical representations thus drives the system slowly towards a Boltzmann equilibrium state [35]. There are two ways out of this dilemma in the case of Vlasov calculations: either one introduces a dedicated Fermionic stabilization [36] or one postpones the critical time for the unwanted transition beyond the intended observation time by using a sufficiently large number of test particles (see Subsection 2.4). We consider here not too long observation times and thus employ the second option. It should nevertheless be noted that the VUU collision term does contain Pauli blocking factors (Subsection 2.3.1) which are sufficient to enforce Fermionic statistics, formally and numerically [37], so that VUU calculations are not affected by a possible loss of Fermi statistics, at least on the time scales considered here. 2.3. Beyond Mean-Field 2.3.1. The VUU equation. An advantage of the semi-classical Vlasov-LDA is that electronelectron scattering effects can easily be included as a Markovian collision term for the phase space distribution f. This has been worked out in the context of nuclear physics applications, where the problem of treating two-body collisions received much attention during the last decade [10]. In the spirit of

216

DOMPS, REINHARD, AND SURAUD

kinetic theories, the Vlasov-LDA equation can thus be complemented by a UehlingUhlenbeck [38] collision term to yield the VUU equation f p f V f + & =I UU (r, p, t) t m r r p

(7)

with the collision term I UU =

|

d 3p 2 d0 d_ |v 12 | [ f 1 f 2(1& f 3 2)(1& f 4 2))&f 3 f 4(1& f 1 2)(1& f 2 2))], (2?) 3 d0 (8)

where v 12 is the relative velocity of the colliding particles 1 and 2 and d_d0 is the differential cross section (depending on the scattering angle 0) evaluated in the center of mass frame of the two colliding electrons. Indices 3 and 4 label the moments of the two particles after an elementary collision and we use the standard abbreviation fi = f(r, p i , t). The collision takes place at one position r=r 1 =r 2 =r 3 =r 4 and the complementing moments p 3 and p 4 are deduced from p 1 and p 2 by conservation of energy, conservation of total momentum, and scattering angle 0. As known from solid-state physics [39], Pauli blocking factors (1&f i 2)(1&f j 2) play an important role for electronic systems. At zero temperature, all the collisions are Pauli blocked and the mean-free path of the electrons becomes infinite, at least where electronelectron collisions are concerned. At low temperature the collision rate can be predicted to increase quadratically with temperature [23]. But if the system becomes highly excited, phase space opens up widely, which strongly activates the collision term. As we shall see below, this typically Fermionic feature plays an important role in controlling the rate of electronelectron scattering in excited metal clusters. 2.3.2. The double counting problem. The use of VUU to describe electron dynamics in metal clusters is not (yet) a standard tool as it is in other fields such as nuclear [10] or plasma physics [40]. This raises various concerns from which the most prominent one is related to a possible double counting of interactions. The collision term is deduced from a second-order correlation diagram and the same diagram is contained formally also in the effective mean field through the energydensity functional for the LDA. This question has been much investigated in the other fields of physics where similar methods are employed. Liquid 3He, for example, constitutes a prototype of a dense Fermion system in which bubble diagrams play a decisive role, in a way similar to the case of the electron gas. It was worked out for that case in [41] that there is no double counting because the collision term (in Markovian approximation) accesses the imaginary part of the correlation diagram while the mean field uses the real part (and that only at zero frequency). Furthermore, it was shown that the collision integral should involve the screened interaction associated to the potential entering the mean field part. A discussion of the analogous problem in the

217

SEMI-CLASSICAL ELECTRON DYNAMICS

case of nuclei can be found in [24] leading basically to the same result. There is thus no danger of double counting if the collision term is treated in Markovian approximation which is the level of approach used here. We have to evaluate, however, the properly screened interaction to be used in the VUU collision term. 2.3.3. Estimate of the cross section. The computation of the VUU collision integral requires the knowledge of the differential cross section d_d0, possibly depending also on the relative velocity v 12 of the scattering electrons. Our present goal is to estimate the gross effects of two-body collisions. We thus simplify the treatment by assuming that the cross-section is isotropic _=_ tot 4?, and by taking its value at Fermi velocity v 12 rv F . This means that we describe all scattering processes by one global number for the effective cross section, a procedure which is commonly used in other applications [10]. The approximation of isotropic scattering is justified as long as we investigate global observables because the various scattering angles average in the course of time. Moreover, the actual cross section turns out to be dominated by s-wave scattering (see next paragraph) such that the angular dependence is very weak anyway. The neglect of velocity dependence is justified for cases of low temperatures because the scattering events are then confined to a narrow vicinity of the Fermi surface. It may become incorrect for more energetic cases where collisions cover a broader range of velocities. Our calculations are then to be considered as being exploratory in that regime of large excitations. This has to be kept in mind for the following discssions of the results. As pointed out above, the effective cross section has to be built from the screened electronelectron interaction V sc which can be estimated as follows. In a first stage, we solve the static KohnSham equations in a large spherical cluster and obtain the self-consistent Coulomb field V. In a second stage, we add a point-like negative charge and record the new Coulomb field V$. The screened potential emerges as V sc =V$&V. We find in actual calculations that this empirical V sc can be nicely fitted by a Yukawa shape er exp(&ra) with a & 3 a 0 , in good agreement with bulk values for sodium [42]. The cross section is then computed from V sc using standard quantum theory of scattering [43]. For electrons close to the Fermi level (v F &20 a 0 fs), the Born approximation fails badly and _ tot has to be calculated as a function of the phase shifts $ l produced by V sc on partial waves with angular momentum l. We have computed these phase shifts up to l=6 (see Table I), which is very sufficient with regard to the short range of the screened two-body potential (k F a &1.5). Taking into account the possible exchange of identical particles during 2

TABLE I Phase Shifts for Electrons with Fermi Velocity in a Screened Coulomb Potential 1

0

1

2

3

4

5

6

$l

0.964

0.346

0.149

0.072

0.029

0.019

0.005

218

DOMPS, REINHARD, AND SURAUD

a collision, and after averaging over spin degrees of freedom, we finally obtain _ tot =101 a 20 . 2.4. Numerical Details 2.4.1. Vlasov dynamics. The Vlasov and VUU equations are solved using the standard test particle method [44]. The phase-space distribution f of the n electrons is projected onto a swarm of & Gaussian functions with moving centers (r i (t), p i (t), i=1, ..., &) f (r, p, t)=

n & : g r(r&r i (t)) g p(p&p i (t)), & i=1

(9)

where g r and g p are normalized Gaussians g r(r&r$)=

1 (2?)

(32)

2

_

3 r

2

e &(r&r$) 2_r ,

g p(p&p$)=

1 (2?)

(32)

2

_

3 p

2

e &(p&p$) 2_p.

(10)

The solution of the Vlasov equation is then mapped into Hamiltonian equations of motion for the test particles (r i , p i ), moving in the mean field V C g r [35]. They read r* i =

pi m

(11)

p* i =&% r V C g r . The symbol C in Eq. (12) stands for the folding operation in coordinate space. These classical equations of motion (12) are solved with a standard leap-frog algorithm [45]. 2.4.2. VUU collision term. The algorithm to evaluate the VUU collision term in Eq. (7) is inherited from nuclear physics (where it was deduced from the intranuclear cascade [46]). The scheme is rather obvious and well documented [10]. The test particle representation simplifies the tedious integration over two-particle phase space by replacing it through a loop over pairs of test particles. The collision rate is governed by the cross section _ tot . It is translated into a scattering distance d phys =- _ tot ? which applies to the scattering of physical particles. The cross section _ needs to be rescaled to the scattering of test particles which yields the effective distance d _ =d _phys - n&. The actual sampling of the collision integral then proceeds in the following steps: 1. Loop over pairs of test particles 1 and 2. Preview the distance of closest approach d min that will emerge within the coming time step t Ä t+$t and select those cases which stay within the effective scattering disk, i.e., for which d min
SEMI-CLASSICAL ELECTRON DYNAMICS

219

2. Scattering of the two momenta (p 1 , p 2 ) Ä (p$1 , p$2 ) in accordance with conservation of energy and total momentum. The remaining freedom in the scattering angles is evaluated stochastically. 3. The final momenta have to be checked with respect to the Pauli blocking factors. But the latter are very difficult to handle technically since they require the distribution f in full 6 dimensional phase space. An appreciable simplification is obtained under the assumption that, for evaluating Pauli blocking factors, f (r, p, t) can be essentially determined by the energy ==p 22m+V(r, t) of the given phase space cell. In the Pauli blocking factor (1&f2), f can then be replaced by its average over an energy shell [36], i.e.,  f (r, p, t) $(=& p 22m&V(r, t)) d 3r d 3p f (r, p, t) Ä n(=)= . 2 2  $(=& p 22m&V(r, t)) d 3r d 3p

(12)

Equation (12) as such holds only when the electron cloud is at rest. When the electrons carry a collective motion with local velocity field u(r, t), we average f in the local rest frame, i.e., we define an energy shell by a fixed value = of (p&mu(r, t)) 22m+V(r, t) in Eq. (12). This provides a reasonable estimate of the average Pauli blocking effect. Finally, the decision for an event to take place is again sampled stochastically proportional to the amount of open phase space factors 1&n(=$i ). 2.5. Observables As mentioned in the Introduction, the most relevant observables for our purposes are dipole response, electron emission, and energy deposit. The number of escaped electrons is computed as [47] N esc(t)=N(t=0)&N(t),

N(t)=

|

d 3r *(r, t),

(13)

V

where N(t) is the number of electrons still residing in a certain analyzing volume V located in the vicinity of the ionic background. Actually, we consider a stripe of 2r s width around the cluster to define V. The question of the impact of the definition of V on the observables we consider has been thoroughly addressed in several papers and we do not repeat this discussion here [31, 47, 48]. The dipole moment D(t) of electrons is evaluated inside the same analyzing volume V [31] as D(t)=

|

d 3r r *(r, t).

(14)

V

The dipole signal D(t) carries all information about the optical response of clusters. This has been discussed at length elsewhere both in the linear and non linear regime, and both at the quantal and semi-classical level [20, 31, 49]. Spectral

220

DOMPS, REINHARD, AND SURAUD

properties can be deduced from the Fourier transformation from the time domain into the frequency domain. But here we will only need the raw signal D(t) as a complementing analyzing tool helping to understand the dynamical evolution. Energy transfer from collective dipole oscillations into intrinsic degrees-offreedom can be described in terms of the thermal energy E th acquired by electrons which we define as E th =

|

d 3p d 3r V

p2 f (r, p, t)& 2m

|

V

j 2(r) 3 d r& 2*

|

(3? 2 ) 23

V

3 2 *(r, t) 53 d 3r. 10m

(15)

The first term corresponds to the total kinetic energy of the electrons. The second term subtracts the collective kinetic energy evaluated from the current

|

j(r, t)= d 3p

p f (r, p, t) m

(16)

and the third term subtracts the local kinetic energy of a Fermi gas with zero temperature (contained inevitably in the Pauli pressure of Fermions). The difference thus characterizes the non-collective contributions to the kinetic energy which represent just the intrinsic heat of the electrons. It is known from nuclear physics that this quantity accounts fairly well for the effects of two-body collisions [50]. 3. RESULTS 3.1. Scaling with Particle Number There are no quantum size effects, as, e.g., shell corrections, in a semi-classical theory like VUU. All features will thus depend smoothly on the electron number n. In fact, all volume properties are independent of n and only small surface effects will add corrective trends B n &13. We have checked the n-dependence of the various observables discussed below and found only very small effects of the system size. We demonstrate that here for one example, namely electron emission following a sudden excitation of the electron cloud. A fast excitation can be realistically modelled by an instantaneous initial shift of the whole electron cloud against the ionic background, a method which we used in several previous investigations [31, 47]. For the sake of consistency, we consider here initial shifts delivering in all the 4 clusters the same excitation energy per electron (with respect to their respective ground-states) E*nt0.44 eV. The time evolution of the number of escaped electrons is shown for the four dif+ + + ferent clusters, Na + 9 , Na 21 , Na 41 , and Na 93 in Fig. 1. The results within each case, VUU or Vlasov, are all very similar in their trends and fairly comparable in numbers whereas the difference between the two cases (upper versus lower panel) is sizeable. We thus can conclude that there is some residual n-dependence which is due to the different ratio of surface to volume in the different clusters. But these differences are unimportant for our purposes of studying principle effects of the collision term. We thus present all further results for the simplest test case Na 9 + for which we performed the most extensive variations of the various excitation parameters.

SEMI-CLASSICAL ELECTRON DYNAMICS

221

FIG. 1. Time evolution of the number of escaped electrons as defined in Eq. (13) for four different Na + n clusters with electron number n as indicated. The lower panel shows results from full VUU propagation and the upper panel from Vlasov-LDA.

3.2. A First Estimate of Relevant Time Scales As we shall see below, the impact of the VUU collision term also depends on the details of the excitation mechanism where several parameters are to be varied: velocity, impact parameter, and charge of the projectile for ionic collisions, frequency, time profile, and fluency for lasers. The time scales of the cluster as such are best explored with the sudden excitation mechanism as sketched in the previous paragraph. This mechanism sets a well defined clock and allows one to control the degree of excitation by the amplitude of the initial shift. A typical result was shown in Fig. 1 and a detailed evaluation of the various time scales has been given in [37, 51]. We summarize here quickly the basic results. There are two one-body mechanisms which lead to a damping of the dipole signal (already at the level of mere TDLDA or Vlasov-LDA): direct electron emission and Landau damping. Both have the same origin in the fact that the dipole excitation is distributed over the various oscillatory eigenmodes of the system. The direct emission is related to those eigenmodes which reside in the electron continuum. Landau damping is related to the bound modes and provides the internal heating of the residual electron cloud. The time scale for direct electron emission is about { esc =13 fs, nearly independent of system size. The time scale for Landau damping

222

DOMPS, REINHARD, AND SURAUD

depends on system size due to quantum size effects. A typical value (for Nr100) is { MF =10 fs. Both processes depend only very weakly on temperature. The two-body collisions are absent in the first stages after excitation due to the still active Pauli blocking. It requires a pre-thermalization through Landau damping to heat the particle distribution and thus to activate the collision term. This already gives a first clue on the impact of dynamical correlations depending on the time scale of the excitation as such. For very short excitation processes, typically below { MF , one expects that the VUU collision term will have less effect on the dynamics whereas longer excitation processes will be much more sensitive to dynamical correlations. We shall see from the practical examples below that this is indeed the case, to some extent. The relaxation time due to two-body collisions { col can be deduced from comparing the Vlasov and VUU results for the approach of the thermal energy E th to equilibrium. Practically, one deduces this time by fitting an exponential law to E (VUU) &E (Vlasov) = th th &1 A(1&exp(&(t&t 0 ){ col )). The inverse relaxation time { col represents the collisional relaxation rate and the offset t 0 stands for the instant when VUU begins to deviate from Vlasov, which stems from a mixing of Landau damping time and the shorter direct emission time scale. One finds that, much unlike { MF , the collisional relaxation time { col strongly depends on excitation energy following a law [37] fs E* &1 & 0.3 eV { col N

.

(17)

The four cases of Fig. 1 have all the same specific excitation energy and thus the results for the relaxation times are compatible with each other, being { col & 8 fs throughout, a value which is also typical for most dynamical situations discussed below. 3.3. Charged-Projectile Induced Collisions Collisions of metal clusters with charged projectiles are typical examples of fast, violent excitations accompanied by plasmon oscillations and electron emission. These phenomena have been investigated experimentally in [4] and theoretically within TDLDA or Vlasov-LDA in [21, 20, 47]. It is thus of importance to check the modifications arising when proceeding to the more general VUU description. In the case of a grazing collision by a projectile, the cluster is exposed to a short electromagnetic pulse V ext(r, t)=

Q proj , |r&r proj(t)|

(18)

where Q proj. is the charge of the projectile and the ionic projectile position r proj(t) follows a Coulomb trajectory inside the electrostatic field of both electrons and ions of the perturbed cluster. Note that in the ionic velocity domain we consider

SEMI-CLASSICAL ELECTRON DYNAMICS

223

here (projectile velocity of order of the electron Fermi velocity of the cluster), the projectile trajectory is almost unperturbed [20]. The projectile motion follows basically a straight line. Feedback effects become important as soon as the projectile velocity is significantly lower than the Fermi velocity. In all the cases presented below the projectile trajectory has been explicitly determined dynamically to eliminate avoidable approximations. 3.3.1. Example of a rapid excitation. As a test case, we consider the cluster Na + 9 excited by the passage of a proton with kinetic energy 10 keV. The corresponding velocity 25 a 0 fs is close to the Fermi velocity of electrons inside the cluster. The impact parameter is chosen equal to 10 a 0 , comparable with the cluster radius R jell & 8 a 0 . This means that the projectile comes very close to the target, providing a fast, violent perturbation of the electron cloud at the instant of closest approach t=10 fs. The orientation is chosen such that the motion of the proton takes place in the x_z plane with the initial velocity parallel to the z-axis and the impact parameter b extending into the x direction. Figure 2 shows the emerging dipole response signals D i (t) along the longitudinal (z) and transverse (x) directions as a function of time. Both the Vlasov and the VUU results are presented for comparison. The number of ejected electrons is also shown in Fig. 3. The first stage of the interaction is an almost instantaneous motion of the electron cloud towards the projectile, with an amplitude higher than 1 a 0 in both the x and z directions. Then, the dipole amplitude quickly falls off, and less than 5 fs after the collision a second

FIG. 2. Dipole signals along the longitudinal z and transverse x directions during a protonNa + 9 collision, as obtained from Vlasov (full line) or VUU (dashed line) computations.

224

DOMPS, REINHARD, AND SURAUD

stage starts consisting in damped plasmon oscillations. We see from Fig. 2 that during the first stage (t<15 fs), the dipole signal is almost unaffected by the addition of the collision term. In other words, the excitation itself as well as the fast initial relaxation reflect pure mean-field dynamics. In particular, the very first excitation, corresponding to the growth of the dipole signal, appears here perfectly

FIG. 3. Numbers of emitted electrons (upper panel) and electron thermal energy E th (lower panel) as a function of time during a protonNa + 9 collision, as obtained from Vlasov (full line) or VUU (dashed line) computations.

SEMI-CLASSICAL ELECTRON DYNAMICS

225

independent of two-body collisions. Soon after that, phase space opens up widely during the cluster-ion interaction due to one-body dissipation (Landau damping, direct electron emission) which activates the collision term Eq. (8) and raises the collision rate per electron up to n coll n & 1 fs &1. But the excitation step is so short here (&1 fs) that this hardly modifies any of the relevant observables at that stage. Furthermore, the fast initial drop in the plasmon amplitude can be partly attributed to simple kinematic factors: the projectile attracts the electron cloud, which enhances the plasmon signal for x>0, but also tends to slow it down for x<0. Clearly, VUU collisions do not interfere with this immediate response mechanism. Another appreciable contribution to the initial damping stems from the electron emission which instantaneously releases part of the available excitation energy. We also see from Fig. 3 that this channel is little affected by two-body collisions. This again is correlated to the very short duration of the escape process (&1 fs). Let us now consider the second stage of the electronic motion, i.e., the relatively slow relaxation for t>15 fs. We see from Fig. 2 that VUU collisions significantly enhance the gradual extinction of plasmon oscillations. This could be expected from the general argument that a collision term is likely to reinforce thermalization. But the plasmon signal is relatively weak here, and it might have been inferred that the electron distribution remains quite close to equilibrium, thereby suppressing twobody collisions through the Pauli blocking factors. Actually, these small amplitude vibrations are performed around equilibrium, but for a hot electron cloud, as it emerges from the first stage of the motion. This appears clearly by looking at the thermal energy E th Eq. (15) as it is plotted in Fig. 3. At t=12 fs, almost 2 eV have already been stored as thermal energy, being shared among the 7 remaining electrons only. The electron distribution is thus already far away from its cold ground state, thereby opening phase space for two-body collisions. After the cluster-proton interaction is over, the thermal energy remains constant within the collisionless Vlasov propagation but continues to increases in the case of VUU where electron electron scattering is taken into account. Thermally generated collisions cause a gradual transfer of energy from mean-field towards disordered, thermal motion. At t=50 fs, VUU yields a hot, thermalized electron cloud. By fitting the final energy distribution n(=) to a FermiDirac shape, we obtain a final temperature of 0.8 eV to be compared to the Fermi energy 3.2 eV. At the same instant, for pure Vlasov propagation, the electron fluid is still oscillating with half the thermal energy. From this first example, we conclude that correlation effects beyond mean-field are inefficient to modify a fast (5 fs) excitation process with respect to the collisionless Vlasov-LDA picture. A quantity such as direct ionization by a projectile, which is governed by a quasi-instantaneous mechanism, therefore remains essentially unaffected. But two-body collisions strongly enhance the subsequent dissipation of the plasmon energy into thermal degrees of freedom of the electron cloud. 3.3.2. Slower collisions. The situation is different for collisions which last several tens of femtoseconds. To illustrate that, we now consider the passage of a relatively slow (v & 2.5 a 0 fs) Ar 8+ ion in the vicinity of the cluster, with an impact parameter

226

DOMPS, REINHARD, AND SURAUD

8+ FIG. 4. Dipole signals in Na + ion with velocity 2.5 a 0 fs and impact 9 during the passage of a Ar parameter b=30 a 0 , as obtained within Vlasov (full line) and VUU (dashed line).

b=30 a 0 . By comparison to the previous example, this leads to comparable magnitudes for the dipole signal and number of emitted electrons, but stretched over a much longer duration. In Fig. 4, D x and D z are plotted as a function of time. The motion of the electron cloud is a quasi-adiabatic polarization and depolarization towards the attracting projectile, with hardly any plasmon oscillations. From Fig. 4 we see that this slow electronic motion is almost unaffected by the addition of two-body collisions. This can be simply explained in terms of Pauli blocking effects. Since the electron cloud experiences an adiabatic motion, its local Fermi equilibrium is hardly disrupted and n(=) remains close to its initial cold Fermi Dirac shape. This, in turn, hinders the development of two-body collisions. At second glance, however, we note that for t>50 fs, the lowering of D x proceeds a bit faster within VUU than within Vlasov. Slight collisional effects also show up in the amount of electrons emitted during the electron-ion interaction, as plotted in Fig. 5. For t>50 fs, the VUU result gradually departs from Vlasov and finally reaches a value which is by 2025 0 higher. The difference represents the contribution

SEMI-CLASSICAL ELECTRON DYNAMICS

227

FIG. 5. Number of electrons emitted from Na + 9 (upper panel) and generated electronic thermal energy (lower panel) during the passage of a Ar 8+ ion with velocity 2.5 a 0 fs and impact parameter b=30 a0 as obtained within Vlasov (full line) and VUU (dashed line).

of thermal emission to N esc , associated to two-body collisions. The VUU term in Eq. (7) plays a statistical role to fill the tail of the FermiDirac distribution and sends electrons into the continuum. These two features appearing for t>50 fs are a consequence of two-body friction and finally show up as a sizable extra deposited thermal energy. In Fig. 5, E th is plotted as a function of time and it exhibits a sensitive warming for t>50 fs. In the Vlasov description, the cluster is brought back

228

DOMPS, REINHARD, AND SURAUD

to a ``cold'' state after depolarization. Two-body collisions, on the other hand, cause a warming of the electrons associated to a higher charge state. The effect remains, however, rather small here (E th 2 eV), because adiabatic motion implies that phase space scarcely opens up for two-body collisions. Still, a combination of a slight opening of phase space together with a long associated time scale seems enough to make a noticeable effect in the thermal excitation energy. From these two examples, we conclude that electronelectron scattering effects remain quite limited during cluster-ion interactions, at the ``high'' projectile velocities considered here. The kinematic regimes investigated here are either so short that two-body collisions hardly step in, or so slow that they hardly develop. Still, on the long run, they produce a visible warming of the electron cloud associated to a faster decay of plasmon oscillations. To probe essentially different situations, we shall now consider excitations by femtosecond laser pulses. They expose the electrons to repeated, non-adiabatic oscillations, and should thus cause a more spectacular development of collisional effects. 3.4. Excitation by Intense Laser Pulses Thanks to recent progress in laser technology, high intensity pulses with duration as short as a few femtoseconds are now available [52]. In cluster science, this opens the door to new experimental studies focusing on multiple ionization and Coulomb explosion of those few particle systems [5]. This ultra-fast laser-induced electronic motion has recently been studied by means of TDLDA [53], and it was shown that the cluster response is highly sensitive to the relation between laser frequency and Mie plasmon frequency. Here we want to probe the modifications brought to this mean-field description by the addition of the VUU collision term. As above, we consider Na + 9 , this time excited by a linearly polarized (z-axis) laser beam. The corresponding external field reads, in dipole approximation V ext(r, t)=E 0 f (t) z sin(|t).

(19)

The pulse envelope f (t) chosen here is a Gaussian with width { pulse =10 fs. Such a short duration furnishes a proper decoupling of ionic and electronic motions, as explained before. The instant of peak maximum is t 0 =30 fs. 3.4.1. Off resonance irradiation. As a first example, we consider a laser pulse with intensity 1.2 10 13Wcm 2 and frequency | las =5 eV, about 2 eV above the plasmon resonance. In Fig. 6, the dipole moment D z and the number of emitted electrons N esc , as obtained within Vlasov and VUU, are plotted as a function of time. The dashed lines indicate the envelope f(t) of the laser pulse. In this far offresonance regime, dipole oscillations die out as soon as the laser excitation switches off. The motion is actually identical to that of a classical harmonic oscillator with spring constant at the plasmon frequency [53]. It is hence no surprise that VUU and Vlasov dipole responses look very similar in Fig. 6. The result for N esc deserves more comments. In contrast to the quasi-instantaneous sucking of electrons by the by-passing proton, the ionization process lasts here several tens of femtoseconds

SEMI-CLASSICAL ELECTRON DYNAMICS

229

FIG. 6. Dipole signal along the direction of laser polarization (D z , upper part) and number of emitted electrons (N esc , lower part) as a function of time during the interaction of a femtosecond laser pulse with 13 Wcm 2, and the pulse envelope Na + 9 . The laser frequency and fluency are | las =5 eV and J=1.2 10 is indicated in dashed line. Both the Vlasov and VUU results are plotted for comparison.

and should be considered as a plasmon enhanced emission. Since both Vlasov and VUU predict very similar dipole amplitudes, they also yield comparable N esc , at least for t<40 fs. In analogy to the Ar 8+ case, we nevertheless notice for the VUU case a slight enhancement of N esc associated to thermal emission. For t>40 fs, this thermal emission dominates and the VUU result gradually departs from the Vlasov curve, the difference reaching some 0.4 electrons at t=100 fs. The subsequent evolution will gradually couple to ionic motion and goes beyond the present treatment. Still, we already see here that two-body collisions can have a sensitive effect on cluster ionization from laser pulses. 3.4.2. On resonance irradiation. To have an essentially different case, we now consider a laser pulse whose frequency matches the plasmon resonance (| las & 2.75 eV). A weaker laser fluency (6_10 11 Wcm 2 ) is then sufficient to cause

230

DOMPS, REINHARD, AND SURAUD

about the same large amplitude oscillations of the valence electrons as before. The resulting D z (t) and N esc(t) are shown in Fig. 7. Due to the resonant laser frequency, the first tens of femtoseconds (t25 fs) are characterized by a fast growth of plasmon oscillations, independent of the two-body collisions. But then the amplitude reaches a maximum around 1 a 0 in VUU, whereas it increases up to

FIG. 7. Dipole signal along the direction of laser polarization (D z , upper part) and number of emitted electrons (N esc , lower part) as a function of time during the interaction of a femtosecond laser 11 2 pulse with Na + 9 . The laser frequency and fluency are | las =2.75 eV and J=6 10 Wcm , and the pulse envelope is indicated by a dashed line. Both Vlasov and VUU results are plotted for comparison.

SEMI-CLASSICAL ELECTRON DYNAMICS

231

FIG. 8. Electron thermal energy E th as a function of time during the interaction of a femtosecond laser pulse with Na + 9 , as obtained within Vlasov-LDA (full line) and VUU (dashed line).

1.7 a 0 in collisionless Vlasov. With growing amplitude, phase space gradually opens up, releasing the collision term and redirecting more and more energy into heat. This can be checked in Fig. 8, where we see the extra thermal energy generated by two body collisions beyond t>25 fs. In fact, the additional energy provided by the laser is instantaneously thermalized by electronelectron collisions, a mechanism which is absent in the mean-field description. The effect of two-body collisions is thus markedly different from what we had at | las =5 eV. In fact, the subtle resonant behaviour investigated here appears much more sensitive to collisional friction than the off-resonance dynamics simulated in the previous case. For t>40 fs, both Vlasov and VUU predict a damping of plasmon oscillations. In analogy to Subsection 3.3, it is however much faster within VUU: a complete vanishing of plasmon vibrations is obtained very soon after laser extinction. This effect is even more pronounced than in the proton case, because of higher thermal energy (8 versus 4 eV) and temperature T (1.3 versus 0.8 eV). Since the collision frequency grows like T 2, collisional damping is sensitively enhanced as compared to the cases of Subsection 3.3. Looking now at N esc(t), we realize that the smaller plasmon amplitude causes a much reduced cluster ionization. At t=60 fs (laser extinction), we obtain N esc =1.8 in Vlasov versus N esc =1 in VUU. Again, we are seeing a collisional effect which was absent in the previous ``brute force'' example. For t>60 fs, thermal emission takes over to bring the VUU result closer to the Vlasov values. The effect of twobody collisions thus appears as a ``retarded'' electron emission, thermal energy

232

DOMPS, REINHARD, AND SURAUD

being first stored during laser-cluster interaction, and later on released as electron evaporation, whereas direct emission has a stronger weight in pure Vlasov propagation. As explained above, the asymptotic value of N esc cannot be extrapolated within the present model dealing with fixed ions. But it is obvious that this retarded electron emission could produce sensitive modifications of the ion dynamics with respect to the collisionless LDA result, showing the importance of collisional effects in that case. 3.4.3. Influence of laser parameters. We now examine briefly the influence of laser parameters around the resonance case. Figure 9 is a comparison of results obtained when varying either the laser frequency (| las Ä 2.6 eV) or fluency (J Ä 1.5 10 11 Wcm 2 ) with respect to the previous resonant example. These two

FIG. 9. Comparison of cluster response to laser irradiation for two sets of laser parameters: (| las = 2.6 eV, J=6.10 11 Wcm 2, left column) and (| las =2.75 eV, J=1.5.10 11 Wcm 2, right column). For each set, we plot the electron dipole moment D z and number of emitted electrons N esc and compare Vlasov to VUU predictions.

SEMI-CLASSICAL ELECTRON DYNAMICS

233

additional sets of laser parameters provide roughly the same amplitudes for the mean-field dynamics and are thus adequate for comparison. The results confirm the trends seen above. Two-body collisions accelerate the decay of plasmon oscillations. They also reduce electron emission during laser irradiation and postpone part of it to later thermal evaporation. As can be seen from Fig. 9, both these features become more pronounced for a frequency exactly matching the plasmon resonance. This is most visible in N esc , which experiences a dramatic reduction by two-body collisions at | las =2.75 eV. This is in agreement with what we have said before in the far off-resonance regime, namely that off-resonant behaviours are less affected by collisional effects.

4. CONCLUSION We have investigated the influence of two-body collisions in the electron dynamics of metal clusters, in various dynamical regimes corresponding to realistic excitation processes. The starting point of these investigations is the Vlasov-LDA, namely the semi-classical counterpart of TDLDA, which is expected to provide a proper account of the gross features of the dynamics in the non-linear domain. The Vlasov-LDA also offers a suitable starting point for extensions to include dynamical correlations by kinetic equations, yielding the VlasovUhlingUhlenbeck (VUU) approach. We have implemented the VUU scheme for metal clusters using the experience gathered in other fields of physics where similar situations are encountered. The basic ingredient is the effective electronelectron cross section which we have computed microscopically and which was then taken approximately as isotropic and velocity independent. In the first stage of the investigation we employed a schematic and generic model for the excitation process which is based on an instantaneous shift of the electron cloud and thus sets a well defined clock. This allows to read off the typical time scales involved in VUU as compared to Vlasov-LDA. It turns out, e.g., that the time-scale for collisional relaxation in VUU depends strongly on the actual temperature of the system. In the further stages, we have investigated two realistic excitation schemes, namely collision with a fast ion and laser irradiation. The time scales found in the schematic excitation model are basically confirmed and we have worked out, in addition, the various effects from the interference with the excitation time. Alltogether, collisional effects are found to depend critically on the dynamical regime as well as on the typical time scales of the excitation process, namely overall duration { pulse and oscillation time 1| ext . The most dramatic modifications to the mean-field behaviour are obtained when the excitation lasts long enough ({ pulse & 10100 fs) and drives the electrons at a sufficiently high pace (| ext & plasmon frequency) to break local equilibrium. Typical examples for that are the excitations by fs laser pulses. The most spectacular effects obtained in that case consist of a sensitive reduction of plasmon oscillations by the electronelectron collisions, associated to the

234

DOMPS, REINHARD, AND SURAUD

appearance of a retarded thermal evaporation. These phenomena are particularly pronounced for laser frequencies close to the plasmon resonance and could become important for simulating the non-BornOppenheimer fragmentation of metal clusters after a fast electron excitation. To conclude, these first investigations have clearly demonstrated that two-body collisions can have a large influence on the electronic and ionic dynamics following violent excitations of clusters. Investigations in even more extreme situations with even higher energy deposit are desirable. Work along this line is in progress.

ACKNOWLEDGMENTS The authors thank FrenchGerman exchange program PROCOPE 95073 and Institut Universitaire de France for financial support during the realization of this work.

REFERENCES 1. 2. 3. 4. 5. 6. 7. 8. 9. 10. 11. 12. 13. 14. 15. 16. 17. 18. 19. 20. 21. 22. 23. 24. 25. 26. 27. 28. 29. 30. 31.

M. Brack, Rev. Mod. Phys. 65 (1993), 677. W. de Heer, Rev. Mod. Phys. 65 (1993), 611. U. Naher et al., Phys. Rep. 285 (1997), 245. F. Chandezon et al., Phys. Rev. Lett. 74 (1995), 3784. B. Lang et al., Z. Phys. D 40 (1997), 1. T. Reiners et al., Phys. Rev. Lett. 74 (1995), 1558. T. Hertel et al., Phys. Rev. Lett. 76 (1996), 535. D. Steinmuller-Nethl et al., Phys. Rev. Lett. 68 (1992), 389. R. H. M. Groenveld et al., Phys. Rev. B 51 (1995), 11433. G. F. Bertsch and S. Das Gupta, Phys. Rep. 160 (1988), 189. E. Suraud, Ch. Gregoire, and B. Tamain, Prog. Part. Nucl. Phys. 23 (1989), 357. T. Reiners et al., Phys. Rev. Lett. 77 (1996), 2440. F. Calvayrac, P. G. Reinhard, and E. Suraud, J. Phys. B 31 (1998), 5023. R. M. Dreizler and E. K. U. Gross, ``Density Functional Theory,'' Springer-Verlag, Berlin, 1990. E. K. U. Gross et al., Top. Curr. Chem. 181 (1996), 81. W. Ekardt, Phys. Rev. Lett. 21 (1984), 1925. F. Calvayrac, P. G. Reinhard, and E. Suraud, J. Phys. B 31 (1998), 1367. K. Yabana and G. F. Bertsch, Phys. Rev. B 54 (1996), 4484. U. Saalmann and R. Schmidt, Z. Phys. D 38 (1996), 153. L. Feret et al., J. Phys. B 29 (1996), 4477. M. Gross and C. Guet, Z. Phys. D 33 (1995), 289. C. Ullrich, P.-G. Reinhard, and E. Suraud, J. Phys. B 31 (1998), 1871. D. Pines and P. Nozieres, ``The Theory of Quantum Liquids,'' AddisonWesley, Reading, MA. K. Gutter, K. Wagner, P.-G. Reinhard, and C. Toepffer, Ann. Phys. (N.Y.) 225 (1993), 339. C. Greiner, K. Wagner, and P.-G. Reinhard, Phys. Rev. C 49 (1994), 1693. G. Vignale, C. A. Ullrich, and S. Conti, Phys. Rev. Lett. 79 (1997), 4878. G. B. Bachelet, D. R. Hamann, and M. Schluter, Phys. Rev. B 26 (1982), 4199. F. Calvayrac, P. G. Reinhard, and E. Suraud, Phys. Rev. B 52 (1995), R17056. P.-G. Reinhard, O. Genzken, and M. Brack, Ann. Phys. (Leipzig) 5 (1996), 1. N. W. Ashcroft, Phys. Lett. 23 (1966), 48. F. Calvayrac, P. G. Reinhard, and E. Suraud, Ann. Phys. (N.Y.) 255 (1997), 125.

SEMI-CLASSICAL ELECTRON DYNAMICS 32. 33. 34. 35. 36. 37. 38. 39. 40. 41. 42. 43. 44. 45. 46. 47. 48. 49. 50. 51. 52. 53.

235

W. Kohn and L. J. Sham, Phys. Rev. 140 (1965), 1133. A. Domps et al., Ann. Phys. (Leipzig) 6 (1997), 455. O. Gunnarson and B. I. Lundqvist, Phys. Rev. B 13 (1976), 4274. A. Domps et al., Ann. Phys. (Leipzig) 6 (1997), 468. A. Domps et al., Ann. Phys. (N.Y.) 260 (1997), 171. A. Domps, P.-G. Reinhard, and E. Suraud, Phys. Rev. Lett. 81 (1998), 5524. E. A. Uehling and G. E. Uhlenbeck, Phys. Rev. 43 (1933), 552. N. W. Ashcroft and N. D. Mermin, ``Solid State Physics,'' Saunders, Philadelphia. S. Ichimaru et al., Phys. Rep. 149 (1987), 91. J. W. Serene and D. Rainer, Phys. Rep. 101 (1983), 221. D. Pines, Phys. Rev. 92 (1953), 626. L. D. Landau and E. Lifschitz, ``Quantum Mechanics,'' Pergamon, London, 1959. W. G. Hoover, ``Computational Statistical Mechanics,'' Elsevier, AmsterdamNew York, 1991. L. Verlet, Phys. Rev. 159 (1967), 89. J. Cugnon, T. Mizutani and J. Vandermeulen, Nucl. Phys. A 352 (1981), 505. C. A. Ullrich, P.-G. Reinhard, and E. Suraud, Euro. Phys. J. D 1 (1997), 303. F. Calvayrac, A. Domps, P. G. Reinhard, E. Suraud, and C. Ullrich, Euro. Phys. J. D 4 (1998), 207. A. Domps, P. G. Reinhard, and E. Suraud, Euro. Phys. J. D 2 (1998), 191. B. Remaud, Ann. Phys. Fr. 21 (1996), 503. P.-G. Reinhard et al., Euro. Phys. J. D, in press. M. Protopapas et al., Rep. Prog. Phys. 60 (1997), 389. C. A. Ullrich, P.-G. Reinhard, and E. Suraud, J. Phys. B 30 (1997), 5043.