text:
# Properties of cluster satellites in hydrodynamical simulations

## 1 Introduction

Clusters of galaxies are the largest virialized structures observed in the universe. Mainly composed of dark matter, they are also large repositories of hot plasma, the intracluster medium (ICM). This makes them very bright X-ray sources, mainly through free-free emission, and enables us to observe them out to reasonably high redshifts.

The local abundance of galaxy clusters and its redshift evolution are powerful tests to discriminate between different cosmological models. In fact, in the framework of hierarchical clustering, clusters are expected to form from the repeated merging of smaller, pre-existing, units. Depending on the amount of matter present in the universe - i.e. on the value of the present matter density parameter $`\mathrm{\Omega }_\mathrm{m}`$ \- we expect clusters to form typically at different redshifts: the more matter present in the universe, the later in time clusters are allowed to assemble.

Observations of the ICM portray the cluster structure and its dynamical stage: relaxed systems may show features consistent with the presence of cooling flows, although evidence is still contradictory (e.g. Lewis, Stocke & Buote 2002; Ettori et al. 2002), while merging events are clearly associated to inhomogeneities in the ICM temperature and morphology (e.g. Markevitch, Vikhlinin & Mazzotta 2001; Mazzotta, Fusco-Femiano & Vikhlinin 2002; Markevitch et al. 2002).

The latest generation of X-ray telescopes, like Chandra and XMM/Newton, have revealed a wealth of details in the morphological and spectral properties of clusters, often suggesting the combined action of several physical processes - one example for all - the role of thermal conduction and magnetic fields in producing the observed structure of cold fronts (Vikhlinin, Markevitch & Murray 2001). These detailed observations are also showing a complexity in cluster structure that was not evident in previous observations, and shows that galaxy clusters are often the site of violent collisions between different subunits.

Although the hierarchical clustering model seems to fit this picture perfectly, there are many uncertainties and open questions on the detailed build-up of matter. Mentioning only the simplest: how often do clusters undergo major mergings? How long does it take for the accreted substructure to relax and virialize in the new gravitational potential? What is the effect of the collision on the internal structure of sub-clumps, and on the overall observables of the system? In order to address these questions one needs to recur to numerical simulations that mimic the formation and evolution of galaxy clusters in a given cosmological model, and analyze the merging history of the forming systems as they grow in time. While a fairly extended literature has been published on the evolution of the dark matter component in the assembly of galaxy clusters (Tormen, Diaferio & Syer 1998; Wechsler et al. 2002; Zhao et al. 2003; Taffoni et al. 2003), surprisingly little work exists which describes the evolution of the ICM counterpart in a cosmological framework. The main reason for this is that any study of the ICM requires performing N-body hydrodynamical simulations; now, producing a sample of simulated clusters with sufficient spatial and mass resolution to follow the details of the ICM evolution is still a demanding task, both for the practical setup (choosing a cosmological simulation or performing several re-simulations, as described below) and for the rather heavy computational requirements. In fact, while there has been a lot of investigation of the physics of merging cluster using simulations, (e.g. Roettiger, Burns & Loken 1996; Ricker & Sarazin 2001; Ritchie & Thomas 2002), these studies mostly used non-cosmological initial conditions in order to address selected and specific theoretical and observational issues, and disentangle them from the convolution of cosmological merging history. In this paper we use the somewhat orthogonal viewpoint of analysing the history and properties of cluster collisions in a self-consistent cosmological environment. This approach guarantees self-consistency of the cluster structural properties (like profiles and shape) and makes it easy to give a statistical meaning to our results.

An important step towards the possibility of realising such a program has been taken by the public release of the Tree-SPH code GADGET (Springel, Yoshida & White 2001), which has soon established itself as one of the reference cosmological codes for the simulations of large-scale structure and of single systems. In this paper we combine the use of GADGET and the technique of single system re-simulation (described below; see also Tormen, Bouchet & White 1997) to produce a sample of 17 high-resolution hydrodynamical clusters.

Whereas in principle simulations should contain as much physics as is reasonable to model, we decided to start from the simplest possible case: that of a non-radiative gas. That means that our simulations model processes like dissipation and shock heating, but do not include radiative cooling, nor more complex processes like star formation or energy feed-back (e.g. Yoshida et al. 2002; Muanwong et al. 2002; Tornatore et al. 2003). This choice is justified by several reasons: first, runs with cooling and extra physics are much more demanding in terms of computational costs, making it much harder to obtain a sample as large as ours at the same resolution. Second, the practical implementation of a simulation becomes more and more uncertain as extra physics is added. In particular, as of today no real agreement has been reached on the effect of, e.g. cooling on the global properties of clusters nor on the numerical convergence of the results (e.g. Balogh et al. 2001). Finally, we deem it necessary to proceed by steps, and the step immediately next to pure N-body simulation is the one we have taken of a non-radiative gas. Proceeding by steps will hopefully allow a better interpretation and understanding of the results and should guarantee that all the physical and numerical effect of the simulation are under control.

The outline of the paper is as follows. In Section 2 we present the simulations and procedures used to create the halo catalogues and merging history trees. In Section 3 we report our results on the evolution of orbital and internal properties of the satellites after merging. In Section 4 we present a few analytical fits modelling the previous results. In Section 5 we apply our findings to a real observation, while in Section 6 we summarize and present our conclusions. In the Appendix we give explicit formulae for the analytic fits presented in the paper, and give some scaling laws useful to translate our results from rescaled to physical units.

## 2 Method

### 2.1 Simulations

#### 2.1.1 Parent simulation

The clusters considered in this paper were obtained using the technique of re-simulating at higher resolution patches of a pre-existing cosmological simulation, as described in Tormen et al. (1997). In this instance the parent simulation is an N-body run with $`512^3`$ particles in a box of side $`479h^1\mathrm{Mpc}`$ (Yoshida, Sheth & Diaferio 2001, see also Jenkins et al. 2001). It assumes a flat universe with a present matter density parameter $`\mathrm{\Omega }_\mathrm{m}=0.3`$ and a contribution to the density due to the cosmological constant $`\mathrm{\Omega }_\mathrm{\Lambda }=0.7`$. The value of the Hubble constant (in units of 100 km/s/Mpc) is $`h=0.7`$. The initial conditions correspond to a cold dark matter power spectrum which is normalized to have at the present epoch a r.m.s. in sphere with radius of $`8h^1\mathrm{Mpc}`$ equal to $`\sigma _8=0.9`$. The particle mass is $`6.8\times 10^{10}h^1\mathrm{M}_{}`$, meaning that cluster-sized haloes are resolved by several thousand particles; the gravitational softening is $`30h^1`$ kpc. From the output of this simulation corresponding to the present time we randomly extracted 10 spherical regions of radius between 5 and $`10h^1\mathrm{Mpc}`$, each containing either a cluster-sized dark matter halo or a pair of haloes.

#### 2.1.2 Generation of new initial conditions

For each of these regions we built new initial conditions using a software package developed by one of us (GT), and named ZIC (for Zoomed Initial Conditions). A summary of the procedure follows. The position of the particles in the initial conditions of the run define a Lagrangian region where the initial density field was resampled using a higher number of particles - on average $`10^6`$ dark matter (DM) high resolution (HR) particles - in order to increase spatial and mass resolution. On top of each HR particle we placed a gas particle with identical velocity; we assumed a baryon density parameter $`\mathrm{\Omega }_\mathrm{b}=0.03`$. The mass resolution of the resimulation ranges from $`2\times 10^9`$ to $`6\times 10^9h^1\mathrm{M}_{}`$ per DM particle. The gravitational softening is given by a $`5h^1`$ kpc cubic spline smoothing for all HR particles.

We drastically reduced the number of particles outside the HR Lagrangian region by interpolating them onto a spherical grid centred on the geometrical centre of the HR region. The angular resolution of the grid was taken between 3 and 5 degrees, corresponding to $`50,000`$ to $`150,000`$ low resolution (LR) particles of varying mass and gravitational softening. Extensive testing has shown that this reduction still guarantees a sufficiently correct description of the tidal field on large scales for most cases.

The HR distribution in the initial conditions contains all the fluctuations of the matter power spectrum realization of the original cosmological run, plus a new and independent realization of high frequency fluctuations from the same spectrum, in order to extend the spectrum up to the Nyquist frequency of the HR distribution.

