text:
# Sub-mm Galaxies in Cosmological Simulations

## 1 INTRODUCTION

The observation of galaxies in the sub-millimeter has become a crucial area of astronomy since the installation of the SCUBA camera on the JCMT (the field is reviewed by Blain et al. 1999a). Because emission from dust in galaxies rises sharply up to a peak at $`100`$ µm, sub-mm observations preferentially detect high-redshift galaxies. While it was known from optical/near-IR observations that the cosmic star formation rate increases steeply with redshift up to at least $`z2`$ (Madau, Pozzetti & Dickinson, 1998), measurement of the extragalactic background with FIRAS (Fixsen et al., 1998) suggested that still more of the high-redshift star formation is hidden in the optical. The first SCUBA surveys (Smail, Ivison & Blain, 1997; Hughes et al., 1998; Barger et al., 1998) immediately showed that a large fraction of this hidden star formation takes place in extremely luminous systems, with star formation rates of $`100M_{\mathrm{}}\text{yr}^1`$.

Many authors have conjectured that the objects dominating the SCUBA surveys are elliptical galaxies or spheroidal components of disk galaxies in the process of formation. The analogy to low-redshift galaxies that are similarly bright in the far-infrared (FIR) suggests that they may be experiencing bursts of star formation due to major mergers. Although the surveys are increasing in depth and quantity, the low resolution, small fields, and large noise of the SCUBA maps compared to optical surveys leads to many systematic uncertainties that are still poorly quantified. Optical counterparts of the SCUBA sources are still difficult to identify, and there are very few detections in the FIR at wavelengths other than 850 µm. As a result, vigorous debates continue about the redshift distribution and star formation rates of the sources, the fraction of their energy supplied by AGN and star formation, their connection to other populations such as the Lyman break galaxies, and the contribution of these sources to the cosmic star formation rate (Lilly et al., 1999; Blain et al., 1999b; Adelberger & Steidel, 2000; Eales et al., 2000).

Comparison of these sub-millimeter observations to current theories of galaxy formation is clearly an important task. Modeling of the SCUBA population has so far been done in a phenomenological manner. For example, Blain et al. (1999b) and Blain et al. (1999c) showed there must be rapid luminosity evolution in the infrared galaxy luminosity function to reproduce the observed counts and sub-mm background. Guiderdoni et al. (1998) found that matching the counts within their semi-analytic model of galaxy formation required the ad hoc addition of a component of starbursting systems whose numbers increase rapidly with redshift. Artificial catalogs of sub-mm sources have also been constructed using Poisson-distributed sources (Eales et al., 2000) or large-scale dark matter simulations (Hughes & Gaztañaga, 2000), to aid in measurement of observational systematics. While useful, such catalogs do not bring one any closer to a physical understanding of these sources.

In this paper, we use numerical N-body/hydrodynamical simulations to study the sub-mm galaxy population. This paper follows several others devoted to the galaxy population in our simulations (Katz, Hernquist, & Weinberg, 1999; Weinberg et al., 1999; Davé et al., 1999b; Weinberg, Hernquist, & Katz, 2000). These papers showed that the optically detected Lyman-break galaxies at $`z3`$ are a natural consequence of CDM models, assuming reasonable amounts of dust extinction. They focused on the dependence of the galaxy abundance on the assumed cosmology and the clustering of these sources. This paper concentrates on the correspondence between the galaxies in our simulation and the bright sub-mm sources detectable with SCUBA. We restrict the study to a single cosmology, the currently favored $`\mathrm{\Lambda }`$-dominated cold dark matter (LCDM) model. We rely primarily on a single simulation of a large volume but only moderate mass resolution, which is well suited to a statistical study of the bright, rare sources. We defer a detailed discussion of fainter sources, the FIR background, and the relation between the sub-mm and UV/optical populations to a future paper using higher resolution, smaller volume simulations that are currently underway (Fardal et al., 2001).

We describe our numerical simulations in § 2, where we also discuss the resulting galaxy population and examine the effects of numerical resolution. In § 3, we describe our recipe for calculating the sub-mm fluxes. Our method ignores the detailed physics of dust absorption and emission, opting instead for simple, parameterized models motivated by observations. In § 4, we constrain these parameters by requiring the simulated galaxy population to match the counts at $`3\text{mJy}`$. We compare the resulting flux and redshift distributions to observations, arguing that the sub-mm population can be reproduced in our simulations for plausible assumptions about dust emission. In § 5, we examine the physical properties of the simulated galaxies that correspond to SCUBA sources. Section 6 discusses our confrontation of the simulations with observations and presents our conclusions.

## 2 SIMULATING THE GALAXY POPULATION

Our simulations are performed with the N-body/hydro code Parallel TreeSPH. This code follows both the gas and the dark matter with discrete particles, using spline kernel interpolation to smooth the gas quantities over a compact set of particles. The code is described elsewhere in detail (Hernquist & Katz, 1989; Katz, Weinberg, & Hernquist, 1996; Davé, Dubinski, & Hernquist, 1997). Simulations with this code have been successful in reproducing the observed clustering pattern and evolutionary trends of galaxies (Katz, Hernquist, & Weinberg, 1999; Weinberg, Hernquist, & Katz, 2000), as well as phenomena such as the Lyman-$`\alpha `$ forest (Hernquist et al., 1996; Davé et al., 1999a), and the abundance of higher-column density systems, including Lyman-limit and damped Lyman-$`\alpha `$ absorbers (e.g. Gardner et al. 1997a,b; 2000).

Star formation in this simulation is determined heuristically using local physical properties of the gas. For star formation to take place within a gas particle, the gas must be cold ($`T<3\times 10^4\text{K}`$), dense ($`\rho _H>0.1\text{cm}^3`$ and overdensity $`\rho _g/\overline{\rho }_g>55.7`$), converging, and Jeans unstable. The star formation timescale is determined by the local cooling and dynamical timescales, and because the star formation rate is an increasing function of gas density, the star formation rate in a simulated galaxy is governed mainly by the rate at which gas condenses and cools onto the central object (see Katz, Weinberg, & Hernquist 1996 for details of the algorithm and Weinberg, Hernquist, & Katz 2000 for further discussion). We identify galaxies in the simulations at the sites of local baryonic density maxima, using the program SKID111http://www-hpcc.astro.washington.edu/TSEGA/tools/skid.html. SKID finds galaxies by first determining the smoothed density field, then moving particles upward along the gradient of the density field. We define the galaxy to be the set of particles that aggregate at a particular density peak after removing particles that do not satisfy a negative energy binding criterion. To be included in a galaxy, gas particles must also have a temperature $`T<3\times 10^4\text{K}`$ and and a density $`\rho >10^3\mathrm{\Omega }_b`$, i.e. we consider only the cold gas and stars, the material that comprises the bulk of the matter in the central, visible regions of galaxies. Adding up the star formation rates for the individual particles in the galaxy gives the total star formation rate for that galaxy. Emission from active galactic nuclei or quasars is not included in the simulation.

Properties of the simulated galaxies are described in Weinberg et al. (1999) and Weinberg, Hernquist, & Katz (2000). Galaxies in the simulations form with a 2-phase gas structure: a central glob of gas at $`10^4\text{K}`$ that forms stars, surrounded by hotter gas at $`10^6\text{K}`$. Galaxies with luminosities comparable to the Lyman break galaxies form in abundance even at $`z>5`$. The clustering of these high-redshift galaxies is highly biased relative to the mass, in accordance with observations (Katz, Hernquist, & Weinberg, 1999). The star formation rates in the simulated galaxies are relatively steady on timescales of 200 Myr.

Comparing simulations of different resolutions and volumes, we find that there is a regime in which the galaxy properties are robust to changing these numerical parameters. The criterion for the low-mass end of this regime is that the galaxy should have $`64`$ cool gas or star particles. At the high-mass end, the volume should be large enough to allow halos of a given size to form in significant numbers. This latter criterion is less easy to quantify (as it is affected not only by Poisson statistics but also by the exclusion of long-wavelength modes in finite boxes), but it can be seen empirically as an upper cutoff in the mass function that depends on the box size. The star formation rates of the galaxies are strongly correlated with the galaxy masses, especially at high redshift. Hence, the star formation rate also has a characteristic resolution limit, though it is not as sharply defined as the mass limit. Because of the restricted range of halo masses, any one simulation represents only a portion of the total star formation.

We refer to the main simulation used in this paper as the L50/144 simulation. This simulation is tuned to study the evolution of large ($`LL_{}/4`$) galaxies down to $`z=0`$. It has $`144^3`$ particles in each of the gas and dark matter components, and models a periodic cube that is $`50h^1`$ comoving Mpc on a side. This gives us a baryonic mass resolution of $`5.4\times 10^{10}M_{\mathrm{}}`$ (based on the 64-particle limit above), and a spatial resolution of $`7.0h^1`$ comoving kpc (equivalent Plummer softening). We assume a $`\mathrm{\Lambda }`$ dominated cold dark matter cosmological model with $`\mathrm{\Omega }_m=0.4`$, $`\mathrm{\Omega }_\mathrm{\Lambda }=0.6`$, $`hH_0/(100\text{km}\text{s}^1\text{Mpc}^1)=0.65`$ and a primeval spectral index $`n=0.93`$. With the tensor mode contribution, normalizing to COBE using CMBFAST (Seljak & Zaldarriaga, 1996; Zaldarriaga, Seljak, & Bertschinger, 1998) implies a normalization $`\sigma _8=0.8`$, which provides a good match to cluster abundances (White, Efstathiou, & Frenk, 1993). We use the Hu & Sugiyama (1996) formulation of the transfer function. We adopt a baryonic density $`\mathrm{\Omega }_b=0.02h^2`$ consistent with the deuterium abundance in high-redshift Lyman-limit systems (Burles & Tytler, 1997, 1998), and the opacity of the high-redshift Lyman-$`\alpha `$ forest (e.g. Rauch et al. 1997).

