Spin and density longitudinal response of quantum dots in time-dependent local-spin-density approximation

The longitudinal dipole response of a quantum dot has been calculated in the far-infrared regime using local spin density functional theory. We have studied the coupling between the collective spin and density modes as a function of the magnetic field. We have found that the spin dipole mode and single particle excitations have a sizeable overlap, and that the magnetoplasmon modes can be excited by the dipole spin operator if the dot is spin polarized. The frequency of the dipole spin edge mode presents an oscillation which is clearly filling factor ($\nu$) related. We have found that the spin dipole mode is especially soft for even $\nu$ values, becoming unstable for magnetic fields in the region $1<\nu \leq 2$. Results for selected number of electrons and confining potentials are discussed. An analytical model which reproduces the main features of the microscopic spectra has been developed.


I. INTRODUCTION
The far-infrared ͑FIR͒ response of quantum dots is a subject of current interest since the experiments carried out by Sikorski and Merkt 1 and by Demel et al. 2 These experiments and subsequent theoretical work ͑see Refs. 3-8 and references therein͒ showed that the excitation spectrum of quantum dots in the FIR region is dominated by the dipole edge magnetoplasmon peak that splits into two different B dispersion branches when a magnetic field B is applied perpendicularly to the dot. These peaks are density ͑charge͒ collective modes excited by the operator D ϭ ͚ iϭ1 N x i . In the case of harmonic confinement by the potential 1 2 m 0 2 r 2 , as a consequence of Kohn's theorem 9 the density mode is not coupled to any other mode, and the dipole operator D excites only two collective states at the energies Ϯ ϭ ͱ 0 2 ϩ 1 4 c 2 Ϯ 1 2 c , where c is the cyclotron frequency. If the confining potential is not harmonic, Kohn's theorem does not hold. On the one hand, the energy of the modes depends on the number of electrons in the dot, and, on the other hand, a richer excitation spectrum appears. Raman spectroscopy has made it possible to observe in the same sample single-particle ͑sp͒, charge, and spindensity excitations, 10,11 whose evolution as a function of B has been studied in recent experiments. This has revealed several interesting features of the sp ͑Ref. 12͒ and of spin collective excitations 13 in quantum dots. Limiting ourselves to the latter, the experiments have determined that the spin mode lies very close in energy to the uncorrelated singleelectron excitations, and that magnetoplasmons can also be detected using spin-dependent probes. In addition, it has been experimentally determined that the spin mode has a much lower energy than the charge mode. These facts constitute the body of experimental results that any theory aim-ing at a quantitative simultaneous description of spin-and charge-density collective modes in quantum dots should reproduce.
The dipole spin response function for unpolarized quantum dots at zero magnetic field was recently addressed by two of us. 14 In the FIR regime, it has been found that the response is dominated by a low-energy collective dipole spin mode excited by the operator D m ϭ ͚ iϭ1 N x i z i , where x i and z i are Cartesian components of the position and spin vectors, and N is the number of electrons in the dot. Similar modes have been described in atomic nuclei, 15 and in alkalimetal clusters. 16,17 The aim of the present work is to extend our previous study to the case of a quantum dot submitted to a perpendicular static magnetic field, which originates a B-dependent spin polarization in the ground state of the dot. We will explicitly show that this not only causes the splitting of the spin dipole mode into two branches, one with negative B dispersion and another with a positive B dispersion, but also its coupling with the dipole density mode mainly excited by the operator D . We shall see that if the confinement is not harmonic and the dot is polarized, that operator also excites the dipole spin mode. Conversely, when the dot is polarized, which is the case if it has an odd number of electrons, or for most cases when B acts on the dot, the spin response is coupled to the density response so that the external operator D m ϭ ͚ i x i z i also excites the density mode. When the system is fully polarized both modes coincide, while at zero polarization they are uncoupled.
To this end, we have self-consistently evaluated the longitudinal response of the dot in the time-dependent localspin-density approximation ͑TDLSDA͒. By longitudinal we mean an external field which is either spin independent, or dependent on the spin component parallel to the magnetic field, i.e., the z component. We present results obtained for selected numbers of electrons and confining potentials. Specifically, we have used a harmonic oscillator potential to describe dots with Nϭ5 and 25, and a disk-confining potential to describe dots with 25 and 210 electrons, for which the FIR response has been determined in detail. 2 The results obtained for the Nϭ5 dot have been presented as preliminary results in Ref. 18. The ground state ͑gs͒ structure of the later two dots in intense magnetic fields has been recently addressed. 19,20 However, to our knowledge, no selfconsistent TDLSDA calculation for a dot as large as N ϭ210 has been carried out before even in the density channel.
To obtain correct collective modes one needs to have a proper description of the ground state these excitations are built on. Several density-functional calculations have addressed this question. [21][22][23][24][25] The LSDA we use in the present work is based on the exchange-correlation energy functional employed in Ref. 21 as an input to construct their current density-functional theory ͑CDFT͒. Within the range of magnetic fields we are interested in, we have checked that both the LSDA and CDFT yield similar results for gs properties other than the current density. Tests of the CDFT against exact and Hartree-Fock ͑HF͒ calculations have been presented in Ref. 21. Tests of unrestricted HF against exact gs energies in the filling factor region 2уу1 were also presented in Ref. 26 for a small number of electrons ͑up to 5͒.
We conclude, from the comparisons presented in the above references, that the TDLSA may yield fairly accurate results for the density and spin response in the range of magnetic fields for which experimental information is currently available. Comparison with these experiments constitute the ultimate test of this essentially parameter-free approximation.
Like Kallin and Halperin 27 and MacDonald, 28 who thoroughly studied the spin and density response of the twodimensional electron gas ͑2DEG͒, we mainly addressed the longitudinal response of quantum-dot configurations corresponding to integer filling factors. In a few cases we have considered configurations which are the finite-size analog of partially filled 2DEG configurations. In these cases, the use of the TDLSDA may be questioned, since these configurations are believed to have very complicated, strongly correlated ground states. Hence one has to regard the corresponding results as qualitative extensions of the ones obtained at integer filling factors. Although qualitative, these results may yield the general trends of the excitation spectrum, and for this reason here we have considered some cases of this kind.
To help understand the microscopic spectra, simpler methods have been developed to reproduce their gross features. One such method is the sum rule approach used to describe multipole density modes in quantum dots. 19 Here we present an analytical model, called vibrating potential model ͑VPM͒, which provides a complementary physical insight into the longitudinal response of quantum dots. The VPM model has been widely used in nuclear physics, 29 where it was developed to describe nuclear collective modes. It has also been applied to the description of simple metal clusters. 30

