text:
# Self-Similar Solutions, Critical Behavior and Convergence to Attractor in Gravitational Collapse

## 1 Introduction

The basic equation of general relativity is the Einstein equations. Since the Einstein equations are simultaneous second-order partially differential equations (PDEs), it is not easy to obtain solutions. It is useful to assume symmetry such as spherical symmetry and staticity in obtaining solutions. However, it is often difficult to obtain inhomogeneous and dynamical solutions. In this context, it is quite useful to assume self-similarity. For example, in spherical symmetry, this assumption reduces the PDEs to a set of ordinarily differential equations (ODEs). Therefore, the assumption of self-similarity is very powerful in finding exact and/or numerical dynamical and inhomogeneous solutions. The application of self-similar solutions is very large including cosmological perturbations, star formation, gravitational collapse, primordial black holes, cosmological voids and cosmic censorship. Among them, the most widely researched model is the spherically symmetric self-similar system with perfect fluids, which was pioneered by by Cahill and Taub (1971) in general relativity and by Penston (1969) and Larson (1969) in Newtonian gravity.

In recent progress in general relativity, self-similar solutions have called attention not only because they are easy to obtain but also because they may play important roles in cosmological situations and/or gravitational collapse. One of the important roles of self-similar solutions is to describe asymptotic behaviors of more general non-self-similar solutions. This role has been shown explicitly in several models such as homogeneous cosmological models. In cosmological context, the suggestion that spherically symmetric fluctuations might naturally evolve from complex initial conditions via the Einstein equations to a self-similar form has been termed the ‘similarity hypothesis’ (Carr 1993, Carr and Coley 1999). The hypothesis has been also suggested to hold in more general situations, including gravitational collapse.

In the present article, we pay attention to a spherically symmetric system with a perfect fluid. For this system, recent numerical simulations of PDEs have revealed two interesting behaviors in gravitational collapse, one is the critical behavior and the other is the convergence to attractor. The critical behavior was observed around the threshold of black hole formation. There is a self-similar solution that sits at the threshold, which is called a critical solution. For near-critical collapse, the solution first approach the critical solution but diverges away from this solution after that. There is observed the scaling law for the formed black hole for supercritical collapse. This is observed only as a result of fine-tuning the parameter which parametrizes initial data. On the other hand, the convergence to attractor is observed without fine-tuning. The density and velocity fields around the center approach those of another self-similar solution in an approach to the formation of singularity at the center, which is called an attractor solution.

The main claim of this review article is that the critical behavior and convergence to attractor in gravitational collapse are both well understood in the context of self-similarity hypothesis and stability analysis in gravitational collapse. The difference between these two behaviors is in the stability of the pertinent self-similar solution. We also briefly review the recent work on how these two behaviors are in Newtonian gravity.

## 2 Self-similarity

In general relativity, self-similarity or homothety is covariantly defined in terms of the homothetic vector $`\xi `$ which satisfies the homothetic Killing equation

$$_\xi g_{\mu \nu }=2g_{\mu \nu }.$$
(1)
We consider a perfect fluid as matter field, the stress-energy tensor of which is given by

$$T^{\mu \nu }=(\rho +P)u^\mu u^\nu +Pg^{\mu \nu }.$$
(2)
In a spherically symmetric spacetime with self-similarity, where the homothetic vector is neither parallel nor orthogonal to the four velocity of the fluid element, the line element is written in terms of the comoving coordinate $`r`$ as (Cahill and Taub 1971)

$$ds^2=e^{\sigma (\xi )}dt^2+e^{\omega (\xi )}dr^2+r^2S^2(\xi )(d\theta ^2+\mathrm{sin}^2\theta d\varphi ^2),$$
(3)
where $`\xi =t/r`$. In this coordinate system, the dimensionless quantity $`h`$ such as $`e^\sigma `$, $`e^\omega `$ and $`S`$ has scale-invariance in the following sense:

$$h(t,r)=h(at,ar)=h(\xi ),$$
(4)
where $`a>0`$ is an arbitrary constant. For a self-similar spacetime, the density field $`\rho `$ is given via the Einstein equations by

$$\rho (t,r)=a^2\rho (at,ar).$$
(5)
It implies that the density field is written as

$$\rho (t,r)=\frac{1}{8\pi r^2}\eta (\xi ).$$
(6)
Therefore, the time evolution is equivalent with the scale transformation. An example of the time evolution of the density field is displayed in Fig. 1.

The above definition of self-similarity can be generalized. One possibility is to discretize the above definition: the dimensionless quantity $`h`$ satisfies

$$h(t,r)=h(at,ar),a=e^{n\mathrm{\Delta }},$$
(7)
where $`n=0,\pm 1,\pm 2,\mathrm{}`$. This kind of self-similarity is referred to as discrete self-similarity. The discrete self-similar solution has been found to be important in the critical behavior observed in the gravitational collapse of a massless scalar field (Choptuik 1993) and gravitational waves (Abrahams and Evans 1993). In this context, the originally defined self-similarity is referred to as continuous self-similarity.

Another possible generalization is to allow the dimensionless quantity $`h`$ to be a function of $`r/t^\alpha `$. The generalization along this direction was done by Carter and Henriksen (1989,1991) and is referred to as kinematic self-similarity. The kinematic self-similarity is covariantly defined by the existence of kinematic self-similarity vector, which satisfies

$$_\xi g_{\mu \nu }=2\delta g_{\mu \nu },_\xi u_\mu =\alpha u_\mu ,$$
(8)
where $`\delta `$ and $`\alpha `$ are arbitrary constants. See Coley (1997), Benoit and Coley (1998) and Maeda et al. (2002a, 2002b, 2003) for more details of kinematic self-similarity and kinematic self-similar solutions. In this context, the originally defined self-similarity is referred to as first-kind self-similarity.

## 3 Self-similar solutions

Hereafter we mainly concentrate on the originally defined self-similarity and refer to it simply as self-similarity unless it is confusing. For a perfect fluid, the Einstein equation requires that the equation of state be given by the following form:

$$P=k\rho ,$$
(9)
where $`k`$ is constant (Cahill and Taub 1971). The fluid is referred to as a dust fluid for $`k=0`$, a radiation fluid for $`k=1/3`$ and a stiff matter for $`k=1`$. There are infinitely many self-similar solutions. These solutions are obtained by solving ODEs. The classification of self-similar solutions of this system has been done in the homothetic approach by Goliath, Nilsen and Uggla (1998a, 1998b) and in the comoving approach by Carr and Coley (2000a, 2000b). For the complementarity of these two approaches, see Carr et al. (2001). For a brief review of self-similar solutions and, in particular, the classifications of them, see Carr and Coley (1999).

In obtaining physical solutions, it will be required that the solution be obtained from regular initial conditions. In particular, the physical solution must have the initial data with regular center. For $`0<k<1`$, the central density parameter $`D_0`$, which is defined by

$$D_04\pi \rho (t,0)t^2,$$
(10)
well parametrizes the regular center of self-similar solutions, where $`t`$ is chosen to be negative. For the regular center of self-similar solutions, the $`D_0`$ is finite and does not depend on $`t`$ because

$$D_0=\underset{r0+}{lim}4\pi \rho (t,r)t^2=\underset{\xi \mathrm{}}{lim}\frac{1}{2}\xi ^2\eta (\xi ).$$
(11)
For $`0<k<1`$, the ODEs have a singular point which is called a sonic point, at which the relative speed $`|V|`$ of the fluid element to the curve $`\xi =\text{const}`$ is equal to the sound speed $`\sqrt{k}`$. The matching condition beyond the sonic point results in nontrivial reduction of self-similar solutions. The reduction of the self-similar solutions was demonstrated by Ori and Piran (1990). If we require regularity, by which we mean the continuity of first derivatives, the possible solutions make a band structure in the space of all solutions. The $`D_0`$ axis is divided into allowed ranges and forbidden ranges (for which the solution terminates at the sonic point). If we further require analyticity (Taylor-series expandability), the possible solutions are discretized. The regular center is also analytic, where it should be noted that the analyticity at the center is defined with respect to the local Cartesian coordinates. However, only a discrete set of values for $`D_0`$ is allowed for the analyticity of the sonic point.

This complicated structure of physical solutions is due to the character of the sonic point as a critical point of the autonomous system, which was studied by Bicknell and Henriksen (1978). The sonic point regarded as a critical point of the autonomous system of ODEs after an appropriate transformation and is considered in the two-dimensional subspace of the whole state space. As a result, sonic points are classified into nondegenerate nodes, degenerate nodes, saddles and foci, among which the foci are unphysical. For every physical sonic point, there are two possible directions along which solutions can cross the sonic point, and there are two analytic solutions corresponding to these two directions. For a nondegenerate node, the two directions are called primary and secondary ones. For the former, there are one-parameter family of regular solutions including one of the analytic solutions, while for the latter, another analytic solution alone is possible. For a degenerate node, the two directions are degenerate and there are one-parameter family of regular solutions including the one analytic solution. For a saddle one, all that are regular at the sonic point are the two analytic solutions alone.

The numerical solutions with analytic initial data were investigated by Ori and Piran (1987, 1990) and Foglizzo and Henriksen (1993). This discrete set of solutions includes the flat Friedmann solution, general relativistic counterpart of the Larson-Penston solution, and that of the Hunter family of solutions. However, it is noted that the whole set of self-similar solutions with analytic initial data is still uncertain because of the numerical difficulty. Hereafter we refer to the general relativistic Larson-Penston solution and Hunter (a) solution as the Ori-Piran solution and Evans-Coleman solution for convenience. The most researched nontrivial self-similar solutions are the Ori-Piran and Evans-Coleman solutions. For the character of the sonic point, the Ori-Piran solution crosses a degenerate node along the secondary direction for $`0<k\begin{array}{c}<\hfill \\ \hfill \end{array}0.036`$, a degenerate node for $`k0.036`$ and a nondegenerate node along the primary direction for $`0.036\begin{array}{c}<\hfill \\ \hfill \end{array}k<1/3`$ (Ori and Piran 1990, Harada 2001). Carr et al. (2001) investigated the detailed nature of the Evans-Coleman solution for $`0<k<1`$. The Evans-Coleman solution exists for $`0<k1`$ (Brady et al. for $`k=1`$). It crosses a saddle for $`0<k\begin{array}{c}<\hfill \\ \hfill \end{array}0.41`$, a nondegenerate node along the secondary direction for $`0.41\begin{array}{c}<\hfill \\ \hfill \end{array}k\begin{array}{c}<\hfill \\ \hfill \end{array}0.89`$, a degenerate node for $`k0.89`$ and a nondegenerate node along the primary direction for $`0.89\begin{array}{c}<\hfill \\ \hfill \end{array}k<1`$. The Ori-Piran solution describes coherent collapse, while the Hunter family of solutions some oscillative collapse. The Hunter family of solutions are well labelled by the number of zeroes in their velocity fields. Ori and Piran (1987, 1990) also pointed out that a naked singularity forms in the Ori-Piran solution for $`0<k\begin{array}{c}<\hfill \\ \hfill \end{array}0.0105`$, which may have some implication to cosmic censorship proposed by Penrose (1969, 1979). (See Harada, Iguchi and Nakao (2002) for recent works in naked singularity physics.) Actually, the general relativistic counterpart of the Hunter family of solutions with more than one oscillation can also have a naked singularity (Ori and Piran 1990, Foglizzo and Henriksen 1993). It is not clear whether nontrivial solutions other than the Evans-Coleman solution exist for every value of $`k(0,1]`$, although numerical studies suggest that the structure of the solutions is common at least for $`0<k9/16`$ (Ori and Piran 1990, Foglizzo and Henriksen 1993). The density profiles for these solutions are shown in Fig. 2.