The ZIC package has been extensively tested and used to generate initial conditions for many resimulations at medium to extremely high resolution (e.g. Tormen et al. 1997; Springel et al. 2001b; Yoshida et al. 2001b; Stoher et al. 2002).

#### 2.1.3 Resimulations

The initial conditions were evolved using the publicly available code GADGET (Springel, Yoshida & White 2001) from a starting redshift $`z_{\mathrm{in}}=3560`$ \- depending on the run - to $`z=0`$; during the runs we took 51 snapshots equally spaced in $`\mathrm{log}(1+z)`$, from $`z=10`$ to $`z=0`$. This choice gives a typical time spacing of $`dt0.5`$ Gyr between successive outputs for $`z<0.8`$.

### 2.2 Halo catalogues

We adopted the spherical overdensity criterion to define collapsed structures in the simulations. For each snapshot of each resimulation we first estimated the local dark matter density at the position of each particle, $`\varrho _{i;\mathrm{DM}}`$, by calculating the distance $`d_{i,10}`$ to the tenth closest neighbour, and assuming $`\varrho _{i;\mathrm{DM}}d_{i,10}^3`$. We then sorted the particles in density and took as centre of the first halo the position of the densest particle. Around this centre we grew spherical shells of matter, recording the total (i.e. DM + gas) mean overdensity inside the sphere as it decreases with increasing radius. We stop the growth and cut the halo when the overdensity first crosses the virial value appropriate for the cosmological model at that redshift; e.g. $`\varrho (<r)/\varrho _\mathrm{b}=323`$ at $`z=0`$. The particles selected in this way belong to the same halo and are used to compute its virial properties (mass, radius, etc.). For the definition of virial overdensity we adopted the model of Eke, Cole & Frenk (1996). We tagged all halo particles as engaged in the list of sorted densities, and selected the centre of the next halo at the position of the densest available particle. We continue in this manner until all particles are screened. We include in our halo catalogue only haloes with at least $`n=10`$ DM particles inside the virial radius. All other particles were considered field particles.

Occasionally it may happen that large haloes host small knots of dark matter which are denser than the main halo itself. In such cases the density maximum would not pick the correct position as putative halo centre. To prevent this kind of errors, whenever a large halo (defined as one with more than 10,000 particles) was found, we obtained a second centre estimate as the converged centre of mass of spheres of decreasing radius (sometimes called the moving centre method), and chose - between the two centres - the one for which the grown halo was more massive.

Another problem originating during the build-up of the halo catalogue is when a newly found halo intersects one or more haloes previously found, meaning that the distance between the centres is less than the sum of the two virial radii. We deal with such cases as follows: whenever a halo $`A`$ is found, we check if it intersects with any other haloes, and undo all intersecting haloes with mass $`M<M_A`$ making their particles available for the build-up of new haloes. We then recalculate the mass of halo $`A`$ considering also the particles disengaged by this procedure. This procedure ensures that - in the final catalogue - particles located inside the intersections are always assigned to the most massive of the intersecting haloes.

At redshift $`z=0`$ this procedure identified 17 massive haloes, containing on average $`N_V200,000`$ dark matter particles within their virial radius, and a similar number of gas particles. The corresponding virial masses are in the range $`M_V=3.1\times 10^{14}`$ to $`1.7\times 10^{15}h^1\mathrm{M}_{}`$.

### 2.3 Merging history tree

For each of these cluster-sized haloes we built a merging history tree using the halo catalogues at all time outputs - separated by redshift intervals $`dz_i`$ \- as follows. Starting with a halo at $`z=0`$, we define its progenitors at the previous output $`z=dz_1`$ to be all haloes containing at least one particle that by $`z=0`$ will belong to the first halo. We call main progenitor at $`z=dz_1`$ the one giving the largest mass contribution to the halo at $`z=0`$. We then repeat the procedure, starting at $`z=dz_1`$ and considering progenitors at $`z=dz_1+dz_2`$, and so on backward in redshift, always following the main progenitor halo. In this merging history tree we term satellites the progenitors which - at any time - merge with the main progenitor. In practice we started from the progenitor list at a given simulation output and selected those haloes for which at least one particles was found in the main progenitor at the next time output.

The actual merging redshift $`z_{\mathrm{mer}}`$ of a satellite with the main progenitor - which we define as the redshift when the satellites first crosses the virial radius of the main halo - is somewhere in-between the redshift $`z_{\mathrm{last}}`$ \- when the satellite is still found as an individual halo in the catalogue - and the redshift of the next simulation snapshot, $`z_{\mathrm{last}}dz`$. We formally assigned $`z_{\mathrm{mer}}`$ to each satellite by choosing a random value uniformly distributed between these two redshifts. The evolution of the satellite orbital and dynamical properties is then studied as a function of the time $`\tau `$ elapsed after its merging: $`\tau t(z)t(z_{\mathrm{mer}})`$, following it through the simulation down to $`z=0`$.

Hereafter we will present results with times in Gyr, while all other quantities will be expressed in the appropriate units of $`h`$ for easier comparison with analogous simulation results (and especially with Tormen, Diaferio & Syer 1998), although remember that $`h=0.7`$.

To get a rough idea of the merging timescale, for a galaxy cluster observed at $`z=0`$ the following relation between lookback-time and redshift holds: $`t_{lb}=1,2,4,6`$ Gyr correspond to $`z=0.075,0.16,0.37,0.65`$ for our cosmological model.

## 3 Time evolution of the satellites

### 3.1 Merging of satellites

In hierarchical clustering the formation of cosmic structures begins at high redshifts with the smallest resolved systems and proceeds to larger ones by repeated merging events. In this picture galaxy clusters are located at the highest level of hierarchy, and gather all the matter accreted through their merging history. However, numerical simulations show that this is only approximately true. The violence of mergers, the variety of orbits and the complexity of the tidal field make the merging history of the collapse a more complex process, so that in general only a fraction $`f_cm_{\mathrm{cap}}/m_V`$ of the satellite mass $`m_V`$ is actually accreted by the main halo; as a consequence, the sum of all progenitor masses exceeds the mass of the final halo by approximately 20 per cent, and the excess mass is finally redistributed in debris within roughly 2 times the halo virial radius (Tormen 1998).

#### 3.1.1 Selection criteria

In our analysis we will consider - without loss of generality - only satellites donating at least half of their virial mass to the main halo (i.e. $`f_c>0.5`$). This criterion ensures that satellites contributing only marginally to the merging history are excluded.

Once we have identified all the satellites that merge onto the main halo progenitor at all redshifts, we can investigate how each satellite evolves by following in time the fate of the $`n_V`$ particles contained inside its virial radius at the last time output before merging happens. That is, we can focus on how substructures evolve in time inside a cluster.

We applied two more selection criteria:

* we only consider satellites with $`n_V100`$ dark matter particles to minimize effect due to numerical resolution;
* we limit our study to satellites merged after $`z_{\mathrm{last}}=0.8`$.

The last selection criterion was applied after we realized that the orbital properties of the satellites, when expressed as function of the time $`\tau `$ after merging, depended somewhat on the merging redshift for $`z_{\mathrm{last}}>\text{ }0.8`$, going in the direction of smaller dynamical timescales for higher values of $`z_{\mathrm{last}}`$; as a simple rescaling of $`\tau `$ did not remove this time dependence, we decided to exclude those data from our study.

Since in our simulations the average mass of the main halo progenitors is $`M_V(z=0.8)3\times 10^{14}h^1\mathrm{M}_{}`$, and by the end it grows to $`M_V(z=0)8\times 10^{14}h^1\mathrm{M}_{}`$, the selection in $`z_{\mathrm{last}}`$ has the side advantage that our results will always refer to mergers onto cluster-size haloes, and so aquire a more physical meaning in any comparison with observations. In total, there are roughly 1,200 satellites satisfying the above requirements, on which we will concentrate for our following investigation.

#### 3.1.2 Definition of satellite mass

As a satellite comes close to the main halo, dynamical interaction within the mutual gravitational field will modify the orbital and internal properties of the satellites. In particular, tidal stripping will start unbinding the external regions of the satellite, as they gradually get heated to the virial temperature of the larger system. This becomes more dramatic as the satellite crosses the virial radius of the main halo at $`z_{\mathrm{mer}}`$.

In order to calculate the satellite’s properties after $`z_{\mathrm{mer}}`$ we need thus to refer to the part of the satellite that - at each time $`\tau `$ \- is still self-bound and resisting against disruption. We calculated the self-bound fraction of a satellite at each $`z<z_{\mathrm{mer}}`$, that is, for $`\tau >0`$, by the following iterative procedure:

* we locate the satellite densest part by the moving centre method and take the resulting position as the satellite centre. We then make a list of all satellite particles inside $`r_V`$, the virial radius of the satellite at $`z_{\mathrm{last}}`$, assuming that all particles outside this radius are physically unbound to the satellite;
* we calculate potential and kinetic energy for each selected particle; velocities are taken in the reference frame of the ten densest particles; for gas particles we add the internal energy to the balance;
* we remove from the list all particles with positive total energy;
* we repeat steps (ii) and (iii) until the self-bound mass converges to some value $`m_{\mathrm{sb}}(\tau )`$ to better than 10 per cent. Obviously the relation holds $`m_{\mathrm{sb}}(\tau )m_V`$. This procedure returns also the converged position and velocity of the bound part of the satellite.

Each satellite donates both dark matter and gas to the main halo. However, while the dark component is non-collisional and in principle can freely stream through the cluster centre, ram pressure and dissipative effects should slow the gas down and bring it to hydrostatic equilibrium. Therefore, as the satellite dives into the main halo, its ICM is stripped out of the dark matter potential well and the two components physically separate. The time it takes for this to happen depends on the mass of the satellite relative to that of the main halo.

Now, as long as dark matter and gas are united, the orbital properties of a satellite should be those of the two components together, while we should follow dark matter and gas as individual systems after they physically separate. To discriminate between the two regimes - united and separated - at each time $`\tau `$ we calculate the distance between the dark matter and gas centres, and assume separate components whenever this distance exceeds $`r_V`$, the virial radius of the satellite at $`z_{\mathrm{last}}`$. This distinction has the practical consequence - for example - that the potential and kinetic energy of each particle is calculated considering dark matter and gas as one system with a common velocity reference frame before separation, but as two independent particle systems with different velocity frames after separation.

### 3.2 Survival of satellites

With our apparatus up and running, we can start asking interesting questions on the mean properties of merged satellites, like how long - on average - a satellite remains gravitationally self-bound before being destroyed and digested by the main halo, what is the rate of this mass loss, and so on. Here and below we will show results for satellites in four different mass ranges: $`m_V/M_V<0.001`$, $`0.001m_V/M_V<0.01`$, $`0.01m_V/M_V<0.1`$ and $`m_V/M_V0.1`$; these intervals will roughly correspond to mergers of - respectively - small galaxies, massive galaxies, small groups, massive groups or small clusters, onto a massive cluster of mass $`M_V`$.

In Fig. 1 we show the evolution of the bound mass fraction for the dark and gas component, and for each mass range. To take into account the fact that only a fraction $`f_c`$ of the satellite’s mass has actually merged, self-bound fractions $`f_{\mathrm{surv}}`$ are re-normalized in units of $`f_cm_V`$. Dots show the median values, while the bands enclose the second and third quartiles of the distributions.

The different behaviour of the two components is evident: dark matter survives longer in smaller satellites, while the opposite holds for the gas. Small, galactic-size objects (upper left panel), maintain a significant fraction of their initial mass bound to themselves for many Gyr; on the other hand their ICM is stripped practically as soon as they cross the virial radius of the main halo. As we move to larger masses (upper right and lower left panels) tidal stripping and dynamical friction become more and more effective in unbinding dark matter from the satellites; at the same time, the deeper potential well of the satellites retains a higher fraction of gas for longer times. For mergers of comparable masses (lower right panel) dark matter and gas share a similar fate as they are unbound from the satellites in 2 to 4 Gyr.

This behavior is further understood by considering Fig.2, where we plot the separation $`\mathrm{\Delta }_{\mathrm{CM}}`$ between the centre of their dark matter and gas distribution, in units of the satellite’s virial radius, as a function of $`\tau `$. The vertical arrows indicate the average time when the dark matter and gas components separate - as previously defined. Median separation times for the first three mass ranges are at $`\tau =`$ 0.27, 0.68 and 1.85 Gyr from small to larger masses. Therefore separation happens sooner for small lumps and later for more massive ones, while the most massive satellites (lower right panel) are able to retain their ICM throughout all their orbit. In all cases the separation reaches a maximum value at around $`\tau =3.5`$ Gyr, which roughly corresponds - as we will see below - to the apocentre time for the satellite’s orbit. The maximum separation is larger for smaller satellites.

Looking back at Fig.1, the same arrows show that, as soon as the ICM separates from the dark matter component, it becomes gravitationally unbound. This is understood as follows: before separation most of the gravitational potential of the satellite is provided by its dark matter component, and the gas sits in equilibrium in its potential well. After separating from its dark counterpart, the gas potential well becomes roughly nine times less deep - the ratio of dark to gas total mass; on the other hand the gas internal energy is the same as before, that is nine times larger relative to the new potential: under such conditions the ICM is abruptly unbound for all satellites. The most massive satellites are those where the satellite gas sticks with its dark matter counterpart for longer; therefore their ICM retains its original properties for a longer time.

### 3.3 Orbital properties of satellites

In this section we will investigate the orbital properties of satellites after they cross the virial radius of the main halo at $`\tau =0`$. We will consider the evolution of orbital distances and velocities for all the particles initially associated to a satellite, and also for the self-bound part of each satellite. Both cases are of interest: the former gives indications on how the satellite matter is redistributed in the new halo, while the latter tells how the core of a satellite moves through the main cluster.

#### 3.3.1 Mean orbital properties

In Figs. 3 and 4 we illustrate the average orbital distance and velocity for all particles belonging to merged satellites, as a function of $`\tau `$. The figures confirm the different behaviour of the DM and gas particles for small satellites: while the DM oscillates in position and velocity, the ICM steadily floats towards smaller and smaller radii, loosing bulk motion. Large satellites show instead a similar trend for DM and gas, as they both slow down and sink towards the centre of the main cluster. Maximum (minimum) orbital velocities are reached at pericentres (apocentres), as expected.

#### 3.3.2 Orbital properties for the self-bound part

The analogue of Figs. 3 and 4 for the satellite self-bound cores is shown in Figs. 5 and 6. Since the figures reflect the behaviour of the self-bound satellites, they illustrate the behaviour of the subset of satellites which - following Fig. 1, have not been disrupted by time $`\tau `$.

The figures show an evolution consistent with the previous results. We observe that the gas follows the same orbit of the dark component for longer and longer times in going from small to massive satellites. In particular, for the smallest satellites (upper left panels of both figures) the gas slows down and separates dynamically already in less than 0.5 Gyr, well before reaching the orbital pericentre. Orbital velocities increase as satellites fly towards the pericentre, then decrease toward the apocentre, and so on. We see that, on average, satellites merged after $`z=0.8`$ have time to complete at most a couple of orbits; for example, in similar size mergers (lower right panel) the satellites can reach their pericentre and bounce back, then finally sink in the core of the main cluster, although by the end they are still not quite at rest in the cluster potential.

As a summary of the last two subsections, we can remark that the dynamical behaviour of dark matter and gas in the satellites is different. The dark matter component can survive and orbit around the cluster for quite a long time, while the diffuse baryonic component is slowed down by the gas already present in the cluster. For the largest satellites these differences are reduced, partially because the dark matter component is quickly dragged toward the cluster centre by dynamical friction, and also because the potential well of the satellites is deep enough to retain part of its ICM.

### 3.4 Internal properties of satellites

In this section we will address the following issue: how do the internal properties of a satellite change after it has merged onto the main cluster progenitor? The question can be further split in two considerations: (i) how long does it take for the satellite to thermalize into the new environment, that is, how many Gyr are required before the total dark matter and gas of a satellite are dispersed and conform to the average internal properties of the main halo? And (ii) how do the self-bound part of a satellite react to its diving into the potential well of the cluster, that is, how will its velocity dispersion, temperature and entropy change?

To answer these questions we looked at the evolution of the dark matter internal velocity dispersion, and of the gas temperature and entropy as a function of $`\tau `$. Again we did this in two ways: in Figs. 78 we consider all satellite particles (bound and unbound) in order to answer the first posed question and see how long it takes for the satellite matter to reach the virial properties of the hosting halo; in Figs. 89 instead we consider only the self-bound satellite particles: this will reveal how the inner part of a satellite - the one that survives stripping - changes its properties in reaction to the merging.