II. LONGITUDINAL RESPONSE
We consider a quantum dot made of N electrons moving in the zϭ0 plane, where they are confined by the circular potential V ϩ (r) in the presence of a constant magnetic field B in the z direction. In the LSDA, the single electron wave functions are given by the solution of the Kohn-Sham ͑KS͒ equations where V H ϭ͐dr ជ Ј(r ជ Ј)/͉r ជ Ϫr ជ Ј͉ is the Hartree potential. V xc ϭ‫ץ‬E xc (,m)/‫͉ץ‬ gs and W xc ϭ‫ץ‬E xc (,m)/‫ץ‬m͉ gs are the variations of the exchange-correlation energy density E xc (,m) in the local approximation taken at the ground state, and (r) and m(r) are the electron and spin magnetization densities. The exchange-correlation energy density E xc has been constructed from the results on the nonpolarized and fully polarized 2DEG ͑Ref. 31͒ in the same way as in Ref. 32, i.e., using the two-dimensional von Barth and Hedin 33 prescription to interpolate between both regimes. We have used effective atomic units (បϭe 2 /⑀ϭmϭ1), where ⑀ is the dielectric constant of the semiconductor, and m is the electron effective mass. In units of the bare electron mass m e one has mϭm*m e . In this system of units, the length unit is the effective Bohr radius a 0 *ϭa 0 ⑀/m*, and the energy unit is the effective Hartree H*ϭHm*/⑀ 2 . For GaAs we have taken ⑀ϭ12.4, m*ϭ0.067, and g*ϭϪ0.44, which yields a 0 *ϭ97.9 Å and H*ϳ11.9 meV ϳ95.6 cm Ϫ1 . In Eq.
͑1͒, c ϭeB/(mc) is the cyclotron frequency and B ϭeប/(2m e c) is the Bohr magneton. The use of the same letter for the effective mass and the spin magnetization, and for the dielectric constant and the single-electron energies should cause no confusion, since neither the mass nor the dielectric constant will explicitly appear in the rest of the work.
As a consequence of circular symmetry, the ␣ 's are eigenstates of the orbital angular momentum l z , i.e., ␣ (r,)ϭu nl (r)e Ϫil , with lϭ0,Ϯ1,Ϯ2, . . . . The gs electron density is given by (r)ϭ͚ ␣ n ␣ ͉u ␣ (r)͉ 2 , while the gs spin magnetization density is expressed in terms of the spin of orbital ␣, ͗ z ͘ ␣ , as m(r)ϭ͚ ␣ n ␣ ͗ z ͘ ␣ ͉u ␣ (r)͉ 2 . The numerical calculations reported in the following have been performed at a small but finite temperature Tр0.1 K, and the KS equations have been solved by integration in r space. The thermal occupation probabilities n ␣ are determined by the normalization condition which fixes the chemical potential . Our iterative method works for weak and strong magnetic fields as well, for which the effective potential is very different. It has proved to be very robust, allowing us to handle without any problem several hundreds of electrons.
Once the gs has been obtained, we determine the induced densities originated by an external field employing linearresponse theory. Following Refs. 34 and 35, we can write the variation ␦ induced in the spin density (ϵ↑,↓) by an external spin-dependent field F, whose nontemporal dependence we denote as Fϭ ͚ f (r ជ )͉͉͗͘: where Ј is the spin-density correlation function. In this limit, the frequency corresponds to the harmonic time dependence of the external field F and of the induced ␦ .
Equation ͑3͒ is a 2ϫ2 matrix equation in the two-component Pauli space. In longitudinal response theory, F is diagonal in this space, and we write its diagonal components as a vector For the operators defined in Sec. I, we then have

