text:
# Warm Dark Matter versus Bumpy Power Spectra

## 1 Introduction

The so-called Cold Dark Matter crisis has led to a vast number of publications trying to solve the problems which seem to be associated with an excess of power on small scales. One possibility to reduce this power is to introduce Warm Dark Matter (i.e. Knebe et al. 2002; Bode, Ostriker & Turok 2001; Avila-Reese et al. 2001; Colin et al. 2000). But another way to decrease power on a certain scale is to introduce a negative feature (“dip”) into an otherwise unperturbed CDM power spectrum.

Several mechanisms have been proposed that could generate such features in the primordial spectrum during the epoch of inflation. Among these are models with Broken Scale Invariance (BSI) (Lesgourgues, Polarski & Starobinsky 1998), and particularly BSI due to phase transitions during inflation (Barriga et al. 2000). Other inflationary models include resonant production of particles (Chung et al. 1999) and non-vacuum initial states as the quantum mechanical origin for inflationary perturbations (Martin, Riazuelo & Sakellariadou 2000). Broad features can also be introduced into the power spectrum by including a running spectral index in slow-roll inflationary models (e.g. Hannestad, Hansen & Villante 2000). However, the parameter range for such broad features is relatively tightly constrained by observational data from the BOOMERanG (de Bernardis et al. 2000) and MAXIMA experiments (Balbi et al. 2000), and some of the features tend to be too broad to be able to yield the scale-dependent effects on the power spectrum we consider here. After inflation, other effects such as pressure-induced oscillations111This is also recorded in the oscillatory nature of the spectrum of the cosmic microwave background fluctuations. in the matter radiation fluid before decoupling also leave an imprint in the overall matter transfer function towards small scales mainly in the form of successive crests and troughs (Eisenstein & Hu 1998). These features will therefore also be part of the cosmological matter power spectrum after matter-radiation equality. Only for extraordinarily large baryon fractions will these features be significant enough to affect the subsequent formation of small and medium scale structures.

Most of the above mechanisms appear capable of producing localized features on any given scale, but to date only large scales have been tested quantitatively. The non-vacuum initial state model by Martin, Riazuelo & Sakellariadou (2000), for instance, can lead to a primordial $`P(k)`$ with a sharply localized peak, where the location $`k_0`$ of that peak is determined by a characteristic built-in scale of their model. On the other hand, a negative feature is possible BSI scenarios which lead to oscillations beyond a certain scale in $`P(k)`$ (Lesgourgues, Polarski & Starobinsky 1998). As $`N`$-body simulations are limited in the $`k`$-range modeled, one might only be able to trace the first dip of those oscillations due to a privileged energy within the inflaton potential that determines the scale $`k_0`$.

Our objective now is to examine the effect of a negative feature in the primordial power spectrum at very small scales that directly affects the formation and evolution of (sub-)galactic halos. In particular we have chosen a negative Gaussian feature on a scale ($`k1.8h\mathrm{Mpc}^1`$). We have not tied our analysis to any of the above models, but we have adopted a more generic approach modeling the dip as a Gaussian in log-space (see below). This differs from earlier studies of ’bumpy power spectra’ (e.g. Lewin & Albrecht, 2001). Such a modification is more similar to the idea of Warm Dark Matter, but the difference lies in the rise of power towards even smaller scales. Where WDM cuts off all power exponentially below a certain threshold, we introduce here a model that has a rather narrow dip in $`P(k)`$ and the detailed shape for the loss of power is completely different, respectively. We have shown in an earlier paper (Knebe, Islam & Silk 2001) that such localized features in the matter power spectrum can obscure the interpretation of the evolution of cluster abundance as an indicator of the cosmological density parameter. In this paper we are now going to show that a similar feature (but on smaller scales) can lead to comparable results as for instance a Warm Dark Matter model. This Dip model can hence be understood as a plausible alternative to WDM without making assumptions about the nature of dark matter and only further studies (and observations) might break the degeneracy between those models found in this analysis.

The structure of this paper closely follows that of the WDM investigation by Knebe et al. (2002). We chose to study the same properties, but this time focusing on the comparison between the $`\mathrm{\Lambda }`$WDM and the new Dip model. We first start with a brief description of the newly introduced Dip model in Section 2, then move on to the simulations themselves in Section 3. In Section 4 we analyze the output from our numerical simulations highlighting the differences in the $`\mathrm{\Lambda }`$WDM and the Dip model. We finally summarize our results and conclude in Section 5.

## 2 The Dip-Model

The $`\mathrm{\Lambda }`$WDM and the fiducial $`\mathrm{\Lambda }`$CDM model used in this paper are the same as those used in Knebe et al. (2002) with the cosmological parameters $`\mathrm{\Omega }_0=1/3`$, $`\lambda _0=2/3`$, $`\sigma _8=0.88`$, $`h=2/3`$, and $`m_{\mathrm{WDM}}=0.5`$keV for $`\mathrm{\Lambda }`$WDM.

For the Dip model we are using the same prescription to introduce a Gaussian feature into an otherwise unperturbed CDM power spectrum as outlined in Knebe, Islam & Silk (2001) and hence the modified power spectrum follows the equation:

$$P_{\mathrm{mod}}(k)=P(k)\left[1A\mathrm{exp}\left\{0.4\left(\frac{\mathrm{log}k\mathrm{log}k_0}{\sigma _{\mathrm{mod}}}\right)^2\right\}\right]$$
(1)
where $`P(k)`$ is the unmodified spectrum and the parameters $`A,\sigma _{\mathrm{mod}}`$, and $`k_0`$ of the dip are given in Table 1.

The power taken away from the standard $`\mathrm{\Lambda }`$CDM model via the dip agrees with the lack of power in the $`\mathrm{\Lambda }`$WDM model to the extend that the integral

$$\sigma ^2=\frac{1}{2\pi ^2}_0^{k=16h\mathrm{Mpc}^1}P(k)k^2𝑑k$$
(2)
is identical for the $`\mathrm{\Lambda }`$WDM and Dip model ($`k=16h\mathrm{Mpc}^1`$ is the point for which the dip meets the unperturbed $`\mathrm{\Lambda }`$CDM power spectrum again). The power in the Dip model only drops in a narrow interval around $`k_0`$ with the width of that interval controlled by the parameter $`\sigma _{\mathrm{mod}}`$. A visualization of the three power spectra used in this study is given in Fig. 1.

Even though the power that has been taken away from the standard $`\mathrm{\Lambda }`$CDM model is similar, both models ($`\mathrm{\Lambda }`$WDM and Dip) are based on different physics. WDM is realized by postulating a finite mass for the dark matter particles in the range of $`m_{\mathrm{WDM}}=0.52`$keV (0.5keV in our case), whereas the Dip model assumes that the CDM picture is still valid but based on a primordial power spectrum that is the result of some non-standard inflationary period as outlined in the Introduction.

## 3 The $`N`$-body Simulations

The simulations of the $`\mathrm{\Lambda }`$WDM and the $`\mathrm{\Lambda }`$CDM model were carried out using the Multiple-Mass ART code (Kravtsov, Klypin & Khokhlov 1997) and the details for these runs can be found in Knebe et al. (2002). The new Dip model was set up using the same multiple-mass technique but evolved with the publicly available Adaptive Mesh Refinement code MLAPM (Knebe, Green & Binney 2001). The number of particles in use for all three runs was 1283 in a box of comoving side length 25$`h^1\mathrm{Mpc}`$.

As both codes (ART and MLAPM) are adaptive mesh refinement codes designed in a very similar fashion, the technical details (i.e. refinement criterion, number of integration steps, etc.) are identical to the ones already outlined in Knebe et al. (2002), and hence we obtained the same force resolution of 3$`h^1\mathrm{kpc}`$ in the new Dip model. A quantitative comparison of these two codes can be found in Knebe, Green & Binney (2001) where it is shown that both produce (nearly) identical results when being run with the same parameter setup. As the realizations of the initial density field were identical too for all simulations, we are able to undertake a one-to-one comparison of the two most massive haloes in all three models.

The halos in each individual model were identified using both the friends-of-friends algorithm (FOF, Davis et al. 1985) and the Bound-Density-Maxima (BDM, Klypin & Holtzman 1997) method.

## 4 Analysis

### 4.1 Power Spectra

Fig. 2 shows $`P(k)`$ as derived from our initial particle distribution at redshift $`z=35`$ along with the non-linear evolution at redshifts $`z=2`$, $`z=1`$ and $`z=0`$.

The Dip is well resolved by our initial conditions and lies inside the $`k`$-range defined by the box size $`B`$=25$`h^1\mathrm{Mpc}`$ ($`k_{\mathrm{min}}=2\pi /B`$) and the particle Nyquist frequency as determined by the total number of particles $`N`$=1283 ($`k_{\mathrm{max}}=\pi /\mathrm{\Delta }x`$ with $`\mathrm{\Delta }x=B^3/N`$).

We can also clearly see that the transfer of power from large to small scales washes out the features initially present in the power spectrum at our starting redshift. The major difference between the WDM and the Dip model in the simulations is that the loss of power is far steeper in the latter but also rises again to the level of the unperturbed $`\mathrm{\Lambda }`$CDM spectrum. In the WDM model the power decreases more gradually and monotonically. We will see that this has an effect on the formation of low-and high-mass halos as well as the amount of substructure to be found in galactic halos. However, the general predictions of both models remain similar.

Unfortunately we are unable to study the formation process of objects smaller than 1010$`h^1\mathrm{M}_{}`$because of the limits imposed by the total number of particles used. Our mass resolution is “only” $`m_p710^8`$$`h^1\mathrm{M}_{}`$ and hence we cannot initially resolve the regime $`k>16h\mathrm{Mpc}^1`$ where the amplitude of the power spectrum in the Dip and the $`\mathrm{\Lambda }`$CDM models join up again.

