text:
# The effect of a cusped matter distribution on the formation of brightest cluster members

## 1 Introduction

Brightest cluster members (BCMs) are giant ellipticals found in the central parts of clusters of galaxies. Some of these galaxies are further classified as cD galaxies because they show in their external parts a projected surface brightness excess over the standard $`r^{1/4}`$ law. Several models have been proposed for the origin of BCMs in general and of cDs in particular. A first class of models relies on the presence of cooling flows in clusters of galaxies. If the central cluster density is sufficiently high, intra-cluster gas can gradually condense and form stars at the bottom of the potential well (Cowie & Binney 1977; Fabian & Nulsen 1977). Such models, however, imply the formation of a lot of new stars, of which there is no good observational evidence. A second class of models relies on the dynamical mechanisms of merging, cannibalism and accretion (Ostriker and Tremaine 1975; Ostriker and Hausman 1977; Hausman and Ostriker 1978; McGlynn and Ostriker 1980 etc.). Schematically the central giant galaxy is formed by the merger of smaller galaxies and its growth continues by further mergings and cannibalism. It also accretes material stripped from the outer parts of the other galaxies, by interactions between galaxies, or between a galaxy and the common background. It has been proposed (Gallagher & Ostriker 1972; Richstone 1976) that such material could be responsible for the luminous envelope of cD galaxies. Merritt (1984) argued that most of the formation happens during cluster collapse.

Most of the mass of clusters or groups of galaxies is not bound to the galaxies that constitute them. Instead it is distributed in a common background encompassing the whole cluster and it includes the high temperature gas emitting X-rays and an unknown fraction of dark matter. This unknown mass distribution must have influenced the formation mechanism of BCMs. Nevertheless, considerable uncertainty still remains as to what fraction of the total cluster mass resides in the dark matter, and how it is distributed.

Observations of galaxy clusters can be used to constrain their mass distribution. The first attempts to obtain the profile of this mass distribution were based on galaxy counts and the projected velocity dispersion of the galaxy system. They rely on the assumption that the cluster has reached virial equilibrium and depend on what is assumed about the possible anisotropy of the system. They lead to mass distributions possessing central cores where the mass density is nearly constant (e.g. Kent & Gunn 1982). More recent studies, however, point to a mass distribution with a cusp in the central parts (Carlberg et al. 1997).

For some galaxy clusters which are very luminous in X-rays, many authors have studied the total mass distribution under the assumption that the hot gas is in hydrostatic equilibrium within the total potential. The correlation between the velocity dispersion of the galaxies and the X-ray temperature (e.g. Wu et al. 1998), as well as the observational evidence that the dynamical properties of clusters have evolved little since $`z0.8`$ (Mushotzky & Scharf 1997; Bahcall, Fan & Chen 1997; Henry 1997; Rosati et al. 1998; Vikhlinin et al. 1998, etc.) indicate that clusters are globally relaxed structures. Further supportive evidence comes from the fact that Crimele et al. (1997) find a relation between the gas and galaxy densities consistent with the assumption of hydrodynamical equilibrium. Thus the Jeans equation together with the hypothesis of hydrostatic equilibrium can be used to describe clusters globally, except for clusters which are still in the process of formation, or of merging with another cluster. In most cases X-ray mass estimates have been obtained with models with flat cores. More recently, however, Tamura et al. (2000) have shown that models of the total mass distribution with cusps give better fits to the data of the cluster A1060 than models with flat cores.

Yet a third way to derive the mass in a cluster is to use the gravitational lensing effect which is produced when the images of background galaxies are distorted by the mass concentration of a galaxy cluster. In the case of strong lensing background galaxies are distorted as giant blue arcs. A massive foreground cluster can also change the shapes (shear effect) and the number density (magnification) of the faint background population. These effects are known as weak lensing. Arcs have been used intensively to map the matter distribution in the central parts of clusters (Fort & Mellier 1994; Mellier 1999 and references therein). In particular considerable progress was achieved from recent spectacular images of arc(let)s obtained with the HST. Many different mass distributions, including profiles with cores or cusps, have been used in the modelling. There seems to be agreement that, if cores exist, their radius should be less than $`50h_{100}^1`$ Kpc (e.g. Mellier 1999 and references therein). Except for that, no clear consensus seem to have been reached on the form of the radial profile in the innermost parts.

Comparison of the results obtained with the various methods described so far lead to considerable controversy. In particular, results from strong lensing are a factor $`24`$ higher that results from X-rays. (e.g. Loeb & Mao 1994, Miralda-Escudé & Babul 1995, Tyson & Fischer 1995; Wu et al. 1998, Mellier 1999). This could be explained by the fact that, although the X-ray gas should be globally in gravitational equilibrium, locally in the center-most parts virial equilibrium may not be achieved. It should, however, be noted that models with smaller core radii, with cusps, or more generally with higher masses in the center-most parts, reduce the differences between the results obtained with the different methods (e.g. Wu et al. 1998, Makino & Asano 1999, Markevitch et al. 1999).

N-body numerical simulations of hierarchical gravitational collapse of an ensemble of cold collisionless particles lead to halo radial profiles with a cusp in the central parts. Navarro, Frenk and White (1996) argued that such profiles can be fitted by a “universal” profile of the form:

$$\rho (r)=\frac{\rho _s}{(r/r_s)(1+r/r_s)^2}$$
(1)
where $`\rho _s`$ and $`r_s`$ are the characteristic density and length respectively. This profile was first obtained in the context of the standard cold dark matter cosmology, but subsequent numerical studies have shown that this profile is independent of halo mass, of the initial density fluctuation and of cosmology (e.g. Cole & Lacey 1996; Navarro Frenk and White 1997; Eke, Navarro & Frenk 1998; Jing 1999). Simulations with higher numerical resolution result in profiles which are steeper in the central parts, i.e. of the form $`\rho r^{1.4}r^{1.5}`$ (Fukushige & Makino 1997, Moore et al. 1998), while cosmologies using collisional dark matter point to even steeper ($`\rho r^2`$) profiles (Yoshida et al. 2000; Moore et al. 2000), at least in the limit of short mean free paths.

For these considerations it seems appropriate to compare the evolution of clusters of galaxies within backgrounds with different mass profiles, in hope that this may lead to clues about which kind of mass distribution is more realistic. Several N-body simulations of the dynamical evolution of large groups or clusters of galaxies with each galaxy being represented by many particles have been reported in the literature. Funato, Makino & Ebisuzaki (1993), Sensui, Funato & Makino (1999) and Garijo, Athanassoula & García-Gómez (1997, hereafter GAG) use model clusters where initially all the mass is in the galaxies. Funato, Makino & Ebisuzaki (1993) argue that the Faber-Jackson relation is a result of the evolution of galaxies driven by interaction with other galaxies and with the tidal field of the parent cluster. Sensui, Funato & Makino (1999) follow the evolution of a cluster of 128 galaxies initially in virial equilibrium with all mass being attached to the galaxies. They find that within a few cluster crossing times half of the total mass has escaped individual galaxies, and that the amount of stripped material is larger for galaxies in the central regions. The density profile of the common halo thus formed is cuspy, with a radial dependence roughly $`\rho r^1`$ in the central region and $`\rho r^4`$ in the outer regions. Bode et al. (1994), GAG and Dubinski (1998) witnessed in their simulations the formation of a massive central galaxy. Bode et al. (1994) modelled poor clusters with 50 galaxies and a varying percentage of the mass in a common background. They find multiple nuclei in the central object in between 10 and 40 per cent of the time. Increasing the percentage of mass in the common background slows the merging rate and this can be sufficient for stalling the formation of the central object if 90 per cent of the mass is in a common background. GAG studied in detail the properties of the central object and found good agreement with those of brightest cluster members. The object formed is roughly spherical, oblate or mildly triaxial. In the latter cases the orientation of the central object correlates well with that of the central group. The triaxiality is in general stronger in the outer parts of the central object. They also find three types of projected density distributions. Objects in the first category have profiles well fitted by an $`r^{1/4}`$ law. In the second category of objects the radial profile follows the $`r^{1/4}`$ law only in the main part of the galaxy and falls bellow it in the outer parts. In the third, most interesting, category the projected surface density of the objects follows the $`r^{1/4}`$ law in the main body of the galaxy, but is systematically above it in the outer parts, as is the case for cD galaxies. The projected velocity dispersion profiles of the central objects also agree with those of brightest cluster galaxies. Dubinski (1998) used an N-body simulation of a cluster of galaxies in a hierarchical cosmological model. He replaced 100 dark matter halos at $`z`$ = 2 with self-consistent disc+bulge+halo galaxy models, thus attaining high resolution. He witnessed formation of the central galaxy through the merging of several galaxies along a filament. The central object was flattened and triaxial, showing alignment with the primordial filament. Its projected surface density was well described by an $`r^{1/4}`$ law and showed no excess in the outer parts, as one would expect from a cD. The kinematics of the object also was in good agreement with the observations.

In this paper we will study the influence of the dark matter distribution on the evolution of the cluster and on the formation and properties of the central galaxy. Our aim is to show that the observable properties of this galaxy can be used to constrain the dark matter distribution in the central parts. The paper is organized as follows. In section 2 we present the simulations and the initial conditions. In section 3 we describe the different evolution of the various models and in section 4 we show how the different initial conditions influence the formation process and the properties of the central giant galaxy. Section 5 discusses the main assumptions and simplifications entering in our simulations. Our results are discussed and summarised in Section 6.

## 2 Simulations

