text:
# Asymmetric Bose-Einstein effect

## 1 Introductory remarks

Recently one observes a renewal of interest in analysing the space-time structure of sources in multiparticle production by means of Bose-Einstein (BE) interference which followed the example of astrophysical investigations. The main motivation of this renewal was the controversy concerning the $`e^+e^{}W^+W^{}`$ process. The existence of interference effects between strings from two $`W`$-s is still debatable. Another reason to analyze the BE effect were the efforts to estimate the shape and the lifetime of the source of particle production (e.g., for the LEP data at the $`Z^0`$ peak).

Investigating such subtle effects became possible when one started to model this effect in Monte Carlo generators. There are several methods of modelling: using the ”afterburner” for which the original MC provides a source,, shifting the momenta or adding weights to generated events,. Another approach was set forward by Andersson and collaborators who rejected the assumption of incoherence, basic for the ”standard” models, and used the symmetrization inside a fragmenting string to model the effect for a single string. In this talk111Presented by K. Fiałkowski we consider the most widely used methods of shifting momenta and weighting events.

In a recent paper we compared the 3-dimensional data for BE effect from LEP with the results of the standard momentum shifting procedure and of the weight method. In both cases differences between the longitudinal and transverse radii appear, but the results disagree with data. Introducing asymmetric weight factors one may describe the data satisfactorily. This is presented in the second section.

In the last section we return to the problem of ”inter-$`W`$” interference. We point out that the BE effect for pairs of pions from different $`W`$-s is in principle unavoidable, and the only problem is a reliable estimate of its magnitude. We mention some preliminary results for the obtained effect using a modified weight method applicable in this case.

## 2 Asymmetry in $`Z^0`$ decays

We discuss only the L3 data which measured the BE ratios using ”uncorrelated background”

$$R_2(p_1,p_2)=\frac{\rho _2}{\rho _2^{mix}}/\frac{\rho _2^{MC}}{\rho _2^{mix,MC}}$$

and three different radii to parametrize the data. The DELPHI data are parametrized with only two radii, and the OPAL data use the like/unlike ratio which requires cutting off the resonance affected regions even in double ratios.

We refer to our paper for the choice of variables and definitions of parameters. In the data the fitted value of the parameter $`\lambda `$ is $`0.41\pm 0.01`$, and the values of radii (in $`fm`$) are:

$$R_L=0.74\pm 0.02_{0.03}^{+0.04},R_{out}=0.53\pm 0.02_{0.06}^{+0.05},R_{side}=0.59\pm 0.01_{0.13}^{+0.03}$$

As shown in the L3 paper, the standard LUBOEI procedure built into the JETSET Monte Carlo generator gives

$$R_L=0.71\pm 0.01,R_{out}=0.58\pm 0.01,R_{side}=0.75\pm 0.01.$$

We confirmed these numbers in our calculations and checked how the results depend on the parameters $`\lambda _{in}`$ and $`R_{in}`$ assumed in the LUBOEI input function. In all cases we get (contrary to the data) $`R_{side}>R_L>R_{out}`$, although the input function was obviously symmetric. No choice of input parameters gives the values of $`R_i`$ compatible with data.

The known problems of LUBOEI procedure,, led to a revival of weight methods, known for quite a long time, but plagued also with many practical problems, some of which have been recently solved. In this method we may repeat the same calculations as done for the LUBOEI procedure. The major features of the results are surprisingly similar: with weight factors depending only on $`Q^2`$ we get different values of fitted $`R_i`$ parameters. Moreover, the hierarchy of parameters is the same: $`R_{side}>R_L>R_{out}`$. This suggests that the assymetry is generated by the jet-like structure of final states and not by any specific features of the procedure modelling the BE effect. Again, no choice of the input parameters allows to describe the data.

One may get more information on the problem of asymmetric BE effect in MC generators using the asymmetric weight method, i.e. introducing weight factors which depend in a different way on $`Q_L=|p_{1L}p_{2L}|`$, $`Q_{side}=|p_{1side}p_{2side}|`$ and $`Q_{out}=|p_{1out}p_{2out}|`$, where the indices denote the momentum components defined in the usual way

$$w_2(Q_L,Q_{out},Q_{side})=\mathrm{exp}([Q_L^2(R_L^{in})^2Q_{out}^2(R_{out}^{in})^2Q_{side}^2(R_{side}^{in})^2]/2)$$

This weight factor reduces approximately to the symmetric weight factor

$$w_2(p_1,p_2)=\mathrm{exp}[(p_1p_2)^2R_{in}^2/2]$$

when $`R_L^{in}=2R_{out}^{in}=R_{side}^{in}=R_{in}`$. The weight attached to the event is given by

$$W(p_1,\mathrm{}p_n)=\underset{P}{}\underset{i}{}w_2(p_ip_{P(i)}).$$

Since for the symmetric weights the resulting fitted values of $`R_{side}`$ are bigger than the values of $`R_L`$ (contrary to the inequality seen in the data), it seemed natural to take the input value of $`R_{side}^{in}`$ smaller than $`R_L^{in}`$. The best set we found is

$$R_L^{in}=0.9fm,R_{out}^{in}=0.3fm,R_{side}^{in}=0.4fm.$$

Then we get the following fitted values of the parameters

$$R_L=0.73fm,R_{out}=0.54fm,R_{side}=0.65fm$$

in agreement with data.

There is a striking difference between the input values of the radii assumed in the weight factors and the resulting best fit values from the double ratio calculated with these weights. Although the hierarchy $`R_L>R_{side}>R_{out}`$ is the same in both cases, the fitted values differ by less than 25%, whereas there is a difference by more than a factor of two between the input values.

Moreover, further decrease of the values of $`R_{out}^{in}`$ and $`R_{side}^{in}`$ hardly affects the resulting double ratio and fitted values of $`R_i`$. This seems to be the inherent property of the JETSET generator, which yields a rather strong suppression of large values of $`Q_i`$ and $`Q^2`$ even without any procedure imitating the BE effect. Apparently this suppression dominates over the weak enhancement of low values of $`Q_i`$ induced by the weight factors with small values of $`R_i`$. For small $`R_i^{in}`$ there is no simple correspondence between the input and output values of radii. This looks analoguous to the effect noted already for a symmetric BE effect described by the LUBOEI procedure. Therefore any direct interpretation of the fit values for BE double ratios in terms of the different radii of the asymmetric source is a rather delicate matter.

## 3 Interference between pions from different W-s

There is no experimental evidence of interference effects between pions from different $`W`$-s. Moreover, some data are shown together with ”MC predictions including inter-W BE effect” which clearly disagree with data. This seems to suggest that the ”inter-W BE effect” does not exist and that the models where BE effect is confined to a single string, are preferred.

We want to stress that this suggestion is false. The BE effect for two pions coming from (incoherent) decays of two $`W`$-s is an inevitable consequence of quantum mechanics and any model neglecting it is not complete. However, for given kinematic configuration this effect may be negligible. The ”MC predictions” for this effect shown with the data are completely unrealistic, as they assume the same effect for pairs of pions coming from the same, and from different $`W`$-s. There is no reason to believe that this is true. We are now working on the realistic estimate of this effect for the weight method. For the decay of two $`W`$-s the assumption of no correlations between the creation point and momentum is obviously wrong. The simplest necessary modification of the standard formula for weight is to take into account the different space location of two sources. We get then

$$W(p_1,\mathrm{}p_n)=\underset{P}{}\underset{i}{}w_2(p_ip_{P(i)})\underset{k}{}w_2(p_kp_{P(k)})W^{}(\mathrm{\Delta }p)$$

where $`\mathrm{\Delta }p=[_i(p_ip_{P(i)})_k(p_kp_{P(k)})]/2.`$ Here two products (sums) extend over particles coming from two $`W`$-s and $`W^{}`$ is the Fourier transform of a distribution of distance between two sources ($`W`$-s). This prescription for weights may lead to results significantly different from the standard ones.

Preliminary results with realistic $`W^{}`$ suggest that even at threshold of the $`WW`$ production process (where the spatial distance between two sources is minimal) the effect for pions from different $`W`$-s is weak, and at energies dominating present data sample it is completely negligible.

Summarizing, we have shown that the weight method of implementing the BE effect into MC generators is quite useful for various applications. It may describe the effects for which other methods fail.

## Acknowledgments

One of us (KF) thanks for the support of FNP (subsydium FNP 1/99 granted to A. Białas).


format: text

text:
# The transverse breathing mode of an elongated Bose-Einstein condensate

(Received November 23, 2001)

## Abstract

We study experimentally the transverse monopole mode of an elongated rubidium condensate. Due to the scaling invariance of the non-linear Schrödinger (Gross-Pitaevski) equation, the oscillation is monochromatic and sinusoidal at short times, even under strong excitation. For ultra-low temperatures, the quality factor $`Q=\omega _0/\gamma _0`$ can exceed 2000, where $`\omega _0`$ and $`\gamma _0`$ are the mode angular frequency and damping rate. This value is much larger than any previously reported for other eigenmodes of a condensate. We also present the temperature variation of $`\omega _0`$ and $`\gamma _0`$.

Collective excitations constitute one of the main sources of information for understanding the physics of many-body systems. Since the discovery of Bose-Einstein condensation (BEC) in weakly interacting bosonic gases Varenna98  much theoretical and experimental effort has been devoted to the study of the low energy excitations of these systems Stringari96 ; Jin96 ; Mewes96 ; Onofrio00 ; Marago00 . At ultra-low temperatures, the measured frequencies of these excitations are found to be in very good agreement with the ones derived from a simple linearization of the Gross-Pitaevskii equation, which describes the evolution of the condensate wave function.

