text:
# Probability distribution of density fluctuations in the non-linear regime

## 1 Introduction

The current model for the formation of large-scale structure is based on the gravitational instability of small initial fluctuations in the energy density. These are assumed to have originated from quantum fluctuations during the inflationary epoch. The fluctuations form a random field, which is fully specified by the $`n`$-point joint probability density functions (hereafter PDFs). The one-point PDF is a useful tool for testing the statistical character of the initial conditions (i.e. whether they are Gaussian or others), since its evolution does not depend on the nature of the dark matter (Trimble 1987). So it would be rather convenient to develop an efficient and accurate procedure for obtaining the present PDFs given that for the initial conditions. Many attempts have been made at solving this problem. Hoffman (1987) used the Zel’dovich (1970) approximation for computing the variance of the density fluctuations in the mildly non-linear regime ($`\sigma ^21`$). However, he did not take into account the fact that evolution and filtering are non-commuting operations. To obtain the variance on comoving scale $`r`$, one should evolve the whole field (with all its structure) and then filter it on scale $`r`$. But in the analytical procedure followed by Hoffman the field is first filtered and then evolved. In the linear regime ($`\sigma 1`$) both procedures lead to the same result. However, for $`\sigma 1`$ the results obtained by these procedures are rather different.

A procedure for obtaining the present PDF given the statistics of the initial conditions, that accounts for non-commutability of evolution and filtering was first described by one of us (Betancort-Rijo 1991). This procedure was the result of taking into account the simple fact that overdense regions with present comoving size $`r`$ correspond to larger comoving regions in the initial field. However, in developing this idea a couple of subtle points were involved that were not properly dealt within that paper.

During more than a decade that has passed since these first attempts at dealing with the problem under consideration (the article by Betancort-Rijo was in fact written in 1988) much attention have been given to it. New methods have been developed for approximating the non-linear evolution (Buchert 1989; Bouchet et al. 1995) and for the computation of statistics of the evolved field. Bernardeau (1994a) developed a procedure for taking into account the effect of filtering on the moments of the PDF in the small variance limit and used it to compute the first few hierarchical amplitudes. Then, assuming that these quantities do not depend on $`\sigma `$, he obtained (Bernardeau & Kofman 1995) the present PDF for any value of $`\sigma `$. However, this assumption, although it holds for some initial conditions, is not universal (as we shall see) and, even when it holds, it is not easy to justify a priori.

Despite all the attention that has been devoted to the present problem, the idea described by Betancort-Rijo (1991), which is at the root of the possibility of developing a simple and direct (it is not necessary to obtain the vertex generation function) procedure for computing the present PDF in the non-linear regime analytically, has not been followed. The exception is an article by Padmanabhan & Subramanian (1993), but there the problem is treated in the same way as was by Betancort-Rijo (1991), thereby having the same drawbacks.

In this paper we shall show how to correct the errors included in previous work and give the correct procedure for obtaining the present PDF. To check the correctness of the procedure we compute the values of $`S_3`$ and $`S_4`$ in the low-$`\sigma `$ limit (that have been computed analytically by Bernardeau 1994a) for several power spectra. Our procedure should render the exact values of all hierarchical amplitudes in this limit, and this is what is found for these two amplitudes. We also study the dependence on $`\sigma `$ of these amplitudes and of the variance in the low limit and for finite values of $`\sigma `$. We then compute the PDF explicitly for several cases for which numerical simulations are available and find them to agree with the numerical results.

In §3 we briefly review the basic idea (fully described in Betancort-Rijo 1991) involved in the present approximation and show how to use it to derive correctly the present PDF given the initial one. In the second section we shall discuss the technique that we use for the treatment of the non-linear evolution. In the last section we present and discuss the results of several explicit calculations.

## 2 Non-linear evolution

The complete Zel’dovich approximation (CZA; Betancort-Rijo & López-Corredoira 2000) is an approximation to the non-linear evolution of density fluctuations depending only on three quantities, $`\lambda _i`$, defined so that that $`1\lambda _i`$ are the eigenvalues of the tensor $`\frac{\stackrel{}{u}}{\stackrel{}{q}}`$, where $`\stackrel{}{u}`$ stands for the peculiar velocity of matter as given by the linear theory in comoving units, i.e. the peculiar velocity $`\stackrel{}{v}=R(t)\stackrel{}{u}`$; $`R(t)`$ being the scale factor of the Universe and $`q`$ stands for the Lagrangian coordinates, which are equal to the initial comoving coordinates (Eulerian). This tensor we call the linear local deformation tensor. The full local deformation tensor is defined in the same manner but using the exact value of $`\stackrel{}{u}`$. This approximation is exact to the second order. An approximation exact to the fourth order may be obtained out of it but this implies introducing additional variables. It takes the form (Betancort-Rijo & López-Corredoira 2000)

$$(1+\delta )^1=\underset{i}{}(1r_i(\lambda )\lambda _i+f_1(\lambda )x_i+f_2(\lambda )u_i),$$
(1)
where $`f_1`$, $`f_2`$, $`r_i`$ are certain functions of the variables $`\lambda _i`$. The variables $`x`$, $`u`$ which are given at any point by certain integrals over the whole field (Betancort-Rijo 2001) but, for statistical purposes we only need their probability distribution for fixed values of the $`\lambda _i`$. In fact, to obtain the moments exactly to the third order and almost exactly to the fourth we can use expression (1) with $`u`$ equal to zero and $`x`$ equal to its average value for fixed values of the $`\lambda _i`$ $`(<x_i>(\lambda ))`$. This is the approximation that in principle we use. It has a dependence on the power spectra through the mean value of $`x`$ (see Betancort-Rijo & López-Corredoira 2000). However, the dependence is so small that may be neglected without practical lost of accuracy. So, what we do in practice is to use simply the CZA (expression (1) with $`x=u=0`$ ), which is strictly exact to the third order for a particular power spectra (such that $`<x_i>(\lambda )=0`$) and practically exact for all of them. The CZA gives the density fluctuation at a point in term of the values of the $`\lambda _i`$ at that point. When dealing with filtered fields, however, the relevant quantity is the mean density within certain region (of any prescribed shape). It may be shown that this mean density fluctuation is also given exactly at least to the third order by expression (1) with the $`\lambda _i`$ representing now the linear eigenvalues of the deformation tensor corresponding to the filtered field; that is, the deformation tensor derived from the peculiar velocity field corresponding to the growing mode of the filtered field. The $`x`$ variables are now defined in a different manner. For non-filtered fields the presence of this variables is due to the different rate of growth of the various contributions to the local tidal field. In the present case, however, these variables have an extra component due to the filtering process. But, despite of the different origin of the $`x`$ variables, expression (1) and the probability distributions for the $`x`$ retain their form at least up to third order. So, setting in this expression $`x`$ equal to its mean value for fixed $`\lambda _i`$ leads also to an approximation exact to third order (for statistical purposes) for filtered fields. There is, however, an important difference between the two mentioned contributions to the $`x`$, because the part due to the tidal field remains small for any power spectra while that due to the averaging process grows indefinitely with increasing power on small scale. This implies, as we shall find, that for power spectra with too much power on small scale the procedure presented here fails even for values of $`\sigma `$ smaller than one. However, even in these cases the results are exact to third order, but many orders contribute for finite but small $`\sigma `$. So, the steeper the power spectra the smaller the value of $`\sigma `$ at which the approximation start failing. This happens, however, only for very steep power spectra not very relevant to cosmology. To deal with these cases we should use expression (1) with the mean values of $`x`$, not simply the CZA.

## 3 PDF of the filtered field

The effect of filtering on the field of density fluctuations have been treated exactly by Bernardeau (1994). In this treatment one may, in principle, obtaining (previously to caustic formation) all moments of the PDF of the field of density fluctuation filtered in some scale as a series of powers of $`\sigma `$ and sometimes (when the corresponding Hedgeworth expansion converges) to reconstruct the PDF. However, even when this works it is complicated. The procedure described here to deal with this problem ,which is basically that in Betancort-Rijo (1991) is more immediate, allowing the direct evaluation of the PDF even for values of $`\sigma `$ well above unity. Furthermore, it provides a clearer visualization of the effect of filtering. We now describe the basic idea of our procedure using the CZA, although the procedure in itself is independent of the approximation used for the evolution.

Consider certain geometrical body, $`B_1`$, centered at some point within the initial field. By assigning to this point the mean value of the density fluctuation within $`B_1`$ and repeating the process for each point, we generate a new field which we call the field filtered within $`B_1`$. We may now obtain the linear peculiar velocity field corresponding to the growing mode of this field. Let $`\lambda _i(B_1)`$ be the proper values of the corresponding linear (linearly calculated) deformation tensor, and let $`B_2`$ be the body which is the result of transforming $`B_1`$ with this deformation tensor. The CZA, in the present context, must be written in the form:

$$1+\delta (B_2)=(\underset{i=1}{\overset{3}{}}(1r_i(\stackrel{}{\lambda }(B_1))\lambda _i(B_1)))^1.$$
(2)
where $`\delta `$ is the present value of the density fluctuations filtered with $`B_2`$ at the point under consideration. In our problem, the filtering body $`B_2`$ is a sphere of radius $`r`$ and $`B_1`$ is an ellipsoid with semiaxis $`r(1r_i(\stackrel{}{\lambda })\lambda _i)^1`$. Hereafter, we write simply $`\lambda _i`$ instead of $`\lambda _i(B_1)`$. It must be noted that in the actual evolution the present sphere does not transform (when evolved back) exactly into an ellipsoid. The real surface has the same quadrupolar component as the mentioned ellipsoid, but contains smaller wiggles. However, for the present problem those wiggles are completely irrelevant.

Equation (2) gives the best one-point expression for $`\delta `$ as given by the CZA approximation. By one-point expression we mean that the value of $`\delta `$ at a given point (the mean value within a $`B_2`$ centered at that point) is expressed only in terms of the $`\lambda _i`$ at that point (the mean value within $`B_1`$). To go beyond this approximation, we need to obtain the values of the $`\lambda _i`$ at different points and integrate over them, or, what is the same, to include the x variables described in the previous section.

It should be noted that in equation (2) the $`\lambda _i`$ depend on $`B_1`$, which, in turn, depends on the $`\lambda _i`$ (for a given $`B_2`$). It is possible to show that, as long as caustics do not form, for every point there is a unique ellipsoid (with certain axis and orientation) such that the local deformation tensor of the field filtered with the ellipsoid transforms that same ellipsoid into a sphere of radius $`r`$ ($`B_2`$ in general). This prescription uniquely determines the values of the $`\lambda _i`$ (hence, $`\delta `$) at any given point. These unique $`\lambda _i`$ (for a given $`B_2`$) we call $`\lambda _i^{}`$. We must now obtain the PDF for the $`\lambda _i^{}`$ at a randomly chosen point in the initial field. With this PDF and equation (2), the derivation of the PDF for $`\delta `$ is a standard problem.