For describing the faintest sub-mm sources and for computing the FIR background, the L50/144 simulation has insufficient resolution. We are conducting a suite of simulations with the same cosmology as above, but with different resolutions and box sizes to test the sensitivity of our results to these parameters. In this paper, we supplement the results of the L50/144 simulation with those of two additional simulations. L11/64 uses $`64^3`$ particles in an $`11.1h^1\text{Mpc}`$ volume evolved to $`z=0`$; this simulation has a mass resolution of $`6.8\times 10^9M_{\mathrm{}}`$ and a spatial resolution of $`3.5h^1`$ comoving kpc. The other simulation, L11/128, uses $`128^3`$ particles in an $`11.1h^1\text{kpc}`$ volume evolved to $`z=3`$; this simulation has a mass resolution of $`8.5\times 10^8M_{\mathrm{}}`$ and a spatial resolution of $`1.75h^1`$ comoving Mpc. When combined, these simulations give consistent results for the number density of galaxies with different star formation rates at high redshift (see Weinberg et al. 1999).

Although we will focus mainly on the bright sources that should be adequately represented in the L50/144 simulation, for some purposes we want to consider the contribution of lower mass systems. To do so, we fit the star formation rate (SFR) distribution across the L144/50, L64/11, and (at $`z>3`$) L128/11 simulations, using a Schechter (1976) function, which turns out to be a reasonable fit:

$`{\displaystyle \frac{dn}{d(\text{SFR})}}`$ $`=`$ $`{\displaystyle \frac{\varphi _{}}{\text{SFR}_{}}}(\text{SFR}/\text{SFR}_{})^{\alpha _s}e^{\text{SFR}/\text{SFR}_{}},`$
$`\text{SFR}_{}`$ $`=`$ $`\text{SFR}_{}^0(1+z)^{\delta _s+\gamma _s\mathrm{ln}(1+z)},`$ (1)
with $`\varphi _{}=6.49\times 10^3h^3\text{Mpc}^3\text{dex}^1`$, $`\alpha _s=1.96`$, $`\text{SFR}_{}^0=33.9M_{\mathrm{}}\text{yr}^1`$, $`\delta _s=2.63`$, and $`\gamma _s=1.04`$. This fit is in good agreement with the L50/144 simulation in its valid range of star formation rates. We choose a lower cutoff on the SFR to correspond to the completeness limit of the L11/64 simulation, which roughly equates to $`\mathrm{𝑆𝐹𝑅}_{\mathrm{𝑚𝑖𝑛}}=1.4(1+z)^{1.6}M_{\mathrm{}}\text{yr}^1`$. We then generate an artificial high-resolution “galaxy” population from these fits. Naturally, this artificial sample can be used only as an indicator of the results expected from a higher-resolution, large-volume simulation. It also cannot be used to study the spatial distribution or merger history of the galaxies. Further discussion of our fitting methods and results with a fuller set of simulations will be presented in Fardal et al. (2001).

The star formation rates in the simulations are the foundation of this paper, so we review here what we already know or suspect about the validity of these rates. For the cosmological model considered here, the simulations reproduce the observed density of Lyman-break galaxies given a reasonable amount of extinction (Weinberg, Hernquist, & Katz, 2000; Davé et al., 1999b). The cosmic star formation rate rises from $`z=0`$ back to $`z2`$–3 and then falls off, though we are still testing the sensitivity of the high-redshift decline to resolution (Weinberg et al., 1999). Comparison to low-redshift observations in the UV (Sullivan et al., 2000) and FIR (Saunders et al., 1990) suggests that the simulated star formation rates at low redshift may be too high by a factor $`2`$, although the observed functions and conversions from luminosity to SFR are subject to much uncertainty. The total amount of stars formed by $`z=0`$ may also be too high by a factor $`3`$ (Weinberg et al., 1999). Tests in Katz, Weinberg, & Hernquist (1996) show that the star formation rates are robust to changes in the algorithm parameters, being determined by gas supply rather than these parameters, although an entirely different algorithm could no doubt change the rates.

In sum, there are various hints that the overall amount of star formation in the simulations may be too high, although the size of the discrepency depends heavily upon the assumed IMF. This discrepancy could indicate observational underestimates of the true star formation or stellar density, numerical errors in the simulations, problems in the cosmological model (e.g., incorrect values of $`\mathrm{\Omega }_b`$ or $`\sigma _8`$), incorrect assumption about the form or constancy of the stellar initial mass function, or an incorrect treatment of the “microphysics” of star formation and feedback. A major goal of comparing simulations to a diverse set of observations is to gain guidance on this issue.

The SCUBA sources may be unusual galaxies, and even simulations that are adequate for explaining typical galaxies may be inadequate for studying rarer objects. This concern is particularly acute because the rarity of the SCUBA sources forces us to use a large simulation with low resolution. The true star formation rates may depend on details of galaxy structure that are unresolved in our simulation. For example, in simulated mergers involving two large spirals (Mihos & Hernquist, 1994a, 1996) or between spirals and less-massive companions (Mihos & Hernquist, 1994b; Hernquist & Mihos, 1995), the peak star formation rate depends strongly on the presence or absence of a bulge. In our simulations the low spatial resolution and large particle noise generally prevent the formation of galactic disks. However, if star formation in the SCUBA sources does not involve nuclear starbursts, then the resolution of the simulation might not be crucial. At any rate, it is not clear a priori whether or not we can accurately obtain the luminosities of the sub-mm galaxies, and the reader should keep these sources of uncertainty in mind.

## 3 EMISSION MODEL

Ideally, to model the emission from our galaxies we would take the distribution of stars and gas, compute the spatially dependent dust opacity by tracking the gas temperature and metallicity, and solve the radiative transfer problem through the distributed opacity. However, our simulations are far from adequate for this purpose: they do not resolve the internal structure of galaxies on the necessary scale (that of individual star-forming clouds), and they do not track the metallicity and dust content of the gas. We will therefore pursue a much simpler approach, working from the global star formation rate in each galaxy to a sub-mm flux using empirically motivated recipes.

The first step is to calculate the bolometric luminosity of each galaxy. In actively star-forming galaxies, the bolometric luminosity is dominated by young stars and is therefore roughly proportional to the star formation rate. We assume a Miller-Scalo (Miller & Scalo, 1979) initial mass function (IMF) extending from 0.1 to $`100M_{\mathrm{}}`$, consistent with the star formation treatment in our code. Using the code STARBURST99 (Leitherer et al., 1999), and assuming a burst length of $`10^8\text{yr}`$ and a metallicity of $`Z=0.02`$, we derive a conversion factor of

$$f_{SF}L_{bol}/(\mathrm{𝑆𝐹𝑅})=5.0\times 10^9L_{\mathrm{}}M_{\mathrm{}}^1\text{yr}.$$
(2)
We emphasize here that the calibration of the luminosity in the parameter $`f_{SF}`$ is quite poorly known. For one thing, there is a dependence on the burst age due to the buildup of intermediate-mass stars. Blain et al. (1999b) use $`f_{SF}=2.2\times 10^9L_{\mathrm{}}M_{\mathrm{}}^1\text{yr}`$ for the same IMF, presumably assuming a burst length of $`10^7\text{yr}`$. (The Miller-Scalo IMF is very sensitive to burst length below $`3\times 10^7\text{yr}`$, because of its steep slope above $`10M_{\mathrm{}}`$.) There is also a strong dependence on the shape of the IMF, particularly at low masses, since the mass is concentrated in low or medium-mass stars but the luminosity is concentrated in high-mass stars. A Salpeter IMF with the same mass range, burst age, and metallicity used in equation (2) gives $`6.5\times 10^9L_{\mathrm{}}M_{\mathrm{}}^1\text{yr}`$. The modification of the Salpeter IMF at the low-mass end suggested by Kroupa (2000) gives only 73% as much mass, or $`f_{SF}=8.9\times 10^9L_{\mathrm{}}M_{\mathrm{}}^1\text{yr}`$. Variations in the IMF have been much discussed, usually with the implication that more rapidly star-forming galaxies may have top-heavy IMFs and hence be more luminous (e.g., Doane & Mathews 1993). Variation in the metallicity can also affect the luminosity by at least several tens of percent.

We adopt the value in equation (2) as our standard conversion factor. However, the true value of $`f_{SF}`$ parameter could easily differ by a factor of two, and it could even vary systematically from one type of galaxy to another.

The next step is to determine the fraction of this luminosity that is absorbed and re-radiated by dust. The particular extinction model we choose is motivated by several observational facts. Extinction estimates for Lyman-break galaxies at $`z3`$, based on the UV spectral index, suggest large (mean of 5–8) and variable amounts of extinction (Adelberger & Steidel, 2000). We also require high average extinctions for our simulations to match the rest UV luminosity function for these galaxies (Weinberg, Hernquist, & Katz, 2000). Furthermore, the bright SCUBA sources tend to be extremely faint in direct optical/UV light (Chapman et al., 2000). All of these facts suggest that star-forming galaxies at high redshift are highly opaque, and that the most luminous ones are the most heavily absorbed, as is the case at low redshift.

The specific model we use assumes that the internal extinction of a galaxy depends upon its star formation rate. Let $`f_{\mathrm{𝑈𝑉}}`$ and $`f_{\mathrm{𝐹𝐼𝑅}}`$ be the fractions of the stellar energy that escape in the UV-optical continuum or are absorbed by dust, respectively. We assume that the average ratio FIR/UV ratio in a galaxy with star formation rate SFR is $`f_{\mathrm{𝐹𝐼𝑅}}/f_{\mathrm{𝑈𝑉}}=\mathrm{𝑆𝐹𝑅}/\mathrm{𝑆𝐹𝑅}_{\mathrm{𝑡ℎ𝑖𝑐𝑘}}`$, where $`\mathrm{𝑆𝐹𝑅}_{\mathrm{𝑡ℎ𝑖𝑐𝑘}}=4M_{\mathrm{}}\text{yr}^1`$ is the SFR at which 50% of the bolometric luminosity is absorbed. We also introduce a random scatter in this ratio of 0.5 dex per galaxy. This model makes the UV emission uncorrelated with the FIR emission at high bolometric luminosity, and it results in a Schechter-like cutoff in the UV luminosity function for Lyman break galaxies. With these opacities, the SFR cutoff in the high-resolution sample corresponds roughly to the population of Lyman-break galaxies seen at $`z3`$ for $`R_{AB}27`$.