͑4͒
The TDLSDA assumes that electrons respond as free particles to the perturbing effective field, which consists of the external plus the induced field arising from the changes produced by the perturbation in the gs mean field. This condition defines the TDLSDA correlation function Ј in terms of the free-particle spin-density correlation function Ј (0) through a Dyson-type integral equation The free-particle spin-correlation function at finite temperature is obtained from the KS sp wave functions, energies, and occupation probabilities:

͑6͒
The label ␣ (␤) refers to a sp level with spin (Ј) and occupation probability n ␣ (n ␤ ). To simplify the analysis of the results, we have added a small but finite imaginary part to the energy . This will make an average of the strength function by transforming the ␦ peaks into Lorentzians of width 2. The kernel K Ј (r ជ ,r ជ Ј) is the residual two-body interaction with ↑ ϭ1 and ↓ ϭϪ1. The last expression is the definition of the K, L, and I functions.
When the system is not polarized, there are only two independent correlation functions. These are and mm , describing, respectively, the density response to D and the spin response to D m . They are given by and the four equations ͑5͒ reduce to two uncoupled equations for and mm , whose kernels are given by 1/r 12 ϩK␦(r 12 ) and I␦(r 12 ), respectively, and the free-particle . This constitutes the paramagnetic limit of the longitudinal response with uncoupled density and spin channels, 14 in which the residual interaction consists of a Coulomb direct plus an exchange-correlation terms in one case, and only of an exchange-correlation term in the other.
When the system is polarized one no longer has ↑↑ ϭ ↓↓ (0) , and there are two more independent correlation functions which produce the density response to D m and the spin response to D , respectively.
Equations ͑5͒ have been solved as a generalized matrix equation in coordinate space after performing an angular decomposition of Ј and K Ј of the kind Only modes with lϭϮ1 couple to the external dipole fields D and D m . This can be readily seen performing the angular integral in Eq. ͑3͒. In practice, we have considered the multipole expansion of the external field, using the dipole vectors For a polarized system having a nonzero magnetization in the gs, the lϭϮ1 modes are not degenerate, and give rise to two excitation branches with ⌬L z ϭϮ1, where L z is the gs orbital angular momentum.
The response functions corresponding to the above dipole fields have been obtained from the lϭϮ1 components of the correlation functions AB (Ϯ1) (r,rЈ;) with A,Bϭ,m as Their imaginary parts are related to the strength functions as S AB ()ϭ(1/)Im͓␣ AB ()͔. Although the excitation energy and strength functions are always positive, it may be easily verified that the following relations formally hold: To check the numerical accuracy of the calculations we have used the f-sum rules for the dipole operators, which can be expressed in terms of gs quantities: 36 where S z is the total spin of the ground state.

