text:
# Large scale ionization of the Radio Arc region by the Quintuplet and the Arches clusters Based on observations with ISO, an ESA project with instruments funded by ESA Member States (especially the PI countries: France, Germany, the Netherlands and the United Kingdom) and with the participation of ISAS and NASA.

## 1 Introduction

With the aim of understanding the thermal balance and ionization of the interstellar medium (ISM) in the Galactic center region (GCR), we have undertaken a systematic study at millimeter and infrared wavelengths of a sample of clouds distributed along the Nuclear Bulge (in notation of Mezger et al. Mezger96 (1996)) and the “Clump 2” complex. These clouds exhibit large column densities of warm H2 (with temperatures of $`150`$ K) which constitute an important fraction of the total column densities of molecular gas (Rodríguez-Fernández et al. RF01 (2001)). In addition, the comparison of fine-structure lines of several ionized atoms with hydrogen radio recombination lines indicate that the GCR ISM should be rather inhomogeneous and that the the ionized gas could be distributed, for instance, in relatively low density ionized bubbles (Martín-Pintado et al. 2000a ; Rodríguez-Fernández et al. in prep.).

The two sources of our sample with more intense fine structure lines (M+0.16–0.10 and M+0.21–0.12, see Martín-Pintado et al. 2000a) are located in the vicinity of the Radio Arc, which is one of the most prominent radiocontinuum features in the GCR. The Radio Arc, located at $`l0\stackrel{}{.}17`$, is perpendicular to the Galactic plane and it is apparently connected to Sgr A by a “bridge” of radiocontinuum emission at positive galactic latitude (see e.g. Güsten Gusten89 (1989)). VLA observations showed that the Radio Arc is composed by long and thin straight filaments that emit non-thermal radiation (Yusef-Zadeh, Morris & Chance YZ84 (1984); Yusef-Zadeh, Morris & Gorkom YZ87 (1987)). These filaments are usually known as the Non-Thermal Filaments (NTFs) and indicate the presence of a strong poloidal component of the magnetic field in the GCR (Yusef-Zadeh, Morris & Chance YZ84 (1984); Yusef-Zadeh & Morris 1987a , 1987b ).

On the other hand, the “bridge” is composed by a set of arched filaments that emit thermal radiation (Thermal or Arched Filaments). In the vicinity of the NTFs and probably interacting with them there are other prominent features that also emit thermal radiocontinuum. These are G0.18–0.04 (the Sickle), which is approximately located where the NTFs cross the galactic plane, and a gun-shaped source located south of the Sickle known as the Pistol Nebula (G0.15–0.05). The thermal nature of the Sickle, the Pistol Nebula and the Thermal filaments has been established by the association of the radiocontinuum emission with radio recombination line emission (Pauls et al. Pauls76 (1976); Paul & Mezger Pauls80 (1980); Yusef-Zadeh, Morris & van Gorkom YZ87 (1987)).