We have performed a total of $`17`$ simulations following the time evolution of a group of $`50`$ galaxies in virial equilibrium within a common halo. In nearly all cases we used a total of $`\mathrm{100\hspace{0.17em}000}`$ particles with the same mass. A fraction of these particles was used to represent the galaxies and the rest were used for the common halo. In one simulation only one third of the total mass was in galaxies and for this case we used $`\mathrm{200\hspace{0.17em}000}`$ particles, while maintaining the total mass constant, in order to have galaxies represented by a reasonable number of particles.

The common halo was represented alternatively as a Plummer or as a Hernquist sphere (Hernquist 1990). The density distribution of a Plummer sphere follows the law:

$$\rho _p(r)=\frac{3M_{ch}}{4\pi a_p^3}\left(\frac{a_p^2}{r^2+a_p^2}\right)^{5/2}.$$
(2)
This distribution is characterised by the parameters $`M_{ch}`$, which represents the total mass in the common halo, and $`a_p`$ which represents the scale-length of the distribution. The Hernquist distribution follows the law:

$$\rho _h(r)=\frac{M_{ch}}{2\pi }\frac{a_h}{r}\frac{1}{(r+a_h)^3}.$$
(3)
Two parameters are also necessary to determine this distribution: the total mass in the common halo $`M_{ch}`$ and the scale-length $`a_h`$. Note that while the Plummer distribution has a central density of $`3M_{ch}/4\pi a_p^3`$, the Hernquist distribution is singular at the origin. Its projected surface density follows closely the $`r^{1/4}`$ law for a large fraction of the radius (Hernquist 1990). It gives a cusped density distribution in the central parts which falls as $`\rho (r)r^1`$ for small radii and as $`\rho (r)r^4`$ at large radii.

For a given mass of the common halo $`M_{ch}`$, the parameters $`a_p`$ and $`a_h`$ were selected as to give a common halo with the same total energy. This gives the following relation between the two parameters:

$$a_h=\frac{16}{9\pi }a_p.$$
(4)
In this way we can compare simulations with the same global properties changing only the mass profile of the common halo.

The positions and velocities of the particles in the common halo are drawn according to the chosen Plummer (cf. eg. Aarseth, Henon & Wielen 1974) or Hernquist model (Hernquist 1990). The positions and velocities of the center of masses of the galaxies are drawn from the same distribution as the common halo, except that we impose a cut-off radius which is a fraction of the half mass radius $`R_{1/2}`$ of the common halo. Thus the galaxies are distributed within a radius selected as

$$R_{cl}=fR_{1/2},$$
(5)
where $`f`$ is a parameter which is used to determine the initial concentration of the galaxy system. Finally, we also vary the ratio between the mass in galaxies and the mass in the common halo. This mass ratio is given by the parameter:

$$r=\frac{N_gM_g}{M_{ch}},$$
(6)
where $`N_g`$ represents the initial number of galaxies and $`M_g`$ the initial mass of each galaxy. We use in all cases $`N_g=50`$ identical galaxies, whose density profiles follow a Plummer law with a scale-length of $`a_{gal}=0.2`$. The mass of each galaxy, $`M_g`$, will depend on the mass ratio selected for each particular initial condition, i.e. on the value of $`r`$.

We repeated one of the simulations using a rigid halo which is not modified during the whole simulation, instead of a live halo composed of particles. We did this in order to study the effect of dynamical friction and the possible differences introduced by a halo capable of responding to the variation of the mass distribution as the central giant galaxy is formed.

In Table 1 we give more detailed information about our simulations. In Column 1 we list the simulation name, coded conveniently so as to describe the initial conditions of each simulation. The first letter P or H is used to note the type of mass distribution, P for Plummer models and H for Hernquist models. This letter is followed by one or two digits indicating the scale-length of the common halo in simulation units. For the case of Hernquist models we indicate the scale-length of the Plummer model with the same total energy. Thus from the name of the simulations we can readily see which of the different simulations are directly comparable. After these digits comes the letter $`f`$ followed by the value of this parameter and finally the letter $`r`$ followed by the value of the mass ratio. In this work we have considered values of $`r`$ between $`1/2`$ and $`20`$, thus covering a range of more than one order of magnitude. The three values of the parameters determining the initial conditions are repeated for clarity in Columns 2, 3 and 4. In Column 2 for the case of the Hernquist simulations we list the value of the scale-length of the Hernquist distribution. In Column 5 we indicate the total number of particles used to represent each galaxy and finally in Column 6 we indicate if the common halo is rigid or a live one evolving freely. In two of the simulations, P10r1 and P10r20, the initial positions of the galaxies were selected anywhere within the common mass distribution without any cut-off. Finally in the Hernquist simulation with low galaxy-to-common-halo mass ratio, namely H10f1r1/2, the total number of particles was doubled, while the total mass was kept constant.

To evolve these systems we used the GRAPE-3AF system of the Marseille Observatory using a version of a treecode specially adapted for this system (Athanassoula et al. 1998). This version, since the tree is descended for a whole group of particles, rather than for every particle separately, is more precise than a standard treecode (Barnes & Hut 1986). For the opening angle we used $`\theta =0.7`$, for the time step $`\delta t=0.0078125`$ and for the softening parameter $`ϵ=0.0625`$. With these values the total energy was conserved to within a few parts in a thousand. The simulations with a Plummer matter distribution were followed for a total of roughly $`60`$ computer units or, equivalently, roughly $`7680`$ time steps. Each simulation required about $`24`$ hours of CPU time on the Marseille five board GRAPE-3AF system. The evolution of simulations with a Hernquist distribution was followed for twice as long, to follow best the evolution of the central object.

In this paper, we use the same units as in GAG, where the simulation units were selected as $`G=1`$, the total mass $`M_T=50`$ and the typical scale-length of a galaxy $`a_{gal}=0.2`$. With the aim of comparing with the observations we set $`M_T=1.5\times 10^{13}M_{}`$ and $`a_{gal}`$ = 6 Kpc. These values give a total of $`1.4\times 10^8`$ years for the time unit, $`210`$ Km/s for the velocity unit and a total simulation time of $`T=8.4`$ and 16.8 Gyr, for Plummer and Hernquist distributions respectively. Note that this selection of conversion factors, albeit reasonable, is not unique and neighbouring values could also be used. This should be taken into account when comparing our simulations with observational data, hence agreements to within a factor of two should be considered to be quite satisfactory.

## 3 Evolution of the simulations

### 3.1 A general impression from snapshots

As in the simulations without common halo (GAG) the most noticeable result of the dynamical evolution is the formation of a giant central galaxy in the central parts of the group. The rate of formation of this galaxy, as well as its final extension, depend on the properties of the common halo and of the galaxy system. This can be seen qualitatively in Figs. 1 to 4, where we show respectively the influence of the scale-length, of the initial extension of the galaxy distribution, of the ratio of the mass in galaxies to the mass in the common halo and the influence of different mass profiles.

Amongst our simulations with Plummer background the one that starts more centrally concentrated is P5f3/2r1 and the one that starts off with the least central concentration is P15f3/2r1. Since the other parameters are the same in the two simulations, the comparison of the two evolutions, given in Fig. 1, shows the effect of the scale-length of the galaxy distribution and of the common halo, which changes by a factor of three between the two cases. The more concentrated systems, with shorter scale-length, have faster evolution rates. In this way in simulation P5f3/2r1 a central galaxy is rapidly formed which controls thereafter the dynamical evolution of the system, merging with the satellite galaxies in a shorter time than in the rest of the simulations. On the other hand, in the simulation P15f3/2r1 no central giant galaxy is formed and practically all the galaxies survive until the end of the simulation.

In Fig. 2 we compare the evolution of simulations P10f3/2r1 and P10f3/4r1 to understand the effect of changing the initial cut-off radius of the galaxy distribution. Since the number of galaxies is initially the same in both cases a decrease of the cut-off radius entails a higher density of the galaxy system. As we can see, if the galaxy distribution is initially denser, the system has a faster evolution rate. Since this is also the effect of a more concentrated distribution, as discussed earlier, we can conclude that it is the density in the central parts that controls the rate of evolution. Thus distributions with higher central densities, coming either from higher central concentrations or higher densities over a larger area, will have a faster evolution than distributions with lower densities in the central area.

Fig. 3 compares the evolution of simulations P10f1r1 and P10f1r20 and shows the effect of the ratio of the mass in individual galaxies to the mass in the common halo. Simulation P10f1r1 has 20 times more mass in the common halo than simulation P10f1r20. The comparison shows clearly that the simulation with the lower fraction of the mass in galaxies (P10f1r1), i.e. with the higher fraction of the mass in the common halo, evolves slower than the simulation with the higher fraction of the mass in galaxies. This in agreement with what was seen in the simulations of Bode et al. (1994) for poor clusters of galaxies, as well as those of Barnes (1985), Bode et al. (1993) and Athanassoula et al. (1997) for compact groups. This may be easily understood (Athanassoula et al. 1997) if we think of the galaxies as a perturbation of the global common halo potential. If the fraction of the mass in galaxies is higher, the perturbation is also stronger, the galaxies will attract each other more strongly and have more chances of colliding and therefore of merging. The opposite will be true if the fraction of the mass in galaxies is small. In the limiting case were the galaxies are considered as test particles in the common halo, they will collide only if their orbits accidentally cross.

Fig. 4 compares the evolution of two simulations, one, H10f1r1, in which the distribution of the common halo and the galaxies follows a Hernquist law, and the other, P10f1r1, in which it follows a Plummer law, the remaining parameters being the same. The difference is quite striking. We note that the evolution is much more rapid for the Hernquist model. Thus, at any given time, for the Plummer model there are more galaxies that have not fallen into the central one than for the Hernquist model.