### 4.2 Dark Matter Density Field

We continue our analysis by inspecting the large-scale structure density field in all three runs. In Fig. 3 we show a slice through the particle distribution for all three models with the particles grey-scaled according to the local overdensity (emphasizing on low density regions).

We can clearly see the ”lumpiness” of $`\mathrm{\Lambda }`$CDM against $`\mathrm{\Lambda }`$WDM and the Dip model, respectively. But when comparing $`\mathrm{\Lambda }`$WDM with the Dip model we notice that the filamentary structure of the Universe appears to be slightly more grainy/evolved in the latter. As pointed out elsewhere (Knebe et al. 2002; Bode, Ostriker & Turok 2000) small size objects tend to form in WDM preferably via filament fragmentation as well as at later times. This phenomenon might not be the case for the Dip model anymore; we expect to find the filaments to be more clumpy initially because the power on very small scales rises again to the level of the ”normal” CDM spectrum. And this can be clearly seen in Fig. 3; the low density dark matter particles appear to be more clustered (within the filaments) as opposed to a more smooth distribution in WDM. The similarity between CDM and the Dip model for low-density regions can be even more emphasized when calculating $`\dot{\rho }`$, the flow of matter. To this extent we used the continuity equation

$$\dot{\rho }=\rho /t=(\rho \stackrel{}{v}),$$
(3)
where the rhs was calculated on a 3003 grid covering the whole computational volume. The function $`\dot{\rho }(\delta )`$ with $`\delta =\rho /\overline{\rho }1`$ has one root in the regime $`\delta <200`$ because matter flows out of low density into high density regions. This defines a ’threshold density’ $`\delta _0`$ via

$$\dot{\rho }(\delta _0)=0.$$
(4)
The redshift evolution of this $`\delta _0`$ is presented in Fig. 4 for all three models. We note that $`\delta _0`$ is at all times larger in $`\mathrm{\Lambda }`$CDM than in $`\mathrm{\Lambda }`$WDM. This can be explained by a faster flow of material out of the voids and into the filaments; the threshold density $`\delta _0`$ where the inflow of matter equals the outflow is shifted towards higher density regions. But the Dip model shows a break at redshift around $`z=1.5`$. At later times it agrees with $`\mathrm{\Lambda }`$CDM whereas it matches the $`\mathrm{\Lambda }`$WDM model at higher redshifts. However, this behavior is difficult to test observationally as we investigated the flow of low density material. An analysis of the two-point correlation function of objects in the mass range $`10^{10}`$$`h^1\mathrm{M}_{}`$$`<M<10^{11}`$$`h^1\mathrm{M}_{}`$ though showed a higher (and matching) clustering amplitude for $`\mathrm{\Lambda }`$WDM and the Dip model when compared to $`\mathrm{\Lambda }`$CDM.

### 4.3 Mass Function of Halos

A quantitative difference between the Dip and the $`\mathrm{\Lambda }`$WDM model with respect to galactic halos becomes apparent when calculating the (cumulative) mass function $`n(>M)`$ as a function of redshift $`z`$. Fig. 5 shows $`n(>M)`$ for all three models at redshifts $`z=5,2,1`$,and $`0`$. In addition to the fact that the Dip and the $`\mathrm{\Lambda }`$WDM model show a similar behavior we also note that, first, at higher redshifts there are fewer high mass objects for the Dip model than for $`\mathrm{\Lambda }`$WDM and second, the curve for the Dip model for $`M<5\times 10^{11}`$$`h^1\mathrm{M}_{}`$ runs parallel at all times (but lower in amplitude) to the fiducial $`\mathrm{\Lambda }`$CDM model.

The former can be explained again by the fact that the drop of power compared to WDM is sharper in the Dip model; the ”excess-lack” of power is on scales that roughly agree with objects of the order $`10^{11}`$$`h^1\mathrm{M}_{}`$, which coincides with the high-mass end of mass function at redshift $`z=5`$. However, during the course of the simulation the high mass end of the (Dip-)mass function converges to the level of the $`\mathrm{\Lambda }`$CDM model.

The latter can be ascribed to the rise of power in the Dip model which seems to counter act the top-down fragmentation of the filaments found in the WDM structure formation scenario (i.e. Knebe et al. 2002; Bode, Ostriker & Turok 2001) in such a way that they tend to be more grainy initially (as already seen in Fig. 3). But those small particle groups do not make it into larger objects as we lack power on the ”connecting” scales; they stay isolated within the filaments.

In Fig. 6 we highlight the number density evolution of objects in the mass range $`10^{10}`$$`h^1\mathrm{M}_{}`$$`<M<10^{11}`$$`h^1\mathrm{M}_{}`$ by plotting the abundance evolution for particle groups within that range out to redshift $`z=5`$. Again, the behavior for WDM and the Dip model are very similar, but we see the trend for a steeper and faster evolution in the latter.

### 4.4 Individual Halos

The remaining analysis will focus on the two most massive halos found in all three runs. To this extent we start with showing the density fields of those two objects. Fig. 7 confirms again that the damping of small scale power in the primordial power spectrum leads to a suppression of substructure in galactic halos. But we can also see that this suppression is far more pronounced in the Dip model than in the WDM structure formation scenario. We attribute this again to the steeper cut-off in the primordial power spectrum.

To quantify these finding we calculated the cumulative circular velocity distribution of satellites within halo #1 and halo #2. The result can be found in Fig. 8 which is accompanied by Table 2 where we give the total numbers of satellites found in each halo at redshifts $`z=0`$ and $`z=1`$, respectively.