The temperature dependence of the damping and frequency of these oscillations has attracted strong attention of both experimental Jin97 ; Stamper98 ; Marago01  and theoretical groups Dalfovo99 ; Guilleumas99 ; Rusch99 ; AlKhawaja00 ; Giorgini00 ; Reidl00 ; Morgan00 ; Kagan01 ; Jackson01 ; Williams01 ; Jackson01b . Indeed, a dilute BEC is a very convenient system for the study of dissipation in a macroscopic quantum object. At present, there still exist significant discrepancies between some available experimental results and the corresponding theoretical predictions (see e.g. Reidl00 ; Jackson01b ). Therefore new measurements on simple and well identified modes are clearly needed to test the various theoretical approaches.

It has been known for some time that in an isotropic trap, the monopole (or breathing) oscillation mode may exhibit very unusual features. For instance, Boltzmann noticed that, within the framework of the equation which bears his name, the breathing mode of a set of particles confined in a 3D isotropic harmonic trap is undamped, irrespective of the ratio between the average collision rate and the trap frequency Cercignani ; Guery99a ; Ghosh00 . More recently it was pointed out that in a 2D isotropic harmonic trap, the breathing mode of identical bosons interacting via a contact potential is also undamped Kagan96 ; Pitaevskii97 . This property, which holds at any temperature, can be related to the existence of a hidden symmetry of the problem described by the two-dimension Lorentz group SO(2,1) Pitaevskii97 .

The purpose of the present work is to investigate experimentally the transverse breathing mode of an elongated condensate. We find that this mode exhibits very unique features, such as an extremely large quality factor $`Q=\omega _0/\gamma _02000`$ for temperatures much lower than the condensation temperature. Here $`\omega _0`$ and $`\gamma _0`$ are the angular frequency and damping rate, respectively. This quality factor, much larger than any previously reported value for an oscillation mode of a gaseous BEC, should provide a sensitive test of the various theoretical models of dissipation at ultra-low temperature. Another very unusual property lies in the fact that the measured frequency $`\omega _0`$ is nearly independent of temperature. Furthermore, we find that the transverse area of the condensate undergoes a pure, isochronal oscillation even for strong excitation. This striking behavior emerges from a scaling invariance of the time-dependent Gross-Pitaevski equation.

Our 87Rb condensate is formed by radio-frequency (rf) evaporation of $`10^9`$ atoms in a Ioffe-Pritchard magnetic trap. The initial temperature of the cloud pre-cooled using optical molasses is 100 $`\mu `$K. We vary the final temperature of the gas and the number of atoms in the condensate by adjusting the final rf, $`\nu _\mathrm{f}`$. The condensation threshold is reached at $`T_c=290`$ nK, with $`N_c=\mathrm{8.4\hspace{0.33em}10}^5`$ atoms calibration . The coldest sample reliably reproducible ($`T40`$ nK) is a quasi-pure condensate with $`N_0=10^5`$ atoms. It is obtained using $`\nu _\mathrm{f}=\nu _0+6`$ kHz, where $`\nu _0`$ is the rf at which the trap is emptied.

The trap has a longitudinal frequency $`\omega _z/(2\pi )=11.8(1)`$ Hz for the atoms prepared in the stretched state $`m=2`$ of the $`5S_{1/2},F=2`$ ground state. The transverse magnetic gradient is $`b_{}^{}=1.33`$ T m-1 and the transverse frequency $`\omega _{}`$ is adjusted by varying the bias field $`B_0`$ at the center of the trap. We choose $`B_0=\mathrm{0.844\hspace{0.33em}10}^4`$ T, which leads to $`\omega _{}/(2\pi )=182.6(4)`$ Hz anisotropy . Both frequencies $`\omega _z`$ and $`\omega _{}`$ are measured by monitoring the oscillation of the center-of-mass of the atom cloud.

Once the condensate is formed, we excite the transverse breathing mode by modifying the bias field from $`B_0`$ to $`B_0^{}`$, which changes the frequency $`\omega _{}`$ to $`\omega _{}^{}>\omega _{}`$. After a time $`\tau 2\pi /\omega _{}^{}`$, we set the bias field back to its original value. We let the cloud oscillate freely in the magnetic trap for an adjustable time $`t`$ and then measure the transverse density profile of the condensate after a period of free expansion. In this pursuit, we suddenly switch off the magnetic field, allow for a 25 ms free-fall, and image the absorption of a resonant laser by the expanded cloud. The imaging beam propagates along the $`z`$ axis. We fit the density profile of the condensate assuming the parabolic shape of the Thomas-Fermi distribution and extract the radii $`R_x`$ and $`R_y`$ in the plane $`z=0`$. The temperature is obtained from a Gaussian fit of the uncondensed part of the distribution.

Within the Thomas-Fermi approximation, the evolution of the trapped condensate after the excitation of the transverse breathing mode is well described by a time-dependent scaling transformation Kagan96 ; Castin96 ; Kagan97 . The spatial density $`\rho (𝐫,t)`$ at time $`t`$ is deduced from the spatial density at time $`0`$ (i.e. just after the transverse frequency is set back to its value $`\omega _{}`$) by:

$$\rho (x,y,z,t)=\frac{1}{\lambda _{}^2\lambda _z}\rho (\frac{x}{\lambda _{}},\frac{y}{\lambda _{}},\frac{z}{\lambda _z},0).$$

The evolution of the scaling parameters $`\lambda _{}`$ and $`\lambda _z`$ is given by Castin96 ; Kagan97 :

$$\ddot{\lambda }_{}=\frac{\omega _{}^2}{\lambda _{}^3\lambda _z}\omega _{}^2\lambda _{}\ddot{\lambda }_z=\frac{\omega _z^2}{\lambda _{}^2\lambda _z^2}\omega _z^2\lambda _z.$$
(1)
We start our analysis with the evolution for relatively short times $`t2\pi /\omega _z`$. In this case one can neglect the variation of $`\lambda _z`$ and integrate the equation for $`\lambda _{}`$:

$$\lambda _{}(t)=\left(\mathrm{cosh}(\alpha )+\mathrm{sinh}(\alpha )\mathrm{sin}(2\omega _{}t\beta )\right)^{1/2},$$

where the parameters $`\alpha `$ and $`\beta `$ depend on the initial conditions. One remarkable consequence of this solution is that the quantity $`A(t)=R_x^2(t)+R_y^2(t)=\lambda _{}^2(t)A(0)`$ always undergoes a sinusoidal oscillation at frequency $`2\omega _{}`$, irrespective of the strength of the excitation TOF . This is confirmed by the results of fig. 1, where we observe a sinusoidal oscillation of $`A(t)`$ with a factor 3 between the maxima and the minima. The frequency of the oscillation is indeed found to be equal to $`2\omega _{}`$, within our experimental uncertainty ($`\pm 0.5`$ %).

At longer times one has, in principle, to take into account the variation of $`\lambda _z`$. We restrict our analysis to small amplitude oscillations and write $`\lambda _i=1+ϵ_i`$, where $`ϵ_{},ϵ_z1`$. The linearized equations of motion

$$\ddot{ϵ}_{}+\omega _{}^2\left(4ϵ_{}+ϵ_z\right)=0\ddot{ϵ}_z+\omega _z^2\left(2ϵ_{}+3ϵ_z\right)=0$$

give rise to two modes Stringari96 ; Mewes96 . One is a fast mode of frequency $`2\omega _{}`$ corresponding to $`ϵ_{}/ϵ_z=\omega _{}/\omega _z`$. The other one is a slow mode of frequency $`\sqrt{5/2}\omega _z`$ with $`ϵ_{}/ϵ_z=1/4`$. Our excitation scheme corresponds to $`\dot{ϵ}_{}0`$ at time 0 while the three other quantities $`ϵ_{}`$, $`ϵ_z`$ and $`\dot{ϵ}_z`$ are approximately zero at initial time. Under these conditions, we mainly excite the fast mode $`2\omega _{}`$ since the relative weight of the slow mode in the evolution of $`\lambda _{}`$ is only $`\omega _z/(\sqrt{40}\omega _{})<0.01`$. For this fast mode, the condensate dynamics consists essentially in a transverse monopole oscillation since $`ϵ_zϵ_{}`$.

We investigate the oscillation and the damping of this mode on a long time scale (up to $`0.7`$s). An example is given in fig. 2 for $`T=40(\pm 20)`$ nK (note that for such cold clouds, corresponding to a quasi-pure condensate, the temperature can only be inferred from the final rf used for the evaporation). The evolution of $`A(t)`$ is well fitted by a damped sinusoidal function, $`A_0+\delta A_0\mathrm{cos}(\omega _0t+\varphi )e^{\gamma _0t}`$, with $`\delta A_0/A_0=0.063(4)`$, $`\omega _0/\mathrm{\hspace{0.17em}2}\pi =366.3(5)`$ Hz and $`\gamma _0=1.2(2)`$ s-1. This corresponds to a quality factor $`Q=\omega _0/\gamma _02000`$ much larger than any previously reported for other eigenmodes of a BEC. For instance, measurements performed with a TOP trap ($`\omega _z=2\sqrt{2}\omega _{}`$) gave $`Q200`$ for the lowest $`m=0`$ mode Jin97 . For a cigar shape trap, measurements of the low frequency mode discussed above ($`\omega =\sqrt{5/2}\omega _z`$) led to $`Q80`$ Stamper98 . A smaller value ($`Q25`$) was measured for the scissors mode Marago01 . The present quality factor is also one order of magnitude larger than what we find for the transverse quadrupole mode under the same experimental conditions Chevy00 .

