text:
# The Kinematics and Metallicity of the M31 Globular Cluster System

## 1 Introduction

The mechanisms involved in galaxy formation remain one of the major unsolved problems in astronomy. Over the past several decades, various groups have proposed a range of models which endeavor to explain the observed properties of galaxies and their globular cluster systems (GCSs). One early formation model from Eggen, Lynden-Bell, & Sandage (1962) argues for a single, large-scale collapse of material to form galactic bodies such as the Milky Way. A principal competing model maintains that formation has occurred through random mergers of fragmented gas clouds over the course of billions of years, implying a hierarchical origin (Searle & Zinn, 1978; Cole et al., 2000). In order to interpret the bimodality seen in the globular cluster metallicity and color distributions of many galaxies (Ashman & Zepf, 1998), additional scenarios have been put forth such as multiphase in situ formation (Forbes, Brodie, & Grillmair, 1997), major mergers (Zepf & Ashman, 1993) and tidal stripping/capture (Côté, Marzke, & West, 1998).

Support for the idea that the Milky Way halo formed in a hierarchical or episodic manner arises from the lack of a significant abundance gradient throughout the entire halo globular cluster population (Parmentier et al., 2000), a notable decrease in the radial component of the velocity ellipsoid of the stellar halo beyond the solar circle (Sommer-Larsen et al., 1997), and the detection of a possible spread in cluster and stellar ages within the outer halo (Sarajedini, Chaboyer, & Demarque, 1997; Rosenberg et al., 1999; Stetson et al., 1999). In addition, some studies have uncovered evidence of kinematic and chemical sub-clustering within the globular cluster population of the Andromeda galaxy which may point to accretion remnants (Ashman & Bird, 1993; Saito & Iye, 2000). However, it is possible that some degree of fragmentation might arise in an infall scenario due to thermal and gravitational instabilities within the collapsing proto-galactic gas cloud. A better understanding of the fragment properties including chronology, chemistry and, in particular, kinematics would help us to identify the origin of sub-clustering and to differentiate between the various formation models.

It has been noted often that globular clusters (GCs) provide us with ideal probes of galaxy structure and formation mechanisms. As the oldest galactic stellar systems, globular clusters hold the key to uncovering the formation record of their host environments. Correlations between chemical content (which points to chronology) and the kinematics of globular clusters allow us to investigate the early formative stages of a galaxy and, to some extent, probe its ongoing evolution.

The M31 globular cluster system provides an obvious and desirable target of observation: (i) its proximity makes it the most accessible large system outside of the Milky Way; (ii) with $`>435`$ confirmed candidate members (Barmby et al., 2000), the M31 GCS is sufficiently populous that it provides a statistically significant sample size; (iii) unlike the Milky Way GCS, the M31 population can be observed over the full extent of the galaxy with less severe effects of line-of-sight contamination; (iv) all of the M31 clusters lie at essentially the same distance, and thus are not subject to such large individual distance uncertainties as are the Milky Way globulars. However, the proximity of M31 introduces certain issues of contamination within the cluster sample, such as a difficulty in distinguishing between GC candidates and star-forming regions in many photometric surveys. Furthermore, faint globular clusters can easily be missed against the bright galaxy background.

Previous studies of the Milky Way and M31 globular cluster systems have revealed that these two populations exhibit some remarkable similarities. Zinn (1985) found that the Galactic GCS demonstrated clear evidence for bimodality in its metallicity distribution, and proposed that it could be considered as two distinct subsystems: a flattened, metal-rich, rapidly rotating, kinematically cold disk system, and a spherical, metal-poor halo with lower net rotation and higher velocity dispersion. Signs of bimodality have been found in the M31 metallicity distribution (Ashman & Bird, 1993) and in metallicity-sensitive color distributions (Barmby et al., 2000). Separating the M31 clusters into metallicity sub-systems along the same lines as the Milky Way GCS, Huchra, Brodie, & Kent (1991) found that the metal-rich clusters within $`R7\mathrm{}`$ form a rapidly rotating ($`100200`$ km s-1) disk-like system, whereas the metal-poor clusters within this radius exhibit no net rotation. Beyond this radius, they were unable to distinguish the metal-rich from metal-poor populations on the basis of kinematics. Later studies have provided additional evidence supporting the view that M31 has a rapidly rotating system of metal-rich globular clusters (Ashman & Bird, 1993; Barmby et al., 2000).

Armandroff (1989) confirmed the presence of a modest metallicity gradient within the system of disk clusters in the Milky Way. In the first extensive spectroscopic survey of the M31 GCS, van den Bergh (1969) found no significant metallicity gradient with position but did report a wide range of abundances at any given radius. More recently, Huchra et al. (1982, 1991) and Sharov (1988) have shown that there is evidence for a weak but measurable metallicity gradient with projected radius in the M31 cluster system. A similar finding was reported by Barmby et al. (2000) based on their large sample of metal-sensitive colors. Huchra, Brodie, & Kent (1991) found a mean M31 GC metallicity of $`[\mathrm{Fe}/\mathrm{H}]=1.2`$, slightly higher than that for the Milky Way.

There remains a lack of persuasive evidence that would allow us to distinguish between halo formation models — or indeed to suggest a clear hybrid of (or alternative to) the classical formation scenarios. It is clear that an extensive database of high-precision data for the kinematic and chemical tracers of a large number of host galaxies is required in order to help solve this currently intractable problem of galaxy formation. Only then can the samples be examined for clear evidence of kinematic anomalies, metallicity gradients and sub-structure which, if present, may point to the formation and enrichment processes which have occurred in the galaxy. With this objective, we present spectroscopic results for $`>200`$ M31 cluster candidates, adding a large number of new observations to the database.

This paper will take the following form: Section 2 contains a description of the target selection and coordinate determinations. In Section 3, we provide a summary of the data acquisition and reduction procedure along with a description of the velocity and metallicity determinations. We incorporate these new results into the broader sample of available M31 cluster data to provide an analysis of the abundance distribution of the globular cluster system of our nearest large spiral galaxy neighbor in Section 4. In Section 5, we examine the kinematics of the cluster system and provide a mass estimate for the galaxy. A summary and discussion of the results of this study are presented in Section 6.

## 2 Target Selection and Coordinates

The list of targets was selected from the catalogue of Battistini et al. (1987), with some additional central targets from unpublished CCD photometry by J. Huchra. Cluster positions accurate to $`0\stackrel{}{\mathrm{.}}2`$ were determined from Automatic Plate Measuring (APM) scans of deep photographic plates of a $`6^{}\times 6^{}`$ area around M31. The accuracy of the APM cluster coordinates represent a marked improvement over existing M31 GCS positions, which have typically been measured to the nearest arcsecond.

Galactocentric coordinates for the GCs were computed relative to an adopted M31 central position of $`\alpha _0=00^h40^m00.1^s`$, $`\delta _0=+40^{}59\mathrm{}43\mathrm{}`$ (B1950):

$`X`$ $`=`$ $`C_1\mathrm{sin}(PA)+C_2\mathrm{cos}(PA)`$ (1)
$`Y`$ $`=`$ $`C_1\mathrm{cos}(PA)+C_2\mathrm{sin}(PA),`$
where $`C_1=[\mathrm{sin}(\alpha \alpha _0)\mathrm{cos}\delta ]`$ and $`C_2=[\mathrm{sin}\delta \mathrm{cos}\delta _0\mathrm{cos}(\alpha \alpha _0)\mathrm{cos}\delta \mathrm{sin}\delta _0]`$. The $`X`$ coordinate represents the position along the major axis of the galaxy increasing towards the north-east along the major axis. The $`Y`$ coordinate is the distance from the major axis above and below the galactic plane, with positive values towards the north-west. The position angle (PA) of the major axis was taken to be $`38^{}`$ (Kent, 1989).

Numerous groups have contributed to the identification of the globular cluster members of M31 (Hubble, 1932; Vetes̆nik, 1962; Sargent et al., 1977; Crampton et al., 1985; Battistini et al., 1980, 1987, 1993). The ultimate result of this combined effort is that most of the clusters are known by more than one name, a fact which has the potential to cause confusion. The letters preceding the identification number in the M31 cluster target names given hereafter indicate the relevant catalogue reference as provided in Table 1.

## 3 The WYFFOS Data

The data for this study were obtained at the William Herschel 4.2 m telescope111The William Herschel Telescope is operated on the island of La Palma by the Isaac Newton Group in the Spanish Observatorio del Roque de los Muchachos of the Instituto de Astrofisica de Canarias. (WHT) in La Palma, Canary Islands, during the nights of November 3 to 6, 1996. The WYFFOS Wide Field Fibre Optic Spectrograph (Bridges, 1998) was used with two gratings to provide a total spectral coverage of $`37005600`$ Å. The H2400B 2400-line grating yielded a dispersion of 0.8 Å/pixel and a spectral resolution of 2.5 Å over the range $`37004500`$ Å. This wavelength range was selected to investigate the CN feature at 3883 Å, the H & K lines of calcium, H$`\delta `$, the CH G-band and the 4000 Å continuum break. The R1200R 1200-line grating observations had a dispersion of 1.5 Å/pixel and a resolution of 5.1 Å over the spectral range $`44005600`$ Å to add absorption features such as H$`\beta `$, the Mgb triplet and two iron lines near 5300 Å.

Six different fibre configurations were selected to target a total of 288 globular cluster candidates over the full extent of the major axis of M31. An overlay of the WYFFOS target fields is shown superimposed on a Digitized Sky Survey image of the galaxy in Figure 1. A record of the observations is provided in Table 2, which includes the coordinates of the field centers, observation dates and exposure times. Regrettably, minor axis fields were not observed during this run due in part to a loss of time from poor weather.

Spectra of twilight flats, argon arc lamps, and velocity template stars were also obtained during the same nights of observation as the targets. The data reduction — including fibre throughput corrections, spectral extraction, sky subtraction and wavelength calibration — was accomplished using the wyf\_red WYFFOS Multi-fiber IRAF222Image Reduction and Analysis Facility, distributed by the National Optical Astronomical Observatories, which is operated by AURA under contract with the NSF. package written by Jim Lewis. Ensemble median RMS values for the wavelength calibration were typically $`0.020.08`$ Å for the H2400B grating images and $`0.080.13`$ Å for the R1200R grating.

Most of the globular clusters in the M31 fields are superimposed on a relatively bright, non-uniform background of galaxy light. For the central WYFFOS fibre configurations, roughly 20 blank-sky fibres were arranged in a regular pattern across the $`1^{}`$ field to characterize the background intensity. The outer fields contained more dedicated sky fibres thanks to a lower density of targets in the area. Spectra obtained through the dedicated background sky fibres were averaged to form a composite sky spectrum which was then subtracted from each target. Attempting to match one or a few nearby sky fibres to each target spectrum (i.e. to associate the local sky level with each cluster) did not prove advantageous, since the resulting losses in signal-to-noise largely outweighed any gains in compensating for background variability even within the central fields. Examples of reduced spectra obtained using both gratings are shown in Figures 2 and 3 with the absorption features of interest labelled.

### 3.1 Radial Velocities