For a field filtered with a given body (the same at every point; usually a sphere of given radius), the PDF for the $`\lambda _i`$ is given by (Doroshkevich 1970):

$$P(\lambda _1,\lambda _2,\lambda _3,\sigma )=\frac{5^{5/2}27}{8\pi \sigma ^6}(\lambda _1\lambda _2)(\lambda _1\lambda _3)(\lambda _2\lambda _3)$$

$$\times \mathrm{exp}\left[\frac{1}{\sigma ^2}\left(3(\lambda _1+\lambda _2+\lambda _3)^2\frac{15}{2}(\lambda _1\lambda _2+\lambda _2\lambda _3+\lambda _1\lambda _3)\right)\right];$$

$$\sigma ^2\delta _L^2=(\lambda _1+\lambda _2+\lambda _3)^2;\lambda _3<\lambda _2<\lambda _1.$$
(3)
$`\sigma ^2`$ is the variance of the initial field of density fluctuations, $`\delta _L`$, filtered with the given body and extrapolated linearly to the present time. The trace of the deformation tensor extrapolated linearly to the present time (is given the tensor used in the Zel’dovich approximation) by $`\delta _L`$. Betancort-Rijo (1991) and Padmanabhan & Subramanian (1993) took the PDF for the $`\lambda _i^{}`$ to be simply proportional to the PDF for the $`\lambda _i`$ (with $`\lambda _i=\lambda _i^{}`$) for the field filtered with the ellipsoid that will contract to a sphere of radius $`r`$; that is, with an ellipsoid with semiaxis $`r(1r_i(\stackrel{}{\lambda ^{}})\lambda _i^{})^1`$. Note that now the filtering body changes from one point to another, so $`\sigma `$ is no longer a constant but a function of the $`\lambda _i^{}`$.

$$\sigma ^2(\stackrel{}{\lambda ^{}},r)=\frac{1}{(2\pi )^3}W(F)^2P(K)d^3\stackrel{}{K};$$

$$W(F)=\frac{3[\mathrm{sin}FF\mathrm{cos}F]}{F^3}$$
(4)
with $`F=\sqrt{_iK_i^2E_i^2}r`$; $`E_i=(1r_i(\stackrel{}{\lambda ^{}})\lambda _i^{})^1`$, $`K=|\stackrel{}{K}|`$; $`P(K)`$ power spectra.

In practice, in the papers mentioned and, to some extent, here, instead of the full equation (4) a more convenient approximation is used. Since $`\sigma `$ depends mainly on the volume of the filtering body (this is exactly so for white-noise power spectra), only the dependence on this quantity is considered. We then have:

$$\sigma ^2(\stackrel{}{\lambda ^{}},r)=\sigma ^2(\stackrel{}{0},r(\underset{i=1}{\overset{3}{}}(1r_i(\stackrel{}{\lambda ^{}})\lambda _i^{})^{1/3}))\sigma ^2(r(\underset{i=1}{\overset{3}{}}(1r_i(\stackrel{}{\lambda ^{}})\lambda _i^{})^{1/3}))=\sigma ^2(r(1+\delta )^{1/3}),$$
(5)
where $`\sigma ^2(r)`$ is the variance of the linear field filtered with a sphere of radius $`r`$. However, if we want to obtain results which are exact to the fourth order,as we shall do in some cases, the small shape dependence of $`\sigma (\stackrel{}{\lambda },r)`$ must be taken into account. In these cases, we use the following excellent approximation to the actual shape dependence of $`\sigma `$:

$$\sigma ^2(\stackrel{}{\lambda },r)=\sigma ^2(r(1+\delta )^{1/3})Z^2(\stackrel{}{\lambda })$$

$$Z(\stackrel{}{\lambda })exp\left[\frac{B}{2}\underset{i<j}{}\left(ln\left(\frac{1\lambda _i}{1\lambda _j}\right)^2\right)\right]$$
(6)
To compute the constant $`B`$ in the exponential we may use the exact form of Z in the case when two of the $`\lambda `$ are equal. In this case and for power law power spectra (with index n) expression (5) reduces to:

$$Z(P,n)=\frac{P^{\frac{2n}{3}}_0^{\mathrm{}}_0^{\mathrm{}}(u^2+P^2s^2)^{\frac{n}{2}}(W(\sqrt{u^2+s^2}))^2s𝑑s𝑑u}{_0^{\mathrm{}}W^2(u)u^{2+n}𝑑u}$$

$$P\frac{1\lambda _3}{1\lambda _1};\lambda _1=\lambda _2\lambda _3$$
(7)
We find: $`B(n=2)`$=0.0486; $`B(n=1)`$=0.0499; $`B(n=0)`$=0.

The approximation used in the mentioned works (which use expression (5)) may be expressed in the form:

$$P^{}(\stackrel{}{\lambda ^{}},\sigma )=CP(\stackrel{}{\lambda },\sigma (r(1+\delta (\stackrel{}{\lambda }))^{1/3}))|_{\stackrel{}{\lambda }=\stackrel{}{\lambda ^{}}},$$
(8)
where $`P^{}`$ is the PDF of the $`\lambda _i^{}`$; $`\delta (\stackrel{}{\lambda })`$ is given by equation (2) (with x=u=o); $`\sigma `$ stands for $`\sigma (r)`$, and $`C`$ is a normalization constant. This constant has to be introduced because, unlike $`P`$, the integral of $`P^{}`$ over all $`\stackrel{}{\lambda ^{}}`$ space is not equal to unity, due to the dependence of $`\sigma `$ on $`\stackrel{}{\lambda }`$.

We now show why equation (8) is incorrect and give the correct one. To this end we consider an approximation in which $`\delta `$ depends only on $`\delta _L(\lambda _1+\lambda _2+\lambda _3)`$. We do this because the relevant fact at this point is basically the same with one and with three variables, but it is much easier to understand in the former case.

For definiteness, let us assume that:

$$1+\delta =\left(1\frac{\delta _L}{3}\right)^3,$$
(9)
an approximation that has been used in Betancort-Rijo (1991).

For Gaussian initial conditions, $`\delta _L`$, which corresponds to a fixed filtering scale, follows a normal distribution. To obtain the PDF for $`\delta _L^{}`$ (the value of $`\delta _L`$ within the sphere of radius $`r/(1\frac{\delta _L}{3})`$ that have contracted at present to the sphere of radius $`r`$) that, by means of equation (9), allows us to obtain the PDF for $`\delta `$, we proceed as follows. To say that at a given point $`\delta _L^{}`$ takes the value $`\delta _L^0`$ means that at this point the linear value of $`\delta `$ filtered with a sphere of radius $`r/\left(1\delta _L/3\right)`$ is just $`\delta _L^0`$. The accumulated probability, $`P^{}(\delta _L^{}\delta _L^0)`$, that at a point chosen at random in the initial field $`\delta _L^{}`$ is larger than or equal to some given value may be obtained immediately. To this end we only need noting that if a sphere of radius $`r/(1\delta _L^{}/3)`$ (with $`\delta _L^{}\delta _L^0`$) have to contract to a sphere of radius r and shell crossing does not happen, the sphere of radius $`r/(1\delta _L^0/3)`$ will have to contract to a radius smaller than $`r`$. So, the value of $`\delta _L`$ within the sphere with radius $`R/(1\delta _L^0/3)`$ must be larger than $`\delta _L^0`$, since if the value of $`\delta _L`$ within this sphere were just $`\delta _L^0`$ it would contract to a sphere of radius r. Thus, the accumulated probability under consideration must be equal to the probability that the value of $`\delta _L`$ within a sphere of radius $`r/(1\delta _L^0/3)`$ be larger than or equal to $`\delta _L^0`$. So, we have

$$P^{}(\delta _L^{}\delta _L^0)=\frac{1}{2}\mathrm{erfc}\left(\frac{\delta _L^0}{\sqrt{2}\sigma (\delta _L^0)}\right)$$

$$\sigma (\delta _L^0)=\sigma \left(r\left(1\frac{\delta _L^0}{3}\right)^1\right)$$
(10)
where erfc stands for the error complementary function. To obtain the probability distribution for the values of $`\delta _L^{}`$ (the value of $`\delta _L`$ within the sphere that will contract to a sphere of radius $`r`$) at a randomly chosen point we only need deriving the cumulative probability with respect to $`\delta _L^0`$.

$$P^{}(\delta _L^{})=\frac{d}{d\delta _L^0}P^{}(\delta _L>\delta _L^0)|_{\delta _L^0=\delta _L^{}};$$

$$P^{}(\delta _L^{})d\delta _L^{}=\frac{1}{\sqrt{2\pi }}e^{\delta _L^2/(2\sigma ^2(\delta _L^{}))}d\left(\frac{\delta _L^{}}{\sigma (\delta _L^{})}\right).$$
(11)
Note that since $`\sigma `$ depends on $`\delta _L^{}`$ we cannot take it out of the differential. So we cannot write:

$$P^{}(\delta _L^{})=CP(\delta _L,\sigma (\delta _L))|_{\delta _L=\delta _L^{}}$$

$$=C\frac{e^{\delta _L^2/(2\sigma ^2(\delta _L^{}))}}{\sqrt{2\pi }\sigma (\delta _L^{})},$$

contrary to the assumption in equation (8). With these considerations in mind, it is easy to obtain the correct distribution for the $`\lambda _i^{}`$. We need only to generalize the above argument to three dimensions. The key argument is that the variables $`\lambda _i^{}/\sigma (\stackrel{}{\lambda ^{}})`$ must follow the distribution followed by the variable $`\lambda _i/\sigma `$ for fixed filtering (when $`\sigma `$ is a constant). So we have:

$$P^{}(\stackrel{}{\lambda ^{}},\sigma )=\left(P(\stackrel{}{\lambda },\sigma (\stackrel{}{\lambda }))\sigma ^3(\stackrel{}{\lambda })\left|\frac{\frac{\lambda _i}{\sigma (\stackrel{}{\lambda })}}{\lambda _j}\right|\right)|_{\stackrel{}{\lambda }=\stackrel{}{\lambda ^{}}};$$

$$\sigma (\stackrel{}{\lambda })\sigma (\stackrel{}{\lambda },r).$$
(12)
The last factor in this expression is the Jacobian of the transformation from the variables $`\lambda _i/\sigma (\lambda _i)`$ to the variables $`\lambda _i`$. Note that the $`P^{}`$ obtained by means of this expression is automatically normalized.

Equation (12) is exact for any initial conditions (provided that the $`\lambda _i`$ determine the evolution), although we shall use it only for Gaussian initial conditions; that is, with $`P`$ given by equation (3).

$`P^{}`$ is the PDF of the $`\lambda _i^{}`$ (the $`\lambda _i`$ associated with the unique ellipsoid that will contract to the sphere of radius r centered at the given point) for a point chosen at random in the initial field; that is, in Lagrangian coordinates. However, we want to obtain the PDF for $`\delta `$ at a point chosen at random in the present “physical” coordinates (i.e. Eulerian coordinates). So, we must first obtain the PDF for the $`\lambda _i^{}`$ in Eulerian coordinates. To this end, we only need to multiply $`P^{}`$ by the Jacobian of the transformation from Lagrangian to Eulerian coordinates. For a non-filtered field this Jacobian is simply