We measure the oscillation frequency and the damping rate over a large temperature range compatible with the detection of a condensate. The results are given in fig. 3, together with the total number of atoms $`N`$. The measurement of $`\gamma _0`$ is made from the decay of the oscillation amplitude between $`0`$ and $`200`$ ms. For the points of lowest and highest temperature in figs. 3b and 3c, we calculate, including mean field and finite size corrections Dalfovo99 , $`T/T_c=0.3`$ and $`0.8`$, respectively.

The frequency $`\omega _0`$ is found to be very close to the Thomas-Fermi prediction $`2\omega _{}`$, and it varies only slightly with temperature (fig. 3b). In particular, we do not observe the 10 % relative decrease predicted in Giorgini00  when $`T/T_c`$ increases from 0 to 0.8.

The damping rate $`\gamma _0`$ varies approximately linearly with temperature in the range $`50200`$ nK. The linear variation of $`\gamma _0`$ for $`k_\mathrm{B}T`$ larger than $`\mathrm{}\omega _0`$ ($`20`$ nK) and than the chemical potential $`\mu `$ ($`60`$ nK for the coldest point in fig. 3) is expected if one assumes either Beliaev-type or Landau-type damping of the monopole excitation (see e.g. Fedichev98 ). In the Beliaev mechanism, the monopole excitation decays into two excitations with a lower energy. For Landau damping, it annihilates in a collision with a thermal excitation to give rise to another excitation. The calculation of the corresponding damping rate is quite complex for inhomogeneous systems and can lead to drastically different results for modes with similar energies, depending on the stochastic or regular nature of the relevant thermal excitations Fedichev98 . In this respect it is worth noting that the frequency of this transverse breathing mode is an integer multiple of the trap oscillation frequency, which is usually not the case for other modes of the condensate.

For a very low temperature ($`T40`$ nK) and an increasing excitation strength, we observe a significant deviation of the long term behavior of $`A(t)`$ with respect to the simple damped sinusoid shown above. Fig. 4 shows the oscillation amplitude for the weak excitation ($``$, $`\tau =75\mu `$s) together with two sets of data obtained for $`\tau _1=300\mu `$s ($`\mathrm{}`$) and $`\tau _2=450\mu `$s ($`\mathrm{}`$). One clearly sees that for strong excitation the amplitude of the oscillations passes through a local minimum around $`t=300`$ ms before increasing again, and decaying finally after 700 ms. We fit our experimental data with the phenomenological function $`A(t)=A_0+\delta A(t)\mathrm{cos}(\omega _0t+\varphi )`$ where the amplitude $`\delta A(t)`$ is now slowly modulated in time, in addition to the exponential decay already discussed: $`\delta A(t)=\delta A_0\left(1\xi \mathrm{sin}^2(\delta \omega t)\right)e^{\gamma t}`$. The best fit parameters are for these two cases $`\gamma _1=1.48(15)`$ s-1, $`\delta \omega _1/(2\pi )=1.08(4)`$ Hz, $`\xi _1=0.5(1)`$, and $`\gamma _2=1.35(15)`$ s-1, $`\delta \omega _2/(2\pi )=0.96(4)`$ Hz, $`\xi _2=0.6(1)`$.

Two scenarios can be invoked to explain this evolution of $`A(t)`$, which is not predicted by (1). A first explanation could consist in assuming that we actually excite two modes whose frequencies $`\omega _0`$ differ by a small amount and that we observe the beating between these modes. However, if we assume that these two modes are excited linearly, the relative amplitude of the beating should not depend on the excitation strength, contrary to the observed. Furthermore, this explanation is not confirmed by a linear analysis of the time dependent Gross–Pitaevskii equation for our experimental conditions CastinPrikopenko . A second, more plausible, scenario involves a non-linear process, in which the monopole excitation may be converted in a reversible way into a sum of other excitations of frequency $`\omega _j`$. This process is efficient if it is resonant: $`\omega _0_j\omega _j`$. In this context, the modulation in fig. 4 can be interpreted as a parametric oscillation.

As a possible candidate for such a parametric (or Beliaev type) oscillation Hodby01 , one may consider the excitation of the pair of modes $`l=2,m=\pm 1`$ Stringari96 ; DGO99 ; Marago00 . The frequency of each mode is $`(\omega _{}^2+\omega _z^2)^{1/2}`$ so that the resonance condition given above is indeed nearly fulfilled. By imaging the condensate along the $`y`$ direction, we observe a weak excitation of this mode. However, no increase in its amplitude is found around $`t=`$300 ms when the amplitude of the breathing mode collapses (cf. fig. 4).

Alternatively one may consider the excitations consisting in a pair of phonon-like excitations propagating along $`z`$ with momentum $`k`$ and $`k`$, each with an energy $`\mathrm{}\omega _0/2\mathrm{}\omega _{}`$ Kagan01 . To test this latter hypothesis, we perform the same experiment for various values of $`\omega _{}`$, so that it spans the whole interval between two successive phonon modes (typically 6 Hz in this energy domain). However, expecting a resonance when $`\omega _0/2`$ coincides with the frequency of a phonon mode, we do not observe any significant change in the evolution of $`A(t)`$.

In any case we note that for the long oscillation times involved here, small non-linear terms in the confining potential (varying as $`x^2z`$ and $`y^2z`$) may play a significant role. Further work investigating this possibility is under progress in our laboratory.

In summary we have presented in this letter a study of the transverse breathing mode of an elongated Bose-Einstein condensate. For low temperatures and short evolution time, the observed isochronal and sinusoidal oscillation is in excellent agreement with the prediction based upon the scaling invariance of the Gross-Pitaevski equation. However, the temperature dependence of the frequency and damping of the oscillation reveals, on the contrary, some unexpected features. First, the observed damping is comparatively ten times smaller than what was previously reported for other modes of a BEC (quality factor of 2000 instead of 100-200). Second, the frequency of the breathing mode is found to be quasi-independent of temperature, in contradiction with recent theoretical predictions. Finally the long time evolution of this mode when strongly excited reveals an unexpected non-linear mixing phenomenon. This could be explored further using the formalism of Castin97 , by looking at a possible dynamical instability of the Gross-Pitaevskii equation in this regime of strong excitation.

###### Acknowledgements.

We thank Y. Castin, Yu. Kagan, L. Pitaevski, and G. Shlyapnikov for useful discussions. K. M. and P. R. acknowledge support by DEPHY and Alexander von Humboldt-Stiftung, respectively. This work was partially supported by CNRS, Collège de France and DRED.


text:
# The Ionized Stellar Wind in Vela X-1 During Eclipse

## 1 Introduction

Massive early type stars show substantial mass loss through a strong wind. Matter is radially accelerated via momentum transfer from the stellar UV radiation field through absorption and scattering in resonance lines of abundant ions (Castor, Abott, and Klein 1975). In general the stellar wind itself is a source of X-rays, which form through reverse shocks generated by high density material in the wind (Lucy and White 1980, Owocki et al. 1988, Feldmeier et al. 1997). Typical observed X-ray luminosities range from a few times 1031 erg s-1 for early type B-stars to a few times 1032 erg s-1 for massive O-stars (Chlebowski et. al 1989, Berghöfer et al. 1995). Recent Chandra and XMM-Newton observations resolved the X-ray emission of several early type stars confirming the existence of collisionally ionized plasmas, but also suggested the possibility that alternative emission mechanisms for the observed X-ray emission need to be investigated. (Schulz et al. 2000, Kahn et al. 2001, Waldron and Cassinelli 2001).

The presence of a neutron star companion leads to a typical Bondi-Hoyle accretion flow (Bondi and Hoyle 1944), where a fraction of the stellar wind gets swept up by the neutron star, which generates a strong X-ray continuum. The X-ray luminosity of such a source is many orders of magnitudes larger than the one expected from shocks in the stellar wind, and X-rays from that source are reprocessed in the wind resulting in discrete emission lines (Stevens and Kallman 1990).

A prototype of such a system is the high-mass X-ray binary (HMXB) Vela X-1 (4U 0900-40). Its early type primary star HD 77581 is a massive B 0.5 Ib type supergiant (Brucato and Kristian 1972) with a mass and radius of 34 $`R_{}`$ and 23 $`M_{}`$ (from Nagase 1989), and an inferred mass loss rate of $`>10^7`$ $`M_{}`$ yr-1 (Sato et al. 1986, Sako et al. 1999). The terminal wind velocity has been determined to be 1700 km s-1 (Dupree et al. 1980). The neutron star exhibits 282 s X-ray pulsations (Forman et al. 1973) and orbits the center of mass of the system at a distance of only about 0.6 stellar radii from the surface of the supergiant (a = 53.4 $`R_{}`$ ) with a period of 8.96 days. This means that the neutron star is deeply embedded within the acceleration zone of the stellar wind.

The X-ray continuum of Vela X-1 has been studied extensively in X-rays. The general shape can be explained by a power law with a high energy cutoff typical for X-ray pulsars (White, Swank, and Holt 1983). Its average X-ray luminosity is $`10^{36}`$ erg s-1 , which is consistent with accretion from a stellar wind without Roche lobe overflow (Conti 1978).