The extinction model here is not unique, and we will examine this issue in more detail in Fardal et al. (2001). However, the essential point of the model for our purposes is that galaxies with large star formation rates are typically highly obscured; this is certainly true for the majority of the SCUBA sources. Our results in this paper thus depend only weakly on the details of the extinction, except perhaps for the estimates of the integrated background. Conversely, the predictions for the bright end of the UV/optical luminosity function are much more sensitive to extinction assumptions.

The next step is to assign an FIR spectral energy distribution (SED) to these galaxies. We adopt the usual emission model with just three parameters: a dust temperature $`T`$; an emissivity index $`\beta `$, such that the emissivity at low frequencies is $`f_\nu \nu ^{2+\beta }`$; and a short-wavelength slope $`\alpha `$, such that $`f_\nu `$ falls no faster than $`\nu ^\alpha `$. In luminous high-redshift galaxies, these parameters are of order $`T=40\text{K}`$, $`\beta =1.5`$, and $`\alpha 1.7`$ (e.g. Blain et al. 1999a). Increasing $`T`$, $`\beta `$, or $`\alpha `$ decreases the flux at 850 µm for galaxies at $`z5`$. These parameters provide reasonable, though not exact, fits to the SEDs of galaxies at low and high redshift (Dunne et al., 2000; Ivison et al., 1998). We frame our discussion of SEDs in terms of these parameters, not because they are truly fundamental parameters of the galactic emission (which in a more realistic model would come from dust at a variety of temperatures) but because they provide a simple and commonly used description.

The dependence of our results on the parameters $`T`$ and $`\beta `$ will be examined below. The short-wavelength slope $`\alpha `$ is unimportant for determining the shape of the spectrum near 850 µm, but it results in an overall shift in the normalization. Hence we keep it fixed at a value $`\alpha =3.0`$, comparable to that seen in Arp 220 (Hughes & Gaztañaga, 2000). For this choice, the quasi-thermal portion of the spectrum is shifted downwards by about 15%. The short-wavelength behavior of the spectrum may be the most variable part in the SEDs of low-redshift, FIR-bright galaxies (e.g. Figure 1 in Hughes & Gaztañaga 2000). Values of $`\alpha =2.2`$ or $`1.7`$, as assumed by Blain et al. (1999b) or Blain et al. (1999a), reduce the quasi-thermal spectrum by factors of 1.4 or 1.9 respectively (assuming $`T=40\text{K}`$, $`\beta =1.5`$).

The following fitting formula may be convenient for comparison with other work. With our choices for $`\alpha `$ and $`f_{SF}`$, the 850 µm flux from a completely opaque galaxy at $`z=2`$ is approximately given by

$$f_{850}=(0.44\text{mJy})\mathrm{exp}[3.7t(10.30t)]f_{\mathrm{𝐹𝐼𝑅}}(\text{SFR}/100M_{\mathrm{}}\text{yr}^1),$$
(3)
where $`t\mathrm{ln}(T/40\text{K})+0.38(\beta 1.5)`$. This formula is accurate to within 20% over the range $`25\text{K}<T<60\text{K}`$, $`0.3<\beta <2.5`$. Similar fits can be constructed at other redshifts, with slightly different coefficients. The effects of changing $`T`$ and $`\beta `$ on the flux are nearly degenerate.

For simplicity, most of our models will use a single SED for all galaxies. However, there seems to be substantial scatter in the observed SEDs at low redshift (Dunne et al., 2000; Adelberger & Steidel, 2000). Using a SCUBA survey of the IR-luminous galaxies in the local universe complemented with IRAS data, Dunne et al. (2000) fitted $`T`$ and $`\beta `$ for each galaxy, finding they were drawn from populations with means and standard deviations of $`T=35.6\pm 4.9\text{K}`$ and $`\beta =1.3\pm 0.2`$. Below we will test the effect of this scatter on our predictions.

Our empirical emission models are heavily based on the properties of IR-luminous galaxies at low redshift. In high-redshift galaxies, the amount of dust will probably be lower due to lower metallicity, while the FIR luminosity will typically be higher; these effects can be expected to raise the temperature. On the other hand, the FIR spectrum is also sensitive to the dust properties through the value of $`\beta `$, in a way that is not quantitatively understood. It is thus difficult to say whether the conclusions drawn at low redshift will also apply at higher redshift.

## 4 COMPARISON TO SUB-MM OBSERVABLES

To compute predictions of sub-mm counts, we take the list of star formation rates by particle group (“galaxy”) at each output redshift. The 27 outputs used from the L50/144 simulation are at redshifts $`z=0`$, 0.125, 0.25… 1.0, 1.25, 1.5…, 4, 4.5, 5… 6, and 7. Conceptually, we take the simulated volume and replicate it over shells centered on these redshifts. Each galaxy enters the observed distribution in an amount proportional to the volume of its redshift shell. The galaxies in any output are partly descended from those in previous outputs, of course, so the different outputs are correlated. We use the parametric models in § 3 to construct catalogs of observed flux, and thus distributions of the sources in flux and redshift. We construct similar catalogs based on the analytic fit to the high-resolution results (eq. 1).

We also construct simulated sub-mm maps based on our L50/144 simulation. Figure 1a shows a map at 1.2mm using the angular resolution of the large millimeter telescope (LMT) currently under construction. Figure 1b is at 850 µm and uses the angular resolution of SCUBA on JCMT. Neither map includes any observational or instrumental noise. These maps use a randomization procedure much like that employed by Croft et al. (2000). Since the sources all have angular extents much smaller than the $`15\mathrm{}`$ SCUBA beam, we treat each galaxy as a point source. The simulation box is replicated many times in comoving space along the line of sight, using the redshift output closest to the required redshift. Each box is given a random offset in the $`x`$, $`y`$, and $`z`$ directions, using the periodic nature of the simulation. The box is given a random parity inversion along each axis, and an axis is chosen randomly to lie along the line of sight. The box is then spun by a random amount around the line of sight. The resulting field size is 20 arcmin on a side, since we include outputs to $`z=7`$ and lose some area due to the spinning part of the randomization algorithm. As discussed below, these maps are incomplete with respect to sources below 1 mJy compared to a higher-resolution simulation.

Figures 1 give a visual sense of how the current SCUBA observations relate to the future observations that will be possible with the LMT. We also use these simulated maps in our analysis of clustering of sub-mm sources discussed in §5 below.

### 4.1 Distribution in flux

Using the catalogs of fluxes from each galaxy in the simulation, we calculate the counts on the sky as a function of flux at 850 µm. The counts that result from using only the L50/144 simulation, and single-temperature models, are shown in Figure 2. This plot demonstrates the large effect that the uncertainty in the SED can have on the predicted counts. The emission models used span most of the range usually considered for the sub-mm sources. Raising either $`T`$ or $`\beta `$ depresses the calculated fluxes by reducing the 850 µm luminosity for a given bolometric luminosity. Also, $`T`$ and $`\beta `$ are nearly degenerate in their effect on the counts, as one can see, for example, by comparing the curves for $`T=40\text{K}`$, $`\beta =1.0`$ and $`T=35\text{K}`$, $`\beta =1.3`$. As discussed above, the short-wavelength slope $`\alpha `$ could shift the curves leftward by tens of percent. The uncertainty in $`f_{SF}`$ could also shift the curves horizontally by as much as a factor of two (in either direction).

Some of the results from SCUBA surveys of galaxies are presented on the plot, with their source marked in the legend. The many different surveys and the large scatter they exhibit deserve some explanation. The solid symbols show the results of blank-field surveys (Barger, Cowie & Sanders 1999; Eales et al. 2000; Borys et al. 2001; Dunlop 2001). The lensing cluster surveys (Blain, Ivison, Kneib & Smail 1999; Chapman et al. 2000a) have two advantages stemming from the magnifying effect of the cluster: they can go deeper in flux, and source confusion is reduced. However, the necessary correction for the cluster lensing amplification introduces some risk of systematic error.

Since this is a cumulative plot, the data points are correlated in any given survey. Also, the error bars generally account for Poisson errors only and thus may be underestimated.

The raw counts are strongly affected by observational biases, because of the poor spatial resolution and S/N currently available. Eales et al. (2000) argue that the combined effect of noise and confusion amplifies the measured source fluxes by about a factor of 1.4 on average in their survey, independent of flux. We have thus shown in gray their survey, and the similar survey of Barger, Cowie & Sanders (1999), corrected down by this factor. However, this argument is controversial, with other workers disagreeing on the constancy (Blain, 2000) or even the sign (Hughes & Gaztañaga, 2000) of the effect. It is also unclear whether the cluster surveys, or the higher flux limit, blank-field surveys (Borys et al., 2001; Dunlop, 2001), are affected by confusion; if not, the disagreement in the measured counts is even greater at the high-flux end.

The sub-mm observations are still clearly in their early stages, and the present results have a large scatter that is not all due to random statistical errors. If we take our results with $`T35\text{K}`$ and $`\beta =1.3`$, roughly the mean parameters found in the local SCUBA survey (Dunne et al., 2000), we appear to have a reasonable fit to the lower envelope of the observations defined in part by the confusion-corrected points. This model would produce an integrated 850 µm background of $`12.9\pm 4.7\text{mJy}\text{arcmin}^2`$, which agrees with the background measurement from FIRAS, $`(12.2\pm 4.7)\text{mJy}\text{arcmin}^2`$ (Fixsen et al., 1998). The background measurements from DIRBE (Hauser et al., 1998) exceed those from FIRAS by about 40% at shorter wavelengths, although the two measurements are formally consistent, so the true background may be slightly higher. This SED is at the “cold” end of the parameters usually considered for sub-mm sources, an issue we return to in § 4.3.

While this model is in agreement with the lowest of the observational points, it is hard to escape the impression in Figure 2 that the true counts are somewhat higher. (Whether this is so depends largely on the validity of the confusion correction introduced by Eales et al. 2000.) It is interesting to see what kinds of models are required for the simulations to match the raw blank-field or the cluster survey counts instead. In Figure 3 we show three attempts to accomplish this match. The first model simply adopts yet smaller SED parameters of $`T=35\text{K}`$, $`\beta =1.0`$. The second model uses a larger $`f_{SF}=8\times 10^9L_{\mathrm{}}M_{\mathrm{}}^1\text{yr}`$, with $`T=35\text{K}`$ and $`\beta =1.3`$, and produces similar results to the first. In both cases the simulated counts may be steeper than the observations.