### 3.2 Global parameters

We will use in the analysis of the simulations the same definitions of the central galaxy and of the merging conditions as in GAG. They allow us to make a more quantitative study of the evolution of the system. The first factor that we can consider is the number of galaxies that survive at a given time. This is plotted as a function of time in Figures 5 and 6. In the upper panel of Fig. 5 we show the results for the simulations with a Hernquist distribution and in the lower panel for the simulations with a Plummer distribution. In this figure we concentrate on simulations with $`a_p=10`$ and $`f=1`$ and vary only the fraction of mass in the galaxies. The evolution shown in Fig. 5 is in good agreement with the qualitative results of Fig. 3. In general a larger percentage of the total mass in galaxies leads to a faster evolution, for the reasons explained in the previous section. For simulations with Plummer distributions the effect is only important in the second half of the evolution when the central galaxy is sufficiently massive. On the other hand for simulations with Hernquist distributions the effect is clear from the start, and is stronger. These differences are due to the different evolutionary histories of the two types of simulations, which will be discussed in section 4.1. Finally we should also note that the fastest of all evolutions is for simulations with Hernquist profiles and the biggest fraction of the mass in galaxies.

Fig. 6 also shows the number of galaxies that have not merged at a given time as a function of this time, but concentrates on Plummer profiles, to see the effect of the three free parameters of the Plummer distribution. The upper panel illustrates the effect of the scale-length of the distribution of the common halo and of the galaxies, the middle one the effect of the cut-off radius and the lower one the effect of the fraction of mass in the galaxies. The results agree well with what is seen qualitatively in Figures 1 to 3. The galaxy number decreases faster in the cases of more concentrated mass distribution, i.e. halos and galaxy distributions with a smaller scale-length. Also configurations with a smaller cut-off radius of the galaxy distribution evolve faster.

The next quantity that we can consider is the evolution of the mean radius of the system of surviving galaxies. It is calculated as the average distance of the centers of the surviving galaxies from the center of the group not taking into account the giant central galaxy. Its time evolution is shown in Fig. 7. As in Fig. 5 we show in this figure the simulations with $`a_p=10`$ for the Hernquist halo distributions (upper panel) and Plummer distributions (lower panel). In all but a couple of cases the mean radius of the system stays constant. This argues that the satellite galaxies that fall into the central one are not only those that have lower binding energies and orbits with a smaller mean radius. Rather all galaxies, independent of their binding energy, are concerned. The case of H10f1r20 is different, since the mean radius of the surviving galaxies increases with time. This is due to the fact that one galaxy is on an outwards going orbit, while the total number of galaxies is rapidly reduced (cf. Fig. 5) so that this one galaxy can influence the statistics.

The mean mass in surviving galaxies, i.e. not taking into account the giant central galaxy, is indicative of whether the mass loss effects due to stripping are important in the dynamical evolution of the groups, or not. The time evolution of this parameter is shown in Fig. 8 for Hernquist systems (upper panel) and for Plummer systems (lower panel). In all the cases the mass in galaxies diminishes with time, which means that stripping is present and that the stripped mass is not all captured by other individual galaxies. Some of it is captured by the giant central galaxy, while another part stays in between the individual galaxies.

### 3.3 The rigid halo case

If we could, in general, use in simulations a rigid halo, rendered by an external forcing, constant in time, instead of a “live” halo formed by a large number of particles, we could save a considerable amount of CPU time. For this reason we will, in this section, compare the effects of a rigid and a live halo on the evolution of a simulation, to see whether the corresponding gain in CPU is possible without introducing a bias in the results. At the beginning of the simulation a live halo will give, on average, the same forces as a rigid halo. This will, however, not be necessarily the case as the simulation evolves, since the mass distribution of the halo may change with time, in response to the evolution of the distribution of the individual galaxies and to the formation of the giant central galaxy. Furthermore, even in cases where the mass distribution in the common halo hardly evolves, we are not assured of the adequacy of the rigid halo, since it will not reproduce the dynamical friction on the satellite galaxies that the live halo exerts naturally. This effect could in principle be estimated by the Chandrasekhar (1943) formula, except for the fact that the common halo does not necessarily have a homogeneous density distribution, while dynamical friction is much more important in strongly centrally concentrated backgrounds than in less concentrated ones (Athanassoula, Makino & Bosma 1997). It is thus necessary to compare a simulation with a live halo to one with a rigid one, to see how big a bias the rigidity of the halo introduces.

To check this effect we have repeated one simulation, namely P10f3/4r1, using a rigid halo instead of the live halo composed of 50000 particles. In Fig. 9 we compare the evolution of the global parameters of the two systems. We can see that the number of surviving galaxies is considerably lower in the case of a live halo (upper panel). The principal reason for this is the fact that the dynamical friction with the live halo brakes the galaxies and, as they move slower, they merge easier with the central galaxy. The radius which contains half the mass of the galaxy system stays in both cases roughly constant with time and the average value is roughly the same in the two simulations (middle panel). The velocity dispersion of the system of galaxies shows in both cases an increase with time, presumably due to the fact that the faster individual galaxies merge less than the slower one. The lower panel of Fig. 9 shows that this increase is stronger for the live halo than for the rigid one. This should be linked to the fact that there are more mergings in the live halo case.

### 3.4 The effect of the cutoff on the galaxy distribution.

In most of our simulations we have included an artificial cutoff in the initial distribution of the galaxies. This partly reflects the fact that the distribution of galaxies in a group or cluster also does not exceed to infinity. It has the added advantage of improving the statistics for a fixed number of particles per simulation. In this section we discuss its effect by comparing the evolution of runs with cut-off with that of runs without cut-off and otherwise the same properties. Fig. 10 compares two such sets of simulations, simulation P10r20 to P10f1r20, and simulation P10r1 to P10f1r1. In this figure we compare two global properties of the system of surviving galaxies: their number (upper panel) and their velocity dispersion (lower panel). The number of surviving galaxies decreases considerably faster with time if we include a cut-off. This is due to the fact that by constraining the volume, while keeping the number of galaxies fixed, we increase the density of the galaxies and this in turn, as discussed in the beginning of this section, leads to a faster evolution of the system.

The lower panel shows that the velocity dispersion is smaller for the cases without cut-off. There can be two possible explanations to this. The first one relies on the fact that it is the galaxies with the lower velocities that will merge, and therefore systems that involve more merging should also have a higher dispersion of velocities. If, however, this was the explanation then we would note differences in the velocity dispersions after the differences in the number of remaining galaxies occurred, which, as can be seen in Fig. 10 is not the case. We will thus explain the difference by the fact that the velocity dispersion of the group decreases with radius. Thus a more extended group has a lower mean dispersion of velocities than a more compact one.

## 4 The central object

The central object is defined as in GAG. The properties of this object are studied using the same methods described in that paper and we refer the reader to that paper for more details. Here we will discuss directly the physical properties of the objects formed in our simulations.

### 4.1 Stripping versus merging

The mass of the central galaxy grows with time by means of two mechanisms: by merging of satellite galaxies and by accreting material tidally stripped from galaxies. We are interested in how the relative importance of these two mechanisms depends on the initial conditions. In Fig. 11 we show for some selected simulations the time evolution of the total mass in the central regions not belonging to the individual galaxies. In the left panels we show the simulations where the initial density distribution follows a Hernquist profile, while in the right panels we show the cases of initial Plummer density distributions. Each horizontal couple of panels corresponds to two simulations with identical or comparable global properties. The solid line corresponds to the total mass as a function of time, while the dotted line corresponds to the part of this mass that comes from stripped material. The difference between the two lines at a given time gives the part of the total mass that comes from mergers with satellite galaxies.

Let us first note that in general the total mass of the central object is larger in simulations starting with Hernquist profiles than in simulations starting with Plummer profiles. The total mass of the central object (at a given time) is also larger in simulations where the total mass that is initially in galaxies is larger, as expected.

As can readily be seen from Fig. 11 the initial distribution of the mass profile (background mass and distribution of galaxies) is determinant not only for the total mass of the central object, but also for the formation process leading to it. In Plummer simulations all the initial mass increase is due to the accretion of stripped material. This is true for a longer period of time in the case of simulation P10f1r1, where the mass in galaxies was initially the same as the mass in the common background. In this simulation also the final fraction of the total mass that comes from tidal stripping is more important than in the others. On the other hand in the Hernquist simulations the merging of satellite galaxies into the central galaxy is present from the start. This should be due to the differences in central concentration. Our results argue that the central concentration of the Hernquist profile is sufficiently high to compensate for the higher relative velocities of the galaxies and to cause a high initial burst of mergings.

Thus, as can be seen from Fig. 11 the relative importance of the stripping is more significant in the Plummer cases than in the Hernquist ones. Conversely the absolute amount of stripped material in the first part of the simulation is larger for the Hernquist cases. This difference is due to the fact that the total mass of the central object in the first part of the simulation is considerably bigger for Hernquist cases.

These differences between simulations with Plummer and with Hernquist profiles can be readily understood from the different evolutionary histories of the two groups of simulations. In the Plummer case the mass distribution in the central parts is nearly flat and as the initial position of the satellite galaxies is selected following this distribution, no initial concentration of galaxies is present. Thus, there is no strong merging of satellite galaxies into a central object during the first steps of the simulations. The satellite galaxies orbiting within the core radius occasionally merge with some of their companions. These mergings, which do not necessarily occur in the center-most parts increase the number of larger galaxies within the core radius, until two or more of these finally merge to give a giant central object.