For $`k=0`$, i.e., a dust fluid, things are somewhat different. Ori and Piran (1990) studied self-similar dust solutions with regular center. There is one-parameter family of solutions. The regularity at the center requires that the central density parameter $`D_0`$ be $`2/3`$. The solutions are parametrized not by the central density parameter but by the central third-order derivative of the density field. However, this regular center is not analytic except for the flat Friedmann solution because the third-order derivative of the density field is not zero. Therefore, if we require that a physical solution has analytic initial data, the only allowed solution is the flat Friedmann solution. There is no sonic point for the dust case. See Carr (2000) for more general class of self-similar dust solutions.

## 4 Critical behavior

The critical behavior in gravitational collapse was numerically discovered by Choptuik (1993) in the spherical system of a massless scalar field. For a one-parameter family of initial data sets for this system, a critical value is found beyond which the evolution results in black hole formation (see Fig. 3(a)). For near-critical values, the critical behavior is observed. The solution first tends to a discrete self-similar solution, which is called a critical solution, and diverges away after that. For near-critical and supercritical collapse, the mass of the formed black hole follows a power law as

$$M_{BH}|pp^{}|^\gamma ,$$
(12)
where the constant $`\gamma `$ is called a critical exponent (see Fig. 3(b)). The critical exponent for this system was estimated to be $`0.37`$. After this discovery, a similar behavior was observed in various systems, such as gravitational waves (Abrahams and Evans 1993), a radiation fluid (Evans and Coleman 1994), a perfect fluid, a stiff matter (Neilsen and Choptuik 2000) and so on. Then, it has been found that the critical solution can be a discrete self-similar solution or continuous self-similar solution, depending on the system. For the fluid system including the radiation fluid system, the critical solution has continuous self-similarity. For the system of a radiation fluid, Evans and Coleman (1994) estimated the critical exponent to be $`0.36`$ and also found that the critical solution obtained as a result of numerical simulations of the PDEs agrees very well with one of the continuous self-similar solutions which are obtained by integrating the ODEs. (That is why we refer to this self-similar solution as the Evans-Coleman solution in this article. We shall avoid the term ‘critical solution’ to refer to this solution because whether this solution is a critical solution may depend on the system, for example, whether we allow nonspherical perturbations or not, as we will mention later.) It has been subsequently found that the value of the critical exponent is not universal but depends on the parameter $`k`$ which parametrizes the equation of state $`P=k\rho `$.

Koike, Hara and Adachi (1995) gave a clear picture for the critical behavior, which is called a renormalization group approach. According to this picture, the critical solution is assumed to be a self-similar solution with analytic initial data with a single unstable mode. We consider the time evolution of the solution with respect to the self-similarity coordinates $`\tau \mathrm{ln}|t|`$ and $`x\mathrm{ln}|\xi |`$. Then, the time evolution is regarded as a flow in the space of functions of $`x`$. The self-similar solution with a single unstable mode, which is a fixed point in the space of functions of $`x`$, has a stable manifold of codimension one. This implies that a one-parameter family of initial data generically has intersection with the stable manifold. The collapse which develops from the initial data set which is the intersection corresponds to the critical collapse, while near-critical collapse develops from the initial data set near the intersection. After a long time, the deviation of the near critical collapse from the critical solution is dominated by the single unstable mode, i.e.,

$$h(\tau ,x)H_{ss}(x)+(pp^{})e^{\kappa \tau }F_{rel}(x),$$
(13)
where $`H_{ss}`$ is the critical solution and $`\kappa `$ and $`F_{rel}`$ are the eigenvalue and eigenfunction of the single unstable mode. Since the apparent horizon appears when the deviation becomes $`O(1)`$, the mass of the black hole is estimated as

$$M=O(r)=O(e^\tau )|pp^{}|^\gamma $$
(14)
where $`\gamma =\kappa ^1`$. The eigenvalue and eigenfunction of the unstable mode are determined by the linearized Einstein equation with the regularity condition both at the center and at the sonic point. Koike, Hara and Adachi calculated the critical exponent to be $`0.3558019`$ by the eigenvalue analysis, which agrees very well with the estimated value $`0.36`$ from the numerical simulation.

Critical phenomena in gravitational collapse have been observed and widely researched in a great variety of systems. The implication to the formation of primordial black holes was also discussed.

One of the interesting systems is the spherical system of a stiff matter. Neilsen and Choptuik (2000) observed the critical behavior in the stiff-matter system and found that the critical solution is a continuous self-similar solution. They also estimated the critical exponent to be $`0.96\pm 0.02`$. Brady et al. (2002) confirmed that the Evans-Coleman solution for $`k=1`$ has a single unstable mode and that the critical exponent derived from the eigenvalue analysis agrees well with that estimated from the result of numerical simulations of the PDEs. This critical behavior for the stiff-matter system seems to be totally different from that for the scalar-field system. On the other hand, a massless scalar field is equivalent with a stiff matter under certain circumstances (Madsen 1988, Brady et al. 2002). Brady et al. (2002) pointed out that the critical solution, which is discrete self-similar, for the system of a massless scalar field cannot be constructed from the stiff matter because the equivalence is broken.

See Gundlach (1998, 1999) for a recent review of critical phenomena in gravitational collapse.

## 5 Stability

The stability of self-similar solutions against spherically symmetric perturbations has called much attention in the context of the critical behavior in gravitational collapse since Koike, Hara and Adachi (1995) showed that the critical behavior is understood in terms of the behavior of the unstable mode of the critical solution.

There are two types of spherically symmetric perturbations, analytic modes and kink modes. First we consider the analytic modes. The perturbation for a dimensionless quantity $`h`$ is written as

$$h(\tau ,x)=h_0(x)+\delta h(\tau ,x).$$
(15)
After the linearization with respect to the perturbation, the perturbation $`\delta h`$ is separated in two factors as

$$\delta h(\tau ,x)=ϵe^{\kappa \tau }F(x).$$
(16)
For $`0<k1`$, when the regularity condition is imposed both at the center and at the sonic point, the $`\kappa `$ must take a discrete set of values. In other words, the allowed values for $`\kappa `$ are determined as eigenvalues. If $`\kappa >0`$, the mode is growing and referred to as an unstable mode. The stability of the Evans-Coleman solution has been investigated by Koike, Hara and Adachi (1995) for $`k=1/3`$ and by Maison (1996), Koike, Hara and Adachi (1999), Harada and Maeda (2001) for more general value of $`k(0,1)`$, and Brady et al. (2002) for $`k=1`$. All these works have shown that the Evans-Coleman solution has a single analytic unstable mode for $`0<k1`$. For the Evans-Coleman solution, which has a single unstable mode, the inverse of the eigenvalue of the unstable mode is interpreted as the value of the critical exponent, as we have seen in Sec. 4. The results are summarized in Fig. 4. The critical exponents directly estimated from the mass-scaling law obtained after the numerical simulations of the PDEs are also plotted.

The stability of other self-similar solutions with analytic initial data were investigated by Koike, Hara and Adachi (1999) and Harada and Maeda (2001). Then, it has been found that the Ori-Piran solution have no unstable mode, while the general relativistic counterpart of the Hunter-family solutions have unstable modes, the number of which agrees with the number of zeroes in the velocity field at least for the solutions examined there. Since the nonuniqueness of unstable modes implies that the codimension of the stable manifold is greater than one, it is expected that the approach to the self-similar solution is realized only as a result of fine-tuning more than one parameters simultaneously. This is the reason why the numerical simulations of the PDEs do not hit the self-similar solution with more than one unstable modes. On the other hand, since the Ori-Piran solution has no unstable mode, it is expected that generic evolution go towards this solution without any fine-tuning. In the next section, we will see that actually this is the case.

Next, for the kink modes, the existence in general relativity was expected by Ori and Piran (1990) and demonstrated by Harada (2001). For the analysis, we consider the initial perturbations which satisfy the following conditions: (For a while, we assume $`t<0`$ and later will mention the case for $`t>0`$.) (1) The initial perturbations vanish inside the sonic point. (2) The density field is continuous. (3) The density gradient field is discontinuous at the sonic point, although it has definite one-sided values inside and outside. See Fig. 5(a) for the initial setting of perturbations. The evolution of the density gradient at the sonic point turns out to be local, which simplifies the analysis. The analysis of the kink mode is possible in full order. Then, it is found that the stability against the kink mode is closely related to the classification of the sonic point, where the instability is defined by the blow up of the density gradient in a finite interval in terms of the time coordinate $`\tau `$.

The result is as follows. All primary-direction nodal-point solutions, degenerate-nodal-point solutions, and saddle-point anti-transsonic solutions are unstable against the kink mode, while all secondary-direction nodal-point solutions, and saddle-point transsonic solutions are stable against the kink mode. The situation is reversed for $`t>0`$ except for the case of degenerate-nodal-point solutions, which are unstable against the kink mode both for $`t>0`$ and $`t<0`$. Figure 5(b) shows the nonlinear evolution of the kink mode for a sonic point which is a nondegenerate node, as a typical example.