The third and fourth altered models allow for random scatter in the galaxy SEDs, as discussed in § 3. In a real galaxy the factors determining dust emission are surely not “random”, but this model reflects our ignorance of what causes the variation in the FIR SEDs of local galaxies. We find that the counts become significantly shallower in this case, improving the agreement with the blank field surveys. A mean temperature of $`T=35\text{K}`$ gives a good fit to the uncorrected counts, while a temperature of $`T=38\text{K}`$ provides a better fit to the confusion-corrected points.

Because of its physical plausibility and its reasonable success in matching the uncorrected blank-field counts, we choose the $`T=35\text{K}`$, $`\beta =1.0`$ model as our fiducial model for the rest of the paper (e.g., for Figures 1 and 7–9). However, to zeroth order, choosing a different $`T`$ and $`\beta `$ combination affects the fluxes from all galaxies in the same way, although there are first order differences depending on redshift. Hence, any such combination that fits the uncorrected counts will contain essentially the same objects at any given flux or counts threshold, and thus imply almost the same physical properties. Also, if one accepts the shift in the fluxes suggested by Eales et al. (2000), the preferred SED parameters change as discussed above; but the implied physical properties for a given flux-limited sample should be almost the same if one simply lowers the flux thresholds by a factor of 1.4.

So far, the calculated counts use only the L50/144 simulation. However, the typical star formation rate at a given baryonic mass is higher at high redshift (when gas accretion rates are higher), so the minimum star formation rate that is resolved in any given simulation rises towards higher $`z`$. With our preferred flux models, this resolution threshold begins to affect the counts at the 2 mJy level by $`z3`$ in the L50/144 simulation. To examine the effect of resolution, we use the artificial high-resolution sample discussed in §2. Since this sample is generated from an analytic fit rather than an actual simulation, it should be viewed as merely illustrative of the effect of higher numerical resolution. However, its properties are well constrained by the behavior in the low and high resolution simulations.

Figure 4 illustrates the effect of this extended high-resolution sample on the simulated counts. For each model, we show the results obtained both with the extended sample and with the L50/144 simulation alone. Extending the sample makes a significant difference below 1 mJy. The models that best agree with the observed counts above 1 mJy no longer agree below that level, and their cumulative surface brightness now exceeds the FIRAS background. Since the background values obtained using only the L50/144 simulation are acceptable, a cutoff in the opacity at $`60M_{\mathrm{}}\text{yr}^1`$ could produce agreement, but this cutoff value is likely too high.

The predicted excess of faint sources seen in Figure 4 may be a genuine problem for our model. However, at low luminosities and $`z<3`$, our analytic fit to the SFR function is constrained mostly by the small L11/64 simulation; hence uncertainty in the fit could be the cause of the discrepancy. The background excess could also be influenced by our simplistic extinction model. In this paper we are mainly interested in the objects above 1–2 mJy, so we postpone a detailed examination of this issue to a future paper (Fardal et al., 2001).

The counts in the L50/144 simulation may be incomplete at the bright end as well, due to the finite volume and consequent exclusion of long-wavelength modes. This probably occurs at a counts level of $`0.01\mathrm{}^2`$, or a flux level of $`10\text{mJy}`$ with our fiducial model, though it is difficult to confirm this without running a simulation in a yet larger volume.

With the best-fitting emission models here, the star formation rates implied for the SCUBA galaxies are large, though smaller than some estimates. For example, at $`z=2.5`$, a galaxy with $`\mathrm{𝑆𝐹𝑅}=100M_{\mathrm{}}\text{yr}^1`$ has a flux of $`1.2\text{mJy}`$ at 850 µm, assuming our fiducial model of $`T=35\text{K}`$ and $`\beta =1.0`$. For comparison, the flux would be only $`0.16\text{mJy}`$ using $`T=40\text{K}`$, $`\beta =1.5`$, $`\alpha =2.2`$, and $`f_{SF}=2.2\times 10^9L_{\mathrm{}}M_{\mathrm{}}^1\text{yr}`$ as in Blain et al. (1999b).

### 4.2 Distribution in redshift

The cumulative redshift distribution of the sources is shown in Figure 5. Here we assume a flux threshold $`f(850\mathrm{\mu m})>4.5\text{mJy}`$, comparable to the observational samples discussed below. The distribution is actually bimodal, as shown by the subtle minimum of the slope at $`z0.5`$. The broad high-redshift portion extends to $`z4`$. Using our fiducial SED model and the artificial high-resolution sample, the median redshift for the 4.5 mJy sample is 2.0; the mean and standard deviation are 2.0 and 1.3. Using only the L50/144 simulation lowers the redshifts, with a median of 1.8, mean of 1.8, and standard deviation of 0.9; the shifts are small because of the high sample threshold. Using a higher $`\beta `$ and $`f_{SF}`$ as indicated in the legend increases the median redshift to 2.4, with a mean of 2.3 and standard deviation of 1.2. If the Eales et al. (2000) confusion correction is applied to all fluxes in these redshift surveys, an appropriate model as discussed above is $`T=35\text{K}`$, $`\beta =1.3`$, and our standard $`f_{SF}`$; but in this case one should use a sample threshold of $`3.0\text{mJy}`$. Since lowering the flux threshold and raising the fluxes from all the objects are equivalent changes, the resulting distribution is nearly identical to the curve with the non-standard $`f_{SF}`$.

The predicted redshift distribution is not strongly affected by our choice of SED among the alternatives we prefer above. The main impact is to change the balance of the low-$`z`$ and high-$`z`$ components. The bimodal distribution arises in our simulation for the following reason. At high $`z`$, the flux from a source of a given luminosity is nearly constant, but the flux rapidly increases as $`z^2`$ for $`z1`$. The mean luminosity of our sources brightens as a power-law in $`(1+z)`$ at low $`z`$, so the characteristic flux reaches a minimum at $`z0.5`$. The low-$`z`$ component is boosted by the low values of $`\beta `$ in our preferred models. Changing the threshold flux also changes the balance of the low-$`z`$ and high-$`z`$ components. Above 10 mJy we find that about half of the sources are at $`z<0.5`$ for our fiducial model. However, at this flux level the counts in our simulations are very incomplete at high $`z`$, due to finite-volume effects. For a threshold of 2 mJy, we find instead a higher median redshift of 2.3 for our fiducial model.

Observational constraints on the redshift distribution of the sub-mm population are still fairly weak. Few sources have spectroscopic redshift measurements, so most of the redshifts are estimates based on the radio-FIR correlation. The derived redshift is dependent on the assumed SED, and is roughly proportional to $`T`$. In many cases the radio emission is undetected, giving only a lower limit to the redshift. The redshift estimates and limits from Smail et al. (2000), Barger, Cowie & Richards (2000), and Eales et al. (2000) all assume a set of SEDs to obtain the allowed redshift ranges. We have combined them here without regard to the slightly different SEDs employed.

Assuming a maximum redshift of $`z_{\mathrm{𝑚𝑎𝑥}}=5`$ for all sources, we have placed the observed sources at their lower redshift limit, their upper redshift limit, or at the middle of the range (the limits are taken from the observational papers and are $`2\sigma `$). The extreme assumptions give the shaded region in Figure 5, while the middle assumption gives the solid line bisecting this region. Most of the difference between these curves is systematic error due to the uncertainty in the SED. With the lower values of $`T`$ and $`\beta `$ that we prefer above, the redshift distribution would be close to the left hand boundary. The mean of the estimated redshifts (i.e., the middle curve, which assumes higher SED parameters than preferred here) is 2.9,222Attempting to incorporate the fact that about half of these sources have only lower redshift limits, Eales et al. (2000) estimate the mean redshift to be $`2.81\pm 0.36`$ from the same surveys, not significantly different from our value. and the median is 3.3. The dispersion in the estimated redshifts of the sources is 0.91. Of course, some of this dispersion could be due to scatter in the SED or the radio emission properties, or to observational error.

Our redshift distribution is in fairly good agreement with these results, with a similar mean (once the effect of the SED is accounted for) and width. We do appear to have an excess of predicted sources at $`z<1`$. This excess could be related to the large amount of low redshift star formation in our simulations. We are currently investigating this issue with a larger set of simulations (Fardal et al., 2001). The brightest observed sources could include sources that are primarily AGN, which we do not include in our modeling here; their inclusion would probably boost the predicted median redshift.

### 4.3 The SED of the sub-mm sources

The SED parameters of our preferred models are consistent with low-redshift observations of FIR-bright galaxies (Dunne et al., 2000). Somewhat higher temperatures have generally been preferred for SCUBA sources in the literature, however, which raises the question of whether our temperatures are too low to be plausible. Smail, Ivison & Blain (1997) preferred a temperature of $`60\text{K}`$ and described $`T=40\text{K}`$ as “very cold”. Other choices have included $`T=47\text{K}`$, $`\beta =1.0`$ (Barger et al., 1998); $`T=50\text{K}`$, $`\beta =1.5`$ (Hughes et al., 1998); and $`T=40\text{K}`$, $`\beta =1.5`$ (Blain et al., 1999a). Changing $`f_{SF}`$ by a reasonable amount would bring our best-fit temperatures into this range, as shown above.

As noted by Eales et al. (2000), the argument that the majority of SCUBA sources are ultraluminous infrared galaxies (ULIRGs) is somewhat circular. If the temperatures are assumed to be 40–60 K, typical of local ULIRGs, rather than 20 K, typical of local normal spirals, one finds that the sources have the bolometric luminosity of a ULIRG, but lower assumed temperatures would imply lower luminosities. Escaping the circularity requires redshifts and multiwavelength observations of sub-mm sources without evidence of AGN activity. There are only a few such sources at present (see discussion in Blain 1999). These do seem to have temperatures of roughly 40–50 K, though the limits are not stringent. In principle, an understanding of the dust mass implied by a given SED combined with an upper limit on the metals in a galaxy’s ISM can give lower limits on the temperature. The temperatures assumed in our best-fit models do not appear to violate these limits.