On the other hand, the dynamical evolution of the simulations with a Hernquist profile is very different. The strong peak of this density profile towards the central parts gives a strong initial concentration of satellite galaxies which merge very fast to create the seed of a giant central object. We refer to this stage of the evolution as an initial burst of mergers. The central object will accrete some stripped material and merge with some other satellite galaxies entering the central parts. After the strong initial merger, the mass increase of the galaxy is milder than in the initial steps of the simulations. As we will see in subsection 4.3 these two different formation stories give galaxies with very different observational properties. Not all the central objects will have profiles that resemble real giant elliptical or cD galaxies.

### 4.2 Volume density profiles.

Figure 12 shows the volume density profile as a function of radius for simulations P10f1r1, P10f1r20, H10f1r1 and H10f1r20 at $`t=8.2`$ Gyrs. Plummer simulations show clearly two distinct regions, the outer one being much steeper than the inner one. Each of these can be relatively well fitted by a straight line in a log-log plane. The slopes of these two lines, as obtained by eyeball fits, are roughly -2 and -5. The outer part is clearly discernible all through the simulation and its slope, which is roughly equal to the slope of the outer part of the Plummer density profile, stays roughly constant all through the simulation. This is not the case for the inner part, which in some earlier times is more difficult to discern than at $`t=8.2`$ Gyrs. Since our simulations stop at $`T=8.4`$ Gyrs, we have no information on the evolution of the profile at later times.

The density profiles coming from Hernquist simulations do not show the clear existence of two distinct regions, although they do steepen at larger radii. For consistency we also made eyeball fits with two straight lines. The inner slopes thus found are steeper than the corresponding ones for the Plummer simulations, while the slopes of the parts further out are roughly the same. The profile of H10f1r20 has a shallower sloper in the outermost parts, which will be discussed at length in the next subsection.

We compared all these profiles to the Hernquist profile (Hernquist 1990) and found that, if we exclude the innermost and outermost parts, the Plummer simulations with little mass in the common halo, as well as all Hernquist simulations, are well fitted by this profile. This was, however, not the case for the Plummer simulations with a high fraction of the mass in the common halo. The physical implications of this result will be discussed further in the next subsection with the help of the projected density profiles.

### 4.3 Projected density profiles.

The most interesting property of the central galaxies formed in these simulations appears when we study their projected density distributions. Remnants of the merging of a pair of progenitor galaxies have a projected density profile which is well fitted, in all but the innermost parts, by the $`r^{1/4}`$ law (see for instance Barnes and Hernquist 1992, or Barnes 1998 and references therein). A similar result was found for multiple mergings in compact groups (Weil & Hernquist 1996, Athanassoula & Vozikis 1999). Furthermore the objects formed in GAG from multiple mergers in large groups or small clusters also showed this property, except for an excess projected density in the outermost parts in some cD-like cases. Our new simulations with a common halo present a larger variety of profiles.

The three-dimensional density profiles of the central object formed in the Plummer models with a considerable fraction of the mass in the background halo can not be well fitted by a Hernquist law, and therefore does not give an $`r^{1/4}`$ in projection. In fact the larger the fraction of mass in the common halo, the worst the agreement is. On the contrary, if the mass in the common halo is very much smaller than the mass in common galaxies (e.g. in simulation P10f1r20), then we recover initial conditions closer to those of the simulations in GAG, and the profiles are better fitted by a Hernquist law. This form of the profiles does not depend either on the initial scale-length of the mass distribution, or on the radius initially containing all galaxies, and depends only on the initial mass fraction in galaxies. These results are best illustrated in Figures 13 and 14.

The mass distribution in the central parts of simulation P10f1r1, where matter is equally distributed between the galaxies and the background, is heavily influenced by the matter in the common halo, and the violent relaxation which would have otherwise led to an $`r^{1/4}`$ projected profile is hampered. This result can be seen in Fig. 13, where we show the projected density profile of the central object formed in simulation P10f1r1 at four different times. In each of the panels we show the result of 10 projections from random viewing angles, thus making clear that this result is projection independent. At no time does this profile follow an $`r^{1/4}`$ law and thus the central object cannot be associated to any giant elliptical or cD galaxy, as are observed in the central parts of galaxy clusters.

The central object formed in simulation P10f1r20 feels only a small influence from the common halo and develops a density distribution which is well fitted by an $`r^{1/4}`$ law. This can be seen in Fig. 14, where we show the projected density profile of the central object at four different steps of the simulation. This effect is again not due to a particular viewing angle. As was already mentioned, this result is very similar to what was found in GAG since in P10f1r20 there is very little influence from the background.

Simulations where the common background is described initially by a Hernquist distribution give very different results. The projected density profile for the central galaxy in simulation H10f1r1 is shown in Fig. 15. As in the case of simulations P10f1r1, only half the mass was initially in the galaxies. We note that the projected density profile of the central object follows an $`r^{1/4}`$ law in the main body of the object and can thus be associated to a real giant elliptical.

Simulations H10f1r5 and H10f1r20 have a much smaller fraction of their mass in the common halo; roughly 17 and 5 per cent respectively. Their projected radial density profiles are shown for six different time steps in Figures 16 and 17. They show clearly that the radial profiles follow well the $`r^{1/4}`$ law over the main body of the object. Here it can not be claimed that this profile is formed as a response to the massive common halo with a Hernquist radial profile, since the mass in the halo is rather small. These central objects were formed in a situation where there is a strong central concentration of galaxies at the start, thus favouring a strong and rapid merger. Thus giant galaxies formed from a rapid multiple merger in the central part of groups or clusters should have an $`r^{1/4}`$ projected density profile over their main body.

Now let us turn to the excess mass over the $`r^{1/4}`$ law in the outermost parts of the central object. This should be associated with the corresponding excess light in the outermost parts of some BCM’s and giant elliptical galaxies. This light excess has, unfortunately, been called “halo”, a term which could bring confusion with the dark halo. To avoid such a confusion, while being coherent with the nomenclature used so far, we will use the term “luminous halo”. Giant ellipticals with such a luminous halo are called cD galaxies. We will now examine which of our simulations present central objects with such a luminous halo and which can therefore be associated with cD galaxies.

Simulation P10f1r20, which is the only one of those with an initial Plummer mass distribution that shows an $`r^{1/4}`$ profile over the main body of the galaxy, has only a very small density excess over the $`r^{1/4}`$ profile in the outermost parts of the central object, which hardly grows with time. It can thus not be easily associated with a cD galaxy.

Central objects formed in simulations with a Hernquist background present such an excess of matter over the $`r^{1/4}`$ profile in the outermost parts. This is small in cases with a large fraction of the mass in the common background, as H10f1r1, and becomes more important for simulations where the fraction of mass in the galaxies is also more important. Thus the central object in simulation H10f1r5 has a clear outer mass excess and that in H10f1r20 an even bigger one (cf. Figures 16 and 17), independent of the projection selected. The existence of this extra mass is clear relatively early on in the simulations, and increases considerably with time. Thus such objects should be associated to cD galaxies, since they form in the central part of the group, follow an $`r^{1/4}`$ projected density law over the main part part of their body and have excess matter over that profile in the outermost parts.

This inverse correlation between the excess mass in the outermost parts and the amount of mass in the common background halo could be due to the fact that a deviation from the $`r^{1/4}`$ is more difficult to form in the strong $`r^{1/4}`$ background. Alternatively it could be due to the fact that the rate of evolution is slower in simulations with a massive common halo, i.e. that at a given time the three simulations are at different stages of their evolution. The total evolution time in these simulations was twice as long as that of Plummer simulations in order to study the possible influence of the dynamical evolution of the group in the distribution of the luminous halo material of the central object. In Figs. 15, 16 and 17 we show that indeed the luminous halo of the central object increases with time, specially in the cases where there is a considerable fraction of mass initially in the galaxies. As we have already pointed out in section 2, the conversion from computer units to physical units is not unique and thus the times given in figures 15, 16 and 17 could be modified if e.g. a different value was adopted for the masses or the scale-length of the initial galaxies. We nevertheless believe that all simulations were evolved sufficiently long for a clear picture to be obtained.

## 5 Limitations of our models

Simulations of the evolution of galaxy clusters, seen the complexity of the objects to be modelled, use a number of simplifications and assumptions, and the simulations presented here are no exception. We will here discuss the simplifications and assumptions specific to our simulations to get a better understanding of how they might have influenced our results.

The initial number of galaxies in all our simulations is 50. Thus they do not apply to rich clusters, but more to large groups, or poor clusters. They can also apply to sub-condensations in rich clusters – provided of course the remaining part of the cluster does not significantly influence the evolution – or subunits which will come together to form a large cluster.

Our initial conditions are not cosmological. They were constructed specifically so as to allow us to study best the physics of the formation of the central object and particularly the effect of different physical parameters on this formation. We have thus created series of initial conditions where one parameter only was varied, so as to best put forward the effect of this specific parameter. We considered in turn the effect of the scale-length of the mass distribution, of the cut-off radius of the galaxy distribution and of the mass fraction in the common halo. Most important, we have considered both matter distributions with a core and matter distributions with a cusp, and have studied the effects of these two very different distributions on the evolution of the group and on the properties of the central galaxy. In other words we have been able to create our initial conditions such as to put forward best the phenomena we wished to study. Of course cosmological simulations and in particular observations will be necessary to tell us whether all the spectrum of initial conditions considered here are possible or realistic. Our approach is similar to that used e.g. in studying the formation of bars in disc galaxies. It is possible to conclude that the growth rate of the bar is smaller in discs which are hot and/or have an important spheroidal component in their center or main parts, without having to worry about the exact form of the radial density profile of the disc or halo component. Once the essential physics of the bar formation is understood, then one can start worrying about whether this or that type of radial density profile is more realistic, follows best from cosmological initial conditions, or is in better agreement with the observations.

