text:
# An Efficient Search for Gravitationally-Lensed Radio Lobes

## 1 Introduction

The 60 known gravitational lenses are powerful astrophysical probes. They are the only galaxies selected on mass rather than luminosity, making them a unique sample with very different properties and biases from traditional luminosity-selected galaxies. All lenses have extraordinarily precise mass measurements (better than 5% including systematics), allowing detailed comparisons of galaxies of constant mass. Lensing has been used to constrain cosmological parameters (see reviews by Schneideret al. 1992; Blandford & Narayan 1992), the structure and evolution of galaxies (e.g., Keetonet al. 1997, Keetonet al. 1998, Kochaneket al. 2000), quasar and AGN host galaxies (e.g., Keetonet al. 2000), and even the interstellar medium and extinction in intermediate-redshift galaxies (e.g., Falcoet al. 1999).

Such astrophysical studies are limited by the sample size. The properties of each individual lens depend on many variables (the source redshift, the lens redshift, and the galaxy mass). With only 60 lenses, it is difficult to subdivide the sample over more than one variable. For example, one cannot separately explore the evolution of high and low mass galaxies, or the average mass or number density of the lenses in each bin. Only 10–20% of lenses are produced by spiral galaxies (e.g. Fukugita & Turner 1991), so hundreds of lensed systems must be found to provide a useful sample of spiral galaxy lenses.

Fortunately, there is no shortage of undiscovered lenses. Statistical estimates suggest that there are $`0.5`$ lenses per square degree for a quasar sample complete to about $`V`$=24 (Kochanek 1996a), or for radio sources with 5 GHz radio flux of $`>`$10 mJy (Kochanek 1996b). Equivalently, 1–10 in 1000 sources with redshifts above unity will be gravitationally lensed depending on the type of source and its flux. The problem is finding them in an efficient manner.

The first surveys for lenses focused on bright quasars (see Blandford & Narayan 1992), with a discovery rate of roughly 1 per 100–200 targets. The high rate is due to the effects of magnification bias (Turner 1990), which is important for all lens samples but critical for surveys of bright quasars. This strategy is limited by the relative rarity of bright quasars, the difficulty of obtaining large numbers of high resolution optical images, and contamination from stars.

Another successful strategy was to look for lensed radio sources. Extra-galactic radio sources are numerous and distant, and confusing low-redshift sources are rare. The three largest radio searches, MG-VLA (Burkeet al. 1993), JVAS (Patnaiket al. 1992; Patnaik 1993) and CLASS (Jacksonet al. 1997), have discovered almost half of the $``$60 known galaxy-mass lenses. Their strategy was to image all radio sources to a fixed flux limit, and then use the morphology of each radio image to identify the lenses. Typically, one lensed radio source is discovered for every 500 survey targets. Radio searches also led to the discovery of an important new class of lenses where the lensed images of an extended background source merge to form an “Einstein ring” of emission around the lens galaxy (e.g., Hewittet al. 1988). Unlike a compact source, an extended background source probes many lines of sight which extend well into the halo of the lens, allowing much better constraints on the dark-matter profile (e.g., Kochanek 1995).

In this paper, we present a much more efficient lens search, which uses combined radio and optical surveys to find lensed radio lobes. The survey selection is outlined in §2, and follow-up observations are described in §3. In §5, we discuss individual candidates for lensing, and in §6, we summarize our conclusions.

## 2 The FIRST-APM Lensed Lobe Search

Most lensed radio sources are being missed by present searches. Powerful radio galaxies are often flanked by two extended radio lobes, whose radio spectrum is markedly steeper than that of the compact core (e.g. FR-type objects; Fanaroff & Riley 1974). Due to their extended nature, the optical depth to lensing of such lobes is significantly higher than for unresolved sources. Indeed, the typical lobe size of a few arcseconds implies that $``$70% of lensed radio sources should be lobes (Kochanek & Lawrence 1990). The JVAS and CLASS surveys select only flat-spectrum radio sources, and thus lobes constitute only $``$20% of their lenses. The smaller MG-VLA search found more lensed lobes ($``$50% of lenses), but still fell short of the $``$70% fraction. Many lensed lobes are simply not recognized, because the background source itself often has complicated structure, and often lacks sharp features which can produce a clear lensing signature. Sensitive, high-resolution radio maps are required to identify lensed structures ($``$1″ scale) in these cases.

Lensed lobes can be found efficiently because they usually lie outside their host galaxy. Radio-galaxy hosts are centered on their compact cores, so radio lobes should not have direct optical counterparts. A galaxy which coincides with a radio lobe is almost certainly a chance alignment, and the relative distances of radio and optical sources strongly favors the galaxy being a foreground lens (Kochanek & Lawrence 1990). Clear examples of this phenomenon were found in the MG-VLA search (MG J1654+1346, Langstonet al. 1989; MG J1549+3047, Leháret al. 1993). Thus, an efficient way to find lenses is to identify radio lobes with superposed optical galaxies.

A useful sample of radio lobes requires a large, low frequency, radio survey with sufficient angular resolution to distinguish lobes from cores. The VLA111 The Very Large Array (VLA) is an NRAO facility. The NRAO is operated by Associated Universities, Inc., under cooperative agreement with the National Science Foundation.  FIRST Survey222see http://sundog.stsci.edu/ (Becker, White, & Helfand 1995), has a flux sensitivity of 1 mJy ($`6\sigma `$) at 20 cm, $`5\mathrm{}`$ resolution (FWHM) and sub-arcsecond positional accuracy. This represents a factor of $``$50 improvement in both sensitivity and resolution over previous large-area surveys. At the time of our initial lens search, the northern FIRST sample covered $``$3000 square degrees (7h$`<`$RA$`<`$18h and 22$`<`$DEC$`<`$58), and comprised $``$3$`\times `$105 radio components.