In this paper we are primarily concerned with X-ray line emission from Vela X-1 during the middle of its eclipse. A first analysis of a moderate resolution X-ray line spectrum of Vela X-1 during eclipse was performed on a set observations with ASCA , which covered an entire eclipse transition from ingress to egress (Nagase et al. 1994). The observations revealed lines from He- and H-like ions species. It was suggested that the He-like transitions are formed by cascades following radiative recombination. Sako et al. 1999 identified recombination lines and radiative recombination continua produced by photoionization in an extended stellar wind. The deduced differential emission measure distribution was consistent with the wind model by Hatchett $`\&`$ McCray (1977). Sako et al. also indictated the existence of fluorescence lines from ion species of lower Z than Fe.

## 2 Chandra Observations and Data Reduction

Vela X-1 was observed with the HETGS on 2000 April 13th (09:57:52 UT) continuously for 28 ks. The HETGS carries 2 different types of transmission gratings, the Medium Energy Gratings (MEG) and the High Energy Gratings (HEG). It allows for high resolution spectroscopy between about 1 and 35 Å with a peak spectral resolution at 12 Å of $`\lambda /\mathrm{\Delta }\lambda 1400`$ and at 1.8 Å of $`\lambda /\mathrm{\Delta }\lambda 180`$ in 1st order HEG. The dispersed spectra were recorded with an array of 6 charged coupled devices (CCDs). We refer to the available Chandra X-ray Center (CXC) docments for more detailed descriptions of the spectroscopic instruments111http://asc.harvard.edu/udocs/docs/docs.html.

We recorded a total of 3200 events in the co-added 1st order MEG and 1700 events in the HEG after standard grade selection. The CXC provided aspect corrected event lists via standard pipeline processing. These data were reprocessed using the latest available data processing products. The wavelength scale has a zero point accuracy of 0.002 Å in MEG and 0.001 Å in HEG 1st order. The current status of the overall wavelength calibration is better than 0.05$`\%`$. For the line analysis of the high resolution spectra we divide first by the instrument effective area and exposure, which is legitimate, because the spectral bin size is larger than the actual response of the instrument.

The image in the 0th order provides a source position as well as a medium resolution X-ray spectrum. The position of HD 77581 is well known to a precision of about 0.3 arcsec (Roeser et al. 1991). After correcting for a 2” aspect drift we find the position of the X-ray source to be within within 0.5” of this position. We also find no other X-ray sources within a radius of 30” around Vela X-1 with a flux above 1.2$`\times 10^{15}`$erg cm-2 s-1 .

The CCD spectrum of that image also allows us to draw comparisons with previous ASCA observations, which are of similar spectral resolution, but cover a larger fraction of the eclipse. This CCD spectrum is only moderately piled up ( 7$`\%`$) and can be used to determine qualitative characteristics of the observed emission as well as the total observed luminosity. The continuum model is the same as used by Sako et al. (1999), which consists of two absorbed power laws from the scattered component and a direct component with indices fixed to 1.7. The choice of this model is only motivated by the fact that we can draw comparisons with specifically the ASCA results. Since the spectrum shows significant line emission and we added several gaussian functions to fit the brighter lines observed in the ASCA spectrum (see table 2 in Sako et al. 1999).

In general we find good agreement with the ASCA results. The most striking difference is the column density of the direct component, which is twice as high and also indicates a much fainter direct component. This is a consequence of the fact that the Chandra observations exposure covers a smaller range in orbital phase around phase zero. We verfied this by re-analysing the ASCA data using the Chandra phase interval. From an X-ray flux of 0.88$`\times 10^{11}`$erg cm-2 s-1 we then deduce 2.1$`\times 10^{33}`$ erg s-1 for the luminosity (0.5-10 keV) assuming a distance of 1.9 kpc (Sadakane et al. 1985).

## 3 The HETGS Spectra

The spectral continuum is well modeled in the zero order spectrum and we use it to subtract the continuum from the summed and background subtracted first order grating data. Figure 1 shows the residual spectrum between 3.1 and 12.5 Å .

### 3.1 Fluorescence Emission

The spectrum in Figure 1 shows a variety of fluorescent lines (green model) from Mg, Si, S, Ar, Ca, and Fe. The line properties, including possible identifications, are listed in Table 1. For identifications we use the adjusted calculations from House (1969). A few are also tabulated in Mewe (1994). The measured lines do not show significant velocity shifts or broadening and can be generally identified with near neutral charge states. For the elements Fe, Ca, and Ar we cannot distinguish between different charge states and here table 1 lists the most likely range based on the wavelengths given in House (1969). In the case of Mg, Si, and S we can identify various states, however for the lowest charge state we can again only identify a range of most likely ionization states.

For each detected line, we calculate an equivalent hydrogen column density required to produce the observed line fluxes using the method described in Sako et al. (1999). Fluorescent yields and photoionization cross sections are adopted from Kaastra and Mewe (1993) and Verner et al. (1996), respectively. Kaastra and Mewe (1993) predict an Fe K$`\beta `$/K$`\alpha `$ flux ratio of $`0.13`$. We fix the width of the K$`\beta `$ line to that of the K$`\alpha `$ line and infer a ratio of $`0.124\pm 0.095`$. Since this consistent with the theoretical value we assume that the iron line flux, and, hence the estimated column density, is not affected by resonant Auger destruction (Ross, Fabian, and Brandt 1996). Resonant Auger destruction occurs in regions, where line optical depths are significant. A fluorescent photon from the K-shell can be resonantely absorbed by ions that have vacancies in the L-shell. However, due to high autoionization yields in atoms of much lower Z than Fe, the L-shell ion will most likely Auger decay. In table 1 the quantity $`Y_{\mathrm{eff}}`$ denotes the fraction of the Fe K$`\alpha `$ column density for each identified charge state of the other elements. For Fe K $`\alpha `$ we derive a column density of 3.3$`\times 10^{22}`$ cm-2. If the $`Y_{\mathrm{eff}}`$s are added up for each element we observe a correlation of the column density with $`Z`$, which could indicate the process of resonant Auger destruction at work. Such an interpretation was also given by Sako et al. (1999). However, for Si such an interpretation is problematic, since the Si VII and Si VIII lines appear exceptionally strong braking that correlation. In this respect the detection of fluorescent lines from L-shell ions suggests a somewhat modified explanation for the observed behavior.