When we apply this criterion to the known solutions, we obtain the following results for the stability against the kink mode. The flat Friedmann (collapsing) solution is unstable for $`0<k<1/3`$ while stable for $`1/3k<1`$. The Ori-Piran solution is stable for $`0<k\begin{array}{c}<\hfill \\ \hfill \end{array}0.036`$ while unstable for $`0.036\begin{array}{c}<\hfill \\ \hfill \end{array}k<1/3`$. The Evans-Coleman solution is stable for $`0<k\begin{array}{c}<\hfill \\ \hfill \end{array}0.89`$ while unstable for $`0.89\begin{array}{c}<\hfill \\ \hfill \end{array}k<1`$. Moreover, all regular but nonanalytic (collapsing) solutions are unstable.

The kink instability of the Evans-Coleman solution for $`0.89\begin{array}{c}<\hfill \\ \hfill \end{array}k<1`$ will affect the critical behavior of a perfect fluid. However, Neilsen and Choptuik (2000) observed the critical behavior even for $`0.89\begin{array}{c}<\hfill \\ \hfill \end{array}k1`$ which is considered to be continuous extrapolation of that for $`0<k\begin{array}{c}<\hfill \\ \hfill \end{array}0.89`$. It has not been explicitly shown but is suggested that the Evans-Coleman solution has a single analytic unstable mode even for $`0.89\begin{array}{c}<\hfill \\ \hfill \end{array}k<1`$ from their result. Moreover, it was also reported that the Evans-Coleman solution has a single analytic unstable mode for $`k=1`$ (Brady et al. 2002). About the reason why the kink instability did not affect the PDE calculations by Neilsen and Choptuik (2000), Harada (2001) speculated several possibilities.

The stability of the Evans-Coleman solution against nonspherical perturbations was investigated by Gundlach (2002). The result is that the Evans-Coleman solution is stable against nonspherical perturbations for $`1/9<k\begin{array}{c}<\hfill \\ \hfill \end{array}0.49`$ while is unstable for $`0<k<1/9`$ or $`0.49\begin{array}{c}<\hfill \\ \hfill \end{array}k<1`$. For the stability of other self-similar solutions with analytic initial data against nonspherical perturbations, little has been known until now.

## 6 Attractor

As seen in the previous section, we have at least one self-similar solution which has no unstable mode. The solution is the Ori-Piran solution for $`0<k\begin{array}{c}<\hfill \\ \hfill \end{array}0.036`$. From the discussion based on the renormalization group, which we have seen in Sec. 4, this solution seems to be an attractor which generic solution approach in the vicinity of the singularity formation. The numerical result obtained by Harada and Maeda (2001) strongly suggests that this is the case. In their numerical simulation of the PDEs, the collapse ended in singularity formation at the center for a certain subset of initial data sets, where the initial data sets were prepared without fine-tuning parameters. Thus obtained numerical solutions were compared with the self-similar solutions obtained by integrating the ODEs. Then, for the numerical solutions obtained by numerical simulations of the PDEs, the good continuous self-similarity is found and the profile of the density field agrees very well with the Ori-Piran solution, as seen in Fig. 6. It was confirmed that this convergence to the self-similar attractor does not depend on the choice of the initial density profile. The above is the results of numerical simulations for $`0<k0.03`$.

However, it is noted that we can consider initial data sets, from which the solution does not approach the Ori-Piran solution. For example, it is clear that the initial data set for the self-similar solutions other than the Ori-Piran solution does not show the approach to the Ori-Piran solution. Furthermore, if the central region is initially homogeneous and sufficiently compact, it can be shown that the solution approach the flat Friedmann solution around the center as the collapse goes towards the singularity formation. It is suggestive that the final fate of all these counterexamples is also described by self-similar solutions. Anyway, the result of numerical simulations strongly suggests that the set of the counterexamples to the convergence to the attractor is zero-measure in the space of initial data sets. These numerical results are well incorporated in the picture of renormalization flow; the critical solution as a saddle and the attractor solution as a sink in the renormalization flow of the function space.

This is one realization of the self-similarity hypothesis in the spherically symmetric gravitational collapse. Moreover, since a naked singularity develops from the analytic initial data in the Ori-Piran solution for $`0<k\begin{array}{c}<\hfill \\ \hfill \end{array}0.0105`$, this numerical study strongly suggests that the violation of cosmic censorship is generic in spherically symmetric spacetime. It should be noted that the numerical study based on the null formulation also obtained qualitatively the same result (Harada 1998).

From these considerations, it is expected that one of self-similar solutions may describe the asymptotic behavior of generic and/or special solutions in more general situations. Actually, several other examples are known to exhibit that a self-similar solution describes the asymptotic behaviors of more general solutions (Carr and Coley 1999). If the self-similarity hypothesis is true in a large variety of systems, not only will the analysis of the singularity and/or the asymptotic behavior of the spacetime be much simplified, but also should this be argued in the context of the structure of generic singularity (cf, Belinsky, Khalatnikov and Lifshitz 1970).

## 7 Newtonian case

Newtonian self-similar solutions have been investigated mainly in an effort to obtain realistic solutions of gravitational collapse leading to the star formation. Penston (1969) and Larson (1969) independently found a self-similar solution, which describes a gravitationally collapsing isothermal gas sphere. Here we call this solution the Larson-Penston solution. Thereafter, Hunter (1977) found a new series of self-similar solutions with analytic initial data and that the set of such solutions is infinite and discrete. Since he named the solutions (a), (b), (c), (d) and so on, here we call these solutions the Hunter (a), (b), (c) and (d) solutions and so on. Shu (1977) found other solution and discussed its astrophysical application. Whitworth and Summers (1985) pointed out the existence of a new family of self-similar solutions with loss of analyticity and found that these solutions make a band structure in the space of all solutions. When Ori and Piran (1990) extended these Newtonian self-similar solutions to general relativity, they showed that these Newtonian self-similar solutions are obtained as the limit of $`k0`$ of general relativistic self-similar solutions. Even for a spherical gas system with the polytropic equation of state, it is possible to extend the definition of self-similar solutions and these solutions were studied by Goldreich and Weber (1980), Yahil (1983) and Suto and Silk (1988). However, here we concentrate on the case of isothermal gas.

There are two kinds of perturbation modes, the one is the analytic modes and the other is the kink modes, as we have seen in general relativistic case. For the former, the stability of self-similar solutions were studied by Hanawa and Nakayama (1997), Hanawa and Matsumoto (2000a, 2000b) and Maeda and Harada (2001) against spherically symmetric and nonspherical perturbations. For spherically symmetric perturbations, it is found that the homogeneous collapse solution and the Larson-Penston solution have no unstable mode, while the Hunter (a), (b), (c), (d) solutions have one, two, three and four unstable modes, respectively. It is expected all Hunter-family solutions have more than one unstable modes except for the Hunter (a) solution.

For the kink-mode perturbation, Ori and Piran (1988) showed that for $`t<0`$ all primary-direction nodal-point solutions, degenerate-nodal-point solutions, and saddle-point anti-transsonic solutions are unstable against the kink mode, while all secondary-direction nodal-point solutions, and saddle-point transsonic solutions are stable against the kink mode. The situation is reversed for $`t>0`$ except for the case of degenerate-nodal-point solutions, which are unstable against the kink mode both for $`t>0`$ and $`t<0`$. Actually, the above result is completely the same as in general relativity. Applying this criterion against the kink mode, it is found that the homogeneous (collapsing) solution is unstable while the Larson-Penston solution and the Hunter (a) solution are stable.

Since the above analysis indicates that the Hunter (a) solution has a single unstable mode, Maeda and Harada (2001) identified that the Hunter (a) solution is a critical solution in Newtonian collapse. The critical exponent calculated from the eigenvalue analysis is $`0.10567`$, which is appropriate value as the $`k0`$ limit of the critical exponent in general relativity. Moreover, the character of the Hunter (a) solution is very similar to that of the Evans-Coleman solution. For example, for both solutions, the central region contracts, the surrounding region expanding. From these observations, Harada and Maeda (2001) and Maeda and Harada (2001) identified that the Evans-Coleman solution is the general relativistic counterpart of the Hunter (a) solution.

Moreover, since the Larson-Penston solution has no unstable mode, it is regarded as an attractor solution. Actually, it has been considered and demonstrated that the Larson-Penston solution describes the generic collapse of an isothermal gas near the center since the discovery of this solution (Penston 1969, Larson 1969, Foster and Chevalier 1993, Tsuribe and Inutsuka 1999).

To find the critical behavior in Newtonian collapse and see the relation between the attractor and critical solutions, Harada, Maeda and Semelin (2003) performed the numerical simulations of the PDEs which govern the spherical system of an isothermal gas. The result is as follows: When a one-parameter of initial data sets is prepared, there is a critical value of the parameter. For the critical value the collapse tends to the Hunter (a) solution. For near-critical subcritical case, the collapse first approaches the Hunter (a) solution, thereafter diverges from that and finally disperses away. For near-critical supercritical case, the collapse first approaches the Hunter (a) solution, thereafter diverges from that and finally approaches the Larson-Penston solution. For the solution to approach the Larson-Penston solution, the parameter fine-tuning is not necessary but the parameter only has to be supercritical. The scaling law appears for the mass of the collapsed core for the near-critical supercritical case, while for the maximum density for the near-critical subcritical case. The critical exponents for these laws were estimated to be $`0.11`$ and $`0.22`$, which agree very well with the expected values, $`0.10567`$ and $`0.21134`$, respectively. Figure 7 shows the trajectory of solutions and the obtained scaling law for the collapsed core mass.

## 8 Summary

The general theory of relativity as well as Newtonian gravity admits self-similar solutions. This is due to the scale-invariance of gravity. The self-similar solutions are important not only because they are dynamical and inhomogeneous solutions easy to obtain but also because they may play important roles in the asymptotic behavior of more general solutions. We have seen that the stability analysis of self-similar solutions gives a unified picture of the critical behavior and convergence to attractor in gravitational collapse, which have been both observed in recent numerical simulations of the PDEs. Since both the critical behavior and convergence to attractor are not only in general relativity but also in Newtonian gravity, these two are considered to be common characters of gravitational collapse physics. In the above sense, gravitational collapse singles out self-similar solutions even if the initial condition given is complicated. Finally, this article is concluded by the schematic picture for gravitational collapse in Fig. 8.


format: text

text:
# Resolving signals in the LISA data

## Acknowledgments