Radial velocities for the M31 globular clusters were calculated by performing standard Fourier cross-correlations of each target spectrum with the high signal-to-noise spectra of template stars of known velocity (Tonry & Davis, 1979). Two spectra were acquired for each target, obtained using each of the R1200R and H2400B gratings. Two template stars, HD12029 (K2 giant) and HD23169 (G2 dwarf), were also observed during the run to yield four cross-correlation combinations for the derivation of velocities for each cluster target. We required that each target yield a minimum of two satisfactory cross-correlations with peak values $`CC0.3`$ to qualify as a successful result (a result of $`CC=1`$ would indicate a perfect correlation). A mean of the cross-correlation results, weighted by the associated Tonry & Davis R-value and corrected for solar motion, was used to calculate heliocentric radial velocities for a total of 202 M31 clusters candidates in the sample. The remainder represent:

1. Confirmed non-clusters (foreground stars, background galaxies, HII regions, open clusters etc.) identified by spectral features grossly atypical of globular clusters or radial velocities that were significantly beyond reasonable limits for the M31 GC population ($`|v_{\mathrm{obs}}v_{\mathrm{sys}}|>3\sigma _v`$ where $`\sigma _v150`$ km s-1). The status of these objects is noted in the final column of Table 3.
2. Target spectra with insufficient signal-to-noise to generate reliable velocities.
3. Seven objects with inaccurate coordinates due to misidentifications or saturation on the APM scans.

The radial velocities obtained from the WYFFOS spectra are provided in Table 3. Over half (109) of the M31 GCs for which velocities were obtained had no previously published spectroscopic data. The majority of the remaining clusters (72) had published velocities with large uncertainties of $`2080`$ km s-1, and the remainder (21) had modest to high-precision velocities with errors $`<20`$ km s-1. The positions of the WYFFOS targets in each of these categories are shown schematically in Figure 4.

Robust fits revealed no significant differences between the velocities determined using either template object with a given grating. A comparison of the velocities calculated for a given target from the two gratings showed only slightly larger differences; the derived velocities found using the R1200R grating data, for example, were at worst $`6`$ km s-1 larger than those determined from the spectra obtained using the H2400B grating. An uncertainty of $`\pm 12`$ km s-1 is adopted here for the M31 clusters in our sample based on the RMS of the velocity deviations for 28 targets observed in multiple fields. On average, this represents an overestimate of the true uncertainty since it incorporates certain internal systematic errors as well as statistical errors.

Figure 5 shows a comparison of those radial velocities in the present study that overlap with other samples from the literature. A linear fit, weighted by the inverse square of the published velocity uncertainties, yields $`v_{\mathrm{publ}}=0.43+1.02v_{\mathrm{obs}}`$ km s-1. Overall, there is good agreement between the WYFFOS velocities derived herein and those of other studies within the adopted errors. Encouragingly, the RMS residuals of independent fits to the separate velocity samples were particularly low for the clusters which overlap with the high-precision echelle studies of Peterson (1989: $`v_{\mathrm{P89}}=5.8+0.96v_{\mathrm{obs}}`$ km s-1 with $`\mathrm{RMS}=5.7`$ km s-1 for 8 objects) and Dubath & Grillmair (1997: $`v_{\mathrm{DG}}=0.58+1.02v_{\mathrm{obs}}`$ km s-1 with $`\mathrm{RMS}=8.2`$ km s-1 for 5 objects). Note that the only high-precision result that deviates significantly from our fit is for B29, a cluster for which Peterson (1989) notes her velocity may be incorrect. This GC was not used in the calculation of the weighted fit or RMS, and it is circled in Figure 5.

There are two other points in Figure 5 which deviate significantly from the fit, and thus deserve mention. The open diamond above the best fit line represents the M31 globular cluster B301-S22: Barmby et al. (2000) obtained $`v_r=30\pm 20`$ km s-1, in poor agreement with the value of $`374\pm 12`$ km s-1 measured in the current study. The Barmby et al. coordinates match well with the APM coordinates for this object as provided in Table 3, so the velocity deviation does not appear to be due to a target misidentification in either study. This cluster was also measured by Federici et al. (1993), who obtained a velocity of $`419\pm 30`$ km s-1, in good agreement with the WYFFOS result. The Barmby et al. sample also had the largest RMS residuals in the comparison fits: $`\mathrm{RMS}=110`$ km s-1 for 16 objects, not including B301-S22. The reasons for this large discrepancy are not known. The open upright triangle which deviates from the fit in Figure 5 represents B109-S170, for which Huchra, Brodie, & Kent (1991) measured a velocity of $`613\pm 24`$ km s-1 (compared to $`372\pm 12`$ km s-1 from the WYFFOS results). Again, the coordinates for this object match well between both studies and the deviation does not appear to be due to a target misidentification.

### 3.2 Absorption and Reddening Corrections

Each spectrum was corrected for the effects of atmospheric absorption and for reddening in both the Milky Way and in M31 itself. Color excess values upon which the reddening calculations are based were determined using the slope ($`S`$) of the continuum between $`40005000`$ Å as tabulated in Crampton et al. (1985), where available:

$$E_{BV}=0.066S+1.17(BV)0.32.$$
(2)
The slope parameter has been shown to be a good estimator of reddening, and is not very sensitive to cluster metallicity (Elson & Walterbos, 1988).

Targets with published $`(BV)`$ colors but with no slope parameters or intrinsic color estimates were de-reddened to a typical color of $`(BV)_0=0.76`$ for M31 globular clusters at large galactocentric radii (Crampton et al., 1985). Clusters with no $`S`$ values (or colors) and those for which the calculated excess values were lower than a minimum foreground color excess of $`E_{BV}=0.10`$ (due to errors in the slope parameter or in the photometry) were assigned this minimum value (Frogel, Persson, & Cohen, 1980; Elson & Walterbos, 1988).

### 3.3 Abundance Determinations

To rank the M31 clusters in the WYFFOS sample by metallicity, we have measured the strengths of various absorption features in their integrated spectra. Line indices were calculated to measure the signal in a wavelength-delineated region centered on each spectral feature relative to the signal in the blue and red continuum zones flanking the feature. In light of the fact that no flux standard star spectra were obtained during the WYFFOS observing run, it was not possible to convert photon counts directly into flux for these targets. Our motivation, however, was to obtain a metallicity ranking of the globular clusters in M31; relative flux calibration was not essential as its absence will not significantly affect narrow features and mean count ratios were used instead. Our line and continuum bandpasses were defined as shown in Table 4 (Burstein et al., 1984; Brodie & Hanes, 1986; Brodie & Huchra, 1990).

Each target spectrum was shifted to zero radial velocity to define consistent line index wavelengths in the rest frame. The velocity dispersion of the M31 globular cluster system is roughly $`150`$ km s-1, a spread which translates into a nominal range of $`4\sigma _v600`$ km s-1 in velocity. This would correspond to a wavelength shift of $`9`$Å at $`\lambda =4500`$ Å, a significant fraction of the width of most of the bandpasses selected for measurement. Thus, metallicity determinations were not pursued for target spectra which had yielded no velocity cross-correlations.

The generalized feature index is defined as

$$I=2.5\mathrm{log}\left[\frac{2\overline{}_I}{\overline{}_{C1}+\overline{}_{C2}}\right]$$
(3)
where $`\overline{}_I`$ , $`\overline{}_{C1}`$ and $`\overline{}_{C2}`$ represent the mean count levels in the feature bandpass and associated continuum bandpasses $`C1`$ and $`C2`$ as given by $`\overline{}=(\lambda _2\lambda _1)^1_{\lambda _1}^{\lambda _2}𝑑\lambda `$ (Brodie & Hanes, 1986; Brodie & Huchra, 1990). The continuum break feature $`\mathrm{\Delta }=2.5\mathrm{log}\left[\overline{}_{C2}/\overline{}_{C1}\right]`$ is a measure of the discontinuity in the continuum due to Fraunhofer line-blanketing shortward of 4000 Å. This feature is somewhat more vulnerable to errors in the reddening and absorption corrections, as it spans a broad range in wavelength and is asymmetric.

Figure 6 shows a comparison of our line indices with published M31 globular cluster \[Fe/H\] values from various sources (Bònoli et al., 1987; Huchra, Brodie, & Kent, 1991; Barmby et al., 2000). The RMS differences between the line indices calculated from spectra observed in more than one field through different fibres were used to estimate the internal errors for each feature, and these are shown in the second column of Table 5. The results of linear fits in the form $`LI=a+b[\mathrm{Fe}/\mathrm{H}]`$ weighted by $`1/\sigma _{[\mathrm{Fe}/\mathrm{H}]}^2`$ are provided in columns 3 and 4 of Table 5. Note that only spectroscopic metallicities for overlapping targets were considered in the calibration and hence the Bònoli et al. (1987) near-infrared photometric metallicities were not used in the fits. In column 5 of Table 5, we provide the RMS residuals of the fit and in column 6 we list the strength of the $`r`$ coefficient, a measure of the linear correlation between the line indices and published metallicities.

To best represent the metal abundances of the WYFFOS targets observed, a line index was selected to contribute to the calibration if it met the following four criteria:

1. The uncertainties associated with the line index did not exceed $`20\%`$ of the full range of observed values. The Fe52 index did not satisfy this condition and was thus rejected from the calibration.
2. The linear correlation coefficient $`r`$ exceeded a minimum level of $`r=0.65`$.
3. Numerical experimentation demonstrated that the feature strength was stable in the presence of significant changes in the adopted value of cluster reddening (e.g. in comparison with any adopted reddenings which differed from unpublished M31 values from Barmby et al. (2000)). The continuum break was found to be susceptible to reddening uncertainties and was thus rejected.
4. Only one line index that characterizes a given spectral feature was incorporated in the estimate. In the case of the magnesium triplet, the indices MgG, Mg2 and Mgb all provided an acceptable measure of essentially the same absorption feature. Of these, the Mgb index demonstrated the highest correlation coefficient and a relatively low RMS uncertainty and was therefore selected for use over the others.

Line indices which met the above criteria include the CH (G-band) feature, Mgb to characterize the strength of the magnesium $`b`$ triplet, and the Fe53 iron line. Final cluster metallicities have been determined from an unweighted mean of the \[Fe/H\] values calculated from the CH(G), Mgb and Fe53 line strengths for 202 of the WYFFOS targets. These metallicities and their associated errors are provided in columns 7 and 8 of Table 3. A plot of the correlation between our \[Fe/H\] values and those in other samples from the literature is presented in Figure 7. The best-fit of the data in Figure 7 has a slope of $`0.94\pm 0.02`$ and the RMS of the fit residuals is $`0.24`$. The median value of the formal errors on our WYFFOS metallicities is $`\pm 0.26`$ dex.

### 3.4 The Best Current Sample

In the subsequent analysis, our results were merged with the velocities and \[Fe/H\] values from other sources in order to produce the best available sample of spectroscopic data. Wherever data from other publications overlapped with our WYFFOS sample and demonstrated smaller associated uncertainties, these data superseded the WYFFOS results; such cases are noted in the final column of Table 3. In this final sample, 191 WYFFOS velocities and 189 metallicities were combined with data from Peterson (1989), Huchra, Brodie, & Kent (1991), Federici et al. (1993), Dubath & Grillmair (1997), and Barmby et al. (2000) to yield a total M31 globular cluster spectroscopic database comprising 321 velocities and 301 metallicities.