Some direct evidence regarding the SEDs of SCUBA sources comes from the 450 µm measurements of Eales et al. (2000), who derived a $`3\sigma `$ upper limit of $`S_{450}/S_{850}<2.4`$ for their 850 µm-selected SCUBA sample. Eales et al. argue that this limit implies surprisingly low dust temperatures ($`T25`$, if $`\beta =2.0`$) for the sub-mm sources. However, we disagree with this inference. First, if we accept their argument that their 850 µm fluxes are biased high by a factor 1.4 from observational confusion, this limit increases to 3.4. 333 Blain, Ivison, Kneib & Smail (1999) report detections of four galaxies at 450 µm, but they do not estimate the mean flux ratio for their sample. Since their derived surface density at 450 µm and 10 mJy is equal to that at 850 µm and 3 mJy, we may loosely infer a flux ratio of $`S_{450}/S_{850}=3.3`$, at least consistent with this limit. Second, Eales et al. (2000) assume a fixed redshift to argue that the temperatures are low. However, these redshifts are ultimately themselves based on an assumed SED, because they are almost entirely based on the radio-sub-mm flux ratio. As discussed by Blain (1999), this ratio depends essentially on the combination $`T/(1+z)`$. The 450/850 µm flux ratio depends on the same combination of parameters. Hence assuming cool dust temperatures alters the implied redshifts, and there is little net effect on the 450/850 µm ratio; the two parameters are degenerate with respect to these two flux ratios.

However, the effects of $`\beta `$ and $`(1+z)/T`$ on the 450/850 µm and radio/sub-mm ratios are not entirely degenerate; $`\beta `$ has a stronger influence on the former ratio than the latter. Hence the low 450 µm fluxes seem to suggest low values of $`\beta `$. Let us make a simple calculation, representing the entire sub-mm population by a single “typical” source. We calibrate the radio luminosity using the results for low-redshift IRAS sources of Yun, Reddy, & Condon (2001). Including their conversion factor of 1.5 between IRAS luminosity and total dust luminosity $`L_{\mathrm{𝑏𝑜𝑙}}^{(\mathrm{𝑑𝑢𝑠𝑡})}`$, the implied FIR-radio correlation is $`S_{1.4GHz}=3.1\times 10^{18}\text{erg}\text{s}^1\text{Hz}^1(L_{\mathrm{𝑏𝑜𝑙}}^{(\mathrm{𝑑𝑢𝑠𝑡})}/L_{\mathrm{}})`$. We also assume a radio continuum slope of $`0.8`$. We require that the 850 µm / 1.4 GHz ratio be 100, typical of the sub-mm source population, and the 450/850 µm ratio be $`<3.4`$. Then we find the limits $`\beta 1.0`$ and $`(1+z)/(T/40K)2.8`$. Admittedly, this calculation is crude and the input parameters are fairly uncertain. However, it does seem to support the low values of $`\beta `$ used in our preferred SEDs.

Unfortunately, further multiwavelength observations of the FIR continuum alone will be just as inadequate in disentangling the temperature and redshift. The key to determining the source SEDs is to obtain the source redshifts by some other means: from optical spectroscopy or photometry, mid-IR spectral features, or molecular lines. At this point, we know very little about the correct choice of SED, which is so important when comparing the simulations to observations. All we can say is that the observations may marginally favor the lower $`\beta `$ values that we have suggested, and do not contradict the lower values of $`T`$.

To summarize our results in this section, there is a galaxy population in our simulations that corresponds to the observed sub-mm sources, but only if the dust temperatures are sufficiently low or the stellar energy output is sufficiently high. To zeroth order the agreement is acceptable, although the counts may be somewhat too steep and there may be too many low-redshift objects compared to current observations. We now go on to describe the physical properties of these galaxies.

## 5 PHYSICAL PROPERTIES OF SUB-MM GALAXIES

In general, our sub-mm sources are quite massive galaxies, as shown in Figure 6. In a flux-limited sample with $`S_{850}>1\text{mJy}`$, the mass in cold gas and stars has a median of $`M_b1.3\times 10^{11}M_{\mathrm{}}`$ at $`z=2.5`$ (assuming $`T=35\text{K}`$ and $`\beta =1.0`$), and the median mass increases with increasing flux threshold ($`3.8\times 10^{11}M_{\mathrm{}}`$ for 3 mJy and $`4.4\times 10^{11}M_{\mathrm{}}`$ for 6 mJy). The distribution is quite broad, however. If we choose instead an SED with $`T=35\text{K}`$ and $`\beta =1.3`$, the star formation rate at a given flux increases by a factor of 1.4; since the star formation rate is roughly proportional to mass, the masses at a given flux threshold also increase approximately by this factor. The stellar mass distribution (not shown here) is yet broader, with a median stellar mass of $`6.6\times 10^{10}M_{\mathrm{}}`$ at 1 mJy, again increasing with flux.

We can define a timescale for star formation in these galaxies by taking this median mass over the median star formation rate of $`116M_{\mathrm{}}\text{yr}^1`$, giving a typical timescale of $`0.60\text{Gyr}`$. For reference, the cosmic time is $`2.4\text{Gyr}`$ at $`z=2.5`$. Thus the high star formation rates of our simulated galaxies are not just due to extremely short bursts. This conclusion is also implied by the correlation at high redshift between star formation rate and mass in our simulations, particularly for the high-mass objects (Weinberg, Hernquist, & Katz, 2000). If discrete bursts were the dominant cause of the high luminosities, they would tend to randomize the mass-luminosity correlation. We note that if the emission model contains random scatter as discussed above, the distributions of the flux-limited samples in Figure 6 become broader, but the general trends remain.

Another way to characterize the size of our sub-mm galaxies is to find their circular velocities $`V_c`$. However, if the galaxies live in groups or clusters, there could be an important distinction between the circular velocity of the galaxy itself and that of the largest virialized structure containing it. To tackle this question, we first use the program FOF444http://www-hpcc.astro.washington.edu/TSEGA/tools/fof.html to group dark matter particles into halos by the friends-of-friends algorithm, using a linking length appropriate to the mean interparticle spacing at the edge of a virialized halo (Kitayama & Suto, 1996). We center a sphere on each sub-mm galaxy and measure $`V_c`$ at a radius of 10 physical $`h^1\text{kpc}`$, the smallest reasonable radius given our spatial resolution of $`2h^1\text{kpc}`$ at this redshift. We also center the sphere on the most bound particle in the FOF halo to which the sub-mm galaxy belongs, and measure $`V_c`$ at the radius where the density matches the virial density.

The circular velocities defined at $`10h^1\text{kpc}`$ and at the halo virial radius are not greatly different, at least at $`z=2.5`$. The median $`V_c`$ of the 1, 3, and 6 mJy samples are 328, 500, and $`569\text{km}\text{s}^1`$ for the galaxies themselves at $`10h^1\text{kpc}`$. They are 282, 403, and $`495\text{km}\text{s}^1`$ for their surrounding virialized halos at the virial radius. It is not surprising that the galaxy circular velocities somewhat exceed those of their surrounding dark matter halos. Each sub-mm galaxy tends to be the dominant baryonic object in its neighborhood, with a mass large enough to increase the circular velocity at $`10h^1\text{kpc}`$. At these redshifts, large virialized structures with velocity dispersions characteristic of clusters (and higher than those of their central galaxies) have yet to form. The high circular velocities, substantially higher than those of $`L_{}`$ galaxies today, confirm that the sub-mm galaxies in this simulation are massive objects residing in massive dark matter halos.

We have computed the time evolution of the L50/144 simulation up to the present day, so we can directly measure the lifetimes of the sub-mm sources. To do this, we define various 850 µm flux-limited samples of simulated galaxies, assuming $`T=35\text{K}`$ and $`\beta =1.0`$. We then construct a merger tree of all of these galaxies down to $`z=0`$, using the same methods as in Murali et al. (2001). Restricting our study to the redshift interval $`1<z<4`$, which includes most of the SCUBA sources (c.f. Figure 5), we find the distribution of lifetimes for sub-mm sources as shown in Figure 7. Our time resolution is limited by the spacing of our discrete simulation outputs, which varies from $`0.1\text{Gyr}`$ at $`z=4`$ to $`0.7\text{Gyr}`$ at $`z=1`$. If a source does not repeat (i.e., appear in more than one output), we assume that it lasts for less than 0.1 Gyr. With higher time resolution, therefore, the time distribution curves would be higher at lifespans $`<0.7`$ Gyr, but the fraction of non-repeating sources is small, so this should not be a large effect. Of course, we are unable to tell whether the emission from a given galaxy turns off briefly between outputs; but if our sub-mm sources were constantly flickering on and off, this would only reduce the apparent lifetimes. The mean age of a 3 mJy source is $`500`$ Myr, although a significant fraction of our sources remain sub-mm sources for more than 1 Gyr. The contribution to the background is even more heavily weighted towards the long-duration sources than is apparent in Figure 7, as they tend to have larger luminosities. For example, considering the sources in the 1 mJy sample over the interval $`1<z<4`$, the sources that repeat at least one time have an average flux of 1.9 mJy, compared to 1.3 mJy for those that do not. Inspecting some of the nonrepeating sources, we find that they too are sources with slowly varying star formation rates, which are just cresting over the sample threshold. Similar conclusions apply at higher flux thresholds as well.

We can also follow the star formation in individual galaxies over cosmic time. While there is some rapid variation in the star formation rates of these galaxies, the amount of such “noise” is typically less than a factor of 2, and the variation in star formation rate is typically dominated by slow evolution scales of about a Gyr. (c.f. Weinberg, Hernquist, & Katz 2000). Our overall conclusion is that our simulated sub-mm sources are forming large numbers of stars in a fairly steady way at high redshifts.

Because they are so massive, the sub-mm galaxies in our simulation are strongly clustered. Over the redshift range $`1<z<4`$, their spatial correlation length is about 4–$`5h^1\text{Mpc}`$, with a correlation power-law exponent $`\gamma 2.0`$. For $`1<z<3`$, the correlation amplitude is enhanced relative to the average galaxy resolved in our simulation by a factor of $`1.7`$. However, this strong clustering will be difficult to detect observationally, because these sources are spread out over a wide range of redshifts. Direct measurement of the redshift-space correlation function requires an accurate redshift survey of the sub-mm sources, something that will not be achieved for some time.