$$\underset{i=1}{\overset{3}{}}(1r_i(\stackrel{}{\lambda })\lambda _i).$$
(13)
However, for filtered fields this expression is not correct. This can be appreciated by noting that the resulting distribution for the $`\lambda _i^{}`$ (the result of multiplying $`P^{}`$ by the Jacobian in question) is not normalized. For non-filtered fields, on the other hand, it is easy to check that the Jacobian given by equation (13) leads to a normalized PDF. The correct Jacobian must exhibit the following properties: it has to be zero at points where $`\delta `$ goes to infinity; it must reduce to expression (13) in the limit of small power on scales below $`r`$; it must render a normalized PDF with vanishing first moment for any power spectra. The following expression meets these requirements and is exact at least to first order

$$J(\lambda )=\underset{i}{}(1(1\gamma )r_i(\stackrel{}{\lambda })\lambda _i),$$
(14)
where $`\gamma `$ is certain spectral constant that tends to zero in the limit of small power on scales below $`r`$. We cannot say a priori whether this expression is exact to higher orders. However,it seems to be exact to second order (within the precision of our computations) and it gives excellent results to all orders. It must be noted that expression (14) is not the value of the Jacobian at a point with the given $`\lambda _i`$. This would be so if the evolution were determined solely by these quantities. However, we know that this it is not strictly true, so expression (14) should represent the mean value of the Jacobian over all points with a given value of the $`\lambda _i`$. To obtain $`\gamma `$ we have considered a fluctuation with spherical symmetry and the evolution of the region such that the mean density fluctuation within a sphere centered at points belonging to it is above some threshold. Using the spherical model with the mean density profile, we obtained the contraction of this region and by comparison with expression (14) we find $`\gamma `$. However, this derivation is lengthy and although it might be interesting in itself, it is not essential to the present purposes, since we may obtain that constant using any relevant analytical result available (i.e., the leading order value of $`S_3`$) and comparing it with the one obtained using expression (14). In either case, we find:

$$\gamma =\frac{1}{6}\frac{d[ln(\sigma ^2(r))]}{d[ln(r)]}$$
(15)
## 4 Results and conclusions

To check the accuracy of the approximation developed in this article, we shall present some explicit calculations and compare them with the results of the corresponding available numerical simulation and with exact results obtained with the perturbation theory. We shall consider Gaussian initial conditions with power spectra of the form $`P(k)=Ck^n`$.

In these cases, equation (5) takes the form:

$$\sigma (\stackrel{}{\lambda },r)=\sigma _0(1+\delta )^{\alpha /2};$$

$$\sigma _0\sigma (r);\alpha \frac{n+3}{3},$$
(16)
with $`\delta `$ given by the CZA

$$(1+\delta )^1=\underset{i}{}(1r_i(\lambda )\lambda _i)$$
(17)
To derive exact results we use for $`r_i(\lambda )`$ its exact expression (Betancort-Rijo & López-Corredoira 2000). To derive result for finite values, however, we find it much more convenient to use the following compact approximation

$$r_i(\lambda )=g(\frac{3}{2}(\lambda _J+\lambda _K),\mathrm{\Omega })$$

$$g(x,\mathrm{\Omega })=\left\{\begin{array}{cc}\frac{3}{x}\left(1\frac{x\mathrm{\Omega }^{2/63}}{1.66}\right)^{0.5533}\hfill & \text{;}x1.584\hfill \\ 1.55\hfill & \text{;}x>1.584\hfill \end{array}\right\}$$
(18)
Using this in expression (12), with $`P(\stackrel{}{\lambda },\sigma )`$ given by (3) and multiplying it by the Eulerian–Lagrangian Jacobian, $`J(\stackrel{}{\lambda },\alpha )`$, (expression (14)), we obtain the PDF for the $`\lambda _i^{}`$ corresponding to a point chosen at random in Eulerian coordinates, which we represent by $`G(\stackrel{}{\lambda ^{}},\sigma _0)`$:

$$G(\stackrel{}{\lambda ^{}},\sigma _0)=\left(\left(P(\stackrel{}{\lambda },\sigma (\stackrel{}{\lambda }))\sigma ^3(\stackrel{}{\lambda })\left|\frac{\frac{\lambda _i}{\sigma (\lambda )}}{\lambda _j}\right|\right)_{\stackrel{}{\lambda }=\stackrel{}{\lambda ^{}}}\right)J(\stackrel{}{\lambda ^{}},\alpha ).$$
(19)
The moments of the PDF of the present value of the density fluctuation top-hat filtered on scale r are then given by:

$$\delta ^n(\sigma _0)=_{\mathrm{}}^{\mathrm{}}𝑑\lambda _1^{}_{\mathrm{}}^{\lambda _1^{}}𝑑\lambda _2^{}_{\mathrm{}}^{\lambda _2^{}}𝑑\lambda _3^{}G(\stackrel{}{\lambda ^{}},\sigma _0)\delta ^n,$$
(20)
where $`\delta `$ is given by equations (17) with $`\stackrel{}{\lambda }=\stackrel{}{\lambda ^{}}`$, provided that each of the factors in (17) is positive and $`\delta `$ is smaller than 170. If this expression leads to a value of $`\delta `$ larger than this, we use the following expression:

$$\delta =170\left|\frac{\delta _L}{1.589}\right|^3.$$

This last prescription represents an attempt to deal with caustic formation in a simple manner, so that equation (23) gives sensible results for any value of $`\sigma `$. However, it is not relevant for the results given here.

The PDF for $`\delta `$, $`P(\delta ,\sigma _0)`$, may be obtained in the very process of computing any of the moments. To this end, we add up all the values of $`G(\stackrel{}{\lambda },\alpha )`$ (multiplied by the corresponding volume of the integration step) that are found during the integration to lie between $`\delta `$ and $`\delta +\mathrm{\Delta }\delta `$ (for a suitably chosen $`\mathrm{\Delta }\delta `$, depending on $`\delta `$) and divide the result by $`\mathrm{\Delta }\delta `$.

The values of $`\delta `$ given by (17) for very under dense regions is not correct this imply $`\delta 0`$. However, this problem disappears if we set $`P(\delta )=0`$ for $`\delta <0.82`$.

The values of $`S_3`$ and $`S_4`$ given by

$$S_3\frac{\delta ^3}{\delta ^2^2};S_4\frac{\delta ^43\delta ^2^2}{\delta ^2^3}$$

in the low-$`\sigma `$ limit, may be computed analytically. The exact value of $`S_3(\sigma =0)=(34/7(n+3))`$ is obtained, for the power spectra under consideration, when expression (17) for $`\delta `$ is used; the value $`4(n+3)`$ is found when the $`r_i`$ are set equal to one. For a non-filtered field we have computed $`S_4(\sigma =0)`$ analytically both for the CZA (equation (17)) and for the Zel’dovich approximation. Note that the non-filtered case may be treated formally as a filtered field with negligible power below the filtering scale, so that we may use expression (30) with $`\sigma `$ independent of $`\lambda `$ and $`\gamma =0`$. For power law power spectra this corresponds to the limit $`n=3`$, but it corresponds more generally to any power spectra having the mentioned properties. The correct values (60712/1323; 279/9; Bernardeau 1994a) are found. For the evolution of these amplitudes and the variance in the small $`\sigma `$ limit we find

$$\delta ^2=\sigma ^2+1.81\sigma ^4+3.8\sigma ^6$$

$$S_3=4.86+9.815\sigma ^2+48\sigma ^4$$

$$S_4=45.89+288.8\sigma ^2+2072\sigma ^4$$
(21)
This results are exact to next leading order for $`S_3`$ (in principle only for $`n=2.7`$ but the $`n`$ dependence is negligible) and almost exact to the same order for $`\delta ^2`$, $`S_4`$. As pointed out in section 2, we could use the statistical properties of the $`x`$, $`u`$ variables to obtain these two quantities exactly to the stated order, but the second coefficient in the expansion for $`\delta ^2`$ will increase by less than a per cent which is about the numerical error for this coefficient. So, it would be an unnecessary complication to use here those additional variables. Furthermore, the purpose of the present procedure is to obtain directly and accurately the PDF for finite values of $`\sigma `$. The question of the maximum order of exactness, although interesting, is only of peripheral relevance to our purposes. These results agree with those exact results given by the perturbation theory by Scoccimarro & Frieman (1995) and the remaining ones are close to those obtained by Fosalba and Gaztañaga (1998) using the spherical model. We see that this model predict well the hierarchical amplitudes but not so well the moments. The above results provide a check to the CZA but not properly to the procedure presented here, since for non-filtered fields expressions (12) and (14) reduce to the standard form. A check to these expressions is provided by the cases of power law power spectra with $`n3`$. For $`n=2`$ we find:

$$\delta ^2=\sigma ^2+0.88\sigma ^4+0.54\sigma ^6$$

$$S_3=3.86+3.12\sigma ^2$$

$$S_4=27.56+58\sigma ^2$$
(22)
Here, to obtain the second coefficient in $`\delta ^2`$ exactly we have considered the small (=0.015) contribution due to the dispersion of the $`x`$ variables. This is not necessary for the amplitudes. It is also important to consider the shape dependence of $`\sigma `$ (see expression (6)), otherwise we would obtain a coefficient 0.1 too large. For $`n=1`$, $`n=1`$ we find respectively:

$$\delta ^2=\sigma ^2+0.223\sigma ^4$$

$$S_3=2.86+0.267\sigma ^2$$

$$S_4=13.892.583\sigma ^2$$

$$\delta ^2=\sigma ^20.173\sigma ^4$$

$$S_3=0.86+0.8\sigma ^2$$

$$S_4=0.55617.8\sigma ^2$$
(23)
For $`n1.5`$, $`S_3`$ and $`S_4`$ are practically independent of $`\sigma `$. For example, for $`n=1.3`$ the change of these quantities up to $`\sigma =1.52`$ is less than a 4%. For smaller values of $`n`$, however, the evolution of this quantities is more apparent. For $`n=2`$ we find (table 1) that $`S_3`$ and $`S_4`$ increase considerably from $`\sigma =0`$ to $`\sigma 0.5`$ and more steeply for higher values of $`\sigma `$. In figure 1 we present the PDF for $`\delta `$ corresponding to $`P(k)k^1`$ and $`\sigma =0.3;0.5`$. Comparing this with figure 7 of Bernardeau & Koffman (1995), we find them to be almost identical. The PDFs of Bernadeau & Koffman were obtained under the assumption that the $`S_p`$ parameters (the $`p`$-th cumulants divided by $`\sigma ^{2(p1)}`$) are independent of $`\sigma `$. We have seen that this is the case for $`n=1`$, at least for $`S_3`$ and $`S_4`$. So, the agreement between both pairs of PDFs is to be expected. The evolution of $`\delta ^2`$ agrees very well with the result by Lokas et al. (1996).