#### 3.4.1 Mean internal properties

As usual we consider four mass bins according to the mass ratio between satellite and main halo. As expected, the satellite’s mean velocity dispersion, temperature and entropy of satellites grows after merging. The dark matter velocity dispersion has different behaviour for small and large satellites. In the upper panels of Fig. 7 (referring to $`m_V/M_V<0.01`$), $`\sigma _v`$ remains roughly constant for $`\tau 4`$ Gyr, then it start growing, without however reaching the virial value. Larger satellites (lower panels) show a more gradual growth and reach values closer to - but still lower than - the mean cluster temperature. Notice that by the end ($`\tau 7`$ Gyr) the average $`\sigma _v`$ of matter in satellites is always lower than that of the mean cluster. This is likely due to the fact that the satellite’s internal velocity dispersion is calculated in a different velocity reference frame (that of the bound part of the satellite) compared to the velocity dispersion of the main halo. Another - more unlikely - explanation is that the results of Fig. 7 do not consider the matter accreted from the field, which may constitute a non-negligible fraction of the total accreted matter. However, we have not verified these hypotheses further.

Considering the temperature, most of the jump in $`T`$ happens in roughly $`2`$ Gyr. After this jump the temperature grows very slowly or stabilizes. Satellites in different mass ranges reach different asymptotic temperatures: compared to the mean temperature $`T_V`$ of the main halo, smaller satellites converge to values of $`T<T_V`$, while the largest satellites reach $`T>T_V`$. This effect is not large, but it does reflect the final spatial distribution of matter in the main cluster: small satellites deposit most of their mass in the outskirts of the cluster because they are stripped of their ICM almost as soon as they cross the cluster virial radius. On the other hand, larger satellites resist stripping and carry part of their mass into the cluster core. Since clusters have a larger temperature in the centre than in the outskirts, the different spatial distribution corresponds to different gas temperatures.

After an initial jump, entropy (solid band in Fig. 8) stabilizes in $`\tau 2`$ Gyr for satellites with $`m_V/M_V<0.01`$ (upper panels), while it increases more gradually for larger satellites. The asymptotic entropy value decreases going from small to large satellites, and is slightly higher (lower) than the mean cluster value for small (large) satellites. This is again expression of the overall different spatial distribution of the gas donated by satellites of different mass, as previously discussed for the temperature.

¿From the overall analysis of Figs. 78 we can conclude that, while the dark matter component of satellites reaches the virial equilibrium of the main halo after a time lag that depends on the satellite mass (determining the dynamical friction time), most of the gas thermalization process happens in the first $`2`$ Gyr after merging, independent of the satellite-to-cluster mass ratio, and of the different history of the gas (which is stripped from small satellites, but carried along with larger ones).

#### 3.4.2 Internal properties for the self-bound part

Let us now turn to the other question we asked at the beginning of this section: what is the behaviour of dark matter and gas for that part of the satellite that remains self-bound? This issue has interesting observational consequencies in the discussion of merging clusters, as the mass of the structures undergoing the merger is estimated either from the velocity dispersion of their galaxies or from their X-ray temperature. It is not a priori guaranteed that these quantities are indeed a measure of the satellite mass in such critical dynamical conditions. The analogue of Figs. 78 for the self-bound part of each satellite is shown in Figs. 8 and 9. Here the behaviour is strikingly different from the one described for Figs. 78. In particular we note that now (i) the starting values of velocity dispersion, temperature and entropy are usually lower than those previously seen; (ii) the evolution of these quantities is also different: they not always grow in time, but often decrease, especially for massive satellites.

We wish to consider these two points in more detail. Assuming the standard scaling laws for velocity dispersion and temperature, $`\sigma _vM^{1/3}`$ and $`TM^{2/3}`$, we expect that satellites with mass ratio $`m_V/M_V=10^4`$, $`10^3`$, $`10^2`$ and $`10^1`$ have velocity dispersion - in units of that of the main cluster - of the order of $`\sigma _v/\sigma _v(\mathrm{main})=0.046`$, $`0.1`$, $`0.215`$ and $`0.464`$ and temperature $`T_V/T_V(\mathrm{main})=0.002`$, $`0.01`$, $`0.0464`$ and $`0.215`$. We stress, however, that these are only rough estimates.

Now, as for point (i), we can see that the values of the velocity dispersion for the self-bound part of satellites (Fig. 9, hatched bands) at merging time, $`\tau 0`$, are indeed consistent with those deduced from the scaling laws. That means that the self-bound satellites are themselves in virial equilibrium at merging time, at least as far as their dark matter component is concerned. The corresponding ICM temperatures (solid bands) are instead systematically hotter than their dark matter counterpart, being a few times larger than the virial estimates. The effect is stronger for smaller satellites. This is telling that the internal ICM structure of satellites is modified by the collision with the main halo already at orbiting distances of the order of the virial radius of the main halo. We stress that this is visible even in the self-bound part of the satellite, which is its central and more shielded part.

Let us compare these results to those of Fig. 7 (mean velocity dispersions and temperatures for all satellite particles: bound and unbound). There at merging time ($`\tau 0`$) both the dark matter and gas components were strongly influenced and heated by the main halo. Mean temperatures in particular were far larger than the virial values, except for the most massive satellites ($`m_V/M_V>0.1`$, lower right panel). The gas behaviour is understood by recalling that small satellites are stripped of their gas content almost immediately after crossing the virial radius of the main halo, as shown in the top panels of Fig. 1.

However, velocity dispersions too were higher than expected from virial relations, although Fig. 1 indicates that at merging time the average self-bound fraction of satellite matter is 80 to 90 per cent. Evidently the missing 10-20 per cent has been heated already to the point of contributing significantly to rising the mean values of Fig. 7 compared to Fig. 9. Entropy shows a similar behaviour: initial values are higher (Fig. 8, hatched band).

Passing on to point (ii), that is the evolution of these quantities in time, we see that the velocity dispersion of the bound part of a satellite decreases in time. This means that, while globally the satellite is heated, its core actually cools down. Looking at the temperature evolution, we find the same behaviour for the most massive satellites (lower panels of Fig. 9). The effect is particularly dramatic for the mass range $`m_V/M_V>0.1`$, where the ICM cools down by almost a factor two; unlike them, smaller satellites are heated also in their self-bound part, though less than they are heated globally. A similar result is finally found for the entropy: it increases in small satellites and decreases in large ones, but it always remains at values lower than the average.

We stress that this result is quite interesting in view of a comparison with observations; it shows that mass estimate of substructures undergoing merging cannot be directly derived from, e.g., their X-ray properties or (assuming that mass traces light) from the velocity dispersion of their galaxies. In Section 5 of the paper we will apply these results to observations of a real merging cluster.

## 4 Modelling the results

In this section we come back to the correlation noticed in the mean behaviour of some of the results presented so far; we will try to model these correlations parametrically, in order to obtain simple analytical fits that could be used to improve the dynamical description of cluster satellites e.g. in semi-analytic models of galaxy formation.

### 4.1 Survived fraction

In Section 3.2 we have investigated the times it takes to unbind dark and gaseous matter from satellites and shed it in the main cluster. Fig. 1 shows a typical exponential decay of the self-bound matter fraction, and it would be nice to be able to model it in more detail. To this end we tried the following form:

$$f_{surv}(\tau )\frac{m_{\mathrm{sb}}(\tau )}{f_cm_V}=\mathrm{exp}\left(\tau /\tau _{\mathrm{dec}}\right),$$
(1)
where we introduce a decay time $`\tau _{\mathrm{dec}}`$ as the only free parameter of this model.

We applied a $`\chi ^2`$ fit to the curves $`f_{\mathrm{surv}}(\tau )`$ of each satellite, and obtained the individual best-fitting value of $`\tau _{\mathrm{dec}}`$. We did this separately for the dark matter and ICM components. In Fig. 10 we show the log-log plot for the result, as a function of the satellite mass ratio $`m_V/M_V`$. The dots show the best-fit $`\tau _{\mathrm{dec}}`$ for each satellite; the symbols represent the median values for dark matter (squares) and gas (circles) and the bands enclose the two central quartiles of the distributions. The behaviour of median curves is well fitted by straight lines (superimposed to the curves and symbols) for all $`m_V/M_V0.5`$, with slope and zero point $`(0.12,0.3)`$ and $`(0.3,0.3)`$ for dark matter and gas, respectively. The scatter in the logarithmic residuals is however quite large, as shown in the figure, being of the order of $`\sigma =0.4`$ for both DM and gas points.