We can predict the angular clustering of sub-mm sources using the simulated maps described at the beginning of §4. The angular correlation function of sources with fluxes above 1 mJy is approximated by $`w(\theta )0.02(\theta /\text{arcmin})^{0.9}`$ over a range 0.04 to 10 arcmin. This angular clustering is about a factor 3 lower than that of the Lyman break galaxies (Giavalisco et al., 1998), which is accurately reproduced in our simulations (Katz, Hernquist, & Weinberg, 1999). The reason for the weaker angular clustering of the sub-mm sources is simply that the stronger three-dimensional clustering is offset by a larger redshift depth. To quantify the angular clustering of the 1 mJy sources to 10%, a survey area of about 5 square degrees will be required. This is a much larger survey than available at present, but it could be easily achieved with the LMT.

The merger trees for our sub-mm sample allow us to follow the sub-mm sources in the L50/144 simulation to the present day. We define a sub-mm descendant here as a galaxy whose merger tree contains a sub-mm source in the redshift range $`1<z<4`$. Of course, these galaxies often contain other material as well. In some cases the particles from a sub-mm source form only a small part of the mass of a central cD-like cluster galaxy.

Naturally, we find that the sub-mm descendants are quite massive, since the sub-mm sources were already quite massive. We can examine the stellar population of these descendants as well. The simulation records each star formation event as it takes place, so we can find the star formation history of any given galaxy. Figure 8 plots the median star formation time for each galaxy in the simulation, versus the total baryonic mass of the galaxy. The sub-mm descendants are shown as the more prominent points. While the sub-mm descendants are usually massive, there are some that are among the less massive galaxies. However, they are almost universally among the oldest galaxies in the simulation, in terms of their stellar population. The vertical lines show the time range in which the middle 80% of the stars formed, for a few randomly chosen galaxies. One can see that the star formation is quite extended, agreeing with our depiction of steady star formation. In the real universe, the oldest galaxies are giant cluster ellipticals, which seem to have formed most of their stars before $`z2`$ (Bower, Lucey & Ellis, 1992; Stanford, Eisenhardt & Dickinson, 1995). Few galaxies in our simulation have a median star formation redshift as high as $`z=2`$, but we only have two galaxy clusters in our simulation and they are only the size of the Virgo cluster. Spheroidal populations, i.e. ellipticals and bulges, may have a smaller median age (Zepf, 1997; Fontana et al., 1999). Qualitatively, however, the large ages shown in the figure suggest an identification with elliptical galaxies.

The descendants of the sub-mm sources are also highly clustered. One way to measure this is by the circular velocities of their surrounding environments. We measure the circular velocities of the virialized dark matter halos in the same way as discussed above. Then, we define “cluster” and “group” halos by the requirements $`V_c>700\text{km}\text{s}^1`$ and $`350\text{km}\text{s}^1<V_c<700\text{km}\text{s}^1`$ respectively. With these definitions, 79% of all galaxies in our L50/144 simulation reside in the field, 18% in groups, and 3.2% in clusters at $`z=0`$. In contrast, 57% of our sub-mm source descendants reside in the field, 35% in groups, and 6.6% in clusters. Another way to show the strong clustering is to examine the fraction of galaxies that are sub-mm descendants as a function of local galaxy number density. We find a trend quite similar to the well-known density-morphology relation (e.g., Postman & Geller 1984), but even more extreme, with almost no descendants below number densities of $`0.1h^3\text{Mpc}^3`$.

It has been suggested that the sub-mm sources are massive ellipticals in the process of formation. One argument for this is that elliptical galaxies have old stellar populations, as shown (at least in clusters) by their tight color-magnitude relation (Bower, Lucey & Ellis, 1992) and the passive evolution of their colors (Stanford, Eisenhardt & Dickinson, 1995). This suggests, although it does not prove, that the galaxies themselves were formed at high redshift. Another argument is the common hypothesis that ellipticals were formed in mergers, coupled with the observation that ULIRGs at low redshift are typically mergers in progress. If the high-redshift sub-mm sources are similar to the low-redshift ULIRGs, this supports the notion that they mark the birth of ellipticals.

In the L50/144 simulation, our spatial resolution is too low to resolve the morphology of our galaxies and distinguish spirals from ellipticals. However, several features of the sub-mm source descendants—high mass, old stellar populations, and strong clustering—do indeed suggest that they correspond to massive ellipticals in the process of formation. Still, we emphasize that in our simulations the sources usually do not owe their high luminosities to very short bursts. Our sub-mm galaxies instead form in a more gradual manner over periods of $`0.1`$$`1\text{Gyr}`$.

## 6 DISCUSSION

The results presented in this paper amount to a physically motivated model for the population of sub-mm galaxies detected in the SCUBA surveys. The underlying basis of this model is the $`\mathrm{\Lambda }`$-dominated cold dark matter cosmological scenario, coupled to our numerical methods for following gravitational evolution, gas dynamics and cooling, and star formation. The essential features of the star formation algorithm are a gas consumption timescale that decreases steadily with increasing gas density, feedback effects that are relatively modest, at least in the massive galaxies that correspond to sub-mm sources, and an IMF similar to that inferred in the solar neighborhood (Miller & Scalo, 1979). Given our numerical simulation of the population of star-forming galaxies, the main uncertainty in our predictions comes from the choice of FIR spectral energy distribution. Matching the observed sub-mm counts requires values of the SED parameters $`T`$ and $`\beta `$ that are fairly close to those inferred for IR-luminous galaxies in the local universe by Dunne et al. (2000), but these values are lower than those used in most models of the SCUBA population. Our default model has $`T=35\text{K}`$, $`\beta =1.0`$, but the values of $`T`$ and $`\beta `$ have nearly degenerate effects, scatter in the SEDs would yield similar counts for higher $`T`$ and/or $`\beta `$ values, and the overall factor for converting star formation rate to bolometric luminosity is itself uncertain at the factor of two level.

The essential features of this model are that the sub-mm population is broadly distributed in redshift, with a median $`z22.5`$, and consists mainly of massive galaxies forming stars fairly steadily over timescales $`10^810^9`$ years at rates of $`100M_{\mathrm{}}\text{yr}^1`$. The descendants of these sub-mm sources are even more massive galaxies, with old stellar populations, found primarily in dense environments. At a qualitative level, these properties support the identification of sub-mm sources with the progenitors of luminous early-type galaxies.

The steady star formation and correspondingly moderate star formation rates are what distinguish this model of the sub-mm population from some of the alternative pictures proposed in the literature (discussed further below). While the typical formation timescales are still a factor of a few shorter than the cosmic time at the same redshift, they are much larger than, e.g., the timescale of $`2\times 10^6`$ yrs suggested for FIR emission by Thronson & Telesco (1986). Typical star formation rates of $`100M_{\mathrm{}}\text{yr}^1`$ require, in turn, relatively low values of $`T`$ and/or $`\beta `$ to yield sufficient flux at sub-mm wavelengths.

The principal successes of this model are its ability to reproduce the bright ($`>1\text{mJy}`$) sub-mm counts and its consistency with existing constraints on the source redshift distribution (Figure 5). Unfortunately, the first cannot be counted as a major success because of the large observational uncertainties and the substantial freedom to change the predicted fluxes by reasonable alternative SED choices. We did not make any adjustments to match the redshift distribution, so success there is reassuring, though the model does predict a potentially significant excess of low-redshift sources. The model also has difficulty matching source counts above $`10\text{mJy}`$, though these could be affected observationally by contributions from AGN and by lack of large-scale power in our finite simulation volume.

Another, potentially more serious failing of this model is its prediction, once we correct for numerical resolution effects, of too many faint sources and an excessively high 850 µm background (see Figure 4). We will save detailed discussion of this issue for a future paper (Fardal et al., 2001), but a few comments are in order here. Part of this problem may simply be due to the small volume of the L11/64 simulation we used to constrain the properties of smaller objects at low redshifts. Overproduction of faint counts and the sub-mm background, and the possible excess of low redshift sources mentioned above, could be connected to other possible shortcomings of the simulated galaxy population that we have noted elsewhere (e.g., Katz, Weinberg, & Hernquist 1996; Weinberg et al. 1999; Aguirre et al. 2001): an excess of low redshift star formation and an overly high fraction of baryons converted into stars. Balogh et al. (2001) have argued that the latter problem is generic to simulations of this sort and can only be solved by appealing to much more vigorous feedback effects from star formation (see Cole et al. 2000 for discussion of similar problems in semi-analytic models). While we agree that feedback could be a solution, there are many other possibilities, including: observational errors in the stellar density (see, e.g., the discussions of luminosity function discrepancies by Cole et al. 2001, Blanton et al. 2001, and Wright 2001); numerical errors causing the simulations to systematically overestimate gas cooling; changes to the cosmological model such as reduced initial power on galactic scales or lower baryon density; or an initial mass function (at least in some galaxies) that contains a large mass fraction of brown dwarfs or differs systematically from the solar neighborhood IMF in some other way. We note also that if feedback alone is to reconcile the simulations with the sub-mm source counts at $`1\text{mJy}`$, then it must be effective in galaxies with baryon masses $`5\times 10^{10}M_{}`$ (see Figure 6), not just in low mass systems.

An obvious risk is that any change to our model that reduces the faint source counts and sub-mm background will also spoil the agreement with bright source counts. A luminosity-dependent SED is possible, and perhaps even probable. A luminosity-dependent IMF could similarly improve the agreement. Both, however, represent unattractively ad hoc solutions at present. Another generic way to make the predicted source counts flatter would be to make star formation episodic, so that more of the bright sources correspond to more common, low mass objects caught during a burst of star formation. In particular, it is possible that our simulations underestimate the importance of merger-induced starbursts because of their limited resolution. In their simulations of galaxy mergers, Hernquist & Mihos (1995) and Mihos & Hernquist (1996) found bursts of star formation enhanced by a factor 10–100 over the star formation in isolated galaxies. The resolution in our cosmological simulation may simply be too poor to obtain such large bursts, or the large particle-induced noise may prevent the buildup of large gas reservoirs. Whereas each galaxy in the Mihos & Hernquist simulations was represented by about $`3\times 10^4`$ gas particles, our galaxies often contain only a few hundred, and they do not resolve galactic disks. Observationally, we know that low-redshift mergers produce enormous bursts of star formation (Sanders & Mirabel, 1996).