III. RESULTS
Figures 1-3 display the dipole strength function of the Nϭ5, 25, and 210 dots for selected B values. Solid lines correspond to the density response to D , and dotted lines to the spin response to D m , that is to S and to S mm . Dashed lines represent the free-particle strength function. For the five-electron dot we have used a parabolic confining potential V ϩ (r)ϭ 1 2 m 0 2 r 2 , with 0 ϭ4.28 meV, and for the other dots we have used a disk confining potential. 19,20 For Nϭ25 and 210, most B values displayed correspond to integer filling factors . It was found 20 that for Nϭ210 and an R disk confining potential, these values follow the law ϭ2cn s /(eB) pertaining to the 2DEG, where n s ϭN/(R 2 ) is the electron surface density. For the Nϭ25 dot, that law yields Bϭ3.29 T as the value at which the system becomes fully polarized. Actually, we have found that the dot is in the maximum density drop ͑MDD͒ state for 3.42 T рBр3.70 T. As the ϭ1 configuration we have taken that corresponding to B 1 ϭ3.56 T, and have defined the other configurations as those corresponding to the value B ϭB 1 /.
The results of Ref. 2 for the Nϭ25 dot seem to indicate that a confining potential of parabolic type might be more adequate ͑see the discussion of the charge mode at the end of this section͒. Consequently, we have also studied it using a harmonic confining potential with 0 ϭ2.78 meV to reproduce the experimental dipole energy at Bϭ0. In this case the system is in the MDD state for 4.46 T рBр4.69 T. We have taken B 1 ϭ4.58 T and have defined the other configurations as indicated before.
The existence of a MDD state for this dot is in agreement with the findings of Ref. 24. It is also worth remarking that even for such a small dot, as defined above also coincides with the number of occupied (n,↑or ↓) bands for values up to ϭ5 -6, which correspond to rather low B intensities. In Fig. 4 we present the strength function corresponding to N ϭ25 with parabolic confinement. Figures 1-4 show that in both channels the response at Bϭ0 is concentrated within a small energy range, with one single peak or with several closely lying fragmented peaks which exhaust most of the f-sum rule. The peak energy is lower in the spin than in the density channel. This is due to the character of the residual interaction, which is attractive in the former channel and repulsive in the latter channel, and shifts the TDLSDA responses from the free-particle response in opposite directions. The residual interaction in the spin channel is weaker than in the density channel, where not only the exchange-correlation term but also the Coulomb direct term contributes. Consequently, the spin response is close in energy to the free response. It is thus difficult to distinguish the collective spin mode from the single-particle spectrum. In large dots it also causes a stronger Landau damping in the spin than in the density channel. These facts were observed and discussed in Ref. 13.
At Bϭ0, as a consequence of Kohn's theorem, if the con- fining potential is harmonic with frequency 0 the excitation energy of the dipole density mode is equal to 0 irrespective of the number of electrons. Otherwise, the excitation energy depends on N; see, for example, Refs. 5 and 14. In the spin channel Kohn's theorem does not hold, and a size dependence appears in the dipole spin mode even for parabolic confining potentials. 14 When B is not zero, the dipole mode in either channel splits into two branches, one with negative B dispersion and another with positive B dispersion. The splitting is due to the breaking of the l-degeneracy of the sp energies by the applied magnetic field. Several phenomena then appear. We first notice that for B values such that the spin of the dot gs is different from zero, the spin and density modes are coupled. This is particularly apparent in the Nϭ210 dot. Indeed, at Bϭ1.71 and 5.14 T the system is almost paramagnetic, having 2S z ϭ2 and 0, respectively ͑see Fig. 4 of Ref.
20͒. As a consequence, the modes are uncoupled, as it can be seen from panel ͑b͒ and especially from panel ͑d͒ of Fig. 3. In contradistinction, at Bϭ3.43 and 7 T we have 2S z ϭ54 and 74: the system has a large spin magnetization in the gs and the spin and density modes are clearly coupled, as displayed in panels ͑c͒ and ͑e͒ of that figure. One can see a distinct peak in the spin response at the energy of the density mode. This effect has been experimentally observed. 13 The strength of this peak increases with S z and when the system is fully polarized, which happens slightly above Bϭ10 T for the Nϭ210 dot, all the strength is transferred from the spin to the density channel. Conversely, the spin mode may be observed in the density channel with some intensity. This effect is hindered because Kohn's theorem prevents it from occurring for parabolic potentials, and for the disk potential it is of order (2S z /N) 2 .
The interplay between charge and spin modes is especially marked in the mixed channel, where the density response to the spin dipole operator D m , and the spin response to the density dipole operator D , are described. This is shown in Fig. 5 for the Nϭ210 dot at ϭ3. One clearly observes two peaks at the energy of the density modes, and another two at the energy of the spin modes. This can be understood casting the mixed response into a sum over ''spin dipole states'' ͉m͘ and another over ''charge dipole states'' ͉͘