Barmby et al. (2000) identify the following targets as having spectral properties which are inconsistent with globular clusters: B55-S116 (stellar), B308 (galaxy), B341-S81 (stellar) and B392-S329 (stellar). We concur with the identification of B308 as a definite non-cluster and with that of B392-S329 as a probable non-cluster, and do not carry these through the remainder of the analysis. We do not find any clear abnormalities in the spectra of the other targets and thus cannot at present reject these as confirmed non-clusters.

In an effort to verify that there are no obvious biases in the metallicity sample studied here, we examine the color distribution of our sample and compare it to that of the overall M31 GCS. The mean estimated intrinsic color of 412 M31 globular clusters with $`(BV)`$ data from Barmby et al. (2000) was $`(BV)_0=0.67\pm 0.04`$, while that for the sample of those with \[Fe/H\] values used in the metallicity analysis presented herein (264 clusters with colors) was determined to be $`(BV)_0=0.65\pm 0.05`$. Based on this and a comparison of the color histograms for the metallicity sample and the overall GCS (Figure 8), we conclude that there are no significant population biases with respect to the color of the metallicity sample presented herein. A Kolmogorov-Smirnov two-sample test on the distributions in Figure 8 confirms what the eye suggests: the color distribution of the metallicity sample is drawn from the same parent population as the overall color sample at a high confidence level (99.5%).

Clearly, the globular clusters at small galactocentric radii will suffer from more severe effects of incompleteness due to obscuration by the disk and bulge of M31. Depending on the color/metallicity distribution as a function of position within the galaxy, this observational effect will conceivably bias the observed metallicity sample. We will return to this point in Section 4.1.

## 4 M31 Cluster Metallicities

The metal abundance properties of M31’s globular clusters are of interest as they bespeak the galaxy’s formative processes and enrichment history. The metallicity distribution for the M31 GCS is presented in Figure 9, along with a similar plot for the Milky Way GCS for comparison. The mean metallicity of the full complement of 301 M31 globular clusters is $`[\mathrm{Fe}/\mathrm{H}]=1.21\pm 0.03`$, comparable to the value of $`1.27\pm 0.05`$ obtained for the Milky Way GCS (based on data supplied in Harris (1996), revised 1999). However, a Kolmogorov-Smirnov test reveals that the Galactic and M31 globular clusters have \[Fe/H\] distributions which are drawn from the same parent population at only the 55% confidence level.

As is frequently the case with binned data, the appearance of the histogram can be ambiguous and potentially misleading, and thus a more robust method of analysis is desirable. In order to examine the shape of the metallicity distribution without relying on binning methods, we turned to the KMM algorithm used by Ashman, Bird, & Zepf (1994). KMM mixture modeling operates under the assumption that the data is independently drawn from a parent population which consists of a mixture of $`N`$ Gaussian distributions. The KMM test has previously been applied to the M31 globular cluster metallicity distribution by Ashman & Bird (1993), albeit with a significantly smaller data sample.

A bimodal distribution was tested against the assumption of a unimodal shape (the null hypothesis) for the GCS metallicities. Since there is no obvious reason to assume that both \[Fe/H\] groups in a bimodal fit would exhibit a common covariance, heteroscedastic fitting was applied such that the covariances for the two distributions were not constrained to be identical. The KMM results revealed that this 2-group fit demonstrated a marginal improvement over a single group fit, although Ashman, Bird, & Zepf (1994) cautioned that the output probability in the heteroscedastic case is more difficult to interpret and may, in and of itself, be suspect. We repeated the KMM test on more than 1500 simulations with bootstrap resampling and found that the bimodal case was preferred over the unimodal distribution at a median confidence level of higher than $`97\%`$. The bimodal test returned Gaussian fits characterized by means at $`[\mathrm{Fe}/\mathrm{H}]=1.44`$ ($`\sigma ^2=0.22`$) for the metal-poor distribution and at $`[\mathrm{Fe}/\mathrm{H}]=0.50`$ ($`\sigma ^2=0.13`$) for the metal-rich clusters. The a posteriori probabilities of group membership returned by KMM assigned 231 clusters to the metal-poor distribution and 70 to the metal-rich population.

In their analysis, Elson & Walterbos (1988) treated the M31 GCS as comprising three distinct metallicity sub-populations: metal-poor, intermediate, and metal-rich. We employed the KMM test to ascertain whether a trimodal distribution was also consistent with the measured \[Fe/H\] values for the M31 cluster system. The trimodal KMM results with bootstrap resampling demonstrated that a 3-group fit was indeed preferred over a single Gaussian distribution at a relatively high confidence level (96%). The output population means obtained for the best 3-group fit were located at $`[\mathrm{Fe}/\mathrm{H}]=1.59`$ ($`\sigma ^2=0.19`$), $`1.18`$ ($`\sigma ^2=0.18`$) and $`0.47`$ ($`\sigma ^2=0.12`$). The numbers of clusters assigned to each group were $`N=`$153, 76 and 72, respectively.

In light of the variation in the magnitude of \[Fe/H\] uncertainties quoted within the WYFFOS sample and others, it is appropriate to test the effects of these uncertainties on the inferred shape of the distribution. There are 229 GCs with metallicity uncertainties lower than 0.45 dex, a value roughly equal to the variance of one of the populations as determined from the bimodal results. We applied the KMM test to this higher-precision sample to test the stability of our earlier conclusions. The bimodal case was again found to represent a statistically significant improvement over the unimodal assumption with peaks similar to those found when we incorporated the complete sample. The input assumption of three peaks, however, did not return adequate fits: the a posteriori probabilities of group membership in the trimodal case allocated only two clusters to the intermediate population.

There is an unavoidable smearing of the metallicity distribution due to metallicity calibration uncertainties, absorption and reddening of bulge/disk clusters and other line-of-sight problems, or errors incurred as a result of combining different datasets. As a further test, simulations were run in which we shuffled the high-precision \[Fe/H\] sample by a random fraction of their associated uncertainties. This returned a result of marginal significance (median p-value $`0.1`$) for the 2-group case, and an insignificant result (p-value $`0.4`$) for the 3-group heteroscedastic trial. Note that the shapes of the 2-group Gaussian fits were once again consistent with those of the original bimodal test using the full sample.

In light of these KMM results, we tentatively reject the possibility of a trimodal GC metallicity distribution for the purposes of the analysis which follows. While it may yet be reasonable to separate the cluster system into more than two populations — say, by culling out the intermediate targets with group membership probabilities between 0.4 and 0.6 in the bimodal case — we elect not to do so here, as it is not clear that these represent a bona fide intermediate population. Based on the available data, we argue that the metallicity distribution of the M31 globular cluster system is not unimodal and better resembles a bimodal distribution similar in general shape to that of the Milky Way (Zinn, 1985). Therefore, “metal-poor” and “metal-rich” populations shall refer to those objects separated into two groups based on the a posteriori probabilities returned from the 2-group heteroscedastic KMM test with the complete available sample. Each population was independently fit using an adaptive kernel smoothing procedure333Adaptive kernel smoothing code kindly provided by Karl Gebhardt. and the results are shown in Figure 10.

The positions of our adopted metallicity peaks ($`[\mathrm{Fe}/\mathrm{H}]=1.44`$ and $`0.50`$) differ only slightly from those quoted by Barmby et al. (2000), who found a bimodal distribution in M31 GC metallicities and metallicity-sensitive colors with peaks at $`[\mathrm{Fe}/\mathrm{H}]=1.4`$ and $`0.6`$. Our peaks are also similar to the values of $`[\mathrm{Fe}/\mathrm{H}]=1.5`$ and $`0.6`$ obtained by Ashman & Bird (1993) using the KMM test on an earlier sample of 144 M31 globular clusters; their study yielded a ratio of metal-rich to metal-poor globular clusters of 49:95, in comparison with a ratio of 70:231 found in the current study.

### 4.1 Spatial Distribution

The projected spatial distributions of the metal-poor and metal-rich cluster populations are shown in Figure 11. The metal-rich population is clearly more centrally concentrated, an effect which has been demonstrated previously with a smaller dataset (Huchra, Brodie, & Kent, 1991). In the Milky Way, the metal-rich clusters have been ascribed as belonging to a thick disk system (Zinn, 1985; Elson & Walterbos, 1988; Armandroff, 1989). Later work has suggested that this is better described as a bulge/bar cluster system (Minniti, 1995; Côté, 1999; Forbes, Brodie & Larsen, 2001). The spatial concentration and distribution of the metal-rich clusters in M31 appear to be consistent with the latter interpretation.

There is a notable deficit of metal-poor clusters in the innermost radial bins shown in Figure 12. The already strong central concentration in the metal-rich population might mask any analogous deficit in this population. A similar lack of central clusters was reported by Elson & Walterbos (1988) within the inner $`15\mathrm{}`$ of the GCS and is at least partially attributable to sampling restrictions towards the middle of the galaxy, as well as to observational biases incurred as a result of obscuration by the bulge and disk. The cluster sample at large galactocentric radii is bound to be more complete at faint magnitudes due to the lower degree of absorption as compared with the inner cluster system.

Examination of Figure 13 reveals a modest lack of clusters fainter than $`V17`$ within the inner $`10\mathrm{}`$ or so of the GCS. We attempt to estimate the incompleteness at small radii in the following manner. Using available $`V`$ photometry (Barmby et al., 2000), we plot the observed luminosity functions (LFs) for the metal-rich and metal-poor clusters beyond a galactocentric radius of $`15\mathrm{}`$. These outer LFs are normalized to the same scale as the distributions within the inner $`15\mathrm{}`$ by summing under the brighter end of the inner and outer histograms ($`V<16.75`$ mag) and multiplying by the ratio of these values. The observed histograms for the inner metal-rich and metal-poor clusters are shown by the solid lines in Figure 14, and the scaled outer histograms are represented by the dotted lines. Note that this method provides only a very crude estimate of the number of missing clusters at $`V16.75`$ as it presumes, perhaps falsely, that the shape of the globular cluster LF does not change within the inner region of the galaxy. Dynamical destruction will have a more severe effect on the inner clusters and will thus change the overall shape of the luminosity function. Evidence for this effect has previously been noted by Kavelaars & Hanes (1997).

From the scaled outer LFs in Figure 14, we infer that there are possibly $`30`$ missing metal-poor clusters at $`V16.75`$, whereas there may be only $`8`$ faint, metal-rich clusters which are missing within the inner $`15\mathrm{}`$ of the galaxy. Therefore, sample incompleteness is a reasonable explanation for the dip in the central bins of the metal-poor radial profile.