One of us (AK) acknowledges support from the National Research Council under the Resident Research Associateship program at the Jet Propulsion Laboratory, California Institute of Technology, and from the Polish State Committee for Scientific Research (KBN) through Grant No. 2P03B 094 17. The research was performed at the Jet Propulsion Laboratory, California Institute of Technology, under contract with the National Aeronautics and Space Administration.


text:
# Are Compact High-Velocity Clouds Extragalactic Objects?

## 1 Introduction

The anomalous velocity clouds of neutral hydrogen known as the high-velocity clouds (HVCs) may represent the continuing infall of matter onto the Local Group (e.g., Oort 1966; Verschuur 1969; Blitz et al. 1999). The distances to the majority of these clouds remains unknown, but in the Local Group scenario the distances would range from a few kpc (for gas currently accreting onto the Galaxy) to beyond Andromeda ($`1`$ Mpc). HVCs range in size from $`0.1100`$ deg2, and it is possible that the sizes reflect their distance from the Galaxy. In particular, some of the HVCs are both compact and isolated from extended emission, and show many similar HI properties to dwarf galaxies. These are referred to as compact HVCs (CHVCs) and may be the pristine building blocks of the Local Group at $`1`$ Mpc.

The survival and composition of these CHVCs against photoionization depends on their volume densities (and therefore distances) and the strength of the extragalactic ionizing photon flux. The presence of an extragalactic ionizing radiation field is inferred from several sets of observations. Maloney (1993) used a deep H I observation of NGC 3198 (van Gorkom 1993), which shows a sharply truncated neutral hydrogen disk at $`N_{\mathrm{HI}}5\times 10^{19}`$ cm-2, to infer that the flux in the energy range $`hv=13.6`$ to $`200`$ eV is approximately $`\varphi _i=10^4`$ phot cm-2 s-1. Estimates from the proximity effect at low redshift (Kulkarni & Fall 1993; Scott et al. 2002), upper limits from sensitive H$`\alpha `$ observations (Madsen et al. 2001; Weymann et al. 2001), and theoretical estimates based on known sources (Haardt & Madau 1996; Shull et al. 1999) are in agreement with a value $`\varphi _i0.51.5\times 10^4`$ phot cm-2 s-1.

In this paper we examine the ionization state of the compact high-velocity clouds at various distances, subject to an extragalactic ionizing radiation field. In the original Local Group model, the CHVCs are confined by dark matter halos, with $`M_{dm}/M_{gas}10`$. In a more recent variant of this model, it is suggested that the CHVCs are condensations in an intragroup medium, bound by external pressure (Blitz, private communication). We therefore consider both constant gas density models and dark matter-dominated models. The H I properties of the CHVCs are briefly reviewed in section 2 and the photoionization model is described in section 3. The results of the models are presented in section 4, and in section 5 we discuss the implications of the results.

## 2 HI Data

The H I properties of the CHVCs are taken from the H I Parkes All-Sky Survey (HIPASS) HVC catalog of Putman et al. (2002). This is a catalog of high-velocity clouds with declination $`\delta <+2^{}`$ and $`|v_{\mathrm{LSR}}|500`$ km s-1, excluding emission with $`|v_{\mathrm{LSR}}|<90`$ km s-1 to avoid confusion with Galactic emission. The spatial resolution of the survey is $`15.5^{}`$ and the velocity resolution is 26 km s-1. The $`5\sigma `$ column density sensitivity of the HVC HIPASS data is $`2\times 10^{18}`$ cm-2, which we adopt as a cutoff. Barnes et al. (2001) give a full description of HIPASS. The original CHVC catalog of Braun & Burton (1999) was based on data of lower resolution and sensitivity, and 50% of the CHVCs in the region of overlap between the two catalogs were reclassified based on the HIPASS data. In order for a high-velocity cloud to be considered a compact and isolated CHVC in the HIPASS catalog, it must have a diameter at 25% of the peak column density less than 2°, and must not be elongated in the direction of any extended emission. The typical parameters of a HIPASS CHVC are an angular size of 0.36 deg2 (to the sensitivity limit), a total H I flux of 19.9 Jy km s-1, and a typical peak column density of $`N_{\mathrm{HI}}10^{19}`$ cm-2. There are a total of 179 CHVCs cataloged from HIPASS and when this is combined with the northern sky data from the Leiden-Dwingeloo survey there are approximately 250 CHVCs in total (de Heij, Braun & Burton 2002). The linewidth distribution shows a prominent peak at a FWHM $`\mathrm{\Delta }V=25\mathrm{km}\mathrm{s}^1`$ (de Heij, Braun & Burton).

For most of the results in this paper, we have used the typical CHVC parameters quoted above. This is a reasonable assumption, since the distributions in column density, angular size, and velocity are fairly sharply peaked. An increase or decrease of less than a factor of three in column density away from the peak value produces a drop by a factor of five or more in the number of clouds111Sensitivity starts to affect the column density distribution for column densities not too far below the peak, and so the drop-off with decreasing column below the peak is not reliably established. However, since clouds with column neutral hydrogen column densities smaller than the peak value will be even more highly ionized than clouds with our typical parameters, the actual number of clouds at the low column density end is irrelevant for the purposes of this paper., and similarly for the cloud areas (cf. Figures 12 and 13 of Putman et al. ). The linewidth distribution of de Heij et al. (2002) is less sharply peaked than the column density or angular size distributions, and is skewed to the low-velocity end: $`70\%`$ of the observed clouds have a FWHM $`\mathrm{\Delta }V=25\mathrm{km}\mathrm{s}^1`$ or less; less than $`10\%`$ have $`\mathrm{\Delta }V`$ in excess of $`42\mathrm{km}\mathrm{s}^1`$. Despite the well-defined typical values, in order to verify that our conclusions are robust to uncertainties in the cloud parameters, we have also included models in which the peak central column density is up to an order of magnitude larger than the typical CHVC value, and the linewidth is doubled to $`\mathrm{\Delta }V=50\mathrm{km}\mathrm{s}^1`$ (see §4.2.2).

## 3 Model

The photoionization models of CHVCs have been calculated with the code described by Maloney (1993), modified to allow for the spherical geometry of the CHVCs. (Modest – tens of percent – departures from sphericity will have no significant impact on the model results.) Spherical geometry raises the ionization rates substantially over plane-parallel models with the same total hydrogen density; integration of the ionization rates over angle is done using 10-point Gaussian quadrature. The photoionization equilibrium of a CHVC exposed to an isotropic background radiation field is calculated iteratively, including the effects of the diffuse radiation. Thermal equilibrium at T $`10^4`$ K is assumed, as is appropriate for gas in the warm neutral phase. An ionizing photon flux $`\varphi _i10^4`$ phot cm-2 s-1 has been assumed; a flux a factor of two smaller would reduce the total gas masses by a factor of $`1.5`$ (see below). A radial grid with 100 depth points (with logarithmic spacing) was used; we have verified that this provides adequate accuracy (less than $`1\%`$ error in the total column density for a fixed neutral column density). In all cases we have assumed that the only ionizing photon flux is due to the extragalactic background; depending on the fraction of ionizing photons that escape from the Milky Way, the Galactic ionizing radiation field could dominate for CHVC distances $`d\stackrel{<}{}100`$ kpc (Maloney & Bland-Hawthorn 1999).

## 4 Results

As noted in the introduction, we have considered both models with uniform gas (total hydrogen) density $`n_\mathrm{H}`$ and models where the gas is confined by dark matter potential wells. We first discuss the constant $`n_\mathrm{H}`$ models, as these contain much of the relevant physics, and then the models that include dark matter.

### 4.1 Uniform density models

If the CHVCs were at local group distances ($`1`$ Mpc), their H I volume densities would be $`n_{\mathrm{HI}}10^4\mathrm{cm}^3`$. We therefore consider total hydrogen densities $`n_\mathrm{H}`$ ranging from $`10^4`$ to $`10^2`$ cm-3; as we will see, the low-density end of this range is untenable for reasonable HVC distances. The model results (total hydrogen columns, gas masses, etc.) are given in Table 1.

In Figure 1, we plot the column density of total hydrogen that is required to produce a neutral hydrogen column $`N_{\mathrm{HI}}=10^{19}\mathrm{cm}^2`$ as a function of the total hydrogen density, when the model CHVCs are exposed to our fiducial ionizing photon flux. The solid curve is for our fiducial ionizing photon flux and the dashed curve is for $`\varphi _i`$ reduced by a factor of two. At the low density end of the range, the ratio NH/$`N_{\mathrm{HI}}`$ is enormous, $`100200`$, in consequence of the small neutral fraction of gas at this density even if exposed only to the extragalactic ionizing background. At the highest densities, this ratio has declined to $`2.53`$, as the gas is substantially neutral over much of the total column.

Since we know both the neutral and total hydrogen column densities as well as the volume densities, we can calculate the gas masses directly from the models. However, in order to compare the results to observations, we need to assign distances to the model clouds. To do this, we adopt the following procedure. For each model, we determine the radius at which the projected H I column drops to the HIPASS sensitivity limit of $`N_{\mathrm{HI}}=2\times 10^{18}\mathrm{cm}^2`$. We then set the distance by requiring that the apparent angular size (using this radius) matches the typical value found for the CHVCs in the HIPASS survey, for which the angular radius $`\mathrm{\Delta }\theta 0.34`$ degrees. This gives the result $`d_{\mathrm{HVC}}=169R_{18.3}`$ kpc, where $`R_{18.3}`$ is the model radius (in kpc) at the $`N_{\mathrm{HI}}=2\times 10^{18}\mathrm{cm}^2`$ level.

In Figure 2 we show both the total gas masses $`M_{gas}`$ and the apparent masses $`M_{\mathrm{HI}}`$, the latter derived from the neutral hydrogen column densities and the apparent (as measured in H I) sizes, as a function of the cloud distance and total gas volume density, for our fiducial ionizing photon flux. As is evident from Table 1, models in which the total gas density is less than $`n_\mathrm{H}3\times 10^3`$ are ruled out, as such clouds are so large that they are required to lie outside of the Local Group, and the gas masses become absurdly large. At the maximum acceptable distance ($`d1`$ Mpc) the total gas mass is an order of magnitude larger than the apparent H I mass.

Figure 3 shows the same quantities as in Figure 2, but we now show the results only for clouds within the Local Group ($`d1`$ Mpc). In addition, this figure also shows the gas masses for $`\varphi _i`$ a factor of two smaller than our assumed value. For the lower photon flux, clouds of a fixed density are at smaller radius compared to the fiducial model. At a fixed distance, the total gas masses are smaller by about a factor of $`1.5`$ for the reduced flux models. For model CHVCs at Local Group distances ($`d\stackrel{>}{}0.5`$ Mpc), the gas masses are $`M_{gas}10^710^8M_{}`$.