A good sample of optical counterparts is provided by the APM333 The Automatic Plate Measuring (APM) machine is a National Astronomy Facility run by the Institute of Astronomy in Cambridge, UK (see http://www.ast.cam.ac.uk/$``$apmcat).  catalog (McMahon & Irwin 1992; Maddoxet al. 1990), based on the Palomar Observatory Sky Survey plates (POSS-I). The POSS-I plates have a detection limit of $`R`$$``$20, and we conservatively estimate that 25–35% of the known lens galaxies (see Keetonet al. 1998) are bright enough to be included in the APM catalog.

Our initial search for lenses, of course, does not fully exploit the potential of this method. The shallow magnitude limit of the APM survey means that our lens sample is no longer mass-selected, given that only the optically brightest third will be included. Moreover, most of the FIRST lobes are much too faint to be followed up at higher resolution, due to their steep spectra and extended structures. Nevertheless, we have searched the available FIRST survey down to its detection limit, to get a sense for how many lensed lobes can be found this way.

### 2.1 Survey Selection

For the FIRST-APM lensed lobe search, we combined these two surveys to find optical counterparts to radio lobes. We classified the FIRST sources by their radio morphology (see Appendix), and identified those components most likely to be lobes. We then identified as lens candidates those radio lobes with APM counterparts. The estimated mean redshift of FIRST sources is $`z`$$``$1 (Cress & Kamionkowski 1998), and the measured redshift distribution of optical surveys such as the APM Galaxy Survey peaks at $`z`$$``$0.12 (Baugh & Efstathiou 1993). Thus it is likely that APM counterparts to FIRST radio lobes are indeed in the foreground, and potentially lenses.

To define a radio-lobe sample, it is necessary to first “collapse” the radio catalog, such that multiple components of a single radio source are associated with each other. The collapse algorithm is described fully in the Appendix. We ran the collapse routine using a $`60\mathrm{}`$ search radius, a maximum component separation of $`90\mathrm{}`$ for triples and multiples (defined as sources with more than three components), and a 10% random association probability cutoff, $`P<0.1`$. For ease of detection and phase calibration in follow-up radio observations, we required all component fluxes to exceed 3 mJy, and at least one component to exceed 10 mJy. To enhance selection of true, FR-type double-lobed objects, we required a symmetry factor of $`S>0.5`$ for triples (where $`S`$ is defined as the ratio of the shorter core-to-lobe distance to the longer), and a configuration angle of $`\theta <30^{}`$ for triples and multiples (see Figure 10 in Appendix). While these requirements will necessarily eliminate some true FR-type radio galaxies, this selection favors easily-identified radio galaxies with the most typical component configuration. The output of the routine is a list, containing positional, flux, and morphological data, of collapsed radio sources containing two, three, or more components, each of which is identified as either a “core” or a “lobe.” Due to the wide array of radio-source morphologies (core-jet, double-lobed, cometary, diffuse, etc.) and the limited resolution of FIRST, there is no simple, obvious manner in which one can correctly and unambiguously collapse all sources. The parameters we have chosen yield a reasonable balance between ensuring an accurate correspondence between collapsed and real sources, and are in good agreement with the conclusions of “by-eye” inspections (D. Helfand, private communication).

The APM-counterpart selection was as follows. For every collapsed FIRST source, we searched for optical counterparts to each lobe component, as well as to the core component (or to the source centroid, if no radio core was identified). Optical objects were taken to be coincident with a radio component if they fell inside an ellipse defined by the deconvolved FWHM major and minor axes of the radio component, as listed in the FIRST catalog. For cases where the catalogued major and/or minor axis was less than $`2\mathrm{}`$, we adopted an effective FWHM value of $`2\mathrm{}`$ for that axis. For sources without a radio core, optical objects that fell within $`\mathrm{}/4`$ from the average component position were considered core counterparts (where $`\mathrm{}`$ is the maximum component separation). Those sources with at least one optical counterpart to a radio lobe, which was not also a counterpart to the radio core, became candidates for lensing. We eliminated candidates in which a single optical counterpart matched to both lobe components of a double radio source, as these usually corresponded to identifications of the host galaxy; for the same reason, we also eliminated cases where a double or multiple radio source was found to have only one optical match, corresponding to a compact FIRST component. To ensure sufficient contrast for lensing signatures, we eliminated all lensed radio lobes larger than $`10\mathrm{}`$ (FWHM), and further eliminated all optical lens candidates whose lobe counterparts fell farther than $`2^{\prime \prime }`$ from the fitted peak of the radio lobe. Each lobe optical counterpart was assigned a letter classification (Table 1) to indicate its likelihood of being a lens (lenses are likely to be red and extended), based on its APM color and morphology.

There were 1037 candidates in our sample, each of which has an APM object within 2′′ of a presumed FIRST radio lobe’s peak. The selection methods applied to the radio lobe and optical counterpart samples are summarized in Table 2. With about one lens candidate for every $``$200 independent FIRST sources, our candidate selection rate using our chosen criteria is roughly twice the discovery rate for earlier radio lens searches, which suggests that many of our candidates are not lensed. Only about 14% of the candidates have optical core detections, compared to the much higher 40% POSS-I identification rate for MG-VLA lobed radio sources (Lehár 1991), partly because of the higher MG-VLA flux cutoff (60 mJy at 6 cm), but also because many lensed-lobe candidates are themselves misidentified radio cores. The lensed-lobe candidate and FIRST optical counterpart distributions are similar, and both have more red, extended objects than the APM catalog in general (see Table 1), probably because both radio samples preferentially select galaxies. There are a large number of blue counterparts to presumed radio lobes, and probably many of these are quasar counterparts to radio cores. However, 62% of these blue counterparts have uncertain optical identifications, because they are close to the blue POSS-I detection limit and are absent in the red catalog.

### 2.2 The VLA Sample

Because of the limited angular resolution of FIRST ($``$5′′ FWHM), higher-resolution radio observations are needed to confirm candidates as lensed. Most radio lobes extend over only a few arcseconds, so FIRST cannot always determine unambiguously which radio components associated with a source are truly radio lobes, and which are compact cores. Moreover, the typical image separation, $`\mathrm{\Delta }\mathrm{\Theta }`$, for galaxy-mass lenses is only $``$1′′, so the effects of lensing cannot be seen directly in the FIRST images. The VLA can achieve sufficient angular resolution to distinguish lobes and resolve gravitational lensing, but only by observing at shorter wavelengths (A-array resolution $``$0$`\stackrel{}{\mathrm{.}}`$4 at 6 cm).

Since radio lobes generally have very steep radio spectra, the reduced source flux at shorter wavelengths drastically reduces the number of lens candidates accessible to short, high-resolution VLA observations. The lensed radio lobe in MG J1549+3047 (Leháret al. 1993) is typical of what we expect to find in our lensed lobe search, in that only the fainter parts of the lobe are multiply imaged by an off-center lens galaxy. The lensed lobe has a FIRST flux density of 877 mJy, but at 6 cm the total flux density is 200 mJy, with the hot-spots peaking at only 16 mJy per beam. The mean surface brightness of the lensed ring is much fainter (typically $``$2 mJy/beam), but can be detected with a signal-to-noise ratio of 3 in a 1.5 minute VLA observation. If our lensed lobes have similar properties, lensing signatures should be detectable in 30 minute observations for FIRST lobes exceeding 200 mJy. A lensed lobe like that in MG J1549+3047, but with a FIRST flux of 100 mJy, would require a 2 hr VLA 6 cm observation for reliable ring detection. Much shorter observations, however, would suffice to rule out cases where we had misidentified a compact core as a lobe.

With these considerations, we selected a small sample of lensed-lobe candidates for VLA confirmation. There were 78 lensed-lobe candidates with FIRST flux $`>100`$ mJy. The limited VLA schedule forced us to choose a subsample of these targets, and this was done by visual inspection. By comparing the FIRST maps to POSS-I fields from the Digitized Sky Survey (DSS444 Based on photographic data of the National Geographic Society – Palomar Geographic Society to the California Institute of Technology. The plates were processed into the present compressed digital form with their permission. The Digitized Sky Survey was produced at the Space Telescope Science Institute under US Government grant NAG W-2166.), two of us (AB and JL) visually ranked each of the 78 candidates by likelihood of lensing, based on optical colors, radio morphology, alignment, and so on. We compared notes on our rankings and eliminated those candidates that seemed the least likely to be lensed (e.g., candidates which had a single APM counterpart corresponding to a fairly compact component between two others).

Our final VLA sample of 33 candidates is listed in Table 3, and, as a visual guide, Figure 1 presents FIRST cutout maps, with optical fields from the DSS. Since the DSS astrometry can be in error by a few arcseconds, we registered the DSS images to the APM catalog, whose absolute positions deviate by $`<`$0$`\stackrel{}{\mathrm{.}}5`$ for small objects. This was done by selecting $``$20 APM objects in each DSS field and finding the least-squares linear coordinate transformation on the DSS images (internal residuals were typically less than $`0\stackrel{}{\mathrm{.}}5`$). The FIRST maps were superposed using the VLA phase calibrator astrometry which routinely yields $``$0$`\stackrel{}{\mathrm{.}}2`$ positional accuracy.

## 3 Observations

Most of the FIRST/APM lensed-lobe candidates are expected to be radio cores which had been misidentified as lobes, or genuine radio lobes with stellar or nonexistent counterparts (APM/POSS artifacts). To eliminate such cases, and to establish remaining gravitational lens candidates, we obtained higher quality radio and optical observations. We used archival data where available, but the bulk of our targets required new observations. The observations we used are summarized in Table 4, and will be referred to by the labels therein.

We obtained VLA observations for the accessible sample (Table 3). Nine of our targets had appropriate archival observations, and we observed the other 24 using the A-array at 6 cm (r98b, see Table 4). The total allocation of 8 hours was divided between the targets based on the lobe flux (see Table 3). Calibrators were observed every $``$15 min to ensure proper phase calibration, and the flux densities were scaled to 3C286 (with an assumed flux density of 7.43 Jy). Calibration and mapping followed standard AIPS555 AIPS (Astronomical Image Processing System) is distributed by the National Radio Astronomy Observatory, which is a facility of the National Science Foundation operated under cooperative agreement by Associated Universities, Inc. procedures, with a single iteration of self-calibration, and a circular convolving beam of $`0\stackrel{}{\mathrm{.}}7`$ (FWHM), to bring out extended structures. The resulting maps have off-source rms noise levels of $``$0.1 mJy/beam, close to the thermal noise limit, except for the brightest targets. As expected, most of our radio lobe candidates turned out to be misidentified radio cores, and sometimes presumed lobe components could be resolved into separate radio sources.

Three of the most promising VLA targets were observed at 18 cm using MERLIN (r98a, see Table 4). MERLIN provides sub-arcsecond resolution at lower frequencies, where the steep-spectrum radio lobes are brightest, but it is not very sensitive to faint sources. J0823+3906, J0828+3552, and J1622+3531 were each observed eight times in 30 minute scans at intervals of 2–6 hours, to maximize baseline sampling. Standard MERLIN phase calibrators were also observed, and the flux densities and polarizations were calibrated using 3C286 and OQ208, respectively. The data were reduced using the standard MERLIN pipeline, with several iterations of phase self-calibration. The resulting maps have off-source rms noise levels of $``$0.2 mJy/beam, and an angular resolution of $``$0$`\stackrel{}{\mathrm{.}}`$25 (FWHM). MERLIN produced an excellent map of J1622+3531, but in most cases, faint radio lobe structures could be seen.

We obtained optical CCD images for the 12 most promising targets. The objective was to confirm whether candidate lenses are galaxies, and to seek optical counterparts to the radio cores. Observations were acquired when the opportunity arose, so a variety of telescopes and filters were used (see Table 4). The data quality was variable, but the integration times were chosen to permit detection of an $`R`$$``$24 mag galaxy. Photometric calibration was achieved, in most cases, by acquiring Landolt (1983) standard star observations during each run, but in a few cases we were forced calibrate either to spectra or to the APM photometry. In all cases, photometric accuracy is no better than 0.1 mag on an absolute scale, although relative photometry in each image will be considerably better. The images were reduced using standard IRAF666 IRAF is distributed by the National Optical Astronomy Observatories, which are operated by the Association of Universities for Research in Astronomy, Inc., under cooperative agreement with the National Science Foundation. procedures. Over half of the observed candidates could be eliminated as either stars or spurious APM objects, and radio cores were detected for all of the remaining lens candidates (see Figure 2). Two of our targets (J0823+3906 and J1605+5439) were observed through several filters, though time considerations did not permit this for most.

Four of our targets (J0823+3906, J0849+2448, J0958+2947, J1622+3531) were also observed spectroscopically, again using a variety of telescopes as the opportunity arose (see Table 4). The spectra were reduced using standard procedures, and were calibrated by observing flux density standards. All but one of the observations yielded unambiguous redshifts for the targets.

## 4 Individual Targets

Here, we discuss in detail the twelve candidates which could not be rejected immediately as misidentified radio cores and received CCD observations (see Table 3 and Figure 2). We will refer to observations using the labels in Table 4, and optical photometry of the counterparts will be compared to passively evolving early- and late-type galaxies at various redshifts (for details, see §5 of Leháret al. 2000).

J1549+3047 is the known lens MG J1549+3047 (Leháret al. 1993). This target was selected by our search method, and serves as a clear example of a successful detection. The APM finds a bright counterpart to the western FIRST component. Here, we present a 45 s excerpt of the 6 cm VLA data from Leháret al. (1993), which gives a lensed ring detection which is comparable to what our other observations should yield. The CCD image clearly detects both the lens galaxy and the core counterpart, and the measured lens redshift ($`z`$=0.111, Leháret al. 1993) is much less than the typical radio source redshift. Moreover, corresponding parts of the lensed ring have matching spectral index and fractional polarization (Leháret al. 1993).

J0816+5003 is lensed. There are three optical counterparts, one on each of the FIRST components, including the central core component (Figure 2). The $`R`$$`<`$16 southern lobe counterpart is stellar, and is too far away from the radio lobe to be of interest to us. However, the northern lobe counterpart, Galaxy 1, lies directly over the northern radio component, which is broken into three VLA components (N1–N3, Figure 3). The unusual structure of the Northern FIRST component is confirmed in archival multifrequency VLA observations (r96: 2 min at 20 cm, 5 min at 6 cm, 8 min at 3.6 cm). We present a lens model in which N1 and N2 are lensed images of a northern lobe hot-spot, and N3 is the northern lobe’s unlensed terminal ridge. The N1–N2 separation of 1$`\stackrel{}{\mathrm{.}}`$9 corresponds to an isothermal $`z`$$``$0.4 lens galaxy with a velocity dispersion of $``$250$`\mathrm{km}\mathrm{s}^1`$. The lens model has been matched in position and shape to Galaxy 1 (deconvolved FWHM of $``$1$`\stackrel{}{\mathrm{.}}`$8, axial ratio $``$0.9, and PA=–64). An alternative view is that N1–N3 compose a separate radio source, with its core at N2 corresponding to Galaxy 1. While we cannot completely exclude this possibility, we favor a lens interpretation for the following reasons: (1) The clear lobe-like appearance of the southern FIRST component leads one to expect a northern radio lobe near Galaxy 1. (2) Component N1 has a linear morphology which is very unusual for a radio lobe, yet fully expected in our lens model. (3) N2, which would be the core of a separate FR-type source, does not coincide with Galaxy 1’s center (N2 is $`0\stackrel{}{\mathrm{.}}4`$ to the North), when the radio and optical images are aligned using the core components. N2 could be an off-set radio jet, but such jets are less common than coincident cores. Furthermore, most radio jets appear on the same side as the brighter and more polarized radio lobe (see Saikia & Salter 1988), unlike the situation for N1–N3. This arrangement, however, is naturally explained by our lens model. (4) N1 and N2 have identical radio spectra which are more similar to the other lobe hot-spots than to the flatter radio Core’s (Figure 4). Radio lobe hot-spots usually have steeper spectra than radio cores, and their spectra often differ within a single radio galaxy. Similar spectra are expected, however, for gravitational lensing, which is achromatic. (5) N1 and N2 have the same polarized emission properties, both in fraction and orientation, which are very different from N3. Normally, there is no correlation between the core and lobe polarizations (see Saikia & Salter 1988), but lensed image polarizations should match. (6) The radio core counterpart is much fainter ($`R`$=21.7) than Galaxy 1 ($`R`$=19.2), as expected if it is more distant, and the optical magnitude of Galaxy 1 is consistent with an $`L_{}`$ early-type galaxy at intermediate redshift (0.3$`<`$$`z`$$`<`$0.5). So, the properties of N1–N3 are very unusual for an FR-type radio source, but typical for a gravitational lens, with a morphology similar to FSC 10214+4724 (Broadhurst & Lehár 1995; Eisenhardtet al. 1996). Optical spectroscopy is planned to confirm this system as lensed, and high-resolution radio observations are scheduled, to seek more arc-like lensed images, and to determine whether N1 is composed of three merged images.

J0823+3906 is a strong candidate, with a galaxy directly on a radio lobe, but there are no clear lensing signatures. The lobe counterpart (Galaxy 1) is clearly an extended galaxy, with a deconvolved size of $``$1$`\stackrel{}{\mathrm{.}}`$2 (FWHM). If we align the optical image so that the optical core coincides exactly with the radio core, the section of the lobe under Galaxy 1 is too faint to allow multiple imaging to be detected in either the VLA or MERLIN maps. Deeper VLA observations (14 min at 20 cm, and 1.7 hr at 6 cm, Figure 5) detect the flat-spectrum radio core ($`\alpha >+0.5`$), and find substantial emission near Galaxy 1, but still show no lensing signatures in the steep-spectrum lobe ($`\alpha `$$``$–1). Lensing would not occur, of course, if the radio source were closer to us than Galaxy 1. Palomar 5 m spectroscopy (Figure 6) could not detect the core, and there were no clear features in the spectrum of Galaxy 1. However, several dips in the continuum suggest $`z`$$``$0.3 as a possibility. V and I band images of J0823+3906 (from o00b), yield $`I`$=19.2 and $`V`$$`I`$=1.6 for Galaxy 1. This is consistent with an $`L<L_{}`$ early-type galaxy at $`z`$$``$0.3 or an $`L`$$``$$`L_{}`$ spiral at $`z`$$``$0.5 . Either way, Galaxy 1 is not likely to be farther than $`z`$$``$0.6 . The distance to the radio source is harder to determine, however, although it is much fainter ($`I`$=21.8) than Galaxy 1. To eliminate lensing for J0823+3906, deeper optical spectroscopy is crucial.

J0849+2448 is not a gravitational lens. There are optical objects close to both radio components, but their $`12\stackrel{}{\mathrm{.}}7`$ angular separation is greater than the radio component separation, and neither of them falls directly on a radio component. Moreover, both are optically compact, and optical spectra (s98) show them to be local M (North) and G (South) stars.

J0958+2947 is not lensed, but includes two quasars at very different redshifts. This source was first noted in the MG-VLA search (MG0959+3047, Lehár 1991). The VLA map resolves the FIRST components into two separate quasars, each with cores and radio lobes, and there are compact, blue optical counterparts to each core. The optical spectra of the cores (Figure 8) show emission lines corresponding to quasars at high redshift ($`z_{north}`$=2.064 and $`z_{south}`$=2.744). J0958+2947 does not show any clear signs of lensing, but the 20 cm VLA map (Figure 7) and the 6 cm map both show a small arc-like extension, N5, about 2$`\mathrm{}`$ northeast of N1. This could either be an unusual structure connected with this quasar, or it could be part of the of the northern lobe associated with S1, which is otherwise not seen. The spectral index ($`\alpha `$$`<`$–0.5) and polarization ($``$5%) of N5 are typical for a lobe hot-spot. However, the high redshift of the Northern QSO makes it an inefficient lens for a source at $`z`$=2.744, so $`\mathrm{\Delta }\mathrm{\Theta }`$$`<`$0$`\stackrel{}{\mathrm{.}}`$4 for an $`L_{}`$ galaxy. Thus, even if N5 is part of the southern radio source, we would not expect to see multiple imaging.

J1320+4639, J1342+3507 (4C+35.30), J1357+4807, and J1425+3027 (MG J1425+3027) are not gravitationally lensed. In each case, deeper optical observations failed to detect any lobe counterpart, so these are likely to be spurious APM detections.

J1605+5439 is not lensed. The southern lobe counterpart has a stellar profile, and there is no evidence of lensing in the radio lobe. Moreover, the photometry (from o00b) of the lobe counterpart ($`I`$=16.8, $`V`$$`I`$=3.2) is inconsistent with a normal galaxy at any redshift. The core counterpart ($`I`$=20.3, $`V`$$`I`$=1.6), however, is consistent with a galaxy at $`z>0.3`$.

J1622+3531 has nearby galaxy in front of a distant lobe, but our observations failed to detect lensing signatures. The VLA and MERLIN (Figure 9) maps show the northern lobe to be broken into two steep-spectrum ($`\alpha `$$``$–1) subcomponents, N1 and N2, and there are optical counterparts to both the northern lobe and the inverted spectrum radio core ($`\alpha `$=+0.5$`\pm `$0.1). Galaxy 1, at $`z`$=0.32, is extended (deconvolved FWHM $``$0$`\stackrel{}{\mathrm{.}}`$9) and at lower redshift than the Quasar, at $`z`$=1.47 (Figure 10). We can roughly reproduce the configuration of N1 and N2, using a simple lens model (Figure 9), and their spectral indices are similar ($`\alpha _{N1}`$$``$–0.4, $`\alpha _{N2}`$$``$–0.6), as expected for lensing. However, the fractional polarizations (Figure 9) of N1 and N2 are inconsistent with lensing. Moreover, the lens must lie between N1 and N2, where no optical counterpart is found. We attempted to model the radio structure using models with two masses at the positions of Galaxy 1 and Galaxy 2, but these invariably produced unobserved bright counter-images near Galaxy 1. Thus, we conclude that the morphology of the northern lobe is probably due to a bent radio jet (e.g. Lonsdale & Barthel 1986). Galaxy 1 must be lensing the more distant radio lobe, but it is probably just distorting the intrinsic shape. At $`z`$=0.32, Galaxy 1 has an expected image separation of $``$0$`\stackrel{}{\mathrm{.}}`$7, so we might find a faint lensed image of the northern hot-spot near Galaxy 1 in deeper radio observations.

J1641+5115 is not a gravitational lens. There are three APM counterparts, one each for the two FIRST components, and one between them. The VLA map shows the western radio component to be the core, and the eastern counterpart falls 2′′ away from the eastern radio source.

## 5 Conclusions

Our initial search for lensed radio lobes using the FIRST survey and the APM catalog has produced two lenses: J1549+3047 is the known lens MG J1549+3047 (Leháret al. 1993); and J0816+5003 is a new lensed radio lobe system that requires only spectroscopic confirmation. Two, more tentative, lens candidates (J0823+3906, J1622+3531), have galaxies superposed on radio lobes, but show no evidence of multiple imaging. One of our targets (J0958+2947) is a pair of quasars in close projection.

Our search method is highly efficient, even assuming that all of our tentative candidates are not lensed. The FIRST/APM search found 2 lensed systems, using only 33 VLA observations. Thus 6% of our VLA targets are lensed, compared to a success rate of only 0.2% for other radio lens searches, and $`<`$1% for bright quasar searches. While these efficiency estimates do not take into account the effort spent on follow-up optical and radio observations, it is clear that combining optical and radio survey information can drastically reduce the required number of targets for detailed observations.

There are many more lensed radio lobes to be found in this way. Since our initial search, the FIRST survey has trebled in size, covering $``$80% of the planned $``$10,000 square degrees. Thus, we should discover 5–10 lensed lobes using the FIRST/APM comparison in the very near future.

The sample can be expanded and improved with fainter optical catalogs. The POSS-I plates detect only $``$30% of known lens galaxies, so the resulting lens sample is far from mass-selected. However, almost all the known gravitational lenses are brighter than $`R`$=23, and should be accessible to optical surveys (see Keetonet al. 1998). We are currently extending our search using the POSS-II survey (Djorgovskiet al. 1998), which should find at least $``$50% more lens galaxies than the APM, leading to $``$50% completeness, and 7–15 expected lensed lobes. Moreover, two thirds of the APM lobe counterparts turned out to be stellar or spurious in CCD images. This suggests that selection based on a higher quality optical survey, such as the Sloan Digital Sky Survey (Gunn & Knapp 1992), could treble the search efficiency.

Even greater gains are possible as radio telescopes improve. We had to reject 92% of our FIRST/APM candidates because the FIRST components were too faint for easy VLA mapping. Thus, as the VLA upgrades to larger bandwidth, and once higher resolution 20 cm observations are possible using the VLA-Pie Town link or the extended A+ array, from scores to hundreds of lensed lobes could be found in the Northern sky.

## 6 Appendix: Collapsing the FIRST catalog

Many radio sources exhibit complex morphology (e.g., double-lobed, FR-type objects), and often include several components in the FIRST catalog. Thus, it is useful to first “collapse” the FIRST catalog into component groups for each source. The collapsing procedure we employ here is based on the circular overdensity algorithm of Buchalter (1999), which has been used to systematically identify and characterize large numbers of multi-component radio sources in the FIRST Survey (e.g. Blantonet al. 1997; Buchalteret al. 1998).

First, we morphologically classify individual FIRST components and collect them into associated groupings. Components classified as side-lobes are removed from the radio catalog and remaining components are classified as “extended” or “compact”, depending on whether the major axis of the fitted Gaussian profile exceeds 2″. For compact (extended) components, the radio flux is taken to be the listed peak (integrated) flux (see Becker, White, & Helfand, 1995). Starting at the coordinates of the first catalogued component in a source, the program looks for other components within a fixed search radius, $`R`$. If another component is found, the average position is determined, and the search continues from the new average position. The process is repeated until no additional components are found within a distance $`R`$ of the last search center.

For each component grouping, an estimated probability of random association is calculated to determine whether the grouping qualifies as a single radio source. Denoting the number of sources found in a given iteration by $`N`$, the largest separation between the original source and any other by $`d_1`$, and the flux of the faintest source in a given iteration by $`F_{\mathrm{min}}`$, the program calculates the probability of obtaining $`N1`$ sources all with fluxes $`F_{\mathrm{min}}`$ and within a distance $`d_1`$ of the original source, based on counting statistics and on the measured two-point correlation function of the FIRST Survey777 The two-point correlation function, $`\xi _2`$, gives only the excess probability of finding pairs separated by a given distance. Inasmuch as the observed clustering is weak, i.e., $`\xi _21`$, we may conservatively take $`\xi _2`$ to be an upper limit for the excess probability of finding configurations of triples, quadruples, and so on, since higher-order $`n`$-point correlation functions, $`\xi _n`$, scale roughly as $`\xi _n\xi _{2}^{}{}_{}{}^{n1}`$ in the quasi-linear regime (Fry 1984). According to Cresset al.(1996) the measured two-point correlation function for FIRST, once the fragmentation of multicomponent sources (such as FR-type objects) is accounted for, extrapolates to a value of unity only at a scale of about 20$`\mathrm{}`$, so that using $`\xi _2`$ as an upper limit to the higher-order $`n`$-point correlation functions should be valid at all but the smallest scales we consider. Even if two components are separated by only the survey resolution limit of $``$5$`\mathrm{}`$, and the two-point function exhibited no turnover at such small scales, we would only find $`\xi _2`$$``$4, thus perhaps slightly overestimating the number of collapsed sources at the very smallest scales.  (Cresset al. 1996). If the observed configuration has a small probability, $`P`$, of random association, compared to a chosen significance cutoff, the components are considered to be associated with a single source. Once a component is identified with a source, it is no longer considered by the collapse routine in any subsequent searches. If, however, the random association probability exceeds the cutoff, such that no association is found, the program initiates a new search centered at the position of the next catalogued FIRST component. The individual components in the “failed” search attempt are still considered as candidates for collapse in later search iterations.

Each source is then characterized by its component properties and a set of parameters which can vary with the number of components $`N`$. For all sources, the average position and aggregate radio flux of all components are tabulated, as well as the flux of the brightest component. For pairs ($`N`$=2) the morphology is characterized as “both extended”, “both compact”, or “one compact/one extended” components. For the case of “both extended” components, the difference in position angle between the two components is calculated as a measure of the alignment.

For triples ($`N`$=3), we identify the “core” component as that which is not an endpoint of the longest leg, $`\mathrm{}`$, of the three triangle legs formed (see Figure 11). Denoting the three legs by $`\mathrm{}_{12}`$, $`\mathrm{}_{23}`$, and $`\mathrm{}_{31}`$, such that $`\mathrm{}=\mathrm{}_{12}>\mathrm{}_{23}>\mathrm{}_{31}`$, we calculate a symmetry factor, $`S`$, given by $`\mathrm{}_{31}/\mathrm{}_{23}`$, as well as a configuration angle, $`\theta `$, defined geometrically by $`\theta =\mathrm{arccos}\frac{\mathrm{}_{12}}{\mathrm{}_{23}+\mathrm{}_{31}}`$. Thus, $`\theta `$ is the angle defined by going from one focus to the other to the endpoint of the semi-minor axis of an ellipse whose foci ($`F_1`$, $`F_2`$) are separated by a distance $`\mathrm{}_{12}`$ and whose major axis is equal to $`\mathrm{}_{23}+\mathrm{}_{31}`$. This configuration angle is computationally simpler than the traditional bending angle for radio sources, with which it should not be confused. Requiring that $`\theta `$ be less than some chosen value, $`\theta _c`$, ensures that the “core” is located within the ellipse defined by $`\theta _c`$. The requirement that $`\mathrm{}_{12}>\mathrm{}_{23}>\mathrm{}_{31}`$ implies that $`\theta <60^{}`$. In addition, the morphology of the “core” component and the “lobe” components are noted for triples.

Associations with more than three FIRST components are termed “multiples.” For these sources, we find the largest separation, $`\mathrm{}`$, between any two components, and take the configuration angle, $`\theta `$, to be the maximum value obtained (in the manner described above) from all triangles in which one side is given by $`\mathrm{}`$. Furthermore, if one or more components lie within a distance $`\mathrm{}/4`$ of the average position, the program identifies as a “core” that component nearest to the center; otherwise, no core is identified.

###### Acknowledgements.

Acknowledgements:  We are grateful to Bob Becker, Emilio Falco, Mike Irwin, Don Schneider, and Ion Yadigaroglu for providing optical spectra and CCD images. It is also a pleasure to thank David Helfand and George Djorgovski for many helpful comments. JL and CSK are supported by the Smithsonian Institution. This work was partially supported by NASA/HST grants GO-7495 and GO-7887, and by NSF/AST-9802732.


format: text

text:
# Dissipative Collapse of a Spherical Cluster of Gas Clouds

## 1 Introduction

The formation of stars and galaxies normally proceeds from the collapse of gas clouds. Hence it is imperative that one investigates the evolution of the size and eccentricity of spherical and spheroidal clusters of gas clouds. Such an investigation for point masses with an isotropic distribution of velocities was made by Chandrasekhar and Elbert (Chandrasekhar and Elbert 1972 (1972)) using the scalar and tensor virial equations. This was followed by Som Sunder and Kochhar (Som Sunder and Kochhar 1985 (1985)), who further examined (Som Sunder and Kochhar 1986 (1986)) the evolution for anisotropic velocity distributions also. Here we examine the dissipative evolution of a spherical cluster of gas clouds with an isotropic velocity distribution, using the scalar virial equation, in an attempt to obtain the time scale, for collapse impeded by virialized mass motions. The gravitational binding energy released during collapse feeds the random motions of the clouds. This kinetic energy is dissipated in cloud-cloud collisions. The evolution of the cluster size is followed using the scalar virial equation. It is supplemented with an evolution equation for the kinetic energy of random motions. From this investigation we get quantitative relations between the collapse time and the collision time under various conditions.

## 2 The model

We consider a spherically symmetric cluster of $`N`$ equal mass clouds of radius $`R_c`$ and mass $`M_c`$ distributed uniformly. The clouds have a one dimensional r.m.s. velocity $`v_{rms}`$. The gas clouds are in pressure balance with an intercloud medium. The mass in the intercloud medium is taken to be too small to affect the dynamics and the contribution from it to the virial as well as the dissipation will be ignored throughout (Mestel and Paris Mestel and Paris 1984 (1984)). The scalar virial equation for the system is (Chandrasekhar and Elbert Chandrasekhar and Elbert 1972 (1972))

$`{\displaystyle \frac{1}{2}}{\displaystyle \frac{d^2I}{dt^2}}=2T+\mathrm{\Omega }`$ (1)
where $`I=4\pi \rho (r)r^4𝑑r`$ is like the moment of inertia of the system about the centre and $`\rho (r)`$ is the density at a distance $`r`$ from the origin. $`T`$ is the kinetic energy associated with mass motions, and $`\mathrm{\Omega }`$ is the self gravitational potential energy of the cluster. The evolution of the system must be consistent with equation 1 and we may use it to investigate the evolution of the gross size of the cluster. For the cluster we get

$`I={\displaystyle \frac{3}{5}}MR^2(1+({\displaystyle \frac{R_c}{R}})^2)`$ (2)
and

$`\mathrm{\Omega }={\displaystyle \frac{3}{5}}{\displaystyle \frac{GM^2}{R}}(15({\displaystyle \frac{R_c}{R}})^2).`$ (3)
Thus for $`R_c`$ much smaller than R using the expressions for I and $`\mathrm{\Omega }`$ corresponding to a homogeneous distribution of matter would not be too erroneous.

The velocity of the clouds may be separated into two parts, a random part and a mean motion. The kinetic energy for the random motions say $`T_R`$ is $`\frac{3}{2}Mv_{rms}^2`$. If we consider only homologous changes in the size of the system with the mean velocity at a point being proportional to the distance of the point from the centre, we get the kinetic energy for the mean motion as $`T_M=\frac{3}{10}M(\frac{dR}{dt})^2`$. The cluster is assumed to be uniform at all times. Introducing

$`I={\displaystyle \frac{3}{5}}MR^2,T_M={\displaystyle \frac{3}{10}}M({\displaystyle \frac{dR}{dt}})^2,\mathrm{\Omega }={\displaystyle \frac{3}{5}}{\displaystyle \frac{GM^2}{R}},`$ (4)
and

$`T_R={\displaystyle \frac{3}{2}}Mv_{rms}^2`$ (5)
the virial equation becomes

$`{\displaystyle \frac{d^2R}{dt^2}}={\displaystyle \frac{5v_{rms}^2}{R}}{\displaystyle \frac{GM}{R^2}}.`$ (6)
The evolution equation for the random kinetic energy will have a growth term due to the feeding of gravitational potential energy into the random motions and a dissipative term due to cloud collisions. In the case of nonuniform collapse the gravitational energy released during contraction goes initially into radial flow, but eventually some of it is converted into the energy of random motions. The possibility that all the released binding energy is converted into random kinetic energy has also been advocated (Odgers and Stewart Odgers and Stewart 1958 (1958)). We consider that all the released potential energy is converted to random kinetic energy. This gives

$`{\displaystyle \frac{dT_R}{dt}}_{growth}={\displaystyle \frac{d\mathrm{\Omega }}{dt}}`$ (7)
For dissipation by collisions the dissipation rate will be proportional to (the number of collisions suffered by a cloud in unit time) x (the energy loss per collision) x (the total number of clouds). The number of collisions per unit time will be proportional to the (geometrical cross section of the cloud) x (the velocity of the cloud) x (the number density of clouds). The dissipation rate is thus proportional to $`Mv_{rms}^3NR_c^2/R^3=`$$`\nu `$$`Mv_{rms}^3/R`$. Here we have expressed $`R_c`$ in terms of the cloud filling factor f as $`R_c=(\frac{f}{N})^{\frac{1}{3}}R`$. Putting all the numerical factors into a parameter $`a`$ we write $`\nu =aN^{\frac{1}{3}}f^{\frac{2}{3}}`$. Here $`\nu `$ is proportional to the number of collisions in a free fall time for virial theorem random motions in the system. As an example, Elmegreen(Elmegreen 1985 (1985)) when considering that the initial speed of a cloud is halved after interacting with a mass equal to itself via collisions, gets $`a=\frac{9\sqrt{6}}{16}`$. We thus get

$`{\displaystyle \frac{dT_R}{dt}}_{loss}=\nu {\displaystyle \frac{Mv_{rms}^3}{R}}`$ (8)
Putting all the terms together we get the evolution equation for the random kinetic energy $`T_R`$ in the form

$`{\displaystyle \frac{dv_{rms}^2}{dt}}={\displaystyle \frac{2GM}{5R^2}}{\displaystyle \frac{dR}{dt}}\nu {\displaystyle \frac{2v_{rms}^3}{3R}}`$ (9)
As the cluster contracts it is the contraction produced by the kinematics of the collapse itself that mainly affects the radius of the clumps rather than any self gravitational effects Lacey(Lacey 1984 (1984)). For clouds which are themselves Jeans’ unstable the collapse time for the clouds will be smaller than the free fall collapse time for the cluster only by a factor $`\sqrt{f}`$. The number of the clouds can vary as there can be cloud disruptions or coalescences also during collisions. We assume a constant $`\nu `$ throughout the evolution.

For $`v_{rms}^2`$$``$$`\frac{GM}{R_0}`$ corresponding to virial equilibrium the collision time is $``$$`1/\nu `$ in units of the free fall time $`t_{ff}`$$``$$`\sqrt{(\frac{R_0^3}{GM})}`$, where $`R_0`$ is the initial radius. Thus since dissipation is due to collisions the dissipation time is also $``$$`1/\nu `$. For the system to be collisional the mean free path of the clouds should be less than twice the radius of the cloud. This yields the condition $`N^{\frac{1}{3}}f^{\frac{2}{3}}>1/6`$. This condition is independent of the size of the cluster since we have assumed a constant filling factor throughout the collapse.

## 3 Results

We have considered evolution of the cluster under two different conditions: (a) the system is initially in virial equilibrium with $`T_R=\mathrm{\Omega }/2`$ (b) the system has zero initial kinetic energy. As case (c) we have considered the evolution of the cluster starting from virial equilibrium and with no gravitational feeding . Equations (6) and (9) were solved for various values of the only free parameter $`\nu `$$``$$`1/t_c`$ where $`t_c`$ is the collision time for virial theorem random motions in units of $`t_{ff}`$. The equations were normalized and integrated with a constant time step. Desired accuracy was achieved by using a time step which was appropriately small. It was checked that with zero dissipation the change in the total energy was less than one percent over 10 units of normalized time.

In Fig. 1 we give R as a function of time for the cases (a) and (b) for two different values of $`\nu `$, $`.2`$ and $`0.5`$. The two curves with faint oscillations are for case (a). The curves with strongly marked oscillations are for case (b) and the curves for no gravitational feeding, case (c) show monotonic collapse. For larger values of $`\nu `$ the oscillations start at smaller radii.

In Fig. 2 we show $`t_{10}`$ the time for collapse to one-tenth size in units of $`t_{ff}`$ as a function of $`t_c(=(aN^{1/3}f^{2/3})^1)`$, The uppermost curve is for case (a), the lowermost for case (b), and the middle curve is for case (c). For large values of $`t_c`$, $`t_{10}`$ depends linearly on $`t_c`$ and a fit to the curves in this region gives, for case (a)

$`t_{10}=0.95+2.0t_c,`$ (10)
and for case (b) we get

$`t_{10}=0.74+0.6t_c.`$ (11)
In the nonlinear region the curves can be represented with very good accuracy by a quadratic in $`t_c`$. For case (a) we get

$`t_{10}=1.25+3.0t_c1.4t_c^2`$ (12)
and for case (b)

$`t_{10}=1.16+0.63t_c0.46t_c^2`$ (13)
Fig.2 also shows the case of no gravitational feeding (case (c)), for comparison with earlier work where the dissipation of virial theorem random motions was investigated by various authors. The cluster is initially in virial equilibrium. For small dissipation i.e. for large collision times $`t_{10}`$ is not very sensitive to $`t_c`$ and increases slowly with $`t_c`$ in a linear fashion as $`t_{10}=2.2+0.2t_c`$. However for large dissipation, there is a fast increase of the collapse time as the collision time increases. This nonlinear portion is fitted by the curve $`t_{10}=1.2+2.1t_c1.38t_c^2`$.

## 4 Discussion

For a nondissipative system starting with zero kinetic energy Chandrasekhar and Elbert (Chandrasekhar and Elbert 1972 (1972)) have shown that the system will execute oscillations. If we put zero dissipation in our system we get the same oscillations. For the finite dissipation cases which we have considered we find these oscillations; however they are of a damped character due to the finite value of the dissipation.

Examining $`R`$ as a function of $`t`$, in all the cases with gravitational feeding we find that collapse takes place in an oscillatory fashion. The oscillations are well marked in case (b). Here there is an initial rapid collapse when $`T_R`$ rapidly builds up and finally reaches a value high enough for the virial equilibrium condition to hold. This takes place when the value of $`\frac{dR}{dt}`$ is quite large and the resulting bounce of the system is quite hard. This causes the oscillations to be strongly marked. As $`\nu `$ increases the bounce is delayed due to the the slower increase of $`T_R`$. In case (a) the bounce can take place quite early on and is quite mild. For large dissipation the bounce gets delayed so much that we do not get it within collapse to one tenth size and the collapse looks apparently monotonic.

It is known that a self gravitating cloud collapsing without dissipation will satisfy the conditions for virial equilibrium when it has collapsed to half radius (Hoyle Hoyle 1949 (1949)). The idea of virialization of mass motions in gravitational collapse and their subsequent dissipation has been considered by many people.(Odgers and Stewart Odgers and Stewart 1958 (1958), Larson Larson 1977 (1977), Mestel Mestel 1977 (1977), Scalo and Pumphrey Scalo and Pumphrey 1982 (1982), Elmegreen Elmegreen 1985 (1985).) These works consider virial theorem random motions and the time scale for their dissipation. Simulations of the collisional dissipation of virial theorem random motions by Scalo and Pumphrey give a dissipation time $``$$`20`$ free fall times. (See also Chieze and Lazareff Chieze and Lazareff 1980 (1980), Hausman Hausman 1981 (1981).) Elmegreen considered a magnetized cloud of clumps and obtained dissipation times of the order $`35`$ cloud crossing times.

We have analyzed collapse with gravitational feeding, of clusters starting from virial equilibrium and also starting with zero kinetic energy (like for example from turn around). For small values of the dissipation like say $`t_c=5t_{ff}`$ we get $`t_{10}`$ about $`11`$ for case (a) and about $`3.8`$ for case (b). In the model of Scalo and Pumphrey (Scalo and Pumphrey 1982 (1982)) the number of clouds get smaller with time due to cloud coalescence. The resultant slower dissipation rate explains the comparatively longer dissipation times they get. As a general statement we may say that for small dissipation, clusters starting from virial equilibrium can take more than twice the dissipation time to collapse, when there is gravitational feeding. With no gravitational feeding we get $`t_{10}>t_{ff}>t_c`$. The increase in the collapse time as the collision time increases is less pronounced in this case and for large collision times lies in the range $`24`$ free fall times. This is consistent with the results obtained by other authors for the dissipation time for virial theorem random motions. As $`t_c`$ goes to zero $`t_{10}`$ rapidly approaches the value $`1.1`$ which is the time for free fall collapse to one tenth size in the units we have chosen. Our results indicate that for large collision times gravitational feeding can significantly increase the time scales for the dissipative collapse of clusters of gas clouds.

For large values of $`t_c`$ i.e. small values of the dissipation we see from fig. 1 that the collapse time depends linearly on $`t_c`$. For case (a) this may be understood as follows. For small dissipation the cluster is collapsing very slowly. Thus growth as well as loss of $`T_R`$ is small and the cluster contraction is governed mainly by the decrease of $`T_R`$ in the virial equation. This decrease is in proportion to $`t_c`$. For case (b) also after an initial rapid collapse, when $`T_R`$ rapidly builds up, the size of the cluster is governed similarly by the decrease of $`T_R`$ in equation (6). For large dissipation however the collapse being rapid, changes in $`\frac{dR}{dt}`$ and $`R`$ are large causing large changes in $`T_R`$. There is thus a strong and nonlinear dependence of the evolution of $`T_R`$ on $`t_c`$. The collapse time also reflects this strong nonlinear dependence.

For the system under consideration the free fall time provides the appropriate time scale for the system to adjust to changes in the virial theorem (Brosche 1970). While noticing that a quasi virial evolution does not hold for dissipation times shorter than the dynamical time, we find distinctive difference in the dependance of $`t_{10}`$ on $`t_c`$ in the two regions $`t_c>t_{ff}`$ and $`t_c<t_{ff}`$.

The consequences of this investigation for the formation of stars and galaxies are being examined at present in a quantitative fashion and will be reported separately.

###### Acknowledgements.

KI thanks Profs. Rajaram Nityananda, and Kandaswamy Subramanian of NCRA, Pune for helpful discussions. She also thanks IUCAA, Pune for an associateship during part of this investigation.


text:
# Eigenvalue Distributions of the QCD Dirac Operator

## Abstract

We compute by Monte Carlo methods the individual distributions of the $`k`$th smallest Dirac operator eigenvalues in QCD, and compare them with recent analytical predictions. We do this for both massless and massive quarks in an SU(3) gauge theory with staggered fermions. Very precise agreement is found in all cases. As a simple by-product we also extract the microscopic spectral density of the Dirac operator in SU(3) gauge theory with dynamical massive fermions for $`N_f=1`$ and 2, and obtain high-accuracy agreement with analytical expressions.

NBI-HE-00-34

CSIT2000-17

NORDITA-2000/61 HE

hep-lat/0007041

In the chiral limit of QCD it is possible to derive exact analytical expressions for eigenvalue distributions and eigenvalue correlations of the Dirac operator. This situation arises if one restricts oneself to the phase of spontaneously broken chiral symmetry (with infinite-volume zero-mass condensate $`\mathrm{\Sigma }`$), and considers the theory in a finite space-time volume $`V`$. Taking the chiral limit by sending quark masses $`m_i`$ to zero in such a way that the combination $`\mu _im_i\mathrm{\Sigma }V`$ is kept fixed as $`V\mathrm{}`$, this leads to an interesting finite-volume scaling regime of QCD. Denoting generically all pseudo-Goldstone masses by $`m_\pi `$, exact scaling is achieved if $`m_\pi V^{1/4}`$ . The only assumption is that chiral symmetry is spontaneously broken from one group $`G`$ to a smaller group $`H`$. The coset $`G/H`$ specifies the universality class .

The most detailed analytical predictions were originally obtained from an intriguing exact relation to universal Random Matrix Theory results . Later, the same results have been derived directly from the effective partition functions  based on the so-called supersymmetric formulation . Very recently, also the replica method  has been shown to yield identical results . Beyond any doubt now, these analytical predictions represent exact statements about the Dirac operator spectrum of QCD in the above finite-volume scaling region. This is an exciting situation from the point of view of lattice gauge theory, which is almost tailored to study such a finite-volume regime. It is one of the few instances where there are exact and non-perturbative analytical predictions that immediately can be compared with lattice Monte Carlo results. There have by now been several Monte Carlo studies of these analytical predictions for Dirac operator eigenvalues in various four-dimensional gauge theories . One peculiarity of these predictions is that they are most easily expressed in gauge field sectors of fixed topological charge $`\nu `$111The sum over topology can be performed , but closed analytical forms are not known except in very simple cases. For instance, in theories with massless (dynamical) fermions only the sector with $`\nu =0`$ contributes.. For fermions sensitive to topology this gives an interesting new possibility, as predictions differ markedly for different topological sectors . But staggered fermions seem oblivious to gauge field topology at the gauge couplings that are presently realistic ; predictions can then only be made with respect to the topologically trivial sector of $`\nu =0`$.

In this letter we shall extend the Monte Carlo analyses of Dirac operator eigenvalues in two directions. First, we shall confront some very recent analytical results concerning the probability distribution of the $`k`$th smallest Dirac operator eigenvalue  with Monte Carlo data. This is important, because there is much more detailed information in these individual eigenvalue distributions than in, for example, the spectral density, which is the sum of all these individual distributions (see below). Second, we shall make the extension of the quenched results of Ref.  to the case of dynamical fermions of non-zero mass. Until now, Monte Carlo analyses of microscopic Dirac operator spectra with dynamical fermions have been limited to a different universality class222An SU(2) gauge theory labeled by the so-called “Dyson index” $`\beta =4`$ for staggered fermions . than the one relevant for QCD  (see also Ref. ). We study gauge group SU(3), which, even with staggered fermions, allows us to compare with the universality class relevant for continuum QCD (the class of Dyson index $`\beta =2`$).

Details of how we compute the lowest-lying Dirac operator eigenvalues have been given in Ref. , and we shall not repeat them here. Instead, let us briefly recall the analytical predictions with which we should compare our Monte Carlo data. We start with the recently derived expression for the distribution of the $`k`$th Dirac operator eigenvalue in the finite-volume scaling region mentioned above. The generator for these individual eigenvalue distributions is the joint probability distribution for the $`k`$ smallest eigenvalues. It turns out that this joint probability distribution can be written very compactly in terms of finite-volume partition functions with additional fermion species. For the universality class relevant for QCD, the explicit expression is :

$`\omega _k(\zeta _1,\mathrm{},\zeta _k;\{\mu \})=Ce^{\zeta _k^2/4}\zeta _k{\displaystyle \underset{i=1}{\overset{k1}{}}}\left(\zeta _i^{2\nu +1}{\displaystyle \underset{j=1}{\overset{N_f}{}}}(\zeta _i^2+\mu _j^2)\right){\displaystyle \underset{i>j}{\overset{k1}{}}}(\zeta _i^2\zeta _j^2)^2{\displaystyle \underset{j=1}{\overset{N_f}{}}}\mu _j^\nu `$ (1)
$`\times {\displaystyle \frac{𝒵_2(\left\{\sqrt{\mu _i^2+\zeta _k^2}\right\},\sqrt{\zeta _k^2\zeta _1^2},\sqrt{\zeta _k^2\zeta _1^2},\mathrm{},\sqrt{\zeta _k^2\zeta _{k1}^2},\sqrt{\zeta _k^2\zeta _{k1}^2},\stackrel{\nu }{\stackrel{}{\zeta _k,\mathrm{},\zeta _k}})}{𝒵_\nu (\{\mu \})}}.`$ (2)
Here the additional $`2(k1)`$ fermion species of masses $`\sqrt{\zeta _k^2\zeta _i^2}`$ are doubly degenerate, as indicated. The partition function in the numerator is evaluated in a sector of fixed topological charge $`\nu =2`$ (independent of the actual topological charge in question). But there is information about $`\nu `$ in the fact that $`\nu `$ additional fermion species (all degenerate in mass $`\zeta _k`$) enter in the partition function. Actually, since we shall be dealing with staggered fermions in the present paper, only the sector with $`\nu =0`$ is relevant. In that case the above formula simplifies:

$`\omega _k(\zeta _1,\mathrm{},\zeta _k;\{\mu \})=Ce^{\zeta _k^2/4}\zeta _k{\displaystyle \underset{i=1}{\overset{k1}{}}}\left(\zeta _i{\displaystyle \underset{j=1}{\overset{N_f}{}}}(\zeta _i^2+\mu _j^2)\right){\displaystyle \underset{i>j}{\overset{k1}{}}}(\zeta _i^2\zeta _j^2)^2`$ (3)
$`\times {\displaystyle \frac{𝒵_2(\left\{\sqrt{\mu _i^2+\zeta _k^2}\right\},\sqrt{\zeta _k^2\zeta _1^2},\sqrt{\zeta _k^2\zeta _1^2},\mathrm{},\sqrt{\zeta _k^2\zeta _{k1}^2},\sqrt{\zeta _k^2\zeta _{k1}^2})}{𝒵_0(\{\mu \})}}.`$ (4)
The fact that in particular the smallest (non-zero) Dirac operator eigenvalue has a distribution determined entirely in terms of the effective partition function was noted earlier , but the present most general expression is actually more compact. The proportionality factor $`C`$ depends on the chosen normalization of the partition functions. We choose a normalization in which

$$𝒵_\nu (\{\mu \})=detA(\{\mu \})/\mathrm{\Delta }(\mu ^2),$$
(5)
where the determinant is taken over the $`N_f\times N_f`$ matrix

$$A(\{\mu \})=\mu _i^{j1}I_{\nu +j1}(\mu _i),$$
(6)
and

$$\mathrm{\Delta }(\mu ^2)=\underset{i>j}{\overset{N_f}{}}(\mu _i^2\mu _j^2).$$
(7)
With this convention the normalization factor is $`C=1/2`$ for all values of $`k`$, $`N_f`$ and $`\nu `$.

For the purpose of comparing with lattice gauge theory data, a more convenient quantity to focus on is the probability distribution of the $`k`$th smallest Dirac operator eigenvalue . One gets it from the joint probability distribution by integrating out the previous $`k1`$ smaller eigenvalues:

$`p_k(\zeta ;\{\mu \})`$ $`=`$ $`{\displaystyle _0^\zeta }𝑑\zeta _1{\displaystyle _{\zeta _1}^\zeta }𝑑\zeta _2\mathrm{}{\displaystyle _{\zeta _{k2}}^\zeta }𝑑\zeta _{k1}\omega _k(\zeta _1,\mathrm{},\zeta _{k1},\zeta ;\{\mu \})`$ (8)
$`=`$ $`{\displaystyle \frac{1}{(k1)!}}{\displaystyle _0^\zeta }𝑑\zeta _1{\displaystyle _0^\zeta }𝑑\zeta _2\mathrm{}{\displaystyle _0^\zeta }𝑑\zeta _{k1}\omega _k(\zeta _1,\mathrm{},\zeta _{k1},\zeta ;\{\mu \}).`$ (9)
The individual eigenvalue distributions by definition sum up to the microscopic spectral density $`\rho _S(\zeta ;\{\mu \})`$:

$$\rho _S(\zeta ;\{\mu \})=\underset{k=1}{\overset{\mathrm{}}{}}p_k(\zeta ;\{\mu \}).$$
(10)
But a much simpler way to get the microscopic spectral density is by direct evaluation. As we will in addition be interested in the case of dynamical fermions with non-zero quark masses $`\mu _i`$, we need the microscopic spectral density for the general massive case of this universality class . Also this quantity can most simply be expressed in terms of the finite-volume partition functions :

$$\rho _S(\zeta ;\{\mu \})=(1)^\nu |\zeta |\underset{j=1}{\overset{N_f}{}}(\zeta ^2+\mu _j^2)\frac{𝒵_\nu (\{\mu \},i\zeta ,i\zeta )}{𝒵_\nu (\{\mu \})},$$
(11)
using the same normalization convention as above. This spectral density is “double-microscopic” in the sense that both eigenvalues $`\lambda `$ and masses $`m_j`$ have been rescaled according to $`\zeta =\lambda \mathrm{\Sigma }V`$ and $`\mu _j=m_j\mathrm{\Sigma }V`$. As the rescaled masses $`\mu _j`$ are sent to infinity, one obtains a series of decoupling relations where the corresponding quark flavors become effectively quenched . Because we precisely want to see the effects of having dynamical fermions, it is thus important to keep the rescaled masses $`\mu _i`$ on roughly the same scale as the eigenvalues, or smaller.

We begin by comparing quenched lattice gauge theory data with the predictions for the individual eigenvalue distributions. We have here large statistics for the 10 smallest Dirac operator eigenvalues on an $`8^4`$ lattice at coupling $`\beta =5.1`$. These parameters and all others used in this paper are summarized in Table 1. In Fig. 1 we show a comparison of these lattice Monte Carlo results with the formula (9) for $`k=1\mathrm{}10`$. $`\mathrm{\Sigma }`$ was determined from a fit of the distribution of the lowest eigenvalue. The formula (9) obviously gets quite tedious to evaluate by numerical integration for large values of $`k`$. The analytic curves for the distributions of the last couple of eigenvalues have actually not been obtained directly from (9), but rather from the eigenvalue distributions of computer generated random matrices (which obviously amounts to the same). The agreement with our lattice Monte Carlo data is nothing less than spectacular for all the eigenvalues we have available. Eventually, on any given finite lattice, $`1/V`$-corrections will take over, and this very precise agreement with the theoretical predictions will no longer hold. But in the limit $`V\mathrm{}`$ the agreement should persist for an infinity of smallest Dirac operator eigenvalues $`\lambda _n`$, rescaled so that $`\zeta _n\lambda _n\mathrm{\Sigma }V`$ remains finite in the infinite-volume limit. Agreement with these lattice data and the microscopic spectral density has already been presented in Ref. , but after just a few oscillations the spectral density quickly becomes an almost constant function, with very little information. For the first time we see here that there is detailed agreement with the analytical predictions for a whole sequence of individual eigenvalues underneath.

We next turn to simulations with dynamical fermions of finite mass. As mentioned above, there has previously been only one Monte Carlo study in this direction , for a different gauge group (and hence different universality class) SU(2). We shall here provide first results for a study with gauge group SU(3), and hence the universality class of continuum QCD. We shall not only compare with the microscopic spectral density, but also focus on individual eigenvalue distributions, just as we did above in the quenched case.

For our dynamical simulations we employed the Hybrid Monte Carlo algorithm with $`N_f=1`$ or 2 species of staggered fermions. In the continuum limit, this would correspond to $`n_f=4N_f`$ continuum fermion flavors. However, in order to have a large physical volume with few lattice sites — our computer resources limit us to the use of systems with $`4^4`$ up to $`6^4`$ lattice sites and the finite size scaling limit we want to consider requires a sufficiently large physical volume, $`V\mathrm{\Lambda }_{QCD}^4`$ — we are forced to work at a large lattice spacing and hence a large gauge coupling. In such a situation the flavor symmetry of the staggered fermions is badly broken, especially in the limit of small masses of interest here. Only the $`N_f^2`$ true pseudo-Goldstone pions from the spontaneous breaking of the $`U(N_f)\times U(N_f)`$ symmetry to the diagonal $`U(N_f)`$ symmetry are light. The flavor symmetry breaking makes the other $`n_f^21N_f^2`$ “would-be” pseudo-Goldstone pions much heavier. As far as the infrared limit is concerned, only $`N_f`$ flavors contribute. This will be verified by our results presented below.

In Fig. 2 we compare the spectral density of the simulation data for different quark masses. $`\mathrm{\Sigma }`$ was determined, for each of the three ensembles, from a fit of the distribution of the lowest eigenvalue. The values agreed within errors, and were averaged to obtain the result in Table 1. This agreement is of course as it should be: $`\mathrm{\Sigma }`$ is the infinite-volume chiral condensate in the chiral limit $`m0`$ (and in particular mass-independent). For each mass we plot the analytical result (11) for that particular mass value (rescaled according to $`\mu =m\mathrm{\Sigma }V`$), and include the curve for the limiting cases of the quenched result corresponding to infinite mass quarks, as well as for one massless quark. We see that the data beautifully interpolate between the two limits and follow the analytical curve almost perfectly. Looking at the individual eigenvalue distributions in Fig. 3 we get a better view on how far the agreement goes. It is clear that only the very first eigenvalues match perfectly and that there is a very tiny discrepancy for the last (of the 6) distributions due to the limited volume, but still the agreement is impressive. We did a single run on a larger volume ($`6^4`$), and indeed the agreement is better, see Fig. 4 B. Finally, data from a run corresponding to two flavors of quarks are presented in Fig. 4 A, with the same degree of agreement.

To conclude, we have found very precise agreement with the analytical predictions for the distributions of Dirac operator eigenvalues in the microscopic finite-volume scaling region. This is seen in both quenched and unquenched simulations, and for the first time there are now high-statistics comparisons between the analytical predictions for the Dirac operator spectrum of massive quarks and lattice Monte Carlo data for gauge group SU(3). Our Monte Carlo calculations have been restricted to staggered fermions, which are insensitive to gauge field topology at the couplings considered here. We have thus compared with analytical predictions for the topological sector of charge $`\nu =0`$. Perhaps the biggest surprise is that the previously seen agreement with analytical predictions for the summed-up spectral density holds at a much more detailed level underneath. The individual smallest Dirac operator eigenvalues are really completely locked up by gauge field topology, and, away from the exact zero modes, by the coset $`G/H`$ of spontaneous chiral symmetry breaking.

Acknowledgments: The work of P.H.D. and K.R. has been partially supported by EU TMR grant no. ERBFMRXCT97-0122, and the work of U.M.H. has been supported in part by DOE contracts DE-FG05-96ER40979 and DE-FG02-97ER41022. P.H.D. and U.M.H. also acknowledge the financial support of NATO Science Collaborative Research Grant no. CRG 971487.


text:
# On minimal eigenvalues of Schrödinger operators on manifolds

## 1. Introduction

In the last years there has been a great interest in the study of optimal properties of eigenvalues of Schrödinger operators of the form $`H=\mathrm{\Delta }+V`$ defined on compact manifolds, when some restrictions are imposed on the potential $`V`$. Some of these problems are related to several physical phenomena such as motion by mean curvature, electrical properties of nanoscale structures, etc (see, for instance, \[A, AHS, EI, EHL, HL, Ke\] and the references therein).

In \[EHL\], the authors considered the case of potentials depending on the curvature $`\kappa `$ and studied the problem of minimizing the first eigenvalue of the operator $`H=d^2/ds^2+\alpha \kappa ^2`$ defined on a closed planar curve with length one. They proved that for $`0<\alpha <1/4`$ the circle is the unique minimizer, while for $`\alpha >1`$ this is no longer the case, leaving open the question of the value of $`\alpha `$ where the transition takes place, and also what happens after this critical value.

More generally, one might consider an operator $`H`$ defined on a compact $`n`$manifold $`(M,g)`$ by $`H=\mathrm{\Delta }+\alpha F(\kappa )`$ and with eigenvalues $`\lambda _0<\lambda _1\mathrm{}`$, and study the problem of determining

$$\mathrm{\Lambda }_j(\alpha )=\underset{\kappa K}{inf}\lambda _j(\kappa ),j=0,1,\mathrm{}$$

where

$$K=\{\kappa 𝒞(M;):\frac{1}{|M|}_M\kappa 𝑑v_g=\kappa _0\}.$$

In particular, we are interested in knowing whether or not there exists a critical value of $`\alpha `$, say $`\alpha _c`$, where the constant potential stops being a global minimizer for the first eigenvalue. In this paper we show that in the one–dimensional case studied in \[EHL\] this critical value is in fact equal to $`1/4`$, and that for $`\alpha `$ larger than $`\alpha _c`$ the infimum is identically equal to $`\pi ^2`$ and is not attained. The first part of this result is a consequence of a more general result which provides an upper bound for $`\alpha _c`$ holding in any dimension. Furthermore, we show that for potentials of the form $`\kappa =\kappa _0+\epsilon q`$ where $`q`$ has zero average, this bound is in fact precise for sufficiently small values of $`\epsilon `$, in the sense that for $`\alpha `$ smaller than the bound, the constant potential $`\kappa _0`$ gives a smaller eigenvalue than $`\kappa `$, while for larger values of $`\alpha `$ this is not always the case.

These results could lead us to expect that results similar to those in one dimension would also hold in higher dimensions, that is, that there would exist a nontrivial interval $`(0,\alpha _c)`$ where the constant potential was the unique minimizer. However, it turns out that for dimensions higher than the first there exist potentials satisfying the given restrictions and which make the principal eigenvalue as close to zero as desired. Thus, we see that in this case the constant potential is never a global minimizer. It remains an open question if it is a local minimizer. A similar statement also holds for higher eigenvalues and for minimizations subject to other types of integral restrictions – see Theorem 4 and the remarks that follow it. The reason for this different behaviour in dimensions higher than the first is directly related to the fact that in this case, given a manifold $`M`$ and a geodesic ball $`B_\delta `$ of radius $`\delta `$ centred at a point $`x_0`$ in $`M`$, the Dirichlet eigenvalues of the Laplacian in $`\mathrm{\Omega }_\delta =MB_\delta `$ converge to those of the Laplacian in $`M`$ as $`\delta `$ approaches zero – for a more precise statement of this property see Section 4 and \[CF\].

Finally, we point out that the results in one dimension enable us to obtain a lower bound for the principal eigenvalue in the case of a general potential (Corollary 3) which, of course, corresponds also to the first eigenvalue of Hill’s equation. Note that one of the motivations behind the study of the minimization of eigenvalues when the potential is subject to integral restrictions was precisely to obtain lower bounds for eigenvalues – see \[Ke\].

## 2. Notation and general local results

Let $`(M,g)`$ be a compact Riemannian $`n`$manifold with metric $`g`$ and let $`\mathrm{\Delta }`$ denote the Laplace–Beltrami operator defined on $`M`$ with eigenvalues $`0=\mu _0<\mu _1\mathrm{}`$ repeated according to their multiplicity. Denote the corresponding orthonormal (with respect to the $`L^2(M)`$ inner product induced by the Riemannian measure $`v_g`$) system of eigenfunctions by $`\left\{v_j\right\}_{j=0}^{\mathrm{}}`$. Consider now the operator defined on $`M`$ by $`H=\mathrm{\Delta }+\alpha F(\kappa )`$, where

$$\frac{1}{|M|}_M\kappa 𝑑v_g=\kappa _0,$$

and $`F:`$ is assumed to be of class $`𝒞^3`$ in a neighbourhood of $`\kappa _0`$.

The main result in this section is then the following

###### Theorem 1.

Assume that $`F^{}(\kappa _0)0`$ and define

$$\alpha ^{}=\frac{\mu _1F^{\prime \prime }(\kappa _0)}{2[F^{}(\kappa _0)]^2}.$$

Then, if $`q`$ is a continuous real valued function with zero average and not identically zero, we have that for $`\kappa =\kappa _0+\epsilon q`$ with sufficiently small $`\epsilon `$ (depending on $`q`$ and $`\alpha `$), the principal eigenvalue $`\lambda _0`$ of $`H`$ satisfies

$$\lambda _0(\kappa )>\alpha F(\kappa _0),\text{ if }\mathrm{\hspace{0.33em}0}<\alpha <\alpha ^{},$$

while for $`\alpha >\alpha ^{}`$ there exist functions $`q`$ as above for which

$$\lambda _0(\kappa )<\alpha F(\kappa _0).$$

###### Proof.

Consider the Schrödinger operator defined on $`M`$ by

$$H_\epsilon =\mathrm{\Delta }+\alpha F(\kappa _0+\epsilon q).$$

Since $`M`$ is compact, the spectrum of $`H_\epsilon `$ is discrete and its first (simple) eigenvalue and the corresponding (normalized) eigenfunction are analytic functions of the (real) parameter $`\epsilon `$ \[Ka\]. We thus expand $`\lambda _0`$ and the corresponding eigenfunction $`u`$ as a power series of $`\epsilon `$ around zero:

$$\begin{array}{c}\lambda _0=\mathrm{}_0+\mathrm{}_1\epsilon +\mathrm{}_2\epsilon ^2+\mathrm{}\hfill \\ u=\varphi _0+\varphi _1\epsilon +\varphi _2\epsilon ^2+\mathrm{}.\hfill \end{array}$$

On the other hand, we also have that $`F(\kappa _0+\epsilon q)=f_0+f_1q\epsilon +f_2q^2\epsilon ^2+o(\epsilon ^2)`$, where

$$f_0=F(\kappa _0),f_1=F^{}(\kappa _0),\text{ and }f_2=\frac{1}{2}F^{\prime \prime }(\kappa _0).$$

Substituting these expressions in the equation giving the eigenvalues for $`H_\epsilon `$ we obtain, equating like powers in $`\epsilon `$,

$$\begin{array}{cc}\epsilon ^0:\hfill & \mathrm{\Delta }\varphi _0+\alpha f_0\varphi _0=\mathrm{}_0\varphi _0\\ \epsilon ^1:\hfill & \mathrm{\Delta }\varphi _1+\alpha f_0\varphi _1+\alpha f_1q\varphi _0=\mathrm{}_0\varphi _1+\mathrm{}_1\varphi _0\\ \epsilon ^2:\hfill & \mathrm{\Delta }\varphi _2+\alpha f_0\varphi _2+\alpha f_1q\varphi _1+\alpha f_2q^2\varphi _0=\mathrm{}_0\varphi _2+\mathrm{}_1\varphi _1+\mathrm{}_2\varphi _0.\end{array}$$

From the first equation it follows that $`\mathrm{}_0=\alpha f_0`$ and that $`\varphi _0`$ is constant, which we take to be one. Substituting this in the equation for $`\epsilon ^1`$ and integrating over $`M`$ gives that $`\mathrm{}_1`$ vanishes and $`\varphi _1`$ satisfies

(2.1)
$$\mathrm{\Delta }\varphi _1=\alpha f_1q.$$

Substituting now this in the last equation gives that $`\varphi _2`$ satisfies

$$\mathrm{\Delta }\varphi _2=\alpha f_2q^2\alpha f_1q\varphi _1+\mathrm{}_2.$$

Again integrating over $`M`$ gives

(2.2)
$$\mathrm{}_2=\frac{\alpha f_2}{|M|}_Mq^2𝑑v_g+\frac{\alpha f_1}{|M|}_Mq\varphi _1𝑑v_g.$$

Taking squares on both sides of (2.1) we get $`\left[\mathrm{\Delta }(\varphi _1)\right]^2=\alpha ^2f_1^2q^2`$. On the other hand, multiplying the same equation by $`\varphi _1`$ and integrating over $`M`$ gives that

$$\alpha f_1_Mq\varphi _1𝑑v_g=_M|\varphi _1|^2𝑑v_g.$$

Substituting these two expressions into (2.2) we finally obtain

$$\mathrm{}_2=\frac{f_2}{\alpha f_1^2|M|}_M(\mathrm{\Delta }\varphi _1)^2𝑑v_g\frac{1}{|M|}_M|\varphi _1|^2𝑑v_g,$$

and it follows from Lemma 2.1 below that $`\mathrm{}_2`$ is always positive for $`\alpha <\alpha ^{}`$.

To give an example of a function $`q`$ for which $`\mathrm{}_2`$ becomes negative when $`\alpha >\alpha ^{}`$ it is sufficient to take $`q=v_1`$. We obtain from (2.1) that in this case

$$\varphi _1=c\frac{\alpha }{\mu _1}f_1v_1,$$

where $`c`$ is an arbitrary constant. Substituting this into the expression for $`\lambda _2`$ yields

$$\mathrm{}_2=\frac{\alpha }{|M|}\left(f_2\frac{\alpha }{\mu _1}f_1^2\right),$$

which is negative for $`\alpha >\alpha ^{}`$. ∎

An obvious consequence of this result is that for all $`F`$ of the form above there exists a value of $`\alpha `$, say $`\alpha ^{}`$ such that for $`\alpha >\alpha ^{}`$ the constant potential is not a minimizer of the first eigenvalue.

In the case where $`F^{}`$ is allowed to vanish, it is also clear that if $`\kappa _0=\kappa _0^{}`$ is a (local) minimizer (resp. maximizer) of $`F`$, it follows that, for positive values of $`\alpha `$, $`\kappa (x)\kappa _0^{}`$ will be a (local) minimizer (resp. maximizer). This is the case, for instance, when $`F(\kappa )=\kappa ^2`$ and $`\kappa _0=0`$, where obviously $`\kappa =0`$ is a global minimizer for all $`\alpha `$.

The result needed to prove that $`\mathrm{}_2>0`$ for $`\alpha <\alpha ^{}`$ is neither new nor difficult, but a specific reference could not be found in the literature and so, for the sake of completeness, we provide a proof here.

###### Lemma 2.1.

The functional

$$I_\alpha (u)=_M\alpha (\mathrm{\Delta }u)^2|u|^2dv_g$$

is nonnegative for $`\alpha 1/\mu _1`$.

###### Proof.

The spectral problem corresponding to $`I_\alpha `$ is

(2.3)
$$\alpha \mathrm{\Delta }^2u+\mathrm{\Delta }u=\gamma u,$$

which has discrete spectrum $`\gamma _0\gamma _1\mathrm{}`$. We will prove that if $`\alpha >1/\mu _1`$ then $`\gamma _j0`$ for all $`j=0,1,\mathrm{}`$. To this end rewrite (2.3) as

$$\mathrm{\Delta }(\alpha \mathrm{\Delta }u+u)=\gamma u.$$

It is not difficult to see that $`u`$ is an eigenfunction if and only if

$$\alpha \mathrm{\Delta }u+u=\beta v_j$$

for some real number $`\beta `$ different from zero. For $`\alpha >1/\mu _1`$ the operator $`\alpha \mathrm{\Delta }+I`$ is invertible and thus this last equation has one and only one solution given by $`u=\beta v_j/(1\alpha \mu _j)`$. Substituting this into (2.3) gives $`\gamma =(\alpha \mu _j1)\mu _j`$ from which the result follows. ∎

## 3. The one–dimensional case

In this section we consider the particular case studied in \[EHL\] with $`F(\kappa )=\kappa ^2`$, and for which

$$\alpha ^{}=\frac{\mu _1}{4\kappa _0^2}.$$

As a consequence of Theorem 1 and the results in \[EHL\] we have the following

###### Theorem 2.

In the one dimensional case and for $`F`$ as above, $`\alpha _c=\alpha ^{}`$. Furthermore, for $`\alpha >\alpha _c`$, $`\mathrm{\Lambda }_0(\alpha )\mu _1/4`$.

###### Proof.

It only remains to show the result for $`\alpha `$ larger than $`\alpha _c`$. Clearly in this case $`\mathrm{\Lambda }_0(\alpha )\mu _1/4`$. Consider now the family of potentials given by

$$\kappa _\delta (s)=\{\begin{array}{cc}\kappa _0/\delta ,\hfill & 0<s<\delta ,\hfill \\ 0,\hfill & \delta <s<\mathrm{}.\hfill \end{array}$$

Note that although $`\kappa _\delta `$ is not continuous on the circle, it can be approximated by continuous functions without affecting our results. For this family of potentials we obtain the functional

$$J_\delta (u)=_0^{\mathrm{}}\left[u^{}\right]^2𝑑s+\frac{\alpha \kappa _0^2}{\delta ^2}_0^\delta u^2𝑑s,$$

where $`u`$ is normalized. We now take $`u(s)=\sqrt{2}\mathrm{sin}(\pi s/\mathrm{})`$ to obtain

$$J_\delta (u)=\frac{\mu _1}{4}+\frac{2\alpha \kappa _0^2}{\delta ^2}_0^\delta \mathrm{sin}^2(\frac{\pi s}{\mathrm{}})𝑑s,$$

and since

$$\underset{\delta 0^+}{lim}\frac{{\displaystyle _0^\delta }\mathrm{sin}^2(\pi s)𝑑s}{\delta ^2}=0,$$

it follows that $`J_\delta `$ can be made to be arbitrarily close to $`\mu _1/4`$. ∎

###### Remark 3.1.

Clearly for $`\alpha >\alpha _c`$ the infimum is not attained, as was conjectured in \[EHL\].

A simple consequence of Theorem 2 is a lower bound for the principal eigenvalue of the Schrödinger operator on the circle.

###### Corollary 3.

Consider the operator $`H=d^2/dx^2+V(x)`$ defined on $`(0,L)`$ with periodic boundary conditions, and define

$$V_m=\underset{x}{inf}V(x)\text{ and }I=\frac{1}{L}_0^L\left[V(x)V_m\right]^{1/2}𝑑x.$$

Then

$$\lambda _0\{\begin{array}{cc}V_m+I^2,\hfill & \text{ if }I\frac{\pi }{L}\hfill \\ V_m+\frac{\pi ^2}{L},\hfill & \text{ if }I>\frac{\pi }{L},\hfill \end{array}$$

with equality for $`I<\pi /L`$ if and only if $`V`$ is constant.

###### Proof.

The first inequality follows directly by writing the eigenvalue problem as $`u^{\prime \prime }+(VV_m)u=(\lambda V_m)u`$ and applying the previous corollary with $`\kappa =(VV_m)/\alpha `$. The second part is a consequence of the fact that for $`\alpha `$ larger than $`\alpha _c`$ the principal eigenvalue must be larger than $`\alpha _c\kappa _0^2`$. ∎

###### Remark 3.2.

It follows from Theorem 2 that the given inequalities are sharp in both cases.

## 4. Higher dimensions

In \[EHL\], the proof of the fact that for $`\alpha `$ smaller than $`\alpha ^{}`$ the constant potential is the unique global minimizer of $`\mathrm{\Lambda }_0`$ relied on a result that is not available in higher dimensions. Namely, while in one dimension we have that

$$_{S^1}\left[(uu_m)^{}\right]^2𝑑s\frac{\mu _1}{4}_{S^1}(uu_m)^2𝑑s,$$

where $`u_m`$ is the minimum of $`u`$ in $`S^1`$, from the results in \[CF\] it is known that there is no similar result in higher dimensions. More precisely, if we impose that a function $`f`$ be zero at a finite number of points of a compact manifold with dimension greater than or equal to two, then there is no relation of the form above with a positive constant on the right–hand side. This suggests that an argument similar to that used in the proof of Theorem 2 can now be used for all positive values of $`\alpha `$, and not just for $`\alpha `$ larger than $`\alpha ^{}`$. This is indeed the case, and we have the following

###### Theorem 4.

Assume that $`F(0)`$ is a global minimum of $`F`$. Then, for $`n`$ greater than one, $`\mathrm{\Lambda }_j(\alpha )\mu _jF(0)`$ for all positive $`\alpha `$ and $`j=0,1,\mathrm{}`$.

###### Proof.

Fix a point $`x_0`$ in $`M`$ and denote by $`B_\delta `$ the geodesic ball centred at $`x_0`$ with radius $`\delta `$. Let now $`\mathrm{\Omega }_\delta =M\overline{B_\delta }`$ and define the potential

$$\kappa _\delta (x)=\{\begin{array}{cc}\frac{\kappa _0}{\left|B_\delta \right|},& \text{ if }xB_\delta (x_0),\hfill \\ 0,& \text{ if }x\mathrm{\Omega }_\delta (x_0)\hfill \end{array}$$

(As before, this is discontinuous but can be approximated by continuous functions without changing the results.) By subtracting $`F(0)`$ on both sides of the equation for the eigenvalues, we can, without loss of generality, take $`F(0)`$ to be zero. We are thus lead to the functional

$$J_\delta (u)=_M|u|^2𝑑v_g+\alpha F\left(\frac{\kappa _0}{|B_\delta |}\right)_{B_\delta }u^2𝑑v_g.$$

Consider now the auxiliary eigenvalue problem defined by

$$\{\begin{array}{cc}\mathrm{\Delta }w=\mu w,& x\mathrm{\Omega }_\delta \hfill \\ w=0,& x\mathrm{\Omega }_\delta ,\hfill \end{array}$$

and denote its eigenvalues by $`0<\mu _0(\delta )<\mu _1(\delta )\mathrm{}`$, with corresponding normalized eigenfunctions $`v_{j\delta }`$. From the results in \[CF\] we have that

$$\underset{\delta 0^+}{lim}\mu _j(\delta )=\mu _j,j=0,\mathrm{}.$$

We now build test functions $`u_{j\delta }`$, $`j=0,\mathrm{}`$ defined by

$$u_{j\delta }(x)=\{\begin{array}{cc}v_{j,\delta }(x),& x\mathrm{\Omega }_\delta \hfill \\ 0,& xB_\delta ,\hfill \end{array}$$

for which

$$J_\delta (u_{j\delta })=_{\mathrm{\Omega }_\delta }|v_{j,\delta }|^2𝑑v_g,$$

and, by the result from \[CF\] mentioned above, this converges to $`\mu _j`$, $`j=0,\mathrm{},`$ as $`\delta `$ goes to zero. Finally, note that for each $`\delta `$ the set $`\left\{u_{j\delta }\right\}_{j=0}^{\mathrm{}}`$ satisfies the necessary orthogonality conditions, since this is the case for $`\left\{v_{j\delta }\right\}_{j=0}^{\mathrm{}}`$

A similar result will also hold in other cases, such as manifolds with boundary with Dirichlet or Neumann boundary conditions, for instance.

## 5. Concluding remarks

As was pointed out in \[EHL\] for the one–dimensional case, it is not difficult to see that for negative $`\alpha `$ the constant potential still maximizes the principal eigenvalue. It is also possible to show that in this case there is no lower bound on this eigenvalue, in the sense that there exist potentials $`\kappa `$ with fixed average $`\kappa _0`$ for which this eigenvalue can be made as large (in absolute value) as desired. It is not completely clear what happens to the supremum of the first eigenvalue for positive values of $`\alpha `$.

Regarding higher dimensions, it was shown that integral restrictions of this and similar type actually impose no restrictions at all as far as minimization is concerned, in the sense that it is possible to approximate the eigenvalues of the Laplacian as much as desired by potentials satisfying the given restrictions. Although we have seen that in this case the constant potential is never a global minimizer for positive $`\alpha `$, the results in Section 2 raise the question of whether or not it is a local minimizer for $`\alpha <\alpha ^{}`$.

We end by remarking that similar results to those in Sections 2 and 3 also hold in the case of manifolds with boundary and Neumann boundary conditions.

### Acknowledgements

This work was carried out while I was visiting the Department of Mathematics of the Royal Institute of Technology in Stockholm, Sweden. I would like to thank the people there and, in particular, Ari Laptev, for their hospitality.