Comparing simple counts of the cluster populations within quadrants centered on the major and minor axes reveals that there is a slight enhancement in the minor axis fields for the metal-rich population: 40 clusters in the East-West sectors versus 30 in the North-South. If one takes into account only the inner $`30\mathrm{}`$ of the GCS, this enhancement largely vanishes for the metal-rich clusters: 28 clusters are found within the major axis sectors compared to 26 along the minor axis. Based on their photometric metallicity estimates, Elson & Walterbos (1988) discovered that their sample of clusters classified as metal-rich defines a flattened system with an ellipticity of $`ϵ>0.4`$ within a major axis radius of $`30\mathrm{}`$. Our cluster counts are not consistent with such a large degree of flattening of the metal-rich system in this region. At present, we cannot provide reasonable limits on the flattening of the metal-poor cluster distribution, since the arrangement of the WYFFOS fields has resulted in a selection bias in favor of clusters in major axis sectors at radii beyond $`30\mathrm{}`$.

### 4.2 Metallicity Gradient

Several previous works have suggested the presence of a weak metallicity gradient in the M31 cluster system out to radii of $`10`$ kpc (Sharov, 1988; Huchra, Brodie, & Kent, 1991). In Figure 15, we plot the metallicity of the M31 globular clusters as a function of galactocentric radius for the full sample and the separate metallicity sub-populations. The large scatter in these plots is due at least in part to the fact that we can only measure projected positions, but nonetheless we find that the mean \[Fe/H\] does indeed decrease with radius.

A sliding bin method was used to fit the mean metallicity in $`12\mathrm{}`$ bins centered on each data point, with results shown in Figure 15. Robust estimates of the mean metallicity in radial bins are consistent with the presence of a modest gradient in \[Fe/H\] within the full sample. Separating the data into sub-populations reveals no clear evidence for a gradient in the metal-rich GCS, although the spread in metallicity does seem to increase towards the center of the galaxy. There is an apparent trend towards decreasing metallicity within both the full sample and the inner metal-poor cluster system, with gradients of roughly $`0.017`$ and $`0.015\mathrm{dex}/\mathrm{arcmin}`$, respectively. Relatively poor sampling of the clusters at large radii makes it difficult to comment decisively on the results beyond a radius of about $`60\mathrm{}`$, although the metal-poor sample appears to level off at a lower mean metallicity beyond this distance.

## 5 GCS Kinematics

The maximum likelihood value for the mean cluster radial velocity is $`v=284\pm 9`$ km s-1, somewhat higher than M31’s systemic velocity of $`300\pm 4`$ km s-1 (RC3). The globular clusters have an overall velocity dispersion of $`\sigma _v=156\pm 6`$ km s-1. The kinematics of metal-rich and metal-poor GCS populations will be discussed later in this section. First, we examine the global kinematic properties of the M31 GCS.

### 5.1 GCS Rotation

A histogram of the cluster line-of-sight velocities is presented in Figure 16. There is strong evidence for bimodality in the velocity distribution, a feature which is attributable to the fact that the cluster system is rotating with respect to our line of sight. Considering, for the moment, a purely rotating disk population of globular clusters in the reference frame of the galaxy:

$`v_{\mathrm{GC}}`$ $`=`$ $`v_{\mathrm{rot}}\mathrm{sin}(i)\mathrm{sin}(\theta +\delta )+v_\sigma `$
$`=`$ $`v^{}+v_\sigma `$
where $`v_{\mathrm{rot}}`$ is the rotation velocity of the population, $`i`$ is the inclination of the rotation axis to the line of sight ($`77.7^{}`$), $`\theta `$ is the angle about the angular momentum axis, $`\delta `$ is an arbitrary phase offset, and $`v_\sigma `$ incorporates the velocity dispersion plus random measurement errors. The combined probability distribution $`P(v_{\mathrm{GC}})`$ for the M31 globular cluster radial velocities is the convolution of $`P(v^{})`$ with $`P(v_\sigma )`$. The probability distribution $`P(v_\sigma )`$ is Gaussian and $`P(\theta )`$ is uniform, thus $`P(v^{})`$ is given by:

$`P(v^{})`$ $`=`$ $`P(\theta ){\displaystyle \frac{d\theta }{dv^{}}}`$ (5)
$`=`$ $`{\displaystyle \frac{1}{2\pi }}{\displaystyle \frac{1}{\sqrt{v_{\mathrm{rot}}^2\mathrm{sin}^2iv_{}^{}{}_{}{}^{2}}}}.`$
This represents a bimodal distribution with maxima at $`\pm v_{\mathrm{rot}}\mathrm{sin}(i)`$ about the systemic velocity of the population. The peaks of the (unshifted) velocity histogram in Figure 16 lie at $`415`$ and $`170`$ km s-1, yielding an estimated rotation amplitude of $`v_{\mathrm{rot}}125`$ km s-1.

The fact that M31 is rotating has been established by various studies (see Hodge, 1992, and references therein). Using a sample of $`150`$ M31 clusters as dynamical tracers, Huchra, Brodie, & Kent (1991) noted a rapid rotation of the handful of inner ($`R<7\mathrm{}`$) clusters in their sample with $`[\mathrm{Fe}/\mathrm{H}]0.8`$ (see also Huchra, Stauffer, & Van Speybroeck (1982)). Huchra, Brodie, & Kent (1991) discovered that the metal-poor population had no rotation signature within the inner galaxy, but shared a moderate rotation of $`60`$ km s-1 with the bulk of the GCS beyond this radius. With our increased sample of 321 velocities, we can now provide a more comprehensive look at the kinematics of the M31 globular cluster system.

Perhaps not surprisingly, there appears to be no significant rotation of the GCS about the major axis of M31 as shown in the left panel of Figure 17. A robust fit to the targets within the inner $`50\mathrm{}`$ along the minor axis (dropping 11 objects at large radii) gives $`v_r=286+0.13Y`$ km s-1, where $`Y`$ is the distance along the minor axis in arcminutes. Such a small deviation from a slope of zero may simply be due to statistical imprecision or variation in the position angle for the galaxy (Hodge & Kennicutt, 1982).

A plot of the observed cluster radial velocities as a function of projected distance along the major axis of M31 is shown in the right panel of Figure 17. A basic $`2\sigma `$ culled robust linear fit to the velocities along the major axis reveals that the globular cluster system indeed has a significant rotational component about the minor axis. This rotation will now be investigated in more detail.

We used an adaptive binning technique to generate the GCS rotation curve shown in Figure 18. Before fitting the rotation curve, the globular cluster coordinates were transformed into the reference frame of the galaxy and the data were folded about each axis. The function $`v=v_{\mathrm{sys}}+v_{\mathrm{rot}}\mathrm{sin}\varphi `$ was fit to the radial velocities of the clusters within radial bins outwards from the galaxy center (Kissler-Patig & Gebhardt, 1998). The angle $`\varphi `$ represents the position angle of the target object in the projected frame of the galaxy, and $`v_{\mathrm{sys}}`$ is the mean velocity of the cluster system. The kinematic properties were determined in a series of annuli centered on each cluster and encompassing a specified number of data points (79 points, with a minimum of 39 at the extrema). With this method, we obtain a rotation amplitude of $`136`$ km s-1 where the curve flattens off beyond $`R=30\mathrm{}`$.

We then employed the robust biweight estimators of location and scale (Beers, Flynn, & Gebhardt, 1990; Hoaglin, Mosteller, & Tukey, 2000) as measures of mean velocity and dispersion in 10 arcminute bins along the major axis of the galaxy. The biweight results and their 90% confidence intervals (based on bootstrap errors) are presented in Figure 19 for the full sample and the two metallicity sub-populations. Based on the appearance of these plots, it seems reasonable to approximate the shape of the curve as consisting of solid-body rotation within the central $`25\mathrm{}`$, with a flattening beyond this radius. One-half of the velocity difference between the flat regions of the curve in Figure 19 yields a rotation amplitude of $`v_{\mathrm{rot}}=138\pm 13`$ km s-1 for the full sample.

Significant differences between the kinematic properties of the metal-rich and metal-poor clusters may indicate that these populations underwent dissimilar formation mechanisms. The maximum likelihood value for the mean velocity of the 70 clusters in the metal-rich sample is $`v=260\pm 18`$ km s-1 with a dispersion of $`\sigma _v=146\pm 12`$ km s-1. For comparison, the larger metal-poor sample of 231 clusters demonstrates a somewhat lower mean velocity of $`v=290\pm 10`$ km s-1 and a slightly higher dispersion of $`\sigma _v=155\pm 7`$ km s-1, although these values are still consistent with those of the higher metallicity sample within the quoted errors. From the fits of the biweight results, we determine a rotation amplitude of $`v_{\mathrm{rot}}=160\pm 19`$ km s-1 for the metal-rich population. This value is somewhat larger than the metal-poor amplitude of $`v_{\mathrm{rot}}=131\pm 13`$ km s-1, although these values are also not inconsistent within their formal errors.

For each cluster, the projected position angle was determined as the angle of maximum positive rotation obtained from the maximum likelihood sinusoid fits (Kissler-Patig & Gebhardt, 1998). Plots of the variation of position angle with radius for the full sample and metallicity sub-populations are presented in Figure 20. The input GC positions already incorporate a de-rotation of $`\mathrm{PA}_{\mathrm{M31}}=38^{}`$ for the galaxy and, all else being equal, one would expect to find that the PA for the cluster system is roughly zero. Despite the sizeable spread at low radius in the full sample, there is a relatively clear increase in PA within the central region of the metal-poor cluster system. This effect is consistent with studies of the central bulge and inner disk region of M31 which demonstrate a similar increase in PA (Hodge & Kennicutt, 1982; Kent, 1989). The M31 GCS shows a notable decrease in the PA of the rotation axis beyond about $`30\mathrm{}`$. Interestingly, the metal-rich clusters seem to be rotating about an axis that is tilted at a slightly lower PA with respect to the major axis throughout the full radial extent of the population.

### 5.2 The Mass of M31

The kinematics of M31’s GCS provide us with a dynamical probe of the underlying mass distribution of the host galaxy (Federici et al., 1990, 1993). The projected mass estimator (PME) was used to place constraints on the underlying mass of the system (Bahcall & Tremaine, 1981; Heisler, Tremaine, & Bahcall, 1985):

$$M_{\mathrm{PME}}=\frac{f_p}{NG}\underset{i=1}{\overset{N}{}}r_iv_i^2,$$
(6)
where $`N`$ is the number of bodies, G is the gravitational constant, $`r`$ is the projected galactocentric radius and $`v`$ is the radial velocity of the body relative to the systemic velocity of the galaxy. The factor $`f_p`$ depends on the assumptions made regarding the distribution of orbits, and is equal to $`32/\pi `$ and $`64/\pi `$ for isotropic and radial orbits, respectively.

Admittedly, this projected mass estimate is hampered by our limited sample size and coverage, as well as a lack of detailed understanding of the cluster orbits. Recent work by Evans & Wilkinson (2000) has shown that the dynamical tracers in the M31 halo exhibit a predominantly isotropic velocity dispersion. Applying the PME to the available cluster data and assuming isotropic orbits within an extended mass distribution, we estimate a mass of $`M_{\mathrm{M31}}=(4.1\pm 0.1)\times 10^{11}M_{}`$ using 319 dynamical targets out to a radius of $`120\mathrm{}`$ from the galaxy center ($`27`$ kpc at a distance of 770 kpc). This result is comparable to previous mass estimates obtained out to similar radii as summarized in Table 6 of Evans & Wilkinson (2000), and to the mass of $`M_{\mathrm{M31}}(3.1\pm 0.5)\times 10^{11}M_{}`$ found by Huchra, Brodie, & Kent (1991) based on their sample of 150 globular clusters. In the (unlikely) case of purely radial orbits, our mass estimate would be increased to twice the value quoted above.

