2 Model description

2.1 Input physics

We use the evolutionary tracks of the Padova group (Bressan et al. 1993; Fagotto et al. 1994a, 1994b, 1994c; Hereafter referred to as the Padova tracks ). The Padova group gives effective temperatures and luminosities as a function of time for many masses ($0.6 - 120\ M_\odot$) and a metallicity range from Z = 0.0001 to Z= 0.05. Their tracks include all stages in stellar evolution from the zero age main sequence to the tip of the RGB and from the zero age horizontal branch to the tip of the EAGB. For the stellar mass loss we use the method described in Bressan et al. (1993) with Reimers law.

For lower masses ($0.08\ M_\odot$ to $0.5\ M_\odot$) we use the calculations from Chabrier & Baraffe (1997). They based their calculation on a new description of the interiour of low mass objects and use non-grey atmospheres. Since their grid in metallicity does not match the grid of the Padova group we have linearly interpolated the values for [M/H] = -1.7, [M/H] = -0.7 and [M/H] = -0.4 while for [M/H] = -2.3 and [M/H] = 0.4 we used the calculations for [M/H] = -2.0 and [M/H] = 0, respectively. The contributions of these low mass stars to the integrated light is very low, therefore the error should be small.

To obtain synthetic colours, the conversion from theoretical quantities to observable quantities is very important. The evolutionary tracks for stars give effective temperatures, bolometric luminosities and gravities at the surface of stars of different masses as a function of time. These values have to be converted to colours in the various bands and to atmospheric indices of the evolutionary states.

For the fluxes and colours we use the theoretical library of model atmosphere spectra for various metallicities of Lejeune et al. (1997, 1998). Lejeune et al. have assembled a coherent library of synthetic stellar atmosphere calculations from Kurucz (see e.g. Kurucz 1979); Fluks et al. (1994); Bessell et al. (1989, 1991). The library has an effective temperature range from $T_{\rm eff} = 2000$ K to 50000 K and covers a broad range of metallicities. To cope with discrepancies between colours derived from model atmosphere spectra and observed colours, they correct (i.e. bend ) the model spectra to give agreement with observed colours for U through K. To fit their metallicity grid to that of the stellar evolutionary tracks, colours and bolometric corrections are interpolated linearly.

For a series of absorption indices, the empirical functions of Worthey et al. (1994) are used. Worthey et al. supply fitting functions that give index strength as a function of $T_{\rm eff}$, $\log g$, and ${\rm [Fe/H]}$. These were obtained from observations of 460 stars, covering a large range in the above parameters.

2.2 Model parameters

Once the input physics database is defined, the only free parameters in an evolutionary synthesis calculation for an SSP are those describing the initial mass function IMF.

The IMF can be expressed as:

$$\phi(m){\rm d}m \propto m^{-(1+x)}{\rm d}m$$ (1)

where the precise exponent x is observationally still somewhat controversal. For globular clusters, Chabrier & Méra (1997) find slopes between 0.5 and 1.5, independent of metallicity. In this paper we use the standard Salpeter IMF with a slope of 1.35. Variations in the slope within the mentioned range have only a small effect on the colours and indices. The low mass cut-off corresponds to the hydrogen-burning limit, which is dependent on the metallicity as shown in Chabrier & Baraffe (1997). It ranges from $0.083\ M_{\odot}$ for [M/H] = -2.0 to $0.075\ M_{\odot}$ for [M/H] = 0.

2.3 Numerical method

 The SSP models presented here are single metallicity single burst models where SF occurs in one timestep, i.e. during the first 10⁷ yr. The exact duration of this burst does not affect the properties of our SSPs at ages of a few Gyr.

Theoretical stellar tracks are supplied for discrete stellar masses only. In real stellar systems, the mass distribution is expected to be continuous. Using only the discrete mass grid of the track libraries would result in severe discontinuities, since all stars of a given mass would move to the red giant branch and die at the same time. We would see many "bumps'' in the luminosity and colour evolution of the stellar population. This effect is large for SSPs where all stars have about the same age. For continuous star formation rates (as in e.g. late spirals) the evolution is much smoother as expected.