In figure 2 the PDFs corresponding to $`\sigma =0.92;1.52`$, $`n=1.3`$ are represented. These should be compared with figure 10 of Bernardeau & Kofman (1995). The agreement is very good, our results being always within the numerical errors of the results given in that article.

It should be noted that the mentioned numerical results assume a cold dark matter (CDM) while we have used a power-law power spectrum with the same local logarithmic derivative. Using the whole CDM power spectrum in our computations would represent only a small complication. If we have not done so, it is because this implies only a slight improvement in the results that is completely unobservable given the precision of the simulations.

The procedure described here for obtaining the PDF for $`\delta `$ is in principle valid only for top-hat filtered fields. However, the basic lines of reasoning followed here are entirely translatable to another smoothing procedures.

For a given smoothing (in the present field), we must determine the smoothing in the initial field rendering the minimum dispersion for the smoothed present field. That is, we must find a smoothing procedure that, when applied to the initial field at a given point, renders a set of $`\lambda _i^{}`$ (extrapolated linearly to the present time) that predict the mean value of $`\delta `$ at the transformed point, given by an expression similar to (17), and so that the dispersion in $`\delta `$ for given $`\lambda _i^{}`$s is a minimum. This is what we have implicitly done here for the spherical top-hat smoothing, and we have found that it corresponds to an initial ellipsoidal top-hat smoothing.

The solution to this problem for Gaussian smoothing is not simple, although one could guess that a triaxial Gaussian, with the axes determined as in the top-hat case, should not be a bad approximation. Within this simple approximation and for a power-law spectrum the PDF for $`\delta `$ must be the same as for the top-hat case with the same $`\sigma `$. For more general power spectra, we cannot decide a priori whether the value of the logarithmic derivative must be taken on the scale where the top-hat has the same volume as the Gaussian or on the scale where it has the same $`\sigma `$.

As an example we discuss one of the cases represented in figure 4 of Kofman et al. (1994), where the Gaussian smoothing radius is $`R_s=6h^1`$Mpc; $`\sigma =0.55`$. If we matched the volumes, the corresponding top-hat radius is 9.38 $`h^1`$Mpc, so that for a CDM power spectrum $`n1`$. The PDF corresponding to a top-hat filtering with $`\sigma =0.55`$; $`n=1`$ is represented in figure 3 and shows very good agreement with the numerical results in the last reference. If we matched the $`\sigma `$’s, the corresponding top-hat smoothing scale imply $`n=1.1`$. The corresponding PDF is almost indistinguishable from the former.

We have seen that to obtain the PDF for $`\delta `$ we must integrate over all values of $`\lambda `$ subject to the condition $`\delta (\lambda )=\delta `$. If we considered, however, an approximation where $`\delta `$ depended on the $`\lambda _i`$ only through $`\delta _L=\lambda _1+\lambda _2+\lambda _3`$, the computation of $`P(\delta )`$ would be greatly simplified. We have already seen that such an approximation cannot be entirely correct, but leads to a simple analytical expression that displays, in a clear manner, the effect of smoothing. Taking for $`\delta (\delta _L)`$ that corresponding to spherical collapse and using the fact that the distribution for $`\delta _L`$ is simply a Gaussian, equation (11) and (14) immediately gives for $`P(\delta )`$:

$$P(\delta )=\frac{1}{\sqrt{2\pi }}\frac{e^{\frac{\delta _L^2(\delta )}{2\sigma ^2(\delta )}}}{\sigma (\delta )(1+\delta )}\left((1+\delta )^{0.602}+\frac{\alpha }{2}\delta _L(\delta )\right)$$

$$\times \left(1(1\frac{\alpha }{2})\frac{\delta _L(\delta )}{1.66}\right)^{1.66};$$

$$\sigma (\delta )=\sigma _0(1+\delta )^{\alpha /2};\alpha \frac{n+3}{3}$$

$$\delta _L(\delta )=1.66(1(1+\delta )^{0.6024}).$$
(24)
Where the last equation represents a very accurate approximation to the spherical model that we propose. However, expression (24) represents only a rough approximation to $`P(\delta )`$. For $`\sigma =1.52`$ the value of $`P(\delta )`$ given by (24) may differ from the real value by a factor greater than three. Furthermore, the value of $`S_3(\sigma =0)`$ given by (24) is $`37/7(n+3)`$, which is 3/7 too large.

We may conclude that for $`\delta (0.9,10)`$ and $`\sigma `$ as large as 1.52, the $`P(\delta )`$ computed as indicated here is, so far as we have been able to check, in excellent agreement with the results of numerical simulations. This agreement is particularly meaningful since our derivation of $`P(\delta )`$ involves neither free parameters nor modeling freedom (except for some minute refinements non-essential for the results given here). Equation (12) is exact and is valid for any field. Equation (14) for the Jacobian of the Eulerian–Lagrangian transformation in a filtered field is also exact (at least to first order in $`\lambda _i`$) and, although derived for a Gaussian field, it seems to be valid for any field.

FIGURE CAPTIONS

Fig. 1: Probability distribution function (PDF) for $`\sigma =0.3`$ and $`\sigma =0.5`$; $`n=1`$.

Fig. 2: Probability distribution function (PDF) for $`\sigma =0.92`$ and $`\sigma =1.52`$; $`n=1.3`$.

Fig. 3: Probability distribution function (PDF) for $`\sigma =0.55`$; $`n=1`$.


format: text

text:
# On black hole masses, radio-loudness and bulge luminosities of Seyfert galaxies

## 1 Introduction

Supermassive black holes (SMBHs), with masses in the range of $`10^6`$ to $`10^9\mathrm{M}_{}`$, have been suggested to exist in the center of quasars and active galactic nuclei (AGNs). Accretion onto these SMBHs may account for the huge power of these energetic objects (Lynden-Bell 1969; Rees 1984). Recently, the masses of central objects in 20 Seyfert galaxies and 17 nearby quasars have been measured with the reverberation mapping technique (Ho 1999; Wandel, Peterson & Malkan 1999; Kaspi et al. 2000), which confirmed the existence of SMBHs in the center of these objects. On the other hand, a lot of observations using gas and stellar dynamics indicated that SMBHs probably also exist in the center of our Galaxy (Ghez et al. 1998; Genzel et al. 1997) and in the nuclei of many normal galaxies (Kormendy & Richstone 1995; Magorrian et al. 1998; Kormendy & Gebhardt 2001). However, the apparent inactive feature of normal galaxies seems to suggest that the central engines of these galaxies may be different from those in quasars and AGNs. It was suspected that advection-dominated accretion flows (ADAFs; see Narayan, Mahadevan & Quataert 1998 and Kato, Fukue & Mineshige 1998 for reviews) with very low accretion rate and low radiative efficiency may exist in the nuclei of normal galaxies (Fabian & Rees 1995; Di Matteo & Fabian 1996), while quasars and most AGNs probably host the standard geometrically thin accretion disks (Shakura & Sunyaev 1973; Novikov & Throne 1973) with higher accretion rate.

One interesting result found recently in the searches of SMBHs in nearby galaxies is the correlation of black hole masses with the properties of galactic bulges. Although with large scatters, the black hole masses seem to correlate with bulge luminosities (Kormendy 1993; Kormendy & Richstone 1995; Magorrian et al. 1998). This also leads to the finding that the black hole mass, $`\mathrm{M}_{\mathrm{BH}}`$, is possibly proportional to the bulge mass, $`\mathrm{M}_{\mathrm{bulge}}`$, though the mass ratio found by different authors was different, in the range of $`0.2\%`$ to $`0.6\%`$ (Kormendy & Richstone 1995; Magorrian et al. 1998; Ho 1999). Laor (1998; 2001) recently found that some nearby quasars and Seyfert galaxies follow nearly the same $`\mathrm{M}_{\mathrm{BH}}`$-bulge luminosity relation as normal galaxies, and suggested a universal nonlinear relation, $`\mathrm{M}_{\mathrm{BH}}\mathrm{M}_{\mathrm{bulge}}^{1.54\pm 0.15}`$, for both normal and active galaxies. This means that the $`\mathrm{M}_{\mathrm{BH}}/\mathrm{M}_{\mathrm{bulge}}`$ ratio is not constant for galaxies with different bulge luminosities. Recently, a significantly tight correlation of $`\mathrm{M}_{\mathrm{BH}}`$ with the bulge velocity dispersion $`\sigma `$ was also found for nearby galaxies (Gebhardt et al. 2000a; Ferrarese & Merritt 2000). More recent studies indicated that 11 Seyfert galaxies with $`\mathrm{M}_{\mathrm{BH}}`$ measured by reverberation mapping follow the same $`\mathrm{M}_{\mathrm{BH}}`$-$`\sigma `$ relation as for normal galaxies (Gebhardt et al. 2000b; Ferrarese et al. 2001), implying another possible universal relation for both normal and active galaxies. The correlations of the SMBH mass with the properties of the galactic bulge strongly suggest a tight connection between the formation and evolution of the SMBH and galactic bulge, though the nature of this connection is still in debate (Haehnelt & Rees 1993; Haiman & Loeb 1998; Silk & Rees 1998; Haehnelt, Natarajan & Rees 1998; Kauffmann & Haehnelt 2000; Adams, Graff & Richstone 2001)

Nonthermal radio emissions of quasars and AGNs are believed to be probably produced by relativistic electrons that are powered by jets (Begelman, Blandford & Rees 1984; Blundel & Beasley 1998). Similar radio emissions have been detected in the nuclei of normal elliptical galaxies (Sadler, Jenkins & Kotanyi 1989) and also in some spiral galaxies (Sadler et al. 1995). Recently, some studies have indicated that the radio power may be directly correlated with the black hole mass. Franceschini, Vercellone & Fabian (1998) found a very tight relation between black hole mass and radio power in a small sample of nearby mostly non-active galaxies. McLure et al. (1999) estimated the black hole masses for a sample of AGNs using $`\mathrm{M}_{\mathrm{BH}}`$-bulge mass relation found by Magorrian et al. (1998) and noted they follow the same correlation with radio power as found by Franceschini et al. (1998). However, Laor (2000) recently argued that the $`\mathrm{M}_{\mathrm{BH}}`$ and radio power of a sample of 87 $`z<0.5`$ Palomar-Green (PG) quasars (Schmidt & Green 1983; Boroson & Green 1992) do not follow the tight correlation suggested by Franceschini et al. (1998), because quasars usually have over 100 times larger radio power than normal galaxies at a given $`\mathrm{M}_{\mathrm{BH}}`$. He suggested that the larger scatters of radio power at a given $`\mathrm{M}_{\mathrm{BH}}`$ may be simply due to the different levels of overall continuum luminosity of different objects. Laor (2000) also noted that the radio-loud quasars seem to host more massive black holes than radio-quiet quasars. Very recently, Ho & Peng (2001) studied the radio-loudness of bright Seyfert 1 galaxies using the nuclear radio and optical luminosities, and suggested that the majority of Seyfert 1 nuclei in their sample are essentially radio loud. Therefore, it is useful and feasible to check if these radio-loud Seyfert nuclei host more massive black holes than radio-quiet ones and if the Seyfert galaxies still follow the same correlation between the radio power and black hole masses of nearby galaxies and quasars (Franceschini et al. 1998; Laor 2000). In this paper we will try to derive the SMBH masses for a number of Seyfert galaxies and study their correlations with radio power and bulge luminosities. Further studies of these correlations will probably provide some important clues to help us understand the nature of radio emissions in active and normal galaxies and the formation and evolution of SMBHs and galaxies.