͑16͒
For a disk confining potential, the matrix element ͗0͉D ͉m͘ is not zero and there is a contribution to S m from the spin modes. For a harmonic confining potential ͗0͉D ͉m͘ is zero, and only the density modes would contribute to S m through the -sum in Eq. ͑16͒.
The B dispersion of the main peaks of the spectrum is reported in Figs. 6-9. In these figures the density modes are represented by dots, and the spin modes by triangles. Solid symbols correspond to integer filling factor values, and the insets show the negative B dispersion branch of the spin mode. As a guide, we have drawn lines starting at the value of the Bϭ0 frequency and following the VPM B dispersion laws ͓see Eqs. ͑24͔͒.
Several features of these figures are worth discussing. Concerning the spin modes, we first see that at low B their energy is much smaller than the energy of the density modes, in agreement with the experimental findings. 13 At higher B the dot is eventually fully polarized and the longitudinal spin and density modes merge, as in the two-dimensional electron gas ͑see Fig. 5 of Ref. 27͒. This is not explicitly shown in the figures. Second, the negative B dispersion branch of the spin mode manifests a clear oscillatory behavior with , similar to that found for the density response, 37 also discussed in Ref. 20: the ''paramagnetic'' evenconfigurations have softer spin modes than the ''ferromagnetic'' oddconfigurations.
Our calculation predicts a spin instability occurring when the energy of the spin mode lower branch goes to zero at a critical B between ϭ1 and 2. This instability also manifests in the static spin polarizability Re͓␣ mm (Ϫ1) (0)͔, which becomes negative at these large B's. This indicates that the gs is no longer an energy minimum, and thus is not stable. It is worth recalling that no collective spin dipole modes are observed in the experiments at these rather high-B values. However, one cannot discard the idea that this might be due to the strong Landau damping existing in this energy region. We want to recall that in this B region, correlations not explicitly taken into account by TDLSDA might be important and could wash out that spin instability.
Finally, we would like to comment on the density dipole mode. For the parabolic confining potential, Figs. 6 and 9 show the well-known result that the density response yields the classical law represented by the first of Eqs. ͑24͒. For the disk confining potential, that law is also fairly obeyed, particularly by the negative B dispersion branch. We have systematically found that the positive B dispersion branch is fragmented, especially for the Nϭ210 dot. Comparing our calculations for the Nϭ25 dot with the experimental results, 2 we conclude that a parabolic potential is better suited than a disk potential to represent the physical situation. Conversely, the confining potential of the Nϭ210 dot is not parabolic, and this is the origin of the second upper branch found in the experiment. The nonharmonicity of the confining potential has been presented as the origin of that branch on the basis of a Hartree plus random-phase-approximation ͑RPA͒ calculation 5 in dots with Nр30 ͑see also Ref. 38͒. Our TDLSDA calculation supports that interpretation.
We can see that the density response of such a large dot displays two satellite branches instead of one branch ͑compare with the results for the small Nϭ25 dot͒, i.e., for TDLSDA seems to produce a high-frequency peak that is more fragmented than in the experiments. We remark that the upper branches disappear at intense B, in agreement with the experimental findings. This gives further support to our explanation and that of Ref. 5 about the origin of these branches. Indeed, at higher magnetic fields one expects that the harmonic cyclotron potential dominates over the other contributions in the KS equation.