The fraction of mass in the common halo is, in all our simulations, considerably smaller than what is believed to be the case in real groups or clusters. Nevertheless, the range explored is sufficient for us to establish how the fraction of mass in the common halo influences the evolution. It also allowed us to establish that the form of the radial profile of the central object in Plummer simulations is incompatible with observed profiles in cases with a high fraction of mass in a common halo. Simulations with a higher fraction of mass in the common halo would necessitate considerably more CPU time for two reasons. First the total number of particles in the simulations would increase considerably if the number of particles per galaxy is kept constant. Also the time it would take for the central object to form would also increase, thus leading to longer CPU times, if the time-step is kept constant. Nevertheless, it would be interesting to pursue in the future simulation with higher percentage of mass in the common halo, particularly in order to test the formation of the luminous halo in the outer parts of the central objects.

Probably the most serious drawback of our simulations is the fact that our galaxy models are very simplistic. Thus we do not distinguish between luminous and dark matter, and have not considered any galaxies with a disc, with anisotropy or with rotation. Thus we can not distinguish between dark and luminous mass in the central object. Simulations with more realistic galaxy models would necessitate a considerably higher number of particles, particularly if disc galaxies are considered. Such simulations would nevertheless be very interesting, and should be envisaged now that increases in computing power make them possible.

We have also neglected the existence of gas, both cold and hot, both in the galaxies and in the common background halo. Thus effects like ram-pressure stripping of individual galaxies are totally neglected, but this should not influence the formation of the central object. It should of course be kept in mind that the fraction of common background halo in our simulations includes all mass encompassing the whole cluster and not linked to individual galaxies, i.e. it includes also the hot gaseous component.

The above summarise the main simplifications and assumptions underlying our simulations. It is obvious that there is room for improvement, but none of our assumptions should alter the main physical results found in this paper.

## 6 Summary and Discussion

We have performed N-body simulations of the dynamical evolution of groups of 50 galaxies embedded in a common background mass distribution. We considered two density distributions, the Plummer and the Hernquist one. The former has a core, while the latter has a cusp. The main goal of this paper was to compare the evolution of the group and the properties of the central giant galaxies formed in these two very different matter distributions.

As in the case of simulations of groups with no common background distribution of matter (GAG), the evolution of the systems within a common halo is driven mainly by the merging instability (Carnevali et al. 1981). In all cases a giant galaxy is formed in the central parts by merging of some secondary galaxies and by accreting material stripped from the rest of the galaxies by tidal forces. The evolution rate of the systems is very sensitive to the initial conditions and is controlled mainly by the density in the central parts. Distributions with higher densities resulting either from higher central concentrations or higher densities over a larger area, have a faster evolution than distributions with lower densities in the central parts. One should nevertheless remember that the dynamical time of denser configurations is shorter.

To study more quantitatively the evolution of the different simulations we have followed the time evolution of some global parameters, like the number of galaxies that have not yet merged, their mean distance from the center of the group, their mean mass and their velocity dispersion. The total number of galaxies decreases with time in all cases. This decrease is faster for configurations with a larger fraction of the total mass in galaxies. This temporal decrease is particularly important for the first $`2`$ Gyrs of the evolution of cusped initial distributions, where there is a strong initial merging. In the Plummer simulations, where the initial dominant effect is the tidal stripping of material, the temporal decrease of the number of galaxies is only important after the first $`2`$ Gyrs. The effect of the initial conditions on the mean distance of the galaxies from the center of the group, on their mean mass and on the velocity dispersion of the galaxy system is less spectacular.

We have also studied the effect of substituting the live common background halo of particles by a rigid halo with the same parameters by comparing two simulations with the same initial conditions, except that the one has a rigid and the other a live halo. The global evolution of the two simulations is quite different. In particular we find that the number of galaxies that have not yet merged decreases much faster in the simulation with the live halo. This can be easily understood since the evolution in a live halo includes the effects of dynamical friction which brakes the galaxy motion and favours the merging with the giant galaxy. Thus, although simulations with rigid common halos require considerably less CPU, they should be avoided, since their results are not reliable.

The different profiles of the mass distribution influence both the evolution of the group and the properties of the giant galaxy formed in the centre. In the case of a Plummer density distribution the individual galaxies undergo pairwise mergers giving rise to some larger galaxies in the core of the cluster. At some later step these galaxies merge and form the giant central galaxy. On the other hand, in the Hernquist case the high central concentration of galaxies in the center leads to a high “‘burst” of mergers in the central parts during the initial steps of the evolution and the giant central galaxy is quickly formed.

Not only the merging histories but also the quantity of stripped material from the satellite galaxies is different in the two cases. In the Plummer cases the fraction of the total mass in the central galaxy that comes from stripping is very important, especially early on in the simulation. In the Hernquist cases the fraction of the mass of the central galaxy that comes from stripped material is not so important, but amounts nevertheless to a considerable mass, since the total mass in the central object is important.

These different evolutionary histories give rise to central objects with very different observable properties. We have studied the projected density distributions of the giant galaxies at different time-steps and for different random projections. The shape of these profiles depends heavily on the initial distribution of the mass. In the cases of Plummer distributions with half the total mass in the common background, the projected density profile of the giant galaxy does not follow the $`r^{1/4}`$ law at any time during the simulation. The radial density profiles of these central objects are not observed in the real giant galaxies in the center of the clusters. Only for the case where most of the mass is initially in the galaxies does the giant galaxy formed within a Plummer distribution follow the $`r^{1/4}`$ law typical of elliptical galaxies.

On the contrary, the galaxies formed within Hernquist distributions do follow $`r^{1/4}`$ laws in the main parts of their bodies in all the cases. Moreover, in the external parts, they show the typical projected density excess of cD galaxies in the case where a considerable fraction of the mass is initially in the galaxies. This excess over the $`r^{1/4}`$ law is not a projection effect and thus reflects the real matter distribution of the central object. Furthermore, this excess gets larger as the dynamical evolution of the group proceeds.

Thus our simulations suggest that central galaxies with $`r^{1/4}`$ profiles can form either if the fraction of mass in the common halo is very small, or if the density distribution is cusped. Since the amount of matter observed in the common background is considerably larger than what would lead to $`r^{1/4}`$ profiles for all density distributions, our simulations argue strongly in favour of cusped matter distributions, at least for those groups and clusters where a BCM has formed.

Acknowledgments. We would like to thank the referee, Yoko Funato, for helpful comments, Albert Bosma for many useful discussions and Jean-Charles Lambert for his help with the administration of the runs. E.A. would also like to thank IGRAP, the INSU/CNRS and the University of Aix-Marseille I for funds to develop the computing facilities used for the calculations in this paper. A.G. and C.G.G. acknowledge financial support by the Dirección de Investigación Científica y Técnica under contract PB97-0411.

References.

Aarseth S.J., Henon M., Wielen R., 1974, A&A, 37, 183

Athanassoula E., Bosma A., Lambert J.C., Makino J., 1998, MNRAS, 293, 369

Athanassoula E., Makino J., Bosma A., 1997, MNRAS, 286, 825

Athanassoula E., Vozikis Ch. L., 1999, in “Galaxy interactions at low and high redshift.”, eds. J. Barnes & D. Sanders, Kluwer Academic Publishers, Dordrecht, p. 145

Bahcall N.A., Fan X., Chen R., 1997 ApJ, 485, L53

Barnes J., 1985, MNRAS, 215, 517

Barnes J., 1998 in “Galaxies: Interactions and induced star formation.” eds. D. Friedli, L. Martinet & D. Pfenniger, Springer-Verlag, Berlin, p. 275

Barnes J. Hernquist L., 1992, ARAA, 30, 705

Barnes J., Hut P., 1986, Nat, 324, 446

Bode P.W., Berrington R.C., Cohn H.N., Lugger Ph.M., 1994, ApJ, 433, 479

Bode P.W., Cohn H.N., Lugger Ph.M., 1993, ApJ, 416, 17

Carlberg R.E., Yee H.K., Ellingson E., Morris S.L., Abraham R., Gravel P., Pritchet C.J., Smecker-Hane T., Hartwick F.D.A., Hesser J.E., Hutchings J.B., Oke J.B., 1997, ApJ, 485, L13

Carnevali P., Cavaliere A., Santangelo P., 1981, ApJ, 249, 449

Chandrasekhar S., 1943, ApJ, 97, 255

Cirimelle G., Nesci R., Trevese D., 1997, ApJ, 475, 11

Cole S., Lacey C., 1996, MNRAS, 281, 716

Cowie L.L., Binney J., 1977, ApJ, 215, 723

Dubinski j., 1998, ApJ, 502, 141

Eke V.R., Navarro J.F., Frenk C.S., 1998, ApJ, 503, 569

Fabian A.C., Nulsen P.E.J., 1977, MNRAS, 180, 479

Fort B., Mellier Y., 1994, A&AR, 5, 239

Fukushige T., Makino J., 1997, ApJ, 477, L9

Funato Y., Makino J., Ebisuzaki T., 1993, PASJ, 45, 289

Gallagher J.S., Ostriker J.P., 1972, AJ, 77, 288

Garijo A., Athanassoula E., García-Gómez C., 1997, A&A, 327, 930 (GAG)