We can immediately translate the above fit into an expression for $`\tau _{\mathrm{dec}}`$ (expressed in Gyr) as a function of the mass ratio $`m_V/M_V`$:

$$\tau _{\mathrm{dec}}(m_V/M_V)=2\left(\frac{m_V}{M_V}\right)^{0.12}$$
(2)
for dark matter and

$$\tau _{\mathrm{dec}}(m_V/M_V)=2\left(\frac{m_V}{M_V}\right)^{0.3}$$
(3)
for gas. These expressions, introduced in the model Eq. (1), give

$`f_{\mathrm{surv};\mathrm{DM}}(\tau )`$ $`=`$ $`\mathrm{exp}\left[{\displaystyle \frac{\tau }{2}}\left({\displaystyle \frac{m_V}{M_V}}\right)^{0.12}\right]`$
$`f_{\mathrm{surv};\mathrm{gas}}(\tau )`$ $`=`$ $`\mathrm{exp}\left[{\displaystyle \frac{\tau }{2}}\left({\displaystyle \frac{m_V}{M_V}}\right)^{0.3}\right]`$ (4)
for the dark matter and ICM components, respectively. The large scatter quoted above should be considered if one is going to use these relations.

### 4.2 Orbital properties

It is interesting to look in more detail at the mean satellite orbits. In particular, a simple study we can do is to consider the mean properties at the first pericentric passage of the satellite, and at the subsequent first apocentric passage. We did this by considering each satellite in turn, and by identifying the first pericentre as the first relative minimum in the satellite orbital distance going from $`\tau =0`$ on. Analogously, the first apocentre was identified as the first relative maximum in orbital distance maximum following in $`\tau `$ the first pericentre. We limit this study to the first peri(apo)centric passages because they are easier to find in an automated way and also because this guarantees a better data completeness (as fewer and fewer satellites reach further orbital stages).

For each satellite we then write down the following quantities: pericentric and apocentric times, distances from the main halo and orbital velocities. In Fig. 11 we show the median values and quartile bands for all these quantities, separately for the dark matter component (square symbols and hatched bands) and for the ICM (circles and solid bands).

Upper panels refer to pericentre quantities, lower to apocentre ones. We can note that - on average - satellites reach their pericentre roughly in 1 Gyr, regardless of their mass. The dark matter component reaches a minimum distance of $`0.3R_V`$, with a slight tendency for large satellites to have more radial and penetrating orbits; these figures agree with the results found by Tormen (1997), in the analysis of much smaller N-body simulations of a scale-free cosmology with $`P(k)k^1`$. At pericentre the satellites reach a maximum speed of the order of twice the circular velocity of the main cluster. As already observed in Fig. 2, in small satellites the ICM is stripped from its dark matter counterpart, and slows down soon after entering the main halo: in fact gaseous (solid bands) pericentric distances are larger (and - correspondingly \- velocities are smaller) for $`m_V/M_V<\text{ }0.01`$.

Apocentres for the dark matter component (hatched bands) are reached in about $`\tau 3`$ Gyr for small ($`m_V/M_V<\text{ }0.01`$) satellites, and gradually decrease to $`\tau 2`$ Gyr for the most massive satellites (bottom left panel). The gas components instead reaches its apocentre in approximately $`2`$ Gyr, regardless of the satellite mass. This again indicates the almost immediate detachment of the ICM in satellites with $`m_V/M_V<\text{ }0.01`$.

The corresponding distances are larger for small satellites, while large satellites stay closer to the centre of the main halo (bottom centre). Apocentric velocities are a decreasing function of satellite mass; this might be counterintuitive, as large satellites have smaller apocentres and one could expect thus larger velocities. Evidently dynamical friction is extremely efficient in slowing the satellites down. The ICM component follows fairly closely the dark matter one, maintaining only slightly smaller apocentric distances and a hint of larger apocentric velocities.

The continuous lines superimposed to the median values of the dark matter component are analytical fits made to model the dark matter results. In the Appendix we write down the expression of these fits for practical use.

## 5 Example Application: the Cluster 1E0657-56

The orbital properties modeled in the previous section can be used as templates to interpret the dynamics of observed merging systems. As an example we choose the cluster 1E0657-56. This is a merging system at $`z=0.296`$, with both X-ray (e.g. Markevitch et al. 2002, 2003, and references therein) and optical (e.g. Barrena et al. 2002) observations. Its average X-ray temperature is of the order of $`15`$ keV, indicating a massive cluster. The Chandra image shows a bullet-like substructure off-centre with respect to the cluster core, probably a galaxy group or small cluster observed in a post-merging phase. The temperature of the bullet is of the order of $`67`$ keV, with a bow-shock feature in front of it. From the jump in temperature and density across the shock front Markevitch et al. (2003) derive a peculiar velocity of $`4500_{800}^{+1100}\mathrm{km}\mathrm{s}^1`$ relative to the main cluster. Assuming that the merger is almost perpendicular to the line of sight, from the X-ray bullet-cluster distance of $`0.4h^1\mathrm{Mpc}`$ (for the cosmology considered in this paper) these authors date the pericentric passage of the bullet at $`0.10.2`$ Gyr ago.

In the optical data of Barrena et al. (2002) the subclump has a l.o.s. galaxy velocity dispersion of $`210\mathrm{km}\mathrm{s}^1`$, and a l.o.s. velocity of $`620\mathrm{km}\mathrm{s}^1`$, while the main cluster has a l.o.s. velocity dispersion of $`1250\mathrm{km}\mathrm{s}^1`$. The distance between the core of the main cluster and that of the bullet is 2.7 arcmin, or $`0.5h^1\mathrm{Mpc}`$. This value is larger than that derived from the X-ray data due to an offset between the galaxies and the ICM. In particular, galaxies in the subclump are $`0.3h^1\mathrm{Mpc}`$ in front of the peak of the X-ray emission, as inferred by Fig.1b of Markevitch et al. (2002). Assuming for the merger plane a 10 degree deviation from the plane of the sky, this gives a three-dimensional velocity of $`3500\mathrm{km}\mathrm{s}^1`$, in rough agreement with the derivation of Markevitch et al. (2002). Barrena et al. (2002) derive a post-merger mass for the bullet of the order of $`10^{13}h^1\mathrm{M}_{}`$, with an uncertainty of a factor of a few. They also suggest that the pre-merger mass of the subcluster could be substantially larger, due to its high X-ray temperature. Note that Markevitch et al. (2003) were able to constrain the dark matter self-interaction cross-section by using four independents methods, based on the observed gas-dark matter offset, the absence of an offset between dark matter and galaxies, the high subcluster velocity and the subcluster survival.

In simulations, dark matter clumps are naturally associated with the galaxy distribution; the distance between the bullet and the cluster derived from the optical data should correspond to the distance between the dark matter structures. Let us interpret the optical data in the light of our simulations. Taking the relations reported in the Appendix, namely eqs.(11) and (12) at face value, from the measured velocity dispersions the dynamical masses of the cluster and bullet should be of the order of $`1.35\times 10^{15}h^1\mathrm{M}_{}`$ and $`1.04\times 10^{13}h^1\mathrm{M}_{}`$, corresponding to a mass ratio $`1:130`$ post merger.

Fig. 8 dates pericentres at $`\tau 1`$ Gyr regardless of mass; taking the Markevitch et al. (2002) datation of 0.1-0.2 Gyr after pericentre, we derive a time after merging of $`\tau 1.11.2`$ Gyr. From Fig. 10 we note that the velocity dispersion of the self-bound part of a subclump (presumably associated with the galaxies associated to it) does not significantly change in the first 1.5 Gyr after merging, regardless of its mass. Therefore the observed ratio $`\sigma _v(\mathrm{bullet})/\sigma _v(\mathrm{main})0.17`$ would place the bullet in a mass intermediate between the second and third panel, i.e. at around 1:100 mass ratio, also at the pre-merger stage.

This points to a first possible model: that of a galaxy group of pre-merger mass of the order of $`10^{13}h^1\mathrm{M}_{}`$ observed at $`\tau 1.11.2`$ Gyr, just after its first pericentric passage through a cluster of mass $`10^{15}h^1\mathrm{Mpc}`$. Fig. 2 shows that subclumps in the mass range around 1:100 (second and third panels) still retain - on average - 70 to 80 per cent of their initial dark matter mass after $`\tau 1.11.2`$ Gyr, while the ICM self-bound mass may range from 0 to 40 per cent. The expected velocity after the first pericentric passage is of the order of $`v_p2V_c3000\mathrm{km}\mathrm{s}^1`$ (Fig. 8), not inconsistent with the speed inferred by X-ray data.