## 2 Estimating the SMBH masses of Seyfert galaxies

Currently the SMBH masses of a few weak AGNs have been well measured by stellar dynamics, ionized gas dynamics and water maser dynamics (for a summary see Table 1s in Ho 1999 and Gebhardt et al. 2000a). Using the reverberation mapping technique, the SMBH masses of 20 Seyfert galaxies and 17 bright quasars have been recently estimated (Ho 1999; Wandel et al. 1999; Kaspi et al. 2000). However, for most Seyfert galaxies, it is difficult to measure the SMBH mass using these methods because of either the large nuclear luminosity or the lack of long-term variability monitoring and precise measurements of characteristic velocity dispersions in the broad emission line region.

Since the tight correlation between the SMBH mass and bulge velocity dispersion was found for both normal and active galaxies (Gebhardt et al. 2000b; Ferrarese et al. 2001), it may be straightforward to estimate the SMBH mass from the measured bulge velocity dispersion. The correlation between the black hole mass and nuclear velocity dispersion for Seyfert galaxies with both measured nuclear velocity dispersions and black hole masses is demonstrated in Figure 1. The SMBH masses of 11 Seyfert galaxies (squares in Figure 1) were estimated by reverberation mapping (Wandel et al. 1999) and 9 (triangles in Figure 1) by the dynamical method (Ho 1999; Gebhardt et al. 2000a). The values of their central velocity dispersions were taken from Nelson & Whittle (1995) and Ferrarese et al. (2001). Figure 1 clearly shows that Seyfert galaxies follow the same $`\mathrm{M}_{\mathrm{BH}}`$-$`\sigma `$ relation as normal galaxies (Gebhardt et al. 2000a; Merritt & Ferrarese 2000). Although measurements of the central velocity dispersions of these Seyfert galaxies were made by different groups, the small scattering around the $`\mathrm{M}_{\mathrm{BH}}`$-$`\sigma `$ relation indicates that these measurements were reliable and the systematic errors may not be important (Ferrarese et al. 2001). In this paper, we adopt the $`\mathrm{M}_{\mathrm{BH}}`$-$`\sigma `$ relation found by Merritt & Ferrarese (2001), namely,

$$\mathrm{M}_{\mathrm{BH}}=1.3\times 10^8\mathrm{M}_{}(\sigma /200\mathrm{km}\mathrm{s}^1)^{4.72},$$
(1)
to derive the SMBH masses for Seyfert galaxies with measured nuclear velocity dispersions. Using a slightly flatter relation found by Gebhardt et al. (2000) does not cause significant changes in our results.

In the next sections we will compare some properties of Seyfert galaxies with those of quasars. Only a few quasars have SMBH masses determined by reverberation mapping. Laor (1998) adopted an empirical relation between the size of the broad line region (BLR) and the bolometric luminosity, $`R_{BLR}L^{1/2}`$ (Kaspi et al. 1996; see also Kaspi et al. 2000 for a slightly steeper relation), and derived the SMBH mass for a number of quasars using the measured $`H\beta `$ velocity dispersion and the continuum luminosity. We have compared the estimated SMBH masses for 9 Seyfert galaxies using Laor’s method (but adopted the BLR velocity $`V=(\sqrt{3}/2)\mathrm{FWHM}(\mathrm{H}\beta )`$, as did Kaspi et al. (2000)) with the those obtained by reverberation mapping, and found that they agree well. This shows that using SMBH masses derived with different methods may not cause serious problems in our present work.

## 3 Sample of Seyfert galaxies

We selected 37 Seyfert galaxies from two well-studied nearby Seyfert samples, the Palomar optical spectroscopic survey of bright ($`B_T12.5\mathrm{mag}`$), northern ($`\delta >0^o`$) galaxies (Ho, Filippenko & Sargent 1995), including 21 Seyfert 1s and 28 Seyfert 2s (Ho et al. 1997a), the most complete and least biased available (Ho, Filippenko & Sargent 1997b; Ho & Ulvestad 2001), and the CfA redshift Survey (Huchra et al. 1983) of galaxies with Zwicky magnitude $`14.5`$, including 33 Seyfert 1s and 15 Seyfert 2s (Huchra & Burg 1992; Osterbrock & Martel 1993). The Seyfert samples from these two surveys seem to complement one another, though the combined sample may not be complete (Ho & Peng 2001).

Our Seyfert sample include 22 Seyfert 1s and 15 Seyfert 2s in Palomar and CfA samples (see Table 1). These Seyfert galaxies have either the measured SMBH masses or the measured central velocity dispersions. Among them, three Seyfert 1s and two Seyfert 2s have dynamical SMBH masses (Ho 1999; Gebhardt et al. 2000a), and ten Seyfert 1s have SMBH masses measured by the reverberation mapping method (Wandel et al. 1999; Ho 1999). Another Nine Seyfert 1s and 13 Seyfert 2s have measured central velocity dispersions but unknown SMBH masses (Nelson & Whittle 1995). The SMBH masses can be estimated by equation (1) using the measured central velocity dispersions. Therefore, our sample consists of 37 Seyfert galaxies with derived SMBH masses (see Table 1).

## 4 Radio properties and black hole masses of Seyfert galaxies

The radio properties of Seyfert galaxies have been investigated in detail using the Very Large Array (VLA) at 3.6 cm by Kukula et al. (1995) for the CfA sample, and at 6 cm and 20 cm recently by Ho & Ulvestad (2001) for the Palomar sample. In this paper, we adopt the radio data for Seyfert 1s from Ho & Peng (2001) and the data for Seyfert 2s from Ho & Ulvestad (2001) and Kukula et al. (1995). The 6 cm data for some sources in the CfA sample (Kukula et al. 1995) were extrapolated from the 3.6cm data assuming $`f_\nu \nu ^{0.5}`$.

Seyfert galaxies have been considered usually as radio-quiet AGNs because most of them have lower radio-loudness, defined as the ratio of the radio to optical luminosities, $`R=L_{6cm}/L_B`$. However, a recently study of Ho & Peng (2001) showed that though the nuclear radio power for Seyfert 1 galaxies on average accounts for about $`75\%`$ of the total radio emission, the nuclear optical luminosity, measured by high-resolution optical images, accounts for merely 0.01% of the integrated light. If the radio-loudness is measured by the nuclear radio and optical luminosities, most Seyfert 1s are in the category of radio-loud AGNs (with radio-loudness larger than 10). In Table 1, we listed the total and nuclear absolute B magnitudes ($`M_B^{\mathrm{tot}}`$ and $`M_B^{\mathrm{nuc}}`$) and nuclear radio power ($`P_{6cm}`$) for Seyfert galaxies. The nuclear radio-loudness, $`R_c`$, is calculated by $`R_c=L_{6cm}^{\mathrm{nuc}}/L_B^{\mathrm{nuc}}`$. The Hubble constant $`H_0=75\mathrm{k}\mathrm{m}\mathrm{s}^1\mathrm{Mpc}^1`$ and the deceleration parameter of $`q_0=0.5`$ were adopted. For Seyfert 2s, the measurements of their nuclear optical magnitudes are not available in the literature.

Figure 2 shows the relation of SMBH masses with the nuclear radio-loudness of Seyfert 1s and the radio power of Seyfert galaxies listed in Table 1. It is clear that the radio-loud Seyfert 1s with larger nuclear radio-loudness seem to host more massive black holes. The similar tendency has been found for PG quasars by Laor (2000). From Figure 2(b) we see a trend that Seyfert galaxies having a larger radio power perhaps host a more massive black hole than those of less radio power. With the same black hole mass, Seyfert galaxies seem to have 100 to 1000 times greater radio power than normal galaxies. Laor (2000) has shown that PG quasars also depart from such a correlation for nearby galaxies, with the radio luminosity of quasars being $`10^4`$ larger at a given $`\mathrm{M}_{\mathrm{BH}}`$. The difference of radio luminosities of quasars, Seyfert galaxies and nearby galaxies may be simply due to different levels of nuclear activity.

## 5 Black hole masses and bulge luminosities

The correlation of black hole mass and galactic bulge luminosity found for nearby galaxies implies a relationship between the black hole and bulge masses (Kormendy & Richstone 1995; Magorrian et al. 1998; Ho 1999). However, it is not clear whether the ratio between black hole and bulge masses remains constant for all kinds of objects. Laor (2001) recently checked the correlation between the estimated black masses and bulge luminosities for 15 nearby quasars and 9 Seyfert galaxies, and found that they probably follow the same $`\mathrm{M}_{\mathrm{BH}}`$-bulge luminosity relation as for 16 nearby galaxies. This suggests a universal nonlinear relation between the estimated black masses and bulge luminosities, $`\mathrm{M}_{\mathrm{BH}}\mathrm{M}_{\mathrm{bulge}}^{1.54\pm 0.15}`$, for both normal and active galaxies. Therefore, $`\mathrm{M}_{\mathrm{BH}}/\mathrm{M}_{\mathrm{bulge}}`$ ratio is not constant for galaxies with different bulge luminosities. However, this needs to be confirmed by larger samples including both normal and active galaxies.

### 5.1 Our sample of Seyfert galaxies