To avoid the discontinuity problem we use a Monte-Carlo -method to calculate the distribution of stars in the HRD at each timestep of our evolutionary synthesis model. This method was developed by Loxen (1992, 1997). For this method, no isochrones with interpolated stellar evolutionary tracks are needed. Instead, at each timestep a large grid in stellar masses and ages is created. The grid is created randomly, hence the Monte-Carlo designation. Each cell in this grid has a size $\Delta m \cdot \Delta t$ in the 2-dimensional mass - time space. The grid ranges in mass from the lower mass to upper mass cutoff. In the case of an SSP the grid ranges in time from zero to the end of the burst, while for a continous SFR this would be from zero to the model age. Each cell represents a pseudo star, which is weighted with the value of the initial mass function (IMF) and the value of the star formation rate (SFR) at the position of the cell. This weight is given by  

$$w = \Delta m \cdot \phi(m) \cdot \Delta t \cdot \psi(t),$$ (2)

where $\phi(m)$ is the value of the IMF at the mass of the pseudo star and $\psi(t)$ is the star formation rate at the time the pseudo star is born. In the case of an SSP, the function $\psi(t)$ is constant over the time of the burst and then equal to zero. The unit of the value w is number of stars , although it can be fractional. For each of these cells, the lifetime of the pseudo star is determined by interpolation between the nearest two stellar tracks supplied in the library. Then the timesteps of the tracks are stretched according to this new lifetime. The cell is then split in two parts with relative weights given by

$$w_1 = \frac{\log(m_2) - \log(m)}{\log(m_2) - \log(m_1)} \cdot w; \quad w_2 = w - w_1,$$ (3)

where m is the mass value of the cell, m₁ is the mass of the track with the next lower mass and m₂ is the mass of the track with the next higher mass.

The only interpolation that is done is to determine the life times and hence the duration of the individual states. No interpolation is done for the luminosities or the effective temperatures which would require a precise definition of equivalent evolutionary stages. In this way, no artificial tracks are created. Especially in mass ranges were there is a strong dependance of the stellar evolution on the mass this is important.

The cells in the grid do not represent individual stars since their weights are not necessarily unity. The weights w are added in a book-keeping list for the individual states in the tracks of the input library to get weights W_i,j. Thereafter, the luminosities of all states on all tracks are summed up, each weighted with the calculated weights:

$$L_{\rm total} = \sum_{i,j} W_{i,j} \cdot L_{i,j},$$ (4)

where $L_{\rm total}$ is the luminosity in some band, W_i,j are the assigned weights and L_i,j are the luminosities of a theoretical star from the library at ith mass and jth state.

Due to the randomly spaced grid, there is a noise on the results, but the finer the grid (i.e. the more cells in the grid), the better is the signal-to-noise ratio. For the models in this paper, we use $200\thinspace000$ masses distributed evenly in $\log$ mass and 1000 ages distributed evenly, which make up a grid of $2 \ 10^8$ cells. Any further increase of the size of the grid does not give significantly better results. For more continuous SFRs, the grid can be coarser. Trivially the grid extends in time such that all stars are born in the burst interval, i.e. in the first 10⁷ years. This implies that each $\Delta t$ from Eq. (2) is equal to around 10⁴ years.

For the same input physics data set, the Monte Carlo method has been tested in detail against the standard synthesis model for SSPs (Fritze - v. Alvensleben 1995). Results for solar metallicity are compared to those of various evolutionary synthesis codes including the isochrone synthesis code of Bruzual & Charlot in Sect. 4.3.

On our Linux Pentium II 300 MHz work stations a single run with a timestep of one Gyr needs a few minutes of cpu time for the evolution of an SSP over a Hubble time.

This method has many advantages. Models at any time can be calculated without the need of any interpolation between stellar evolutionary tracks avoiding the problem of defining equivalent evolutionary stages, i.e. only the time intervals of the individual tracks are interpolated (but not the temperatures or gravities and therefore the colours). It is also possible to use an arbitrary star formation history. In the future, this code will be extended to be able to work with arbitrary metallicities not covered by the stellar tracks available. Since this is an evolutionary synthesis code rather than a population synthesis code, it can also be extended to calculate the chemical enrichment of galaxies.