In the recent years, the origin of the ionization of these thermal features has been a subject of great interest. It was first thought that the Sickle and the Thermal Filaments were the surfaces of molecular clouds ionized by collisions between the neutrals and the ions that spiral around the magnetic field lines associated with the NTFs (Yusef-Zadeh & Morris 1987c ; Morris & Yusef-Zadeh Morris89 (1989); Serabyn & Güsten Serabyn91 (1991); Serabyn & Morris Serabyn94 (1994)). However, with the discovery of the outstanding stellar cluster known as the \`\`Quintuplet" or AGFL 2004 (Nagata et al. Nagata90 (1990); Okuda et al. Okuda90 (1990); Glass et al. Glass90 (1990)) the effect of ultra-violet (UV) radiation on the ionization of the Sickle was revised. Nowadays, it is commonly assumed that this cluster is the main ionization source of the Sickle H ii region (Timmermann et al. Timmermann96 (1996); Simpson et al. Simpson97 (1997)). Furthermore, there is increasing evidence that the Pistol Nebula is also ionized by the Quintuplet instead of by the Pistol star, located in the center of curvature of the nebula (Moneti et al. Moneti99 (1999); Figer et al. Figer99 (1999)).

On the other hand, the Object 17 in the infrared survey undertaken by Nagata et al. (Nagata93 (1993)) was found to be a cluster of young stars (Nagata et al. Nagata95 (1995); Cotera et al. Cotera96 (1996)). This cluster, also known as G0.121+0.017 or the Arches cluster have enough OB stars to be the ionization source of the Thermal Filaments (Cotera et al. Cotera96 (1996); Serabyn et al. Serabyn98 (1998)). However, to explain the homogeneity in the excitation observed in the fine structure line emission, Colgan et al. (Colgan96 (1996)) proposed that more hot stars should be distributed rather uniformly through the region of the Arched Filaments. Nonetheless, the recent analysis of fine structure lines observed by ISO (Cotera et al. 2000a ) and the recombination lines VLA observations of Lang et al. (Lang01 (2001)), suggest that the Arches cluster alone could account for the ionization of the Arched Filaments.

Recently, mid-infrared images (7-25 $`\mu `$m) of the GCR taken by the Midcourse Space Experiment (MSX) satellite have shown a new morphological feature in the Radio Arc region (Egan et al. Egan98 (1998) ). It is a ring of warm dust which is also present, although it is rather weak, in the 20 cm radio image of Yusef-Zadeh & Morris (1987c ). Although it is not yet clear its actual morphology (bubble, ring or cylinder), in the following we will call it the Radio Arc Bubble (RAB). The RAB is also seen in the IRAS HIRES maps (see Simpson et al. Simpson97 (1997)) where it appears as a ring of warm dust emitting mainly at 25 and 60 $`\mu `$m. The PAHs emission detected by ISO in the RAB shows the presence of carbon-rich material in the bubble (Levine et al. Levine99 (1999)). The part of the Sickle which is perpendicular to the NTFs just lie in the north-west edge of the RAB, suggesting that the Sickle is part of this large structure of warm dust. The geometrical center of the RAB is located just on the most intense NTF at $`(l,b)(0\stackrel{}{.}16,0\stackrel{}{.}11)`$, north-west of the non-thermal maximum (G0.16-0.15).

Figer et al. (Figer99 (1999)) suggested that winds and radiation pressure from O-type stars in the Quintuplet cluster or even a supernova explosion cloud have blown the dust bubble. They proposed that the offset of the Quintuplet cluster from the center of the RAB could be due to the relative motion of the stellar cluster with respect to the ISM or to large-scale density fluctuations in the medium where the bubble is expanding. However, the nearly circular shape of the RAB seen in the radio images of Yusef-Zadeh & Morris (1987c ) rules out the hypothesis of a density gradient.

The purpose of this paper is to test the scenario of a very inhomogeneous medium —proposed by Martín-Pintado et al. (2000a ) to explain the ionization of M+0.12–0.12 and M+0.16–0.10— by studying the large scale ionization structure in the region of the Radio Arc. We have analyzed selected fine structure lines of the observations taken with the Infrared Space Observatory (ISO; Kessler et al. Kessler96 (1996)) towards the Radio Arc. Previous works were restricted to the Sickle and the Pistol Nebula. We study the Sickle in the broader context of the RAB and investigate the possible connection between the ionization in the Radio Arc area and the existence of the RAB. The paper is organized as follows: Section 2 defines the sample of sources selected from the ISO Data Archive (IDA) and Sect. 3 comments on the data reduction and results. The analysis of the results is presented in Sect. 4 and, with the help of photoionization models, in Sect. 5. The results are discussed in Sect. 6 and 7. Finally, Sect. 8 summarizes the main conclusions.

## 2 Observational data

To study the large scale ionization structure of the Radio Arc region and its possible connection to the ring or bubble of warm dust (hereafter RAB) we have used Short Wavelength Spectrometer (SWS; de Graauw et al. deGraauw96 (1996)) and Long Wavelength Spectrometer (LWS; Clegg et al. Clegg96 (1996)) ISO data of our sources M+0.16-0.10 and M+0.21-0.12. We have also retrieved other LWS and SWS observations towards this region from the ISO Data Archive (IDA). We have analyzed the Ne iii 15.6 $`\mu `$m, Ne ii 12.8 $`\mu `$m and the S iii 33.5 and 18.7 $`\mu `$m fine structure lines observed with the SWS. In the LWS observations we have analyzed the fine structure lines of N iii at 57 $`\mu `$m, N ii at 122 $`\mu `$m and O iii at 52 and 88 $`\mu `$m. Table 1 lists IDA identification numbers (TDTs), coordinates and observing modes (AOTs), for all the selected sources. Figure 1 shows the position of the data points overlaid on the 18.2-25.1 $`\mu `$m band MSX image in which the RAB is clearly visible.

Positions M1 (M+0.16-0.10) and M2 (M+0.21-0.12) are the two clouds from our sample located in the Radio Arc region (Martín-Pintado et al. 2000a ; Rodríguez-Fernández et al. RF01 (2001)). M1 lies close to the geometrical center of the RAB over the most intense NTF. Positions B1 to B3 lie in the southern part of the RAB (see Levine et al. 1999) while B4 and B5 are located in the northeast edge of the RAB, close to M2. Q1 and Q2 are the Quintuplet sources GCS3 II and GCS4 (see Chiar et al. Chiar00 (2000)), S1 and S2 are two observations toward the Sickle and P toward the Pistol Nebula. R is a raster map of 4$`\times `$9 points (4 rows, 9 columns), with a spacing of 50$`^{^{\prime \prime }}`$ and rotated 50° north from west. In the following, individual points in the raster will be named as R$`xy`$, i.e., R11, R12,…

## 3 Data reduction and results

The observations consisted of SWS full grating spectra (SWS01), SWS grating spectra of selected lines (SWS02), and LWS full grating spectra (LWS01). The telescope aperture at the wavelengths of the lines of interest are listed in Table 2 and Table 3. M1 and M2 have been processed through the Off Line Processing (OLP) version 7.0. The other observations have been automatically reprocessed when retrieval from the IDA with OLP versions 9.1 to 9.5. Further reduction and analysis have been carried out with ISAP 2.0.

### 3.1 SWS AOTs 01 and 02

As a first stage, we have removed bad data points. Afterwards, for the AOT 01 observations, we have shifted the different detectors to a common level and we have averaged across scan directions and detectors. For the AOT 02 observations we have averaged across scan directions before shifting the different detectors to a common level and averaging across detectors. Baselines of order 1 were removed from the spectra. The fluxes of the lines with their corresponding rms errors as derived from Gaussian fits are listed in Table 2. Errors in the line fluxes due to calibration uncertainties are smaller than 25$`\%`$ for the Ne iii and Ne ii lines and less than 20$`\%`$ and 30$`\%`$ for the S iii 18.7 and 33.4 $`\mu `$m, respectively (Leech et al. 2001).

### 3.2 LWS AOT 01

For all the LWS01 observations, detector SW2 showed memory effects at wavelengths shorter than 55 $`\mu `$m. To analyze the O iii 52 $`\mu `$m line we have treated the two scan directions independently, averaging only across scans. For both scan directions the fluxes are in agreement within 15 $`\%`$. The fluxes of the O iii 52 $`\mu `$m lines listed in Table 3 are the average of the fluxes derived independently for the two scans directions. For the lines at other wavelengths we have shifted the different scans to a common level before averaging across the two scan directions and across scans. We have subtracted from the spectra baselines of order 1. The line fluxes and rms errors as derived from the Gaussian fits are listed in Table 3.

The O iii 88 $`\mu `$m line is in the overlap region of detectors SW5 and LW1. We have analyzed both detectors separately and the measured fluxes agree within 15$`\%`$. The fluxes listed in Table 3 are the average values of both determinations. Before fitting Gaussian profiles to the O iii 88 $`\mu `$m (detector LW1) and the N ii 122 $`\mu `$m (detector LW3) lines we have defringed the spectra. The line flux calibration uncertainties are expected to be smaller than 30$`\%`$ (Swinyard et al. 1996).

The raster map, R, has been reduced using the same procedure described above and Fig. 2 shows the emission maps for all the lines. The maximum of the emission is found towards the Sickle for all but the N ii 122 $`\mu `$m line, which peaks in the western part of the map where the intensity of the N iii line also increases. The emission, however, does not completely follow the Sickle. The part of the Sickle parallel to the NTFs that points toward the center of the RAB is not clearly seen in any of the emission maps. The N iii map also shows a local maximum at the position of the Pistol Nebula.

### 3.3 Comparison with previous observations

Some observations presented in this paper have also been previously carried out with other instruments. One can compare the LWS O iii 52 and 88 $`\mu `$m and N iii 57 $`\mu `$m lines fluxes in the Sickle region with those measured by Simpson et al. (Simpson97 (1997)). Their points P2, P5, P3, P4, P8 and P6 are included in the LWS raster map. Assuming extended homogeneous emission and taking into account the differences in the beams sizes of both observations, one finds that for the sources located in the Sickle (P2, P5, P3 and P4) the fluxes measured by Simpson et al. are 60$`\%`$ higher than those obtained with the LWS. On the other hand, for P8 and P6 they are rather similar. These results indicate that the emission is indeed quite homogeneous in the area of their P8 and P6 positions (approximately corresponding to R34 and R28 respectively) but not in the region of the Sickle. However, even in the Sickle region, the clumpiness of the fine structure line emission in the LWS beam is relatively small and the approximation of extended and homogeneous emission is rather good (see also Martín-Pintado et al. 2000a). A similar conclusion is derived if one compares our O iii 88 $`\mu `$m and N iii 57 $`\mu `$m data with the observations of the Sickle by Timmermann et al. (Timmermann96 (1996)). Taking into account that they had a beam of 22$`^{^{\prime \prime }}`$ and assuming homogeneous emission, the expected flux in the LWS beam ($`80^{^{\prime \prime }}`$) is larger than that observed for some points and lower for other points being on average 17 $`\%`$ larger than those measured with the LWS. Given the different beam sizes and the flux calibration uncertainties, we conclude that the LWS data are in good agreement with previous data and that the fine structure lines observed with the LWS show a relatively uniform emission.

## 4 Analysis

### 4.1 Extinction correction

Before undertaking further analysis of the data it is necessary to correct them for the dust extinction effects. For the typical average visual extinction of $`25`$ mag measured towards the GCR (Catchpole et al. 1990; Schultheis et al. 1999) even at mid and far IR wavelengths the extinction corrections are non negligible. We can use different methods to estimate the extinction affecting our data. One possibility is to use the ratio of the S iii lines at 18.7 and 33.5 $`\mu `$m to derive lower limits to the extinction for the sources where these two lines are available. This is because the S iii 18.7/33.5 ratio cannot be lower than the theoretical limit of 0.5 for the case of very low density. If the measured ratio is lower than the predicted lower limit, this indicates a higher extinction for the S iii 18.7 $`\mu `$m than for the S iii 33.5 $`\mu `$m line. From our S iii data, we derived a lower limit to the visual extinction assuming extended emission for all the positions of $`20`$ mag at most.

A direct determination for the overall extinction in the area of the Sickle has been recently made by Cotera et al. (2000b ) both by analyzing the stellar content and by comparing the radiocontinuum emission with the Br$`\gamma `$ emission. They found an extinction at 2.2 $`\mu `$m of $`2.53.2`$ mag. These values imply a visual extinction of $`2532`$ mag, in agreement with the large scale studies of Catchpole et al. and Schulteis et al. and consistent with our lower limits. In the following we will assume that our data is affected by an average visual extinction of 30 mag. To extrapolate the visual extinction to the wavelengths of the SWS lines we have used the extinction law derived by Lutz (Lutz99 (1999)) towards Sgr A. This extinction law is almost equivalent to that of Draine (Draine89 (1989)) at the wavelengths of the lines presented in this paper. To derive the extinction for the LWS lines we have used the following extinction law:

$$\mathrm{A}_\lambda =0.014A_\mathrm{V}(30/\lambda )^{1.5}$$
(1)
where $`\lambda `$ is the wavelength in microns, $`\mathrm{A}_\lambda `$ is the extinction at that wavelength, and $`A_\mathrm{V}`$ is the visual extinction.

### 4.2 Electron densities

Table 4 lists the extinction corrected O iii 52 to 88 $`\mu `$m and S iii 18.7 to 33.5 $`\mu `$m line intensity ratios (hereafter O iii 52/88 and S iii 18/33 ratios, respectively) and Fig. 3b shows the map of the extinction corrected O iii 52/88 ratio for the raster. The most prominent feature in this map is a ridge that goes northeast-southwest at the position of the Pistol Nebula. The O iii 52/88 $`\mu `$m ratio in the map ranges between 0.78 and 1.43. Table 4 also list the electron densities ($`n_\mathrm{e}`$) derived from these ratios assuming that the lines are excited by collisions with electrons and that the sources are extended (see Rubin et al. Rubin94 (1994)). The derived O iii 52/88 ratio implies $`n_\mathrm{e}`$ between 101.8 and 102.6 cm-3. From the S iii 18/33 ratio we can only derive an upper limit to $`n_\mathrm{e}`$ of $``$100 cm-3 for most sources. For the Pistol Nebula and the Quintuplet sources (P, Q1 and Q2) the derived densities from the S iii lines are as high as $`10^{3.23.8}`$ cm-3. These densities are higher than those derived for the same region from the O iii lines (with a larger beam), suggesting the presence of small scale structures with high density in the Sickle and the Pistol Nebula as previously mentioned in Sect. 3.3.

### 4.3 Ionization structure

Table 4 lists the extinction corrected N iii 57 to N ii 122 $`\mu `$m and Ne iii 15.5 to Ne ii 12.8 $`\mu `$m ratios (hereafter N iii/N ii and Ne iii/Ne ii ratios, respectively). The extinction corrected map of the N iii/N ii ratio corresponding to the raster is shown in Fig. 3a. The circles in Fig. 3a indicate the location of the Quintuplet, the triangle indicates the location of the Pistol Nebula and the square indicates the position of the Arches cluster. The N iii/N ii ratio ranges from 0.33 to 2.52. For the range of $`n_\mathrm{e}`$ derived from the O iii 52/88 ratio, the ratio of emissivities in the two nitrogen lines ($`ϵ_{57}/ϵ_{122}`$) only changes from $`6.3`$ to $`7.7`$ (Rubin et al. 1994). Thus, in spite of the relatively low critical density of the N ii line, the N iii/N ii ratio basically traces changes in the ionization fraction of N iii relative to N ii and not density variations. Although the N iii/N ii ratio shows a local maximum at the location of the Sickle, the absolute maximum is found towards the Pistol Nebula, close to the Quintuplet (although it is elongated towards the southwest). The N iii/N ii ratio shows a general trend to decrease as the distance to the Quintuplet cluster increases. This is illustrated in Fig. 4, where we show all the N iii/N ii ratios as a function of their projected distance to the Quintuplet cluster. The N iii/N ii ratio shows a considerable dispersion which could be due to the difference between projected and actual distances and the possible minor contribution of other ionizing sources. It is remarkable that the degree of ionization toward sources like M1, M2, B4 and B5 located in the RAB, to the south and east of the Quintuplet, is similar to the ionization toward sources located in other directions but at similar distances from the cluster. This clearly points to the Quintuplet as the main ionization source of the RAB and not only of the Sickle and the Pistol Nebula.

The northwest part of the N iii/N ii map shows another gradient that does not point to the Quintuplet but towards other ionization source located in the direction of the Arches cluster. This cluster is believed to be the source of the ionization of the Arched Filaments (Colgan et al. 1996; Cotera et al. 2000a ; Lang et al. Lang01 (2001)). In fact, the western-most part of the raster map covers the northern part of the E1 arched filament (for notation see for instance Lang et al. Lang01 (2001)). Figure 4 also shows this effect. The N iii/N ii ratio measured for the raster points at distances larger than 0$`\stackrel{}{.}`$05 show a systematic trend. For a given column $`y`$ of the raster, the N iii/N ii ratio decreases from the points in the first row of the raster (R1$`y`$) to the points in the fourth column (R4$`y`$). The effect is even more clear when one plots the N iii/N ii ratio derived for the first four columns of the raster as a function of their projected distance to the Arches cluster (see Fig. 5).

As for the N iii/N ii ratio, the maximum of the Ne iii/Ne ii ratio is also found in the Quintuplet and Pistol Nebula sources. Figure 6 displays the Ne iii/Ne ii ratio as a function of projected distance to the Quintuplet. Although the number of available measurements is less than those of the nitrogen lines, the general behavior of the Ne iii/Ne ii ratio is similar to that of the N iii/N ii ratio, decreasing as distance to the Quintuplet increases. However, the gradient in the Ne iii/Ne ii ratio with distance is smaller than that exhibited by the N iii/N ii ratio ($`3`$ instead of $`10`$). The critical densities of both the Ne iii and the Ne ii lines are $`\stackrel{>}{}`$$`\mathrm{2\hspace{0.17em}10}^5`$ cm-3. Thus, for the range of $`n_\mathrm{e}`$ derived for all the sources, the emissivities of both lines do not depend on $`n_\mathrm{e}`$, and like the N iii/N ii ratio, the Ne iii/Ne ii ratio traces the ionization structure.

## 5 Origin of the ionization

To study the large scale ionization structure we have compared the trend observed in the N iii/N ii and Ne iii/Ne ii ratios with the predictions of photoionization models considering the Quintuplet and the Arches clusters as the main ionizing sources. We have used CLOUDY version 90.04 (Ferland Ferland96 (1996)) and the MICE interface developed by H. Spoon at the MPE to compute models for a grid of physical parameters. Namely, effective temperatures, $`T_{\mathrm{eff}}`$, and Lyman continuum photons fluxes, $`Q`$(H), of the ionizing sources. To define the shape of the incident continuum radiation we have used the models of stellar atmospheres developed by Kurucz (1994) and the Schaerer & de Koter (Schaerer97 (1997)) “CoStar” models, both for solar abundances. We have assumed nebular N and Ne abundances of 3 $`10^4`$ and 2.8 $`10^4`$, respectively, as derived by Rubin et al. (1995) for H ii regions in the GCR. These abundances are a factor of 2-3 higher than those measured for Galactic disk H ii regions. We have checked with the Kurucz atmospheres that increasing the stellar metallicities in a factor of 3 changes the predicted lines ratios in less than 15$`\%`$. Therefore, enhanced metal abundances in the atmospheres does not affect significantly the results presented in the following sections. As a typical density of the ionized material we have used $`10^{2.2}`$ cm-3 as derived from the O iii lines ratio.

### 5.1 The Quintuplet cluster

Since the Quintuplet cluster seems to dominate the ionization for most of the observed positions, we first considered this cluster as the ionizing source. We have used CLOUDY to estimate the size of the ionized region for the case in which the ionizing source is surrounded by material with an uniform $`n_\mathrm{e}`$ of 102.2 cm-3. We assumed $`Q`$(H)=1050.9 s-1 as derived by Figer et al. (1999) for the Quintuplet cluster and a $`T_{\mathrm{eff}}`$ of 36300 K (B2 model of Schaerer & de Koter). The predictions of CLOUDY for the size of the ionized region around the Quintuplet is of $`\stackrel{<}{}`$6 pc, much smaller than what is observed from the fine structure line ratios (which is $`18`$ pc assuming a distance to the Galactic center of 8.5 kpc). One can, however, explain the presence of Ne iii and N iii at such large distances from the Quintuplet with the $`T_{\mathrm{eff}}`$ and $`Q`$(H) derived for the Quintuplet if the material between the Quintuplet an the observed positions has a density of $``$ 10 cm-3. This is consistent with the ring, or almost empty bubble, morphology observed in the MSX images.

We have therefore used a simple model in which the ionized lines observed in our beam arise from clouds of density 102.2 cm-3 located at a distance to the Quintuplet cluster equal to their projected distance on the sky (we assume a distance the GCR of 8.5 kpc). Under these conditions, one can characterize the effect of the radiation from a ionizing source on a cloud located at a distance $`D`$ using an ionization parameter ($`U`$) defined as:

$$U=\frac{Q(\mathrm{H})}{4\pi D^2n_ec}$$
(2)
where $`c`$ is the velocity of light.

Following the models used by Simpson et al. (1997) for the Sickle and the Pistol Nebula, we have first computed the N iii/N ii and Ne iii/Ne ii ratios as a function of the distance of the clouds to the cluster (i.e. as a function of $`U`$) for $`T_{\mathrm{eff}}`$ in the range of 34000-38000 K and $`Q`$(H) in the range of 1049.5-50.5 s-1 with Kurucz (1994) atmospheres. We found that $`Q`$(H)=1049.5 s-1 produces too small N iii/N ii ratios at large distances. On the opposite, $`Q`$(H)=1050.5 s-1 can explain the derived ratios if $`T_{\mathrm{eff}}`$ is in the range of 35000-36000 K (dot-dashed lines in Fig. 4). These values are in good agreement with those of Simpson et al. (Simpson97 (1997)) and indicate that our results extend to larger distances the results of Simpson et al. derived for the Sickle (see also Philipp et al. in prep). However, this model with the Kurucz atmospheres predicts Ne iii/Ne ii ratios which are a factor of $`10`$ smaller than the measured values. This is likely related to the well known Ne iii/Ne ii problem associated to the Kurucz atmospheres, which predict Ne iii/Ne ii ratios that are always much smaller than the measured ones. More recent model atmospheres give larger Ne iii/Ne ii ratios, much closer to the observed values (see e.g. Sellmaier et al. 1996).

We have tried to explain both the Ne iii/Ne ii and the N iii/N ii ratios consistently with similar parameters for the ionizing source using CLOUDY models with $`Q`$(H)=1050.9 s-1 using the stellar atmospheres of Schaerer & de Koter. Figure 4 shows, as dashed lines, that the decrease of the N iii/N ii ratio with the distance is well explained with $`T_{\mathrm{eff}}`$=32600-33300 K (B3 and A2 CoStar models). The discrepancies between the observational results and the model are within a factor of $`1.5`$. This is quite reasonable in view of the simple model we use and the uncertainties introduced by its approximations. These are, mainly, the assumptions that projected distances are actual distances to the cluster and that there is not material in between the cluster and the observed sources (however, the model takes into account the attenuation of the radiation within the clouds). In addition, the model uses a single value of the $`n_\mathrm{e}`$ for all the positions. Figure 6 shows, also with dashed lines, that the predictions of this model for the Ne iii/Ne ii ratio are in better agreement with the ratios derived from the observations than those derived using the Kurucz atmospheres. However, the high $`Q`$(H) value used to explain the N iii/N ii ratios predicts curves for the Ne iii/Ne ii ratio as a function of distance that are rather flat. Moreover, the minimum $`T_{\mathrm{eff}}`$ of $`34500`$ K required in order to explain the Ne iii/Ne ii ratio is somewhat higher than the $`T_{\mathrm{eff}}`$ ($`33000`$ K) required to explain the N iii/N ii ratio.

For completeness, Fig. 4 and 6 also show, as dotted lines, the predictions for model calculations using CoStar atmospheres but with a slightly lower $`Q`$(H) of 1050.5 s-1. Figure 4 shows that we can reproduce the observed behavior of the N iii/N ii ratio for $`T_{\mathrm{eff}}`$ in the range of 33300 and 35500 K. On the other hand, Fig. 6 shows that we also have a rather good agreement between the predictions and the measured Ne iii/Ne ii ratio for $`T_{\mathrm{eff}}`$$`3450035500`$ K, which is similar to the $`T_{\mathrm{eff}}`$ needed to explain the N iii/N ii ratio. In summary, the N iii/N ii and Ne iii/Ne ii ratios observed with ISO in the RAB region can be explained by assuming that the fine structure lines arise from clouds located at different distances from the Quintuplet cluster, which would be the main ionizing source. The required Lyman continuum photons flux of $`10^{50.550.9}`$ s-1 and effective temperatures of $`35000`$ K are in good agreement with the estimations based on the stellar content of the cluster (see Figer et al. 1999).

### 5.2 The Arches cluster

As mentioned in section 4.3 the observed gradient in the N iii/N ii ratio in the northwest region of the raster cannot be explained by the effect of the Quintuplet. In this section we will discuss the influence of the Arches cluster in the ionization of the Radio Arc material. Figure 5 shows the N iii/N ii ratio derived for the positions of the raster map which are closer to the Arches cluster as a function of their projected distance to this cluster. Figure 5 also displays, as long-dashed lines, the results of CLOUDY models computed with CoStar atmospheres, $`Q`$(H)=1051.4 s-1 and $`n_\mathrm{e}`$=102.2 cm-3. The observed gradient in the N iii/N ii ratio is pretty well explained with a $`T_{\mathrm{eff}}`$ in the range of 30400-32600 K. These parameters are in agreement with those estimated for the Arches cluster (Morris et al. in prep see also Lang et al. 2001). Therefore, we conclude that the observed N iii/N ii gradient to the southwest is likely produced by the Arches cluster and that the influence of the radiation from this cluster should be considered up to distances of 7-13 pc.

### 5.3 The combined effect of both clusters

We now consider the combined effects of the two clusters on the ionization of the whole region surrounding them. As discussed in Sect. 5.1, the N iii/N ii and Ne iii/Ne ii ratios are directly related to the ionization parameter $`U`$. Since the effective temperature of the radiation derived for both clusters are rather similar (see above), we have considered the most simple model to account for the combined effect of the two clusters.

As a first approximation one can estimate a “total ionization parameter” ($`U_{\mathrm{tot}}`$) as the sum of two different ionization parameters, one due to the Quintuplet ($`U_\mathrm{Q}`$) and other due to the Arches cluster ($`U_\mathrm{A}`$). Thus, for a cloud at distances $`D_\mathrm{Q}`$ and $`D_\mathrm{A}`$ from the Quintuplet and the Arches cluster, respectively, $`U_{\mathrm{tot}}`$ will be given by:

$$U_{\mathrm{tot}}=U_\mathrm{Q}+U_\mathrm{A}=\frac{Q(\mathrm{H})_\mathrm{Q}}{4\pi D_\mathrm{Q}^2n_ec}+\frac{Q(\mathrm{H})_\mathrm{A}}{4\pi D_\mathrm{A}^2n_ec}$$
(3)
We have plotted in the upper panel of Fig. 7 some contours of equal $`U_{\mathrm{tot}}`$ assuming a constant density of $`n_\mathrm{e}`$=102.2 cm-3, a $`Q(\mathrm{H})_\mathrm{Q}`$ of $`10^{50.9}`$ s-1 and a $`Q(\mathrm{H})_\mathrm{A}`$ of $`10^{51.4}`$ s-1. The agreement of the iso-$`U_{\mathrm{tot}}`$ curves with the contour map of the N iii/N ii ratio (shown with thick lines) is very good taking into account the simplicity of the model, which for instance, does not consider any shielding of the radiation. It reproduces the observed gradients towards the Arches and Quintuplet clusters. However, the model does not account for the local maxima found in the Sickle. The lower N iii/N ii ratios derived towards points located southward (in galactic latitude) of the Quintuplet which lie on the same iso-$`U_{\mathrm{tot}}`$ curve than others located northward of the cluster should be due to the fact that $`U_{\mathrm{tot}}`$ is more homogeneous in the LWS beam in the north of the cluster than in the south.

Furthermore, this simple model also reproduces the observed overall distribution of the warm dust. This is also illustrated in the upper panel of Fig. 7, where the iso-$`U_{\mathrm{tot}}`$ contours are superposed on the 25 $`\mu `$m band image taken by MSX. One can see that the morphology of the warm dust emission approximately follows the shape of the iso-$`U_{\mathrm{tot}}`$ curves. The warm dust in the area of the Thermal Filaments is under the strong influence of the Arches cluster while the warm dust emission in the Sickle follows the shape of the iso-$`U_{\mathrm{tot}}`$ curves that are dominated by the Quintuplet. The elliptical rather than circular distribution of the warm dust around the clusters, in particular around the Arches cluster, can be explained by the combined effect of both clusters. Even at large distances from the clusters ($`\stackrel{>}{}`$15 pc) the distribution of the warm dust also seems to follow the shape of the iso-$`U_{\mathrm{tot}}`$ curves, see for instance in Fig. 7 the 25 $`\mu `$m dust emission in $`(l,b)(0\stackrel{}{.}2,0\stackrel{}{.}15)`$ or in $`(l,b)(0\stackrel{}{.}05,0\stackrel{}{.}12)`$. The agreement between the morphology of the warm dust emission with the iso-$`U_{\mathrm{tot}}`$ curves strongly suggest a same origin for both the warm dust and the ionized gas (see next section).

## 6 The Radio Arc bubble of warm dust

The origin of the ring or bubble of warm dust it is not clear but the ISO data presented in this paper shows that the ionization of the gas is mainly produced by the Quintuplet instead of any source in the center of the ring. One unlikely possibility is that the hot dust observed in the RAB is not heated by the UV radiation from the clusters that explain the ionization structure (see Sect. 5.3). In this case, one would expect not to see a correlation between the warm dust and the ionization. We can check this possibility by comparing for instance the O iii 88 $`\mu `$m line intensity with the continuum emission from warm dust for different positions in the RAB. The intensity of this fine structure line for points “inside the ring” (M1, R19, R29; i.e. where the emission from warm dust is relatively weak, see Fig 1) is lower than that for points like M2, B4 or B5 located on the edge of the ring . For instance, the average O iii 88 $`\mu `$m line fluxes for the points on the ring is $`\mathrm{12\hspace{0.17em}10}^{18}`$ W cm-2 while for those points inside the ring is $`\mathrm{7\hspace{0.17em}10}^{18}`$ W cm-2. From the IRAS HIRES 25 $`\mu `$m map (with a resolution of $`66^{^{\prime \prime }}\times \mathrm{\hspace{0.17em}35}^{^{\prime \prime }}`$) one obtains that the averaged intensity of the dust emission at this wavelength for the points on the edge of the ring, M2, B4 and B5, is $`1500`$ MJy sr-1 and $`750`$ MJy sr-1 for the points inside the ring (M1, R19, R29). Furthermore, towards the Sickle positions (R16-R46) the average O iii 88 $`\mu `$m line flux is $`\mathrm{36\hspace{0.17em}10}^{18}`$ W cm-2 while the 25 $`\mu `$m intensity is $`5000`$ MJy sr-1. Therefore, it seems that the ratio of the luminosity of the fine structure lines (ionized gas column density) to the warm dust emission is approximately constant for all positions in the RAB area. This constant line-to-continuum ratio shows that the warm dust should also be heated by the UV radiation from the cluster. The apparently lack of warm dust in the central part of the ring is just an effect of column density but not of UV radiation since the N iii/N ii ratio is higher for points like M1 inside the RAB, than for points at the edge of the RAB like M2.

The smaller column density of material towards the center of the ring is consistent with different geometries like a ring (or a cylinder with its axe close to the line of sight) and a shell or almost empty bubble. For a shell geometry with constant density and internal and external radius of 0$`\stackrel{}{.}`$05 and 0$`\stackrel{}{.}`$075, respectively, the expected column density in the line of sight to the edge of the shell is a factor of $`2`$ larger than that in the line of sight towards the center of the shell. This is the typical ratio of the O iii 88 $`\mu `$m line or the 25 $`\mu `$m dust emission between M2, B4, B5 at the edge of the shell and M1, R19, R29 at the center. Although this does not demonstrate that the observed structure has a shell morphology instead of a ring one, the data are consistent with a shell morphology (bubble) of dust heated, and gas ionized, by the UV radiation of the Quintuplet and the Arches cluster.

The geometrical center of the bubble is approximately located at $`(l,b)=(0\stackrel{}{.}16,0\stackrel{}{.}11)`$, over the NTFs. We have searched in the SIMBAD database for sources in the central area of the RAB without success. The closest sources to the geometrical center of the bubble are the near-IR  sources 48 and 49 of Nagata et al. (Nagata93 (1993)) and our source M1 but no link between them and the RAB can be established. Somewhat further from the center there are other IR sources detected by Nagata et al. (1993) like Objects 43-46 but they are late type M-stars candidates. However, close to this sources, there is also a luminous X-ray source: 1E 1742.9-2849 (represented as a filled circle in the lower panel of Fig. 7), with an X-ray luminosity of $`\mathrm{8\hspace{0.17em}10}^{34}`$ erg s-1 (Predehl & Trümper Predehl94 (1994), also recently detected by Chandra, Yusef-Zadeh priv. communication). It is $`0\stackrel{}{.}04`$ ($`6`$ pc) far from the center of the shell and its position suggests that it could be responsible for the elongation of the maximum in the N iii/N ii ratio toward the southwest instead of being centered in the Quintuplet and the Pistol Nebula as expected from the models. There is another X-ray source inside the shell at a similar distance from its center: 1E 1743.1-2852 (triangle in the lower panel of Fig. 7).

It is remarkable that the diffuse hard X-rays (2-10 keV) continuum emission observed by the ASCA satellite towards the GCR (Koyama et al. Koyama96 (1996)) shows extended emission which covers an important fraction of the RAB. This would indicate that the RAB might be filled with hot gas that gives rise to the extended X-ray emission. One interesting prediction for the proposed scenario of a bubble filled with hot gas is that the morphology of the 6.4 keV line of neutral Fe, which traces gas irradiated by hard X-rays (that is X-ray dominated regions -XDRs) should trace the walls of the bubble. The lower panel of Fig. 7 shows the Fe 6.4 keV line emission map of Koyama et al. (1996) overlaid on the MSX image. The Fe 6.4 keV line emission shows two maxima in the Sgr A neighborhood. One is located at $`(l,b)(0\stackrel{}{.}04,0\stackrel{}{.}01)`$, contouring the hard X-rays emission of the SgrA region while the other maximum is located at $`(l,b)(0\stackrel{}{.}13,0\stackrel{}{.}11)`$ within the RAB. The Fe 6.4 keV emission in the RAB peaks in between 1E 1742.9-2849 and 1E 1743.1-2852, although closer to the last source. It clearly fills the bubble covering the region between the NTFs at $`l0\stackrel{}{.}18`$ and the edge of the bubble at $`l0\stackrel{}{.}08`$.

In summary, the RAB seems to be filled with hot gas emitting hard X-rays continuum which create a large XDR that is revealed by the Fe 6.4 keV line emission. We speculate with the possibility that the X-ray sources might be related with the energetic events that created the bubbles of warm dust. Relative motions of the ISM and the X-ray source would explain the off-center position of the X-ray sources in the RAB. Considering the shape, the size of the bubble and the presence of the X-rays Fe lines, a supernova explosion is the most likely origin for the RAB.

## 7 General implications

As suggested by Martín-Pintado et al. 2000a, the picture arising from the comparison of the observed ionization and the model calculations is that the ISM in the Radio Arc region is highly inhomogeneous as indicated by the presence of bubbles filled with relatively low density ionized gas and hot gas emitting in X-rays. This structure allows the radiation from hot stars to reach large distances giving rise to extended low density ionized regions. The Quintuplet and the Arches cluster are the dominant sources of ionization in the Radio Arc region. Alternative ionization mechanisms (for instance magnetic effects) would play a minor role and should be restricted to smaller regions such as the Sickle.

This inhomogeneous structure seems to be a typical characteristic of the GCR ISM. In fact, this kind of structures has been invoked to explain the non-detection of hydrogen recombination lines and the detection of fine-structure lines from ions like Ne ii, S iii or O iii in a sample of clouds distributed all along the GCR (Martín-Pintado et al. 2000a ; Rodríguez-Fernández et al. in prep.). On the other hand, inhomogeneities like shells and filaments are known to exist in other regions of the GCR like those revealed by the CO survey of Oka et al. (1998) or those found by Martín-Pintado et al. (MP99 (1999)) in the envelope of Sgr B2 which, like for the RAB, also seem to be correlated with the Fe 6.4 keV line emission (see Martín-Pintado et al. 2000b ).

## 8 Conclusions

We have presented ISO observations of selected fine structure lines in the Radio Arc in the Galactic center region. The analyzed lines are N iii 57 $`\mu `$m, N ii 122 $`\mu `$m, O iii 52 and 88 $`\mu `$m, Ne iii 15.6 $`\mu `$m, Ne ii 12.8 $`\mu `$m and S iii 18.7 and 33.5 $`\mu `$m. The main results derived from these observations can be summarized as follows:

* The electronic densities derived from the O iii 52 to 88 $`\mu `$m ratio range between $`10^{1.8}`$ and $`10^{2.6}`$ cm-3. The electronic densities derived from the S iii 18.7 to 33.5 $`\mu `$m ratio in the smaller beam of the SWS are $`10^2`$ cm-3 for all the sources but those in the Quintuplet and Pistol Nebula which reach values of $`10^{3.5}`$ cm-3.
* After correcting for extinction, the N iii/N ii and the Ne iii/Ne ii ratios vary from position to position from 0.3 to 2.5 and from 0.05 to 0.30, respectively. Both ratios exhibit a clear trend to decrease with increasing distance to the Quintuplet cluster. For some observed positions, the N iii/N ii ratio also shows a gradient that points towards the Arches cluster.
* Photoionization models confirm that the observed ratios and their trend with the distance to the clusters are consistent with the derived parameters for the clusters, that is $`T_{\mathrm{eff}}`$=32000–35000 K, and a flux of Lyman continuum photons of $`10^{51.4}`$ s-1 and $`10^{50.9}`$ s-1 for the Arches and the Quintuplet cluster, respectively. Thus, the combined effects of both the Quintuplet and the Arches cluster completely dominate the ionization in the Radio Arc region. They create a large ionized region with a size of more than $`30\times 30`$ pc2. Alternative ionization mechanisms (for instance magnetic effects or more hot stars) seems to play a minor role at large scales.
* A simple model of a total ionization parameter expected for the two clusters can explain not only the ionization structure but also the distribution of warm dust and, in particular, the elliptical rather than circular symmetry of some features like the Thermal Filaments. The warm dust seems to be correlated with the ionized gas and therefore also seems to be heated by the two clusters that ionize the gas.
* The effect of the radiation from the Quintuplet stars over large distances is at least in part due to the presence of a shell or an almost empty bubble (Radio Arc Bubble, RAB) of warm dust and ionized gas.
* The RAB is filled with hot gas emitting hard X-rays which creates a large X-ray dominated region (XDR) observed in the Fe 6.4 keV line emission.
* The origin of the bubble it is not clear, but it could be related to a couple X-ray sources located $`6`$ pc off the shell center. Taking into account the size and shape of the bubble and the presence of the Fe X-rays lines, a supernova origin is certainly possible.

###### Acknowledgements.

NJR-F acknowledges Consejería de Educación de la Comunidad de Madrid for a predoctoral fellowship. The authors acknowledge support by the Ministerio de Ciencia y Tecnología under grants 1FD97-1442 and PNE 014-2000-C. NJR-F thanks H. Spoon for his introduction to CLOUDY and MICE. MICE, SWS and the ISO Spectrometer Data Center at MPE are supported by DLR (DARA) under grants 50 QI 86108 and 50 QI 94023. The ISO Spectral Analysis Package (ISAP) is a joint development by the LWS and SWS Instrument Teams and Data Centers. Contributing institutes are CESR, IAS, IPAC, MPE, RAL and SRON. This research has made use of the SIMBAD database, operated at CDS, Strasbourg, France. This research made use of data products from the Midcourse Space Experiment. Processing of the data was funded by the Ballistic Missile Defense Organization with additional support from NASA Office of Space Science. The data were accessed by services provided by the NASA/IPAC Infrared Science Archive.


format: text

text:
# Finite size effects on the collective oscillations of a trapped Bose gas

## I Introduction

After the experimental realization of Bose-Einstein condensation (BEC) in atomic trapped alkali gases  the study of their dynamical properties has become an important subject of experimental and theoretical research (see  and references therein). The measurements of the collective frequencies have confirmed the validity of the mean field approach  in describing the dynamical behaviour of a trapped Bose-Einstein condensed gas at zero temperature, properly accounting for the important role played by two-body interaction and by Bose statistics. The high precision of these frequency measurements is stimulating new perspectives in the study of finer effects like the temperature dependence of the collective frequencies, the occurrence of beyond mean field corrections , etc. In this perspective it is important to have a full analytical control of the values of the collective frequencies predicted by the mean field scheme at $`T=0`$, including finite size corrections which are neglected in the Thomas-Fermi limit. The purpose of the present work is to calculate these corrections in the framework of the mean field Gross-Pitaevskii theory, employing a sum rule approach . In general sum rules provide only an upper bound to the frequency of the lowest modes. In this paper we will show that, due to the harmonic nature of the confinement, the first corrections to the Thomas-Fermi limit can be calculated with logarithmic accuracy using this method. We will show that the corrections are fixed by the kinetic energy term which, in these trapped Bose-Einstein condensed systems, plays the role of an effective surface tension .

## II Sum rules and collective modes

The usefulness of the sum rule technique in studying the dynamical behaviour of trapped Bose-Einstein condensates has been already discussed in several papers (see, for example, ). This method, typical of many-body physics, is based on the calculation of the energy-weighted moments $`m_p=𝑑\omega S_{}(\omega )(\mathrm{}\omega )^p`$ of the dynamic structure function associated with a given operator $``$. The dynamic structure function, at zero temperature, takes the form

$$S_{}(\omega )=\underset{n}{}|n||0|^2\delta (\mathrm{}\omega \mathrm{}\omega _{n0}),$$
(1)
where $`\mathrm{}\omega _{n0}=\mathrm{}\omega _n\mathrm{}\omega _0`$ is the excitation energy of the eigenstate $`|n`$ of the Hamiltonian, whose ground state is $`|0`$. At $`T=0`$ the quantity $`S_{}(\omega )`$ is related to the dynamic response function of the system relative to $``$ through the expression

$$\chi _{}(\omega )=_{\mathrm{}}^{\mathrm{}}𝑑\omega \left[\frac{S_{}(\omega ^{})}{\omega \omega ^{}+i\eta }\frac{S_{}(\omega ^{})}{\omega +\omega ^{}+i\eta }\right].$$
(2)
In this work we make the following choice for the Hamiltonian

$$H=\frac{1}{2M}\underset{i=1}{\overset{N}{}}p_i^2+\frac{M}{2}\underset{i=1}{\overset{N}{}}\left(\omega _x^2x_i^2+\omega _y^2y_i^2+\omega _z^2z_i^2\right)+g\underset{i<j}{}\delta (𝐫_i𝐫_j),$$
(3)
where the parameter $`g=4\pi \mathrm{}^2a/M`$ is the strength of the two-body contact interaction, fixed by the $`s`$-wave scattering length $`a`$. Ratios between different sum rules provide rigorous upper bounds

$$\mathrm{}\omega _{\mathrm{lowest}}\sqrt{m_{p+2}/m_p}$$
(4)
to the energy of the lowest state excited by the operator $``$. These bounds coincide with $`\omega _{\mathrm{lowest}}`$ only if the strength of the operator $``$, defined by Eq. (1), is exhausted by a single mode.

By using the closure relationship involving the eigenstates $`|n`$, the energy-weighted and the cubic moments can be easily rewritten in the form of sum rules as follows

$`m_1`$ $`={\displaystyle \frac{1}{2}}[^{},[H,]]`$ (5)
$`m_3`$ $`={\displaystyle \frac{1}{2}}[[^{},H],[H,[H,]]],`$ (6)
where the averages are taken over the ground state of the system which, in the mean field scheme, is provided by the stationary solution of the Gross-Pitaevskii equation .

Let us first consider the simplest situation where the harmonic confinement is isotropic, i.e. $`\omega _x=\omega _y=\omega _z=\omega _{\mathrm{ho}}`$. In this case the calculation of the above moments is straightforward (see Appendix A) if the operator $``$ coincides with the monopole

$$=\underset{i=1}{\overset{N}{}}r_i^2$$
(7)
where $`r_i^2=x_i^2+y_i^2+z_i^2`$, or quadrupole

$$𝒬=\underset{i=1}{\overset{N}{}}x_iy_i$$
(8)
operators, yielding the following results for the ratio $`\sqrt{m_3/m_1}`$

$`\omega _{_{3,1}}`$ $`=\sqrt{5}\omega _{\mathrm{ho}}\left[1{\displaystyle \frac{E_{\mathrm{kin}}}{5E_{\mathrm{ho}}}}\right]^{1/2}`$ (9)
$`\omega _{𝒬_{3,1}}`$ $`=\sqrt{2}\omega _{\mathrm{ho}}\left[1+{\displaystyle \frac{E_{\mathrm{kin}}}{E_{\mathrm{ho}}}}\right]^{1/2}.`$ (10)
In Eqs. (9)-(10) $`E_{\mathrm{kin}}=_{i=1}^Np_i^2/2M`$ and $`E_{\mathrm{ho}}=M\omega _{\mathrm{ho}}^2_{i=1}^Nr_i^2/2`$ are the expectation values of the kinetic and harmonic oscillator energies respectively, $`E_{\mathrm{kin}}=_{i=1}^N(p_{x_i}^2+p_{y_i}^2)/2M=2E_{\mathrm{kin}}/3`$ and $`E_{\mathrm{ho}}=M\omega _{\mathrm{ho}}^2_{i=1}^N(x_i^2+y_i^2)/2=2E_{\mathrm{ho}}/3`$ are their corresponding radial contributions. In deriving Eq. (9) we have also used the virial identity

$$2E_{\mathrm{kin}}2E_{\mathrm{ho}}+3E_{\mathrm{int}}=0$$
(11)
relating in an exact way the various contributions to the ground state energy

$$E=E_{\mathrm{kin}}+E_{\mathrm{ho}}+E_{\mathrm{int}}$$
(12)
where $`E_{\mathrm{int}}=g_{i<j}\delta (𝐫_i𝐫_j)`$ is the average value of the mean field interaction energy.

Another possibility for calculating the monopole and the quadrupole frequencies makes use of the ratio $`\sqrt{m_1/m_1}`$, which involves the inverse energy weighted moment. In general the following inequality holds

$$\sqrt{m_1/m_1}\sqrt{m_3/m_1}$$
(13)
showing, by comparison with (4), that the lowest excitation frequency is in general better approximated by the ratio $`\sqrt{m_1/m_1}`$ than by $`\sqrt{m_3/m_1}`$. Differently from $`m_1`$ and $`m_3`$, the moment $`m_1`$ cannot be written in terms of commutators, but is related to the static response function via the following general relationship (compressibility sum rule )

$$\chi _{}(0)=2\underset{n}{}\frac{|n||0|^2}{\mathrm{}\omega _{n0}}=2m_1,$$
(14)
which provides a natural rule for calculating $`m_1`$ through the explicit determination of the polarization $`+^{}=\gamma \chi _{}`$ induced in the system by an external static field interacting through the Hamiltonian $`\gamma (+^{})`$.

We present in Appendix B the details of the derivation and report here only the final results :

$`\omega _{_{1,1}}`$ $`=2\omega _{\mathrm{ho}}\left[r^2/\left(r^2{\displaystyle \frac{N}{2}}{\displaystyle \frac{}{N}}r^2\right)\right]^{1/2}`$ (15)
$`\omega _{𝒬_{1,1}}`$ $`=\sqrt{2}\omega _{\mathrm{ho}}\underset{\epsilon 0}{lim}\left({\displaystyle \frac{\delta }{\epsilon }}\right)^{1/2},`$ (16)
where $`\epsilon =(\omega _x^2\omega _y^2)/(\omega _x^2+\omega _y^2)`$ and $`\delta =y^2x^2/y^2+x^2`$ are, respectively, the deformation of the confining potential and of the condensate in the $`x`$-$`y`$ plane. Equations (15) and (16) show that the monopole frequency $`\omega _{_{1,1}}`$ is determined by the $`N`$-dependence of the square radius $`r^2=_{i=1}^Nr_i^2/N`$, while the quadrupole by the dependence of the deformation $`\delta `$ of the cloud on the deformation $`\epsilon `$ of the trap. In (16) we have explicitly taken the limit of a symmetric trap.

It is worth pointing out that the results (9)-(10) and (15)-(16) are valid in general in the framework of the mean field Gross-Pitaevskii theory and apply, in particular, to both repulsive and attractive forces . In particular they reproduce the correct values for the frequencies in both the relevant limits of the non interacting gas, where $`E_{\mathrm{kin}}=E_{\mathrm{ho}}`$, and of the Thomas-Fermi approximation where the kinetic energy is negligible with respect to $`E_{\mathrm{ho}}`$ . In the former case both the monopole and quadrupole frequencies approach the value $`2\omega _{\mathrm{ho}}`$, while in the latter they tend to the analytic hydrodynamic values $`\sqrt{5}\omega _{\mathrm{ho}}`$ and $`\sqrt{2}\omega _{\mathrm{ho}}`$ for the monopole and quadrupole oscillations respectively, as predicted in Ref. . This feature allows us to conclude that in the hydrodynamic regime a single mode exhausts the sum rules $`m_3`$, $`m_1`$ and $`m_1`$. In the following section we will show that, with logarithmic accuracy, the same is true also if one includes the leading finite size corrections to the Thomas-Fermi limit.

## III Finite-size corrections in spherically symmetric systems

In the Thomas-Fermi regime the asymptotic behaviour of the kinetic energy per particle of a harmonically isotropic Bose gas at $`T=0`$ is given by

$$\frac{E_{\mathrm{kin}}}{N}=\frac{5}{2}\frac{\mathrm{}^2}{MR^2}\mathrm{log}\left(\frac{R}{Ca_{\mathrm{ho}}}\right)$$
(17)
where

$$R=(15\frac{Na}{a_{\mathrm{ho}}})^{1/5}a_{\mathrm{ho}}$$
(18)
is the Thomas-Fermi radius of the condensate, with $`a_{\mathrm{ho}}=\sqrt{\mathrm{}/M\omega _{\mathrm{ho}}}`$, and $`C1.3`$ is a numerical factor. Such a behaviour arises from a careful analysis of the shape of the order parameter near the border of the condensate, where the Thomas-Fermi approximation fails and the wave function of the condensate should be calculated with proper accuracy. By using the virial identity (11), and the usual definition $`\mu =E/N`$ for the chemical potential one can easily show that the first corrections to $`E_{\mathrm{ho}}/N`$ and $`E_{\mathrm{int}}/N`$ beyond the Thomas-Fermi approximation are given by:

$`{\displaystyle \frac{E_{\mathrm{ho}}}{N}}`$ $`={\displaystyle \frac{3}{7}}\mu _{\mathrm{TF}}+{\displaystyle \frac{\mathrm{}^2}{MR^2}}\left[\mathrm{log}\left({\displaystyle \frac{R}{Ca_{\mathrm{ho}}}}\right)+{\displaystyle \frac{3}{8}}\right]`$ (19)
$`{\displaystyle \frac{E_{\mathrm{int}}}{N}}`$ $`={\displaystyle \frac{2}{7}}\mu _{\mathrm{TF}}{\displaystyle \frac{\mathrm{}^2}{MR^2}}\left[\mathrm{log}\left({\displaystyle \frac{R}{Ca_{\mathrm{ho}}}}\right){\displaystyle \frac{1}{4}}\right],`$ (20)
where

$$\mu _{\mathrm{TF}}=\frac{1}{2}M\omega _{\mathrm{ho}}^2R^2$$
(21)
is the Thomas-Fermi chemical potential. Note that in deriving the above equations we have consistently taken into account all the terms of order $`1/R^2`$, which renormalize the constant $`C`$ inside the $`\mathrm{log}`$ argument. Eqs. (19)-(20) were also derived in Ref.  by employing a boundary layer theory.

From Eqs. (17) and (19) one finds that the first corrections to the ratio $`E_{\mathrm{kin}}/E_{\mathrm{ho}}`$ can be written as

$$\frac{E_{\mathrm{kin}}}{E_{\mathrm{ho}}}=\frac{35}{3}\left(\frac{a_{\mathrm{ho}}}{R}\right)^4\times \mathrm{log}\left(\frac{R}{Ca_{\mathrm{ho}}}\right)$$
(22)
and becomes smaller and smaller in the large $`R`$, Thomas-Fermi limit. Equations (19) and (17) also reveal that the kinetic energy plays the role of an effective surface energy, providing the first correction to the total energy (12) in the Thomas-Fermi regime.

By using the above results it is immediate to determine the finite-size corrections to the hydrodynamic value of the frequencies predicted by Eqs. (9) and (10). It is interesting also to evaluate these corrections starting from Eqs. (15) and (16). Appendix C is devoted to the technical details of this calculation, which shows that ratios $`\omega _{3,1}`$ and $`\omega _{1,1}`$ do coincide for both the monopole and the quadrupole modes, apart from a different value of the coefficient inside the $`\mathrm{log}`$ argument. The results for the collective frequencies, including the first corrections due to finite size effects, take the form

$`\underset{Ra_{\mathrm{ho}}}{lim}\omega _{}`$ $`=\sqrt{5}\omega _{\mathrm{ho}}\left[1{\displaystyle \frac{7}{6}}\left({\displaystyle \frac{a_{\mathrm{ho}}}{R}}\right)^4\mathrm{log}\left({\displaystyle \frac{R}{C_{}a_{\mathrm{ho}}}}\right)\right]`$ (23)
$`\underset{Ra_{\mathrm{ho}}}{lim}\omega _𝒬`$ $`=\sqrt{2}\omega _{\mathrm{ho}}\left[1+{\displaystyle \frac{35}{6}}\left({\displaystyle \frac{a_{\mathrm{ho}}}{R}}\right)^4\mathrm{log}\left({\displaystyle \frac{R}{C_𝒬a_{\mathrm{ho}}}}\right)\right],`$ (24)
with $`C_{}=Ce^{1/8}`$ and $`C_𝒬=Ce^{3/5}`$ for $`\omega _{1,1}`$, while $`C_{}=C_𝒬=C`$ if one uses the estimate $`\omega _{3,1}`$. Notice that our results differ from the ones found in Ref. , where corrections behaving like $`(a_{\mathrm{ho}}/R)^4`$ (without the logarithmic factor) were obtained by employing a different procedure based on a perturbative solution for the excited state wave function.

The physical role played by the kinetic energy and the $`R`$ dependence of the finite size corrections can be understood in terms of the classical picture developed in Ref.  for the surface modes of a Bose condensed gas confined by a linear potential. In this case the resulting dispersion relation takes the form $`\omega ^2=Fq/m`$ plus $`q^4\mathrm{log}q`$ corrections originating from the kinetic energy term, where $`q`$ is the wave vector of the surface wave. In this equation $`F=m\omega _{ho}^2R`$ is the harmonic force evaluated at the Thomas-Fermi radius $`R`$. The correspondence with our result (24) for the surface quadrupole mode is achieved by imposing the discreteness condition $`q1/R`$ to the wave vector.

It is also worth noticing the different sign of the finite size corrections entering the monopole and quadrupole frequencies. Furthermore the correction for the quadrupole is about $`5`$ times larger than for the monopole.

It is finally interesting to compare the finite size corrections predicted by Eqs. (23) and (24) with the ones due to beyond mean field effects . These are fixed by the gas parameter

$$a^3n(0)=\frac{15^{2/5}}{8\pi }\left(N^{1/6}\frac{a}{a_{\mathrm{ho}}}\right)^{12/5}$$
(25)
calculated at the central density $`n(0)`$ of the gas. These effects are not included in the mean field Gross-Pitaevskii theory and take into account the changes in the equation of state originating from quantum correlation effects according to the theory of Lee-Yang-Huang . The dependence of the gas parameter (25) on the relevant physical variables ($`N`$, $`a`$ and $`a_{\mathrm{ho}}`$) is quite different with respect to the one exhibited by the ratio (22). Furthermore in Ref.  it was found that only the monopole oscillation is affected by beyond mean field effects. This should introduce enough flexibility to distinguish, in the study of the collective frequencies, between finite size and beyond mean field effects. First experimental efforts in this direction were made on 85Rb where the value of the scattering length can be tuned, thereby amplifying the effects of quantum correlations .

## IV Results for deformed traps

Until now our analysis was restricted to the case of spherical traps. It is interesting, especially from the experimental viewpoint, to extend the above results to the case of deformed external potentials. We will focus on harmonic traps having axial symmetry ($`\omega _x=\omega _y=\omega _{}`$ and $`\omega _z=\lambda \omega _{}`$, the parameter $`\lambda `$ characterizing the shape of the trap). Having proven in the previous section that the accuracy of the sum rule approach is logarithmic, we will present in the following only the corrections to the estimate $`\sqrt{m_3/m_1}`$, whose calculation is straightforward.

We are interested in particular in the study of the collective excitations described by the ($`\mathrm{}=2`$, $`m=2`$) quadrupole operator (8) and by the operator

$$=\underset{i=1}{\overset{N}{}}\left(\alpha r_i^2+\beta z_i^2\right),$$
(26)
which accounts for the coupling between the $`\mathrm{}=2`$, $`m=0`$ and the the $`\mathrm{}=0`$ excitations. The parameters $`\alpha `$ and $`\beta `$ fix the relative weight of the two modes. Moreover one can show that by treating them as variational parameters within a sum rule calculation, one recovers the correct dispersion for the two decoupled modes in Thomas-Fermi regime, first obtained in Ref.  by solving the linearized hydrodynamic equations of superfluids.

Let us first discuss the $`m=2`$ mode. By using the results (A3) and (A4) reported in Appendix A it is straightforward to find the result

$$\underset{R_{}a_{}}{lim}\omega _𝒬=\sqrt{2}\omega _{}\left[1+\frac{7}{4}\frac{E_{\mathrm{kin}}}{N\mu _{\mathrm{TF}}}\right],$$
(27)
where the Thomas-Fermi chemical potential is the natural generalization of Eq. (21)

$$\mu _{\mathrm{TF}}=\frac{1}{2}M\omega _{}^2R_{}^2,$$
(28)
with the radial size of the condensate given by

$$R_{}=(15\lambda \frac{Na}{a_{}})^{2/5}a_{}$$
(29)
and $`a_{}=\sqrt{\mathrm{}/M\omega _{}}`$.

For the $`m=0`$ mode, by using Eqs. (A5) and (A6), we find

$$\omega ^2(m=0)=\frac{1}{NM}\frac{\alpha ^2(10E_{\mathrm{ho}}2E_{\mathrm{kin}})+2\alpha \beta (10E_{\mathrm{ho}z}2E_{\mathrm{kin}z})+\beta ^2(6E_{\mathrm{ho}z}+2E_{\mathrm{kin}z})}{\alpha ^2r^2+2\alpha \beta z^2+\beta ^2z^2},$$
(30)
for the ratio $`(\mathrm{}\omega )^2=m_3/m_1`$, where $`E_{\mathrm{kin}z}=_{i=1}^Np_{z_i}^2/2M`$ and $`E_{\mathrm{ho}z}=M\omega _z^2_{i=1}^Nz_{i}^{}{}_{}{}^{2}/2`$ are the axial contributions to $`E_{\mathrm{kin}}`$ and $`E_{\mathrm{ho}}`$ respectively. By imposing that the ratio (30) be stationary with respect to the variation of $`\alpha `$ (or $`\beta `$) one can determine the decoupled eigenfrequencies. We present here only the explicit results in the simpler cases of strongly deformed traps. In the limit $`\lambda 1`$ (cigar) one finds that the lowest of the two frequencies of Eq. (30) reads

$$\underset{R_{}a_{}}{lim}\omega (m=0)=\sqrt{\frac{5}{2}}\omega _z\left[1+\frac{7}{20}\frac{E_{\mathrm{kin}}}{N\mu _{\mathrm{TF}}}\right].$$
(31)
Conversely when $`\lambda 1`$ (disk) one gets

$$\underset{R_{}a_{}}{lim}\omega (m=0)=\sqrt{\frac{10}{3}}\omega _{}\left[1+\frac{14}{15}\frac{E_{\mathrm{kin}z}}{N\mu _{\mathrm{TF}}}\right].$$
(32)
According to Eqs. (27), (31) and (32), in order to make quantitative estimates of the finite size corrections one needs an explicit expression for $`E_{\mathrm{kin}}`$ and $`E_{\mathrm{kin}z}`$. We have calculated these quantities by employing the method developed in Ref. . In Appendix D we report the results, which explicitly depend on the deformation of the trap.

## V Conclusions

In this paper we have used a sum rule approach to calculate the frequency of the monopole and quadrupole modes of dilute Bose gas confined by a harmonic trap at $`T=0`$. We have shown that with logarithmic accuracy the sum rules are exhausted by a single mode not only in Thomas-Fermi regime but also if one includes the first corrections due to the finite size of the could. We have provided an explicit expression for the frequency shift of the monopole and quadrupole modes due to finite size effects both for spherical and deformed traps.

## ACKNOWLEDGMENTS

Useful discussions with L.P. Pitaevskii are acknowledged.

## A Results for the energy weighted and the cubic sum rules

The $`m_1`$ and $`m_3`$ sum rules can be calculated by carrying out the commutators of Eqs. (5)-(6) with the mean field Hamiltonian (3). The results for $`m_1`$ are instead derived in Appendix B.

For the monopole operator (7) one finds:

$`m_1^{}`$ $`={\displaystyle \frac{2\mathrm{}^2}{M}}Nr^2`$ (A1)
$`m_3^{}`$ $`={\displaystyle \frac{2\mathrm{}^4}{M^2}}\left(4E_{\mathrm{kin}}+4E_{\mathrm{ho}}+9E_{\mathrm{int}}\right).`$ (A2)
For the quadrupole operator (8) the results read:

$`m_1^𝒬`$ $`={\displaystyle \frac{\mathrm{}^2}{2M}}Nx^2+y^2`$ (A3)
$`m_3^𝒬`$ $`={\displaystyle \frac{2\mathrm{}^4}{M^2}}\left(E_{\mathrm{kin}}+E_{\mathrm{ho}}\right).`$ (A4)
For the operator (26) we find:

$`m_1^{}`$ $`={\displaystyle \frac{2\mathrm{}^2}{M}}N\left(\alpha ^2r^2+2\alpha \beta z^2+\beta ^2z^2\right)`$ (A5)
$`m_3^{}`$ $`={\displaystyle \frac{2\mathrm{}^2}{M^2}}\left(\alpha ^2(10E_{\mathrm{ho}}2E_{\mathrm{kin}})+2\alpha \beta (10E_{\mathrm{ho}z}2E_{\mathrm{kin}z})+\beta ^2(6E_{\mathrm{ho}z}+2E_{\mathrm{kin}z})\right).`$ (A6)
## B The inverse energy weighted sum rule

The aim of this appendix is to determine the inverse energy weighted sum rule relative to the monopole and quadrupole operators (7) and (8) when the Hamiltonian of the system is given by Eq. (3), with $`\omega _x=\omega _y=\omega _z=\omega _{\mathrm{ho}}`$. The procedure is based on the determination of the static response function through Eq. (14).

Monopole: in this case the problem is equivalent to finding the static response when the system is perturbed by an external field proportional to the monopole operator (7), which can be described by the Hamiltonian:

$$H_{\mathrm{pert}}=\gamma \underset{i=1}{\overset{N}{}}r_i^2,$$
(B1)
where $`\gamma `$ is the strength of the perturbation. Notice that adding such a perturbation to the Hamiltonian (3) has simply the effect of renormalizing the frequency of the external confinement, which now becomes

$$\omega _\gamma =\sqrt{\omega _{\mathrm{ho}}^2+\frac{2\gamma }{M}}.$$
(B2)
As a consequence, the average square radius will now depend on $`\gamma `$.

By definition, the static polarizability is given by

$$\chi _{}(\omega =0)=N\left(\frac{r^2_\gamma }{\gamma }\right)_{\gamma =0}.$$
(B3)
In order to evaluate (B3), one should exploit the $`\gamma `$ dependence of $`r^2_\gamma `$ in an explicit way. To this purpose it is convenient to write the stationary Gross-Pitaevskii equation in a suitable form, by rescaling the lengths, the energies and the order parameter as suggested in Ref. :

$`𝐫`$ $`=a_\gamma 𝐫_1`$ (B4)
$`E`$ $`=\mathrm{}\omega _\gamma E_1`$ (B5)
$`\psi (𝐫)`$ $`=\sqrt{{\displaystyle \frac{N}{a_\gamma ^3}}}\psi _1(𝐫_1),`$ (B6)
with $`\psi _1(𝐫_1)`$ normalized to $`1`$. According to such notation the time-independent GP equation for $`\psi _1(𝐫_1)`$ reads

$$[_1^2+r_1^2+u_1|\psi (𝐫_1)|^2]\psi _1(𝐫_1)=2\mu _1\psi _1(𝐫_1),$$
(B7)
where $`u_1=8\pi Na/a_\gamma `$, $`a_\gamma =\sqrt{\mathrm{}/M\omega _\gamma }`$ and $`\mu =\mathrm{}\omega _\gamma \mu _1`$. Consequently one can write the square radius of the system as

$$r^2_\gamma =a_\gamma ^2r_1^2_\gamma =a_\gamma ^2𝑑𝐫_1r_1^2|\psi _1(𝐫_1)|^2.$$
(B8)
Note that $`r_1^2_\gamma `$ exhibits a dependence on $`\gamma `$, since the Gross-Pitaevskii equation, in the form (B7), depends on $`\gamma `$ through the mean field term proportional to $`u_1`$. The derivative of $`r^2_\gamma `$ with respect to $`\gamma `$ can then be expressed in terms of the derivative with respect to $`u_1`$ or $`N`$.

By using Eqs. (B8) and (B3) one finally finds

$`\chi _{}(\omega =0)`$ $`=N\left[r_1^2_\gamma \left({\displaystyle \frac{a_\gamma ^2}{\gamma }}\right)+a_\gamma ^2\left({\displaystyle \frac{r_1^2_\gamma }{\gamma }}\right)\right]_{\gamma =0}`$ (B9)
$`=N\left[{\displaystyle \frac{a_\gamma ^2}{M\omega _{\mathrm{ho}}^2}}\left(r_1^2_\gamma {\displaystyle \frac{N}{2}}{\displaystyle \frac{r_1^2_\gamma }{N}}\right)\right]_{\gamma =0},`$ (B10)
which, combined with Eq. (14), gives the relevant result

$$m_1^{}=\frac{N}{2M\omega _{\mathrm{ho}}^2}\left(r^2\frac{N}{2}\frac{r^2}{N}\right).$$
(B11)
Quadrupole: in this case the idea is that in the presence of a harmonic potential deformed in the $`x`$-$`y`$ plane ($`\omega _x\omega _y`$) the quadrupole operator (8) can be written also in the form

$$𝒬=\frac{1}{iM\mathrm{}}\frac{[H,L_z]}{(\omega _x^2\omega _y^2)},$$
(B12)
where

$$L_z=\underset{i=1}{\overset{N}{}}(x_ip_{y_i}y_ip_{x_i})$$
(B13)
is the third component of the angular momentum operator. By using the relationship (14) combined with the expression (B12) for the quadrupole operator, one easily shows that $`m_1^𝒬`$ can be written as

$$m_1^𝒬=\frac{1}{2}\frac{[L_z,𝒬]}{i\mathrm{}M(\omega _x^2\omega _y^2)}.$$
(B14)
From Eq. (B14), taking the spherical limit, one immediately finds the result

$$m_1^𝒬=\frac{N}{4M\omega _{\mathrm{ho}}^2}\underset{\epsilon 0}{lim}x^2+y^2\left(\frac{\delta }{\epsilon }\right),$$
(B15)
where the deformation of the could $`\delta `$ and the one of the trap $`\epsilon `$ have been introduced in Sec. III.

## C Comparison between $`𝝎_{\mathrm{𝟑}\mathbf{,}\mathrm{𝟏}}`$ and $`𝝎_{\mathrm{𝟏}\mathbf{,}\mathbf{}\mathrm{𝟏}}`$

We report here the full derivation of the finite size corrections relative to both monopole and the quadrupole modes in the isotropic case.

Monopole: It is immediate to derive the corrections to $`\omega _{3,1}`$ by using the Eq. (9) combined with the Eq. (22). One finds the result

$$\omega _{_{3,1}}=\sqrt{5}\omega _{\mathrm{ho}}\left[1\frac{7}{6}\left(\frac{a_{\mathrm{ho}}}{R}\right)^4\mathrm{log}\left(\frac{R}{Ca_{\mathrm{ho}}}\right)\right].$$
(C1)
It is quite simple also to calculate the finite size corrections to $`\omega _{1,1}`$ starting from Eq. (15). One has to use the identity

$$r^2=\frac{2}{M\omega _{\mathrm{ho}}^2}\frac{E_{\mathrm{ho}}}{N},$$
(C2)
with the harmonic oscillator energy given by Eq. (19). Note that the $`N`$-dependence of the mean square radius $`r^2`$ originates from both the Thomas-Fermi chemical potential (21) and from the correction behaving like $`1/R^2\mathrm{log}(R/a_{\mathrm{ho}})`$, where $`R`$ is given in Eq. (18). By carrying out carefully the calculation and taking into account all the terms of order $`1/R^4`$ , one finally finds

$$\omega _{_{1,1}}=\sqrt{5}\omega _{\mathrm{ho}}\left[1\frac{7}{6}\left(\frac{a_{\mathrm{ho}}}{R}\right)^4\mathrm{log}\left(\frac{Re^{1/8}}{Ca_{\mathrm{ho}}}\right)\right],$$
(C3)
which coincides with the result (C1) for $`\omega _{3,1}`$ within logarithmic accuracy. Note that the inequality (13) is correctly fulfilled by Eqs. (C1) and (C3), the relative difference between the two estimates becoming smaller than $`10^4`$ for $`R>6a_{\mathrm{ho}}`$.

Quadrupole: the calculation of the corrections starting from the ratio $`\sqrt{m_3/m_1}`$ is immediate if one combines Eqs. (10) and (22). This yields the result

$$\omega _{𝒬_{3,1}}=\sqrt{2}\omega _{\mathrm{ho}}\left[1+\frac{35}{6}\left(\frac{a_{\mathrm{ho}}}{R}\right)^4\mathrm{log}\left(\frac{R}{Ca_{\mathrm{ho}}}\right)\right].$$
(C4)
The calculation of the corrections to $`\omega _{1,1}`$ is less trivial since it requires the explicit knowledge of the dependence of the condensate deformation $`\delta `$ on $`\epsilon `$. To this purpose we determine the finite size corrections to the square radii $`x^2=2E_{\mathrm{ho}x}/NM\omega _x^2`$ and $`y^2=2E_{\mathrm{ho}y}/NM\omega _y^2`$ by using the virial relationship (11) written for the $`x`$ direction

$$E_{\mathrm{kin}x}E_{\mathrm{ho}x}+\frac{1}{2}E_{\mathrm{int}}=0,$$
(C5)
the analogous equations holding for $`y`$ and $`z`$. One then finds

$`y^2x^2`$ $`={\displaystyle \frac{1}{NM}}{\displaystyle \frac{\omega _x^2+\omega _y^2}{\omega _x^2\omega _y^2}}\left[E_{\mathrm{int}}(\epsilon )\epsilon +(1+\epsilon )E_{\mathrm{kin}y}(\epsilon )(1\epsilon )E_{\mathrm{kin}x}(\epsilon )\right]`$ (C6)
$`y^2+x^2`$ $`={\displaystyle \frac{1}{NM}}{\displaystyle \frac{\omega _x^2+\omega _y^2}{\omega _x^2\omega _y^2}}\left[E_{\mathrm{int}}(\epsilon )+(1+\epsilon )E_{\mathrm{kin}y}(\epsilon )+(1\epsilon )E_{\mathrm{kin}x}(\epsilon )\right],`$ (C7)
where we have explicitly indicated the dependence of $`E_{\mathrm{int}}`$, $`E_{\mathrm{kin}x}`$ and $`E_{\mathrm{kin}y}`$ on the trap deformation. By expanding the above equations with respect to $`\epsilon `$ one can easily determine the limit entering Eq. (16) getting, to the lowest order,

$$\omega _{𝒬_{1,1}}=\sqrt{2}\omega _{\mathrm{ho}}\underset{\epsilon 0}{lim}\left(\frac{\delta }{\epsilon }\right)^{1/2}=\sqrt{2}\omega _{\mathrm{ho}}\left[1+\frac{7}{4}\frac{E_{\mathrm{kin}x}^{}(0)E_{\mathrm{kin}y}^{}(0)}{N\mu _{\mathrm{TF}}}\right],$$
(C8)
with $`E_{\mathrm{kin}x}^{}(0)=(E_{\mathrm{kin}x}/\epsilon )_{\epsilon =0}`$ and analogously for $`y`$. In order to estimate the difference $`E_{\mathrm{kin}x}^{}(0)E_{\mathrm{kin}y}^{}(0)`$ entering the last equation, we can use the results reported in the next Appendix D for the different contributions to the kinetic energy in Thomas-Fermi regime. By setting $`\lambda =\omega _y/\omega _x`$ and using the relationship $`\epsilon (1\lambda )`$ holding for small $`\epsilon `$, one finds that

$$E_{\mathrm{kin}x}^{}(0)E_{\mathrm{kin}y}^{}(0)=\frac{2}{3}E_{\mathrm{kin}}N\frac{\mathrm{}^2}{MR^2},$$
(C9)
which combined with Eq. (C8) finally gives

$$\omega _{𝒬_{1,1}}=\sqrt{2}\omega _{\mathrm{ho}}\left[1+\frac{35}{6}\left(\frac{a_{\mathrm{ho}}}{R}\right)^4\mathrm{log}\left(\frac{R}{Ca_{\mathrm{ho}}e^{3/5}}\right)\right].$$
(C10)
It is worth pointing out that, with logarithmic accuracy, Eq. (C10) coincides with the result (C4) derived for $`\omega _{3,1}`$.

## D Kinetic energy in Thomas-Fermi regime

By explicitly calculating the various contributions $`E_{\mathrm{kin}}`$ and $`E_{\mathrm{kin}z}`$ of the total kinetic energy using the method developed in Ref.  we find the following results:

Disk shaped trap, $`\lambda >1`$

$`{\displaystyle \frac{E_{\mathrm{kin},}}{N}}`$ $`={\displaystyle \frac{\mathrm{}^2}{2MR_{}^2}}\left[{\displaystyle \frac{5}{2}}\mathrm{log}\left({\displaystyle \frac{2}{\lambda ^{2/3}}}{\displaystyle \frac{R_{}}{\delta _0}}\right){\displaystyle \frac{5}{2}}\left(1+{\displaystyle \frac{1}{3}}{\displaystyle \frac{1}{\lambda ^21}}\right){\displaystyle \frac{\mathrm{arctan}\left(\sqrt{\lambda ^21}\right)}{\sqrt{\lambda ^21}}}+{\displaystyle \frac{5}{18}}{\displaystyle \frac{1+2\lambda ^2}{\lambda ^21}}\right]`$ (D1)
$`{\displaystyle \frac{E_{\mathrm{kin},\mathrm{z}}}{N}}`$ $`={\displaystyle \frac{\mathrm{}^2}{2MR_{}^2}}\left[{\displaystyle \frac{5}{4}}\lambda ^2\mathrm{log}\left({\displaystyle \frac{2}{\lambda ^{2/3}}}{\displaystyle \frac{R_{}}{\delta _0}}\right)+{\displaystyle \frac{5}{6}}\left(1+{\displaystyle \frac{1}{\lambda ^21}}\right){\displaystyle \frac{\mathrm{arctan}\left(\sqrt{\lambda ^21}\right)}{\sqrt{\lambda ^21}}}{\displaystyle \frac{5}{18}}{\displaystyle \frac{1+2\lambda ^2}{\lambda ^21}}{\displaystyle \frac{5}{18}}{\displaystyle \frac{5}{9}}\lambda ^2\right]`$ (D2)
Cigar shaped trap, $`\lambda <1`$

$`{\displaystyle \frac{E_{\mathrm{kin},}}{N}}`$ $`={\displaystyle \frac{\mathrm{}^2}{2MR_{}^2}}\left[{\displaystyle \frac{5}{2}}\mathrm{log}\left({\displaystyle \frac{2}{\lambda ^{2/3}}}{\displaystyle \frac{R_{}}{\delta _0}}\right){\displaystyle \frac{5}{4}}\left(1{\displaystyle \frac{1}{3}}{\displaystyle \frac{1}{1\lambda ^2}}\right){\displaystyle \frac{1}{\sqrt{1\lambda ^2}}}\mathrm{log}\left({\displaystyle \frac{1+\sqrt{1\lambda ^2}}{1\sqrt{1\lambda ^2}}}\right){\displaystyle \frac{5}{18}}{\displaystyle \frac{1+2\lambda ^2}{1\lambda ^2}}\right]`$ (D3)
$`{\displaystyle \frac{E_{\mathrm{kin},\mathrm{z}}}{N}}`$ $`={\displaystyle \frac{\mathrm{}^2}{2MR_{}^2}}\left[{\displaystyle \frac{5}{4}}\lambda ^2\mathrm{log}\left({\displaystyle \frac{2}{\lambda ^{2/3}}}{\displaystyle \frac{R_{}}{\delta _0}}\right)+{\displaystyle \frac{5}{12}}\left(1{\displaystyle \frac{1}{1\lambda ^2}}\right){\displaystyle \frac{1}{\sqrt{1\lambda ^2}}}\mathrm{log}\left({\displaystyle \frac{1+\sqrt{1\lambda ^2}}{1\sqrt{1\lambda ^2}}}\right)+{\displaystyle \frac{5}{18}}{\displaystyle \frac{1+2\lambda ^2}{1\lambda ^2}}{\displaystyle \frac{5}{18}}{\displaystyle \frac{5}{9}}\lambda ^2\right].`$ (D4)
In Eqs. (D1)-(D4) $`\delta _0`$ represents the surface thickness at $`z=0`$. Its dependence on the relevant parameters of the system is

$$\delta _0=\frac{1.776}{2^{1/3}}(15\lambda \frac{Na}{a_{}})^{1/15}a_{}.$$
(D5)
The above results are consistent with the results of Ref.  for the total kinetic energy $`E_{\mathrm{kin}}=E_{\mathrm{kin}}+E_{\mathrm{kin}z}`$, which reduces to (17) in the isotropic case $`\lambda =1`$.


text:
# 1 Introduction

## 1 Introduction

We present combined results from the ALEPH, DELPHI, L3 and OPAL Collaborations on searches for the neutral Higgs bosons decaying into “invisible” particles such as neutralinos  or majorons . The results are obtained by combining the data collected in the year 2000 at centre-of-mass energies between 202 and 209 GeV with earlier data collected at lower energies . The new data represent an integrated luminosity of approximately 868 $`\mathrm{pb}^1`$ in total.

Unless explicitly specified, all cross-sections, branching ratios and many other physics quantities which are used in this combination of data, are calculated within HZHA .

Each experiment has generated Monte Carlo event samples for the Higgs boson signal and the various background processes, typically, at 202, 204, 206, 208 and 210 GeV energies. Cross-sections, branching ratios, distributions of the reconstructed mass and other discriminating variables relevant to the combination have been interpolated to energies which correspond to the data sets. In this procedure special care has been taken to the regions of kinematic cutoff where the signal and background distributions vary rapidly. It has been established that the interpolation procedure does not add significantly to the final systematic uncertainties.

The statistical procedure adopted for the combination of the data and the precise definition of the confidence levels $`CL_b,CL_{s+b},CL_s`$ by which the search results are expressed, follow our usual definitions.

## 2 Combined search for ‘invisible’ Higgs boson decays

At LEP the SM Higgs boson is expected to be produced mainly via the Higgs-strahlung process $`\mathrm{e}^+\mathrm{e}^{}`$$``$ HZ, while contributions from the WW$``$ H fusion channel, $`\mathrm{e}^+\mathrm{e}^{}`$$``$ H$`\nu _e\overline{\nu _e}`$, are typically below 10%. The Higgs boson may not make is presence obvious. For instance, in supersymmetric theories, depending upon the parameters, the decay of a Higgs boson into neutralinos might dominate. The Higgs boson could then be invisible at LEP. Majoron models can also produce dominantly invisible decay modes. However, if the Higgs boson is produced through the Higgs-strahlung process, the Z can be detected, and the presence of the Higgs boson inferred. This production processes is assumed here.

The four LEP collaborations performed searches for acoplanar jets (H$``$invisible)(Z$``$$`\mathrm{q}\overline{\mathrm{q}}`$) or leptons (H$``$invisible)(Z$``$ll).

The analysis procedures of the four LEP experiments producing the inputs for the present combination are described in individual documents ; we merely summarise the results in Table 1. For DELPHI, the channel (Z$``$ll) includes the $`\tau \tau `$ decay mode of the Z boson, which is not used by the other experiments.

The large spread in the numbers of selected candidates reflects substantial differences in the selection methods and optimisation procedures. For instance, ALEPH and L3 use a sliding analysis technique. Only the relevant candidates for a Higgs mass hypothesis are reported. For table 1 the Higgs boson mass was taken as 110 $`\mathrm{GeV}/c^2`$.

The test-statistic as a function of the mass $`m_\mathrm{H}`$, computed for the observed results, is shown in Figure 1. In the presence of a signal it should have a minimum near the true Higgs boson mass. A negative value would indicate a preference for the signal hypothesis and the more negative the value the more significant the signal. The full-line curve representing the observation is in agreement with the dashed line representing the background hypothesis, and deviates from the dotted curves which represent the most likely signal + background situation. In fact, in the region 112 to 115 $`\mathrm{GeV}/c^2`$there is a deficit of events correspondig to around 1.5 standard deviations from the background expectation, which means that the observed limit will be somewhat stronger than expected.

The compatibility with background of the result is given by $`1CL_b`$, which is plotted as a function of $`m_\mathrm{H}`$ in Figure 2. The dotted line shows the expectation in the presence of a signal; its crossing with the $`5\sigma `$ line at 109.5 $`\mathrm{GeV}/c^2`$ indicates the range of sensitivity of the data to a discovery. As expected, there is no suggestion of any signal.

A 95% confidence level lower limit on the Higgs boson mass may be set by identifying the mass region where $`CL_s<0.05`$, as shown in Figure 3. The $`CL_s`$ curve shown gives the limit on $`m_H`$ assuming a 100% branching ratio into invisible modes. The median limit expected in the absence of a signal is 113.5 $`\mathrm{GeV}/c^2`$ and the limit observed by combining the LEP data is 114.4 $`\mathrm{GeV}/c^2`$. The inclusion of systematic errors, which have been neglected, is expected to reduce this lower bound by approximately 100 $`\mathrm{MeV}/c^2`$.

The upper limit on the rate of $`Hinvisible`$ s a function of $`m_H`$ is shown in Figure 4. The scale is:

$$\xi ^2=\frac{\sigma _{HZ}}{\sigma _{HZ}^{StandardModel}}Br(Hinvisible)$$
(1)
as a fraction of the rate expected from a Standard Model H decaying 100% invisibly. The same information is shown in figure 5 but as a limit on the cross-section at 206 GeV. There is a somewhat stronger exclusion than expected around the $`Z`$ region, due to a deficit of candidates.

## 3 Summary

The LEP working group for Higgs boson searches has updated its previous combined limit for the mass of a Higgs boson decaying invisibly, including the data collected in the year 2000 at energies between 200 and 209 GeV, for a total integrated luminosity of approximately 868 $`\mathrm{pb}^1`$. In the absence of a statistically significant excess in the data, a new lower bound of 114.4 $`\mathrm{GeV}/c^2`$ has been obtained at the 95% confidence level and assuming the Standard Model production cross-section and that the Higgs boson exclusively decays invisibly.

ALL THE RESULTS QUOTED IN THIS NOTE ARE PRELIMINARY.