There is clear evidence that the Dip model produces even fewer satellites than WDM. Even though there are roughly the same number of satellites at redshift $`z=1`$ in all three models, their number only marginally increases in the Dip model (especially for halo #2). We previously interpreted the suppression of substructure in WDM by balancing the accretion and destruction rate of satellites in the halo (Knebe et al. 2002); as such, the even higher destruction rate in the Dip model can best be explained by the presence even more loosely bound satellites. This again can be ascribed to the faster drop-off in power and can therefore be controlled by the amplitude $`A`$ and the width $`\sigma _{\mathrm{mod}}`$ of the Dip. However, tuning the Dip parameters to bring the number of satellites into agreement with the WDM model will automatically have an effect on, for instance, the mass function as given in Fig. 5. The lower the amplitude the closer it should get to the unperturbed $`\mathrm{\Lambda }`$CDM model. Hence the Dip model provides a way of suppressing substructure but otherwise reproducing the standard CDM structure formation scenario. The even milder increase in number of satellites for halo #2 in the Dip model (cf. Table 2 and Fig. 8) is readily explained by the fact that that particular halo underwent a major merger near redshift $`z=0.2`$. The tidal effects of such a merger tend to destroy the loosely bound satellites in that model more successfully and hence we are left with even fewer than in the $`\mathrm{\Lambda }`$WDM model.

## 5 Summary and Discussion

We investigated the differences of two models that both suppress small scale fluctuations in the primordial power spectrum, namely a warm dark matter and the newly introduced Dip model. While it may be argued that the $`\mathrm{\Lambda }`$WDM model is better constrained in terms of the allowed masses and thus free streaming length of the potential particle candidates, the Dip model considered here makes no assumption about the nature of the dark matter component. This happens at the expense of introducing two additional parameters (the amplitude $`A`$ and width $`\sigma _{\mathrm{mod}}`$ of the dip), although these are determined entirely by inflationary physics and thus before decoupling (on the basis of the models considered in the Introduction). Using both the $`\mathrm{\Lambda }`$WDM and the Dip model as the input to $`N`$-body simulations we showed that their predictions are fairly similar with a few notable differences.

The dip in our model was chosen to mimic the loss of power in a WDM Universe in a sense that the integrated power (out to the point where the dip reaches the unperturbed $`\mathrm{\Lambda }`$CDM level again) in both models is identical. Therefore all differences to be found in the structure formation process are due to the way in which power on small scales is suppressed and the detailed shape of the power spectrum, respectively. The analysis presented in this paper showed that the way the loss of power is realized is only of minor importance for objects within the mass range under investigation. However, there are some quantitative differences, mainly in the shape of the cumulative mass function at the low mass end and the flow of low density material.

For the Dip model the low density material appears to be more in line with that expected for CDM (cf. Fig. 4). This agreement may reflect the rise in power to that of CDM on very small scales. To quantify this effect more accurately, the resolution has to be improved significantly, but on theoretical grounds this behavior is expected. In filaments the ‘background overdensity’ contributed by large wave length modes tends to be lower than in highly clustered regions and consequently the onset of non-linear evolution of the smallest scales will occur later than in cluster environments. In this case traces of the rise in the linear power spectrum at very small scales may be more likely to be observable here than in the clustered regions at the present day. Observations of structures residing in filaments could therefore potentially help to constrain the power spectrum on these very small scales.

Despite the similarities in filamentary structure of the Dip and CDM model compared to the WDM model, the suppression of low-mass halos is even more pronounced. This even stronger suppression of substructure (in galactic halos) can easily be explained by the amplitude of the dip and an adequate choice for $`A`$ will definitely eliminate the deviation from the WDM model. But in that case one will also automatically alter the formation times and sites for low mass halos in general bringing it more into agreement with CDM. The Dip model hence provides a way of fine tuning the number of satellites without deviating too far from the standard CDM paradigm (cf. Fig. 4 and Fig. 5) otherwise.

The other difference between WDM and the Dip model is the shape of the cumulative mass function and the evolution and formation times of halos. Here we showed that even though small mass objects (as well as high mass objects) are suppressed, the shape of $`n(>M)`$ follows more closely that of the $`\mathrm{\Lambda }`$CDM model; the power law behavior at the low mass end is preserved. We interpreted this result with respect to the rise of power in the Dip model which counter acts the top-down fragmentation of filaments. The filaments are more grainy/evolved initially due to the relative excess of power on very small scales compared to the WDM model (cf. Fig. 3). Moreover, the lack of high-mass halos at high redshift in the Dip model might allow us to observationally distinguish between WDM and the Dip model.

The resolution of our simulations did not allow for a more detailed comparison of the models on very small (mass) scales where the Dip model has the same power as the $`\mathrm{\Lambda }`$CDM model again. This we postpone to a future investigation.

To conclude, we remark that the Dip model has proven to be a viable alternative to the WDM structure formation scenario and only future observations as well as better resolved numerical simulations will be able to disentangle the degeneracy between those two models. The only way of breaking this degeneracy might lie within the filamentary structures of the Universe and the low-mass end of the cumulative mass function. Our results suggested a simple power-law behavior for these scales (similar to the $`\mathrm{\Lambda }`$CDM model) whereas the shape of the mass function in the WDM case is more complicated. But to strengthen these findings we clearly need more detailed and much better resolved simulations in the future.

## Acknowledgments

We benefitted from valuable discussions with Brad Gibson. The simulations presented in this paper were carried out on the Beowulf cluster at the Centre for Astrophysics & Supercomputing, Swinburne University and the Oxford Supercomputer Centre. We are grateful to Anatoly Klypin and Andrey Kravtsov for kindly providing a copy of the ART code that was used for running the $`\mathrm{\Lambda }`$CDM and $`\mathrm{\Lambda }`$WDM simulation, respectively. We acknowledge the support of the Swinburne University Research Development Grants Scheme.


format: text

text:
# The effect of non–gravitational gas heating in groups and clusters of galaxies

## 1 Introduction

The $`X`$–ray emission from clusters of galaxies offers a unique means for studying the physics of cosmic baryons and their connection with the processes of galaxy formation and evolution. The first attempt to model the ICM in the framework of the hierarchical clustering scenario assumed its thermodynamical properties to be entirely determined by gravitational processes, such as adiabatic compression during collapse and shock heating by supersonic gas accretion (Kaiser 1986). Since gravity does not have characteristic scales, for an Einstein–de-Sitter cosmology this model predicts rich massive clusters and poor groups to appear as scaled versions of each other. Under the assumptions of emissivity dominated by free–free bremsstrahlung and of hydrostatic equilibrium of the gas, this model predicts $`L_XT^2(1+z)^{3/2}`$ for the shape and evolution of the relation between $`X`$–ray luminosity and gas temperature. Furthermore, if we define the gas entropy as $`S=T/n_e^{2/3}`$ ($`n_e`$: electron number density; e.g. Eke et al. 1998), then the self–similar ICM has $`ST(1+z)^2`$.

This simple model was quickly recognized to fail at accounting for several observational facts: (a) the $`L_X`$$`T`$ relation for nearby clusters is steeper than predicted, with $`L_XT^3`$ for $`T>\text{ }2`$ keV clusters (e.g., David et al. 1993, White, Jones & Forman 1997, Allen & Fabian 1998, Markevitch 1998, Arnaud & Evrard 1999), with a possible further steepening at the group scale, $`T<\text{ }1`$ keV (Ponman et al. 1996, Helsdon & Ponman 2000); (b) no evidence of evolution for its amplitude has been detected out to $`z>\text{ }1`$ (e.g., Mushotzky & Scharf 1997, Donahue et al. 1999, Fairley et al. 2000, Della Ceca et al. 2000, Borgani et al. 2001b, Holden et al. 2002); (c) the gas density profile in central regions of cooler groups is relatively softer than in clusters and, correspondingly, the entropy is higher than predicted by self–similar scaling (e.g., Ponman, Cannon & Navarro 1999, Lloyd–Davis et al. 2000, Finoguenov et al. 2002); (d) colder clusters contain a relatively smaller amount of gas than very hot ones (e.g. Neumann & Arnaud 2001).

A common interpretation for such observational facts requires non–gravitational energy input, which should have taken place during the past history of the ICM. This extra heating would place the gas on a higher adiabat and sustain the gas during the gravitational collapse of the cluster dark matter (DM) halo, preventing it from reaching high density in central regions and suppressing the $`X`$–ray emission. For a fixed specific amount of heating energy per gas particle, the effect is larger for the smaller systems, while being negligible for more massive, hotter systems. Kaiser (1991) suggested for the first time that gas pre–heating and subsequent adiabatic collapse turns into a differential steepening of the $`L_X`$$`T`$ relation at the cluster scales. Evrard & Henry (1991), Navarro, Frenk & White (1995) and Bower (1997) assumed the gas in the cluster core to have a minimum entropy level, which was established at some pre–heating epoch.

Cavaliere, Menci & Tozzi (1998, 1999) were able to predict the correct $`L_X`$$`T`$ slope, both at the cluster and at the group scale, by assuming a pre–heated gas at a temperature $`T0.5`$ keV, which is shocked with different strengths when falling into groups and clusters. A pre–heating producing an isentropic gas distribution has been assumed by Balogh, Babul & Patton (1999) and Tozzi & Norman (2001, TN01 hereafter). In particular, TN01 worked out a series of predictions of observables properties of the ICM with isentropic pre–heating, including the amount of energy feedback required to produce the correct $`L_X`$$`T`$ relation and entropy threshold at the group scale. Brighenti & Mathews (2001) discussed the effect of heating the gas either when it has still to collapse within clusters and groups (external heating) or when it is already accreted (internal heating). They concluded that internal heating is more efficient at accounting at the same time for both the slope of the $`L_X`$$`T`$ relation and the excess entropy in poor clusters.

While most of such analysis agree that an extra heating energy of $`0.5`$–1 keV per gas particle is required, yet no general consensus has been reached about its astrophysical origin. Supernovae (SNe) have been advocated by several authors as a possible source for pre–heating (e.g., Wu, Fabian & Nulsen 1998) and metal enriching (e.g. Loewenstein & Mushotzky 1996, Renzini 1997, Finoguenov & Ponman 1999, Pipino et al. 2002) the ICM. Although some analyses showed that SNe can actually provide an adequate energy budget (Menci & Cavaliere 2000, Lloyd–Davis et al. 2000), other authors claimed that this requires a very high, possibly unrealistic, efficiency for the thermalization of the released energy (Wu, Fabian & Nulsen 2000, Kravtsov & Yepes 2000, Bower et al. 2001).

In this case, alternative sources for ICM heating are required, like SN from a primordial star population, the so–called Pop III stars (e.g., Loewenstein 2001), and nuclear galactic activity associated with QSOs (e.g., Valageas & Silk 1999, Mc Namara et al. 2000, Yamada & Fujita 2001, Nath & Roychowdhury 2002). As for Pop III stars, one expects them not to significantly heat the diffuse medium, in order not to destroy Ly–$`\alpha `$ absorbers at lower redshift, $`z<\text{ }3`$, and not to overpollute the IGM with metals.

To date, only few attempts have been pursued to address in detail the effects of non gravitational heating on $`X`$–ray observable quantities with numerical hydrodynamical simulations (Navarro, Frenk & White 1995, Lewis et al. 2000, Borgani et al. 2001a, Valdarnini 2002). Bialek, Evrard & Mohr (2001) have realized a large set of moderate resolution simulations of clusters and groups, setting different entropy floors at very high redshift. They found that agreement with observational constraints on $`X`$–ray cluster scaling properties requires an entropy floor of about 100 keV cm2, which, at their heating redshift $`z_h21`$, correspond to an extra energy of about 0.2 keV/part.

As an alternative to extra–heating, radiative cooling has been also suggested as a mechanism to break ICM self–similarity. In this case, $`X`$–ray luminosity is suppressed as a consequence of the reduced amount of gas in the hot diffuse phase (Pearce et al. 2000, Muanwong et al. 2001). In addition, low–entropy gas in central cluster regions is selectively removed from the hot phase as a consequence of its short cooling time, thus leaving only high entropy gas (Voit & Bryan 2001). However, cooling in itself is known to convert into a cold stellar phase a large fraction of gas, $`>30\%`$ (e.g., Lewis et al. 2000, Davé et al. 2001), much larger than the $`<\text{ }10\%`$ baryon fraction locked into stars indicated by observations. This underlines the need for extra heating to prevent overcooling of the ICM (e.g. Suginohara & Ostriker 1998, Prunet & Blanchard 1999, Balogh et al. 2001).

The main aim of this paper is to study in detail the effect of non–gravitational energy injection into the ICM and how much of this energy is needed to correctly reproduce different $`X`$–ray observable properties of groups and clusters of galaxies. In a previous paper (Borgani, Governato, Wadsley et al. 2001a, Paper I hereafter) we concentrated on the effect of non–gravitational heating on the entropy pattern of the ICM. As a main conclusion, we showed that the entropy excess in the central regions of poor clusters (Ponman et al. 1999) requires an extra heating of about 1 keV/particle. This paper presents an extended analysis of a larger set of simulations, which include several schemes of ICM extra heating. We will also show results from a simulation of a poor group, which includes gas cooling and a star formation algorithm to remove dense cold gas particles from the SPH computation. This simulation will be used to discuss the rôle that cooling plays in determining the ICM thermodynamics.

The plan of the paper is as follows. In Section 2 we give a short description of the Tree+SPH code Gasoline and provide details about the simulations. Section 3 describes the different recipes for non–gravitational heating that we used. Section 4 is devoted to the presentation of the results. In Section 5 we discuss the effect of cooling, while a summary of our main results and our conclusions are presented in Section 6.

## 2 Numerical Method

### 2.1 The Code

We use Gasoline, a new N-body/Smoothed Particle Hydrodynamics code. Gasoline performs gravity and hydrodynamics operations using a binary tree-structure with periodic boundary conditions for cosmology. It uses concurrent multiple timesteps for increased throughput and runs in parallel on a wide variety of architectures using, e.g., pthreads, shared memory or the Message Passing Interface (MPI). The code handles gravity, pressure gradients, hydrodynamic shocks, radiative cooling, photoionizing UV background and was extended for this work to incorporate external heating sources. For a full description of Gasoline, including extended tests of the code, see Wadsley, Stadel & Quinn (2002, preparation). For a description of the gravity code on which Gasoline is based, see Stadel (2001). The external heating SN described in Section 3 is treated within the code as an additional component added locally to the standard SPH compressive heating (“PdV”) term.

As a demonstration of the effectiveness of Gasoline for cluster simulations, we evolve the initial conditions from the “Santa Barbara Cluster Comparison” (Frenk et al. 2000). We resolve the resulting luminous $`X`$–ray cluster with $`15000`$ gas particles within the virial radius ($`2.7`$ Mpc) at the current epoch with $`64^3`$ gas and $`64^3`$ dark matter particles in the whole box. This resolution was used by the majority of the cluster comparison SPH codes including the widely used HYDRA code (Couchman 1991). The test evolved the gas adiabatically without additional heat sources or sinks. The gas profiles for Gasoline with curves for HYDRA and the average result from the cluster comparison are shown in Figure 1. The innermost radial bin plotted in each case is indicative of the resolution. The variation seen is similar to that found between cluster comparison code results, even among different SPH codes, due to differing implementations and the interpretation of the initial condition data which were given to the coders at $`256^3`$ resolution. The central temperature in particular is sensitive to shock treatment and, most of all, minor variations in the timing of small subclump mergers (for an extended discussion see Wadsley et al. 2002, in preparation). Except for these small variations, the Gasoline results are in excellent agreement with the other codes.

### 2.2 The Simulations

We simulated at high resolution four halos taken from a DM–only cosmological simulation within a 100 Mpc box of a $`\mathrm{\Lambda }`$CDM model with $`\mathrm{\Omega }_m=0.3`$, $`\mathrm{\Omega }_\mathrm{\Lambda }=0.7`$, $`\sigma _8`$=1 and Hubble constant $`H_0=70`$ km s-1 Mpc-1. The baryon density parameter was chosen to be $`\mathrm{\Omega }_b=0.019h^2`$, as inferred from the deuterium abundance detected in QSO absorption lines (e.g., Burles & Tytler 1998), which corresponds to a baryon fraction $`f_{\mathrm{bar}}0.13`$.

From the simulated box we selected four halos with masses corresponding to virial temperatures in the range $`0.5<\text{ }T<\text{ }8`$ keV. With this choice, we are sure to encompass the temperature interval where the $`L_X`$$`T`$ relation is observed to steepen substantially and where the effect of energy injection should be more important. The sizes of the four simulated structures are of the same order as those of the Coma cluster (e.g. Geller, Diaferio & Kurtz 1999), the Virgo cluster (Bingelli, Tammann & Sandage 1987), the Fornax group (e.g. Drinkwater et al. 2000) and of a smaller Hickson–like group (e.g. Pildis, Bregman & Evrard 1995). In the following we indicate with these names the corresponding simulated structures.

Using the so–called renormalization technique (Katz & White 1983), the simulation within the periodic box is centered around each re–simulated halo, the mass distribution and the perturbation spectrum are resampled at higher resolution within the Lagrangian region that contains the halo of interest. The high–resolution region extends typically for a few comoving Mpc around the halo. The outer regions are resampled at decreasing resolution. As the four halos span almost two decades in mass, the numerical resolutions differ accordingly, so as to keep comparable the number of particles within the virial radius and the force resolution in units of the virial radius. Typically 0.5–1.5 million particles are used for each run, with about half of them to simulate the high resolution central regions. A summary of the characteristics of the simulations are reported in Table 1. Note that halo particle number refers to the number of particles within the virial radius, and not just in the high resolution region, which typically contains a few times more particles. In order to investigate the effect of resolution, the Virgo, Fornax and Hickson runs have been realized with different mass and force resolution. The high–resolution version of each run has a force softening which is at most 1% of the virial radius, while the mass resolution is chosen so that at least $`5\times 10^4`$ gas particles fall within that radius by $`z=0`$. For all the runs we used $`\eta =0.2`$ and a force accuracy parameter of the tree algorithm of $`\theta =0.7`$ ($`\theta =0.5`$ for $`z>2`$), as suggested by Moore et al. (1998) and later verified by Power et al. (2002), for the parameters regulating time accuracy of the integration and the accuracy of the tree algorithm.

We show in Figure 2 the gas density maps for the high–resolution regions of the simulated structures at $`z=0`$. The quite high resolution achieved allows us to resolve numerous substructures. Some of them retain part of their individual gas content, which survives for a few crossing times before being stripped, as well as discontinuities associated with bow shocks (see (see Paper I for a discussion of these features as seen from the entropy maps of simulations). It is tempting to associate such structures to the small–scale features revealed by Chandra observations, which in some cases are associated with bow shocks or cold fronts from merging structures (e.g., Markevitch et al. 2000, Mazzotta et al. 2001, Ettori & Fabian 2001, and references therein). At $`z=0`$ the Coma–like halo is undergoing several merging events with massive sub–groups, at different stages of advancement, much like what is seen in the XMM observations of the real Coma cluster (Briel et al. 2001).

### 2.3 Definitions of observables

The bolometric luminosity for a set of $`N_{gas}`$ gas particles is defined as

$$L_X=\frac{m_{gas}}{\mu m_p}\underset{i=1}{\overset{N_{gas}}{}}\frac{\rho _i}{\mu m_p}\mathrm{\Lambda }_c(T_i),$$
(1)
where $`m_p`$ is the proton mass, $`m_{gas}`$ is the mass of a gas particle (see Table 1), $`\rho _i`$ and $`T_i`$ are the density and temperature at the position of the $`i`$-th particle, respectively, and $`\mu =0.6`$ is the mean molecular weight for a primordial gas composition with 76% mass provided by hydrogen. Assuming that the main contribution to the $`X`$–ray emission is from free–free bremsstrahlung, the cooling function turns out to be $`\mathrm{\Lambda }_c(T)=1.2\times 10^{24}(kT/\mathrm{keV})^{1/2}`$ erg cm3 s-1 (e.g., Navarro et al. 1995).

Under the assumption of isothermal gas following the spherically symmetric DM distribution, it can be shown that the bolometric $`X`$–ray luminosity is

$$L_X=\left(\frac{f_{gas}}{3\mu m_p}\right)^2\mathrm{\Delta }_{vir}M_{vir}\overline{\rho }F(c)\mathrm{\Lambda }_c(T),$$
(2)
(e.g., Eke et al. 1998) where $`f_{gas}`$ is the gas fraction within the cluster, $`\mathrm{\Delta }_{vir}`$ the overdensity corresponding to virialization ($`178`$ for $`\mathrm{\Omega }_m=1`$; Eke, Cole & Frenk 1996), $`\overline{\rho }`$ the average cosmic density and $`M_{vir}`$ the virial cluster mass (i.e., the mass enclosed within the radius $`R_{vir}`$ within which the average density is $`\mathrm{\Delta }_{vir}\overline{\rho }`$). Under the above assumptions, mass and temperature are related according to

$`k_BT`$ $`=`$ $`{\displaystyle \frac{1.38}{\beta }}\left({\displaystyle \frac{M_{\mathrm{vir}}}{10^{15}h^1M_{}}}\right)^{2/3}`$ (3)
$`\times `$ $`\left[\mathrm{\Omega }_m\mathrm{\Delta }_{vir}(z)\right]^{1/3}(1+z)\mathrm{keV},`$
(e.g., Eke et al. 1996). The $`\beta `$ parameter is defined as the ratio of the specific kinetic energy of the collisionless DM to the specific thermal energy of the gas,

$$\beta =\frac{\mu m_p\sigma _v^2}{k_BT},$$
(4)
being $`\sigma _v`$ the one–dimensional velocity dispersion of DM particles. For a Navarro-Frenk-White (1997, NFW) density profile, $`F(c)`$ is a function of the concentration parameter $`c`$ only:

$$F(c)=c^3\frac{1(1+c)^3}{[\mathrm{ln}(1+c)c/(1+c)]^2}.$$
(5)
For the sake of comparison between the prediction of eq.(2) with results from simulations, we compute the concentration parameter as a function of the halo mass for our choice of cosmological parameters following the prescription of Eke, Navarro & Steinmetz (2001).

Besides the bremsstrahlung emission, line emission gives a non–negligible contribution at temperatures $`T<2`$ keV. To account for it and correctly compare our results to observations, we use the Raymond–Smith (1977) code in the XSPEC package to estimate the emissivity of a plasma with metallicity $`Z=0.3Z_{}`$.

We define the mass–weighted temperature of the cluster as

$$T_{mw}=\frac{1}{N_{gas}}\underset{i}{}T_i,$$
(6)
where the sum is over all the particles falling within the cluster virial radius. A more useful quantity for comparison with observations is the emission–weighted temperature, which is defined as

$$T_{ew}=\frac{\underset{i}{}\rho _iT_i^{3/2}}{_i\rho _iT_i^{1/2}},$$
(7)
This definition strictly holds only for bremsstrahlung emissivity. We verified that final values of $`T_{ew}`$ are essentially unchanged if we account for the contribution from metal lines.

Finally, we define the gas entropy carried by the $`i`$-th particle as

$$s_i=\frac{kT_i}{n_{e,i}^{2/3}}\mathrm{keV}\mathrm{cm}^2,$$
(8)
where $`n_{e,i}`$ is the electron number density associated with that particle.

In Table 2 we give the values of bolometric luminosity, mass and emission–weighted temperatures for the simulated structures and for the different pre–heating schemes (see below).

### 2.4 Effects of resolution

It is generally regarded that a few thousands of DM particles within the virial radius are required to describe global quantities, such as total virial mass $`M_{vir}`$ and velocity dispersion $`\sigma _v`$, and as many gas particles are also sufficient to correctly estimate the ICM temperature. However, a resolution higher by about one order of magnitude is required to correctly describe the $`X`$–ray luminosity: since this quantity depends on the square of the local gas density, its correct computation requires the details of cluster structure in the innermost regions and gas clumping to be accurately resolved (e.g., Anninos & Norman 1996, Navarro et al. 1995, Eke et al. 1998, Bryan & Norman 1998, Lewis et al. 2000).

Resolution is even more an issue if one wants to properly resolve the internal structure of the ICM, including shock patterns of the infalling gas, which carry information on the physics of the diffuse medium. The higher the resolution, the larger the number of substructures surviving within the virial radius of the main cluster halo. The largest sub–halos are able to keep a fraction of their gas and can carry it down to the center of the main halo, possibly shocking the gas and raising the central entropy. The $`X`$–ray luminosity profiles also depend on the underlying matter profile, for which about a hundred thousand particles are necessary to avoid underestimating its central density (Moore et al. 1998).

In Figure 3 we compare the relation between $`X`$–ray luminosity and mass–weighted temperature, $`T_{mw}`$, for the HR and LR runs with gravitational heating (GH) only. Increasing the resolution has a non–negligible effect on the estimated $`X`$–ray luminosity, while, as expected, has only a marginal effect on the mass–weighted temperature. The $`L_X`$ value for the HR run of the Virgo cluster is $`25\%`$ higher than for the LR run, the difference increasing to $`50\%`$ for the Fornax group and to $`100\%`$ for the Hickson group. We also verified that decreasing the softening by a factor two for the Fornax and Hickson runs, while keeping the mass resolution fixed, increases $`L_X`$ by about 10%, thus showing that the adopted spatial resolution is adequate to resolve all the structures which are responsible for the $`X`$–ray emission. In fact, at the highest resolution of our runs, the number of massive subhalos (i.e., with circular velocity larger than one–tenth of that of the main halo) is likely to have converged (Ghigna et al. 2000). We expect spatial resolution to be even less of an issue when considering simulations with extra heating (see below). Indeed, in this case the gas is put on a higher adiabat and, therefore, does not follow the small–scale potential wells which are characterized by a virial temperature of a few tens of keV.

We stress that the HR runs predict a $`L_X`$$`T`$ relation which agrees well, both in slope and normalization, with the prediction of the scaling model of eq.(2) for a pure bremsstrahlung cooling function. We consider this agreement as a convincing indication that the HR runs provide a correct description of the gas distribution, which requires using at least $`5\times 10^4`$ particles within the virial radius. This result is in agreement with that from resolution studies involving the collisionless component only (e.g. Moore et al. 1998), which established the minimum number of particles required to correctly model the central profile of dark matter halos. Fig. 3 also highlights a common problem of simulations of large cosmological volumes at fixed mass resolution: as resolution becomes worse for smaller masses, $`L_X`$ becomes underestimated. For instance, taking the same mass and spatial resolution for the Virgo, Fornax and Hickson runs would produce spurious steepening of the $`L_X`$$`T`$ relation. Finally, while a pure bremsstrahlung emissivity is a good approximation at $`T>\text{ }2`$ keV, the effect of line emission is shown to become important for smaller/cooler systems. The net effect of properly accounting for them is that of flattening the $`L_X`$$`T`$ relation, thus further increasing the discrepancy with respect to observations at the scale of poor clusters and groups.

In Figure 4 we show the effect of the resolution on the profiles of gas density and mean entropy per particle. There is an overall good convergence of these quantities, at the radii encompassing at least 100 gas particles (tick curves), with slight deviations at somewhat larger scales for the Hickson run. This number of particles corresponds to about three times the number of neighbors over which the SPH smoothing length is estimated. Therefore, the radius containing this many particles indicates the smallest scale which is resolved by the SPH kernel with the given softening. On smaller scales the gas density tends to be somewhat underestimated and so is the $`X`$–ray luminosity, while two–body heating causes some flattening of the entropy profile (e.g., Steinmetz & White 1997, Yoshida et al. 2001). As mass and spatial resolutions are increased, gas and entropy profiles steepen in the innermost regions. The absence of an entropy core agrees with predictions from analytical models based on gravitational gas heating and spherical accretion: gas residing in the central cluster regions has been accreted at the very beginning and, therefore, never significantly shocked, thus preserving its initial low entropy.

Besides extra heating, a further possibility to increase the entropy level in central cluster regions could be by anisotropic accretion of gas which has been previously gravitationally shocked by large–scale filaments. In order to verify this possibility, we identify at $`z=0`$ gas particles in the central region of a halo and trace them back at higher redshift. As an example, we show in Figure 5 the density–entropy scatter plot for the gas particles falling at $`z=0`$ within $`0.1R_{vir}`$ in the HR Virgo simulation, and their corresponding distribution at $`z=1`$. In the left panel these particles define a darker strip, which correspond to the cluster core, while the lighter strip at lower entropy witnesses the presence of a merged subclump which still preserves its identity (see also Fig. 2). At $`z=1`$ these same particles define several individual subclumps, which by the present time have merged together at the center of the main halo. Such groups at $`z=1`$ correspond to lower entropy structures, in line with the expectation that the gas within smaller structures undergoes weaker shocks. Gas possibly accreted along filaments would correspond at $`z=1`$ to structures with density contrast $`\delta 10`$ containing shocked gas with $`S>\text{ }50`$ keV cm2. Actually, there is no trace of such gas particles in the right panel of Fig. 5. This shows that gas accreted along filaments does not penetrate efficiently into the cluster central region and contribute to increasing the entropy there.

## 3 Beyond gravitational heating

### 3.1 Entropy floor

Our first scheme for non–gravitational heating is based on setting a minimum entropy value at some pre–collapse redshift (e.g. Navarro et al. 1995, Bialek et al. 2001, TN01). For gas with local electron number density $`n_e`$ and temperature $`T`$, expressed in keV, at redshift $`z`$, we define the entropy as

$$S=\frac{T}{n_e^{2/3}}=\left[\frac{f_{\mathrm{bar}}}{m_p}\frac{1+X}{2}\overline{\rho }(z_h)(1+\delta _\mathrm{g})\right]^{2/3}T\mathrm{keV}\mathrm{cm}^2,$$
(9)
where $`\overline{\rho }(z)=\overline{\rho }_0(1+z)^3`$ is the average cosmic matter density at redshift $`z`$, $`\delta _\mathrm{g}`$ the gas overdensity and $`X`$ the hydrogen mass fraction. We choose two values for this entropy floor, $`S_{\mathrm{fl}}=50`$ and 100 keV cm2, that bracket the observed values in small groups and clusters (Ponman et al. 1999, Lloyd–Davis et al. 2000).

We assume $`z_h=3`$ for the reference heating redshift, since it is close to the epoch at which sources of heating, like SN or AGNs, are expected to reach their maximum activity (see below). At $`z=3`$, we select all the gas particles with overdensity $`\delta _\mathrm{g}>5`$, so that they correspond to structures which have already undergone turnaround. After assuming a minimum floor entropy, $`S_{\mathrm{fl}}`$, each gas particle having $`s_i<S_{\mathrm{fl}}`$ is assigned an extra thermal energy, so as to bring its entropy to the floor value. We estimate the amount of energy injected in the ICM in these pre–heating schemes by selecting at $`z=0`$ all the gas particles within the virial radius and tracing them back to $`z=3`$. We find that taking $`S_{\mathrm{fl}}=50`$ keV cm2 amounts to give an average extra heating energy of $`E_h=\frac{3}{2}T_h1.4`$ keV/part for particles that end up within the virial radius of the Fornax and Hickson groups at $`z=0`$, and $`E_h0.9`$ keV/part for the Virgo cluster and $`E_h0.8`$ keV/part for the Coma cluster. Such values are twice as large for $`S_{\mathrm{fl}}=100`$ keV cm2.

In the semi–analytical model by TN01, the gas is assumed not to have undergone any significant gravitational shock heating before $`z_h`$, so that eq.(9) actually gives the extra energy to be provided to the gas particles to bring them to the appropriate adiabat. The hierarchical clustering scenario, instead, predicts that a significant amount of non–linear structures already exist at $`z_h=3`$, with gas particles raised by gravitational shocks to an entropy level already higher than the floor to be set. In our heating scheme, the thermal energy of such particles is left unchanged.

This heating scheme can be defined as an external one, in the sense that it is also targeted at particles not belonging to collapsed structures. In principle this may violate constraints on the entropy of the diffuse high–$`z`$ IGM, as inferred from observations of Ly–$`\alpha `$ absorption systems. If we identify a Ly–$`\alpha `$ system as a structure with typical gas overdensity $`\delta _\mathrm{g}30`$, it has at $`z=3`$ a temperature $`T2.5\times 10^3`$ keV, as determined by the condition of ionization equilibrium with the UV–background (e.g. Haardt & Madau 2001). This corresponds to a typical entropy $`S5`$–10 keV cm-2, thus about an order of magnitude smaller than that relevant for ICM heating. Therefore, an external heating such as the one imposed in our simulations would destroy Ly–$`\alpha `$ absorbers. In this sense, setting an entropy floor must just be considered as a approach guided by ICM phenomenology rather than motivated by general astrophysical arguments.

Even in this simple heating scenario, a further parameter is represented by the epoch at which the entropy floor is created. As the heating redshift is changed, two opposite effects compete in determining the energy budget required to create a given entropy floor. On the one hand, the higher $`z_h`$, the smoother the particle distribution, the smaller the number of particles at high $`\delta _\mathrm{g}`$, which require a lot of heating to increase their entropy to the desired level. On the other hand, the higher $`z_h`$, the higher the overall cosmic mean density $`\overline{\rho }(z)`$, the smaller the amount of heating energy needed at fixed $`\delta _g`$. We check the effect of changing $`z_h`$ by also running low–resolution simulations of the Fornax group with $`z_h=1`$, 2 and 5. We will discuss in the following the effect of changing $`z_h`$ on the final results. The heating energies corresponding to the different $`z_h`$ are shown in the top panel of Figure 13.

### 3.2 Energy feedback from supernovae

This pre–heating scheme is based on computing the star–formation rate within clusters as predicted by a semi–analytic model of galaxy formation, a technique first introduced by White & Frenk (1991), Kauffmann, White & Guiderdoni (1993) and Cole et al. (1994), and subsequently adopted and developed by several authors (e.g. Somerville & Primack 1999, Wu, Fabian & Nulsen 2000, Cole et al. 2000). Here we use a variation of the scheme described by Menci & Cavaliere (2000), and we refer to this paper for a detailed description of the method.

In our approach, the merging history of dark–matter halos, having the same mass of the simulated ones, is predicted according to the Extended Press–Schechter theory (EPST, e.g. Lacey & Cole 1993). The processes of gas cooling, star formation and stellar feedback within galaxy–sized structures are described by means of a suitable parametrization of the corresponding star formation rate and reheated mass, which is ejected back into the hot diffuse medium. The free parameters of the model are chosen so as to reproduce observed properties of the local galaxy population (Tully–Fisher relation, B– and K–band luminosity functions, disk–sizes).

In Figure 6 we show the integrated star formation history $`\dot{m}_{}(z,M_0)`$ of all the condensations which are incorporated into a structure of total mass $`M_0`$ by the present time. The curves in Fig. 6 correspond to $`M_0`$ values, which span the whole range from galaxies to rich clusters. As $`M_0`$ grows the local star formation rate (SFR) decreases, while the value at $`z>\text{ }2`$ increases and its peak is attained at higher redshift. This is expected in the hierarchical clustering picture, since more massive structures are originated from more biased regions, where initially the star formation rate was higher. The large consumption of cold baryons in the progenitors of $`M_0`$ at large $`z`$ results in a smaller amount of available star-forming gas at small $`z`$, yielding the sharply declining $`\dot{m}_{}`$ at $`z0`$.

After assuming a Salpeter initial mass function (Salpeter 1955), the SFRs are then used to derive the SN heating appropriate to each system that we simulate, based on its virial mass. During the evolution of each halo, this energy is shared among all the gas particles having $`\delta _g>50`$. This overdensity threshold, which roughly corresponds to the density contrast at the virial radius, guarantees that gas heating takes place inside virialized regions. We verified that the final results do not change if only gas particles with $`\delta _g>500`$ are heated. Under the extreme assumption that all the energy released by SN is thermalized into the ICM (i.e. $`\eta =1`$ for the SN efficiency) this scheme dumps a total amount of about 0.35 keV per gas particle. For the Hickson group we also run two more simulations with $`\eta =0.1`$ and $`\eta =4`$. In the first case, the energy budget is so small that the final results are essentially indistinguishable from the GH run (see Table 2) and, therefore, we will not comment further on this low–efficiency case. As for the second case, it corresponds to an amount of extra energy similar to that of the S–50 entropy floor. This allows us to check whether the final results depend on the way a fixed amount of energy is dumped into the diffuse medium. Such a large amount of energy can be interpreted as associated with any heating source whose evolution follows that of the SFR. This is likely the case for AGNs, whose emissivity per unit volume has been suggested to evolve according to the SFR (e.g. Cavaliere & Padovani 1989, Boyle & Terlevich 1998, Franceschini et al. 1999).

A word of caution should be mentioned here about our implementation of SN heating. We assume that all SNe release a fixed amount of energy, $`E_{SN}=10^{51}`$ ergs. In the so–called hypernovae scenario, energies larger by at least one order of magnitude may be expected (e.g. Loewenstein 2001).

Furthermore, in our scenario only Type II SN are included, while Type Ia SN may in principle give a significant contribution to the ICM energetics (e.g., Pipino et al. 2002). Also, most of our simulations assume no radiative losses in SN explosions, so that the released energy is all thermalized in the ICM, which is a reasonable approximation if most of the SN explode in a hot and rarefied medium. Finally, the feedback used in the semi–analytical modeling of galaxy formation is not the same as that adopted for heating gas particles in the simulations. For all these reasons, we tend to consider our SN heating scheme as a realistic but approximate recipe to estimate the amount of energy expected from a stellar population, rather than a self–consistent approach to including the effect of SN explosions on the ICM energetics.

## 4 Results

### 4.1 Gas–density and temperature profiles

A unique signature of the lack of self–similarity in galaxy clusters is provided by the different profiles of gas density observed for clusters of different temperatures (e.g. Ponman et al. 1999, Lloyd–Davis et al. 2000), with colder systems having shallower internal profiles. In the upper panels of Figure 7 we show how extra heating can affect the gas profiles. In the absence of any extra heating, all the simulated structures display the same gas density profile. This result agrees with the expectation of the self–similar scaling paradigm and confirms that non–linear gravitational processes, like accretion shocks, do not introduce any characteristic scale. The only structure slightly deviating from self–similarity in the GH runs is the Coma cluster. Looking at its evolution, it turns out that this structure underwent a very recent major merger event, which redistributed the gas in the central regions. Therefore, its slightly softer gas profile is just the signature of this recent merger.

Self–similarity is broken as additional heating is included. SN feedback has some effect only on the gas profile for the Hickson group, our less massive halo. This is in line with the expectation that ICM thermodynamics ought to be changed only when the non–gravitational heating is comparable to the gravitational one. In fact, SN heating provides $`1/3`$ keV per gas particle, which is similar to the mass–weighted temperature of the Hickson group. As for the S-50 runs, the entropy floor corresponds to a large enough heating that only the Coma cluster remains almost unaffected. This result is in line with that based on the shape of entropy profiles presented in Paper I.

A relevant quantity which is directly related to the gas profile is the baryon fraction, $`f_{\mathrm{bar}}`$, within a given radius. The measurement of the baryon fraction in clusters, in combination with constraints on $`\mathrm{\Omega }_b`$ from primordial nucleosynthesis, is considered as one of the fundamental tests for the density parameter $`\mathrm{\Omega }_m`$ (e.g., White et al. 1993, Evrard 1997, Ettori 2002). We show in Figure 8 the radial dependence of $`f_{\mathrm{bar}}`$ for the Virgo and Hickson runs with the different heating schemes. If gas is only subject to gravitational heating, the cosmic baryons fraction is recovered at about 0.3–0.4 $`R_{\mathrm{vir}}`$, which is comparable to the effective scale where $`f_{\mathrm{gas}}`$ is determined from $`X`$–ray data. As for the Virgo cluster, the S-50 scenario would lead only to a $`10`$–15% underestimate of the correct value at about half virial radius. The situation is quite different for the Hickson group, whose gas fraction is diminished out to the virial radius by heating from the entropy floors. This provides a warning about the determination of $`f_{\mathrm{gas}}`$ and the inference of $`\mathrm{\Omega }_m`$ based on poor clusters and groups. The heating from the SN scheme with 400% efficiency (SN-4) produces a $`f_{\mathrm{gas}}`$ profile quite similar to the S-100 scheme within the virial radius. Quite interestingly, the S-100 scheme involves an extra–heating energy, $`E_h2.8`$ keV, which is almost twice as large than for the former scheme. This indicates that a redshift-modulated heating is more efficient, in terms of the required energy budget, than an impulsive pre–collapse heating.

As for the three-dimensional temperature profiles (lower panels of Fig. 7), the gas looks almost isothermal, in the absence of any extra heating, out to $`0.2R_{\mathrm{vir}}`$, while a fairly steep gradient appears at larger radii. Dips in the $`T`$–profiles indicate the presence of merging structures, whose lower virial temperature brings colder gas in the main body of the simulated system, before being stripped by ram pressure and thermalized. Such structures are progressively washed out in the pre–heated runs, as a consequence of the smaller amount of gas that subgroups are able to keep in such cases. The extra heating also steepens the temperature profiles in the internal cluster regions, with a more pronounced gradient for smaller systems. An increase of gas temperature in the central regions is the consequence of suppressing gas density while keeping unchanged the external gas pressure.

In the last few years, observational data have reached good enough precision to allow for spatially–resolved spectroscopic observations of the ICM for fair number of galaxy clusters. Yet, a general agreement about the temperature profiles in the external cluster regions, not affected by cooling, has still to be reached (e.g., Markevitch et al. 1998, White 2000, Irwin & Bregman 2001). Based on Beppo–SAX observations, De Grandi & Molendi (2002) analyzed projected temperature profiles for a set of 21 clusters. They found the gas to be isothermal out to $`0.2R_{180}`$, with a fairly steep negative gradient at larger radii. In Figure 9 we compare result from our simulations to the average profile by De Grandi & Molendi (2002). Since their clusters are all quite hot ($`T>\text{ }4`$ keV), we use only results from the Coma run. The simulation results are obtained by averaging over the temperature profiles projected along three orthogonal directions. The simulation produces a fairly steep gradient at $`R>\text{ }0.2R_{180}`$, quite similar to the observed one. Forthcoming data from Chandra and Newton–XMM (e.g. Arnaud et al. 2001, Allen, Schmidt & Fabian 2001) will certainly refine the determination of temperature profiles in central cluster regions, where cooling is also expected to play a significant role. There is no doubt that such observational advances will challenge the ability of numerical simulations to treat the relevant physical processes.

### 4.2 The mass–temperature relation

The $`M`$$`T`$ relation represents a key ingredient when comparing observational data on the cluster $`X`$–ray temperature function with the mass function predicted by cosmological models. Quite recently, the different results obtained by different authors on the power–spectrum normalization, $`\sigma _8`$, from cluster XTF and XLF have lead to a discussion on the correct $`M`$$`T`$ relation to be used (e.g. Pierpaoli, Scott & White 2001, Borgani et al. 2001b, Seljak 2002, Viana, Nichol & Liddle 2002, Ikebe et al. 2002). In past years, different analyses have shown that eq. (3) provides a good fit to results of simulations including only gravitational heating, with values of $`\beta `$ ranging within the interval $`0.9<\text{ }\beta <\text{ }1.3`$ (e.g. Evrard, Metzler & Navarro 1996, E96 hereafter, Bryan & Norman 1998, Frenk et al. 1998). In particular, E96 pointed out that temperature measurements provide a quite accurate determination of the mass at radii where the mean cluster density is 500–$`2500\rho _{\mathrm{crit}}`$.

Besides investigating the effect of non–gravitational heating on the $`M`$$`T`$ relation, the excellent resolution of our simulations allows to better resolve the structure of cluster cores, which is relevant for determinations of the emission–weighted temperature. Results on the $`M_{\mathrm{vir}}`$$`T_{mw}`$ relation from our simulations (left panel of Figure 10) closely follow the relation of hydrostatic equilibrium of eq.(3) with a very small scatter and almost no dependence on non–gravitational heating. Although extra heating adds internal energy to the gas particles, they quickly get back in hydrostatic equilibrium, with their temperature determined by the depth of the gravitational potential well which is established by the dynamically dominant DM component.

As for observational data, the $`M`$$`T`$ relation has been measured in the last few years by different groups using ASCA (Horner et al. 1999, Nevalainen et al. 2000, Finoguenov et al. 2001, F01 hereafter) and Beppo–SAX data (Ettori, De Grandi & Molendi 2002). Such analyses consistently find that: (a) $`MT^{3/2}`$ for $`T>\text{ }4`$ keV, but with a normalization significantly lower than that found by E96 from simulations (cf. Ettori et al. 2002); (b) a steeper slope for colder systems, possibly interpreted as an effect of pre–heating. We compare in the right panel of Figure 10 data on the $`M_{500}`$$`T_{ew}`$ relation by F01, which include data on systems down to $`T_{ew}<\text{ }1`$ keV, to results from our simulations and to the relation found by E96. The somewhat lower normalization with respect to the E96 results on the group scale is likely to be due to our improved resolution. However the effect of extra heating is at most marginal and not sufficient to reconcile simulations to data, thus indicating that some other physical process should be at work in establishing the $`M`$$`T`$ scaling.

Finoguenov et al. (2001, F01 hereafter) suggest that the difference between observed and simulated $`M`$$`T`$ relation is due to the combined effect of pre–heating and the effect of formation redshift on the temperature of the system. However, our simulations show that pre–heating has a minor impact on this relation. As for the effect of formation redshift, it may introduce a bias in the definition of the observational data set: observations could tend to select fairly relaxed systems, which formed at higher redshift and, therefore, are characterized by a somewhat higher temperature at a fixed mass (e.g. Kitayama & Suto 1996, Voit & Donahue 1998). Our simulations define an $`M`$$`T`$ relation with a small scatter, thus suggesting that differences in the formation epoch or differences in the current dynamical status among systems should have a small effect. A larger set of simulated clusters would be required to properly address this point. Ettori et al. (2002) detect a segregation in the $`M`$$`T`$ relation for cooling–flow and non cooling–flow clusters, the latter being characterized by a larger scatter.

TN01 showed that their ICM model, which incorporates the effects of pre–heating and cooling, reproduces the observed $`M`$$`T`$ relation. They also find that the predicted relation is weakly sensitive to the value of the pre–collapse entropy floor. This suggests that cooling should be responsible for the lower normalization of the relation, through the steepening of the temperature profiles in cluster central regions. The effect of cooling on the $`M`$$`T`$ relation will be further discussed in Section 5 below.

### 4.3 The luminosity–temperature relation

The observed relation between bolometric luminosity and temperature is considered a standard argument against the self–similar behavior of the ICM. Bremsstrahlung emissivity predicts $`L_XM\rho _{\mathrm{gas}}T^{1/2}`$. Therefore, as long as clusters of different mass are scaled versions of each other, then the $`M`$$`T`$ scaling from hydrostatic equilibrium gives $`L_XT_X^2(1+z)^{3/2}`$ or, equivalently $`L_XM^{4/3}(1+z)^{7/2}`$ for $`\mathrm{\Omega }_m=1`$ (Kaiser 1986, see Eke et al. 1998, for an extension to low–$`\mathrm{\Omega }_m`$ cosmologies). As we also discussed in the introduction, this prediction is at variance with respect to observational evidence of a steeper relation, $`L_XT^3`$ for $`T>\text{ }2`$ keV and, possibly, even steeper for colder systems. This result is also in line with the observed slope of the $`L_X`$$`M`$ relation, $`L_XM^\alpha `$ with $`\alpha 1.8\pm 0.1`$ (Reiprich & Böhringer 2002).

The first determinations of the $`L_X`$$`T`$ relation for clusters showed that it has a quite large scatter (e.g. David et al. 1993, White, Jones & Forman 1997). A significant part of this has been recognized to be the effect of cooling: central spikes associated with cooling regions provide a large fraction of total $`X`$–ray luminosity, so that differences in the cooling structure among clusters of similar temperature induce a spread in the corresponding $`L_X`$ values. After correcting for this effect, different authors (Allen & Fabian 1998, Markevitch 1998, Arnaud & Evrard 1999) were able to calibrate a much tighter $`L_X`$$`T`$ relation. Allen & Fabian (1998) also noticed that the $`L_X`$$`T`$ relation for the hottest clusters tends to flatten to the self–similar scaling prediction, possibly suggesting that the effect of extra heating becomes negligible for such systems.

Figure 11 shows the effect of extra heating on the $`L_X`$$`T`$ relation of nearby clusters. The GH runs are clearly at variance with respect to data over the whole sampled range of temperatures. Heating with the SN-1 recipe provides some steepening of the $`L_X`$$`T`$ relation at the group scale, although not enough to reach agreement with data points. A larger suppression of $`L_X`$ is achieved with the S-50 and S-100 heating schemes and, for the Hickson group, with the SN-4 run. Quite interestingly, the SN-4 heating scheme requires the same heating energy, $`E_h1.4`$ keV/part as the S-50 entropy floor, but provides a significantly smaller luminosity, as a consequence of the lower gas density in the central region of the Hickson group (see also Fig. 8). Therefore, a better efficiency is reached in this case by gradually dumping energy within the virialized regions of the ICM, rather than imposing an impulsive pre–collapse heating on the whole turn–around region.

In a similar way, Figure 12 shows a similar comparison with observations for distant, $`0.5<\text{ }z<\text{ }1.3`$, clusters. Independent analyses have confirmed that data on $`L_X`$$`T`$ relation of distant clusters are consistent with a lack of evolution (e.g. Mushotzky & Scharf 1997, Donahue et al. 1999, Della Ceca et al. 2000, Borgani et al. 2001b, Stanford et al. 2002, Holden et al. 2002). Since no data are available on $`T<\text{ }1`$ keV groups in the above $`z`$–range, we do not include in this comparison results for the Hickson runs. Although temperature determinations for distant clusters are prone to larger errorbars, the results are quite in line with what shown in Fig. 11: the $`L_X`$$`T`$ relation of distant clusters requires pre–heating of the ICM with $`E_h>\text{ }1`$ keV/part.

### 4.4 The effect of changing the epoch of heating

Imposing an entropy floor at some pre–collapse redshift is a useful approximation, in that it allows the characterization of the ICM evolution by a single quantity which is conserved in adiabatic processes. As already mentioned, the heating energy required to establish a given entropy floor is a non–trivial function of the heating redshift, $`z_h`$, through the evolution of $`\overline{\rho }`$ and $`\delta _{\mathrm{gas}}`$ (see eq.9). In general, one may also expect the resulting global ICM properties to depend on $`z_h`$. In their simulations, Bialek et al. (2001) heated at a very high redshift, $`z_h20`$, and found results on the required $`S_{\mathrm{fl}}`$ to reproduce the observed $`L_X`$$`T`$ relation which are generally consistent with ours. This may suggest that the choice for $`z_h`$ is not relevant. However, a close comparison between their runs and ours is not straightforward: mass and dynamical resolutions are higher in our runs and, as discussed in Section 2, this is likely to affect observable quantities, such as $`L_X`$, which are sensitive to the details of gas clumpiness.

As shown in Figure 13, the required $`E_h`$ for the S-50 Fornax run changes by about a factor 2, $`z_h=3`$ being the less efficient choice in terms of energy budget. Both $`T_{ew}`$ and central gas entropy at $`0.1R_{\mathrm{vir}}`$ are quite stable, with variations of at most $`10\%`$ and $`20\%`$ respectively. As expected, the $`X`$–ray luminosity is the most sensitive quantity to changes of $`z_h`$; it decreases by about a factor two as the heating epoch is pushed back in time from $`z_h=1`$ to 5. Such relatively modest variations suggest that different choices for $`z_h`$ induces only moderate changes in the value of $`S_{\mathrm{fl}}`$ required to reproduce the observed $`L_X`$$`T`$ relation.

## 5 The effect of cooling and star formation

The simulations discussed so far have been performed with the aim of understanding in detail the effect of non–gravitational heating on ICM observable quantities. However, excluding cooling in our simulations is clearly a simplification, and its inclusion is likely to change some results significantly. Cooling times can be significantly shorter than the Hubble time in central regions of clusters, thus causing a significant fraction of the gas there to cool down and drop out of the hot diffuse, $`X`$–ray emitting phase (see Fabian 1994, for a review on cooling processes in clusters).

Cooling has been suggested as an alternative to non–gravitational heating for breaking ICM self–similarity. As gas undergoes cooling in central cluster regions, $`X`$–ray luminosity is suppressed as a consequence of the reduced amount of diffuse hot gas (Pearce et al. 2000), thus possibly producing a steepening of the $`L_X`$$`T`$ relation (Muanwong et al. 2001, Wu & Xue 2002). The decreased pressure support causes the recently shocked, higher entropy external gas to flow in, thus causing a net increase of entropy for the hot gas left in central cluster regions (Bryan 2000). Recently, Voit & Bryan (2001) have argued that only gas with entropy in excess of 100–200 keV cm2 has long enough cooling time to remain in the diffuse phase. Accordingly, the entropy excess in central regions of poor groups should be interpreted as an effect of cooling, rather than of extra heating.

In order to check this effect, we run a high–resolution simulation of the Hickson group, with the SN-1 heating scheme, but also including gas cooling. The cooling function adopted assumes gas of primordial composition, with zero metallicity, and follows the recipe of Wadsley & Bond (1997). We also include heating by UV background by adopting the latest version of the model by Haardt & Madau (2001). We follow the recipe by Katz, Hernquist & Weinberg (1992) to convert dense cold gas particles into collisionless stars. Furthermore, when computing properties of the diffuse ICM, we exclude gas particles with temperatures $`T<3\times 10^4`$ K and overdensities $`\delta >500`$ (see also Kay et al. 2002). These particles are assumed to belong to the cold medium, although not yet removed from the gas phase by the star formation (SF) algorithm. The effect of including cooling and SF is that of collecting gas into dense knots, while suppressing gas density in the diffuse filamentary structures (see Figure 14). The resulting density profile in central regions turns out to be steeper, $`R^\alpha `$, with $`\alpha 2`$ than that for a CDM only run (see also Treu & Koopmans 2002).

In Figure 15 we show the resulting entropy profile, as compared to the GH and SN-1 runs, with no cooling. We detect some increase of gas entropy inside $`R_{vir}`$, due to the infall of shocked external gas. However, the inclusion of cooling and SF do not produce any entropy plateau at $`100`$ keV cm2. Although the entropy increase goes in the right direction, it may be marginal in accounting for the observed entropy excess in central group regions (e.g. Ponman et al. 1999; see also Paper I). A higher entropy could be possibly attained by assuming a non–negligible metallicity in the cooling function, which would cause an increasing efficiency of gas cooling.

Besides suppressing the amount of hot $`X`$–ray emitting gas, cooling also causes the formation of high–density concentration of star particles in the innermost cluster region. Such concentrations deepens the potential well and causes temperature there to significantly increase. Therefore, while the $`X`$–ray luminosity decreases, the emission–weighted temperature increases. The effect of cooling in the $`M`$$`T_{ew}`$ and $`L_X`$$`T`$ relations of the Hickson run is shown with the filled square in Figs. 9–11. Such changes actually go in the right direction of improving $`L_X`$$`T`$ relation (see Figure 11) and possibly solve the discrepancy with the observed $`M_{500}`$$`T_{ew}`$ relation (see also Thomas et al. 2002; Voit et al. 2002).

However, cooling in itself is a runaway process, due to the quick increase of cooling efficiency with local gas density. As a consequence, if not counteracted by some feedback process, a too large fraction of gas can cool down to a collisionless phase (e.g., Prunet & Blanchard 1999, Lewis et al. 2000, Davé et al. 2001, Balogh et al. 2001). We see in our simulation that while the overall baryonic fraction inside $`R_{vir}`$ at $`z=0`$ coincide with the cosmic value, 37% of the gas is locked in the stellar phase. This large fraction is far in excess of the $`<\text{ }10\%`$ fraction indicated by observations (e.g., Balogh et al. 2001). Therefore, our SN feedback scheme falls short in preventing the cooling runaway.

Addressing in detail the issue of preventing overcooling through a suitable feedback mechanism is beyond the scope of this paper. Due to the biased nature of halo formation in hierarchical models, the first massive halos to collapse into a protocluster region will eventually merge to form its core. These massive halos host the progenitors of giant ellipticals in today clusters (Governato et al. 2001, Springel et al. 2001) and will likely contribute to most of the energy released in the ICM. Whatever the astrophysical source for this energy is, it should act in a self–regulated way: cooling switches on the feedback mechanism and the subsequent energy release inhibits too much gas to leave the hot phase. Besides SN, also AGN activity has been suggested as a possible way of preventing overcooling, with energy released with a duty cycle of a few $`10^7`$ years (e.g. Binney & Tabor 1995, Ciotti & Ostriker 2001, Quilis et al. 2001, Böhringer et al. 2002). If AGN are responsible for the required feedback, then one should expect to see some signature of their recent activity in the central regions of those clusters whose cooling has been recently switched off. A succesfull and self–consistent inclusion of cooling and feedback in numerical simulations has still to be realized, the main difficulty lying in the treatment of physical processes which involve a large dynamical range. The physics describing AGN activity or SN explosions involve small scales that are out of the reach of current and foreseeable simulations of galaxy clusters, thus leaving as the only possibility the inclusion of such process as external recipes.

## 6 Conclusions

We presented results from high resolution Tree+SPH simulations of galaxy clusters and groups aimed at studying the effect of non–gravitational heating of the intra–cluster medium (ICM) on observable $`X`$–ray properties of galaxy systems. Four simulated structures have been chosen so as to encompass the temperature range from a poor Hickson–like group to a rich Coma–like cluster. Two different kinds of non–gravitational gas heating have been chosen: (i) impulsive heating to set an entropy floor at some redshift ($`z_h=3`$ as a reference epoch); (ii) gradual heating, which follows the star–formation rate inside the cluster region, as predicted by semi–analytical modelling of galaxy formation, so as to mimic the effect of SN energy feedback. Our main results can be summarized as follows.

In order for simulations to reproduce analytical predictions for the $`L_X`$$`T`$ relation of gas sitting in hydrostatic equilibrium within a NFW potential well, the structure of the gas distribution needs to be resolved with at least 5–$`10\times 10^4`$ gas particles within the virial radius, with a spatial resolution of the order of 1–2% of $`R_{vir}`$.

The observed slope and amplitude of the $`L_{bol}`$$`T_{ew}`$ relation are well reproduced by setting an entropy floor of 50–100 keV cm2. Assuming negligible radiative losses, SN feedback provides $`0.35`$ keV/part and, as such, has a too small effect on the $`L_X`$$`T`$ relation. A four times larger energy, possibly consistent with that obtainable from AGN, is required to reproduce the $`L_X`$$`T`$ relation at the scale of a Hickson group.

In the absence of cooling, the $`M_{500}`$$`T_{ew}`$ relation is found to be higher by $`40\%`$ than that observed, over the whole sampled $`T`$–range, almost independent of the presence and amount of non–gravitational heating. Therefore, non–gravitational heating in itself is not able to account at the same time for the $`L_{bol}`$$`T_{ew}`$ and the $`M_{500}`$$`T_{ew}`$ relations.

The effect of including cooling in the simulation of a Hickson group with SN heating is to suppress $`L_X`$ by about 40% and increase $`T_{ew}`$ by about 30%. Therefore, the effect of cooling goes in the right direction of decreasing the amplitude of the $`L_X`$$`T`$ and $`M_{500}`$$`T_{ew}`$ relations at the scale of groups. However, the runaway nature of cooling is not efficiently counteracted by our SN feedback scheme, so that 37% of the gas drops out of the hot $`X`$–ray emitting phase. Such a large amount of cold gas is in contrast with observational data on the stellar content of galaxy systems.

The results presented in this paper demonstrate that no unique simple recipe exists to explain $`X`$–ray observable properties of galaxy systems. Supernova heating provides a sizeable, but insufficient, fraction of the non–gravitational heating energy required. Further feedback energy from AGN may be required. Yet, the mechanisms to thermalize the energy associated with QSO activity in the diffuse medium has to be understood in detail.

While much further work is needed, our study suggests that a combination of different heating sources (SN and AGN) and cooling will hopefully be able to reproduce both the correct $`L_X`$$`T_{ew}`$ and $`M_{500}`$$`T_{ew}`$ relations, and to avoide gas overcooling. A self–consistent treatment of such complex physics is mandatory if the numerical description of the ICM has to keep pace with the progress in the observational description of the hot baryons in clusters.

## Acknowledgements

Simulations were run at the CINECA Supercomputing Center, with CPU time provided by a grant of the Italian National Council for Astronomy and Astrophysics (CNAA). We thank the referee, Jerry Ostriker, for useful comments. We thank Fabrizio Brighenti and Mark Voit for reading the paper and Sabrina De Grandi and Alexis Finoguenov for providing us with the files of the data points shown in Figs. 9 and 10, respectively. We acknowledge useful discussions with Richard Bower, Stefano Ettori and Colin Norman. FG is supported in part by the Brooks Fellowship.