The observed separation between the galaxy and X-ray density peaks of the bullet - estimated from the superimposed X-ray and optical images \- is of the order of $`0.11h^1\mathrm{Mpc}`$ (proper), or $`0.06`$ times the virial radius of the main cluster, which is derived from the cluster mass and the same scaling laws Eq.(11) and (12). Comparing this value with the one found in simulations (Fig. 3) would place the mass ratio of the bullet to main cluster more in the range $`1:100`$ to $`1:10`$ (third panel) than in the range $`1:1000`$ to $`1:100`$ (second panel). However, the most puzzling observation is the temperature of the bullet. From Fig. 10 we see that the expected X-ray temperature of the self-bound ICM should be at most 20 per cent of the cluster temperature, or $`3`$ keV, while the observations suggest a ratio of the order of $`1:2`$. We can reconcile this discrepancy if we accept the idea that the X-ray measurements really refer to the total (bound and unbound) ICM, which from Fig. 8 is shown to be much higher, at or above the observed 7 keV.

Let us consider the alternative picture by Barrena et al. (2002), of a more massive subclump of initial velocity dispersion of $`700\mathrm{km}\mathrm{s}^1`$, corresponding to a virial temperature of 7 keV. In such a case, the initial mass ratio is of order 0.2, and the post-merger mass would be of order of 4 per cent. Now, a subclump this big takes of the order of 4 Gyr to be disrupted (Fig. 2), and the high speed of the bullet would imply that we are observing it at its second pericentric passage (Fig. 5), when the separation between DM and gas is still negligible (Fig. 3). However, from Fig. 10 we see that by this time the DM velocity dispersion is unacceptably high to match the observed $`210\mathrm{km}\mathrm{s}^1`$; equally, no ICM has survived to this time, and the X-ray temperature of the ICM originally in the subclump has reached the virial value of the cluster (Fig. 8).

In the light of these data we suggest that the most likely interpretation is the first one, and that the high X-ray temperature measured for the bullet is due to its ongoing thermalization in the cluster potential.

## 6 Summary and Conclusions

This paper studies the properties of the ICM in high-resolution hydrodynamical simulations of galaxy clusters. Using a fully cosmological scenario we focus on the build-up of baryonic and dark matter through hierarchical clustering. This decomposition allowed us to see the ICM of the final clusters as the evolution of the ICM of its progenitors, which undergo a number of dynamical and thermal changes during and after their merging in the potential well of the main system. In a way, we tried to answer questions related to the nature versus nurture of the ICM. While the analysis presented in

this paper certainly does not exhaust all questions nor gives all answers, we think it constitutes a successful attempt to statistically characterise a few properties of ICM evolution in a cosmological environment, which may be useful both to understand theoretical issues and to interpret observational X-ray data.

The main results of the present work may be summarised as follows:

* satellites merging onto a cluster undergo a fate that is well described by one free parameter, namely the mass ratio $`m_V/M_V`$ of the satellite to the main object. The survival time is a decreasing function of $`m_V/M_V`$ for the DM component, but an increasing function of $`m_V/M_V`$ for the satellite ICM. For the largest satellites (almost equal mass ratio) the DM and ICM have similar survival times.
* The dark matter and ICM of satellites exhibit quite different behaviours for different mass ratios. In small satellites DM and gas decouple quickly after they cross the virial radius of the main cluster: while on average the ICM is promptly stripped and unbound by ram-pressure effects, the DM freely moves and oscillates in the main halo; this causes a clear spatial separation between the two satellite components. In particular, the behaviour of the DM component is quite consistent with that found by Tormen, Diaferio & Syer (1998) and Taffoni et al. (2003), showing that the presence of a diffuse gas has at most a marginal effect on the global satellite properties. This might be expected, as the gravitational contribution of the ICM is of the order of the cosmological baryonic fraction.
* The mean internal properties of satellites (DM velocity dispersion and ICM temperature) evolve differently from those of the self-bound part: globally satellites are heated by their encounter with the main halo, whereas their cores (self-bound parts) can have different behaviours.

On average, the bulk of the ICM initially associated to satellites thermalizes in the new cluster potential in roughly $`2`$ Gyr, while the DM component takes much longer to adapt to the velocity dispersion of the main cluster.

Considering the self-bound part of satellites, the DM cores always cool down for all mass ratios, while the core ICM temperature heats up for small satellites and cools down for massive ones. This demonstrates that the galaxy velocity dispersion and the X-ray temperature of substructure cannot be easily linked to its actual mass, unless a more thorough analysis of the kinematics and dynamics of the system is performed.
* The ICM entropy shows a trend similar to that of the ICM temperature: entropy always increases globally, but decreases in the cores of massive satellites.
* We have modeled the mean orbital properties of satellites at their first pericentre and apocentre, using simple analytic fits. Applications of these results to optical and X-ray observations of cluster 1E0657-56 indicates that this system is more likely a merger between a small group and a massive cluster, observed after its first pericentric passage.

As a final remark, we remind that our results are based on a set of adiabatic hydrodynamic simulations. SPH simulations which include extra physics, like cooling, star formation and supernova feedback (e.g. Lewis et al. 2000; Yoshida et al. 2002; Tornatore et al. 2003; Borgani et al. 2003) show that these processes can modify the global properties of clusters, and this in turn might affect qualitatively the results presented here. However, the present uncertainties on exactly how to model these processes in simulations, and the possible influence of other physical mechanisms not yet considered (e.g. thermal conduction and magnetic fields), makes it impossible to estimate here the impact of these effects.

## Acknowledgments.

This work has been partially supported by Italian MIUR (Grant 2001, prot. 2001028932, “Clusters and groups of galaxies: the interplay of dark and baryonic matter”), CNR and ASI. We are grateful to Massimo Ramella for clarifying discussions.

## Appendix A Analytic fits and scaling laws

Here we briefly summarize the analytic fits to the results of Fig. 11, as a function of the mass ratio $`xm_V/M_V`$, where $`m_V`$ and $`M_V`$ are the pre-merger virial masses of the satellite and of the main cluster progenitor, respectively. In order to make these quantities physical they must be converted using the scaling laws also given in Sec. A.2 below.

### A.1 Analytic fits

* pericentric time:

$$\tau _p(x)=0.9[Gyr]$$
(5)
* pericentric distance:

$$d_p(x)=0.2(x+0.02)^{0.15}[R_V]$$
(6)
* pericentric velocity:

$$v_p(x)=1.8[V_c]$$
(7)
* apocentric time:

$$\tau _a(x)=1.6(x+0.02)^{0.17}[Gyr]$$
(8)
* apocentric distance:

$$d_p(x)=0.33(x+0.02)^{0.3}[R_V]$$
(9)
* apocentric velocity:

$$v_p(x)=0.33x^{0.11}[V_c]$$
(10)

### A.2 Scaling laws

The results presented above are given in rescaled units \[distances in units of $`R_V(z)`$ and velocities in units of $`V_c(z)`$\] in order to superimpose results coming from clusters of different mass and at different redshifts. In order to apply these models to real observations, we need to provide suitable conversions to physical or observational quantities. We also give the mean conversion relation between mass and internal velocity dispersion. The virial masses, radii, circular velocities and one-dimensional velocity dispersion of our simulated clusters - at redshift $`z<\text{ }0.8`$ \- are well related by the following relations, which hold for $`10^{14}h^1\mathrm{M}_{}M_V1.3\times 10^{15}h^1\mathrm{M}_{}`$:

$`\left({\displaystyle \frac{V_c(M_V)}{1000\mathrm{k}\mathrm{m}/\mathrm{s}}}\right)`$ $`=`$ $`1.410\left({\displaystyle \frac{M_V}{10^{15}\mathrm{h}^1\mathrm{M}_{}}}\right)^{0.353},`$
$`\left({\displaystyle \frac{\sigma _v(M_V)}{1000\mathrm{k}\mathrm{m}/\mathrm{s}}}\right)`$ $`=`$ $`0.982\left({\displaystyle \frac{M_V}{10^{15}\mathrm{h}^1\mathrm{M}_{}}}\right)^{0.367},`$
$`\left({\displaystyle \frac{R_v(M_V)}{\mathrm{h}^1\mathrm{Mpc}}}\right)`$ $`=`$ $`2.164\left({\displaystyle \frac{M_V}{10^{15}\mathrm{h}^1\mathrm{M}_{}}}\right)^{0.294},`$ (11)
where virial radii are in comoving $`h^1\mathrm{Mpc}`$. The corresponding observables in physical coordinates, for systems at redshift $`z`$, are:

$`V_c(\mathrm{phys})`$ $`=`$ $`V_c(M_V)\sqrt{1+z},`$
$`\sigma _v(\mathrm{phys})`$ $`=`$ $`\sigma _v(M_V)\sqrt{1+z},`$
$`R_v(\mathrm{phys})`$ $`=`$ $`R_v(M_V)/(1+z).`$ (12)


format: text

text:
# Hall effect and resistivity in underdoped cuprates

## Abstract

The behaviour of the Hall ratio $`R_H(T)`$ as a function of temperature is one of the most intriguing normal state properties of cuprate superconductors. One feature of all the data is a maximum of $`R_H(T)`$ in the normal state that broadens and shifts to temperatures well above $`T_c`$ with decreasing doping. We show that a model of preformed pairs-bipolarons provides a selfconsistent quantitative description of $`R_H(T)`$ together with in-plane resistivity and uniform magnetic susceptibility for a wide range of doping.

\]

The theory of high temperature superconductivity remains the biggest challenge in condensed matter physics today. One way of thinking is that this phenomenon is of purely electronic origin and phonons are irrelevant PWA ; KIVE ; ZHAN ; PINE . Other authors (see, for example SANF ; emi ; dev ; tru ; kab ; dzh ) explore an alternative view, namely that the extension of the BCS theory towards the strong interaction between electrons and ion vibrations describes the phenomenon. High temperature superconductivity could exist in the crossover region of the electron-phonon interaction strength from the BCS to bipolaronic superconductivity as was argued before the discovery ale0 . In the strong coupling regime, $`\lambda 1,`$ pairing takes place in real space (i.e. individual pairing)  due to a polaron collapse of the Fermi energy aleran . At first sight, polaronic carriers have a mass too large to be mobile; however it has been shown that the inclusion of the on-site Coulomb repulsion leads to the favoured binding of intersite oxygen holes ALEXAND ; CATL . The intersite bipolarons can then tunnel with an effective mass of about 10 electron masses ALEXAND ; CATL ; alekor ; tru .