On the other hand, there are reasons to suspect that bursts and mergers do not play a dominant role in the global history of star formation, even if they are important in some objects. From this same L50/144 simulation, Murali et al. (2001) conclude that smooth accretion dominates over merging by at least 3:1 in mass during the assembly of galaxies. While the mass limit of the simulation is quite high, making it impossible to distinguish between low-mass mergers and accretion, they argue from the mass spectrum of the merging objects that true accretion is still dominant. The effect of individual mergers may also be less dramatic at high redshifts than at the present day. Cosmological timescales at high $`z`$ are small (e.g., the time from $`z=3`$ to $`z=2`$ is only 1.1 Gyr), and every merger takes some finite amount of time $`0.1`$–1 Gyr. Moreover, the merger rate of galaxies is larger at high redshift. Thus the bursts induced by different merger events may overlap for long stretches at high redshifts. From studying the star formation histories alone, it may be unclear whether the star formation is steady because mergers are unimportant, or because they are ubiquitous. In any case, higher-resolution simulations and detailed observations of high-redshift galaxies should clarify the role of mergers in producing high-luminosity sources.

There have been a number of models for the sub-millimeter sources in the literature, mostly phenomenological in nature. For example, Blain et al. (1999b) constructed sets of models based on pure luminosity evolution of the IRAS low-redshift luminosity function. They constrained the average SED and the low-redshift evolution by comparing the 60 µm and 175 µm counts from low-redshift sources. They then found several fitting functions for the typical luminosity, with “anvil” and “peak” shapes, that adequately represented the sub-mm counts. In a follow-up paper, Blain et al. (1999c) used a Press-Schechter formalism for merging objects to describe the sub-mm sources. To enable a fit to the counts, similarly peaked fitting functions were used for the amount of stars formed in an average merger. Unfortunately, the power spectra and dark halo masses used in this formalism appear somewhat disconnected from those expected within the CDM model, and the physical nature of the sub-mm sources in this model is somewhat obscure. Both of these papers put the sub-mm sources at somewhat higher redshifts (median $`z`$ of 2.5–4.5) than in our work, and derive much higher star formation rates (both in individual objects and in total). This is partly because of the different SED, and partly from the smaller assumed energy output from star formation.

Recently, Shu, Mao, & Mo (2001) have constructed a simpler analytic model for the sub-mm sources, in a sense performing the usual semi-analytical calculation backwards. They start from the observed distributions of the sizes and star formation rates of Lyman break galaxies, and use the Schmidt law for galactic disks to derive the gas masses, circular velocities, and star formation timescales. Within this model, sub-mm sources are simply defined to be Lyman break galaxies with $`\mathrm{𝑆𝐹𝑅}>1000M_{\mathrm{}}\text{yr}^1`$ (corresponding to the bright end of the observed SCUBA sources), and are assumed to have a top hat redshift distribution from $`2.5<z<3.5`$. The brightest sub-mm galaxies in their model are quite massive and highly correlated ($`r_07h^1\text{Mpc}`$), with even more highly clustered descendants corresponding to giant cluster ellipticals. The source timescales are a few tenths of a Hubble time, probably consistent with our results. Despite the crude approach, the basic picture they derive is fairly consistent with that drawn in this paper.

Taken at face value, our model makes several testable predictions. Probably the most robust — because it is directly connected to the long timescale and moderate rate of star formation — is that the FIR SED parameters of SCUBA sources should be ones that produce a relatively large 850 µm flux for a given star formation rate. Rest-frame UV and optical observations of sub-mm sources may also yield constraints on the star formation rates and timescales in these systems. Improved constraints on the redshift distribution of sub-mm sources can also test the model more stringently and give clearer guidance to the origin of possible discrepancies. Finally, our model predicts that the sub-mm galaxies are strongly clustered: their redshift-space clustering should exceed that of typical Lyman-break galaxies, though their angular clustering will be lower because of the large redshift range. Testing these predictions will require substantial improvements in the multi-wavelength and spectroscopic identification of sub-mm sources. These improvements will not come easily, but we have gone from knowing essentially nothing about the sub-mm source population to knowing quite a lot in the space of a few years, and we can expect similar improvements in the years to come.

We thank Simon Lilly, Amy Barger, David Hughes, Stephen Eales, Andrew Blain, and Scott Chapman for helpful discussions. We also thank Jeff Gardner for the use of his program SO. This work was supported by NASA Astrophysical Theory Grants NAG5-3922, NAG5-3820, NAG-4064, and NAG5-3111. By NASA Long-Term Space Astrophysics Grant NAG5-3525, and by the NSF under grants ASC93-18185, ACI96-19019, and AST-9802568. The simulations were performed at the San Diego Supercomputer Center, NCSA, and the NASA/Ames Research Center.


format: text

text:
# Preparation of the COROT space mission: fundamental stellar parameters determinations from photometric and spectroscopic analyses for COROT potential target candidates

## 1 Introduction

In order to prepare the target selection for the asteroseismology space missions MONS and COROT, a photometric and spectroscopic survey is planned of bright main target candidates as well as fainter targets to be observed. The data obtained from different sites will be used for an automatic classification to derive effective temperature (Teff), gravity (log g), and metallicity (\[Fe/H\]) of stars.

During the COROT-MONS Workshop, existing tools which must be applied for these tasks, and which should be adapted and further developed to give the maximal support to the space missions through groundbased support were presented. In this contribution, we only focus on the COROT (COnvection and ROTation) mission and the results we obtained from different methods: a detailed spectroscopic study, and 3 methods based on various photometric calibrations. One of this method (the BaSeL library) can be easily automatized for future applications for the COROT targets stars, so we put particular emphasis on this method in the context of the present workshop.

In this contribution, we first present the sample (Section 2) and the different methods (Sect. 3) used to assess the validity of the BaSeL library for the present purpose. In Sect. 4 we discuss the results derived from the BaSeL models, with comparison with the other methods presented in Sect. 3. Finally in Sect. 5 we present future development of the BaSeL models for the COROT target stars, and Sect. 6 draws the general conclusions.

## 2 The sample: 9 F-type stars

We use a sample of 9 nearby ($``$ 45 pc) F-type stars for which a detailed spectroscopic analysis has been performed. These stars are potential candidates of the COROT space experiment main programme (see for instance Catala et al. 1995, Michel et al. 1998, and Baglin et al. 1998 for details), and has been observed at the 193cm telescope at Observatoire de Haute-Provence (OHP, France) as part of the target selection process. High signal to noise ratio (S/N $``$ 150) spectra have been obtained using Elodie echelle spectrograph and analysed by comparison with theoretical spectra (Lignières et al. 1999). This set of new spectroscopic data together with both Johnson and Strömgren photometry from the literature provides a unique test of the predicted results of the BaSeL library. In addition, other calibration methods, namely the Templogg programme (based on Napiwotzki et al. 1993, and Künzli et al. 1997) and the Marsakov & Shevelev (1995) catalogue, provide further comparisons.

The photometric data in the Johnson (B$``$V and U$``$B) and Strömgren (b$``$y, m1$`=`$(v$``$b)$``$(b$``$y), c1$`=`$(u$``$v)$``$(v$``$b)) systems of our working sample of F-type stars are presented in Table 1, along with cross-identifications, and Hipparcos parallaxes ($`\pi `$). The issue of reddening is discussed in detail in Lastennet et al. (2001a).

## 3 Description of the methods

In this section, we briefly present the methods we used to derive Teff, log g and \[Fe/H\]: the BaSeL library, two photometric calibrations and a detailed spectroscopic analysis.

### 3.1 The BaSeL library

The Basel Stellar Library (BaSeL models) is constituted of the merging of various synthetic stellar spectra libraries, with the purpose of giving the most comprehensive coverage of stellar parameters. This is a library of theoretical spectra corrected to provide synthetic colours consistent with empirical colour-temperature calibrations at all wavelengths from the near-UV to the far-IR (see Lejeune et al. 1998 and references therein for a complete description). It has been corrected for systematic deviations detected in respect to UBVRIJHKLM photometry at solar metallicity, and can then be considered as the state-of-the-art knowledge of the broad band content of stellar spectra. These model spectra cover a large range of fundamental parameters (2000 $``$ Teff $``$ 50,000 K, $``$5 $``$ \[Fe/H\] $``$ 1 and $``$1.02 $``$ log g $``$ 5.5) and their photometric calibrations are regularly updated and extended to an even larger set of parameters (see e.g. Westera et al. 1999).
The BaSeL library spectra have been calibrated directly for standard dwarf and giant sequences at solar abundances and using UBVRIJHKLM broad-band photometry, and are hence expected to provide excellent results in these photometric bands (see for instance Lastennet et al., 2001b for a recent applications to AGB stars in near-infrared JHKLM photometry). Since they are based on synthetic spectra, they can in principle be used in many other photometric systems taken either individually or simultaneously, and this is another major advantage - used in this work - of these models.

Thus, because the BaSeL library has only been calibrated in the UBVRIJHKLM colours, the best parameters derivations come from this system. Whilst parameters derivations in equivalent bandwidth systems such as Washington seem to be as good with the BaSeL libraries than with empirical methods (Lejeune 1997), parameter derivations from the BaSeL library in narrow band systems such as Strömgren photometry should be of poorer quality. However, even in this non optimal case, good results were obtained for individual stars of detached eclipsing binaries (Lastennet et al. 1999).

In order to derive simultaneously the effective temperature (Teff), the metallicity (\[Fe/H\]), and the surface gravity (log g) of each star, we minimize a $`\chi ^2`$ -functional (see Lastennet et al. 2001a for all details).

#### 3.1.1 The ”Templogg” method

We have also run the ”Templogg” program which is designed to determine effective temperature and log g from either Strömgren or Geneva photometry. For Strömgren photometry it uses a Fortran program written by E. Fresno which relies upon the grids of Moon & Dworetsky (1985) in the Teff-log g parameter space relevant to this paper, with the improvements by Napiwotzki et al. (1993). For Geneva photometry it uses a Fortran program written by M. Kunzli (see North et al. 1994). The program chooses among eight different regions in the HR diagram for selecting the best calibration within the Strömgren system and three different regions for the Geneva system. The results from this method are gathered in Table 3.

#### 3.1.2 The catalogue of Marsakov & Shevelev (1995)