## 6 Discussion and Conclusions

We have presented a sample of new spectroscopic observations of M31 globular clusters, contributing $`200`$ high-precision velocities and \[Fe/H\] values for what represents a large fraction of the dynamical tracers of the galaxy. A comprehensive analysis of the kinematics and abundance properties of M31’s GCS demonstrates some interesting similarities and differences when compared with the Milky Way’s globular cluster population.

The Galactic and M31 globular cluster systems have comparable mean metallicities of $`[\mathrm{Fe}/\mathrm{H}]1.2`$ to $`1.3`$, although a Kolmogorov-Smirnov test provides no conclusive evidence that the two distributions are drawn from the same parent population. A KMM mixture-modeling test reveals that the M31 GCS metallicity distribution can convincingly be separated into two groups, a bimodality which is analogous to that found in the Milky Way GCS. The peaks in the M31 metallicity distribution lie at $`[\mathrm{Fe}/\mathrm{H}]=1.44`$ and $`0.50`$. A KMM test on the Milky Way GCS metallicities from the McMaster catalogue (Harris, 1996) returns peaks at $`[\mathrm{Fe}/\mathrm{H}]=1.59`$ and $`0.56`$, indicating that the metal-poor Galactic globulars have, on average, a somewhat lower metal abundance than their M31 counterparts. Of the 301 M31 globular clusters with available metallicities, KMM assigns 70 to the metal-rich peak. This represents a somewhat smaller fraction of high-metallicity clusters than is found in the Milky Way system, for which 45 of 145 are identified with the metal-rich population.

M31’s metal-rich globular clusters appear to constitute a distinct kinematic subsystem that is quite spatially concentrated, consistent with membership in a bulge population akin to the inner metal-rich clusters of the Milky Way (Minniti, 1995). Furthermore, the velocity dispersion of the metal-rich cluster system, $`\sigma _v=146\pm 12`$ km s-1, is indistinguishable from the mean bulge dispersion of $`\sigma _v(\mathrm{bulge})=146\pm 6`$ km s-1 as determined from the kinematics of faint planetary nebulae (Lawrie, 1983; van den Bergh, 1991). The mean metallicity of this population is quite similar to their Milky Way counterparts and demonstrates no obvious abundance gradient with galactocentric radius. Contrary to observations in the Milky Way, however, we find no clear signs of flattening in the metal-rich component of the M31 GCS despite our relative oversampling of major axis fields in the current spectroscopic sample. We find that the rotation axis of the metal-rich cluster population is tilted at a slightly lower position angle: it is offset by about $`510^{}`$ with respect to that of the minor axis of M31 and the bulk of its GCS. This observation may point to a disparate origin for these clusters. Further searches for potential cluster candidates hidden by the bulge of M31 will be helpful in providing a clearer understanding of the kinematics and nature of this inner system.

In M31, there are more than three times the number of metal-poor clusters as there are metal-rich. It is not clear if it is possible (or even necessary) to further separate this larger, more spatially distended metal-poor population into sub-components. If the metal-poor cluster distribution does indeed incorporate a thick disk component in addition to halo clusters, this would account for the relatively high net rotation amplitude and large velocity dispersion observed therein. While grouping the clusters based on metallicity remains an instructive exercise, \[Fe/H\] alone is not an adequate discriminator of membership in the spatially- or kinematically-defined components of the galaxy. This is primarily a result of the substantial spread and measurable gradient in \[Fe/H\] within the metal-poor cluster system. The task of assigning globular clusters to their appropriate component (e.g. bulge, thin or thick disk, inner or outer halo) remains somewhat simpler within the Milky Way GCS, since here we have three-dimensional spatial information rather than projected positions.

We find evidence for a radial metallicity gradient ($`0.015\mathrm{dex}/\mathrm{arcmin}`$) in the metal-poor population of M31 out to about $`60\mathrm{}`$ from the galaxy center. Despite a notable lack of cluster data beyond this radius, it does seem that mean metallicity levels off, although the scatter remains large. This gradient in the metal-poor clusters is not inconsistent with a single collapse scenario, yet we recognize that a metallicity gradient, in and of itself, is not sufficient to confirm an early dissipational collapse. Such a gradient could be mimicked by a hierarchical halo formation if the metallicity of the accreting fragments correlates with their mass and mean density (Freeman, 1996). The search for substructure in position/velocity/metallicity parameter space is clearly of interest to investigate the accretion history of the galaxy. An analysis of potential sub-clustering within the M31 GCS will follow in a subsequent paper.