The possibility of real-space pairing, as opposed to the Cooper pairing, has been the subject of much discussion. Experimental ZHAO ; mic ; ita ; TIM ; ega ; LANZ  and theoretical phil ; MOTT ; all ; gor  evidence for an exceptionally strong electron-phonon interaction in all novel superconductors is now so overwhelming that even some advocates of non-phononic mechanisms kiv  accept this fact. Nevertheless, the same authors KIVE ; kiv  dismiss any real-space pairing, suggesting a collective pairing (i.e the Cooper pairs in the momentum space) at some temperature $`T^{}>T_c`$ but without phase coherence. The existence of noncoherent Cooper pairs might be a plausible idea for the crossover region of the coupling strength, as was proposed still earlier by Dzhumanov dzh . However, apart from this, Refs.KIVE ; kiv  argue that the phase coherence and superconducting critical temperature $`T_c`$ are determined by the superfluid density, which is proportional to doping $`x`$, rather than to the density of normal state holes, which is ($`1+x)`$ in their scenario. On the experimental side, the scenario is not compatible with a great number of thermodynamic, magnetic, and kinetic measurements, which show that only holes doped into a parent insulator are carriers in both the normal and superconducting state of underdoped cuprates. On theoretical grounds, this preformed Cooper-pair (or phase-fluctuation) scenario contradicts a theorem leg , which proves that the number of supercarriers (at $`T=0`$) and normal-state carriers should be the same in any clean translation-invariant superfluid. A periodic crystal-field potential does not affect this conclusion because the coherence length is larger than the lattice constant in cuprates. Objections against real-space pairing also contradict a parameter-free estimate of the renormalised Fermi energy $`ϵ_F`$aleF , that yields $`ϵ_F`$ less than the normal state charge pseudogap $`\mathrm{\Delta }/2`$ in underdoped cuprates. The condition for real-space pairing, $`ϵ_F\pi \mathrm{\Delta },`$ is well satisfied if one admits that the bipolaron binding energy, $`\mathrm{\Delta }`$, is twice the pseudogap.

Bipolarons in cuprates could be formed by the Fröhlich interaction of holes with optical phonons ALEXAND , and by molecular phonons (Jahn-Teller interactions kab ) because of molecular-like crystal structure of these materials. Mott and Alexandrov proposed a simple model MOTT  of the cuprates based on bipolarons. In this model, all the holes (polarons) are bound into small intersite bipolarons at any temperature. Above $`T_c`$ this Bose gas is non-degenerate and below $`T_c`$ phase coherence of the preformed bosons sets in, followed by superfluidity of the charged carriers. Of course, there are also thermally excited single polarons in the model. There is much evidence for the crossover regime at $`T^{}\mathrm{\Delta }/2`$ and normal state charge and spin pseudogaps in the cuprates mih . Many experimental observations have been satisfactorily explained using this particular approach including the in-plane BRAT ; in  and out-of-plane resistivity AKM ; out , magnetic susceptibility AKM ; MULLE , tunneling spectroscopy NDROV , isotope effect ALEXA ; ZHAO , upper critical field and specific heat anomaly ZAV . ARPES measurements indicate the presence of a pseudogap as well. They also indicate an angular dependent narrow peak and a featureless background. In the polaronic model, the ARPES spectrum can be numerically explained if one considers a charge transfer Mott-insulator and the single polaron spectral function DENT .

Like many other properties, the Hall ratio in high-$`T_c`$ cuprates shows a non-Fermi-liquid behaviour tak ; bat ; car . A Fermi-liquid approach may describe $`R_H(T)`$ for $`TT^{}`$ only if vertex corrections are included. However, the advocates of this approach kon  admit that it is inappropriate for $`T<T^{}.`$ On the other hand, the bipolaron model described an enhanced magnitude, doping and $`1/T`$ temperature dependence of the Hall ratio BRAT ; ALEXAND , but the maximum of $`R_H(T)`$ in underdoped cuprates well above $`T_c`$ tak ; bat ; car  was not addressed. Also, a nonlinear temperature dependence of the $`in`$-plane resistivity below $`T^{}`$ remains one of the unsolved problems. In this paper, we give an explanation of these long-standing problems from the standpoint of the bipolaron model.

Thermally excited phonons and (bi)polarons are well decoupled in the strong-coupling regime of electron-phonon interaction SANF , so that the conventional Boltzmann kinetics for renormalised carries is applied. Here we use a ‘minimum’ bipolaron model, which includes a singlet bipolaron band and a spin 1/2 polaron band separated by $`T^{}`$, and the $`\tau `$approximation ANSE  in an electric field $`𝐄=\stackrel{}{}\varphi `$ and in a weak magnetic field $`𝐁`$ $`𝐄.`$ The bipolaron and single-polaron non-equilibrium distributions are found as

$$f(𝐤)=f_0(E)+\tau \frac{f_0}{E}𝐯\left\{𝐅+\mathrm{\Theta }𝐧\times 𝐅\right\},$$
(1)
where $`𝐯=E/𝐤,`$ $`𝐅=\stackrel{}{}(\mu 2e\varphi )`$ and $`f_0(E)=[y^1\mathrm{exp}(E/T)1]^1`$ for bipolarons with the energy $`E=k^2/(2m_b)`$, and $`𝐅=\stackrel{}{}(\mu /2e\varphi )`$ and $`f_0(E)=\{y^{1/2}\mathrm{exp}[(E+T^{})/T]+1\}^1`$, $`E=k^2/(2m_p)`$ for thermally excited polarons. Here $`m_b2m_p`$ and $`m_p`$ are the bipolaron and polaron masses of $`quasitwo`$-dimensional carriers, $`y=\mathrm{exp}(\mu /T),`$ $`\mu `$ is the chemical potential, $`\mathrm{}=c=k_B=1`$, and $`𝐧=𝐁/B`$ is a unit vector in the direction of the magnetic field. Eq.(1) is used to calculate the electrical resistivity and the Hall ratio as

$`\rho `$ $`=`$ $`{\displaystyle \frac{m_b}{4e^2\tau _bn_b(1+An_p/n_b)}},`$ (2)
$`R_H`$ $`=`$ $`{\displaystyle \frac{1+2A^2n_p/n_b}{2en_b(1+An_p/n_b)^2}},`$ (3)
where $`A=\tau _pm_b/(4\tau _bm_p)`$. The atomic densities of carriers are found as

$`n_b={\displaystyle \frac{m_bT}{2\pi }}|\mathrm{ln}(1y)|,`$ (4)
$`n_p={\displaystyle \frac{m_pT}{\pi }}\mathrm{ln}\left[1+y^{1/2}\mathrm{exp}\left(T^{}/T\right)\right].`$ (5)
and the chemical potential is determined by doping $`x`$ using $`2n_b+n_p=xn_L`$, where $`n_L`$ is the number of carriers localised by disorder. Here we take the lattice constant $`a=1`$. Polarons are not degenerate. Their number remains small compared with twice the number of bipolarons, $`n_p/(2n_b)<0.2`$, in the relevant temperature range $`TT^{}`$, so that

$$y1\mathrm{exp}(T_0/T),$$
(6)
where $`T_0=\pi (xn_L)/m_bT_c`$ is about the superconducting critical temperature of the (quasi)two-dimensional Bose gas SANF . Because of this reason, experimental $`T_c`$ was taken as $`T_0`$ for our fits.

Then using Eqs.(2,3) we obtain

$$R_H(T)=R_{H0}\frac{1+2A^2y^{1/2}(T/T_c)\mathrm{exp}\left(T^{}/T\right)}{[1+A(T/T_c)y^{1/2}\mathrm{exp}\left(T^{}/T\right)]^2},$$
(7)
where $`R_{H0}=[e(xn_L)]^1`$. In the following, we assume that the number of localised carriers depends only weakly on temperature in underdoped cuprates because their average ionisation energy is sufficiently large, so that $`R_{H0}`$ is temperature independent if $`TT^{}`$. As proposed in Ref.BRAT  the scattering rate ($`T^2`$) is due to inelastic collisions of itinerent carriers with those localised by disorder. Here we also take into account the scattering off optical phonons ANSE , so that $`\tau ^1=aT^2+b\mathrm{exp}(\omega /T)`$, if the temperature is low compared with the characteristic phonon energy $`\omega `$. The relaxation times of each type of carriers scales with their charge $`e^{}`$ and mass as $`\tau _{p,b}m_{p,b}^{3/2}(e^{})^2`$, so that $`A=(m_b/m_p)^{5/2}6`$ in the Born approximation for any short-range scattering potential. As a result we obtain the in-plane resistivity as

$$\rho (T)=\rho _0\frac{(T/T_1)^2+\mathrm{exp}\left(\omega /T\right)}{[1+A(T/T_c)y^{1/2}\mathrm{exp}\left(T^{}/T\right)]},$$
(8)
where $`\rho _0=bm_b/[2e^2(xn_L)]`$ and $`T_1=(b/a)^{1/2}`$ are temperature independent. Finally, we easily obtain the uniform magnetic susceptibility due to nondegenerate spin 1/2 polarons as AKM

$$\chi (T)=By^{1/2}\mathrm{exp}\left(T^{}/T\right)+\chi _0,$$
(9)
where $`B=(\mu _B^2m_p/\pi )`$, and $`\chi _0`$ is the magnetic susceptibility of the parent Mott insulator.

The present model fits the Hall ratio, $`R_H(T)`$, the in-plane resistivity, $`\rho (T)`$, and the magnetic susceptibility $`\chi (T)`$ of underdoped $`YBa_2Cu_3O_{7\delta }`$ with a selfconsistent set of parameters (see Fig. 1-2 and the Table). The ratio of polaron and bipolaron mobilities $`A=7`$ used in all fits is close to the above estimate, and $`\chi _01.5\times 10^4emu/mole`$ is very close to the susceptibility of a slightly doped insulator coop . The comprehensive analysis by Mihailovic et al. mih  yileds $`T^{}`$ in the range from $`200K`$ to $`1000K`$ depending on doping, and $`\omega `$ should be about $`500K`$ or so from the optical data by Timusk et al. TIM .

| $`\delta `$ | $`T_c`$ | $`\rho _0`$ | $`R_{H0}`$ | $`10^4B`$ | $`10^4\chi _0`$ | $`T^{}`$ | $`\omega `$ | $`T_1`$ |
| --- | --- | --- | --- | --- | --- | --- | --- | --- |
|  | K | $`m\mathrm{\Omega }cm`$ | $`\frac{10^9m^3}{C}`$ | $`\frac{emu}{mole}`$ | $`\frac{emu}{mole}`$ | K | K | K |
| 0.05 | 90.7 | 1.8 | 0.45 |  |  | 144 | 447 | 332 |
| 0.12 | 93.7 |  |  | 2.6 | 2.1 | 155 |  |  |
| 0.19 | 87 | 3.4 | 0.63 | 4.5 | 1.6 | 180 | 477 | 454 |
| 0.23 | 80.6 | 5.7 | 0.74 |  |  | 210 | 525 | 586 |
| 0.26 | 78 |  |  | 5.4 | 1.5 | 259 |  |  |
| 0.28 | 68.6 | 8.9 | 0.81 |  |  | 259 | 594 | 786 |
| 0.38 | 61.9 |  |  | 7.2 | 1.4 | 348 |  |  |
| 0.39 | 58.1 | 17.8 | 0.96 |  |  | 344 | 747 | 1088 |
| 0.51 | 55 |  |  | 9.1 | 1.3 | 494 |  |  |

As shown in Fig.1 and Fig.2, the model describes remarkably well the experimental data with the parameters in this range (Table), in particular the unusual maximum of the Hall ratio, well above $`T_c`$ (Fig.1), and the non-linear temperature dependence of the in-plane resistivity. The maximum of $`R_H(T)`$ is due to the contribution of thermally excited polarons into transport, and the temperature dependence of the in-plane resistivity below $`T^{}`$ is due to this contribution and the combination of the carrier-carrier and carrier-phonon scattering. It is also quite remarkable that the characteristic phonon frequency from the resistivity fit (Table) decreases with doping as observed in the neutron scattering experiments and the pseudogap $`T^{}`$ shows the doping behaviour as observed in other experiments mih . The temperature dependences of $`R_H(T),`$ $`\rho (T)`$ and $`\chi (T)`$ in underdoped $`La_{2x}Sr_xCuO_4`$ and in other underdoped cuprates are very similar to $`YBa_2Cu_3O_{7\delta }`$, when the temperature is rescaled by the pseudogap. It should be noted that adding a triplet bipolaron band MOTT ; mic  could improve the fit further (in particular, above room temperature) however increasing the number of fitting parameters.

To conclude: we applied the multi-polaron approach based on the extension of the BCS theory to the strong-coupling regime SANF  to describe peculiar normal state kinetics of underdoped cuprates. The low energy physics in this regime is that of small bipolarons and thermally excited polarons. Using this approach, we have explained the temperature dependence of the Hall ratio, the in-plane resistivity and the bulk magnetic susceptibility of underdoped cuprates. A direct measurement of the double elementary charge $`2e`$ on carriers in the normal state could be decisive. In 1993, Mott and Alexandrov NEV  discussed the thermal conductivity $`\kappa `$; the contribution from the carriers provided by the Wiedemann-Franz ratio depends strongly on the elementary charge as $`(e^{})^2`$ and should be significantly suppressed in the case of $`e^{}=2e`$. Recently, a new way to determine the Lorenz number has been applied by Zhang et al. zha , based on the thermal Hall conductivity. As a result, the Lorenz number has been directly measured in $`YBa_2Cu_3O_{6.95}`$. Remarkably, the measured value of $`L`$ just above $`T_c`$ is the same as predicted by the bipolaron model, $`L0.15L_e`$ ( $`L_e`$ is the conventional Lorenz number). A breakdown of the Wiedemann-Franz law has been also explained in the framework of the bipolaron model lor .

This work was supported by the Leverhulme Trust (grant F/00261/H) and by the Royal Society (grant ref: 15042). We are grateful to J.R. Cooper for helpful discussion of his experiments, and to A.J. Leggett for elucidating his theorem.