TABLE 1
FLOURESCENT LINE ENERGIES AND FLUXES.     Ion type1 $`\lambda _{exp}^1`$ $`\lambda _{meas}`$ Flux2 Y$`{}_{}{}^{3}{}_{eff}{}^{}`$     Å Å     Fe K $`\beta `$ 1.756 $`1.748\pm 0.020`$ $`1.29\pm 0.88`$ 0.950    K $`\alpha `$ 1.937 $`1.937\pm 0.001`$ $`10.40\pm 1.08`$ 1.000    II - XI       Ca II-VII 3.359-3.352 $`3.356\pm 0.005`$ $`0.77\pm 0.04`$ 1.015   Ar VI-IX 4.186-4.178 $`4.184\pm 0.010`$ $`0.82\pm 0.04`$ 0.460   S IX 5.320 $`5.322\pm 0.018`$ $`0.43\pm 0.04`$ 0.046    IV-VIII 5.370-5.356 $`5.365\pm 0.018`$ $`1.26\pm 0.06`$ 0.147   Si IX 6.947 $`6.937\pm 0.004`$ $`0.33\pm 0.09`$ 0.017    VIII 7.007 $`7.008\pm 0.005`$ $`0.56\pm 0.14`$ 0.036    VII 7.063 $`7.059\pm 0.004`$ $`0.75\pm 0.14`$ 0.054    VI 7.117 $`7.117\pm 0.005`$ $`0.52\pm 0.14`$ 0.042    II-VI 7.126-7.122 $`7.124\pm 0.005`$ $`0.74\pm 0.13`$ 0.048   Mg V 9.814 $`9.827\pm 0.011`$ $`0.15\pm 0.07`$ 0.0097    II-IV 9.890-9.883 $`9.893\pm 0.008`$ $`0.17\pm 0.08`$ 0.0115      1from House 1969
2in units of $`10^5`$ photons cm-2s-1
3fraction of column density with respect to Fe K (uncertainties are $`15\%`$, in case of Mg $`45\%`$

### 3.2 Radiative Recombination Continua (RRC)

In photoionized plasmas, free electrons are captured by ions, which result in a continuous emission feature above the recombination edge (for a review, see Liedahl et al. 2001 and references therein). Since the plasma is overionized relative to the local electron temperature, the widths of these continua are narrow and appear as line-like features in the spectrum. RRCs from Ne X and Ne IX at 9.10 Å and 10.37 Å, respectively, are clearly observed. Although some of the features are rather weak, they are consistent with kT = $`10\pm 2\mathrm{eV}`$ (for Ne X) or a temperature of $`T1.2\times 10^5\mathrm{K}`$. The Mg XI RRC is blended with the Si fluorescence lines, and a possible RRC of O VIII is observed at $`14.25`$ Å. The detection of such narrow RRCs clearly indicates that there is a warm photoionized medium associated with the wind.

For Ne X at 12.13 Å and 9.10 Å , the measured Ly$`\alpha `$/RRC flux ratio is $`2.7\pm 0.58`$ after correcting for interstellar column density, and may be higher if there is additional absorption above the Galactic value. Liedahl and Paerels (1996) predict for the optically thin limit a ratio of 1.3 for a purely recombining plasma. Since this is inconsistent with the observed ratio, we assume that the assumption of a purely recombining plasma is not correct here. A more recent calculation by Kinkhabwala et al. (2001) also includes photoexitations in a more self-consistent manner and predicts a ratio of Ly$`\alpha `$/RRC $`=2.7`$ at $`kT=10\mathrm{eV}`$ which agrees well with our observations. The spectrum also exhibits very weak, but detectable amounts of line flux from Fe L ions. This is consistent with the prediction by Kallman et al. (1996) that recombination emission from K-shell ions dominates the spectrum in an overionized plasma.

### 3.3 H- and He-like Lines

The spectrum shown in Figure 2 also shows very strong emission lines (blue model) from H- and He-like S, Si, Mg, Ne, and possibly O. The measured line fluxes range between 0.2 and $`3\times 10^5\mathrm{photons}\mathrm{cm}^2\mathrm{s}^1`$, and appear unresolved with no significant systematic line shifts. In the case of Ne X, four lines from the Lyman series are detected, with Ly$`\alpha `$:$`\beta `$:$`\gamma `$:$`\delta `$ line ratios of $`1.0`$:$`0.32`$:$`0.24`$:$`0.10`$ (with uncertainties of the order of 25$`\%`$). The fluxes of higher lines in the Lyman series are larger than those expected from either a purely recombining plasma or a purely collisionally ionized plasma, and implies the presence of photoexcitation. This is consistent with the high Ly$`\alpha `$/RRC flux ratio found in the previous section. For a photoionized plasma including the effects of photoexcitation self-consistently, the theoretical ratios for Ne X at $`kT=10\mathrm{eV}`$ are $`1.0:0.35:0.16:0.08`$ at a column density of $`N_{\mathrm{NeX}}=3\times 10^{17}\mathrm{cm}^2`$ and a turbulent velocity of $`300\mathrm{km}\mathrm{s}^1`$ (Kinkhabwala et al. 2001). The Lyman series for Mg XII is difficult to assess since the Ly$`\beta `$ and $`\gamma `$ appear blended.

The He-like lines of Si XIII and Mg XI are shown in Figure 2 (middle and right). In both cases the resonance and forbidden line fluxes have roughly equal strengths, which is not expected in a predominantly low density recombining plasma, but consistent with the fact that effects from photoexcitation have to be considered. It would therefore be incorrect to interprete these triplets using the calculations by Porquet and Dubau (2000), which in fact would predict a hot and dense collisional-dominated plasma. Such an interpretation would be in contrast with the results above, which indicate the presence of cool and warm photoionized media through the detection of narrow RRC and moderate line strength. Here Wojdowski et al. (2001) suggest that resonance scattering is the most probable probable process to enhance the resonance line in the He-like triplets.

We have yet to include contributions from the companion wind itself. In this case we have to consider collisionally ionized plasmas. Using the Astrophysical Plasma Emission Code (APEC, Smith et al. 2001), we can calculate line ratios between the Si XIII resonance line and Fe L lines at representative temperatures. Temperatures in O-stars range between $`5\times 10^6`$ K (as measured in $`\zeta Ori`$, (Waldron and Cassinelli 2001) and in an extreme case 5$`\times 10^7`$ K (as measured in $`\theta `$ Ori C, Schulz et al. 2000) from highly resolved X-ray lines. For now we assume that this is also a representative range of possible temperatures in the winds of B0.5 Ib supergiants. In the high temperature case we can exclude a significant contribution, because here Fe XXIV lines should dominate the spectrum, which we do not observe. Similarly, in the low temperature case Fe XVII lines should dominate the spectrum, which we do not observe either. For the intermediate range we estimate, under the assumption that the sporadic and weak detections of Fe L lines indeed reflect the existence of collisionally ionized emission, a contribution of a few percent at most.

## 4 Conclusions

The X-ray spectrum of Vela X-1 during eclipse shows a large variety of emission features including lines from H-like and He-like ions, RRCs, a population of very weak Fe L lines, and a number of fluorescence emission lines from mid-$`Z`$ elements and iron indicating a complex structure of the X-ray emitting plasma.

From the fit of the shape of the Ne X RRC at 9.10 Å we deduce a temperature of $`1.2\times 10^5\mathrm{K}`$ for the plasma. Its flux ratio with the Ne X Lyman $`\alpha `$ line is consistent with an optically thin photoionized plasma. The spectrum shows clear evidence of photoexcitation processes, for which resonance scattering from the strong eclipsed X-ray continuum may be a likely source. More detailed calculations on this issue are required.

Fluorescence lines from various charge states were detected and partially identified. They include strong emission from M-shell ions, but also from various L-shell ions. Fluorescence yields are generally extremely low in mid-$`Z`$ elements even at moderate optical depths. The correlation of the inferred column density with Z implies that autoionization dominates the fluorescence yields in these ions and the fluorescent plasma is not optically thin. If the fluorescent plasma coexists with an optically thin photoionized plasma or in a region outside the warm photoionized region, its has to be in the form of dense cool clumps. Both scenarios may be the case, which is not inconsistent with current stellar wind models. For example, Feldmeier (1995) pointed out that the wind can be devided into two distinct regions, an active inner one where shell to shell collisons occur and an outer, quiescent region with ”old” material.

The detection of fluorescence from L-shell ions is still mysterious. Given the relative low fluorescent yields for these ions, resonant Auger destruction would have almost completely destroyed these lines. Further investigations of the possible mechanisms and physical conditions that prevent the lines from being destroyed are required.

The authors want to thank the MIT HETG team and the Chandra X-ray Center for their support, and Ali Kinkhabwala for providing results from calculations of photoionized plasmas in progress. This research is funded by contracts SV-61010 and NAS8-39073.


text:
# Helioseismology as a New Constraint on SUSY Dark Matter

## 1 Introduction

The dynamical behaviour of various astronomical objects, from galaxies to galaxy clusters and to large-scale structure, in the observed universe can only be understood if the dominant component of the mean matter density is dark, amounting to $`\mathrm{\Omega }_m=0.3\pm 0.1.`$ Constraints from primordial nucleosynthesis of the light elements provide a compelling measure of the mean baryon abundance, $`\mathrm{\Omega }_b=0.05\pm 0.02.`$ The bulk of the dark matter is consequently non-baryonic, and the existence of particles that interact with ordinary matter on the scale of the weak force, so-called Weakly Interacting Massive Particles (WIMPs), provides one of the best-motivated candidates, arising from the lightest stable particle predicted by super-symmetry theories (SUSY), for resolving this problem.

One can show (Lee & Weinberg 1977) that if such a stable particle exists, its relic abundance is $`\mathrm{\Omega }_xh^23\times 10^{27}cm^3s^1/\sigma _av`$, where $`\sigma _av`$ is the thermally averaged product of annihilation cross-section and relative velocity. Alternatively, if a new particle with weak-scale interactions exists, then its annihilation cross-section can be estimated to be $`\sigma _av10^{25}cm^3s^1`$ (Jungman et.al. 1996). Indeed, this constitutes one of the strongest arguments for considering the WIMP to be one of the best particle candidates for the non-baryonic component of the dark matter.

These particles are predicted by super-symmetry theories (SUSY) to be thermally produced at the early stages of the universe (Ellis 2001). The best WIMP candidates are neutralinos, but other types of WIMPs (axions, axinos, gravitinos and Wimpzillas) are not excluded as an alternative to the constituents of the dark matter component (Roszkowski 2001). Nor are charged super-heavy particles with strong interactions (Albuquerque et.al. 2000).

It is usually assumed that the WIMPs are a thermal relic of the big bang. The resulting thermic relic particle should have a mass inferior to $`340TeV`$, the unitary limit (Griest and Kamionkowski 1990). However, for SUSY models, this bound is typically 2 orders of magnitude stronger (Jungman et.al. 1996). The neutralino is the best candidate to be produced thermally, while CDM axinos, gravitinos and Wimpzillas are the best candidates to be produced by a non-thermal mechanism. For example, in the early Universe, vacuum fluctuations during or after inflation can produce a class of very weak super-massive particles with mass in the range of $`10^{12}GeV`$ up to $`10^{13}GeV`$ (Kuzmin & Tkachev 1999).

Many experiments around the world are currently engaged in searching for dark matter by assuming that it is dominated by the light neutralinos, as predicted by the minimal supersymmetric extension of the Standard Model (Ellis 2001). In particular, the DAMA collaboration (Bernabei et al. 2000) has reported evidence of an annual modulation of recoil energy, which was interpreted as being due to scattering of some cold dark matter particles with masses between 50 and 100 GeV, and spin-independent cross-sections on a proton of between $`10^{42}cm^2`$ and $`10^{41}cm^2`$.

Neutralinos might not constitute all the cold dark matter, but could be complemented by other particles such as axions and superheavy relics. Neutralinos are neutral Majorana particles, the lightest SUSY particles that are massive and stable (Jungman et.al. 1996). Neutralino properties as dark matter candidates are quite dependent on the neutralino type, such as the bino, the wino and the two neutral higgsinos (Rowszkowski 2001). The neutralino cross-section for elastic scattering on a nucleus is expected to be typically very small, roughly $`10^{42}\mathrm{cm}^2`$. This is because the elastic cross-section is related to the cross-section of neutralino annihilation in the early Universe which must be only a fraction of the weak interaction strength in order to have $`\mathrm{\Omega }_xh^21`$ (Bergström 2000). In this paper we adopt the neutralino as our WIMP test particle, unless stated otherwise.

The Sun traps WIMPs in the course of its lifetime, as we discuss below. If the Sun were to contain even a minute mass fraction of WIMPs, there could be a significant influence on its central thermal structure. Helioseismology provides a means for an independent test of the validity of this idea and can be used to discriminate between possible candidates for the Dark Matter (Kaplan et al. 1991). Indeed, helioseismology has become a mature discipline through the detection of more than three thousand solar acoustic modes with which it has been possible to constrain the internal structure of the Sun (Gough 1996). The large quantity of seismic data obtained by the ground-based network during the last 10 years, as well as the three seismic instruments on board the SOHO mission during the last 5 years, contributed to improve the accuracy of the frequency determination to a level of 1 part in $`10^5`$ (Bertello et al. 2000; Garcia et al. 2001). This has allowed us to obtain a solar model for which the equilibrium thermodynamical structure accurately reproduces the observed seismic data, the so-called seismic solar model.

In this paper, we discuss how WIMPs can modify the evolution of the Sun and ultimately the present internal structure, namely in the nuclear region. In particular, we are interested in WIMPs that interact with baryonic matter through a weak interaction with a spin-independent scattering cross-section with values in the range of the current experiments for the detection of dark matter, such as the DAMA collaboration, and other ongoing or future experiments. It follows from our analysis that the possible WIMP candidates proposed by the DAMA experiment are in disagreement with our current models.

The next section contains a description of the present status of the standard solar model, and the current view of the different difficulties in modelling the observed Sun. In Section 3, we present a theoretical description of WIMPs and how they interact with baryonic matter. In particular, we discuss qualitatively the type of WIMPs capable of modifying the evolution of the Sun. In Section 4, we discuss solar models that evolve in halos of dark matter, and we discuss how those WIMPs can modify the evolution of the Sun. We conclude with a discussion of our new results, stressing their pertinence to the current status of the dark matter problem and Cosmology.

## 2 THE PRESENT SUN: THE SOLAR STANDARD MODEL

The Sun is a unique star for research because its proximity allows a superb quality of solar data, enabling precision measures of its luminosity, mass, radius and chemical composition. The Sun therefore naturally becomes a privileged target for testing stellar evolution theory. In recent years, different groups around the world have produced solar models in the framework of classical stellar evolution, taking into account the best known physics as well as all the available observational seismic data. This has led to determination of a well-established model for the Sun, the so-called standard solar model (Turck-Chiéze & Lopes 1993; Christensen-Dalsgaard et al. 1996; Brun, Turck-Chièze & Zahn 1999; Provost, Berthomieu & Morel 2000; Turck-Chiéze, Nghiem, Couvidat & Turcotte 2001, Turck-Chiz̀e S. et al. 2001a; Bahcall, Pinsonneault, Basu 2001), for which the acoustic modes are in very good agreement with observation. Furthermore, this model has established considerable consensus among the different research groups, concerning the predictions of the solar neutrino fluxes, and has unambiguously helped define the difference between the theoretical predictions and the experimental results.

The combined effort between helioseismology and stellar astrophysics has contributed to strongly constrain the internal structure of the present Sun. Indeed, among the different seismic diagnostics that the solar model is tested against, the radial profile of the square of the sound speed is known with a precision as high as $`0.3\%`$, and in particular in the nuclear region this precision goes up to $`0.2\%`$ (see Fig 1). This progress has been achieved through a detailed description of the microscopic physics, such as the nuclear reaction rates, the equation of state, the coefficient of opacities and the gravitational settling of chemical elements. More recently, the macroscopic mixing induced by differential rotation, seems to present a determining role in the evolution of the star. These physical processes take place in the thin shear layer, located between the radiative interior and the convection zone, the so-called tachocline (Elliott & Gough 1999).

The standard stellar evolution of the Sun assumes that the star is in hydrostatic equilibrium, is spherically symmetric and that the effects of rotation and magnetic fields are negligible. The present solar structure and evolution are computed starting from an initially homogeneous star of a given composition. The solar models used in this paper have been computed by using the CESAM code (Morel 1997). The main physical inputs are as follows. We used the solar composition of Grevesse, Noels & Sauval (1996), the OPAL96 tables of Iglesias & Rogers(1996), the equation of state of Rogers, Stenson & Iglesias (1996), the nuclear reaction rates of Adelberger et al. (1998), the Mitler (1977) prescription for the treatment of the nuclear screening rates, and the microscopic diffusion coefficients suggested by Michaud & Proffit (1993). The present structure of the Sun is obtained by evolving a initial star from the pre-main sequence, around 0.05 Gyr from the ZAMS, until its present age, 4.6 Gyr. The present solar model is obtained by choosing the initial Helium content, as well as the mixing length parameter of convection that best predicts the present solar luminosity and radius. The helioseismolgy requires that the present structure of the Sun, including the global quantities, should be determined with a very high precision. In particular, the calibration of the solar radius is done with a precision of $`10^5`$.

## 3 THE EVOLUTION OF THE SUN IN THE PRESENCE OF DARK MATTER

In order to explain the rotation of spiral galaxies, it is necessary to assume that these are immersed in a halo of dark matter. In the present model, we assume the halo of the Milky Way to be an isothermal gas of WIMPs, the distribution of their velocities being Maxwellian with a mean velocity $`\overline{v}240Kms^1`$. Furthermore, we assume that the local dark matter density near the Sun is about $`0.3GeVcm^3`$ (Jungman et.al. 1996). Typically, the local density by number is of the order of $`\mathrm{3\hspace{0.33em}10}^2cm^3`$ for WIMPs with $`m_\mathrm{x}10\mathrm{G}\mathrm{e}\mathrm{V}`$, and $`\mathrm{3\hspace{0.33em}10}^3cm^3`$ for WIMPs of $`m_\mathrm{x}100\mathrm{GeV}`$.

### 3.1 The WIMP scattering cross-section

Any WIMP that crosses the interior of the Sun may interact with a nucleus and lose enough energy to be trapped by the star’s gravitational potential well. The WIMP interactions with the baryonic matter will depend upon their mass, $`m_\mathrm{x}`$ and on their scattering cross-section on the nucleons, $`\sigma _p`$. The elastic scattering cross-section of relic WIMPs scattering off a nucleus in the solar interior depends on the individual cross-sections of the WIMPs scattering off of quarks and gluons. For non-relativistic Majorana particles, such as the neutralino, these can be divided into two separate types. The coherent part described by an effective scalar coupling between the WIMP and the nucleus is proportional to the number of nucleons in the nucleus. The incoherent component of WIMP-nucleus cross section results from an axial current interaction of a WIMP with the constituent quarks, and couples the spin of the WIMP to the total spin of the nucleus. These two types of cross-sections have been computed exactly (Goodman & Witten 1985, Jungman et.al. 1996, and references therein), and in the case of scalar interaction it strongly depends on the particle physics model used. At this stage, it is worth noticing that in the case of scalar (or coherent) interactions, the total scattering cross-section is proportional to the mass of the nucleus, $`\sigma ^{sc}\mu _i^2`$, where $`\mu _i=m_xm_i/(m_x+m_i)`$ is the reduced mass, and in the case of axial (or incoherent) interactions, it can be shown that the cross-section, $`\sigma ^{ax}\mu _p^2`$, depends on the spin of the nucleus as $`J(J+1)`$. In the case of scalar interactions, because of the constructive interference between all nucleons in the nucleus, the cross-section rises rapidly with the atomic number and, as a consequence, WIMPs scatter much more efficiently on the heavier nuclei. Inside the Sun, elements heavier than hydrogen contribute substantially to the energy transfer and to the capture rate. Nevertheless, the Sun’s evolution on the main sequence is strongly dominated by the production of helium rather than heavier chemical elements. In such a scenario, the helium nucleus will be the dominant source of scalar scattering, mainly due to its high abundance, given that heavier nuclei such as iron, nitrogen, carbon and oxygen are sub-dominant species. However, in more detailed studies on the capture of WIMPs by the Sun, the contribution of heavier elements seems to be more significant than was previously thought (Gould 1987). In fact, the scattering on these species will increase the total number of WIMPs accreted by the Sun leading to a more significant change on the structure of the solar core. Consequently, by choosing to consider the scattering only for helium, we make a conservative approach to study the impact of WIMPs on stellar evolution and Helioseismology.

The essential particle physics parameters which enter the dark matter problem are the mass of the WIMP $`m_x`$ and the elastic scattering cross-section, $`\sigma _p`$, of the various nuclei which constitute the solar material, and $`\sigma _av`$ the thermal average of the annihilation cross-section times the relative velocity of the WIMPs at freeze-out. The annihilation cross section can be calculated for each model (Krauss et.al. 1985) however, since we lack certain theoretical guidance, we decided to concentrate on a typical simple case for WIMP-nucleus scalar scattering, as previously discussed. We relate the cross-section of the WIMP on nuclei to the WIMP-proton cross section $`\sigma _p`$ taken as a free parameter.

### 3.2 The Capture and Annihilation Rates of WIMPs in the Sun

The abundance of WIMPs in the Sun depends on the WIMPs accumulated in the Sun by capture from the Galactic halo and is depleted by annihilation. If $`N_x`$ is the number of WIMPs in the Sun, then the differential equation governing the time evolution of $`N_x`$ is (Gould et.al 1987)

$`{\displaystyle \frac{dN_x}{dt}}=\mathrm{\Gamma }_cC_aN_x^2,`$ (1)
where $`\mathrm{\Gamma }_c`$ is the rate of accretion of WIMPs onto the Sun. The determination of $`\mathrm{\Gamma }_c`$ depends on the nature of the particle. If the halo density of WIMPs remains constant in time, $`\mathrm{\Gamma }_c`$ is time-independent. The second term accounts for depletion of WIMPs, and is twice the annihilation rate in the Sun, $`\mathrm{\Gamma }_a=1/2C_aN_x^2`$. The quantity $`C_a`$ depends on the WIMP total annihilation cross section times the relative velocity in the non-relativistic limit, $`\sigma _av`$, and the distribution of WIMPs in the Sun. The total annihilation cross-section of the new particle should possibly have a value that is consistent with the cosmological density of dark matter as specified, $`\sigma _av\mathrm{3\hspace{0.17em}10}^{26}cm^3s^1\mathrm{\Omega }_x^1h^2`$. It follows that, not knowing the initial content of relic WIMPs, $`\mathrm{\Omega }_x`$, we will consider different scenarios, such that $`\sigma _av<10^{26}cm^3s^1`$ (Jungman et.al. 1996).

Solving the previous first-order differential equation, the total number of WIMPs at time $`t`$, is given by

$`N_x(t)=\mathrm{\Gamma }_c\tau \mathrm{tanh}\left({\displaystyle \frac{t}{\tau }}\right),`$ (2)
and the annihilation rate is given by

$`\mathrm{\Gamma }_a(t)={\displaystyle \frac{1}{2}}\mathrm{\Gamma }_c\mathrm{tanh}^2\left({\displaystyle \frac{t}{\tau }}\right),`$ (3)
where $`\tau =1/\sqrt{\mathrm{\Gamma }_cC_a}`$ is the time scale for capture and annihilation of WIMPs to equilibrate. As the age of the Sun is much larger than the equilibrium time-scale, $`t_{}/\tau 1`$, equilibrium is reached very rapidly, then $`\mathrm{\Gamma }_a=1/2\mathrm{\Gamma }_c`$. It follows that

$`{\displaystyle \frac{t_{}}{\tau }}=322\left({\displaystyle \frac{\mathrm{\Gamma }_c}{s^1}}\right)^{1/2}\left({\displaystyle \frac{\sigma _av}{cm^3s^1}}\right)^{1/2}\left({\displaystyle \frac{m_x}{10GeV}}\right)^{3/4},`$ (4)
where $`t_{}`$ is the present age of the Sun (Jungman et.al. 1996).

The rate of accretion of WIMPs in the Sun was first calculated by Press and Spergel (1985). Even though the detailed computation of the accretion can be quite elaborate, the basic idea is simple (Srednicki et.al. 1987); if the WIMP elastically scatters from a nucleus with a velocity smaller than the escape velocity, $`v_{esc}=\sqrt{2GM_{}/R_{}}`$, then the WIMPs will be captured. The WIMPs have a typical velocity of the order of $`300km/s`$, which is smaller than the escape velocity at the surface of the Sun $`618km/s`$, consequently, they are trapped quite efficiently in the solar interior. If we consider that all the WIMPs that hit the Sun lose enough energy to be trapped, the trapping rate is given by the product of the surface area of the Sun, about $`\mathrm{6.1\hspace{0.33em}10}^{22}cm^2`$, and the flux of WIMPs, of about $`n_\mathrm{x}\overline{v}`$. It reads $`\mathrm{\Gamma }_\mathrm{c}=4\pi R_{}^2n_\mathrm{x}\overline{v}`$. Gravitational focussing effects enhance the previous trapping rate by $`2GM_{}/\overline{v}^2`$. This occurs because of the enlargement of the Sun’s effective area due to the gravitational well and the requirement that the impact parameter be less than $`R_{}`$. This yields a trapping rate of $`\mathrm{\Gamma }_\mathrm{c}=\pi (2GM_{}R_{})n_\mathrm{x}/\overline{v}`$. Bouquet and Salati (1989) have generalized and refined the capture rate for main-sequence stars (Salati & Silk 1989)

$`\mathrm{\Gamma }_c=10^{32}\mathrm{s}^1\left({\displaystyle \frac{\rho _\mathrm{x}}{1M_{}\mathrm{pc}^3}}\right)\left({\displaystyle \frac{m_p}{m_x}}\right)\left({\displaystyle \frac{300\mathrm{kms}^1}{\overline{v}_x}}\right)`$
$`\left[1+0.16\left({\displaystyle \frac{\overline{v}_x}{300\mathrm{kms}^1}}\right)^2\left({\displaystyle \frac{M_{}}{M}}\right)\left({\displaystyle \frac{R}{R_{}}}\right)\right]`$
$`\left({\displaystyle \frac{M}{M_{}}}\right)\left({\displaystyle \frac{R}{R_{}}}\right)\mathrm{min}[1,{\displaystyle \frac{\sigma _p}{\sigma _{}}}{\displaystyle \frac{M}{M_{}}}\left({\displaystyle \frac{R_{}}{R}}\right)^2].`$ (5)
The total number of WIMPs accreted, $`N_x`$, can be computed taking into account that $`t_{}/\tau >>1`$, and is

$`N_x\mathrm{\Gamma }_c\tau .`$ (6)
The total number of WIMPs at present strongly depends on the ratio $`\sqrt{\sigma _av/\mathrm{\Gamma }_c}`$. If we consider $`\sigma _av`$ to be constant, the capture is weak and proportional to the scattering cross-section and $`N_x`$ increases with $`\mathrm{\Gamma }_c`$, for the case of low scattering cross-sections. On the contrary, for large scattering cross-sections, any WIMP which enters the star is captured, the captured flux saturates at the level of the incoming flux and the system reaches an equilibrium. Note that evaporation is unimportant for WIMPs in the mass range considered here, namely $`m_x>10`$ GeV (Gould 1987).

This concentration of WIMPs in the solar core has a marginal effect on the evolution of the Sun. Luckily, the precision presently obtained by seismic diagnostics allows to determine in certain cases the effect of WIMPs on the solar core. In fact, a more accurate expression for the capture should take into account the way that WIMPs scatter from the nucleus through scalar interactions and/or axial-vector interactions, as well as second-order effects (Jungman et.al. 1996). However, we choose to focus on a simplified picture in this first approach to the problem.

### 3.3 The energy transport by WIMPs in the solar interior

The energy transport by WIMPs in the solar interior is governed by three natural length scales: the mean free path of the WIMP, $`l_x`$, the inverse of the logarithmic temperature gradient, $`\mathrm{ln}T`$, and a typical geometric dimension of the system, such as the WIMP scale height, $`r_x`$. Bouquet and Salati (1989) showed that the WIMP distribution is approximatively gaussian and therefore

$`n_x(r)=n_0\mathrm{exp}\left[{\displaystyle \frac{r^2}{r_x^2}}\right]`$ (7)
where $`n_0=N_x/\pi \sqrt{\pi }r_x^3`$, and $`r_x^2=3T_c/(2\pi m_xG\rho )`$ which typically is $`r_x10^2R_{}m_{100}^{1/2}`$ with $`m_{100}`$ the WIMP mass in units $`100`$ GeV.

There are two extreme regimes for the energy transport by WIMPs,

characterized by the values of the Knudsen number, $`K`$, defined as the ratio $`l_x/r_x`$. The Knudsen number strongly depends on the scattering cross-sections of the WIMPs on nuclei, and typically goes like $`K30m_{100}^{3/2}/\sigma _{36}`$ where $`\sigma _{36}`$ is the scattering cross section measured in units $`10^{36}`$ cm<sup>2</sup>.

In the case of small cross-sections, the WIMP mean free path is larger than the scale length $`r_x`$, $`K>1`$, and two successive collisions are widely separated. This is the so-called large Knudsen number regime (or Knudsen regime). The transport of energy is non-local. Thus the WIMPs undergo few interactions on each orbit and therefore are not in thermal equilibrium with the nuclei. In such a case, it is convenient to define the average temperature of the WIMP core, $`T_x`$, as being the temperature of the star at $`r_x`$. At the center of the star, the WIMP temperature (or the averaged kinetic energy) is lower than the temperature of the nuclei, and in such conditions the net effect of the WIMP-nucleus collisions is to transfer energy from nuclear matter to WIMPs. Conversely, in the outskirts of the WIMP core, in certain cases the WIMPs give back the energy gained in the center of the star. This process is a very efficient mechanism at evacuating the energy from the nuclear region. In such a regime no precise analytic result exists. Nevertheless, an approximative solution has been found by Spergel and Press (1985). The luminosity carried by WIMPs can then be evaluated by

$`L_{\mathrm{lp}}(r)=32\sqrt{2\pi }{\displaystyle \frac{\sigma _p\mu _p}{m_x+m_p}}{\displaystyle _0^r}n_p(r)\overline{T}(r)n_x(r)r^2𝑑r`$ (8)
where $`\overline{T}`$ reads

$`\overline{T}(r)=\sqrt{{\displaystyle \frac{m_pT_x+m_xT(r)}{m_xm_p}}}\left[T(r)T_x\right].`$ (9)
In the opposite limit, $`K<1`$, the transport of energy is local and the energy transport proceeds by conduction, allowing for a much simpler treatment (Gilliland et.al. 1986). In this regime, the WIMP mean free path is much smaller than the dimension of the region where they are trapped. The WIMPs interact sufficiently often to be in local thermal equilibrium with the nuclei and therefore it is possible to define a temperature of the WIMP, which is equal to the temperature of the nuclei.

$`L_{\mathrm{sp}}(r)=4\pi r^2n(r)\sqrt{{\displaystyle \frac{T(r)}{m_x}}}l_{sp}(r)T(r)`$ (10)
where $`l_{sp}(r)`$ is the local mean free path of the WIMP, it reads

$`l_{\mathrm{sp}}(r)=\left[{\displaystyle \frac{\rho (r)}{m_u}}{\displaystyle \sigma _i\frac{X_i(r)}{A_i}}\right]^1.`$ (11)
The range of WIMP masses allowed by the cosmological model is quite large, even if we restrain ourselves to the WIMPs that are produced in a thermal scenario. Consequently, depending also on their scattering cross-section, the WIMPs can transport energy, not only in the conductive region and in the Knudsen regime, but also in an intermediate case. We have chosen to define the intermediate case as the one where the transport of energy is determined by the interpolation formula (Gilliland et.al. 1986)

$`L_x(r)={\displaystyle \frac{L_{lp}(r)L_{sp}(r)}{L_{lp}(r)+L_{sp}(r)}}.`$ (12)
This formula appropriately converges to the right approximation in each of the two regimes, accordingly to the value of $`K`$. The errors introduced are within the precision of our computation.

## 4 DISCUSSION AND CONCLUSIONS

The thermodynamical structure in the interior of the Sun, namely in the nuclear region, is presently known with a precision of much less than a few per cent. This level of accuracy in constraining the solar interior has been achieved by a systematic study of the differences between the acoustic spectrum obtained from helioseismology experiments and the theoretical spectrum. Presently, this difference is less than $`10\mu Hz`$ for almost all of the 3000 modes that probe the interior of the Sun (Gough 1996). This level of precision in the description of the solar core allows us to discriminate physical processes that could not be discussed otherwise, in particular those that present a peculiar behaviour, as seems to be the case for WIMPs trapped in the solar core. In Fig. 1, we compare the square of the sound speed, $`c_s^2`$, of different solar models evolving within the presence of a halo of WIMPs and the solar standard model (Brun, Turck-Chièze & Morel 1998). The changes induced by the presence of WIMPs are concentrated in the inner core within 10% of the solar radius, typically seen in the profiles of the temperature, density and molecular weight.

Indeed, the WIMPs are thermalized within the solar core and are on Keplerian orbits around the solar center, interacting through elastic scattering with the solar nuclei, such as helium, and thereby providing an alternative mechanism of energy transport other than radiation. The result is a nearly flat temperature distribution, leading to an isothermal core. Consequently, the central temperature is reduced. This reduction of temperature has two main consequences: since central pressure support must be maintained, due to the hydrostatic equilibrium, the central density is increased in the solar models with WIMPs, and since less hydrogen is burnt at the centre of the Sun, the central helium abundance and the central molecular weight are smaller than in standard solar models. The increase of the central density and hydrogen partially offset the effect of lowering the central temperature in the central production of energy. In fact, this is the reason why minor changes are required to the initial helium abundance and the mixing-length parameter in order to produce a solar model of the Sun with the observed luminosity and solar radius. This readily leads to a balance between the temperature, and the molecular weight in the core, leading to the peculiar profile of the square of the sound speed, $`c_s^2T/\mu `$. This seems to be the case for most of the solar models within WIMP halos. The balance between the temperature and molecular weight is critical for the profile in the center of the star, leading to some of the profiles presented in Fig. 1. In the same figure, we display the inversion of the sound speed obtained from the data of the three seismic experiments on board the SOHO satellite. It follows from our analysis that the presence of WIMPs in the solar core produces changes in the solar sound speed of the same order of magnitude as the difference between the sound speed of the standard solar model and the inverted sound speed. It is important to remark that the inversion of the sound speed still presents some uncertainty in the central region due to the lack of seismic data, mainly due to the small number of acoustic modes that reach the nuclear region. Furthermore, the inversions are not very reliable at the surface, above $`98\%`$ of the solar radius, due to a poor description of the interaction of acoustic waves with the radiation field and the turbulent convection, namely, in the superadiabatic region (Lopes & Gough 2001). However, we can establish with certainty that the difference between the inverted sound speed and the theoretical sound speed is known with a precision of $`0.3\%`$ in the solar interior, within $`95\%`$ of the solar radius. The evolution of the Sun in a halo of WIMPs will increase the evacuation of energy from the solar core. The WIMPs will work as a ’cold bridge’ between the core and the more external layers of the Sun. The magnitude of the effect is proportional to the total number of WIMPs concentrated in the solar core. Nevertheless, even if some systematic effect is present in the inversion of the sound speed, the presence of WIMPs in the solar core leads to a quite different nucleosynthesis history from the solar standard model case, and from that to a peculiar radial profile of the sound speed. In this way, the effect of WIMPs in the solar core can be inferred on the basis of seismic diagnostics, such as the inversion of the square of the sound speed, among other possible techniques. The proposed method constitutes a new way to disentangle the contribution of different non-baryonic particles to the dark matter.

The luminosity in the core of the Sun is presently known with a precision of one part in $`10^3`$. In the coming years, it is very likely that the new seismic data available from the SOHO experiments, will allow us to obtain a seismic model of the Sun with an accuracy of $`10^5`$. In such conditions, the Sun can and should be used as an excellent probe for dark matter in our own galaxy. In Figs. 2 and 3, we compute the ratio of the WIMP luminosity against the Sun’s luminosity produced in the inner core of $`5\%`$ of the solar radius. A significative region of the $`\sigma _{scat}\sigma _av`$ plot shows changes in the solar luminosity of the order of $`10^3`$. This order of magnitude on the luminosity produced in the solar core can be tested through seismological data. In particular, we are interested in the lighter WIMPs, $`m_x<100GeV`$, and the smallest scattering cross-section, $`10^{45}cm^2`$ up to $`10^{40}cm^2`$. This range of parameters are presently being tested by the DAMA and CDMS experiments, among others, and are also well within the range of the parameters of future helioseismological experiments (see Fig. 3). It is interesting to note, that the X-ray Quantum Calorimeter (XQC) experiment (D. McCammon et.al. 2001) may exclude scattering cross sections bigger than about $`10^{29}`$ cm<sup>2</sup> for the mass range considered here, relevant for strongly interacting massive particles recently hypothesised as an alternative dark matter candidate to WIMPs and that penetrate neither subterranean laboratories nor the solar core (Wandelt et.al. 2000), hence being complementary to helioseismology in potentially excluding overlapping regions in parameter space.

If we believe in the simple model presented in this paper, DAMA candidate WIMPs with masses of 60 GeV and annihilation cross-section of the order of $`10^{32}cm^3/s`$ cannot exist (see Fig. 1), otherwise their effect in the solar core should already have been identified by seismic diagnostics. However, we stress that in order to determine with certainty the range of masses and cross-sections for WIMPs that is in disagreement with the helioseismological results, a more careful analysis of the different regimes of energy transport by WIMPs should be done. Furthermore, the microscopic physics of this region of the star is not fully established, as there remains uncertainty in certain nuclear reaction rates, such as the p+p reaction and the dynamical screening in other nuclear interactions such as $`{}_{}{}^{3}He+^4He`$ and $`{}_{}{}^{4}He+^4He`$.

An important diagnostic of the core can be obtained from the seismology of gravity modes. Indeed, it is the low-degree internal gravity modes that are the most sensitive to the conditions in the core, the region where substantial deviations from the so-called standard solar model might occur. Solar models evolving in the presence of dark matter have a g-mode period spacing that is drastically different from that of other solar models, mainly due to the peculiar distribution of density on the nuclear region that occurs as a consequence of the energy transport enhanced by the WIMPs. Therefore, g-mode observations could provide a sensitive test to the radial distribution of density in the solar core, and ultimately to the presence of dark matter in the solar neighborhood. The observation of gravity modes by SOHO seismic experiments, such as GOLF (Turck-Chiéze et al. 2001b), could ultimately strongly constrain the physics of the solar core.

We have in this paper focused on a slightly simplified version of WIMP energy transport. There are several improvements which can be included in the analysis (Lopes et al. 2001), e.g. in the Knudsen regime the deviation from isotropy leads to a radius dependent luminosity suppression (Gould & Raffelt 1990b); in the conductive regime, the mass dependence of the thermal conductivity gives an additional factor $`2\sqrt{m_{100}}`$ for scattering on helium (Gould & Raffelt 1990a). Also inclusion of scattering off other light elements might affect the results slightly. Finally, the transition region between non-local and conductive energy transport could be expressed in terms of the Knudsen number (Dearborn et al. 1991; Kaplan et al. 1991) instead of the intuitive interpolation formula used here, which would be useful for a detailed investigation of the transition region.

In conclusion, we have identified the range of WIMP masses, scattering cross-sections and annihilation cross-sections which reduce the luminosity in the core of the Sun. These effects on the solar structure are now within the range of effects capable of being probed by the diagnostic capabilities of helioseismology. We did not concern ourselves with particular particle physics models, but considered a generic case, for WIMPs which have a large range of masses, scattering cross-sections and annihilation cross-sections. The effect of WIMPs in the core of the Sun is of the same order of magnitude as the microscopic physics and dynamical processes that are now being discussed in the framework of stellar evolution in the light of the most recent results of helioseismology. Studies of the sun and of cosmology may have much to gain from each other.

## Acknowledgments

The authors thanks P. Morel for using the CESAM code and S. Turck-Chièze and R. Garcia for stimulating discussions and helioseismic analysis. The authors would like also to thank the referee D. Jungman for his valuable suggestions that lead to improve the original manuscript. SHH is supported by a Marie Curie Fellowship of the European Community under the contract HPMFCT-2000-00607. IPL is grateful for support by a grant from Fundação para a Ciência e Tecnologia.