IV. VIBRATING POTENTIAL MODEL
The intuitive idea behind the TDLSDA is that a small amplitude time-dependent variation of the mean field around the static equilibrium configuration produces an oscillation in the electronic density, which causes a small-amplitude collective motion of the system. This motion is self sustained if the induced density is precisely that needed to generate the oscillating potential. The vibrating potential model naturally arises when one considers the first iteration of the ͑perturbed mean field͒ ↔ ͑induced density͒ self-consistent scheme. In homogenous systems, where translational invariance determines the shape of the induced density, this yields the exact solution. In finite systems ͑nuclei, metal clusters, dots͒, the model provides a useful approximation.
Using the general method described in the Appendix, we consider the following VPM Hamiltonian: FIG. 9. Same as Fig. 7, using a parabolic confining potential with 0 ϭ2.78 meV instead of a disk confining potential.
with a time-oscillating potential ␦V͑r ជ ,t ͒ϭϪ where 0 2 ϭn s /R. The time dependence of ␦V(r ជ ,t) is in the ͗ • • • ͘ spatial foldings with the densities induced by a time-dependent external field; see the Appendix. To obtain the static part H 0 of the Hamiltonian ͓Eq. ͑17͔͒ we assume an exact cancellation between the Hartree and external potentials, and have taken a parabolic approximation for the exchange-correlation potential at Bϭ0. V is the exchangecorrelation constant introduced in the Appendix.
This VPM Hamiltonian can now be solved analytically within the RPA by finding the operators O ϩ solution of the equations of motion, We have used the methods illustrated in Ref. 36 to compute the commutators with the Hamiltonian as where H 0 is the static Hamiltonian, and ␦V(O) the variation arising from the induced densities. It can be shown that the solutions to Eq. ͑19͒ are given by The corresponding frequencies are and it is easy to verify that The states ͉ Ϯ ,m ͘ϵO Ϯ ,mϩ ͉0͘ are orthonormal, and carry an orbital angular momentum L 0 Ϯ1 and a spin S 0 , where L 0 and S 0 are the orbital and spin angular momenta of the ground state, respectively.
The charge dipole and spin dipole strengths are distributed among the above states as follows: It is a simple matter to check that the above solutions exhaust the sum rules equations ͑15͒. This vibrating potential model reproduces the gross features of the full self-consistent calculation. Its parabolic form guarantees that Kohn's theorem is fulfilled, and as a consequence the response to D is shared by just two peaks, which have the same strength N/4 . Accordingly, the spin dipole modes ͉ Ϯ m ͘ are not excited by the dipole operator D , and the corresponding matrix elements vanish. Another consequence of Kohn's theorem is that within the VPM only the density mode contributes to the mixed response. The model predicts that, in the mixed channel, the charge dipole modes ͉ Ϯ ͘ are excited with the same strength S z 2 / N by the spin dipole operator D m .
In the spin channel, the spin dipole operator D m excite both the charge ͉ Ϯ ͘ and spin ͉ Ϯ m ͘ modes. These peaks have strengths ͉͗0͉D m ͉ Ϯ m ͉͘ 2 ϭN͓1Ϫ(2S z /N) 2 ͔/4 and ͉͗0͉D m ͉ Ϯ m ͉͘ 2 ϭS z 2 / N, which are the same for both ⌬L z ϭϮ1 branches. Finally, when S z ϭ0, the density and spin modes are decoupled, and when the system is fully polarized, i.e., 2S z ϭN, all the strength is transferred to the density channel.
Besides this qualitative agreement, the strengths given by the first and last lines of Eq. ͑26͒ agree well with the result of the microscopic calculation. Also the ratio (2S z /N) 2 between the strength of the density peaks excited by D m and D is reproduced.
The second of Eqs. ͑24͒ can be used to determine the energy of the spin dipole mode Ϯ m , as the first one has been often used in the case of the charge dipole mode, if we fix the parameters entering that equation. We first take xc 2 ϭ2n s /R 2 . This estimate is obtained by identifying the kinetic energy per particle in the exchange-correlation harmonic-oscillator potential with that of the 2DEG n s /2, and approximating the mean-square radius of the dot ͗r 2 ͘ by R 2 /2ϭr s 2 N/2, where r s 2 ϭ1/(n s ). One then obtains The value of r s 2 V is related to the spin susceptibility of the two-dimensional electron gas 0 /ϭr s 2 V ϩ1 ͑Ref. 31͒. Equation ͑27͒ yields values for Ϯ m which agree with the TDLSDA ones. In particular, we have checked that the N dependence of the Bϭ0 calculations reported in Ref. 14 is well reproduced taking r s 2 V ӍϪ0.3 at r s Ӎ1. Note that the N and r s dependence of energy of the spin mode is different from that of the density mode, which is given by when we take for the frequency 0 of the confining potential the estimate 0 ϭr s Ϫ3/2 N Ϫ1/4 obtained from the disk potential in the r/RӶ1 limit.
In Fig. 10 we report the energies given by Eq. ͑27͒ for a dot of Nϭ200 together with the experimental spin dipole mode. 13 A value of r s ϭ0.65 has been used which yields r s 2 V ϭϪ0. 24. From the figure one sees that the gross trends of the experimental B dispersions are reproduced by the simple model.
The VPM works better for large dots than for small dots. This is not surprising in view of the approximations leading to the VPM. It can be seen by comparing with the TDLSDA results in Figs. 6-9 after being advised that, for the Nϭ5 dot, the dotted line representing ͉ ϩ m ͘ does not pass near the main peak energies of the positive B dispersion branch, but near to the satellite peak energies. Equation ͑27͒ shows that, within the VPM, the spin instability at Bϭ0 mentioned in Sec. III occurs when r s 2 V ϩ1/2р0. In the 2DEG this happens 31 when r s 2 V ϩ1р0. The difference between these conditions is due to finite-size effects, which are crucial in determining the values of r s at which the instability appears in quantum dots. Whereas in the bulk the spin instability sets in at r s Ӎ37, Eq. ͑27͒ shows that in dots it happens at r s ϳ3. This value is well within the range of those found in Ref. 32 for the onset of a spindensity-wave instability in small magic dots, and a factor of 2 larger than the typical values obtained in that reference for open shell dots.