Hausman M.A., Ostriker J.P., 1978, ApJ, 224, 320

Henry J.P., 1997, ApJ, 489, L1

Hernquist L., 1990, ApJ, 356, 359

Jing Y.P., 1999, ApJ, L69

Kent S.M., Gunn J.E., 1982, AJ, 87, 945

Loeb A., Mao S., 1994, ApJ, 435, L109

Makino N., Asano K., 1999, ApJ, 512, 9

Markevitch M., Vikhilinin A., Forman W.R., Sarazin C.L., 1999, ApJ, 545, 553

McGlynn T.A., Ostriker J.P., 1980, ApJ, 241, 915

Mellier Y., 1999, ARAA, 37, 127

Merritt D., 1984, ApJ, 276, 26

Miraldá-Escudé J., Babul A., 1995, ApJ, 449, 18

Moore B., Gelato S., Jenkins A., Pearce F.R., Quilis V., 2000, ApJ, 535, L21

Moore B., Governato F., Quinn T., Stadel J., Lake G., 1998, ApJ, 499, L5

Mushotzky R.F., Scharf C.A., 1997, ApJ, 482, L13

Navarro J.F., Frenk C.S., White S.D.M., 1996, ApJ, 462, 563

Navarro J.F., Frenk C.S., White S.D.M., 1997, ApJ, 490, 493

Ostriker J.P., Hausman M.A., 1977, ApJ, 217, L125

Ostriker J.P., Tremaine S.D., 1975, ApJ, 202, L113

Richstone D.O., 1976, ApJ, 204, 642

Rosati P., Ceca R.D., Norman C., Giacconi R., 1998, ApJ, 492, L21

Sensui T., Funato Y., Makino J., 1999, PASJ, 51, 1

Tyson J.A., Fischer Ph., 1995, ApJ, 446, L55

Vikhlinin A., McNamara B.R., Forman W., Jones C., Quintana H., Hornstrup A., 1998, ApJ, 498, L21

Weil M.L., Hernquist L., 1996, ApJ, 460, 101

Wu X.P., Chiueh T., Fang L.Z., Xue Y.J., 1998, MNRAS, 301, 861

Yoshida N., Springel V., White S.D.M., Tormen G., 2000, ApJ, 535, L103


format: text

text:
# Charged Particle Fluctuations as a Signal of the Dynamics of Heavy Ion Processes

## 1 Introduction

In the analysis of the hadronic multi-particle production (for a recent review see ) a key observation has been the local compensation of charge. The charge fluctuations connected to the soft hadronic part of the reactions were found to involve only a restricted rapidity range. This observation limited the applicability of statistical models to rather local fluctuations (see e.g.  where the charge fluctuation between the forward and backward hemisphere was discussed).

In heavy ion scattering charge flow measurements should be analogously decisive. It is a central question of an unbiased analysis whether the charges are distributed just randomly or whether there is some of the dynamics left influencing the flow of quantum numbers. This is not an impractical conjecture. In heavy ion experiments the charge distribution of particle contained in a central box with a given rapidity range $`[Y_{\mathrm{max}.},+Y_{\mathrm{max}.}]`$ can be measured and the dispersion of this distribution:

$$<\delta Q^2>=<(Q<Q>)^2>$$
(1)
can be obtained to sufficient accuracy. For sufficiently large gaps this quantity also reflects the long range charge flow.

It was proposed to use this quantity to distinguish between particles emerging from a equilibrized quark-gluon gas or from a equilibrized hadron gas. In a hadron gas each particle species in the box is taken essentially poissonian. In a central region where the average charge flow can be ignored, one obtains a simple relation for particles like pions with charges $`0`$ and $`\pm 1`$ :

$$<\delta Q^2>=<N_{\mathrm{charged}}>.$$
(2)
It is argued in the cited papers that this relation would change in a quark gluon gas to:

$$<\delta Q^2>=\underset{i}{}q_i^2<N_i>=0.19<N_{\mathrm{charged}}>$$
(3)
where $`q_i`$ are the charges of the various quark species and where again a central region is considered. The coefficient on the right was calculated for a two flavor plasma in a thermodynamical consideration which predicts various quark and gluon contributions with suitable assumptions. A largely empirical final charged multiplicity $`N_{\mathrm{charged}}=\frac{2}{3}(N_{\mathrm{glue}}+1.2N_{\mathrm{quark}}+1.2N_{\mathrm{antiquark}})`$ was used.

It should be pointed out that the estimate is not trivial. In the considered $`m_{\mathrm{quark}}=0`$ theory the observable $`Q`$ is not infrared safe and makes no sense. A way to make it well defined is to consider the quantity:

$$Q_{\mathrm{quark}}\stackrel{~}{Q}_{\mathrm{quark}}=Q_{\mathrm{quark}}<Q_{\mathrm{quark}}>$$
(4)
which avoids the influence of extra sea quarks. Such a correction is also needed to have an observable which can be expected to survive hadronization. Numerically the effect is not very big and the problems can be ignored if only a rough description is sufficient.

There are a number of sources of systematic errors in the above comparison. The result strongly depends on what one chooses as primordial and secondary particles. Considering these uncertainties we follow Fiałkowski’s conclusion that a clear cut distinction between the hadron- and the quark gluon gas is rather unlikely. This does not eliminate the interest in the dispersion. The hadron gas model is anyhow no optimal reference point to compare with.

In the next section we discuss various possible measures to observe such fluctuations. We favor the dispersion of the charge transfer. Using an explicit Dual Parton model calculation we observe a clear distinction between models with local compensation of charge and equilibrium approaches. In section 3 a simple interpretation of the dispersion in terms of quark lines is outlined. This suggests to compare the dispersion to the particle density as it is done in section 4.

## 2 Various Measures for Fluctuations in the Charge Distribution

For the analysis of the charge structure several quantities were discussed in the recent literature. It was proposed to look at the particles within a suitable kinematic region and to measure just the mean standard deviation of the ratio $`R`$ of positive to negative particles:

$$<\delta R^2>=<\left(\frac{N_+}{N_{}}<\frac{N_+}{N_{}}>\right)^2>$$
(5)
or the quantity $`F`$:

$$<\delta F^2>=<\left(\frac{Q}{N_{\mathrm{charged}}}<\frac{Q}{N_{\mathrm{charged}}}>\right)^2>$$
(6)
where $`Q=N_+N_{}`$.

We consider them not attractive. The quantities are not suitable for small intervals, as there are actually undefined in a certain region. They are less clean than the dispersion, $`<\delta Q>`$, as they are not exclusively connected to the flavor structure and as they mix up charge and density fluctuations.

For large nuclei at high energies this is not a problem as the density fluctuations are small and all three quantities are connected by the following relations:

$$<N_{\mathrm{charged}}><\delta R^2>=4<N_{\mathrm{charged}}><\delta F^2>=4\frac{<\delta Q^2>}{<N_{\mathrm{charged}}>}.$$
(7)
To examine the range where these relations hold, all three quantities were calculated in the Dual Parton model implementation DPMJET . For the most central $`5\%`$ Pb-Pb scattering at LHC energies ($`\sqrt{s}=6000`$ A GeV) there is indeed a perfect agreement between all three quantities as shown in figure 1. For the most central $`5\%`$ S-S scattering at RHIC energies ($`\sqrt{s}=200`$ A GeV) the agreement is no longer as good. For the minimum bias S-S scattering at RHIC energies the agreement is lost and the first two expressions behave rather erratic. The same erratic behavior is seen for the proton-proton case which is shown figure 2. As any conclusion will strongly depend on a comparison of central processes with minimum bias events, we advocate to stick to the dispersion of the net charge distribution $`<\delta Q^2>`$.

We observed no significant difference between rapidity and pseudo-rapidity boxes.

## 3 A Simple Relation between the Quark Line Structure and Fluctuations in the Charge Flow

To visualize the meaning of charge flow measurements it is helpful to introduce a general factorization hypothesis. It is not exact but it is widely expected to hold to good accuracy. It postulates that the flavor structure of an arbitrary amplitude can be described simply by an overall factor, in which the contribution from individual quark lines factorize111 The hypothesis is based on the exchange degeneracy of octet and singlet Regge trajectories effectively changing the $`SU(N_{\mathrm{effective}})`$ flavor symmetry to an $`U(N_{\mathrm{effective}})`$ symmetry in which this relation is exactly valid. One of the important corrections to the hypothesis originates in the special behaviour of the masses of the lowest lying mesons of the trajectories, which is especially significant in the pseudo-scalar sector between the $`\pi _0`$ and the $`\eta `$ meson. That one of the two neutral states is sometimes suppressed, introduces a mild anticorrelation between neigbouring flavors, which can be ignored for our consideration centered at long range charge transfers. If a higher accuracy is desired the hypothesis can be restricted to primordial particles generated at a high “temperature” and “secondary” charges produced during the decay of primordial particles can be considered extra. With simple assumptions their contribution can be related to the corresponding particle spectrum. If all charged particles were secondaries the dispersion of the charge transfer across an arbitrary rapidity boundary would be given by the Quigg-Thomas relation $`<\delta Q^2(y)>=\sigma \frac{1}{2}\rho _{\mathrm{charged}}(y)`$ where $`\sigma =1`$ if widening and narrowing effects balance. .

The hypothesis can be used to obtain the following generalization of the Quigg-Thomas relation. This generalization states that the correlation of the charges $`Q(y_1)`$ and $`Q(y_2)`$, which are exchanged during the scattering process across two separate kinematic boundaries, is just:

$$<\left\{Q(y_1)<Q(y_1)>\right\}\left\{Q(y_2)<Q(y_2)>\right\}>=n_{\mathrm{common}\mathrm{lines}}<(q<q>)^2>.$$
(8)
where $`n_{\mathrm{common}\mathrm{lines}}`$ counts the number of quark lines intersecting both borders and $`q`$ is the charge of the quark on such a line. Depending on the flavor distribution average values $`<(q<q>)^2=0.22\mathrm{}0.25`$ are obtained.

Most observables of charge fluctuations depend on this basic correlation. The fluctuation of the charges within a $`[Y_{\mathrm{max}.},+Y_{\mathrm{max}.}]`$ box discussed above contains a combination of three such correlations:

$$<\delta Q^2>=<\delta Q(y_1)^2>+<\delta Q(y_2)^2>\mathrm{\hspace{0.25em}2}<\delta Q(y_1)\delta Q(y_2)>.$$
(9)
Using (8) the dispersion of the charges in a box subtracts to:

$$<\delta Q[\mathrm{box}]^2>=n_{\mathrm{lines}\mathrm{entering}\mathrm{box}}<(q<q>)^2$$
(10)
where $`n_{\mathrm{lines}\mathrm{entering}\mathrm{box}}`$ is the number of quark lines entering the box.

## 4 Calculation of the Dispersion of the Charge Distribution within a Box

Let us consider the prediction of a thermodynamic model in more detail. In the thermodynamic limit with an infinite reservoir outside and a finite number of quarks inside all quark lines will connect to the outside as shown in figure 3.

The dispersion of the charge transfer is therefore proportional to the total number of particles inside. In the hadron gas all particles contain two independent quarks each contributing with roughly $`1/4`$ yielding the estimate of (2). For the quark gas only one quark originates in the thermalization and is taken to be responsible for charge transfer yielding a considerable smaller result. Obviously there are several refinements to this simple picture.

Let us consider the limit of a tiny box. Looking only at the first order one trivially obtains:

$$<\delta Q^2>/<N_{\mathrm{charged}}>=1$$
(11)
which corresponds to the hadron gas value.

If the box size increases to one or two units of rapidity on each side this ratio will typically decrease, as most models contain a short range component in the charge fluctuations usually attributed to secondary interactions. One particular short range fluctuation might be caused by the hadronization of partons of the quark gluon gas. The quark antiquark pair needed for the hadronization is assumed to be short range so that for a box of a certain size one just obtains the charge dispersion of the original partons. This is responsible for the reduction of the fluctuation discussed above. The decreasing is however not very distinctive. Common to many models are secondary interactions which involve decay processes and comover interaction. In hadron hadron scattering processes such correlations are known to play a significant role and there is no reason not to expect the same for the heavy ion case.

After a box size passed the short range the decisive region starts. In all global statistical models the ratio will have to reach now a flat value. Only for very large rapidity ranges charge conservation will force the ratio to drop $`\left(1y_{\mathrm{max}.}/Y_{\mathrm{kin}.\mathrm{max}.}\right)`$. This is different in string models as it is illustrated in figures 1 and 2 . The model calculation shows with its rapid fall off a manifestly different behavior. It is a direct consequence of the local compensation of charge contained in string models. The effect is illustrated in Figure 4 in which only quark lines which contribute to the charge flow and which intersect the boundary are shown.

The local compensation of charge allows now only a contribution of lines originating around the boundaries. If the distance is larger that the range of charge compensation the dispersion will no longer increase with the box size. The total contribution will now be just proportional to the density of the particles around the boundaries

$$<\delta Q^2>\rho _{\mathrm{charged}}(y_{\mathrm{max}.}).$$
(12)
This resulting scaling is illustrated in a comparison between both quantities shown in Figure 5 for RHIC and LHC energies.

The agreement is comparable to the proton proton case shown in Figure 6.

The proportionality is expected to hold for a gap with roughly $`\frac{1}{2}\delta y>1`$ . For smaller boxes some of the quark lines seen in the density do not contribute as they intersect both boundaries. For large rapidity sizes there is a minor increase from the leading charge flow $`Q_L`$ originating in the incoming particles. In a more careful consideration one can subtract this contribution

$$<\delta Q^2>_{\mathrm{leading}\mathrm{charge}\mathrm{migration}}=<Q_L>(1<Q_L>)$$
(13)
and concentrate truly on the fluctuation.

The prediction for the proportionality factor for the case of mere short range fluctuations would be roughly a factor one (see footnote 1). In string models primordial particles are responsible for a longer range charge transfer coming from the contributions of the quark resp. diquark fragmentation chains. Taking everything together one obtains

$$<\delta Q^2>=\underset{\mathrm{left}+\mathrm{right}}{}\{n_{\mathrm{strings}}<(q<q>)^2>+f_{\mathrm{secondary}}\sigma \frac{1}{2}\rho _{\mathrm{charged}}(y)\}$$
(14)
where $`n_{\mathrm{strings}}`$ is the number of strings where $`f_{\mathrm{secondary}}`$ is the fraction of secondary to primary particles and where the width of the local fluctuations $`\sigma `$ is roughly unity. In our Dual Parton Model calculation we observe factor of roughly $`1.2`$ between the density and the dispersion. It means that most of the fluctuation originate in secondary interaction and that the effective larger coefficient of the first term which significantly rises the factor plays only a lesser role.

## Conclusion

In the paper we demonstrated that the dispersion of the charge distribution in a central box of varying extend in rapidity is an extremely powerful measure. Within the string model calculation the dispersion seen in relation to the spectra shows no difference between simple proton proton scattering and central lead lead scattering even though both quantities change roughly by a factor of 400.

The dispersion allows to clearly distinguish between conventional string models and thermal models. As most expected changes in the dynamic are somehow connected to the onset of thermalization they will introduce more fluctuations. Even if the truth should lie somewhere in between string models and quark gluon plasma models it is therefore quite reasonable to hope that the position of the transition can be determined in a quantitative way.

## Acknowledgments

F.W. Bopp acknowledges partial support from the INTAS grant 97-31696.


text:
# Highlights of Cern Workshop on Charm Production in A+A Collisions

## 1 Introduction

Basic aim of the Workshop  organised at Cern by U. Heinz and C. Lourenco in December 1999 was the discussion of two questions concerning open charm production in nucleus–nucleus (A+A) collisions at Cern SPS energies:

* do we want to measure it?
* can we measure it?

The workshop should be seen as an important element of an ongoing debate on future of heavy ion physics at Cern SPS. Two proposals were formulated for the experimental studies of A+A collisions in year 2000 and beyond:

* a measurement of open charm production at the top SPS energies (158–200 A$``$GeV) by NA60  and NA49  and
* a study of onset of deconfiment at low SPS energies (20, 30 and 80 A$``$GeV) by NA49 .

The Workshop allowed for a review and critical discussion of theoretical and experimental aspects of open charm proposals.

During the last 25 years a pQCD–based framework for a treatment of open (charm mesons and hyperons) and hidden (charmonia) charm production was formulated and it dominated the interpretation of the experimental data. It is based on the assumption that the perturbative QCD apparatus is applicable for open charm production in hadronic and nuclear collisions and that charmonium creation is due to binding interaction between perturbatively produced $`c`$$`\overline{c}`$ quarks. Since recent analysis of $`J/\psi `$ and dimuon production in Pb+Pb collisions put the standard framework to question, a new and vivid discussion has started on the subject. Several non–standard models for open and hidden charm production were proposed, the most controversial among them are statistical models. During the Workshop data and models were discussed mainly within the standard framework. There are, however, many excellent reviews  written by experts in this subject. Thus in my brief report I will focus on the new developments instead of repeating details of a well known standard approach. For completeness I will mention recent results which were not presented during the Workshop.

In section 2 of this report $`J/\psi `$ production and its relation to open charm are discussed. A summary of various approaches to open charm production is given in section 3. Finally a first attempts to measure an open charm yield in A+A collisions at SPS and future plans for continuation of this effort are presented.

## 2 $`J/\psi `$ Production

The Workshop was opened by Satz  who discussed $`J/\psi `$ production and its relation to open charm in nuclear collisions. The $`J/\psi `$ has been measured in A+A collisions over the last 15 years by NA38 and NA50 Collaborations. This experimental program was mainly motivated by a hypothesis of Matsui and Satz  formulated within standard framework. They suggested that $`J/\psi `$ may serve as a probe of the state of matter created in the early stage of the collision. Recent results and their new interpretations indicate that Matsui–Satz approach may be incorrect.

Two years ago Gorenstein and myself  observed that the $`J/\psi `$ yield in nuclear collisions is proportional to the pion multiplicity. We further showed  that the data may be explained assuming statistical production of $`J/\psi `$ at hadronization. Recently several models in which $`J/\psi `$ creation is due to coalescence of $`c`$$`\overline{c}`$ quarks or $`D`$$`\overline{D}`$ mesons at different stages of collision were formulated. A brief summary of these approaches is given below.

### 2.1 Standard Approach

In the standard approach reviewd by Satz  charmonium production is considered as a three stage process: the creation of a $`c`$$`\overline{c}`$ pair, formation of a bound $`c`$$`\overline{c}`$ state, and its subsequent interaction with the surrounding matter. The first process is calculated within pQCD, whereas non–perturbative QCD is needed to describe the last two stages. The interaction of the $`c`$$`\overline{c}`$ bound state with matter results in suppression of the final $`J/\psi `$ yield in comparison with its initial number. In order to comapre the model to the data, the assumption is made that the initial $`J/\psi `$ yield is proportional to multiplicity of Drell–Yan pairs. Within this framework the significant suppression of $`J/\psi `$ production relative to Drell–Yan which is observed when going from peripheral to central Pb+Pb interactions at 158 A$``$GeV  is attributed to the formation of a QGP in the latter collisions .