Marsakov & Shevelev (1995) (hereafter \[MS95\]) have computed effective temperatures and surface gravities using Moon’s (1985) method, which is also based on the interpolation of the grids presented in Moon & Dworetsky (1985). According to Moon (1985), the standard deviations of the derived parameters are Teff$`=`$$`\pm `$ 100 K and log g $`=`$ $`\pm `$0.06. The metallicities of Marsakov & Shevelev (1995) are obtained with the equation of Carlberg et al. (1985) which relates \[Fe/H\] with the colour excess $`\delta `$m1 and the Strömgren $`\beta `$ index. All the \[MS95\] results relevant for our sample are given in Tab. 3.

### 3.2 Spectroscopic analysis

The fundamental stellar parameters have been derived from a detailed analysis of spectra with high signal to noise ratio (S/N $``$ 150) obtained at OHP with the Elodie echelle spectrograph (spectra ranging from 3906 Å to 6811 Å at a resolution of $`\lambda /\mathrm{\Delta }\lambda `$ $`=`$ 42000). After reduction, the observed spectra were compared with theoretical ones constructed from a combination of Kurucz atmospheric models (ATLAS9 - Mixing Length Theory of convection with $`l/H_\mathrm{P}=0.5`$ and without overshooting111This choice of parameters different than those used by Kurucz is fully justified in van ’t Veer-Menneret & Mégessier (1996), the VALD-2 atomic database (Kupka et al. 1999), and the SYNTH radiative transfer codes (Piskunov 1992) and BALMER9 (Kurucz 1993)). In addition, the Least-Squares Deconvolution method (Donati et al. 1997) provided accurate determination of the projected rotational velocities listed in Tab. 2. While we confirm the low v sin$`i`$ values of HD 43587, HD 43318, HD 45067, HD 49933 and the high v sin$`i`$ of HD 46304 provided by the SIMBAD database, we find slightly largest values (by $``$25 to 40%) for HD 162917, HD 171834 and HD 171834. We also derive a new rotational velocity determination for HD 49434: v sin$`i`$ $`=`$ 79 km s-1.
Details about these determinations are given in the two next subsections and the results are summarized in Tables 2 and 3.

#### 3.2.1 Determination of Teff from the H$`\alpha `$ line

The effective temperature can be determined by taking advantage of the sensitivity of H$`\alpha `$ line wings. Detailed studies (e.g. van ’t Veer-Menneret & Mégessier 1996), have shown that the H$`\alpha `$ line is independent of the surface gravity (for non-supergiant stars) for effective temperatures ranging from 5000 K to $``$8500 K, and depends only slightly on the metallicity. Teff is therefore obtained for each star of the sample by fitting the observed H$`\alpha `$ line with synthetic spectra computed from a grid of solar metallicity atmospheric models separated by 250K.

#### 3.2.2 Determination of the surface gravity and metallicity

Within the temperature range that we found, we noticed that Fe I absorption lines depend only on temperature and metallicity, being practically independent of the surface gravity, while Fe II lines are sensitive to the temperature, metallicity and gravity. Consequently, the temperature being known from the H$`\alpha `$ line, the metallicity along with the microturbulence velocity are determined first by fitting a set of weak and strong Fe I lines. Then, the gravity is obtained by fitting Fe II lines. The spectral region near 6130 Å proved to be suitable for this analysis. However, the line broadening induced by rotation tends to mix neighbouring lines and prevent the analysis of individual Fe lines for high values of v sin$`i`$: HD 49434, HD 46304, HD 171834 and HD 164259 (all these stars have v sin$`i`$$``$60 km s-1).

### 3.3 Other determinations in the literature

To be as complete as possible, we looked for other determinations available in the literature and the SIMBAD database. One of the most comprehensive sources for our purpose is the fifth Edition of the catalogue of Cayrel de Strobel et al. (1997), which includes \[Fe/H\] determinations and atmospheric parameters ($`T_{\mathrm{eff}}`$, log g) obtained from high-resolution spectroscopic observations and detailed analyses, most of them carried out with the help of model atmospheres. However, the stars of our sample are not included in this catalogue. Since the catalogue comprises the literature (700 bibliographical references) up to December 1995, we only looked for more recent references. To the best of our knowledge, the catalogue of metallicities of Zakhozhaj & Shaparenko (1996) (hereafter \[ZS96\]) is the only one which contains useful information for our purpose. These metallicities are obtained from photometric UBV data and are available for two stars of our sample: HD 43587 and HD 164259 (see Table 3).

## 4 Results

Fig. 1, 2 and 3 show the Teff-\[Fe/H\]-log g results by using the following synthetic colour combination: B$``$V, U$``$B, and b$``$y. The results are also summarized in Table 3. Other photometric combinations are presented and discussed in detailed elsewhere (Lastennet et al. 2001a).
A general comparison with the other methods shows that the BaSeL solutions are very satisfactory for the 3 fundamental parameters: effective temperatures are in excellent agreement, and \[Fe/H\] and gravities show good agreement. The improvement using the B$``$V, U$``$B, b$``$y combination is clear when compared with the original results using the B$``$V, U$``$B, b$``$y, m1, c1 combination presented in Figs. 1 and 2 of Lastennet et al . (2001a) because there is no longer any systematic trend towards lower temperatures and because the determinations of log g and \[Fe/H\] are also much better, without systematic deviations.
There is only one exception (HD 46304): its Teff is in perfect agreement with its spectroscopic value, but the predicted log g and \[Fe/H\] are still low. While the BaSeL contours are consistent with the Templogg gravity only at the 2-$`\sigma `$ level, the metallicity predicted from the BaSeL models is poor in comparison to the result of the Templogg calibration ($`\mathrm{\Delta }`$\[Fe/H\]$``$0.7). What could explain this persistent difference ? It is worth noticing that this star has a large v sin$`i`$ (200 km s-1, the largest in our working sample), and it is well known that high rotational velocities modify the colours. Since the expected colour effect due to rotation is typically a few hundredths of a magnitude222These values are highly dependent of spectral type, age, and chemical composition, see for instance Maeder (1971) and Zorec (1992). in B$``$V and increases with v sin$`i`$, this is probably part of the reason why the predictions of the BaSeL models disagree with the Templogg method for this star. However this explanation is not satisfying because rotation is not taken into account in the Templogg method.
In conclusion, except in the case discussed before, reliable and simultaneous estimates of the 3 atmospheric parameters can be derived for F-type stars from only the three synthetic BaSeL colours, B$``$V, U$``$B and b$``$y. This is a very useful criterion for further applications.

A comment on the impressively good accuracy of the BaSeL determinations for the Johnson-Strömgren B$``$V, U$``$B, b$``$y combination (see Figs. 1, 2 and 3) is in order.
The Strömgren b$``$y index provides a reliable measure of the continuum and, therefore, a good temperature index, as discussed in detail in Lastennet et al. (2001a). This is particularly interesting, given the known difficulties of matching empirical and theoretical B$``$V-Teff scales to within better than about 0.03 mag (see, e.g., Sekiguchi & Fukugita, 2000)333 The b$``$y index is less vulnerable than B$``$V to the secondary effects of surface gravity and metallicity, even at the low spectral resolution given by the model spectra. . It is also well known that if one does not ask UBV data to provide the temperature in the first place (e.g., by using an independent source, such as spectroscopy, spectral classification, or else suitable other photometry, such as Strömgren b$``$y), the sensitivities of both U$``$B and B$``$V can be used to full advantage for determining log g and \[Fe/H\]. Moreover, although these sensitivities change with temperature, they are near or even at their maxima in the F-dwarf star domain (e.g., Buser & Kurucz 1992). This means that the derived values of \[Fe/H\] and log g presented in Figs. 1, 2 and 3 are as reliable as they can possibly be, given the uncertainties in the colours.

## 5 Future development of the BaSel models: application for COROT target stars

Usually, the BaSeL models provide colours in selected photometric systems, given a set of effective temperature, surface gravity and metallicity (Teff, log g, \[Fe/H\]). The new ”BaSeL interactive server”, the web version of the BaSeL models hosted by the Coimbra Observatory since the end of 2000 (http://tangerine.astro.mat.uc.pt/BaSeL/) already provides such information. We propose to develop in the near future a basic interactive tool based on the inverse method: an automatic tool based on the method of Lastennet et al. (1999) to complete the facilities of the ”BaSeL interactive server”. This method would be applied to the $``$ 1000 remaining potential targets of the COROT exploratory programme to provide their fundamental stellar parameters as discussed by Lastennet, Lejeune & Cuisinier (2001).

## 6 Conclusion

Several methods of determination of the fundamental stellar parameters Teff, log g and \[Fe/H\] are compared for nine single F stars COROT potential targets. Particular attention has been paid to the simultaneous predictions of the BaSeL models in two photometric systems, Johnson and Strömgren.

We presented the best combination to determine stellar parameters within the BaSeL library, using a combination of two Johnson colours, U$``$B and B$``$V, and a Strömgren colour, b$``$y.
This BaSeL combination is the best because on the one hand the b$``$y synthetic index gives reliable and accurate estimates of the effective temperature and, on the other hand, B$``$V and U$``$B give good estimates of \[Fe/H\] and the surface gravity.

We also note that the agreement between Templogg and BaSeL for the hottest stars of the sample could be especially useful in view of the well known difficulty of spectroscopic determinations for fast rotating stars.
As far as the astrophysical applications are concerned, the BaSeL synthetic colours are of particular interest in evolutionary synthesis and colour-magnitude diagram studies of stellar populations, such as open clusters and young associations, where (i) F-type stars are highly common and, of course, (ii) a vast abundance of data are available in B$``$V, U$``$B and b$``$y colours. Concerning the determination of the atmospheric parameters of the COROT potential targets, a result of the present analysis is that all the methods presented give consistent solutions. In the context of the further $``$1000 potential targets of the COROT exploratory programme, it will be of first interest to compare the results of the BaSeL models with those of other automated spectral analysis methods (e.g. Katz et al. 1998, and Bailer-Jones 2000).

### Acknowledgments

We thank J.C. Bouret, C. Catala and D. Katz for their participation in the spectral caracterization of the stars. This research has made use of the SIMBAD database operated at CDS, Strasbourg (France).