V. SUMMARY AND OUTLOOK
In this work we have thoroughly discussed the longitudinal dipole response of quantum dots. We have shown that the TDLSDA is able to reproduce the main features of the experimental results. In particular, we have found that the density and spin modes are clearly coupled in the spin channel if the system is partially polarized, and that the frequency of the spin edge magnetoplasmon presents an oscillatory behavior as a function of , being especially soft for even filling factors.
Our numerical scheme has allowed us to study large-sized dots whose spectrum has been experimentally examined in detail, instead of relying on extrapolation of the results obtained for small-sized dots. This is crucial to identifying the behavior of physical quantities like excitation frequencies.
TDLSDA can be easily applied to other multipole spin and density excitations. This is relevant in view that recent experiments have been able to identify the spin monopole and quadrupole modes, and likely charge modes different from dipole. 13 Work to extend the present study to other multipolarities, and to study the spin transverse channel in large dots along the line described in Ref. 39, is underway.

ACKNOWLEDGMENTS
This work was performed under Grant Nos. PB95-1249 and PB95-0492 from CICYT, Spain, and GRQ94-1022 from Generalitat of Catalunya. A.E. acknowledges support from the Dirección General de Enseñanza Superior ͑Spain͒.

APPENDIX
In this appendix we provide a general derivation of the VPM Hamiltonian leading to Eq. ͑17͒. We start from the density and magnetization variations ͓␦͑ r ជ ,t ͒,␦m͑ r ជ ,t ͔͒ϭ␣͑ t ٌ͓͒ x 0 ͑ r ជ ͒,ٌ x m 0 ͑ r ជ ͔͒ ϩ␤͑t ٌ͓͒ x m 0 ͑ r ជ ͒,ٌ x 0 ͑ r ជ ͔͒, ͑A1͒ obtained from the static ground-state density 0 and spin magnetization m 0 through the unitary transformations e ␣(t)͚ i ٌ x i and e ␤(t)͚ i ٌ x i z i , respectively, in the limit of small deformations. In Eq. ͑A1͒, ␣(t) and ␤(t) give the amplitude of the oscillations.
Assuming the forms V xc ϭV xc () and W xc ϭmV / 0 for the exchange-correlation potentials, which amounts to expanding E(,m) around mϭ0 up to m 2 order and identifying V with 0 ‫ץ‬ 2 E xc ( 0 ,m)/‫ץ‬m 2 ͉ mϭ0 , one obtains the variation in the one body potential of Eq. ͑1͒ induced by the density variations ͑A1͒: where with V 0 ϭV xc (ϭ 0 ). We have further assumed that m 0 ϭ(2S z /N) 0 in the gs. Using the results it is then possible to write the variations in the one-body potential of Eq. ͑1͒ in the separable form

͑A5͒
Following Ref. 30, one could now express the various responses to an oscillating potential of form ͑A5͒ in terms of the independent particle response functions , Ј 0 through RPA-type equations. However, our aim here is to develop an analytical model which allows us to understand in the numerical results of Sec. III a simple way. To this end, in Eqs. ͑A3͒ we take a harmonic oscillator for the one-particle potential V 0 ϭ 1 2 xc 2 r 2 to simulate the short-range effects, and a step function for the electronic density 0 entering the long-range contribution ٌ͐ x 0 (r ជ )/͉r ជ Ϫr ជ Ј͉dr ជ Ј. We also assume Q (r ជ i )ӍϪ2V x i /͗r 2 ͘.
Equation ͑17͒ of Sec. IV is then obtained.