A comprehensive investigation of the kinematics of the dynamical tracers in M31 would not be complete without the addition of planetary nebulae data. Fortunately, a study of a large sample of M31 is currently underway (Halliday et al., in preparation. Moreover, these data can be combined with the globular cluster sample to better probe the underlying mass distribution of the galaxy (Wilkinson et al., in preparation). An analysis of M31 cluster ages is also in progress (T. Bridges, private communication) which should provide a better sense of the chronology of GC formation within the galaxy. Despite its long and rich history of investigation, there is little doubt that there are missing globulars in M31’s database which are awaiting discovery and identification. Continuing searches and surveys for M31 cluster candidates (Lee et al., 2001) promise to help remedy this situation and to further fuel this “booming industry”. It is clear that detailed studies of the kinematics, metallicities and ages of the globular cluster systems of a large sample of galaxies can provide us with the clues necessary for a more complete understanding of the mechanisms which form and evolve galaxies of all kinds.

## 7 Acknowledgments

This work was supported in part by operating grants to DAH from the Natural Sciences and Engineering Research Council of Canada. JPB and JPH acknowledge funding support from the National Science Foundation. The authors wish to extend their thanks to Jim Lewis, Karl Gebhardt and Tim Beers for making their software available, and to Pauline Barmby for providing unpublished reddenings. The authors are also very grateful to the anonymous referee for the helpful comments and to the kind staff at WHT for their support during the observations.


format: text

text:
# Quantum Interference between Impurities: Creating Novel Many-Body States in 𝑠-wave Superconductors

## Abstract

We demonstrate that quantum interference of electronic waves that are scattered by multiple magnetic impurities in an $`s`$-wave superconductor gives rise to novel bound states. We predict that by varying the inter-impurity distance or the relative angle between the impurity spins, the states’ quantum numbers, as well as their distinct frequency and spatial dependencies, can be altered. Finally, we show that the superconductor can be driven through multiple local crossovers in which its spin polarization, $`s_z`$, changes between $`s_z=0,1/2`$ and $`1`$.

Over the last two years, several beautiful experiments have studied quantum interference of electronic waves that are scattered by multiple impurities Man00 ; Der02 ; Chi02 ; Hol01 . In a groundbreaking experiment, Manoharan et al. Man00  used a corral of magnetic impurities on the surface of a metallic host to demonstrate that quantum interference can lead to the focusing of electronic waves into a quantum image. Moreover, using scanning tunneling spectroscopy (STS), Derro et al. Der02  were the first to observe four resonance states in the local density of states (DOS) of the one-dimensional chains in YBa2Cu3O6+x. These states were interpreted as arising from quantum interference of electronic waves scattered by two magnetic impurities Morr01 . Quantum interference effects were also studied in optical quantum corrals by Chicanne et al. Chi02  and between impurities located on quantum dots by Holleitner et al. Hol01 . Some first theoretical work theory  has focused on impurity geometries in metallic systems similar to the one studied by Manoharan et al. In contrast, quantum interference in strongly correlated electron systems, such as superconductors (with the exception of NbSe2 Flatte ), charge- and spin-density-wave systems, or even semi-conductors, have not yet been addressed. However, the study of interference effects in theses systems involving spin impurities is not only of great fundamental interest, but might also possess important applications in the field of spin electronics Wolf01  and quantum information technology Kane98 .

In order to describe the properties of complex impurity structures such as quantum corrals, it is first necessary to understand interference effects associated with the presence of few impurities. In this Letter we therefore consider two impurities embedded in a general $`s`$-wave superconductor (SC). The presence of two magnetic impurities allows for a coupling of the bound states associated with a single impurity Shiba68 , and gives rise to the emergence of novel many-body states. We show that the nature of these novel states, i.e., their quantum numbers, can be altered by varying the distance between the two impurities, $`\mathrm{\Delta }r`$, or the relative angle between the directions of their spin moments, $`\alpha `$. Moreover, we demonstrate that these changes are accompanied by local crossovers in which the spin polarization of the superconductor changes between $`s_z=0,1/2`$, and $`1`$. We predict that the interplay between the states’ quantum numbers and the inter-impurity distance determines the distinct frequency and spatial dependence of the two-impurity bound states. Finally, we discuss the implications of our work for systems with a larger number of impurities.

Starting point for our calculations is the $`\widehat{T}`$-matrix formalism Shiba68  which we generalized to treat the case of $`N`$ impurities of spin $`S`$ with non-magnetic and magnetic scattering potentials Morr01 ; Balatsky . In the following, we focus on the case $`N=2`$, and, following Ref.Shiba68 , treat the impurity spins as classical, static variables, corresponding to the limit $`\beta _0=JS/2=const.`$ and $`S\mathrm{}`$. In a fully gaped s-wave SC, this approximation is well justified since no Kondo-effect occurs for sufficiently small coupling between the impurities and the delocalized electrons. Within this approach, any interaction between the impurities is only important to the extent that it determines the angle, $`\alpha `$, between the direction of the impurity spins. Within the Nambu-formalism and in Matsubara frequency space the electronic Greens function in the presence of $`N`$ impurities is given by

$`\widehat{G}(r,r^{},\omega _n)`$ $`=`$ $`\widehat{G}_0(r,r^{},\omega _n)`$ (1)
$`+{\displaystyle \underset{i,j=1}{\overset{N}{}}}\widehat{G}_0(r,r_i,\omega _n)\widehat{T}(r_i,r_j,\omega _n)\widehat{G}_0(r_j,r^{},\omega _n),`$
where the $`\widehat{T}`$-matrix is obtained from the Bethe-Salpeter equation

$`\widehat{T}(r_i,r_j,\omega _n)`$ $`=`$ $`\widehat{V}_{r_i}\delta _{r_i,r_j}`$ (2)
$`+\widehat{V}_{r_i}{\displaystyle \underset{l=1}{\overset{N}{}}}\widehat{G}_0(r_i,r_l,\omega _n)\widehat{T}(r_l,r_j,\omega _n).`$
In the case of two impurities

$`\widehat{V}_{r_1}={\displaystyle \frac{1}{2}}\left(U_1\sigma _0+J_1S\sigma _3\right)\tau _3;`$
$`\widehat{V}_{r_2}={\displaystyle \frac{1}{2}}\left(U_2\sigma _0+J_2S\sigma _3\mathrm{cos}\alpha +J_2S\sigma _1\mathrm{sin}\alpha \right)\tau _3;`$
$`\widehat{G}_0^1(𝐤,i\omega _n)=\left[i\omega _n\tau _0ϵ_𝐤\tau _3\right]\sigma _0+\mathrm{\Delta }_𝐤\tau _2\sigma _2.`$ (3)
Here, $`\widehat{V}_{r_1,r_2}`$ is the scattering matrix for the impurities located at $`𝐫_1`$ and $`𝐫_2`$, respectively. Without loss of generality, we take the spin of impurity $`1`$ to be parallel to the $`\widehat{z}`$-axis, while that of impurity $`2`$ is rotated from the $`\widehat{z}`$-axis into the $`zx`$-plane by an angle $`\alpha `$. $`\widehat{G}_0(𝐤,i\omega _n)`$ is the Greens function of the unperturbed (clean) system in momentum space, and $`\sigma _i`$, $`\tau _i`$ are the Pauli-matrices in spin and Nambu-space, respectively. $`U_i`$ and $`J_i`$ are the potential and magnetic scattering strengths of the impurities. We consider a two-dimensional (2D) electronic system whose normal state dispersion is given by $`ϵ_𝐤=k^2/2m\mu `$ ($`\mathrm{}=1`$), where $`\mu =k_F^2/2m`$ is the chemical potential, and $`k_F=\pi /2`$ is the Fermi wave-vector (we set the lattice constant $`a_0=1`$). The results and conclusions presented below are qualitatively robust against changes in the form of $`ϵ_𝐤`$, the dimensionality of the $`s`$-wave SC, or the size of the momentum-independent SC gap, $`\mathrm{\Delta }_\mathrm{𝟎}`$. For definiteness we set $`\mu =370`$ meV and $`m^1/\mathrm{\Delta }_0=15`$, but quantitatively similar results are obtained for $`m^1/\mathrm{\Delta }_0=30`$. The DOS, $`N(𝐫,\omega )`$, presented below is obtained from a numerical computation of Eqs.(1)-(3) with $`N(𝐫,\omega )=A_{11}+A_{22}`$ and $`A_{ii}(𝐫,\omega )=\mathrm{Im}\widehat{G}_{ii}(𝐫,\omega +i\delta )/\pi `$.

For a single magnetic impurity in an s-wave SC, the $`\widehat{T}`$-matrix possesses poles at frequencies $`\omega _{res}^{(1,2)}`$, reflecting the presence of two bound states. The spectroscopic evidence for these bound states are two peaks in the DOS, as shown in Fig. 1, where we present the DOS obtained from Eqs.(1)-(3) at the impurity site; for comparison, we also plot the DOS of the clean system. These results are in general agreement with those of STS experiments Yaz97 , which provides further support for the validity of the $`\widehat{T}`$-matrix approach. Assuming for definiteness that the impurity spin $`𝐒\widehat{𝐳}`$ and $`J>0`$, we find that the bound state at $`\omega _{res}^{(1)}<0`$ ($`\omega _{res}^{(2)}>0`$), which we denoted by $`|p,`$ ($`|h,`$), is particle-like (hole-like) with spin along the $`z`$-direction ($`+z`$-direction).

We next consider two magnetic impurities with parallel spins, $`U_i=0`$ and $`J_i=J`$. For $`\mathrm{\Delta }r=\mathrm{}`$, the two sets of bound states given by $`|p,,i`$ and $`|h,,i`$ ($`i=1,2`$) are degenerate. However, for $`\mathrm{\Delta }r<\mathrm{}`$, the probability that an electron scattered by one of the impurities is also scattered by the second one is non-zero. Hence, quantum interference of electronic waves that are scattered by both impurities leads to the formation of novel even and odd (or bonding and anti-bonding) states, $`|p,_{e,o}=(|p,,1\pm |p,,2)/\sqrt{2}`$, and similarly for the hole-like states. This picture is confirmed by the numerically computed DOS shown in Fig. 1 for two impurities located at $`𝐫_1=(0,0)`$ and $`𝐫_2=(2,0)`$, and $`\beta _0=300`$ meV (the DOS shown is that on one of the impurity sites). As expected, the DOS exhibits four mid-gap peaks with peak $`(1),(2)`$ corresponding to the particle-like states $`|p,_{e,o}`$ and peak $`(3),(4)`$ to the hole-like states $`|h,_{e,o}`$.

To determine which peaks in the DOS correspond to the even and odd states, we plot in Fig. 2a the spatial dependence of the particle-like states $`(1)`$ and $`(2)`$ along the $`\widehat{x}`$axis with $`𝐑=(r,0)`$ (the location of the impurities at $`r=0`$ and $`r=2`$ are indicated by arrows).

Since the DOS of the odd states vanishes by symmetry at the midpoint between the two impurities, i.e., at $`r=1`$, peak $`(2)`$ and $`(1)`$ correspond to the odd and even particle-like states, respectively. Note that their spatial dependence is remarkably different: while the odd state exhibits oscillations well beyond the two impurity region, the even state is primarily confined to the region between the two impurities. This qualitative difference is associated with the $`(k_Fr)`$-oscillations of the $`|p,,i`$-states. Since $`k_F=\pi /2`$ and $`\mathrm{\Delta }r=2`$, the wave-functions of $`|p,,1`$ and $`|p,,2`$ are shifted by a phase $`\mathrm{\Delta }\varphi =k_F\mathrm{\Delta }r=\pi `$ outside the two-impurity region ($`r>2`$ in Fig. 2a) and their spatial oscillations are consequently out-of-phase. Thus, $`|p,,1`$ and $`|p,,2`$ interfere destructively for $`|p,_e`$, and only weak spatial oscillations are observable in the DOS for $`r>2`$. In contrast, $`|p,_o`$ shows constructive interference of $`|p,,i`$ ($`i=1,2`$) and its spatial oscillations are enhanced for $`r>2`$. By changing the inter-impurity distance to $`\mathrm{\Delta }r=4`$ , with $`\mathrm{\Delta }\varphi =2\pi `$, the interference pattern between the even and odd states is exchanged, and $`|p,_e`$ ($`|p,_o`$) now exhibits strong (weak) oscillations beyond the two impurity region, $`r>4`$, as shown in Fig. 2b.

In Fig. 3 we present the bound state energies of the states $`|p,_{e,o}`$ as a function of $`\mathrm{\Delta }r`$ for $`\alpha =0`$ (the frequencies of the $`|h,_{e,o}`$-states are obtained via $`\omega _{res}\omega _{res}`$).

The amplitude of the oscillations in $`\omega _{res}`$ decays as $`\mathrm{\Delta }r^{d/2}e^{\mathrm{\Delta }r/\xi }`$ in d-dimensions, but since $`\xi >v_F/\mathrm{\Delta }_020`$, the exponential decay is barely perceivable in Fig. 3. The oscillations’ spatial period, $`1/k_F`$, directly reflects that of the ($`k_Fr`$)-oscillations in the wave-functions of $`|p,_{e,o}`$, in contrast, to the ($`2k_Fr`$)-oscillations of the DOS. For certain $`\mathrm{\Delta }r`$, the splitting between the even and odd states vanishes and the two sets of bound states are decoupled. This decoupling arises whenever the bound state wave-functions of one impurity possess a node at the position of the other impurity. The same type of oscillations were obtained in Ref.Flatte  for the case of NbSe2.

At small $`\mathrm{\Delta }r`$, the bound state energy of $`|p,_e`$ crosses zero, and the state becomes hole-like, $`|h,_e`$ (at the same time, the state $`|h,_e`$ transforms into $`|p,_e`$). We find that this zero-crossing of $`\omega _{res}`$ is accompanied by a crossover in the spin-polarization of the superconducting system which at $`T=0`$ is given by

$$s_z=\frac{1}{2}d^2r_{\mathrm{}}^0𝑑\omega \left[A_{11}(𝐫,\omega )A_{22}(𝐫,\omega )\right].$$
(4)
This crossover is similar to the one predicted to occur when the scattering strength, $`\beta _0`$, of a single magnetic impurity in an $`s`$-wave SC exceeds a critical value, $`\beta _c`$ Balatsky ; Sak70  (for the band parameters chosen, we obtain $`\beta _c460`$ meV). At this point, the impurity breaks a Cooper-pair and forms a bound state with one of its electrons. Specifically, for $`𝐒\widehat{z}`$ and $`J>0`$, the spin polarization changes from $`s_z=0`$ for $`\beta <\beta _c`$, to $`s_z=1/2`$ for $`\beta >\beta _c`$. Similarly, for $`\beta _0>\beta _c/2`$, the system undergoes a crossover at $`\mathrm{\Delta }r_c`$, and for $`\mathrm{\Delta }r<\mathrm{\Delta }r_c`$ one electron of the broken-up Cooper-pair forms a single bound state with both impurities. For $`\mathrm{\Delta }r0`$, the DOS reduces to that of a single magnetic impurity with scattering strength $`2\beta _0`$. Accordingly, $`|p,_o`$ moves towards the particle-hole continuum and vanishes for $`\mathrm{\Delta }r0`$. Finally, a comparison of Fig. 2 and Fig. 3a shows that, as expected, the spatially more confined bound state possesses a larger $`|\omega _{res}|`$ than the spatially more extended one.

As $`\beta _0`$ approaches $`\beta _c`$ from below, the number of crossovers increases, as shown in Fig. 3b for $`\beta _0=400\mathrm{meV}0.87\beta _c`$. Due to the oscillatory behavior of $`\omega _{res}`$, the bound state energy of the $``$-state crosses zero for several value of $`\mathrm{\Delta }r_c`$. As a result, the spin-polarization oscillates between $`s_z=0`$ and $`1/2`$ and the electronic system can be tuned through multiple crossovers by varying $`\mathrm{\Delta }r`$. The same tuning could also be achieved by keeping $`\mathrm{\Delta }r`$ constant and altering $`k_F`$ through changes in the doping level using recently developed field-effect transistor geometries Sch00 .

For $`\beta _0\beta _c`$ the superconductor exhibits a different crossover in which its spin polarization changes from $`s_z=1`$ to $`s_z=1/2`$. For $`\mathrm{\Delta }r=\mathrm{}`$, each impurity breaks one cooper pair and the spin polarization of the superconducting system is $`s_z=1`$. As $`\mathrm{\Delta }r`$ decreases, one of the bound state energies crosses zero at least once, as shown in Fig. 3c where we present $`\omega _{res}`$ for the $``$-states and $`\beta _0=600\mathrm{meV}>\beta _c`$. Accordingly, $`s_z`$ changes from $`1`$ to $`1/2`$. Since for $`\mathrm{\Delta }r0`$, the odd bound state vanishes by symmetry, the spin polarization reaches $`s_z=1/2`$ for any value of $`\beta _0\beta _c`$.

Changes in $`s_z`$ are also reflected in the spatially resolved spin polarization, $`s_z(r)=_\mathrm{\Delta }^0𝑑\omega (A_{11}A_{22})`$, which due to the limited frequency integration is experimentally more easily accessible. In Fig. 3d, we plot $`s_z(r)`$ along $`𝐑=(r,0)`$ for two impurities located at $`r=0`$ and $`r=2`$ and $`s_z=0,1/2,`$ and $`1`$, corresponding to $`\beta _0=300,400`$ and $`800`$ meV, respectively. For $`s_z=0`$, the spin polarization near the impurities is negative, as expected for $`𝐒\widehat{z}`$ and $`J>0`$. For $`s_z=1/2`$, both impurities form a single bound state with an electron in the $`|_o`$-state (see Figs. 3b). Thus, $`s_z(r)`$ is substantially increased at the impurity sites, but remains practically unchanged at $`r=1`$. In contrast, for $`s_z=1`$, the electron from the second broken Cooper-pair joining the two-impurity bound state is in the $`|_e`$-state, and consequently, $`s_z(r)`$ increases primarily around $`r=1`$. Note, that for two impurities separated by $`\mathrm{\Delta }r=4`$, the first electron to join the two-impurity bound state is in the $`|_e`$-state, while the second one is in the $`|_o`$-state, with corresponding changes in $`s_z(r)`$.

The superconducting system can also be tuned through a crossover by changing the angle, $`\alpha `$, between the two impurity spins, as shown in Fig. 4, where we plot $`\omega _{res}`$ for all four bound states as a function of $`\alpha `$ (the impurities are located at $`𝐫_1=(0,0)`$ and $`𝐫_2=(1,0)`$). Since we choose $`\beta _0=400`$ meV $`>\beta _c/2`$ we have $`s_z=1/2`$ for $`\alpha =0`$, corresponding to the vertical dashed line in Fig. 3b. As $`\alpha `$ increases from zero, the frequencies of the even bound states move towards $`\omega =0`$ which they cross zero at $`\alpha 0.27\pi `$. Simultaneously the spin polarization changes from $`s_z=1/2`$ to $`s_z=0`$. The frequency separation between the even and odd states of a given spin direction decreases with increasing $`\alpha `$ and vanishes at $`\alpha =\pi `$. This is expected since for antiparallel impurity spins ($`\alpha =\pi `$), the bound states for impurity 1 ($`|p,,1`$ and $`|h,,1`$) and impurity 2 ($`|p,,2`$ and $`|h,,2`$) possess different quantum numbers; thus they cannot be coupled and remain degenerate. However, since the bound states of one impurity are subjected to the repulsive potential of the second impurity, their resonance frequencies are larger than those of a single impurity with the same $`\beta _0`$ (indicated by the arrows on the right). This repulsion leads to the disappearance of all bound states for $`\mathrm{\Delta }r0`$.

Finally, a non-zero $`U`$ transfers spectral weight between the particle- and hole-like states and increases $`\beta _c`$, but does not affect our above conclusions. Moreover, a self-consistent approach that allows for a gap suppression near the magnetic impurity does not change the qualitative features of the DOS discussed above Balatsky , in agreement with experiment Yaz97 .

The results presented above suggest that a superconducting system with $`N`$ impurities for which $`\beta _0>\beta _c/N`$ can be tuned through multiple crossovers with spin polarizations ranging from $`s_z=0`$ to $`N/2`$, depending on $`\beta _0`$, the inter-impurity distances, and the angles between the spin moments. Work is currently under way to study these crossovers in more complex impurity geometries, such as quantum corrals, as well as the extensions to other host materials, such as unconventional SC, charge-density-wave systems, or semi-conductors Morr02b .

In summary, we show that quantum interference of electronic waves scattered by two magnetic impurities in an $`s`$-wave SC gives rise to novel bound states. We predict that by varying the inter-impurity distance or the angle between the impurity spins, the states’ quantum numbers can be altered, and the SC can be driven through multiple local crossovers in which its spin polarization changes between $`s_z=0,1/2`$ and $`1`$.

We would like to thank A. Balatsky, J.C. Davis, A. de Lozanne, M. Randeria and A. Yazdani for stimulating discussions.


text:
# Vortex Excitations in 𝛼 Cluster Condensed Nuclei

## Abstract

We point out an intriguing possibility that there exist topological soliton-like excitations in nuclei with $`\alpha `$ cluster condensation. The excitations are vortices in superfluid, a property possessed by such Bose condensed nuclei. We show that the vortex-like states possess quantum numbers $`J^P=n^{(1)^n}`$ and energies $``$ several MeV, measured from the ground state energies of the $`n\alpha `$ condensed states. We point out actual candidates, for instance, $`2^+(3\text{MeV})`$ for 8Be, $`3^{}(9.6\text{MeV})`$ for 12C, etc, although the topological nature as a vortex is obscure in such light nuclei.

preprint: Nisho-02/2

Topological excitations arise in many branches of physics and are fascinating objects in nature. Especially, vortex excitations are well known in superfluid such as He liquid, superconductor, or in quantum Hall liquid. These excitations are characterized as topological solitons in Bose condensed matter with infinite degree of freedom. Namely, they are coherent states of bosons, e.g. He, Cooper pairs of electrons or bosonized electrons. Even in a finite system of atomic gas involving the order of $`10^5`$ atoms, such vortex excitations have been observed. Now, an interesting question arises; whether or not vortex state may exist in a smaller system whose number of dynamical degrees of freedom is such as $`10^4,10^3,10^2,10`$ or $`3,4`$. We expect that vortex excitations or vortex-like excitations are still present even in such a small system although the topological properties ( topological stability, vorticity, etc. ) of the excitations become obscure gradually as the number of degrees of freedom in the system decreases. But some specific properties characterizing the vortex solitons are preserved even in such states with small dynamical degrees of freedom. In this letter we point out that vortex-like excitations may arise in nuclei with $`\alpha `$ cluster condensation.

It is well known that many nuclear states involve cluster structures in themselves such as $`\alpha `$, C, O etc. Among them, the states with the $`\alpha `$ cluster condensation have recently been argued to exist around threshold energies at which the nuclei are dissociated into each cluster. The excited states with larger radii than that of the ground state of 12C or 16O have been shown to involve possibly $`\alpha `$ clusters interacting weakly with each other. For example, a detail calculation has suggested that the second $`0^+`$ state in 12C is such a state with $`3\alpha `$ clusters with its rms radius $`4.3`$ fm. Thus, the nuclear density is much smaller than that of the ground state whose radius is about $`2.5`$ fm. The small density of the nucleus implies that the $`\alpha `$ clusters interact weakly with each other so that they may condense into an identical state. Similarly, the calculation have shown possible existence of $`n\alpha `$ clusters condensed nuclei with larger mass numbers $`4n`$ ( $`n3`$ ).

Based on the possible existence of the $`\alpha `$ cluster condensed states in nuclei, we wish to see the vortex-like excitations in the nuclei. Our concern is whether or not their energies are accessible for observation in such small system. Furthermore, we wish to clarify what quantum state do the clusters occupy in such a vortex-like nuclear state of the finite system, although topological excitations can only arise in infinite system, rigorously speaking.

Topological states are characterized by nontrivial phases of the condensates and they are not eigenstates of the particle number. The number of the particles in the states is fluctuating. Nuclei are, however, the eigenstates of the number of nucleons. On the other hand, the number of clusters in some nuclear states is not necessarily definite. The wave functions of the states involve non-cluster components as well as cluster components. That is, the number of the clusters is fluctuating generally in such nuclear states. Therefore, it is reasonable to expect the presence of such nuclear states involving the vortex-like excitation.

Hereafter, we are only concerned with symmetric nuclei ( $`N=Z`$ ). Some excited states of such nuclei are supposed to be $`\alpha `$ condensed states. We use a model of a complex boson field $`\varphi `$ in order to describe condensed states and vortex excitations. The boson field represents quanta of $`\alpha `$ clusters which interacts with each other through a contact interaction in the model. The Hamiltonian is supposed to be given by

$$H=d^3x\left(\frac{|\stackrel{}{}\varphi |^2}{2M}+\frac{\lambda (|\varphi |^2\rho )^2}{2}\right)+V_c$$
(1)
where $`M`$ ( $`\rho `$ ) denotes mass of $`\alpha `$ particle ( average number density of $`\alpha `$ clusters in Bose condensed nuclei ). The term $`|\varphi |^4`$ represents the contact interaction between $`\alpha `$ clusters, i.e. $`V(xy)=\lambda \delta ^3(xy)`$. The strength $`\lambda `$ can be extracted from Gaussian potentials and is given by $`1.7\times 10^4/(\text{MeV})^2`$. The term $`V_c`$ represents Coulomb interaction: $`V_c=e^2(2|\varphi (x)|^2)\frac{1}{2|xy|}(2|\varphi (y)|^2)d^3xd^3y3(2ne)^2/5R_0`$ where the last term represents the Coulomb energy of spherical nuclei with uniform charge distributions. We are only concerned with the difference between the energy of the ground state with $`\alpha `$ condensation and the energies of the excited states with vortices.

We have assumed in the treatment of the $`\alpha `$ particles that the $`\alpha `$ particles interact very weakly in the $`\alpha `$ condensed nuclei and that their internal structures can be neglected; they are weakly interacting point particles. The assumption seems to be reasonable because the proposed $`\alpha `$ condensed nuclei have much smaller nuclear density ( $`1/5\times \text{standard nuclear density}`$) so that the momentum of each cluster is not so large as the internal structure of the cluster to affect seriously the properties of the condensed states. Furthermore, we assume hereafter that the properties of the finite nuclei can be extracted from those of infinite nuclear matter. Actually, the Hamiltonian describes the infinite $`\alpha `$ clusters matter. In order to describe precisely the finite $`\alpha `$ condensed nuclei we should use more realistic potentials involving attractive potentials. In our model only repulsive potential is included and this is sufficient for the formation of the $`\alpha `$ condensation and the vortex in the infinite nuclear matter. Obviously, attractive potentials for the formation of the $`\alpha `$ condensed finite nuclei are needed. Thus, more realistic treatment by which our analysis is now in progress, is to use the Gaussian potentials in . Such a treatment leads to vortices as non-topological solitons.

We comment that Coleman’s Q ball is also non-topological solitons. His model Hamiltonian involves also attractive potentials of fields. Especially, they are point attractive interactions, while ours are nonlocal potentials such as Gaussian potentials. As far as we know, such nonlocal potentials have not been used for the formation of nontopological solitons. But they are inevitably needed for the formation of such solitons in finite nuclei.

The ground state of the Hamiltonian is given by $`<\varphi >\varphi _0=e^{i\theta _0}\sqrt{\rho }`$ with a constant $`\theta _0`$. Obviously, this represents an infinitely large condensed state of the $`\alpha `$ clusters. In order to describe finite nuclei we simply cut the finite matter out of this infinite matter. Thus, the mass number $`A`$ of the nuclei under consideration is given by $`A=4n=d^3x|\varphi _0|^2=4\pi R_0^3\rho /3`$ where $`R_0`$ is the radius of the $`n\alpha `$ condensed nucleus with no vortex excitations. Here, the form of the nucleus is assumed to be spherical.

Vortex excitation is characterized by non trivial phase factor $`e^{i\theta }`$ in $`\varphi `$ with angle $`\theta `$ around $`z`$ axis. Namely, this factor induces a current $`\stackrel{}{}\theta `$ of $`\alpha `$ clusters around the center of the vortex; the center is along $`z`$ axis. In order to find such vortex excitations, by taking $`\varphi =\sqrt{\rho }f(r)e^{i\theta }`$, we solve an equation derived from the Hamiltonian,

$$\frac{1}{2M}\stackrel{}{}^2\varphi +\lambda (|\varphi |^2\rho )\varphi =\rho e^{i\theta }\{\frac{1}{2M}\left(f(r)^{\prime \prime }+\frac{f(r)^{}}{r}\frac{1}{r^2}\right)+\lambda \rho (f(r)^21)f(r)\}=0,$$
(2)
where $`r`$ denotes radial coordinate in x-y plane and the prime denotes a derivative in $`r`$. We have assumed axial symmetric solutions which has no dependence on the coordinate $`z`$. The Coulomb interaction is taken into account perturbatively; after finding the solutions we calculate its Coulomb energy.

In order to avoid a singularity at $`r=0`$, $`f(r)`$ has to satisfy the boundary condition, $`f(r=0)=0`$. $`f(r)`$ should also satisfy the boundary condition, $`f(r=\mathrm{})=1`$. The latter condition is, however, not necessarily needed for the finite nuclei although the condition is needed for the infinite system in order to for the vortex energy per unit length being finite. It has to be replaced by the condition $`f(r=\mathrm{})=0`$ when we take into account a realistic potential of the $`\alpha `$ clusters. But since the spatial extension ( coherent length $`=1/(\sqrt{2\lambda M\rho })0.4`$ fm ) of the vortex is much less than the radius of the $`\alpha `$ condensed nucleus $`4`$ fm, the condition $`f(r=\mathrm{})=1`$ should be satisfied effectively. Here, we adopt this condition for simplicity. ( In real finite system $`f`$ starts with $`0`$ at $`r=0`$ and goes to $`1`$ as $`r`$ goes beyond the coherent length. Finally, it goes to zero as $`r\mathrm{}`$. Thus in our approximation we neglect surface effects of the nuclei although it may be important in small nuclei.) In our treatment the nucleus with the vortex excitation is supposed to be an axial symmetric ball; the density of $`\alpha `$ cluster does not depend on both $`z`$ and $`\theta `$, but depends on $`r`$ in the spherical ball. Although it is very crude approximation, this simple treatment reveals significant properties of the states with vortices.

The vortex solutions of the equation (2) can be obtained numerically. The energies of those solutions are obtained by inserting the solutions into the Hamiltonian,

$`E_v`$ $`=`$ $`{\displaystyle \frac{2\pi \rho }{M\sqrt{M\lambda \rho }}}K,K={\displaystyle _0^{S_v}}𝑑ss\sqrt{S_v^2s^2}\left((f^{})^2+f^2/s^2+(f^21)^2\right)`$ (3)
$``$ $`6.2\text{MeV}\sqrt{{\displaystyle \frac{n}{3}}}\left({\displaystyle \frac{4\text{fm}}{R_0}}\right)^{3/2}(S_v/7.1)(\mathrm{log}(S_v/7.1)+1)`$ (4)
with $`S_v=\sqrt{\lambda M\rho }R_v=3.8\sqrt{n}(R_v/R_0)\sqrt{4\text{fm}/R_0}`$, where we have not included their Coulomb energies in this formula. $`R_v`$ is the radius of the nucleus with the vortex and is defined such that $`n=|\varphi |^2d^3x\frac{0.05}{\sqrt{n}}(R_0/4\text{fm})^{3/2}_0^{S_v}s𝑑s\sqrt{S_v^2s^2}f[s]^2`$. $`E_v`$ represents the energy of the state with the vortex excitation measured from the energy of $`\alpha `$ condensed state $`<\varphi >=\varphi _0`$ without the vortex. We note that since the radius $`R_v`$ ( and $`R_0`$ ) may becomes large such as $`R_vA^{1/3}`$ because of the constancy of the nuclear density even in the Bose condensed nuclei, the energy of the excited state with the vortex goes such as $`E_vKA^{1/3}=(4n)^{1/3}R_v`$. Namely, when the density $`\rho `$ does not depend on $`n`$, $`E_v`$ approaches the energy of the vortex in the infinite system where the vortex energy is proportional to its length $`R_v`$.

We do not still take account of the contribution of the Coulomb energy. But, we find that the difference between the Coulomb energy of the vortex state and the state $`<\varphi >=\varphi _0`$ is much less than $`1`$ MeV. This is because the radius $`R_v`$ of the vortex state is nearly equal to $`R_0`$ and so the Coulomb energies of the both states are almost identical. $`R_v`$ becomes slightly larger than $`R_0`$ owing to the lack of the $`\alpha `$ clusters around the center of the vortex; $`R_v/R_01.05`$

Up to now, we use the terminology of “nuclear state with vortex”. Since the number of the alpha particles in nuclei discussed below is at most 4 or 5, topological nature is very obscure. But there are states possessing some properties of the vortex solitons. Thus, we should use a terminology of “vortex-like excitation”.

12C has been argued to possess the condensed states of 3 $`\alpha `$ clusters. Since its radius is $`4`$ fm, the state with the vortex-like excitation has the energy $`6`$ MeV measured from the energy of the condensed state. Similarly, the state with the vortex in 16O has the energy $`9`$ MeV. We expect in general that the energies of the states with the vortex-like excitations increase such as $`E_v6(n/3)^{1/3}`$ MeV with the number $`n`$ of $`\alpha `$ clusters, since $`E_vR_v(4n)^{1/3}`$. It is important to stress that although our numerical results come from very rough approximations, they indicate the energies of such vortex-like states sufficiently small for experimental detection ( probably, several MeV $`10`$ MeV for 12C, 16O, or 20Ne smaller than the above estimation because attractive interactions should be included in more realistic treatments ).

We should mention that since these topological excitations are absolutely stable in the infinite system, they may be rather stable even in the finite system. Namely, they ( $`A=4n`$ ) would decay mainly into the lighter vortex-like states ( $`A=4(n1)`$ ) by emitting an $`\alpha `$ particle so as to preserve the vorticity, although the vorticity is obscure in the finite system. If this decay mode is inhibited energetically, the other decay mode would be the decay with $`\gamma `$ emission. This corresponds to the disappearance of the vortex due to the finiteness of the system; the vortex gets away from the condensed state without passing infinite energy barrier. This decay mode is relatively rare. Thus the decay width of the vortex-like nuclei would be small when $`\alpha `$ emission is forbidden. This tendency can be seen more prominently in heavier nuclei, since the number of alpha particles increases and the vorticity becomes clearer in such heavier nuclei.

In order to see what quantum states each cluster occupies in the vortex-like nuclei, we calculate the angular momentum of the vortex soliton, especially, the $`z`$ component of the angular momentum. It is easy to find

$$J_z=d^3x(xP_yyP_x)=d^3x|\varphi |^2=n$$
(5)
with no use of the explicit form of the solution, where the momentum $`P_i`$ of the field is given by $`P_i=\frac{1}{2i}(\varphi ^{}_i\varphi _i\varphi ^{}\varphi )=|\varphi |^2_i\theta `$. We also find $`J_x=J_y=0`$ with the use of the assumption of the axial symmetry of the solution. The result implies that all clusters in the finite nucleus take an identical state with angular momentum $`l=1`$ and rotate around $`z`$ axis, i.e. $`l_z=+1`$. We notice the configulation of the vortex field $`\varphi e^{i\theta }`$, which indicates that the states with $`l_z=1`$ compose the vortex-like state. Therefore, we expect that the states with the vortex-like excitation are such states with $`J^P=3^{}`$ for 12C, $`J^P=4^+`$ for 16O, etc. Generally, $`J^P=n^{(1)^n}`$ for the nucleus $`A=4n`$. All of the clusters occupy a p state with $`l_z=1`$ so that their wave functions vanish at the center axis ( $`r=0`$ ) of the nucleus. This lack of clusters at the center axis of the nucleus is consistent with the form of the vortex solution, i.e. $`f(r=0)=0`$. In this sense the states have the similar property to that of the vortex soliton.

Now, we wish to find some candidates for the vortex-like nuclear states. The states possess quantum number $`J^P=n^{(1)^n}`$ with $`J_z=n`$ and energies $``$ several MeV measured from the $`\alpha `$ cluster condensed nuclear states. First, we take 8Be for which 2 alpha clusters play. The state $`2^+(3\text{MeV})`$ is a candidate for the nuclear vortex-like state. The state is well known as a molecular state with large nuclear radius and decays into two alphas. This state is the lightest vortex-like nuclear state, although its character as a vortex is very obscure. Second, the state $`3^{}(9.6\text{MeV})`$ is a candidate for the vortex-like nuclear state of 12C. The state is also well known theoretically as a molecular state of equilateral triangle. Each cluster occupies a p state with $`l_z=1`$. This molecular state is not so rigid that each cluster can move easily from its original position at the corner of the equilateral triangle. The molecular state is soft. Therefore, we may regard the state as a vortex-like state. Furthermore, its main decay mode is the emission of gamma ray with its width being small. This fact also supports our interpretation. Finally, we point out a candidate for the nuclear vortex-like state of 16O. This is the state $`4^+(16\text{MeV})`$, which is not yet fully understood. The state will be understood as a soft molecular state of equilateral square. This is a our prediction based on the idea of the nuclear vortex-like states.

Finally we mention that the vortex state with higher vorticity $`\varphi e^{im\theta }`$ ( $`m>1`$ ) may exist with much higher energy as well as the topological state possessing two vortices with each vorticity of $`e^{i\theta }`$. These are states of $`\alpha `$ clusters occupying higher angular momentum ( $`l>1`$ ). These vortex states would have various decay modes satisfying the conservation of the vorticity. We also point out that our formulation is based on hypothetical $`\alpha `$ cluster condensed nuclei. Although they have been argued to arise around the threshold energy, our results hold even if they arise in any regions of energies. It is difficult, in general, to confirm the existence of such Bose condensed nuclei. But the observation of the topological nuclear states above the condensed states result in the confirmation of the superfluidity of the Bose condensed nuclei.

To summarize, we found the following picture of the $`\alpha `$ cluster condensed finite nuclei with vortex-like excitations. Their excitation energies are several MeV$`10`$ MeV for 12C $``$ 20Na measured from the $`\alpha `$ condensed states without vortices. Since all $`\alpha `$ clusters occupy an p state with $`l_z=1`$ in the vortex-like nuclei, the states are characterized by $`J^P=n^{(1)^n}`$ for the nuclei with the mass number $`A=4n`$. Their main decay modes are expected such that a vortex-like state with $`J^P=n^{(1)^n}`$ decays predominantly into a vortex-like state with $`J^P=(n1)^{(1)^{n1}}`$ with emission of an $`\alpha `$ particle. But the decay width is small because the $`\alpha `$ particle emitted is dominantly in a s state.

We thank Y. Akaishi, O. Morimatsu and Y. Kanada-En’yo in KEK for fruitful discussions and H. Nemura for useful information of $`\alpha `$ particles. We also thank N. Itagaki for useful imformations of nuclear molecular states in 12C and 16O.


text:
# Electron-electron interactions in a weakly screened two-dimensional electron system

## ACKNOWLEDGMENTS

The authors wish to acknowledge M.I. Dykman, D. Herman, and H. Mathur for helpful conversations. This work was supported in part by NSF grant DMR-0071622.