There is an additional constraint that we can apply to the CHVC models. The linewidth distribution of CHVCs is sharply peaked at a FWHM $`\mathrm{\Delta }V25\mathrm{km}\mathrm{s}^1`$ (de Heij, Braun, & Burton 2002). We have therefore calculated the expected line FWHM for the models (i.e., using the velocity dispersion for a self-gravitating system), both with and without a dark matter component (in the former case the total mass has simply been scaled up by the assumed ratio of dark to baryonic mass); $`\mathrm{\Delta }V`$ has been calculated assuming a uniform spherical mass distribution. In Figure 4 we show the predicted line FWHM as a function of distance, for dark matter to gas mass ratios of ten and zero. In the former case, the models are consistent with the observations for a distance $`d330`$ kpc; with no dark matter the clouds must lie at approximately 1.5 Mpc. For $`\varphi _i`$ reduced by a factor of two, these distances increase somewhat, to $`d430`$ kpc and 2.1 Mpc for $`M_{dm}/M_{gas}=10`$ and no dark matter, respectively.

While these models rule out the possibility that the CHVCs sample a population of cosmological objects (i.e., gas in dark matter-dominated potential wells) at distances characteristic of the Local Group ($`d1`$ Mpc), they suggest that they could represent such objects at distances of $`300400`$ kpc. However, as we will see in the next section, models of CHVCs in realistic dark matter potential wells require that they lie much closer to the Galaxy, with an upper limit to the distance of about 200 kpc. This distance scale rules out models in which the CHVCs represent continuing infall onto the Local Group.

### 4.2 Dark matter halo models

Motivated by the order-of-magnitude agreement between the properties of model clouds with $`M_{dm}/M_{gas}10`$ and the observations for characteristic distances $`d350`$ kpc, we have constructed models of CHVCs using realistic dark matter potential wells. Numerical simulations of halo formation in cold dark matter cosmologies predict a specific form for the density profile, as initially found by Navarro, Frenk, & White (1996). However, there is considerable disagreement as to whether this profile describes the halos of real galaxies, particularly dwarf and low-surface brightness galaxies (e.g., Flores & Primack 1994; Moore 1994; Burkert 1995). A number of authors have argued that these latter objects have nearly constant-density cores rather than the cusped power-law NFW density profile, but this point is controversial. We have therefore adopted two different models for the halo profile. The first is the the Navarro, Frenk, & White (1996) halo density profile, but modified to allow for the possible presence of a density core:

$$\rho _{dm}(r)=\rho _{crit}\frac{\delta _cr_s}{(r+r_o)(1+r/r_s)^2}$$
(1)
where $`\rho _{crit}=3H^2/8\pi G`$ is the closure density, the scale length $`r_s`$ is related to the virial radius of the halo by the concentration parameter $`c=R_{vir}/r_s`$, $`r_o`$ is the core radius, and $`\delta _c`$ is a characteristic density contrast, determined by the requirement that the halo mass

$$M_h(R_{vir})=\mathrm{\Delta }\rho _{crit}(4\pi /3)R_{vir}^3,$$
(2)
where $`\mathrm{\Delta }`$ is the overdensity parameter. For a flat ($`\mathrm{\Omega }+\mathrm{\Lambda }=1`$) universe, $`\mathrm{\Delta }`$ is given by (e.g., Eke, Navarro, & Frenk 1998) $`\mathrm{\Delta }=178\mathrm{\Omega }^{0.45}111`$ at $`z=0`$ for the $`\mathrm{\Omega }_{dm}=0.35`$ cosmology adopted here. For the relevant mass range of halos, $`c10`$.

In the limit of no core, the density contrast is given by the standard result (e.g., Eke, Navarro, & Steinmetz 2001)

$$\delta _c=\frac{\mathrm{\Delta }}{3}\frac{c^3}{\mathrm{ln}(1+c)c/(1+c)}.$$
(3)
Since we have no physical basis for choosing the core radius (this point is discussed further below) we have kept the concentration parameter fixed (for a given halo mass) at the no-core value, and calculated $`\delta _c`$ from the requirement that the integral of equation (1) to $`R_{vir}`$ equal the mass given by equation (2). The NFW halo parameters ($`c`$, $`R_{vir}`$) for a given mass have been calculated using the results of Mo, Mao & White (1998) and Eke et al. (2001).

The second halo mass profile is that suggested by Burkert (1995) based on observations of dwarf galaxies,

$$\rho _{dm}(r)=\frac{\rho _or_o^3}{(r+r_o)(r^2+r_o^2)}$$
(4)
where $`\rho _o`$ is the core density and $`r_o`$ is again the core radius. Like the NFW profile, the density of the Burkert halo falls as $`r^3`$ at large radius. Burkert found that the core density and radius were tightly correlated for the observed galaxies, so that the dark matter profiles are described by only one free parameter, which can be taken to be $`r_o`$. Burkert also concluded that the halo virial radius $`R_{vir}3.4r_o`$, so that these halos are much less centrally concentrated than the NFW halos.

We have calculated the gas density profile within these halos, including the effects of the gas self-gravity. The dark matter density profile is fixed, however: we do not include the response of the halo to the gas. (This would only strengthen our conclusions on the viability of models with dark matter.) In most cases, the gas velocity dispersion is taken to be $`\sigma _g=10.6\mathrm{km}\mathrm{s}^1`$, which produces a line FWHM equal to the typical observed value $`\mathrm{\Delta }V=25\mathrm{km}\mathrm{s}^1`$, and is the expected value for gas at $`T10^4`$ K, characteristic of the warm neutral medium. For NFW halos, in the case of no core and negligible gas self-gravity, the gas density profile is given by

$$\rho _g(\stackrel{~}{r})=\rho _g(0)e^{C_h\left[\mathrm{ln}(1+\stackrel{~}{r})/\stackrel{~}{r}1\right]}$$
(5)
where $`\stackrel{~}{r}=r/r_s`$ and

$$C_h=\frac{4\pi G\rho _{crit}\delta _cr_s^2}{\sigma _g^2}=7.34\times 10^3\frac{\delta _cr_s^2}{\sigma _g^2}$$
(6)
for $`r_s`$ in kpc and $`\sigma _g`$ in km s-1. For Burkert halos the gas profile is not analytic even in the absence of self-gravity.

We first discuss the results of models with NFW halos, and then the Burkert halo models.

#### 4.2.1 NFW halo models

The uniform density models that are consistent with a dark matter to baryon ratio of about 10 have total hydrogen volume densities of approximately $`3\times 10^3\mathrm{cm}^3`$ and total gas masses of about $`M_{gas}1.4\times 10^7M_{}`$. We have therefore calculated the gas distribution and the neutral fraction for a halo with a mass $`M_{dm}=1.4\times 10^8M_{}`$. The virial radius of this halo is 13 kpc and the concentration parameter $`c=12.44`$. We have also allowed for tidal truncation of the halo due to the Galaxy’s tidal field; this introduces a cutoff at a radius of 10.2 kpc.

In Figure 5 we show the gas density profile for the model that satisfies the requirement that the peak $`N_{\mathrm{HI}}=10^{19}\mathrm{cm}^2`$. This profile (which is characteristic of all the halo models) immediately reveals the problem that plagues all of these models. Most of the gas mass in the halo is at large radius, where the density is low and therefore the neutral fraction is very small. Hence, in order to match the observed neutral column density, the density at the center of the halo must be raised substantially over that in the uniform density model. As Figure 6 (which plots the cumulative dark matter and gas masses as a function of radius) demonstrates, this large central gas density results in a gas to dark matter mass ratio of approximately unity within the halo, which is obviously untenable. This conclusion is unaltered by reducing the ionizing flux by a factor of two or by allowing the dark matter profile to have a core. Although this latter modification results in gas density profiles that are more nearly constant over a larger range in radius, and therefore produce the same total hydrogen column density for a smaller central density, the slower drop-off of $`n_\mathrm{H}`$ with $`r`$ due to the core means that these models always have a higher value of $`M_{gas}/M_{dm}`$ at the outer boundary than the no-core models. (In fact, the models with cores invariably require larger total hydrogen columns than the models without cores, which makes the problem even worse.)

The only way to produce reasonable baryon to dark matter ratios in these models is to increase the mass of the halo. In order to produce a model in which the ratio $`M_{gas}/M_{dm}`$ is less than about 0.1222Assuming that galaxy clusters provide a fair sample of the universe, the baryon to dark matter mass ratio is estimated to be 0.13 for $`H_o=70`$ km s-1 Mpc-1 (Mohr, Mathiesen, & Evrard 1999). at the outer boundary prior to tidal truncation, the halo mass must be $`M_{dm}\stackrel{>}{}5\times 10^8M_{}`$. The mass increase is smaller than the factor of ten one might expect from Figure 6 because higher mass halos are both larger than lower mass ones and denser at the same physical radius, and therefore the observed neutral hydrogen column can be produced with a smaller total gas mass than is required in the lower mass halo.

In addition to the mass problem, to which we return below, these NFW dark matter-dominated models suffer from another, severe problem for the cosmological hypothesis: the physical size of the neutral hydrogen distribution is so small that the model CHVCs are restricted to lie less than $`d200`$ kpc from the Galaxy. This is illustrated in Figure 7, which shows the projected H I column as a function of impact parameter for a variety of halo masses, core radii, and ionizing fluxes. In no case does the model CHVC radius at the HIPASS threshold, $`R_{18.3}`$, significantly exceed 1 kpc. As discussed in §4.1, this constrains the clouds to distances of no more than about 200 kpc from the Galaxy. In fact, except for the lowest mass ($`M_{dm}=1.4\times 10^8M_{}`$) model shown in Figure 7, which, as discussed above, has an unacceptably large gas to dark matter mass ratio, all of the models have $`R_{18.3}\stackrel{<}{}0.5`$ kpc unless we include a core in the halo density distribution. As we noted earlier, we have no physical basis for choosing the halo core parameters, and have included a core simply to see whether, given the uncertainties in the actual halo dark matter profiles, this could alleviate the problems arising from the small distance scale mandated by the small physical size of the model CHVCs. As is evident from Figure 7, it does not. The reason is simply that if the halo is massive enough that the observed neutral hydrogen column density can be reproduced with a reasonable baryon to dark matter mass ratio, the velocity dispersion of the baryonic component (which is fixed by the observed linewidth) is only a fraction of the velocity dispersion characterizing the dark matter potential well, with the result that the gas is confined to the core of the halo. Increasing the size of the core beyond $`r_or_s`$ does not lead to any significant further increase in $`R_{18.3}`$, as is seen in the right-hand panel of Figure 7.