We use the 37 Seyfert galaxies with derived SMBH masses in our sample to check the $`\mathrm{M}_{\mathrm{BH}}/\mathrm{M}_{\mathrm{bulge}}`$ ratio. The absolute total B magnitudes ($`\mathrm{M}_B^{\mathrm{tot}}`$) for our Seyfert sample (see Table 1) were taken from Ho et al. (1997b) and Whittle (1992) (adapted to $`H_0=75\mathrm{k}\mathrm{m}\mathrm{s}^1\mathrm{Mpc}^1`$). The absolute bulge B magnitudes ($`\mathrm{M}_B^{\mathrm{bulge}}`$) were obtained based on the relation between $`\mathrm{M}_B^{\mathrm{bulge}}`$ and $`\mathrm{M}_B^{\mathrm{tot}}`$ (Simien & de Vaucouleurs 1986) and the Hubble stages (defined in de Vaucouleurs, de Vaucouleurs & Corwin 1976) of the host galaxies. The B-V colors were taken from Verón-Cettty & Verón (2000). In order to compare our result with those obtained for quasars and nearby galaxies, we convert $`M_B^{\mathrm{bulge}}`$ of Seyfert galaxies to the absolute bulge V magnitude, $`M_V^{\mathrm{bulge}}`$ (adapted to $`H_0=80\mathrm{k}\mathrm{m}\mathrm{s}^1\mathrm{Mpc}^1`$) by assuming that B-V color in the bulge is nearly the same as the total B-V color of Seyfert galaxy. In Figure 3 we show the relation between the SMBH masses and bulge V-band absolute magnitudes for Seyfert galaxies and the comparison of this relation with those for quasars and nearby galaxies. The SMBH masses of 15 quasars and the absolute V magnitudes of their inner host were taken from Laor (1998, 2001).111Because we used BLR velocity dispersion as $`V=(\sqrt{3}/2)\mathrm{FWHM}(\mathrm{H}\beta )`$, the SMBH masses of 15 quasars are smaller by a factor of 0.75 than those given by Laor (1998).  The SMBH masses determined by stellar and gas dynamics for 33 nearby galaxies were taken from Kormendy & Gebhardt (2001) (we omitted M 81, a Seyfert 1 already in our Seyfert sample, and 3 objects with black hole masses measured by maser dynamics). Their absolute bulge V magnitudes were derived from B-band absolute bulge magnitudes by adopting a standard bulge color, B-V = 0.94 (Worthey 1994). It is evident that both Seyfert and nearby galaxies, as well as quasars, seem to follow the same $`\mathrm{M}_V^{\mathrm{bulge}}`$-$`\mathrm{M}_{\mathrm{BH}}`$ relation. Note that unlike others, nearby galaxies NGC 4342, NGC 4486B and M 32 deviate significantly from the $`\mathrm{M}_V^{\mathrm{bulge}}`$-$`\mathrm{M}_{\mathrm{BH}}`$ relation. These offset galaxies are usually fainter because the outer regions of them may have been stripped away in the tidal interactions with more massive companions (Faber 1973). These three galaxies, together with two Seyfert 1 galaxies, NGC 1275 with peculiar Hubble type and NGC 4395 with only the upper limit of measured SMBH mass, are not included in our statistical studies. The measurement errors of $`\mathrm{M}_{\mathrm{BH}}/\mathrm{M}_{}`$ of our Seyfert sample were listed in Table 1, and those of normal galaxies were adopted from Kormendy & Gebhardt (2001). For nearby quasars, such errors are not available in Laor (1998) and assumed to be 60%. The measurement errors of $`\mathrm{M}_V^{\mathrm{bulge}}`$ were seldom mentioned in the literature and thus are difficult to estimate. These errors for Seyfert galaxies were allocated to be 0.5, 0.75 and 1.0 mag respectively based on their quality assessment factors (Whittle 1992). Those of normal galaxies and nearby quasars were adopted to be 0.5 and 0.75 mag, respectively. Taking into account these ‘typical’ errors of $`\mathrm{M}_V^{\mathrm{bulge}}`$ and $`\mathrm{M}_{\mathrm{BH}}`$, the least $`\chi ^2`$ fit for 35 Seyfert galaxies, 15 quasars and 30 normal galaxies gives:

$$M_V^{bulge}=11.01\pm 0.78(1.22\pm 0.10)\mathrm{log}(\mathrm{M}_{\mathrm{BH}}/\mathrm{M}_{}).$$
(2)
The Spearman rank-order correlation coefficient is $`r_S=0.76`$, which has a probability of $`P_r=3.3\times 10^{16}`$ occuring by chance. Considering the errors in both parameters, we adopted a bootstrap method to estimate the uncertainty of this correlation, and obtained the mean correlation coefficient $`r_S=0.63\pm 0.05`$. Using the standard relation

$$M_V^{bulge}=4.832.5\mathrm{log}L_{bulge}/L_{},$$
(3)
and the relation between the bulge mass and luminosity (Magorrian et al. 1998),

$$\mathrm{log}(\mathrm{M}_{\mathrm{bulge}}/\mathrm{M}_{})=1.11+1.18\mathrm{log}(L_{bulge}/L_{}),$$
(4)
we then can get from eq. (2),

$$\mathrm{log}(\mathrm{M}_{\mathrm{BH}}/\mathrm{M}_{\mathrm{bulge}})=11.06\pm 1.11+(0.74\pm 0.14)$$

$$\mathrm{log}(\mathrm{M}_{\mathrm{bulge}}/\mathrm{M}_{}).$$
(5)
This gives $`\mathrm{M}_{\mathrm{BH}}\mathrm{M}_{\mathrm{bulge}}^{1.74\pm 0.14}`$. The fitting for a sample of 35 Seyfert galaxies and 30 normal galaxies alone gives almost identical results. Laor (2001) recently found $`\mathrm{M}_{\mathrm{BH}}\mathrm{M}_{\mathrm{bulge}}^{1.54\pm 0.15}`$ for a sample of objects including 9 Seyfert galaxies and 15 quasars. Our result shows that this nonlinear relation is more evident from a larger sample including more Seyfert galaxies. Eq. (5) also shows that $`\mathrm{M}_{\mathrm{BH}}/\mathrm{M}_{\mathrm{bulge}}`$ is about 0.02 when $`\mathrm{M}_{\mathrm{BH}}=10^6\mathrm{M}_{}`$, and is 0.3% when $`\mathrm{M}_{\mathrm{BH}}=10^9\mathrm{M}_{}`$. This clearly indicates that galaxies with more massive black holes have a larger $`\mathrm{M}_{\mathrm{BH}}/\mathrm{M}_{\mathrm{bulge}}`$ ratio.

### 5.2 Seyfert sample in reverberation mapping studies

We noted that the above universal $`\mathrm{M}_\mathrm{V}^{\mathrm{bulge}}`$-$`\mathrm{M}_{\mathrm{BH}}`$ relation (eq. 2) is contrary to the result obtained by Wandel (1999), who found that Seyfert 1 galaxies significantly depart from the $`\mathrm{L}_{\mathrm{bulge}}`$-$`\mathrm{M}_{\mathrm{BH}}`$ relations of quasars and nearby galaxies ( Laor 1998; Magorrian et al. 1998). However, Mclure & Dunlop (2001) used the decomposed bulge luminosities of the same sample of Seyfert galaxies as in Wandel (1999) and found no evidence for a different $`\mathrm{L}_{\mathrm{bulge}}`$-$`\mathrm{M}_{\mathrm{BH}}`$ relation from quasars. We noted that Wandel (1999) took the absolute bulge B magnitude, $`\mathrm{M}_\mathrm{B}^{\mathrm{bulge}}`$, for Seyfert galaxies from Whittle (1992)(adopting $`H_0=50\mathrm{k}\mathrm{m}\mathrm{s}^1\mathrm{Mpc}^1`$) and calculated the bulge luminosity from $`M_B^{bulge}`$ with the standard expression for $`\mathrm{M}_\mathrm{V}^{\mathrm{bulge}}`$ (see eq. (3)) by assuming $`\mathrm{M}_\mathrm{B}^{\mathrm{bulge}}\mathrm{M}_\mathrm{V}^{\mathrm{bulge}}`$ for Seyfert 1s. We carefully checked these points and re-derived the absolute bulge V magnitude of 17 Seyfert 1 galaxies with SMBH masses measured by reverberation mapping (Wandel et al. 1999), using $`H_0=80\mathrm{k}\mathrm{m}\mathrm{s}^1\mathrm{Mpc}^1`$ and assuming the bulge B-V color being nearly the same as the B-V color of the whole galaxy. The $`\mathrm{M}_\mathrm{V}^{\mathrm{bulge}}`$-$`\mathrm{M}_{\mathrm{BH}}`$ relation of 17 Seyfert 1 galaxies, as compared with those for 15 quasars (Laor 1998, 2001) and 30 nearby galaxies with $`\mathrm{M}_{\mathrm{BH}}`$ measured by stellar dynamics (Gebhardt et al. 2000a), is shown in Figure 4. Allocating the ‘typical’ uncertainties of $`\mathrm{M}_\mathrm{V}^{\mathrm{bulge}}`$ and $`\mathrm{M}_{\mathrm{BH}}`$ as we did above, the least $`\chi ^2`$ fit for 62 objects gives,

$$M_V^{bulge}=10.00\pm 0.96(1.34\pm 0.12)\mathrm{log}(\mathrm{M}_{\mathrm{BH}}/\mathrm{M}_{}).$$
(6)
Considering the errors in both parameters, we estimated the mean correlation coefficient as $`r_S=0.63\pm 0.06`$. This result is consistent with that found above (see eq. (2)) and also indicates a nonlinear $`\mathrm{M}_{\mathrm{BH}}`$-$`\mathrm{M}_{\mathrm{bulge}}`$ relation.

### 5.3 Nelson & Whittle’s sample of Seyfert galaxies

In this section we explore the Seyfert sample in Nelson & Whittle (1995), where the measurements of nuclear velocity dispersions of about 70 Seyfert galaxies were reported. After excluding several LINERs and normal galaxies, we got 33 Seyfert 1s and 32 Seyfert 2s. We estimated the SMBH masses of these Seyfert galaxies using the $`\mathrm{M}_{\mathrm{BH}}`$-$`\sigma `$ relation (eq. (1)). The total radio power at 5 GHz of them was calculated from the 1.4 GHz data in Nelson & Whittle (1995) by assuming $`f_\nu \nu ^{0.5}`$ and $`H_0=80\mathrm{k}\mathrm{m}\mathrm{s}^1\mathrm{Mpc}^1`$. The radio-loudness was calculated by the 5 GHz radio luminosity and the B-band optical luminosity. The V-band bulge absolute magnitude was estimated from the B-band bulge absolute magnitude (adapted to $`H_0=80\mathrm{k}\mathrm{m}\mathrm{s}^1\mathrm{Mpc}^1`$) in Nelson & Whittle (1995) by assuming the B-V color of the bulge being the same as the total B-V color (taken from Verón-Cetty & Verón (2000)) of Seyfert galaxy.

The radio-loudness and total radio power against the SMBH masses for 29 Seyfert 1s and 25 Seyfert 2s with available radio data in the sample of Nelson & Whittle (1995) are plotted in Fig. 5. Now we see that most Seyfert galaxies seem to be radio quiet ($`R<10`$) when we adopted the total radio and optical luminosities to calculate the radio-loudness. There is a weak tendency that Seyfert galaxies with larger radio-loudness have larger SMBH masses. A comparison with the result in Figure 2(a) indicates that the nuclear radio-loudness may be more fundamental and reflect the nature of central engine of Seyfert galaxies. From Figure 5(b) we see that there is a strong correlation between the total radio power and the SMBH mass, though the scatters are large. This confirms again the previous result that AGNs with larger radio power may host more massive SMBHs (Franceschini et al. 1998; McLure et al. 1999). However, as is shown in Figure 2, Seyfert galaxies depart significantly from the tight relation between the radio power and SMBH mass found for normal galaxies by Franceschini et al. (1998).