The data on $`J/\psi `$ production in proton–nucleus (p+A) interactions serve as a first test of the model. Aichelin  pointed out that recent precise results on $`x_F`$ dependence of the charmonium yield in p+A collisions  are difficult to understand within the standard framework.

### 2.2 Statistical Approach

The statistical approach, presented by Gorenstein and myself , assumes that $`J/\psi `$ mesons are created at hadronization according to the available hadronic phase–space. Thus, within this model, the $`J/\psi `$ yield is independent of the open charm yield. Moreover, the $`J/\psi `$ production is, in good approximation, insensitive to the state of matter created at the early stage. This is because charmonium is created only at the hadronization stage. The model offers natural explanation of the proportionality of the $`J/\psi `$ and pion yields and the magnitude of the multiplicity of $`J/\psi `$ mesons in hadronic and nuclear collisions .

It should be noted that a consistent picture of hadron production within this model is still missing. As an example one can consider $`\psi ^{}`$ production. The statistical approach used for $`J/\psi `$ description works for $`\psi ^{}`$ production in central A+A collisions . It does not give, however, a natural explanation of a substantial increase of the $`\psi ^{}/J/\psi `$ ratio measured in p+p and p+A interactions.

### 2.3 Microscopical Coalescence Model

Microscopical coalescence model , presented by Lévai , assumes that $`J/\psi `$ mesons are formed at hadronization as a result of coalescence of $`c`$$`\overline{c}`$ quarks created in the earlier stages of the collision. By introduction of microscopic coalescence factors and accounting for quark number conservation one can relate $`J/\psi `$ and open charm production. Starting from the measured multiplicity of $`J/\psi `$ mesons in central Pb+Pb collisions at SPS, one predicts mean multiplicity of $`c`$$`\overline{c}`$ pairs, $`n_{c\overline{c}}3`$.

### 2.4 Statistical Coalescence Model

Also in the statistical coalescence model, recently introduced by Braun–Munzinger and Stachel , charmonium states are produced at hadronization as a coalescence of earlier created $`c`$$`\overline{c}`$ quarks. It is further assumed that the number of $`c`$$`\overline{c}`$ quarks is given by pQCD and that the redistribution of these quarks among hidden and open charm hadrons follows the maximum entropy principle. The model reproduces the measured $`J/\psi `$ yield in central Pb+Pb collisions at 158 A$``$GeV using parameters fitted to the light hadron sector.

This model, however, does not give a natural explanation for the observed proportionality of $`J/\psi `$ and pion yields in nuclear collisions.

### 2.5 Coalescence in QGP

A possible contribution from the coalescence of $`c`$$`\overline{c}`$ quarks originating from different elementary interactions was calculated by Thews, Schroedter and Rafelski . The coalescence process is assumed to take place in the QGP phase in parallel with the expected disintegration of charmonium states. Starting from the pQCD charm yield they calculated that the coalescence contribution can lead to enhanced production of $`J/\psi `$ mesons in A+A collisions at RHIC and LHC energies. This prediction should be confronted with the expectation of $`J/\psi `$ suppression derived within the standard approach .

### 2.6 Secondary Hadronic Production

Charmonium production resulting from the interaction of $`D`$ and $`\overline{D}`$ mesons in hadronic matter was considered by Braun–Munzinger and Redlich . The contribution from this process to the total $`J/\psi `$ yield was found to be small even at LHC energies. This is not the case for $`\psi ^{}`$ production, where the secondary production may substantially exceed the primary yield.

Presented models predict very different relations between hidden and open charm production. In the statistical model the $`J/\psi `$ and open charm yields are independent. Initial charmonium production is proportional to $`n_{c\overline{c}}`$ in the standard approach. This relation is however modified by the suppression processes. In the coalescence models $`J/\psi `$ multiplicity increases approximately as $`n_{c\overline{c}}^2`$.

It is therefore clear that data on open charm in A+A collisions are necessary in order to understand mechanism of $`J/\psi `$ production.

## 3 Open Charm Production

Open charm production is not measured in A+A collisions. The presentation of Capelli  was, however, devoted to the question whether open charm yield can be indirectly estimated from the current data. NA50 observes enhancement of dimuon pairs in the mass region between $`\varphi `$ and $`J/\psi `$ peaks (intermediate mass region, IMR) in central A+A collisions at SPS . A possible interpretation of this effect, favoured by NA50, is that it is due to enhanced open charm production. The enhancement factor is relative to the pQCD prediction. Note that the possible enhamcement already questions the validity of the standard approach .

The basic ideas concerning open charm production presented during the Workshop are summarised below.

### 3.1 Standard Approach

In standard model, revied by Vogt , one assumes validity of pQCD calculations for open charm production in hadronic and nuclear interactions. A comparison with the open charm measurements in hadronic interactions serves as a first test of the model. In fact the data can be described assuming reasonable values of input parameters (charm quark mass, $`\mathrm{\Lambda }_{QCD}`$, renormalization and factorisation scales). One should note, however, that uncertainties in this parameters lead to the variation of the calculated open charm multiplicity by a factor of about 100 (Frixione et al. ).

The model with the parameters fitted to the hadronic data can be used to predict open charm multiplicity in p+A and A+A collisions. Unfortunately present p+A data are too sparse for a precise test of the model .

In good approximation the following $`A`$–dependence is expected for central A+A collisions:

$$n_{c\overline{c}}_{AA}=n_{c\overline{c}}_{NN}A^{4/3}.$$
(1)
This results in an estimate of the mean multiplicity of $`D^0+\overline{D}^0`$ mesons in central Pb+Pb collisions at 158 A$``$GeV: $`D^0+\overline{D}^0210^1`$ .

### 3.2 Secondary Production

The pQCD based calculations on open charm production from secondary collisions occurring in the expanding QGP were reported by Braun–Munzinger . This additional contribution to the initial open charm yield was found to be negligible at Cern SPS energies.

### 3.3 Statistical Approach

The statistical model, introduced by Gorenstein and myself , assumes that charm quarks and antiquarks are created according to the early stage partonic phase space. This non–standard assumption was motivated by a success the statistical model of the early stage of A+A collisions in description of strangeness and entropy production. Within the model the expected number of $`c`$$`\overline{c}`$ pairs is about 8 in central Pb+Pb collisions at 158 A$``$GeV and the corresponding number of $`D^0+\overline{D}^0`$ mesons is about 6.

For central collisions of large enough nuclei ($`A>30`$) the grand canonical approximation can be used and consequently one expects that open charm multiplicity increases as $`A^1`$.

The pQCD based and statistical estimates of open charm yield in central Pb+Pb collisions at 158 A$``$GeV differ by a factor of about 30. The predicted $`A`$–dependence in these two approaches is also different. Only direct measurements of open charm production in A+A collisions can uniquely distinguish between these spectacularly different models. In particular open charm measurements should give a unique opportunity to set limits to the applicability of pQCD and statistical models of strong interactions.

## 4 Open Charm Measurements

### 4.1 Present

NA50 (Capelli ) reported increased production of dimuon pairs in the intermediate mass region over the standard sources in A+A collisions. The enhancement factor is about 3 for central Pb+Pb collisions at 158 A$``$GeV. The observed effect can be attributed to enhanced production of open charm.

However other interpretations are also possible. Wang argued  that the data can be reproduced assuming that only the momentum spectrum of charm hadrons is modified in A+A collisions with respect to p+p and p+A interactions. Thermal radiation of dimuons was also shown to reproduce the experimental results on dimuon enhancement .

Finally one can argue  that the background subtraction procedure used by NA50 may lead to biased results in the case of central Pb+Pb collisions.

NA49  made first attempt to estimate an upper limit of mean multiplicity of $`D`$ and $`\overline{D}`$ mesons in central Pb+Pb collisions at 158 A$``$GeV by a direct measurement. In this case invariant mass distribution of identified kaons and pions is studied . Using the currently analysed number of events ($`410^5`$) the estimated upper limit is on the level of the statistical model prediction.

### 4.2 Future

NA49  is planning to significantly increase (up to $`510^6`$ events) the statistics of central Pb+Pb collisions at 158 A$``$GeV during this year’s Pb–run. Statistical resolution of this new data should allow to exclude models predicting $`D^0+\overline{D}^0>1`$.

Shahoian  reported that one of the main goals of the recently proposed and approved experiment NA60  is to study the origin of the dimuon enhancement in the IMR observed by NA50. Main components of NA60 experiment are:

* the present NA50 detector with the muon arm as a basic instrument,
* a new pixel vertex spectrometer and
* a new beam scope.

This set–up will allow for a precise measurements of the track impact parameter at the interaction point. It will also result in significant improvement of the momentum resolution in comparison to that achieved by NA50. Consequently NA60 will be able to distinguish between prompt dimuons (Drell–Yan and thermal contributions) and dimuons from decays of charm hadrons. This should lead to clarification of the origin of dimuon enhancement. NA60 was approved for p–beam in 2001. The decision concerning requested ion beams in 2002 and 2003 will be made in September this year.

Systematic charm measurements in nuclear collisions are key elements in the upcoming RHIC and LHC programs. Indirect and direct measurements of open charm should be possible in PHENIX at RHIC (Averbeck ) and in ALICE at LHC (Safarik).

I would like thank Marco van Leeuwen for comments to the manuscript.

## References