#### 4.2.2 Burkert halo models

In Figure 8 we show results for models with a Burkert halo of mass $`M_{dm}=10^8M_{}`$. For these models we have also explored a larger region of parameters for the CHVCs, considering neutral hydrogen columns up to ten times larger than the typical value seen in the HIPASS survey and linewidths up to a factor of two larger than the typical value. The left-hand panel assumes the standard velocity dispersion of $`\sigma _g=10.6`$ km s-1. From left to right, we plot the neutral hydrogen column as a function of impact parameter for peak neutral hydrogen column densities $`N_{\mathrm{HI}}=10^{19}\mathrm{cm}^2`$ (solid line), $`N_{\mathrm{HI}}=2.5\times 10^{19}\mathrm{cm}^2`$ (long-dashed line), $`N_{\mathrm{HI}}=5.0\times 10^{19}\mathrm{cm}^2`$ (dotted line), $`N_{\mathrm{HI}}=7.5\times 10^{19}\mathrm{cm}^2`$ (short-dashed line) and $`N_{\mathrm{HI}}=10^{20}\mathrm{cm}^2`$ (solid line). This Burkert halo model does not suffer from the unacceptably high baryon to dark matter mass ratio problem that afflicts the NFW halos of comparable mass, because the Burkert halo is physically much smaller (virial radius of 1.8 kpc compared to 11.7 kpc for the same mass NFW model) and therefore its mean density is much higher than the corresponding NFW halo. However, this small size scale becomes a serious problem just as it does for the NFW halos. Even for a peak column $`N_{\mathrm{HI}}=10^{20}\mathrm{cm}^2`$, the projected size at the HIPASS sensitivity threshold barely reaches 1 kpc, again implying that the halos are no more than 200 kpc distant. The right-hand panel shows the same models, only for a doubled velocity dispersion, so that the line FWHM $`\mathrm{\Delta }V=50\mathrm{km}\mathrm{s}^1`$. Even with this assumption, the physical size of the clouds remains small: raising the peak column density to $`N_{\mathrm{HI}}=5.0\times 10^{19}`$ \- five times the typical HIPASS value - only increases the cloud radius to about 1.6 kpc, which would push the clouds out to no more than 270 kpc. (Lower mass halos - e.g., $`10^7M_{}`$ \- with Burkert profiles are so small that they can lie no further than $`d100`$ kpc from the Milky Way, independent of total column density.)

In Figure 9 we present similar results, but for a halo with $`M_{dm}=10^9M_{}`$. In the left-hand panel, with a velocity dispersion corresponding to the typical observed value, the results scarcely differ from those for the $`10^8M_{}`$ halo, and the cloud size at the HIPASS threshold again never exceeds 1 kpc. The results for the models with twice the linewidth are also very similar to those for the lower mass Burkert halo. The chief difference is that, since the $`10^9M_{}`$ halo is physically larger (by a factor of 2.7) than the $`10^8M_{}`$ halo, the highest column density models are larger in the higher mass model, since the column density cutoff is not set by the physical edge of the halo as it is in the $`10^8M_{}`$ model.

Hence the CHVC models with Burkert halos suffer from the same distance problem as the NFW halo models: the sizes at the HIPASS sensitivity threshold are simply too small to allow them to be placed at distances much larger than $`d200`$ kpc from the Milky Way. Our results are in excellent agreement with those of Kepner, Babul & Spergel (1997) and Kepner et al. (1999), who calculated the chemical, thermal, and ionization equilibrium of gas in dark matter minihalos with Burkert density profiles. They found that the characteristic radius at the $`N_{\mathrm{HI}}10^{18}\mathrm{cm}^2`$ level never exceeded $`r1`$ kpc, even when the central hydrogen column density reached values much larger than any we consider here.

The combination of these two factors (the large mass that is required in order to have plausible gas to dark matter ratios for NFW halos, and the close distances required by the small physical size of the clouds for either NFW or Burkert halos) is fatal for any model in which the CHVCs represent a population of cosmological objects. In Figure 10 we show the halo mass distribution for a cosmological simulation of the Local Group (Moore et al. 2001). The $`2000`$ halos were created by analyzing the formation of a binary pair of massive halos in a hierarchical universe dominated by cold dark matter. The masses, separation and relative velocities of the binary halos are close to those observed for the Milky Way and Andromeda. The binary system was chosen from a large cosmological simulation such that a nearby massive cluster similar to Virgo was present. The simulation is described in detail in Moore et al. (2001). The filled circles show all of the halos (aside from those representing the Milky Way and Andromeda galaxies) within the simulated Local Group volume and the open circles show the halos that lie within 200 kpc of the Galaxy. As we showed above, the minimum acceptable mass for CHVC halos is $`M_{dm}5\times 10^8M_{}`$. From this mass limit up to $`M_{dm}=2\times 10^9M_{}`$, there are 16 halos in the simulation333For a lower mass limit of $`M_{dm}10^8M_{}`$, as allowed for Burkert profile halos, the number of halos would be a few times larger. However, except for extreme assumed CHVC parameters (i.e., peak column densities and linewidths much larger than the typical values), the CHVC models with Burkert halos will be physically even smaller than those with NFW halos. Hence these objects would have to lie closer than 200 kpc, and so the total number predicted would be about the same as for the NFW-halo models. To explain the CHVCs as cosmological, dark matter-dominated objects, we need approximately 250. Hence unless the hierarchical, cold dark matter-dominated cosmology grossly under-predicts the number of halos in this mass range, our results do not allow CHVCs to represent such a population.

### 4.3 CHVC HI distribution

There have now been several detailed studies of the HI structure of CHVCs. Wakker & Schwarz (1991) investigated two CHVCs and were the first to notice the core/halo structure of the clouds. Braun & Burton (2000) imaged six CHVCs at 1 resolution and found compact cores with linewidths as low as 2 - 10 km s-1, as well as the large diffuse halo modeled here with linewidths of $`25`$ km s-1. The warm ($`10^4`$ K) halos were imaged at 3 resolution with Arecibo by Burton, Braun & Chengalur (2001; hereafter BBC).444 Of the ten CHVCs imaged by BBC at Arecibo, 8 are covered by the northern extension of the HIPASS survey. We have examined the HIPASS data for these clouds to see whether the HIPASS results have been seriously compromised by resolution effects. With the exception of the faintest cloud in the sample (for which the Arecibo peak column is about 4.5 times the HIPASS value: however, the column density peaks are due to probably unresolved, non-centered structure, and the entire western half of the cloud has a typical column extremely close to the HIPASS value), the Arecibo and HIPASS peak column densities differ by at most a factor of two. Hence there is no reason to think that the HIPASS observations systematically underestimate the cloud column densities significantly enough to affect the results of this paper. They argued that the column density distribution at the edges of the clouds (for $`N_{\mathrm{HI}}`$ $`<10^{18.5}`$ cm-2) drops off as an exponential with radius, indicating a spherical exponential distribution of neutral hydrogen volume density as a function of radius. (Figures 11 and 12 of BBC demonstrate the quality of the exponential fits.) Kinematically the cores and halos of the CHVCs do not appear to be related. BBC derived distances of typically a few hundred kpc for the clouds in their sample, based on the following assumptions:

* The neutral gas density distribution is reasonably well-described by an exponential with radius, so that the central density can be inferred from the peak column density and the fit scalelength;
* The gas neutral fraction is close to unity in the center of the CHVC, and the temperature $`T10^4`$ K;
* The total gas pressure at the center is $`P/k100`$ cm-3 K, which Braun & Burton (2000) argued is the characteristic pressure at which the cold cores can co-exist with the warm neutral medium.

With these assumptions, the angular scale length of the exponential can be converted to a physical scale length, thereby yielding a cloud distance.

In this section we show that the characteristic exponential $`N_{\mathrm{HI}}`$ profile arises as a simple consequence of the physics of photoionization, and in the absence of fine-tuning of conditions, is very unlikely to provide reliable distance estimates.

In Figure 11 we plot the projected neutral hydrogen column density as a function of impact parameter (normalized to the cloud radius) for five different uniform density models, with total densities $`n_\mathrm{H}`$ between $`10^4`$ and $`10^2\mathrm{cm}^2`$. For all of these models the $`N_{\mathrm{HI}}`$ profile resembles an exponential over most of the cloud radius, with the exception of a core which becomes more pronounced with increasing density, and of a cutoff as the cloud outer boundary is approached. The exponential appearance of the H I column density distribution indicates that the neutral hydrogen volume density distribution is also approximately an exponential, at least over a substantial fraction of the cloud volume. This is simply a consequence of the photoionization physics, and reflects the attenuation of the ionizing radiation (and therefore the increase in the neutral fraction) with depth into the clouds; the total gas volume densities in all these models are uniform.

BBC fit $`N_{\mathrm{HI}}`$ profiles of the form

$$N_{\mathrm{HI}}(p)=2n_oh\left(\frac{p}{h}\right)K_1(p/h)$$
(7)
to their observations, where $`p`$ is the impact parameter, $`n_o`$ is the central density, $`h`$ is the scale-length, and $`K_1`$ is a modified Bessel function. This is the form expected for an exponential volume density distribution of infinite extent (van der Kruit & Searle 1981). This is modified for a cloud of finite extent, however, and so instead we use the finite-cloud profile

$$N_{\mathrm{HI}}(p)=2n_oh\left(\frac{p}{h}\right)\left[e^{R/h}\left(R^2/p^21\right)^{1/2}+k_1(p/h,R/p)\right]$$
(8)
where $`R`$ is the cloud radius and $`k_1`$ is an incomplete Bessel function,

$$k_1(\rho ,x)=\rho _1^xe^{\rho t}\left(t^21\right)^{1/2}𝑑t$$
(9)
which goes to the usual $`K_1(\rho )`$ as $`x\mathrm{}`$.