The relation between V-band bulge luminosities and SMBH masses for 33 Seyfert 1s and 32 Seyfert 2s in the sample of Nelson & Whittle (1995) is shown in Figure 6. It again indicates a significant universal relation between V-band bulge luminosities and SMBH masses for both Seyfert galaxies and quasars, as well as nearby galaxies. Adopting the ‘typical’ uncertainties of $`\mathrm{M}_V^{\mathrm{bulge}}`$ and $`\mathrm{M}_{\mathrm{BH}}`$ as we did above, the least $`\chi ^2`$ fit for 65 Seyfert galaxies, 15 quasars and 30 nearby galaxies gives,

$$M_V^{bulge}=11.33\pm 0.66(1.19\pm 0.09)\mathrm{log}(\mathrm{M}_{\mathrm{BH}}/\mathrm{M}_{}),$$
(7)
The simple Spearman rank-order correlation coefficient is $`r_S=0.76`$, which has a probability of $`P_r=1.1\times 10^{21}`$ occuring by chance. Considering the errors in both parameters, we estimated the mean correlation coefficient as $`r_S=0.61\pm 0.04`$. This result is also almost identical to what we found for the sample including 35 Seyfert galaxies. Therefore, the nonlinear $`\mathrm{M}_{\mathrm{BH}}`$-$`\mathrm{M}_{\mathrm{bulge}}`$ relation is confirmed by the substantially enlarged Seyfert sample.

### 5.4 Effects of larger systematic errors

A significant correlation between black hole mass and bulge luminosity was obtained above by considering the ‘typical’ errors of both parameters. These errors are mainly due to the measurement uncertainties of bulge magnitude and black hole mass. However, the systematic errors are probably substantially larger for both bulge magnitude and black hole mass. For example, if $`\mathrm{M}_{\mathrm{BH}}`$ is determined by reverberation mapping, the systematic errors caused by the unknown BLR geometry, inclination may be as large as a factor of 3 or more (Krolik 2001). Using the $`\mathrm{M}_{\mathrm{BH}}`$-$`\sigma `$ relation to estimate $`\mathrm{M}_{\mathrm{BH}}`$ also has larger systematic errors due to the uncertainty of the slope of this relation. In addition, in deriving bulge magnitude from galaxy magnitude for Seyfert galaxies we adopted a statistical relation given by Simien & de Vaucouleurs (1986). The systematic errors caused by applying this relation are probably quite large and are difficult to estimate quantitatively.

In order to estimate the effects of possible larger systematic errors of black hole mass and bulge magnitude on our result, we adopted a bootstrap method by adding the systematic errors into the uncertainties of both parameters. Assuming the possible systematic errors for bulge absolute magnitude being 1 mag and those for black hole mass being 90% for all 35 Seyfert galaxies, 15 quasars and 30 normal galaxies, we obtained the average correlation coefficient using the bootstrap approach as being $`r_S=0.40\pm 0.09`$, which indicates a moderately significant correlation. The minimum chi-square fit considering both measurement and systematic errors in both parameter gives, $`\mathrm{M}_V^{\mathrm{bulge}}=10.49\pm 1.77(1.29\pm 0.22)\mathrm{log}(\mathrm{M}_{\mathrm{BH}}/\mathrm{M}_{})`$. Comparing with the result we obtained in eq. (2), we found the slope of this relation does not change very much, though the uncertainty of the slope is doubled when the larger systematic errors are considered. The correlation coefficient also substantially decreases and its uncertainty increases accordingly. However, these changes are quite limited and have no significant effects on the result we have obtained.

## 6 Discussions

Our results obtained for a larger sample of Seyfert galaxies show that AGNs with a larger nuclear radio-loudness seem to have more massive black holes. This conclusion was obtained recently by Laor (2000) for nearby quasars and was supported by our study for Seyfert galaxies. Our results also strengthen the argument made by Ho & Peng (2001) that the majority of Seyfert 1s are essentially radio-loud AGNs. The total radio power of Seyfert galaxies increases with the black hole mass. At a given $`\mathrm{M}_{\mathrm{BH}}`$, quasars and Seyfert galaxies seem to have greater radio power than that of nearby galaxies. These difference may simply be due to the different level of nuclear activity for different kinds of galaxies. A recent study using the quasars from the FIRST Bright Quasar Survey found evidence for the dependence of radio-luminosity on accretion rate and SMBH mass (Lacy et al. 2001). This may help us to understand the origin of scatters in the relation between the radio power and SMBH mass. If we describe the nuclear activity of galaxies using the accretion rate $`\dot{M}`$, the difference of radio power at a given $`\mathrm{M}_{\mathrm{BH}}`$ may show that quasars and Seyfert galaxies have larger $`\dot{M}`$ than nearby galaxies. This seems also consistent with the picture that the accretion process in these systems may be different (Fabian & Rees 1995; Di Matteo & Fabian 1996). Most likely ADAFs with very low accretion rate exist in the nuclei of nearby galaxies, while quasars and most AGNs probably host standard geometrically thin accretion disks with higher accretion rate. However, the radio emissions from radio-quiet AGNs and nearby galaxies are not well understood at present. Whether they are from ADAFs (Narayan et al. 1998) or from weak jets (Falcke & Biermann 1996, 1999) still remains uncertain. For radio-loud AGNs, the radio emissions are thought to be mainly from the jet and are probably related to the magnetic fields or black hole spin (Blandford & Payne 1982; Blandford & Znajek 1977). If the radio emissions correlate with the black hole mass and accretion rate, we need to explain the possible relations of these parameters with magnetic fields and black hole spin. However, no satisfactory theory can provide clear physics about these relations at present.

Using the sample including 37 Seyfert galaxies in two well-defined bright Seyfert samples, we studied the correlations of black hole masses with V-band bulge absolute magnitudes and bulge masses for active and normal galaxies. We find that the correlation between the black hole masses and the V-band absolute bulge magnitudes for Seyfert galaxies is almost consistent with those found for quasars and nearby galaxies. The combined sample of 37 Seyfert galaxies, 15 quasars and 30 nearby galaxies seem to follow a universal nonlinear relation between the black hole and bulge masses, $`\mathrm{M}_{\mathrm{BH}}\mathrm{M}_{\mathrm{bulge}}^{1.74\pm 0.14}`$, which is slightly steeper than that found recently by Laor (2001) for a sample including 9 Seyfert galaxies. Our results support the suggestion that the ratio of $`\mathrm{M}_{\mathrm{BH}}`$ and $`\mathrm{M}_{\mathrm{bulge}}`$ is not constant and galaxies with larger bulge luminosities or more massive SMBH probably have a larger $`\mathrm{M}_{\mathrm{BH}}/\mathrm{M}_{\mathrm{bulge}}`$ ratio. Including 65 Seyfert galaxies in the larger sample of Nelson & Whittle (1995) led to almost the same result. In fact, the nonlinear relation between $`\mathrm{M}_{\mathrm{BH}}`$ and $`\mathrm{M}_{\mathrm{bulge}}`$ has been predicted by some theories and supported by some recent studies. For example, using a collapse model for black hole and bulge formation, Adams et al. (2001) predicted $`\mathrm{M}_{\mathrm{BH}}/\mathrm{M}_{\mathrm{bulge}}\sigma `$, which gives $`\mathrm{M}_{\mathrm{BH}}/\mathrm{M}_{\mathrm{bulge}}\mathrm{M}_{\mathrm{BH}}^{1/4}`$ if $`\mathrm{M}_{\mathrm{BH}}\sigma ^4`$. This states clearly that the ratio of black hole and bulge masses is not constant and can be larger for galaxies with more massive SMBHs. A similar result can be obtained by exploring the model of Wang, Biermann & Wandel (2000) who derived $`\mathrm{M}_{\mathrm{BH}}/\mathrm{M}_{\mathrm{bulge}}\sigma ^{1.4}`$. In addition, recent studies on narrow line Seyfert 1s showed that their mean $`\mathrm{M}_{\mathrm{BH}}/\mathrm{M}_{\mathrm{bulge}}`$ ratio is significantly smaller than that for normal Seyfert galaxies (Mathur, Kuraszkiewicz, & Czerny 2001), which is also consistent with our results because narrow-line Seyfert 1s probably have smaller $`\mathrm{M}_{\mathrm{BH}}`$ than normal Seyfert galaxies (Boller, Brandt & Fink 1996).

Finally we would like to mention that the tight correlation between the black hole mass and bulge luminosity is not surprising and is actually expected from some existing well-known relations. From the $`\mathrm{M}_{\mathrm{BH}}`$-$`\sigma `$ relation and Faber-Jackson relation between the bulge luminosity and $`\sigma `$ (Faber & Jackson 1976), we certainly expect a relation between black hole mass and bulge luminosity. The detailed studies of this relation and the $`\mathrm{M}_{\mathrm{BH}}`$-$`\mathrm{M}_{\mathrm{bulge}}`$ relation are necessary because these relations are probably more fundamental than the Faber-Jackson relation and are more closely related to the physics of black hole and galaxy formation. We noticed that our result on the $`\mathrm{M}_{\mathrm{BH}}`$-$`\mathrm{M}_{\mathrm{bulge}}`$ relation implies a rather flatter $`\mathrm{L}_{\mathrm{bulge}}`$-$`\sigma `$ relation than the Faber-Jackson relation. This seems to be supported by the finding of Nelson & Whittle (1996). However, the slope of the $`\mathrm{L}_{\mathrm{bulge}}`$-$`\sigma `$ relation depends sensitively on the statistical method, sample selection and uncertainties of both parameters. In addition, in deriving the $`\mathrm{M}_{\mathrm{BH}}/\mathrm{M}_{\mathrm{bulge}}`$ ratio, we adopted the same bulge mass-to-light ratio for Seyfert galaxies as for nearby galaxies obtained by Magorrian et al. (1998). Such a ratio may be smaller for Seyfert galaxies (Whittle 1992). Introducing a smaller bulge mass-to-light ratio may affect our results. We expect that further high quality observations on Seyfert galaxies, normal galaxies and quasars using the Hubble Space Telescope and larger ground-based telescopes could diminish the uncertainties in measuring the galactic bulge properties and central black hole masses of these objects. These efforts will undoubtedly help us to understand better the physics of formation and evolution of SMBHs and galaxies.

###### Acknowledgements.

We are very grateful to the referee, Todd Boroson, for helpful comments and suggestions. We thank Peter Biermann, Xinwu Cao, Jiansheng Chen, Zugan Deng, Jun Ma, Tinggui Wang, Hong Wu, Xiaoyang Xia and Suijian Xue for stimulating discussions. This work was partially supported by the Pandeng Project, the National Natural Science Foundation, and the National Key Basic Research Science Foundation (NKBRSF G19990752) in China.


text:
# Vortex waistlines and vortex gas

## Abstract

We investigate the manner in which a linear potential arises from fluctuations due to vortices linked with the Wilson loop. In particular, the customary naive picture is critically reexamined.

A large body of work over the last few years  supports the picture of long thick vortices as the most important long range degrees of freedom for confinement.

In this note we wish to examine more closely the manner in which a linear potential arises from vortex fluctuations. One customarily assumes a picture of randomly distributed vortices of certain thickness and basically arbitrary length. Furthermore, in the presence of a Wilson loop, one assumes that the free energy cost and hence the probability for a vortex to pierce anywhere inside the area enclosed by the loop is constant for all loop sizes. Summing then independently over all vortices linking with the loop one gets area law behavior for the loop expectation:

$`W[C]`$ $``$ $`1+(1)k|A_C|+{\displaystyle \frac{((1)k|A_C|)^2}{2!}}+\mathrm{}`$ (1)
$`=`$ $`\mathrm{exp}(k|A_C|).`$
How accurate is this picture?

We consider explicitly only the case of $`SU(2)`$ for simplicity. We work on a lattice $`\mathrm{\Lambda }`$ of linear extension $`L_\mu `$ in spacetime direction $`\mu =1,\mathrm{},d`$. Let $`𝒪[𝒱_{\mu \nu }]`$ denote the operator that flips the sign of the coupling (introduces a $`Z(2)`$ ‘twist’) on a coclosed set of plaquettes $`𝒱_{\mu \nu }`$ winding around the periodic lattice perpendicularly to the $`\mu \nu `$directions. The twist introduces a discontinuous gauge transformation with multivaluedness in $`Z(2)`$, i.e. a vortex trapped in the periodic lattice. The expectation of this operator defines the vortex free energy

$$\mathrm{exp}(F_v[\mu \nu ])=𝒪[𝒱_{\mu \nu }].$$
(2)
The expectation depends only on the directions in which $`𝒱`$ winds through the lattice, not its exact shape or location, nor the number mod 2 (mod $`N`$) of homologous sets $`𝒱`$ carrying a twist. This expresses the mod 2 conservation of flux, a simple consequence of which is

$$𝒪[𝒱_{\mu \nu }]𝒪[𝒱_{\mu \nu }^{}]=1.$$
(3)
Our argument in the following assumes that the vortex free energy (2) behaves as follows. For sufficiently large $`|A_{\mu \nu }|`$ and dimension $`d4`$

$$F_v[\mu \nu ](\underset{\lambda \mu \nu }{}L_\lambda )\mathrm{exp}(\rho (\beta )|A_{\mu \nu }|),$$
(4)
where $`|A_{\mu \nu }|=L_\mu L_\nu `$.

This is the optimal behavior under exponential transverse spreading of the flux introduced by the twist on $`𝒱`$, with $`\rho `$ approaching, at least asymptotically, the exact linear potential string tension. It is dictated by physical reasoning, and explicitly seen in strong coupling expansion. Recently, it has become possible to exhibit this behavior at large $`\beta `$ in extensive numerical simulations .

The $`Z(2)`$ Fourier transform of (2):

$$\mathrm{exp}(F_{\mathrm{el}})=\frac{1}{2}(\mathrm{\hspace{0.17em}1}𝒪[𝒱_{\mu \nu }])$$
(5)
gives the dual (w.r.t. the center) color-electric free energy. The mod 2 conservation of the magnetic flux is now expressed by the projection property

$`{\displaystyle \frac{1}{2}}(\mathrm{\hspace{0.17em}1}𝒪[𝒱_{\mu \nu }]){\displaystyle \frac{1}{2}}(\mathrm{\hspace{0.17em}1}𝒪[𝒱_{\mu \nu }^{}])`$
$`=(\mathrm{\hspace{0.17em}1}𝒪[𝒱_{\mu \nu }])/2.`$ (6)
There are several identities involving the Wilson loop operator and the operator $`𝒪[𝒱]`$ that can serve as the starting point for deriving rigorous inequalities relating the Wilson loop to the flux order parameters (2), (5) . A simple such relation, which is easily obtained and the only one needed here, is :

$$W[C]=\frac{1}{2}\mathrm{tr}U[C]\frac{1}{2}\left(\mathrm{\hspace{0.17em}1}𝒪[𝒱^{}]\right).$$
(7)
In (7), $`W[C]=\frac{1}{2}\mathrm{tr}U[C]`$ is the Wilson loop expectation for a loop $`C`$, and $`𝒱^{}`$ is a coclosed set of plaquettes linking with $`C`$.

Now by the property (6), one may insert multiple $`\frac{1}{2}(\mathrm{\hspace{0.17em}1}𝒪[𝒱])`$ factors in (7) corresponding to a collection of coclosed sets $`\{𝒱_i\}`$ linking with the loop $`C`$. Imagine enclosing each such $`𝒱_i`$ in a ‘vortex container’ , i.e. a sublattice $`\mathrm{\Lambda }_i\mathrm{\Lambda }`$ containing $`𝒱_i`$ and wrapping around $`C`$ (figure 1).

Furthermore, imagine integrating over the bond variables in the interior of each container, keeping bond variables on its boundary $`\mathrm{\Lambda }_i`$ fixed, and define

$$f_{\mathrm{\Lambda }_i}(U_{\mathrm{\Lambda }_i})=\frac{1}{2}\left(\mathrm{\hspace{0.17em}1}\frac{z_{\mathrm{\Lambda }_i}^{()}(U_{\mathrm{\Lambda }_i})}{z_{\mathrm{\Lambda }_i}^{(+)}(U_{\mathrm{\Lambda }_i})}\right).$$
(8)
In (8), $`z_{\mathrm{\Lambda }_i}^{()}(U_{\mathrm{\Lambda }_i})`$ denotes the partition function for (the interior of) $`\mathrm{\Lambda }_i`$ with $`()`$ or without $`(+)`$ twist on $`𝒱_i`$. Note that $`f_{\mathrm{\Lambda }_i}(U_{\mathrm{\Lambda }_i})`$ is nothing but the electric free energy now defined on a lattice $`\mathrm{\Lambda }_i`$ (the vortex container) with fixed (instead of periodic) b.c. in the directions transverse to $`𝒱_i`$. By these steps one then obtains from (7) the exact relation

$$W[C]=\frac{1}{2}\mathrm{tr}U[C]\underset{i}{}f_{\mathrm{\Lambda }_i}(U_{\mathrm{\Lambda }_i}),$$
(9)
which, in turn, immediately yields the bound:

$$W[C]\underset{i}{}\underset{U_{\mathrm{\Lambda }_i}}{\mathrm{max}}|f_{\mathrm{\Lambda }_i}(U_{\mathrm{\Lambda }_i})|.$$
(10)
In (10) the maximum is taken over all values of the bond variables on the boundary $`\mathrm{\Lambda }_i`$.

The bound given by (10) leads to area law according to the naive notion of independent vortices linking with the loop provided they grow long enough to be able to pierce at any point of the area spanning the loop. For this to happen one must assume, according to (4), that they grow thick enough to reach the regime:

$$\mathrm{ln}(z_{\mathrm{\Lambda }_i}^{()}(U_{\mathrm{\Lambda }_i})/z_{\mathrm{\Lambda }_i}^{(+)}(U_{\mathrm{\Lambda }_i}))|𝒱_i|\mathrm{exp}(\rho d_i^2)$$
(11)
for all $`U_{\mathrm{\Lambda }_i}`$. Here $`d_i`$ denotes the size of $`\mathrm{\Lambda }_i`$ in each of the two directions transverse to the enclosed set $`𝒱_i`$; its longitudinal size is given by $`|𝒱_i|`$. It follows that to keep the free energy cost of each vortex less than a fixed value $`f`$, we need:

$$d_i^2\frac{1}{\rho }\mathrm{ln}\left(|𝒱_i|/f\right).$$
(12)
Then also $`\mathrm{max}f_{\mathrm{\Lambda }_i}(U_{\mathrm{\Lambda }_i})<\frac{1}{2}(1e^f)<\mathrm{const}`$. But, $`|𝒱_i|`$ (in $`d`$ dimensions) is of order $`R^{d2}`$ in order to have linkage through points away from the perimeter of a rectangular loop of side lengths $`T`$ and $`R`$, $`T>R`$. Such a loop can then accommodate $`RT/\mathrm{ln}R`$ containers wrapped around it. Thus we obtain a confining but not quite purely linear potential

$$V(R)\text{const}R/\mathrm{ln}R.$$
(13)
By the same reasoning we now see that the familiar argument for area law summarized in (1) must also fail in the same manner. It ignores the actual free energy requirements for having vortices of sufficient length link with a large loop: the type of discussion just given above for (10) applies to each term in the summation in (1). Thus, for one vortex to link anywhere with a large loop of side lengths $`T`$ and $`R`$ ($`TR`$), at bounded fixed cost, a vortex cross section area of order $`\mathrm{ln}R`$ is required. This leads at best to (13), not (1).

The discrepancy comes from the treatment of vortices as localized and independent. For each vortex enclosed in a vortex container of fixed, but sufficiently large, width $`d`$, the vortex bulk free energy cost inside is correctly estimated from (11). In addition, however, the vortex is surrounded by the pure gauge long tail that encodes its nontrivial topology and flux quantization. The tail incurs no additional action cost, but is of infinite range, and communicates the presence of a vortex inside the container to other vortices, or other topological obstructions, outside enforcing flux conservation mod N. This long distance interaction allows a system of vortex excitations to adjust the amount of flux spreading, i.e adjust the thickness of vortices to minimize the free energy of the system. The thickness of vortex cores then is not fixed, but is adjusted relative to their length as required by the presence of other vortices and/or other obstructions (e.g. Wilson loop legs) sensitive to the presence of topological $`Z(N)`$ flux.

Thus, in the presence of the Wilson loop source, the optimal configurations for the system are not those of multiple isolated linked vortices, each of some fixed free energy cost, hence length $`|𝒱|R^{(d2)}`$ and fixed width $`d^2\mathrm{ln}RR,T`$. It is more advantageous, in terms of free energy cost, for multiple linking vortices to thicken and merge, the total topological flux being conserved mod N. Since the Wilson loop operator is affected by the topological flux through it only mod N, this should optimize the expectation.

The exact expression (9) holds for any number of factors in the product inside the expectation. In view of the above discussion, optimization requires combining containers into ones as thick as possible by integrating over the boundary fields of neighboring containers. With $`TR`$, this amounts to taking containers in the product in (9), (10) having transverse area $`R^2`$, and longitudinal extension $`(\text{const}R)^{(d2)}`$. (10) now gives

$$V(R)\text{const}R\text{const}(\mathrm{ln}R+\text{const})/R$$
(14)
replacing (13). For loops with $`T\stackrel{}{\stackrel{>}{}}R`$, strict linear potential arises essentially from thick vortex fluctuations spanning the entire loop area.

Our discussion shows that an effective picture of the long distance confining fluctuations as isolated, independent vortices winding over long distances - i.e., as some kind of an approximately dilute or weakly interacting vortex gas - is not generally applicable. In general vortices cannot be considered isolated, and a definite number of vortices, specified more precisely than mod N, cannot necessarily be unambiguously assigned to every configuration. This accords with the picture of ‘condensation’ of vortex flux and ‘percolated’ vortices over distances $``$ 1 fm, as suggested by the simulations. The effective description of this vacuum is discussed in