In Figure 12 we show two examples of fits of equation (8) to the projected neutral hydrogen column density distributions for the CHVC models described in §4.1. These are not formal fits to the distributions, and are intended solely to show that these uniform density models are very well described over a large range in radius (in fact, everywhere except for the cores) by an exponential $`n_{\mathrm{HI}}`$ distribution. Although BBC argued based on observations of dwarf galaxies that there is a characteristic value of $`h`$, that conclusion does not hold for these models, as $`h`$ depends nontrivially on the volume density. For the pair of models shown in Figure 12, which differ in $`n_\mathrm{H}`$ by a factor of two, the fitted values of $`h`$ differ by nearly a factor of four.

The BBC distance method generally yields erroneous results when applied to these models, for two reasons: (1) many of the CHVC models have substantial ionized fractions even in their cores, so that the central neutral density - even if determined with reasonable accuracy from the exponential fit - can be substantially below the actual total gas density, and (2) the central pressure is generally not equal to $`P/k100`$ cm-3 K. 555For the $`n_\mathrm{H}=0.01\mathrm{cm}^3`$ model, for which $`P/K`$ does equal 100 at cloud center, the BBC distance method gives a result within about $`25\%`$ of the value assigned to the cloud ($`D=67`$ kpc, Table 1) in order to make its angular size match the typical HIPASS size, provided that the true central column density is used, rather than the fitted column density, since in these relatively high-density models the presence of a core means that the exponential fit overestimates the central column density - cf. Figures 11 and 12. It is not at all obvious that the value of $`P/k100`$ quoted by Braun & Burton (2000) for the co-existence of cold gas with the warm neutral medium is generally applicable. This conclusion was based on a phase diagram calculated by Wolfire et al. (priv. comm.) for “Local Group conditions”, which are not otherwise specified. In the calculations of Wolfire et al. (1995a,b), which presumably use similar physics, gas heating is dominated by grain photoelectric heating, and the only source of ionizing photons is the extragalactic soft X-ray background, which produces an ionizing flux about an order of magnitude smaller than assumed here. Furthermore, the resulting phase diagram is rather sensitive to the adopted gas parameters (see Figures 1 and 6 in Wolfire et al. 1995b, in which the range of pressures over which a two-phase medium is allowed can extend over three orders of magnitude). In fact, in the phase diagram presented by Braun & Burton (their Figure 13), the minimum pressure $`P_{min}`$ for which the cold phase exists is approximately five times larger ($`P_{min}/k300`$ cm-3 K) for a CHVC warm gas column density $`N_{\mathrm{HI}}=10^{19}\mathrm{cm}^2`$, as appropriate for the HIPASS CHVCs discussed here, than it is for $`N_{\mathrm{HI}}=10^{20}\mathrm{cm}^2`$ (for which $`P_{min}/k`$ is actually about 60), as adopted by Braun & Burton for their sample. This would reduce the inferred distances by a factor of three. It is also not yet clearly established that the presence of cold gas is a generic property of CHVCs: there are only six clouds in the Braun & Burton (2000) sample, and while all of these exhibit cold cores, a much larger sample is clearly needed. The Magellanic Stream shows no evidence for cold gas (Mebold et al. 1991).

Interestingly, the models in which the gas is confined by a dark matter potential are not well described by profiles of the form (7) or (8). We show one such NFW model in Figure 13. This is the $`M_{dm}=5\times 10^8M_{}`$, no-core model shown as the solid line in the center panel of Figure 7. Inspection of that figure shows that the shape of the neutral hydrogen column density distribution in these dark matter-dominated models is generic. As is evident from Figure 13, these models cannot be described by a single component of the form (7) or (8). This is because the neutral hydrogen column density always exhibits a core with an extended wing, in consequence of the gas density profiles in these halos (cf. Figure 5). Hence the profile fit to the core of the $`N_{\mathrm{HI}}`$ distribution (we have used equation (7) for these models, since the models are actually much larger than the plotted region, and hence the cutoff is unimportant here), shown by the dashed line in Figure 13, grossly under-predicts $`N_{\mathrm{HI}}`$ outside $`r0.5`$ kpc, while the fit to the column density distribution at larger radius (dotted line) falls an order of magnitude short of $`N_{\mathrm{HI}}`$ at $`r=0`$. The results for Burkert profile halos are similar; in Figure 14 we show the neutral column density distribution for the largest column density model ($`N_{\mathrm{HI}}(0)=1.0\times 10^{20}`$ cm-2) shown in the lefthand panel of Figure 8 ($`M_{dm}=10^8M_{}`$, $`\sigma _g=10.6`$ km s-1). As for the NFW halo model, there is no reasonable fit with a single component of the form (7) or (8); in fact, the core-halo structure is even more pronounced than in the NFW halo. As before, fits to the core of the $`N_{\mathrm{HI}}`$ distribution (here we had to use equation ) severely underestimate $`N_{\mathrm{HI}}`$ beyond $`r0.8`$ kpc, while fits to $`N_{\mathrm{HI}}`$ at large radius underestimate the column density at small radii by an order of magnitude.

## 5 Summary and Discussion

Models in which the compact high-velocity clouds (CHVCs) lie at distances of hundreds of kpc to $``$1 Mpc from the Galaxy have such low volume densities of hydrogen that exposure to the extragalactic ionizing background radiation alone will largely ionize them. For models of uniform density clouds (§4.1), we find that the total hydrogen column needed to produce the typical neutral hydrogen column of $`N_{\mathrm{HI}}`$$`=10^{19}\mathrm{cm}^2`$ ranges between 3 and 200 times $`N_{\mathrm{HI}}`$ (Figure 1) as the total hydrogen density varies between $`10^2`$ and $`10^4\mathrm{cm}^3`$. The total gas masses of these clouds are therefore much larger than the apparent neutral gas masses. To match the typical angular size seen for CHVCs in the HIPASS survey, low-density clouds must be at large distances and have correspondingly large masses. Models in which the cloud hydrogen densities are less than $`n_\mathrm{H}23\times 10^3\mathrm{cm}^3`$ are ruled out, as such clouds would have to lie outside the Local Group (cf. Figures 2 and 3).

At Local Group distances, $`d0.71`$ Mpc, the gas mass alone of the individual CHVCs would be $`M_{gas}10^8M_{}`$, and the mass of the observed population would be $`M4\times 10^{10}M_{}`$. Such CHVCs could have dark matter/gas mass ratios of only $`23`$ without producing line widths larger than the observed value. Models in which the CHVCs are at smaller distances, $`d350`$ kpc, with gas masses $`M_{gas}10^7M_{}`$, appear to be more reasonable, and could be dynamically compatible with plausible ($`10`$ to 1) dark matter to gas mass ratios.

A closer examination of this scenario, however, using realistic dark matter halos, shows that it runs into severe difficulties (§4.2). This is because: (1) the sizes of the detectable H I clouds are too small ($`R\stackrel{<}{}1`$ kpc), which requires that the CHVCs are at distances $`d\stackrel{<}{}200`$ kpc (this is true whether we use NFW profiles or Burkert profiles for the halo dark matter density distribution); and (2) the halo masses required for reasonable ratios of $`M_{dm}/M_{gas}`$ for NFW halos are too large ($`M_{dm}10^9M_{}`$). Hence the predicted number of such objects in a hierarchical CDM structure formation scenario is more than an order of magnitude smaller than the observed number of CHVCs, and CHVCs cannot represent a population of infalling, cosmological objects.

We have also examined the suggestion that the approximately exponential $`N_{\mathrm{HI}}`$ profiles seen in high resolution studies of a sample of CHVCs could provide an estimate of the cloud distances (§4.3). Such profiles are a natural consequence of an external source of ionizing photons for clouds with uniform total density. The derived scale lengths are not constant and depend on the cloud density. We argue that this technique cannot be used to derive distances to the CHVCs in a model-independent way, as the scale lengths are not fixed and the alternative assumption of a constant pressure for the CHVC population is debatable. The CHVC $`N_{\mathrm{HI}}`$ profiles can be explained by photoionization even if the clouds are near to the Galaxy ($`d\stackrel{<}{}100`$ kpc).

Our models have assumed that the gas is smoothly distributed, and therefore do not include substantial clumping. This is arguably a plausible assumption for the CHVC envelopes that we have modeled in this paper, as the gas in these envelopes presumably is in the warm ($`T10^4`$ K, as is indicated by the observed linewidths) neutral phase, and is supported by the generally smooth distributions and reflection symmetry seen in high-resolution studies of the envelopes (BBC). The only way in which this assumption can be substantially in error, and the gas highly clumped, is if the clumps are in the cold neutral phase, with temperatures of order 100 K. In this case the observed linewidths must reflect bulk motion rather than thermal velocities. Since the (thermal) linewidths of individual clumps would be $`\mathrm{\Delta }V2\mathrm{km}\mathrm{s}^1`$, a large number of clumps would need to be present along any given sightline in order to explain the smoothness of the velocity profiles. Such clumps would also have to be substantially smaller than the spatial resolution of the high resolution imaging studies, which frequently do find (and resolve) small numbers of cores with narrow linewidths (e.g., BBC). If pressures $`P/k100\mathrm{cm}^3`$ K are necessary to drive the gas into the cold neutral phase (BBC), then the hydrogen densities must be at least $`n_\mathrm{H}1\mathrm{cm}^3`$. It is rather difficult to see how to make such a density structure compatible with the observations. If the clumps are self-gravitating, then they must have column densities that are much larger than the beam-averaged column densities measured observationally, which requires that the areal filling factor is small even though the filling factor in velocity space must be large to reproduce the observed line profiles. If they are not self-gravitating, we are left with the question of what confines them (or generates them, if they are transient objects). Hence we consider it unlikely that the CHVC envelopes that we have considered in this paper actually consist of dense, cold clumps rather than relatively smooth warm neutral gas.

What possibilities does this leave us with for the nature of CHVCs? Clouds without dark matter located within the confines of the Local Group are not massive enough to be self-gravitating, and must be either transient or confined by some external medium. Either of these scenarios would seem to argue for the CHVCs being associated with the Galaxy rather than the Local Group, as both the relatively short dynamical timescale ($`t10^8`$ years) and the presence of a confining medium are much more easily understood if the CHVCs are connected to the flux of mass and energy through the Galactic halo. This suggestion is supported by the overall similarities of the CHVCs to the remainder of the HVC population (Putman et al. 2002). The total and neutral hydrogen gas masses of the entire observed CHVC population ($`250`$ objects) would be $`M_{gas}7\times 10^8M_{}`$ and $`M_{\mathrm{HI}}7\times 10^7M_{}`$, respectively, at $`d200`$ kpc, and scale nearly as $`d^2`$ for smaller distances.

We are grateful to several colleagues, especially Nick Gnedin and Jessica Rosenberg, for helpful discussions. Nick kindly provided the software (written by Oleg Gnedin) for calculating the properties of NFW halos. PRM is supported by the National Science Foundation under grant AST 99-00871 and by NASA through HST grant AR-08747.02-A. MEP is supported by NASA through Hubble Fellowship grant HST-HF-01132.01-A awarded by the Space Telescope Science Institute, which is operated by the Association of Universities for Research in Astronomy, Inc., for NASA, under contract NAS 5-26555.


text:
# Lattice dynamics and electron-phonon coupling in transition metal diborides

## Abstract

The phonon density-of-states of transition metal diborides TMB2 with TM = Ti, V, Ta, Nb and Y has been measured using the technique of inelastic neutron scattering. The experimental data are compared with ab initio density functional calculations whereby an excellent agreement is registered. The calculations thus can be used to obtain electron-phonon spectral functions within the isotropic limit. A comparison to similar data for MgB2 and AlB2 which were subject of prior publications as well as parameters important for the superconducting properties are part of the discussion.

preprint: APS/123-QED

The discovery of superconductivity in MgB2 at the unexpected high temperature of 39 K has renewed the interest in diboride compounds. The well known AlB2 structure ($`P6/mmm`$) is formed with quite a number of transition metals (TM) Buzea . However, a clear proof of superconductivity could not be established for any of these candidates. Obviously the influence of lattice faults on the transport properties is important. A unique feature for MgB2 is the depopulation of B 2$`p\sigma `$ bands due to a particular Mg-B interaction Kortus . It has been shown by band structure calculations for the TM-diborides that $`E_F`$ moves up into the TM 4$`d`$-states with the consequence that these compounds behave more like normal metals Vajeeston . Most important for the transport properties is the electron-phonon interaction in these systems. Very recently results from point contact spectroscopy have been reported for TM = Zr, Nb, and Ta. The authors conclude on an only modest electron-phonon interaction parameter for these samples Naidyuk . A general investigation of the phonon spectra of TM-diborides has so far not be performed. In the present paper we confront measurements of the generalized phonon density of states for selected compounds with TM = Ta, V, Nb, Ti, and Y with ab initio density functional theory (DFT) calculations. Results for electron-phonon spectral functions are derived and placed in context to previous investigations of superconducting MgB2.

All diborides were prepared from stoichiometric mixtures of TM elements and amorphous 11B. We have chosen the less absorbing 11B in view of the intended inelastic neutron scattering (INS) experiments. The mixed powders were compacted to pellets and arc melted. The final polycrystalline samples showed metallic brightness and were found to be single phase by $`x`$-ray diffraction. All of our samples failed to show superconductivity above 4 K although transition temperatures for the pure elements of Nb and Ta with 9.2 K and 4.39 K are known. Our INS experiments were performed on the IN6 time-of-flight spectrometer at the HFR in Grenoble, France, with an incident neutron energy of 4.75 meV at 300 K in the upscattering mode. A high chopper speed of 201 rms and focusing in the inelastic region was used to improve the resolution. The generalized phonon density of states (GDOS) has been calculated from the recorded intensities integrated over a scattering region from 14 to 114. For the data evaluation we have applied multi-phonon corrections in a self-consistent procedure. The GDOS implies a weighting of vibrational modes of by $`\sigma /m`$ (scattering cross-section over the mass, see Tab. 1). For comparison the scattering power of 11B is 0.525 barn/amu.

Fig. 1 shows a comparison of GDOS spectra for all investigated TM-diborides. Each spectrum is normalized to 1. The most significant common feature is a gap around 40 to 50 meV. When analyzing the respective intensities we find that the high-frequency band is essentially composed of boron vibrations. The gap is thus a consequence of the large mass difference, which leads to a decoupling of transition metal and boron vibrations. A clearly different behavior is observed for MgB2 and AlB2 where due to the smaller mass of the metal atoms a hybridization of metal and boron modes occurs. In the low-frequency band of the transition metal diboride spectra we observe distinct shifts which however do not simply follow a $`\sqrt{m}`$ relation. This indicates that there are appreciable differences in the B-TM coupling. This conclusion is confirmed by the concomitant variations of interatomic distances (Tab. 1). Somewhat more surprising is the fact that significant frequency shifts are equally registered for the peak at the high-frequency end of the spectra where the lowest frequency is found in YB2. The E2g and B1g modes, which are important for superconductivity in MgB2, lie within this region where neighboring boron atoms move against each other. Strong electron-phonon coupling is expected to cause a significant broadening. This broadening is visible in the spectrum of MgB2 but not in those of the other compounds. The low frequency of the peak in YB2 alone can, therefore, not be taken as a signature of strong electron-phonon coupling. A simpler explanation is offered by comparatively weak B-B-interactions. As can be seen from Tab. 1 the B-B distance is exceptionally large for this compound. If any, then the $`G(\omega )`$ for NbB2 shows signatures of a stronger electron-phonon coupling, a compound for which some authors report superconductivity.

First principals density-functional calculations give parameter free insight into the electronic structure of all of these compounds. They in addition yield the lattice dynamical properties. Corresponding results are shown in Fig. 2. For the calculations we used the mixed basis pseudopotential method, which is described in some detail in Ref. Heid ; Bohnen . Norm-conserving pseudopotentials for V, Ta, Nb, Y and Ti were constructed according to the description of Bachelet-Hamann-Schlüter Bachelet ; Ho , whereas for boron a Vanderbilt-type potential was created Vanderbilt . The mixed basis scheme uses a combination of local functions together with plane waves for a representation of the valence states. The good agreement between calculated and experimental structural parameters (see Tab. 1) as well as between calculated and experimental phonon GDOS (Figs. 1 and 2) proves the validity of our description.

A very unique feature of electron-phonon coupling in MgB2 is the downshift of the in-plane E2g mode well below the out-of-plane B1g mode. Such an inversion of the usual sequence of mode frequencies is not encountered in any of the other diborides as shown in Tab. 2, where we compile the calculated $`\mathrm{\Gamma }`$-point frequencies of the boron vibrations.

The present calculations prove capable of describing the lattice dynamics of various diborides with quite different electronic transport properties. For the TM diborides, the isotropic character of the d-bands has the consequence that these compounds behave electronically very much like normal 3D metals Vajeeston . Although the band structures show strong similarities, the Fermi level $`E_F`$ moves up and down within the region of the transition metal $`d`$-bands. This variable band filling leads to appreciable differences in the various interatomic distances. The E2g frequencies exhibit an expected inverse correlation with the in-plane lattice constant, which indicates that the B-B covalent bond is electronically rather similar for all these compounds. On the contrast, the out-of-plane boron mode B1g is very insensitive to changes in the lattice constants.

Somewhat exceptional in our series of transition metal diborides is only YB2. Peculiarities in the electronic structure lead to a drastic increase of both the a and c lattice constants. As a consequence of the increased B-B bond length, the longitudinal force constant of the B-B bond is strongly reduced, and becomes even smaller than the transverse coupling. This is the reason for the observed softening of the high-frequency B-modes in YB2.

In order to get a more quantitative picture of the electron-phonon coupling within the series of diborides we have calculated the isotropic Eliashberg function $`\alpha ^2`$F($`\omega `$) which describes the phonon mediated pairing interaction (Fig. 3). Details of the computations are outlined in Ref. Heid . Again MgB2 proves to be singular. The extremely strong coupling at high-phonon frequencies is known to be connected to the E2g boron in-plane mode. It is found as a general trend that the region of strong coupling for TM diborides shifts down to lower energies as found for conventional superconductors thus leaving hope only for modest Tc’s. The existence of 2D- and 3D-type Fermi sheets derived mainly from boron $`\sigma `$ and $`\pi `$ bands, respectively, has been shown to be important for a quantitative understanding of the transport properties of MgB2 Tsuda . The 2D Fermi surface is connected to the hole doping of boron $`2p`$ bands and does not seem to be important for transition metal diborides.

We have also calculated the relevant parameters for the description of the superconducting properties within the Eliashberg formalism. The results are compiled in Tab. 3. A common feature unfavorable for superconductivity found in all transition metal diborides is the comparatively low $`\omega _{log}`$. This quantity represents the effective average frequency of the coupling modes and sets the energy scale for the pairing interaction. Its small value as compared to MgB2 indicates that the pairing interaction is mainly mediated by the TM vibrations and not by the boron modes. The isotropic coupling constant $`\lambda `$ shows a variation from weak to intermediate coupling strength. In Tab. 3, we have also included values for the superconducting transition temperature as obtained by solving the linearized isotropic gap equations. The Tc values depend strongly on the screening properties expressed by $`\mu ^{}`$ and drop appreciably for $`\mu ^{}>0`$. The values shown for a typical metallic screening of $`\mu ^{}=`$0.13 thus should only be taken as indicators of the general trend. Furthermore, as mentioned before, for MgB2 it is necessary to go beyond the isotropic limit and take into account the multigap structure introduced by the particular Fermi surface geometry for a proper quantitative description of the pairing state. Anisotropy is expected to play a minor role for the coupling in the low-frequency region of the transition metal vibrations that are dominated by the 3D-metal bands.

The values of $`\lambda `$ and $`\omega _{log}`$ are most favorable for MgB2. The next promising candidates according to the present calculations are TaB2 and NbB2, for which some authors find superconductivity at 9.5 K and 6.4 K, respectively Heid . Differences in stoichiometry might be able to explain different experimental results. As far as a comparison to recent point-contact measurements is possible we can state that the maxima in our $`G(\omega )`$ correspond well to the peaks found in the second derivative of the $`IV`$ characteristic. There is also agreement in a strongly reduced electron-phonon coupling strength for transition metal diborides with respect to MgB2. We, however, want to mention that the reported very small numerical values for $`\lambda _{PC}`$ should not be compared directly to our values since point contact experiments do not measure the McMillan $`\lambda `$.

To summarize, we have presented measurements of the generalized phonon density-of-states for various transition metal diborides. The observed trends in the lattice dynamics can be well ascribed by first principles calculations and do not exhibit indications of strong electron-phonon interaction. This view is supported by theoretical calculations of the electron-phonon coupling. While they do not exclude the possibility of superconductivity with low Tc mediated by TM vibrations, they underline the outstanding properties of MgB2 among the class of diborides.