text:
# Measuring the metric: a parametrized post-Friedmanian approach to the cosmic dark energy problem

## I INTRODUCTION

Modern cosmology is in a somewhat equivocal state of affairs. The good news is that a recent avalanche of high-quality data are well fit by an emerging “standard model” whose roughly ten free parameters are being constrained with increasing precision . The bad news is that this emerging model is more complicated than anticipated. There is not one kind of invisible substance, but three: cold dark matter (CDM) to explain clustering, dark energy to close the Universe and a smidgeon of massive neutrinos  that may well be too small to be cosmologically important. Moreover, problems involving small-scale clustering have triggered increasingly complicated models for the dark matter — self-interacting CDM , annihilating CDM , warm dark matter , fuzzy dark matter  and fluid dark matter , to mention a few — and a large number of dark energy models (e.g. ) have appeared where Einstein’s single parameter $`\mathrm{\Lambda }`$ is replaced by a “quintessence” field that varies temporally and perhaps spatially.

This perceived profusion of bells and whistles has caused unease among some cosmologists  and prompted concern that these complicated dark matter flavors constitute a modern form of epicycles. There have been numerous suggestions that, just as in the days of Ptolemy, the apparent complications can be eliminated by modifying the laws of gravity  For instance, there are scalar-tensor theories that can reproduce the observed accelerating Universe without invoking a cosmological constant .

It will undoubtedly take time to settle these issues. In the interim, however, it is desirable to quote cosmological measurements in a language that is fairly theory-independent. This is the topic of the present paper, delimited to gravity in the linear regime, on large scales.

### A Parametrizing our ignorance with two functions

There is general consensus that space-time is well described by a metric, at least on non-microscopic scales, so competing theories of gravity differ in their predictions for how the metric evolves with time and responds to the presence of matter. General Relativity (GR) makes specific predictions in the form of the Einstein Field Equations. The parametrized post-Newtonian approximation  has been highly successful in describing possible departures from GR in the solar system, and the parametrized post-Keplerian approximation  has done the same for for binary orbits, notably binary pulsars. Although GR has so far passed all experimental tests with flying colors, these tests tend to probe only the present epoch. Broad classes of theories have been shown to evolve into GR at late times, so substantial departures from the Einstein field equations are still possible in cosmology.

Assuming merely the cosmological principle (that our position is not special), it follows directly from the observed near-isotropy of the cosmic microwave background (CMB) that the Universe is approximately isotropic and homogeneous. This means that the large-scale metric consists of small perturbations on a Friedman-Robertson-Walker (FRW) background :

$$ds^2=dt^2a(t)^2\left[\frac{dr}{1kr^2}r^2d\theta ^2\mathrm{sin}^2\theta d\varphi ^2\right]$$
(1)
for some function $`a(t)`$ that describes the expansion history of the Universe. $`k=1`$, $`0`$ or $`1`$ for an open, flat and closed Universe, respectively. If GR is correct, then $`a`$ evolves according to the Friedman equation

$$H^2\left(\frac{\dot{a}}{a}\right)^2=\frac{8\pi G}{3}\rho ,$$
(2)
where we define the mean matter density $`\rho (z)`$ to include a curvature contribution $`\rho _k\frac{3kc^2}{8\pi Ga^2}`$. Other theories make different predictions. As summarized in Figure 1 and reviewed in the next section, the function $`H(z)`$ is directly measurable is a variety of ways that make no assumptions about dark matter properties or gravitational field equations, but rely on spacetime geometry alone.

Let us now turn to the evolution of perturbations. The evolution equations are by definition linear to first order, regardless what underlying theory of gravity is being linearized. Assuming that this theory is well approximated by partial differential equations on large scales, these decouple into ordinary differential equations for each Fourier mode, since we are perturbing around a homogeneous metric and derivatives become local in Fourier space. In short, the key observable predicted by theory is the linear growth factor $`g(z,k)`$, conventionally defined so that a density perturbation $`\widehat{\delta }(𝐤)`$ with wave vector $`𝐤`$ evolves with redshift as

$$\delta (z,𝐤)=\frac{g(z,k)}{(1+z)}\delta (0;k).$$
(3)
A wide variety of upcoming observations will probe $`g(z,k)`$, including galaxy and Lyman $`\alpha `$ forest clustering, the cosmic microwave background, gravitational lensing and type Ia supernovae. When multiple matter components are present, there will of course be one growth factor for each one, whose evolution couple. For photons and neutrinos, the full phase-space distribution is relevant to the dynamics. On large scales, a broad class of models make the simple prediction that $`g(z,k)`$ is independent of $`k`$.

### B Is matter the matter?

It is obviously premature to guess as to whether upcoming data will favor dark matter or modified gravity. However, it is worth noting that the distinction between these two cases is in fact somewhat blurry. As pointed out by Eddington , we can always choose to define matter as that which equals $`(8\pi G)^1`$ times the Einstein tensor, in which case the Einstein field equations $`G_{\mu \nu }=8\pi GT_{\mu \nu }`$ become satisfied by definition. This “matter” $`T_{\mu \nu }`$ would have very strange properties for a generic metric, so the predictive power of GR arises from the fact that observed matter has simple, formalizable properties.

For many modifications of GR, the so defined $`T_{\mu \nu }`$ is in fact rather regular. Indeed, in the familiar case of adding a cosmological constant, it is largely a matter of taste whether to call it modified gravity or dark energy, corresponding merely to whether we insert it on the left or right hand side of the Einstein equations. In the currently popular class of GR generalizations known as scalar-tensor theories , of which Brans-Dicke theory is a special case, we can play the same trick: the Friedman equation (2) remains valid in the conventional (Jordan) frame  if we augment the density $`\rho `$ with the extra “dark energy” term

$$\rho _x=\left(\frac{1}{\zeta F}1\right)\rho _m+\frac{1}{8\pi G_{\mathrm{eff}}F}\left[U+\frac{Z\dot{\mathrm{\Phi }}^2}{2}3H\dot{F}\right],$$
(4)
where the theory is specified by the free functions $`F(\mathrm{\Phi })`$, $`U(\mathrm{\Phi })`$ and $`Z(\mathrm{\Phi })`$ of the dilaton field $`\mathrm{\Phi }`$, $`\rho _m`$ is the density of conventional matter, and

$$\zeta \frac{G_{\mathrm{eff}}}{G_{\mathrm{bare}}}=\frac{1}{F}\left(\frac{2ZF+4(F^{})^2}{2ZF+3(F^{})^2}\right)$$
(5)
is the factor by which $`G_{\mathrm{eff}}`$, the effective small-scale Newton’s constant measured in laboratory experiments today, differs from the “bare” one $`G_{\mathrm{bare}}`$ which is truly constant.

As another example, consider a flat toy model where equation (2) is replaced by

$$H^2=\frac{8\pi G}{3}\frac{\rho ^2}{\rho _{}}$$
(6)
for some constant $`\rho _{}`$, i.e., a model where $`\rho H`$ rather than $`\rho H^2`$. An observer assuming that the Friedman equation equation (2) is valid would then conclude that there is a “dark energy” density $`\rho _x`$ given by

$$\rho _x\frac{3H^2}{8\pi G}\rho =\frac{\rho ^2}{\rho _{}}\rho =\frac{\rho _0^2}{\rho _{}}(1+z)^6\rho _0(1+z)^3,$$
(7)
a quantity which goes negative for small $`\rho `$. Here $`\rho _0`$ denotes the current density of baryonic and dark matter.

The modified behavior of linear perturbations in non-GR theories can clearly be attributed to new “matter” as well. Indeed, it has recently been shown  that the linearized properties of very general types of dark matter can be conveniently characterized by an effective sound speed and an effective viscosity.

In light of this theory/matter ambiguity, we propose to use the Friedman equation as a definition of an effective mean density

$$\rho _{\mathrm{eff}}(z)\frac{3H(z)^2}{8\pi G_{\mathrm{eff}}}.$$
(8)
Thus $`\rho _{\mathrm{eff}}(z)=\rho (z)`$, the true density, if the Friedmann equation is correct. As illustrated in Figure 1, expressing observational results in terms of $`\rho _{\mathrm{eff}}(z)`$ rather than $`H(z)`$ has the advantage of directly visualizing the dark matter problem. The contribution to $`\rho (z)`$ from vacuum density is of course constant. Curvature density scales as $`(1+z)^2`$, matter density as $`(1+z)^3`$, radiation density as $`(1+z)^4`$, and a “quintessence”. component with constant equation of state $`w`$ as $`(1+z)^{3(1+w)}`$.

A large body of work has been performed during the last few years on constraints and forecasts for more specific models. Cast in our present notation, a particularly useful and popular parametrization of $`\rho (z)`$ has been

$$\rho (z)=\rho _m^0(1+z)^3+\rho _\kappa ^0(1+z)^2+\rho _X^0(1+z)^{3(1+w)}$$
(9)
for $`z10^3`$, where $`\rho _m^0`$, $`\rho _\kappa ^0`$ and $`\rho _X^0`$ are the present densities of matter, curvature and dark energy, respectively, and the dark energy equation of state $`w`$ may in turn vary with redshift. In this paper, we treat $`\rho _{\mathrm{eff}}(z)`$ rather than $`w(z)`$ as the free function for two reasons:

1. It can be estimated from the data without assumptions about gravity (the Einstein field equations) or about the dark energy (which may be out of equilibrium and not have any simple equation of state uniquely determining its pressure from its density).
2. It is more directly related to observed data since, as will be elaborated below, it is given by the first rather than the second derivative of a measured function.

A similar choice was made in  as well.

Much progress has also been made on a reconstruction program, aiming to reconstruct the microphysics of specific theories from data. The present work complements this approach, since such reconstructions of, e.g., the quintessence potential  or the functions specifying scalar-tensor gravity  can be performed with $`\rho _{\mathrm{eff}}(z)`$ and $`g(z,k)`$ as a starting point.

The rest of this paper is focused on the 0th order function $`\rho _{\mathrm{eff}}(z)`$, with the 1st order function $`g(z,k)`$ studied in a companion paper . In the next section, we review techniques for measuring $`\rho _{\mathrm{eff}}(z)`$, compute the constraints from current SN 1a data and make forecasts for the proposed SNAP satellite. We discuss our conclusions in Section III.

## II Measuring $`\rho _{\mathrm{eff}}(z)`$

The bars at the top of Figure 1 show the redshift ranges over which various types of observations are sensitive to $`\rho _{\mathrm{eff}}(z)`$. We will first review them briefly in turn, then return the SN 1a in greater detail.

### A Space-time geometric tests

All cosmological tests based on luminosity, angular size and age (see e.g. ) can be described as noisy measurements of some quantities $`x_n`$ at redshifts $`z_n`$, $`n=1,\mathrm{},N`$, modeled as

$$x_n=a\mathrm{ln}d(z)+b_n+\epsilon _n,$$
(10)
where $`a`$ and $`b`$ are constants independent of $`\rho _{\mathrm{eff}}(z)`$, the function $`d`$ incorporates the effects of cosmology and $`\epsilon _n`$ is a random term with zero mean ($`\epsilon _n=0`$) including all sources of measurement error. As briefly review below, these observables all probe different weighted averages of the quantity $`\rho _{\mathrm{eff}}(z)^{1/2}`$.

For luminosity tests like SN Ia, $`x_n`$ is the observed magnitude of the $`n^{th}`$ object and $`d`$ is the dimensionless luminosity distance :

$$d_{lum}=(1+z)\frac{S(\kappa \eta )}{\kappa },\eta (z)=H_0_0^z\frac{dz^{}}{H(z^{})},$$
(11)
where the $`\kappa \sqrt{|1\mathrm{\Omega }_0|}`$ is the spatial curvature, and $`\mathrm{\Omega }_08\pi G\rho _{\mathrm{eff}}(0)/3H^2`$ is the fraction of critical density contributed by $`\rho _{\mathrm{eff}}`$ today. $`S(x)\mathrm{sinh}x`$, $`x`$ or $`\mathrm{sin}x`$ when $`k=1`$, $`0`$ or $`1`$, respectively. From the definition of magnitudes, $`a=5/\mathrm{ln}10`$ and $`b_n=25+5\mathrm{lg}[H_0^1c/1\mathrm{M}\mathrm{p}\mathrm{c}]+M_n`$, where $`M_n`$ is the (hopefully known) absolute magnitude of the $`n^{th}`$ object. The errors $`\epsilon _n`$ include errors in extinction correction and intrinsic scatter in the “standard candle” luminosity.

For tests involving the observed angular sizes $`\theta _n`$ of objects at redshifts $`z_1,\mathrm{},z_N`$ with (hopefully known) linear sizes $`D_1,\mathrm{},D_N`$, we can take $`d`$ to be the dimensionless angular size distance : $`d_{ang}=d_{lum}/(1+z)^2`$. This gives $`x_n\mathrm{ln}\theta _n`$ (angles measured in radians), $`a=1`$ and $`b_n\mathrm{ln}[H_0D_n/c]`$. For such tests (e.g. ), $`\epsilon _n`$ includes scatter in the “standard yardstick” size.

For tests involving estimates $`t_n`$ of the age of the Universe at redshifts $`z_n`$, we define $`x_n\mathrm{ln}H_0t_n`$. Setting $`a=1`$ and $`b_n=0`$, this gives

$$d_{age}=H_0_z^{\mathrm{}}\frac{dz^{}}{(1+z^{})H(z^{})}.$$
(12)
Equation (10) implicitly assumes that the the measurement errors $`\epsilon _n`$ for $`x_n`$ are additive. This is a popular approximation in many cases, notably for luminosity distance where the actual measurement is a magnitude proportional to $`\mathrm{ln}d`$, and tends to be accurate if the relative error (on flux, angular size, age, etc.) is small ($`1`$). In other cases, a more appropriate error model should be employed.

Ideally, we would like to measure the function $`\rho _{\mathrm{eff}}(z)`$ in ways that are sensitive only to the global spacetime geometry and do not depend on $`g(z,k)`$. The measurements that come closest to this ideal are arguably those involving the brightness, angular size and age of distant objects, in situations where the calibration of the “standard candle”, “standard ruler” or “standard clock” depends mainly on microphysics. However, great care must be taken when constraining competing theories of gravity etc. (as opposed to more mundane dark energy), since slight time variations in microphysical constants could affect the standard calibrators and masquerade as modifications to $`\rho _{\mathrm{eff}}(z)`$. For instance, time-variation in a non-gravitational “constant” (say the fine-structure constant “alpha”) could alter the intrinsic energetics and time-evolution of a supernova explosion, making it a bad standard candle. If such an effect were present but neglected, then we would draw an incorrect conclusion about the luminosity distance to the supernova and consequently also about the cosmic geometry.

Terrestrial constraints are are now so strong that time-variation of the fine structure constant $`\alpha `$ is unlikely to cause such misinterpretations: variations in $`\alpha `$ are limited to around a part per million over cosmological timescales  (c.f. ). A more serious concern in this context is time variations in the effective gravitational constant, expected in many scalar-tensor theories according to equation (5). Here the observational constraints are weaker, permitting variations of a few percent over cosmological time scales , which could among other things cause distant supernovae to be systematically brighter or fainter than nearby ones.

As Figure 1 indicates, the function $`\rho _{\mathrm{eff}}(z)`$ is probed by the growth of clustering as well, via the power spectrum growth in the Ly$`\alpha `$ forest, weak lensing and galaxy surveys as well as via the abundance evolution of galaxy clusters and lensing systems. However, these probes all depend on the function $`g(z,k)`$ as well, i.e., on the next order of perturbation theory. The Alcock-Paczynski test elegantly eliminates the dependence on $`g(z,k)`$ caused by number density evolution, but involves certain bias-related issues .

Although the CMB power spectrum $`C_ł`$ depends on the growth $`g`$, the effect is “vertical” rather than “horizontal”. In other words, as long as some detectable acoustic peak structure is visibly preserved, it should be possible to extract the angular scale corresponding to the horizon size at recombination fairly independently of $`g`$. With additional knowledge of the dark matter density, this provides a clean measurement of $`d_{lum}`$ at $`z10^3`$.

If GR is correct so that $`g`$ and the vertical structure of the acoustic peak heights can be computed from first principles, then the CMB is a sensitive probe of the high-redshift asymptote in of the $`\rho _{\mathrm{eff}}`$ curve in Figure 2, where the Universe becomes matter dominated. Writing $`\rho _m(z)=\rho _m^0(1+z)^3`$, the MAP and Planck satellites are expected to measure $`\rho _m^0`$ to accuracies of $`10\%`$ and $`1\%`$, respectively. However, in contrast to, e.g., SN 1a, the CMB is mainly sensitive to the physical matter density $`\rho _m`$, not to the plotted quantity $`\rho _m/h^2`$, since this is what matters for the acoustic oscillations at $`z>10^3`$. The ultimate accuracy with which we can measure $`\rho _m^0/h^2`$ (equivalently $`\mathrm{\Omega }_m`$) is therefore likely to be limited not by CMB issues but by our ability to measure $`h`$, either directly or by combining various complementary measurements .

### B $`\rho _{\mathrm{eff}}(z)`$ from existing SN 1a data

Of the various low-redshift probes of $`\rho _{\mathrm{eff}}(z)`$ that we discussed, type 1a supernovae are arguably the most sensitive at the present time. We will therefore compute constraints on $`\rho _{\mathrm{eff}}(z)`$ from current SN 1a data and make forecasts for the accuracy attainable with future measurements.

Figure 2 shows the measured luminosity distance for 92 SN 1a, combining the high redshift sample of 42 from  with the 50 reported in  with the MLCS method. Apart from systematic and calibration issues, our desired function $`\rho _{\mathrm{eff}}(z)`$ is simply given by the derivative of the curve around which they scatter: $`\rho _{\mathrm{eff}}(z)=3H_0^2/8\pi G_{\mathrm{eff}}\eta ^{}(z)^2h^2/\eta ^{}(z)^2`$.*** From here on, we will limit our discussion to the currently favored case of flat space, i.e., $`k=0`$, where $`\eta (z)=d_{lum}(z)/(1+z)`$ is directly measured for each supernova. In the general case, $`\eta (z)=k^1\mathrm{sin}^1[kd_{lum}(z)/(1+z)]`$ or $`\eta (z)=k^1\mathrm{sinh}^1[kd_{lum}(z)/(1+z)]`$ depending on the sign of the curvature, introducing an extra degeneracy with $`k`$ that has been exhaustively treated elsewhere . Since CMB experiments are measuring the acoustic peak locations to increasing precision, and these are mainly sensitive to curvature (given by $`\rho _{\mathrm{eff}}(0)/h^2`$), degeneracies between different evolution scenarios for $`\rho _{\mathrm{eff}}(z)`$ will be more important than the curvature degeneracy in the post-MAP era.  Since the slope is seen to fall, the cosmic density clearly increases at high redshift. If standard GR is correct ($`\rho _{\mathrm{eff}}=\rho `$), then this fact that $`d\rho _{\mathrm{eff}}/dz>0`$ corresponds to the weak energy condition  that the density measured by any observer is non-negative, and is equivalent to the pressure-to-density ratio constraint $`p/\rho 1`$. A more detailed reconstruction of $`\rho (z)`$ is shown in Figure 3, and was computed by fitting a homogeneous quartic polynomial (dashed curve in Figure 2) to $`\eta (z)`$. Using the notation of equation (10), the fit was done by minimizing

$$\chi ^2=\underset{n=1}{\overset{n}{}}\left(\frac{x_na\mathrm{ln}d(z_n)b_n}{\sigma _n}\right)^2,$$
(13)
where $`\sigma _n`$ are the quoted 1 sigma magnitude errors . This gives the fitting parameters (polynomial coefficients) as linear combinations of the data $`x_n`$. The 1 sigma errors in Figure 3 are obtained by propagating the errors on these fitting parameters.Equation (13) assumes that the magnitude errors on the different supernovae are uncorrelated. It is straightforward to generalize the procedure to the case of correlated errors, where the $`N\times N`$ magnitude covariance matrix is non-diagonal given a physical model for the nature of these correlations. A perfectly correlated term, corresponding to a common offset in all magnitudes, would have no effect on the results since it is only the derivative that matters in Figure 2, but redshift-dependent offsets can of course bias the results. Indeed, a deeper understanding of possible systematic errors of this type is the main challenges for taking full advantage of the statistical power of future experiments like SNAP.  The figure clearly illustrates the familiar fact that dark energy/modified gravity is required: the slope of the $`\rho _{\mathrm{eff}}(z)`$-curve is too shallow at low redshift to be explainable in terms of ordinary matter $`(\rho [1+z]^3`$), spatial curvature $`(\rho [1+z]^2`$) or some combination of the two.

Since $`\rho _{\mathrm{eff}}(z)h^2/\eta ^{}(z)^2`$, the SN 1a give us a direct measurement of $`\rho _{\mathrm{eff}}(z)/h^2`$, where $`h`$ is the Hubble parameter today in units of 100 km/s/Mpc. Changing $`h`$ thus shifts $`\rho _{\mathrm{eff}}`$ vertically on a logarithmic plot like Figure 2 without affecting its problematic slope.

To assess the sensitivity of our results to method details, we show alternative analyses in figures 4 and 5. Here we replaced the polynominal fitting function by piecewise quadratic and piecewise linear interpolating functions for $`\eta (z)`$, respectively. For the quadratic case, we enforced that the function have continuous derivative, thereby making $`\rho _{\mathrm{eff}}`$ a continuous function. The results are seen to be quite robust, in all cases agreeing well with the “concordance” model and indicating a shallow slope at low redshifts incompatible with just matter and curvature alone. A very nice non-parametric SN 1a analysis was recently performed by Wang & Garnavich , and it is interesting to compare the two studies since they were performed concurrently and independently. The quantitative agreement is good taking into account that  treat the dark energy density $`\rho (z)\rho _m^0(1+z)^3`$ rather than $`\rho _{\mathrm{eff}}(z)`$ as the free function and marginalize over $`\rho _m`$ imposing the weak energy condition. The SN 1a data thus speak loud and clear: the dark energy puzzle is robust to theoretical assumptions and method details.

### C With upcoming SN 1a data

The proposed SNAP satelliteMore information on the SNAP mission is available at
$`http://snap.lbl.gov`$ would measure of order $`N=2000`$ SN 1a out to redshifts around 1.2 or higher . There has been substantial debate in the literature about the accuracy with which the equation of state can be reconstructed . We will now compute how accurately SNAP can measure $`\rho _{\mathrm{eff}}(z)`$. It is straightforward to do this numerically in a “black box” fashion, and figures 3-5 show the result of doing this with simulations in the same way as we did it for the currently existing data. However, because of the strong interest and substantial resources focused on designing such future probes, it is also worthwhile to gain intuition on how experimental specifications translate into model constraints. We will therefore derive an analytic result which has the advantage of explicitly showing how the results scale when changing survey specifications, the smoothing scale, etc.

#### 1 An analytic error formula

Defining $`f(z)H(0)/H(z)`$, we have $`\eta (z)=_0^zf(z^{})𝑑z^{}`$ and wish to know how accurately $`f`$ can be recovered from noisy measurements of $`\eta `$. In other words, we simply want to compute the derivative of a function measured with sparse sampling and noise. The errors on $`\rho _{\mathrm{eff}}(z)1/f(z)^2`$ then follow trivially from those on $`f(z)`$.

The accuracy with which $`f(z)`$ can be measured is readily computed with the Fisher information matrix formalism  with infinitely many parameters to be estimated (the value of $`f`$ at each $`z`$). The amount of information about these parameters in the SN 1a data is given by

$$𝐅(z^{},z^{\prime \prime })=\frac{a^2N}{(\mathrm{\Delta }m)^2}_0^{\mathrm{}}g(z)w_z^{}(z)w_{z^{\prime \prime }}(z)𝑑z,$$
(14)
where $`a5/\mathrm{ln}10`$, $`\mathrm{\Delta }m`$ is the rms SN 1a magnitude error, $`N`$ is the number of supernovae, $`g(z)`$ is their redshift distribution (normalized to integrate to unity), and

$$w_z^{}(z)\frac{\mathrm{ln}\eta (z)}{f(z^{})}=\frac{\theta (zz^{})}{\eta (z)}.$$
(15)
since $`\eta ^{}(z)=f(z)`$. Here $`\theta `$ denotes the Heaviside step function ($`\theta (z)=1`$ for $`z0`$, vanishing otherwise). Equation (14) thus gives

$`𝐅(z,z^{})`$ $`=`$ $`{\displaystyle \theta (z^{\prime \prime }z)\theta (z^{\prime \prime }z^{})h(z^{\prime \prime })^1𝑑z^{\prime \prime }}`$ (16)
$`=`$ $`{\displaystyle _{\mathrm{max}(z,z^{})}^{\mathrm{}}}h(z^{\prime \prime })^1𝑑z^{\prime \prime },\text{where}`$ (17)
$`h(z)`$ $``$ $`{\displaystyle \frac{(\mathrm{\Delta }m)^2}{a^2N}}{\displaystyle \frac{\eta (z)^2}{g(z)}}.`$ (18)
The covariance matrix $`𝐂`$ giving the best attainable error bars on our estimated parameters is the inverse of this infinite-dimensional matrix, $`𝐂=𝐅^1`$. Fortunately, this inverse can be computed analytically, giving

$$𝐂(z,z^{})=h^{}(z)\delta ^{}(zz^{})h(z)\delta ^{\prime \prime }(zz^{}),$$
(19)
where $`\delta `$ denotes the Dirac delta function. This follows from the identity

$$𝐂(z,z^{\prime \prime })𝐅(z^{\prime \prime },z^{})𝑑z^{\prime \prime }=\delta (zz^{}),$$
(20)
which is proven by substituting equations (16) and (19), integrating by parts twice and using that $`\theta ^{}(x)=\delta (x)`$. In practice we obviously want to smooth our estimated $`f`$-function ($`\stackrel{~}{f}(z)`$, say) to reduce noise. Using weighted averages $`\stackrel{~}{f}_i\psi _i(z)\stackrel{~}{f}(z)𝑑z`$, the covariance matrix $`𝐂_{ij}`$ for these averages is

$$𝐂(z,z^{})\psi _i(z)\psi _j(z^{})𝑑z𝑑z^{}=h(z)\psi _i^{}(z)\psi _j^{}(z)𝑑z,$$
(21)
employing equation (19) and integrating by parts. Once these errors on $`f(z)`$ have been computed, it is of course trivial to map them into errors on $`\rho _{\mathrm{eff}}(z)`$. When the relative errors are small $`(\mathrm{\Delta }ff)`$, we have simply $`\mathrm{\Delta }\rho _{\mathrm{eff}}/\rho _{\mathrm{eff}}\mathrm{\Delta }f/f`$.

#### 2 Smoothing functions

Using Gaussian smoothing functions

$$\psi _i(z)\frac{1}{(2\pi )^{1/2}\mathrm{\Delta }z}e^{\frac{1}{2}\left(\frac{zz_i}{\mathrm{\Delta }z}\right)^2},$$
(22)
the error bars $`\mathrm{\Delta }f_i𝐂_{ii}^{1/2}`$ are given by

$$\mathrm{\Delta }f_i\frac{\mathrm{ln}10}{10\pi ^{1/4}}\frac{\mathrm{\Delta }m\eta (z_i)}{(\mathrm{\Delta }z)^{3/2}N^{1/2}g(z_i)^{1/2}}$$
(23)
for the case where the smoothing scale $`\mathrm{\Delta }z`$ is much smaller than the scale on which $`\eta `$ and $`g`$ vary substantially. In the same approximation, the dimensionless correlation coefficients between the different averages are

$$r_{ij}(2n^2)e^{n^2},\text{where}n\left(\frac{z_iz_j}{\mathrm{\Delta }_z}\right)$$
(24)
is the number of smoothing lengths by which the two redshifts are separated, so neighboring measurements are uncorrelated if they are separated by $`\sqrt{2}\mathrm{\Delta }z`$. Another interesting case for our discussion below is that of “triangular” smoothing functions

$$\psi _i(z)=\{\begin{array}{cc}\frac{\mathrm{\Delta }z|zz_i|}{\mathrm{\Delta }z}\hfill & \text{if }|zz_i|<\mathrm{\Delta }z\text{,}\hfill \\ \multicolumn{2}{c}{}\\ 0\hfill & \text{otherwise,}\hfill \end{array}$$
(25)
which gives

$$\mathrm{\Delta }f_i\frac{\sqrt{2}\mathrm{ln}10}{5}\frac{\mathrm{\Delta }m\eta (z_i)}{(\mathrm{\Delta }z)^{3/2}N^{1/2}g(z_i)^{1/2}}$$
(26)
in the same approximation. Here there are strong anticorrelations $`r=2^{1/2}0.7`$ between neighboring bins, and all longer-range correlations vanish. “Boxchar” smoothing where we simply average in bins, $`\psi (z)=\theta (\mathrm{\Delta }z/2|zz_i|)/\mathrm{\Delta }z`$, may seem like a natural choice. This is a disaster, however, since the resulting $`\delta `$-functions in $`\psi _i^{}`$ blow up when squared. This problem is easy to understand physically. In a toy example where we have many SN 1a equispaced in redshift and estimate $`\eta ^{}(z)`$ by simply differencing neighboring supernovae, averaging a segment of these derivative estimates will reduce to simply differencing the first and the last, placing all the statistical weight on merely two objects.

As we will return to below, our black-box numerical method can also understood in terms of smoothing functions of a particular form.

The exact result of equation (21) and the approximation of equation (23) illustrate a number of issues relevant to measuring $`\rho _{\mathrm{eff}}(z)`$ in practice, as discusses in the following subsections.

#### 3 Redshift resolution

Errors scale as $`(\mathrm{\Delta }z)^{3/2}`$ rather than $`(\mathrm{\Delta }z)^{1/2}`$ because of the derivative nature of what we are measuring, which means that very large number of SN 1a are needed to probe the small-scale structure of $`\rho _{\mathrm{eff}}(z)`$. It also means that care must be taken in interpreting error forecasts to avoid coming away with a misleadingly pessimistic impression. The superiority of SNAP over existing data is partly hidden in figures 4 and figures 5 since $`\mathrm{\Delta }z`$ is halved for SNAP, roughly tripling error bars. The error bars from the numerical SNAP simulations in figures 4 and figures 5 agree well with the analytic ones of equation (25) with $`\mathrm{\Delta }z=0.1`$ as well as those of , but are substantially larger than those from the polynomial fit shown in Figure 3. This is because the quartic polynomial enforces more smoothness, effectively increasing $`\mathrm{\Delta }z`$. The non-parametric error bars in Figure 3 from equation (23) using Gaussian smoothing are a factor of $`2^{3/2}\pi ^{1/4}4`$ smaller than those from equation (26) using triangle smoothing, which is related to the fact that the triangle smoothing functions have standard deviations $`\sqrt{6}`$ times narrower than the Gaussian ones.

Which error bars are most appropriate depends on what the data are to be used for. Constraining models where $`\rho _{\mathrm{eff}}`$ may change abruptly with $`z`$ requires high redshift resolution and associated large error bars like in fugures 4 and 5. Most published models, however, predict rather smooth functions $`\rho _{\mathrm{eff}}(z)`$, so the small error bars from a low-order polynomial fit like in Figure 3 give the most accurate representation of how well SNAP could distinguish them from one another.

#### 4 High-$`z`$ sensitivity

Since errors scale as $`\eta (z)`$, one does worse at large $`z`$. However, since $`\eta (z)`$ flattens out around $`z1`$ and asymptotes to a constant as $`z\mathrm{}`$, this is not a fundamental limitation to probing very large redshifts — the challenge is simply to find suitable high-$`z`$ standard candles or other reference objects to prevent $`g(z)`$ from going to zero in the denominator.

#### 5 Optimal estimation

Equation (21) also indicates what the optimal estimator is for non-parametric measurement of $`\rho _{\mathrm{eff}}(z)`$. It is easy to verify that when the smoothing scale $`\mathrm{\Delta }z`$ is much smaller than the scale $`h/|h^{}|`$ on which $`h`$ varies substantially, weighting the actual data with the function $`\psi _i^{}(z)h(z)`$ and normalizing appropriately will recover the minimal Fisher matrix error bars. In practice, $`g(z)`$ is of course not a smooth function but a sum of $`\delta `$-functions at the observed SN 1a redshifts, so an additional correction for the non-even spacing between adjacent SN 1a is necessary.

The black box of $`\chi ^2`$-minimization can also be understood in terms of smoothing functions and equation (21). Since $`\chi ^2`$ depends quadratically on the free parameters, the best fit parameters and hence also the recovered function $`\stackrel{~}{f}(z)`$ depends linearly on the observed $`\eta `$-values for each supernova, once the redshifts are given. In other words, we can write

$$\stackrel{~}{f}_i=\underset{j=1}{\overset{N}{}}\psi _i^{}(z_j)\eta _j$$
(27)
for some function $`\psi _i`$, where $`\eta _j`$ is the observed $`\eta `$-value for the $`j^{th}`$ supernova. In an attempt to demystify the numerical results somewhat, we have plotted these weight functions $`\psi _i^{}`$ in Figure 6 for estimating $`f(z)`$ with SNAP at $`z0.8`$. We see that for both the quadratic spline and the linear interpolation cases, there is substantial ringing. Moreover, the measurement of $`\eta ^{}(0.8)`$ is seen to involve supernovae over a very broad range of redshifts, not merely near $`z=0.8`$. In contrast, the triangle smoothing function would be more local, giving $`\psi ^{}(z)=1`$ for $`z=0.70.8`$, $`\psi ^{}(z)=1`$ for $`z=0.80.9`$ and zero elsewhere, simply estimating the derivative $`\eta ^{}(0.8)`$ by subtracting measurements to the left from ones to the right. The Gaussian smoothing function behaves similarly, merely with less abrupt behavior. Just as sophisticated data analysis techniques developed by numerous groups have enhanced the science return of recent CMB missions, there may well be room for methodological improvement here. For instance, what is a good parametrization of $`f`$ that leads to compact, easy-to-interpret weight functions?

## III DISCUSSION

We have argued that the effective density evolution $`\rho _{\mathrm{eff}}(z)`$ and the linear growth factor $`g(z,k)`$ constitute the natural meeting point for theory to confront observation in terms of linear cosmology, since they can be measured with very few assumptions about unseen matter and the underlying theory of gravity. This generalizes studies of the equation of state $`w`$.

We analyzed current SN 1a data in this framework, and showed that the statistical constraints on $`\rho _{\mathrm{eff}}(z)`$ are already fairly strong for $`z\mathrm{}<0.5`$. This alternative approach confirms the results of previous SN 1a analyses , and provides a clean and model-independent case for either dark energy or modified gravity, since $`\rho _{\mathrm{eff}}(z)`$ has the wrong logarithmic slope at low redshift to be explained by either matter (slope $`3`$) or spatial curvature (slope $`2`$).

We estimated the accuracy with which the proposed SNAP satellite should be able to measure $`\rho _{\mathrm{eff}}(z)`$, both with a numerical simulation involving polynomial fitting and by deriving a non-parametric analytic formula, and found that precision density measurements should be attainable in about a dozen independent redshift bins out to $`z1.2`$. Since errors grow as $`(\mathrm{\Delta }z)^{3/2}`$ as the bin width $`\mathrm{\Delta }z`$ is reduced, the large numbers of SN 1a of SNAP are absolutely necessary to resolve the small-scale structure of $`\rho _{\mathrm{eff}}(z)`$. The reason that $`w(z)`$ is so much harder to constrain accurately is that it involves a second derivative of the data, causing error bars to blow up as $`(\mathrm{\Delta }z)^{5/2}`$.

We focused on the 0th order function $`\rho _{\mathrm{eff}}(z)`$, with the 1st order function $`g(z,k)`$ studied in detail in a companion paper . To place our discussion in its observational context, Figure 7 shows the rough ranges of scale and redshift over which various observations are likely to probe $`g(k,z)`$ over the next few years. Regardless of what detailed assumptions are made, it is clear from figures 1 and 7 that upcoming measurements of the CMB, the Ly$`\alpha `$ Forest, weak lensing, cluster abundance and galaxy clustering have the potential to probe both functions back to high redshift over several orders of magnitude in scale, allowing numerous consistency checks between different types of measurements. Measuring $`g(z)`$ for pressureless matter could provide an independent determination of $`H(z)`$ .

It may be possible to extract further constraints on dark energy/modified gravity from going to higher order in perturbation theory or studying the nonlinear regime. Such nonlinear issues were indeed the motivation for some recent suggestions for dark matter and modified gravity . Although many proposed quintessence fields do not cluster on small scales, early work in this direction  has discussed how temporal and perhaps spatial variations in Newton’s constant may be an observable signature.

In conclusion, the evidence for dark energy/modified gravity has emerged as one of the key puzzles of modern cosmology. Fortunately, there are at least half a dozen observational approaches that can directly probe the time-evolution of the effective mean density and clustering, so by going after this problem with the full arsenal, there is real hope that the cosmology community can resolve it over the next decade.

The author wishes to thank Arthur Kosowsky and Jim Peebles for provocative questions triggering this work, and Bernd Brügman, Angélica de Oliveira-Costa, Daniel Eisenstein, Gilles Esposito-Farese, Wayne Hu, Dragan Huterer, Arthur Kosowsky, David Polarski, Alexei Starobinski, Paul Steinhardt, Mike Turner and Robert Wald for helpful comments.

Support for this work was provided by NSF grant AST00-71213, NASA grant NAG5-9194 and the University of Pennsylvania Research Foundation.


format: text

text:
# The Mass-Temperature Relation of 22 Nearby Clusters

## 1 Introduction

In the past ten years, the ROSAT and ASCA observations of galaxy clusters have revealed that there exist tight correlations between the total gravitating mass of clusters ($`M_{tot}`$) and the X-ray luminosity ($`L_X`$) and the temperature ($`T_X`$) of the intra-cluster medium (ICM) (e.g. David et al., 1993; White, Jones & Forman 1997; Markevitch 1998; Wu, Xue & Fang, 1999; Horner, Mushotzky & Scharf 1999; etc.). Presence of such relations has been also motivated by both the theoretical studies and the $`N`$-body numerical simulations. A close comparison of the theoretical/simulated results and the observed data can therefore provide valuable constraints on the prevailing cosmological models and even on the nature of dark matter which dominates the mass distribution and the dynamical evolution of clusters. However, in addition to its dependence on the hydrostatic equilibrium hypothesis, an accurate determination of $`M_{tot}`$ often suffers from instrumental limitations, especially in those faint and distant clusters. This compares with the ICM temperature measurements which can be relatively easily performed through the X-ray spectroscopy in most clusters. Consequently, a well-established $`M_{\mathrm{tot}}T_X`$ relation can instead be used as a powerful cluster mass estimator. In particular, based on the $`M_{\mathrm{tot}}T_X`$ curve, one is able to link the observable quantities with the theoretical Press-Schechter mass function, which appears to be crucial for the determinations of the cosmological parameters using the cluster abundances at different redshifts.

Within the framework of pure gravitational infall (e.g., Lilje 1992), the cluster mass and the gas temperature should scale simply as $`M_{tot}T_X^\alpha `$, where the index $`\alpha (=3/2)`$ is almost unrelated to the different settings of the cosmological parameters. This self-similarity has been justified by cosmological $`N`$-body simulations (Evrard et al. 1996; Bryan & Norman 1998; Eke et al. 1998), and has also been tested by numerous observations. For example, Hjorth et al. (1998) have inferred the cluster mass from the weak and strong lensing for 8 clusters and confirmed such a tight mass-temperature relation. Neumann et al. (1999) analyzed the ROSAT PSPC and HRI data of 26 clusters and estimated the cluster mass by fitting the X-ray brightness profiles with the isothermal $`\beta `$-model, coupled with the hydrostatic equilibrium assumption. For a sub-sample of 15 targets whose temperature are available, they achieved a result consistent with $`M_{tot}T_X^{3/2}(1+z)^{3/2}`$.

However, it should be noticed that in the above observational analyses, only high-temperature ($`kT>45`$ keV) clusters are involved. As has been realized in many recent studies (e.g. Mohr, Mathiesen & Evrard 1999; Ponman, Cannon & Navarro 1999), the influence of energy feed-back into ICM in the low-temperature clusters/groups can become more significant than that in the high-temperature systems. Consequently, the self-similarity may break at the low-temperature end. Indeed, this has been confirmed in Nevalainen et al. (2000) using a small sample of 9 clusters (6 clusters at $`4`$ keV and 3 groups at $`1`$ keV). From the ASCA and the ROSAT temperature profiles, the authors have found that for the whole sample $`M_{tot}`$ scales as $`T_X^{1.79\pm 0.14}`$, which is apparently steeper than the theoretical prediction, whereas $`M_{tot}T_X^{3/2}`$ still holds if the low-temperature groups are excluded.

In the scenario of preheated ICM, such a similarity break can develop naturally in the low-temperature clusters/groups. A robust constraint on the break temperature ($`T_{\mathrm{break}}`$) can be obtained in principle by the observations of the $`L_X`$-$`T_X`$ or $`M_{vir}`$-$`T_X`$ relations, which will in turn yield valuable information about the heating processes in the early phase of cluster formation. In the recent analysis of Nevalainen et al. (2000), only a lower limit of 1 keV could be set to $`T_{break}`$ due to the temperature gap between $`14.5`$ keV in their sample. On the other hand, Ponman et al. (1999) evaluated the entropy of ICM in different systems and argued that $`T_{break}`$ could be as high as $`3`$ keV. Yet, this seems to be inconsistent with the finding of Horner et al. (1999) that the similarity of $`M_{\mathrm{tot}}T_X`$ may continue until $`T2`$ keV.

In this paper, we report our study of the $`M_{\mathrm{tot}}T_X`$ relation using a combined sample of 22 galaxy clusters observed with ROSAT and ASCA. The radial X-ray surface brightness profiles are derived and fitted with different models. Based on the best-fit parameters, the total gravitational masses scaled to the virial radius and the ICM entropies at the cluster central regions are calculated under the assumption of isothermality. Throughout this paper, we employ the MEKAL model (Mewe, Kaastra & Liedahl 1995) to calculate the thermal emission from the plasma with a mean metallic abundance of 0.4 solar, and set the Hubble constant to $`50h_{50}`$ km s-1 Mpc-1.

## 2 The Sample and Data Reduction

Our sample (Table 1) consists of 22 galaxy clusters, whose redshifts range from $`0.01`$ to $`0.09`$. At this distance range, the effect of cosmic evolution can be ignored. ROSAT PSPC data archived at the Institute of High Energy Physics (Beijing) are analyzed. In general, the observations are chosen in such a way that targets are located near the center of FOV, which makes it easier to study their diffuse X-ray emission as a whole. All the targets are sufficiently bright with their 2–10 keV flux being at least several $`10^{11}\mathrm{ergs}\mathrm{s}^1`$ cm-2. Thus, with sufficient exposure times, high signal-to-noise ratio can be achieved. We selected targets that all demonstrate relatively good circular symmetry in their X-ray morphology, and contain no distinct substructures. Such systems may have been undisturbed for rather a long time, for which the assumption of spherical symmetry is more reasonable.

We have searched in the literature for the emission weighted ASCA temperatures for our targets. Essentially, we choose the temperatures obtained when the central cluster regions are excluded. This is to avoid the complication due to the additional cool emission component (e.g. Buote 2000). When the error bars of the temperatures are asymmetric, we adopt the values of the larger sides. Among the targets, A262 has the lowest temperature at $`kT=2.15\pm 0.06`$ keV, and A2029 is the hottest one at about 10 keV. It appears that the sample almost covers the whole temperature range of clusters. The galaxy groups with even lower temperatures at $`1`$ keV are not included in this work.

## 3 The Models

We have extracted the radial X-ray surface brightness profiles for the targets by directly using their 0.5–2 keV count rates (Figure LABEL:fig:f1). Two approaches are adopted to analytically approximate the observed $`S_X(r)`$: the $`\chi ^2`$ fit via the conventional $`\beta `$ model and the $`\chi ^2`$ fit via the theoretically expected $`S_X(r)`$ tracing a priori given dark matter density profile. For the latter, we use the cusped NFW profile and the empirical Burkert profile with a definite core radius.

### 3.1 $`\beta `$ Model

If both galaxies and ICM are assumed to be isothermal and trace the common gravitational potential in a cluster, the ICM density $`n_g`$ can be represented by the so-called $`\beta `$ model when the galaxy distribution is assumed to follow the King profile (Cavaliere & Fusco-Femiano 1976):

$$n_g(R)=n_{g0}\left[1+(R/R_c)^2\right]^{3\beta /2}.$$
(1)
This density profile can also be obtained by inversion of the observed $`S_x(r)`$ in form of

$$S_X(r)=S_0\left[1+\left(\frac{r}{R_c}\right)^2\right]^{3\beta +\frac{1}{2}},$$
(2)
where $`R_c`$ is the core radius and $`\beta `$ is the slope parameter. Here, we use $`R`$ to represent the three-dimensional radius and $`r`$ for the two-dimensional radius, both measured from the cluster center. The total cluster mass within radius $`R`$ can be derived from the hydrostatic equilibrium assumption as

$$M_{\mathrm{tot}}(R)=\frac{kT_XR}{\mu m_\mathrm{p}G}\left\{\frac{\mathrm{dln}n_g(R)}{\mathrm{dln}R}+\frac{\mathrm{dln}T_X(R)}{\mathrm{dln}R}\right\}.$$
(3)
### 3.2 NFW profile

The second model we apply is the NFW profile, suggested by the high resolution $`N`$body simulations (Navarro, Frenk, & White 1995):

$$\rho (R)=\frac{\delta _c\rho _{crit}}{(R/R_s)(1+R/R_s)^2},$$
(4)
where $`\rho `$ is the halo density, $`R_s`$ is the scale radius, $`\rho _{crit}`$ is the critical density of the universe, and $`\delta _c`$ is the characteristic density contract, respectively. In 1998, Makino, Sasaki & Suto derived an analytic form of gas distribution from Eq(4) by assuming (a) spherical symmetry, (b) hydrostatic equilibrium, and (c) that the ICM is an ideal gas:

$$\rho _g(R)=\rho _{g0}\mathrm{exp}(27b/2)(1+R/R_s)^{27b/(2R/R_s)},$$
(5)
where $`b`$ is a dimensionless constant. With Eq(5), we can straightforwardly calculate the ICM X-ray emission and predict the corresponding X-ray surface brightness profile $`S_X(r)`$. Fitting this theoretically expected $`S_X(r)`$ to the X-ray observed one, we will be able to check the goodness of fit by the NFW profile and simultaneously fix the free parameters $`R_s`$ and $`\delta _c`$. Finally, we obtain the total cluster mass through

$$M_{\mathrm{tot}}(R)=4\pi \delta _c\rho _{\mathrm{crit}}R_s^3\left[\mathrm{ln}\left(1+\frac{R}{R_s}\right)\frac{R}{R+R_s}\right].$$
(6)
### 3.3 Burkert profile

Many recent observations, mainly based on the measurements of the rotation curves of dwarf galaxies, have shown that galactic halos may have flatter cores rather than the cusped density profiles such as the NFW profile. Indeed, such flat cores may develop if the dark matter particles are weakly self-interacting (e.g., Spergel & Steinhardt 1999; Burkert 2000). An empirical density form for such halos is proposed by Burkert (1995):

$$\rho (R)=\frac{\rho _0R_0^3}{(R+R_0)(R^2+R_0^2)},$$
(7)
in which $`\rho _0`$ and $`R_0`$ are the central density and core radius, respectively. Although Eq(7) was originally suggested for dwarf galaxies, it should also work for clusters, if the dark halos are self-similar. Indeed, the recent simulations of Yoshida et al. (2000) have shown that Eq(7) is applicable to clusters, provided that the halo particles are indeed self-interacting.

Following Makino et al. (1998), we derive the radial gas distribution as

$$n_g(R)=n_{g0}\mathrm{exp}\left\{k_0\mathrm{ln}\left[\frac{(R+R_0)^{\frac{2}{R}+\frac{2}{R_0}}(R^2+R_0^2)^{\frac{1}{R}\frac{1}{R_0}}}{R_0^{4/R}}\right]2k_0\mathrm{arctg}\left(\frac{R}{R_0}\right)\left(\frac{1}{R}+\frac{1}{R_0}\right)\right\},$$
(8)
where $`k_0`$ is a dimensionless constant. Then, in a similar way to the NFW profile, we can predict the resulting X-ray surface brightness profile as a result of thermal emission and compare it with X-ray observations. This would allow us to check the goodness of fit and work out the free parameters in the Burkert profile and hence, the total cluster mass.

## 4 Results

### 4.1 Goodness of fittings

We find that the goodness of fit of $`S_X(r)`$ is almost indistinguishable between the $`\beta `$ model and the NFW model (Table 2), but is apparently worse for the Burkert model for which $`\chi _r^2`$ are larger by a factor of 2-3. We demonstrate one example (A2199) in Fig. LABEL:fig:f2 by plotting the residuals of the fittings for the Burkert profile. It is apparent that the model provides a significant overestimate of $`S_X(r)`$ at about $`2^{}`$ and an underestimate at the adjacent regions. The situation is quite common in our fittings, implying that the Burkert model tends to give a core which is too flat. Consequently, the failure of the Burkert profile in the fittings of $`S_X(r)`$ essentially excludes this empirical model in our mass estimates of clusters below.

### 4.2 The Break of Mass- and Entropy-Temperature Relations

Based on the best-fits with the $`\beta `$ model and the NFW model, we calculate the total cluster mass within the virial radius $`r_{200}`$, which is defined as the radius where the mean cluster density is $`200`$ times of the critical density of the universe. As can be seen in Table 2, the values of $`r_{200}`$ are typically in the range 1.5–3.5 Mpc which are apparently larger than the observable radii of the clusters. Here we estimate $`r_{200}`$ for each cluster using the extrapolated mass profile. The uncertainties in mass are calculated by taking into account the corresponding measurement errors in both the temperature and the surface brightness profile. It is found that the resulted mass errors are always dominated by the errors of the surface brightness profile. This last point is crucial because the estimates of virial radius and total cluster mass depend on gas temperature through equation (3). An apparent correlation between $`M_{\mathrm{tot}}`$ and $`T_X`$ may arise if the errors of the gas temperature are relatively large. The mean error level for our cluster sample is about 20%, which is consistent with the predictions in the numerical studies. The resulting mass and temperature relation is plotted in Figure LABEL:fig:f3.

For the entire sample, we find that the best-fit mass-temperature relations are

$`M_{\mathrm{tot}}=5.75\pm 0.03\times 10^{13}h_{50}^1\left({\displaystyle \frac{T_X}{\mathrm{keV}}}\right)^{1.60\pm 0.04}M_{},`$ $`(\beta );`$ (9)
$`M_{\mathrm{tot}}=3.43\pm \times 10^{13}h_{50}^1\left({\displaystyle \frac{T_X}{\mathrm{keV}}}\right)^{1.81\pm 0.04}M_{},`$ $`(\mathrm{NFW}).`$ (10)
Both of the relations are steeper than the theoretical prediction of $`MT^{3/2}`$. Since such steepness has been suggested as a result of the pre-heating of ICM, which could be more remarkable in the low-temperature systems, we have excluded the low-temperature clusters one by one to examine beyond which break temperature the similarity is restored. We find that when those clusters at $`kT<3.5`$ keV are excluded, the $`M_{\mathrm{tot}}T_X`$ scaling become self-similar again with $`\alpha =1.54_{0.05}^{+0.06}`$ and $`1.59\pm 0.04`$ for $`\beta `$ and NFW model, respectively. This result agrees well with the one found by Ponmann et al. (1998), who discovered that in ICM at $`<3`$ keV the gas entropy is significantly higher than achievable through gravitational collapse alone. Recently, Tozzi & Normal (2000) used a generalized model to simulate the X-ray properties of the clusters. Their $`M_{\mathrm{tot}}T_X`$ in the low density, cold dark matter universe with an initial entropy $`K=0.4\pm 0.2\times 10^{34}`$ erg cm-2 g-5/3 is again in good accordance with our results, showing that the emission weighted temperature of the $`kT<2`$ keV clusters are about 25% higher than what is predicted by the self-similar scaling.

In order to confirm our result of $`T_{\mathrm{break}}`$ obtained via the $`M_{\mathrm{tot}}T_X`$ , we calculate the ICM entropy at $`0.1r_{200}`$ and illustrate the entropy-temperature distribution in Figure LABEL:fig:f4, where we use the entropy variable $`S=T_X/n_g^{2/3}`$. It appears that there is a clear tendency of the entropy increase for clusters at $`kT_X>4`$ keV, suggesting again that $`T_{\mathrm{break}}`$ should be between 3 and 4 keV.

## 5 Summaries

We have analyzed the archival data of ROSAT PSPC pointing observations of 22 galaxy clusters with gas temperatures ranging from 2 to 10 keV. The radial surface brightness profiles are extracted from the 0.5–2 keV count rates, and are fitted with the $`\beta `$ model, the NFW model, and the Burkert model, respectively.

We have found that the goodness of fit is indistinguishable between the $`\beta `$ model and the NFW model, but is significantly worse for the Burkert model. For the latter, the simple reason is that the Burkert model tends to present too flat a core. If the Burkert model can be used to describe the gravitational potential wells in both the galactic halos and cluster halos within the framework of the self-interacting scenario, we may conclude that the cross-sections of dark matter particles may vary with halo sizes. Yoshida et al. (2000) have recently arrived at the similar conclusion by analyzing their $`N`$-body simulations.

Based on the best-fit parameters of the $`\beta `$ model and the NFW model, we have calculated the cluster mass within the virial radius ($`r_{200}`$). An examination of the $`M_{\mathrm{tot}}T_X`$ shows that the $`M_{\mathrm{tot}}T_X`$ is steeper than the self-similarity prediction unless the clusters at $`<3.5`$ keV are excluded. A further study of the ICM entropy as a function of the temperature shows that an entropy floor is likely to present in the clusters cooler than $`4`$ keV. Therefore, we conclude that the break temperature of the self-similarity is at about 3–4 keV. We have noticed that based on a more complicated sample with spatially resolved temperature profiles, Finoguenov, Reiprich & Böhringer (2000) have arrived at a similar conclusion.

Finally, our $`M_{\mathrm{tot}}T_X`$ may be moderately affected if the intracluster gas deviates from the simple isothermality. Indeed, previous studies have arrived at conflicting results regarding the radial temperature gradients in clusters. For examples, Irwin & Bregman (2000) and Tamura et al. (2001) claimed an essentially flat temperature profile, while Markevitch et al. (1998) and Kaastra et al. (2000) reported a noticeable decline of gas temperature at large radii. However, according to the recent analysis of Finoguenov, Reiprich & Böhringer (2000), it is unlikely that the $`M_{\mathrm{tot}}T_X`$ can be significantly altered even if the temperature gradient is taken into account.

This research has made use of ROSAT archival data obtained from the Institute of High Energy Physics in Beijing, and is supported by the National Science Foundation of China, under Grant No. 19803003 and 19725311, and the Ministry of Science and Technology of China, under Grant NKBRSF G19990754.

Figure Captions


text:
# Classical information capacity of superdense coding

## A

Proof of Eq. (62). Suppose $`\rho _{AB}`$ is disentangled, then the density matrix may be written in the form $`\rho _{AB}=_ip_i\omega _A^i\omega _B^i`$, with $`_ip_i=1`$ and $`p_i>0`$, where the reduced density matrices $`\omega ^i`$ are all pure states. By the convexity of the expression $`S(\rho _B)S(\rho _{AB})`$ , we have,

$`S(\rho _B)S(\rho _{AB})`$ $``$ $`{\displaystyle \underset{i}{}}p_iS(\omega _B^i){\displaystyle \underset{i}{}}p_iS(\omega _A^i\omega _B^i)`$ (A1)
$`=`$ $`0`$ (A2)
and hence $`S(\rho _{AB})S(\rho _B)`$. Similarly for $`S(\rho _{AB})S(\rho _A)`$.


text:
# 1 Introduction

## 1 Introduction

Conventionally deep inelastic lepton-proton scattering is described in terms of scale-dependent parton distributions, $`a(x,\mu ^2)`$, where $`a=xg`$ or $`xq`$. These distributions correspond to the density of partons in the proton with longitudinal momentum fraction $`x`$, integrated over transverse momentum up to $`k_t=\mu `$. They satisfy DGLAP evolution in $`\mu ^2`$. The kinematic region $`k_t<\mu `$ gives the leading $`\mathrm{ln}\mu ^2`$ approximation to deep inelastic scattering.

For less inclusive processes it is, however, necessary to consider distributions unintegrated over the transverse momentum $`k_t`$ of the parton. The unintegrated distributions have the advantage that they exactly correspond to the quantity which enters the Feynman diagrams and therefore allow for the true kinematics of the process even at leading order (LO). These distributions $`f_a(x,k_t^2,\mu ^2)`$ depend on two hard scales111This property is hidden in the conventional distributions as $`k_t`$ is integrated up to the scale $`\mu `$.: $`k_t`$ and the scale $`\mu `$ of the probe. The scale $`\mu `$ plays a dual role. On the one hand it acts as the factorization scale, while on the other hand it controls the angular ordering of the partons emitted in the evolution .

Clearly it is desirable to also include $`\mathrm{ln}(1/x)`$ BFKL-type contributions in the evolution. Recall that both DGLAP and BFKL evolution are essentially equivalent to ordered evolution in the angles of the emitted partons222At LO we have strong ordering of the emission angles, $`\mathrm{}\theta _i\theta _{i+1}\mathrm{}`$; on the other hand if, at one step of the evolution $`\theta _i\theta _{i+1}`$, then this contribution is included inside the NLO splitting function.. In the DGLAP collinear approximation the angle increases due to the growth of $`k_t`$, while in BFKL the angle $`(\theta k_t/k_{\mathrm{}})`$ grows due to the decreasing longitudinal momentum fraction as we proceed along the emission chain from the proton. The factorization scale $`\mu `$ separates the partons associated with emission from different parts of the process, that is from the beam and target protons (in $`pp`$ collisions) and from the hard subprocess. For example it separates emissions from the beam (with polar angle $`\theta \stackrel{<}{}90^{}`$) from those from the target (with $`\theta \stackrel{>}{}90^{}`$), and from the intermediate partons from the hard subprocess. This separation was proved in Ref.  and originates from the destructive interference of the different emission amplitudes in the angular boundary regions. If the longitudinal momentum fraction is fixed by the hard subprocess, then the limits on the angles can be expressed in terms of a factorization scale $`\mu `$ which corresponds to the upper limit333The $`t`$-channel parton may have $`k_t`$ up to $`\mu /z`$, characteristic of BFKL effects, whereas for LO DGLAP the $`s`$ and $`t`$-channel partons are both limited by $`k_t<\mu `$. Of course, some $`k_t>\mu `$ contribution will arise from the NLO splitting functions. on the allowed values of the ($`s`$-channel) parton $`k_t`$.

Since the parton distributions depend on two scales we potentially have to deal with complicated (CCFM ) evolution equations for the $`f_a(x,k_t^2,\mu ^2)`$ functions. Of course it is possible to work with two-scale distributions, but this is much more complicated  and up to now has only proved practical with Monte Carlo generators . However, the evolution process is essentially controlled by one quantity, the emission angle, and on this basis we may expect to be able to obtain the distributions $`f_a(x,k_t^2,\mu ^2)`$ from single-scale evolution equations. Therefore it should be possible to follow an analytic approach where the physical assumptions are much more evident and where, in principle, NLO corrections can be included. Moreover, in practice, it is much easier to use the same unintegrated distributions to describe different hard processes and to perform global parton analyses.

The outline of this paper is as follows. The key observation is that the $`\mu `$ dependence of the unintegrated distributions enters at the last step of the evolution, and so we may use single-scale evolution equations. The procedure is first described in Section 2 in the case of pure DGLAP evolution, and then extended to include $`\mathrm{ln}(1/x)`$ effects in Section 3. In the latter case we use the solution of a single-scale equation which unifies DGLAP and BFKL evolution , and perform a final evolution step which brings in the dependence on the second scale. Ref.  also generated the two-scale unintegrated gluon from the same unified evolution equation, but with a different procedure444In this work we impose the angular ordering constraint in both the BFKL and DGLAP terms, and as a result do not have an exact equality between the integral up to $`\mu ^2`$ of the unintegrated distributions and the value of the integrated distribution. Ref.  takes the opposite approach; that is, exact equality with the integrated distribution is imposed and as a result angular ordering of the BFKL contribution is not complete. The difference is a NLO effect.. The unintegrated gluons obtained using the procedures described in Sections 2, 3 and Ref.  are compared in Section 4. In Section 5 we describe how the structure function $`F_2`$ is calculated from the unintegrated parton distributions, and in Section 6 we discuss the relationship between the unintegrated and integrated distributions. Finally in Section 7 we give our conclusions.

## 2 Unintegrated DGLAP partons

It is informative to review how unintegrated distributions $`f_a(x,k_t^2,\mu ^2)`$ may be calculated from the conventional (integrated) parton densities, $`a(x,\mu ^2)`$, in the case of pure DGLAP evolution. The procedure was explained in Ref. . We start from the DGLAP equation555For the $`ggg`$ splitting we have to insert a factor $`z^{}`$ in front of $`P_{gg}(z^{})`$ in the last integral of (1) to account for the identity of the produced gluons.

$$\frac{a(x,\mu ^2)}{\mathrm{ln}\mu ^2}=\frac{\alpha _S}{2\pi }\left[_x^{1\mathrm{\Delta }}P_{aa^{}}(z)a^{}(\frac{x}{z},\mu ^2)𝑑za(x,\mu ^2)\underset{a^{}}{}_0^{1\mathrm{\Delta }}P_{a^{}a}(z^{})𝑑z^{}\right]$$
(1)
where in the first term a sum over all possible parent partons $`a^{}`$ is implied. This first term on the right-hand-side describes the number of partons $`\delta a`$ emitted in the interval $`\mu ^2<k_t^2<\mu ^2+\delta \mu ^2`$. Such emission clearly changes the transverse momentum $`k_t`$ of the evolving parton. If we were to neglect the virtual contribution in (1), then the unintegrated parton density would be given simply by

$`f_a(x,k_t^2)`$ $`=`$ $`{\displaystyle \frac{a(x,\mu ^2)}{\mathrm{ln}\mu ^2}}|_{\mu ^2=k_t^2}`$ (2)
$`=`$ $`{\displaystyle \frac{\alpha _S(k_t^2)}{2\pi }}{\displaystyle _x^{1\mathrm{\Delta }}}P_{aa^{}}(z)a^{}({\displaystyle \frac{x}{z}},k_t^2)𝑑z.`$
The virtual contribution in (1) does not change the parton $`k_t`$ and may be resummed to give the survival probability $`T_a`$ that parton $`a`$ with transverse momentum $`k_t`$ remains untouched in the evolution up to the factorization scale. The survival probability is given by

$$T_a(k_t,\mu )=\mathrm{exp}\left(_{k_t^2}^{\mu ^2}\frac{\alpha _S(k_t^2)}{2\pi }\frac{dk_t^2}{k_t^2}\underset{a^{}}{}_0^{1\mathrm{\Delta }}P_{a^{}a}(z^{})𝑑z^{}\right),$$
(3)
à la Sudakov form factor. Thus the probability to find parton $`a`$ with transverse momentum $`k_t`$ (which initiates a hard subprocess with factorization scale $`\mu `$) is

$$f_a(x,k_t^2,\mu ^2)=T_a(k_t,\mu )\left[\frac{\alpha _S(k_t^2)}{2\pi }_x^{1\mathrm{\Delta }}P_{aa^{}}(z)a^{}(\frac{x}{z},k_t^2)𝑑z\right].$$
(4)
It is at this last step of the evolution that the unintegrated distribution becomes dependent on the two scales, $`k_t^2`$ and $`\mu ^2`$.

We now have to take care to specify the value of the infrared cut-off $`\mathrm{\Delta }`$, which is introduced to protect the $`1/(1z)`$ singularity in the splitting functions arising from soft gluon emission. In the original DGLAP equation, (1), which describes the evolution of the integrated distributions, this singularity is cancelled between the real emission and virtual contributions. However after the resummation of the virtual terms, the real soft gluon emission must be accounted for explicitly since it changes the $`k_t`$ of the parton. Thus we have to find the physically appropriate choice of the cut-off $`\mathrm{\Delta }`$ to provide the angular ordering of the gluon emissions666Although the splitting functions $`P_{gq}`$ and $`P_{qg}`$ are not singular at $`z=1`$ it is natural to use the same prescription for both the quark and the gluon distributions..

In Ref.  the cut-off was taken to be $`\mathrm{\Delta }=k_t/\mu `$. As a consequence the two-scale unintegrated distributions $`f_a(x,k_t^2,\mu ^2)`$ of  vanish for $`k_t>\mu `$, in accordance with the DGLAP strong ordering in $`k_t`$. However we can do better and impose the more correct angular ordering in the last step of the evolution. It was shown in Refs.  that this leads to a constraint on the scale $`\mu `$, namely

$$\mathrm{\Theta }(\theta \theta ^{})\mu >zk_t/(1z).$$
(5)
Thus the maximum allowed value of the integration variable $`z`$ is

$$z_{max}=\frac{\mu }{\mu +k_t}$$
(6)
and the corresponding cut-off $`\mathrm{\Delta }=k_t/(\mu +k_t)`$. Of course the same $`\mathrm{\Delta }`$ must be used both in the real emission integral in (4) and in the survival probability $`T`$ in (3).777In equation (3), $`\mathrm{\Delta }=k_t^{}/(\mu +k_t^{})`$ is the appropriate cut-off for $`z^{}`$. In fact we shall see that the imposition of angular ordering at the last step of the evolution leads to physically reasonable parton $`k_t`$ distributions which extend smoothly into the domain $`k_t>\mu `$.

## 3 Inclusion of $`\mathrm{𝐥𝐧}\mathbf{(}\mathrm{𝟏}\mathbf{/}𝒙\mathbf{)}`$ effects

We wish to generalize the above method to include the leading $`\mathrm{ln}(1/x)`$ contributions. Clearly there can be different forms of the ‘unified’ evolution equation summing up the leading DGLAP and BFKL logarithms, where the ambiguity is at the subleading level. The aim is to find a good prescription which is not too complicated, but which can account for all the physically relevant kinematic effects just at LO level. In other words we seek an equation which sums up the major part of the subleading corrections in a LO framework.

Let us consider, for the moment, just the gluon distribution. Recall that the unintegrated distribution $`f(x,k_t^2,\mu ^2)`$ depends on two scales. As in the pure DGLAP case of Section 2 we wish to work in terms of a single-scale evolution equation, and then to restore the scale $`\mu `$, and the full kinematics, at the last step of the evolution. This is illustrated schematically in Figure 1. For an analysis which incorporates BFKL effects, the appropriate single-scale distribution is the auxiliary function

$$h(x,k_t^2)=\frac{(xg(x,k_t^2))}{\mathrm{ln}k_t^2}.$$
(7)
Note that $`h(x,k_t^2)`$ is precisely the function which satisfies the BFKL equation in the low $`x`$ limit.

Both BFKL and DGLAP evolution correspond to angular ordering of the emission angles, and are single-scale equations. At LO, strong angular ordering automatically comes either from strong ordering in $`x`$ ($`z1`$) for BFKL or from strong $`k_t`$ ordering $`(k_t^2k_t^2)`$ for DGLAP. For the unintegrated gluon, we face a problem when $`k_t\mu `$ in the DGLAP framework, and similarly we have a problem when $`z1`$ for BFKL. Following the procedure of Ref. , we first neglect the subleading $`k_t\mu `$ and $`z1`$ effects to obtain and solve a unified BFKL/DGLAP equation for $`h(x,k_t^2)`$. Then, in the last step of the evolution, we take account of the precise kinematics so that the $`k_t\mu `$ and $`z1`$ domains are treated correctly.

The unified equation for $`h_g`$, which closely follows that presented in Ref. , takes the form

$`h_g(x,\mu ^2)`$ $`=`$ $`h_g^0(x,\mu ^2)+{\displaystyle \frac{\alpha _S(\mu ^2)}{2\pi }}{\displaystyle _0^1}dz{\displaystyle _{k_0^2}^{\mu ^2}}{\displaystyle \frac{dk_t^2}{k_t^2}}[\mathrm{\Theta }(zx)\overline{P}(z)h_g({\displaystyle \frac{x}{z}},k_t^2)`$ (8)
$`zP_{gg}(z)h_g(x,k_t^2)+\mathrm{\Theta }(zx)P_{gq}(z){\displaystyle }h_q({\displaystyle \frac{x}{z}},k_t^2)P_{gq}(z){\displaystyle }h_q(x,k_t^2)]`$
$`+{\displaystyle \frac{\alpha _S(\mu ^2)}{2\pi }}\mathrm{\hspace{0.25em}2}N_C{\displaystyle _x^1}{\displaystyle \frac{dz}{z}}{\displaystyle _{k_0^2}}{\displaystyle \frac{dq^2}{q^2}}\mathrm{\Theta }(k_t^2k_0^2)[\mathrm{\Theta }(\mu ^2zq^2){\displaystyle \frac{\mu ^2}{k_t^2}}h_g({\displaystyle \frac{x}{z}},k_t^2)`$
$`\mathrm{\Theta }(\mu ^2q^2)h_g({\displaystyle \frac{x}{z}},\mu ^2)],`$
where $`𝒌_t^{}=𝒌_t+(1z)𝒒`$, see Fig. 2.

We have introduced single-scale unintegrated quark auxiliary functions $`h_q(x,\mu ^2)`$ on the same footing as $`h_g`$,

$$h_q(x,\mu ^2)=\frac{(xq(x,\mu ^2))}{\mathrm{ln}\mu ^2},$$
(9)
and in (8) we sum over all $`2n_F`$ active flavours $`q`$ of quarks and antiquarks with $`m_q<\mu `$. The $`h_q`$ distributions satisfy the equation

$`h_q(x,\mu ^2)`$ $`=`$ $`h_q^0(x,\mu ^2)+{\displaystyle \frac{\alpha _S(\mu ^2)}{2\pi }}{\displaystyle _0^1}dz{\displaystyle _{k_0^2}^{\mu ^2}}{\displaystyle \frac{dk_t^2}{k_t^2}}\{\mathrm{\Theta }(zx)[P_{qg}(z)h_g({\displaystyle \frac{x}{z}},k_t^2)`$ (10)
$`+P_{qq}(z)h_q({\displaystyle \frac{x}{z}},k_t^2)]P_{qq}(z)h_q(x,k_t^2)\},`$
which is the usual DGLAP equation for quark evolution written in terms of $`h_q`$ of (9). The last term of (8) is the BFKL contribution which sums up all the leading $`(\alpha _S\mathrm{ln}1/x)^n`$ terms, while the remaining terms on the right-hand-side describe the conventional DGLAP evolution of the gluon distribution with respect to scale $`\mu ^2`$. To avoid double counting, we have excluded the singular part of the $`P_{gg}`$ splitting function in the real emission DGLAP term and used

$$\overline{P}(z)=P_{gg}(z)\frac{2N_C}{z}.$$
(11)
The $`2N_C/z`$ term is already included in the BFKL contribution to (8).

The driving terms, $`h^0`$, which describe the low $`k_t^2<k_0^2`$ domain are given by

$`h_g^0(x,\mu ^2)`$ $`=`$ $`{\displaystyle \frac{\alpha _S(\mu ^2)}{2\pi }}{\displaystyle _0^1}dz\{\mathrm{\Theta }(zx)[P_{gg}(z){\displaystyle \frac{x}{z}}g({\displaystyle \frac{x}{z}},k_0^2)`$
$`+P_{gq}(z){\displaystyle }{\displaystyle \frac{x}{z}}q({\displaystyle \frac{x}{z}},k_0^2)]zP_{gg}(z)xg(x,k_0^2)P_{gq}(z){\displaystyle }xq(x,k_0^2)\},`$
and similarly for $`h_q^0`$. The integrated input distributions $`a(x,k_0^2)`$ are not known and, as usual, must be determined from the data or from some non-perturbative QCD model.

We emphasize that in the unified BFKL/DGLAP equation we choose the scale $`\mu `$ to be $`k_t`$ for the DGLAP contribution, to be consistent with the BFKL term (which is independent of $`\mu `$ at LO). Recall that angular ordering leads to a redefinition of the scale

$$\mathrm{\Theta }(\theta \theta ^{})\mu >zk_t/(1z).$$

If this modified scale were to be adopted then it would be impossible to obtain a simple one-scale unified evolution equation throughout the whole $`x,k_t^2`$ domain. On the other hand, within the LO framework, we may omit the $`z`$ dependence in the scale, so that the DGLAP part of the evolution becomes ordered in transverse momenta. We stress again that it is sufficient to implement the precise constraints coming from angular ordering (relevant to the $`z0`$ and 1 domains for the BFKL and DGLAP terms respectively) at the last step of the evolution.

An advantage of the single-scale unified BFKL/DGLAP equation is that it is straightforward to incorporate a major (kinematical) part of the subleading order $`\mathrm{ln}(1/x)`$ (BFKL) effects888The large NLO BFKL corrections, which have recently been computed , appeared to have put the application of the BFKL framework into question. However a major part of the corrections is kinematic in origin and, when summed to all orders  using the theta function $`\mathrm{\Theta }(\mu ^2zq^2)`$, brings the BFKL approach back under control, see also . by imposing a consistency condition to ensure that the virtuality of the exchanged gluon is dominated by its transverse momentum squared . This is achieved by the inclusion of the theta function $`\mathrm{\Theta }(\mu ^2zq^2)`$ in the real emission contribution shown in the last term of (8). Note that other subleading effects arising from using the complete DGLAP splitting function and running $`\alpha _S`$ are automatically included in the unified equation.

We see that the evolution equations for the auxiliary distributions $`h_a(x,k_t^2)`$ depend on the single scale $`k_t^2`$. The dependence of the unintegrated distributions $`f_a(x,k_t^2,\mu ^2)`$ on the second scale $`\mu `$ will enter when we consider the last step of the evolution. It is sufficient to ensure that the final emitted parton explicitly satisfies the requirements of angular ordering. The angular ordering conditions of the previous steps of the evolution are automatically ensured at LO by virtue of either the strong ordering in $`k_t`$ (in the DGLAP part) or the strong ordering in $`z`$ (in the BFKL part).

To ensure angular ordering in the last step of the evolution, we note that $`z`$ is limited by (5). This condition implies $`z<\mu /(\mu +k_t)`$, and so we take this as the upper limit of the $`z`$ integration in (13) and (14) below. Thus the number of gluons produced at the last step (with transverse momentum $`k_t`$ which initiate a hard subprocess with factorization scale $`\mu `$) is999The low $`k_t<k_0`$ domain of the integrals in (13), (14) should be treated as the driving terms $`h^0`$ in (3).

$`f_g(x,k_t^2,\mu ^2)`$ $`=`$ $`T_g(k_t,\mu ){\displaystyle \frac{\alpha _S(k_t^2)}{2\pi }}\{{\displaystyle _x^{\mu /(\mu +k_t)}}dz{\displaystyle ^{k_t^2}}{\displaystyle \frac{dk_t^2}{k_t^2}}[\overline{P}(z)h_g({\displaystyle \frac{x}{z}},k_t^2)`$ (13)
$`+P_{gq}(z){\displaystyle }h_q({\displaystyle \frac{x}{z}},k_t^2)]+\mathrm{\hspace{0.25em}2}N_C{\displaystyle }_x^{\mu /(\mu +k_t)}{\displaystyle \frac{dz}{z}}{\displaystyle }{\displaystyle \frac{d^2q}{\pi q^2}}`$
$`[{\displaystyle \frac{k_t^2}{k_t^2}}h_g({\displaystyle \frac{x}{z}},k_t^2)\mathrm{\Theta }(k_t^2q^2)h_g({\displaystyle \frac{x}{z}},k_t^2)]\},`$
and the number of a particular quark species $`q`$, produced at the last step, is

$`f_q(x,k_t^2,\mu ^2)`$ $`=`$ $`T_q(k_t,\mu ){\displaystyle \frac{\alpha _S(k_t^2)}{2\pi }}{\displaystyle _x^{\mu /(\mu +k_t)}}dz{\displaystyle ^{k_t^2}}{\displaystyle \frac{dk_t^2}{k_t^2}}[P_{qg}(z)h_g({\displaystyle \frac{x}{z}},k_t^2)`$ (14)
$`+P_{qq}(z)h_q({\displaystyle \frac{x}{z}},k_t^2)].`$
This last step of the evolution is shown schematically in Fig. 1. The unintegrated distributions $`f_a(x,k_t^2,\mu ^2)`$ represent the number of partons $`\delta a`$ emitted in the $`\delta \mathrm{ln}k_t^2`$ interval from $`k_t^2`$ to $`k_t^2+\delta k_t^2`$, and include the factors $`T_a(k_t,\mu )`$ of (3), which give the probabilities that partons $`a=g,q`$ with transverse momentum $`k_t`$ remain untouched in the DGLAP evolution up to the factorization (probe) scale $`\mu `$. These survival probabilities $`T_g`$ and $`T_q`$ resum the virtual DGLAP contributions occurring in (8) and (10) respectively.

It is important to note that the function $`h_g(x,k_t^2)`$ already includes the leading $`\mathrm{ln}(1/x)`$ virtual corrections, which have the effect of reggeizing the exchanged gluon; thus $`k_t`$ is the total momentum transferred via the Regge gluon trajectory. Moreover, the distributions $`f_a(x,k_t^2,\mu ^2)`$, evolved in the final step from the auxiliary functions $`h_a(x,k_t^2)`$, also incorporate the virtual DGLAP contributions, via the survival probabilities $`T_a(k_t,\mu )`$. The final expressions, (13) and (14), are thus more symmetric in that all the LO virtual corrections are included. That is, from a Feynman diagram viewpoint, the function $`f_a(x,k_t^2,\mu ^2)`$ corresponds to the propagator of a $`t`$-channel parton of transverse momentum $`k_t`$, initiating a hard subprocess at scale $`\mu `$, in which all the LO virtual corrections to the parton distribution have been taken into account.

## 4 The unintegrated gluon distribution

We have described how to obtain the (two-scale) unintegrated parton distributions $`f(x,k_t^2,\mu ^2)`$ from the solution $`h(x,k_t^2)`$ of a one-scale equation which unifies DGLAP and BFKL evolution. The link is (13), (14), which represent the last step of the evolution. Only at this stage does the scale $`\mu `$ of the subprocess, initiated by $`f(x,k_t^2,\mu ^2)`$, enter.

It is informative to compare the unintegrated gluon distribution $`f_g(x,k_t^2,\mu ^2)`$ with the behaviour of the auxiliary function $`h_g(x,k_t^2)`$. Consider first the pure BFKL limit, in which DGLAP evolution is neglected, that is $`\alpha _S\mathrm{ln}\mu ^21`$, but $`\alpha _S\mathrm{ln}1/x\stackrel{>}{}1`$. Then, in (13), the survival probability $`T_g(k_t,\mu )=1`$ and there are no DGLAP contributions. Thus we obtain

$$f_g(x,k_t^2,\mu ^2)=h_g(x,k_t^2)=\frac{(xg(x,k_t^2))}{\mathrm{ln}k_t^2}.$$
(15)
This is an expression which is frequently used at low $`x`$.

When $`x`$ is sufficiently large, the derivative

$$h_a=\frac{(a(x,k_t^2))}{\mathrm{ln}k_t^2}$$
(16)
becomes negative, even for $`k_t\stackrel{<}{}\mu `$. The reason is that the negative virtual DGLAP term exceeds the real emission DGLAP contribution, which is suppressed by the large lower limit $`z>x`$ in (8), (10), (13) and (14). After the virtual contributions are resummed into the $`T_a`$ factors the unintegrated parton distributions remain positive everywhere.

We already stated in Section 2 that imposing angular ordering on the last step of “DGLAP” evolution caused the distribution to extend smoothly into the $`k_t>\mu `$ domain. The gluons also populate this domain due to BFKL diffusion in $`\mathrm{ln}k_t^2`$. This raises the question of how much enhancement we find in this domain due to the inclusion of the BFKL contributions. To investigate this point, we compute the unintegrated gluon distribution $`f_g(x,k_t^2,\mu ^2)`$ in two different ways.

* In the first approach, the function $`h_g(x,k^2)`$ of Ref.  is used as the auxiliary function to drive the last-step evolution in (13) and (14). This incorporates essentially maximal BFKL effects.101010In comparison with (10), the evolution of $`h_q`$ in Ref.  contains some resummation of the BFKL-like leading $`(\alpha _S\mathrm{ln}1/x)^n`$ terms.
* In the second case, we use unintegrated “DGLAP” partons to drive the two-scale unintegrated $`f_a`$ via (4). This approach is essentially pure DGLAP, but with the crucial modification that the cut-off in (3) and (4) is motivated by *angular ordering*. Here we use the MRST99  set of partons as input.

Ideally, we should refit to the deep inelastic and related scattering data and perform global *unintegrated* parton analyses111111In Section 5, we describe the theoretical calculation of $`F_2`$, the most important observable for constraining parton sets, from the two-scale unintegrated $`f_g`$ and $`f_q`$. in terms of $`f_a(x,k_t^2,\mu ^2)`$, since now we have added an extra last step to the evolution. However, such global analyses are beyond the scope of the present paper.

The unintegrated gluon distribution $`f_g(x,k_t^2,\mu ^2)`$, obtained by the two alternative procedures (a) and (b), is shown in Fig. 3 for $`\mu =10\mathrm{GeV}`$. We can compare the continuous curves of approach (a) directly with the previous unintegrated gluon distribution (dashed curves) calculated in Ref. , which also used the same auxiliary function $`h_g(x,k_t^2)`$ of  as input. For $`k_t<\mu `$ the new results tend to lie above the previous determination , while for $`k_t>\mu `$ they are increasingly smaller as $`k_t`$ increases. In the present work, the decrease at large $`k_t`$ arises from the restriction $`z<z_{\mathrm{max}}=\mu /(\mu +k_t)`$, whereas the earlier calculation, based on eq. (23) of Ref. , omitted the limit $`z<z_{\mathrm{max}}`$ on the BFKL contribution, which causes $`f_g`$ to increase in the large $`k_t`$ domain at small values of $`x`$.

Finally we compare the DGLAP-like unintegrated gluon, obtained in approach (b), with that of approach (a) which embodied BFKL evolution. That is we compare the dotted with the continuous curves of Fig. 3. It is interesting to note that, with the more precise angular-ordered (CCFM ) cut-off from (5), $`\mathrm{\Delta }=k_t/(\mu +k_t)`$, the unintegrated gluon distribution (b) generated from (4) using the pure DGLAP MRST99  partons turns out to be very similar to the result (a) obtained using the auxiliary function $`h_g`$. This coincidence can be explained by the facts (i) that both the DGLAP partons and the analysis of Ref.  (which yielded $`h_g`$) fit the deep inelastic data well, and (ii) that both these input are used with the same angular-ordered constraint, (5). We draw the conclusion that the role of angular ordering in the last step of evolution is particularly important, even more so than BFKL effects in the HERA domain.

## 5 $`𝑭_\mathrm{𝟐}`$ calculated from the unintegrated partons

To check the reliability of our unintegrated parton distributions, and to demonstrate how to use these distributions in calculations of observables, we compute the deep inelastic structure function $`F_2`$. We wish to treat the unintegrated gluons and unintegrated quarks on an equal footing as input to the subprocess cross sections, so we explicitly separate gluon and (direct) quark contributions to $`F_2`$.

The gluon contributes to $`F_2`$ via the quark box and crossed-box diagrams of Fig. 4. These generate, via the $`gq\overline{q}`$ splitting, a sea quark contribution $`S_q`$ to $`F_2`$ of the form

$$F_2^{gq\overline{q}}(x,Q^2)=\underset{q}{}e_q^2S_q(x,Q^2),$$
(17)
with

$`S_q(x,Q^2)={\displaystyle \frac{Q^2}{4\pi ^2}}{\displaystyle \frac{dk_t^2}{k_t^4}_0^1𝑑\beta d^2\kappa _t\alpha _S(\mu ^2)f_g(\frac{x}{z},k_t^2,\mu ^2)\mathrm{\Theta }\left(1\frac{x}{z}\right)}`$
(18)
$`\left\{\left[\beta ^2+(1\beta )^2\right]\left({\displaystyle \frac{𝜿_t}{D_1}}{\displaystyle \frac{𝜿_t𝒌_t}{D_2}}\right)^2+\left[m_q^2+4Q^2\beta ^2(1\beta )^2\right]\left({\displaystyle \frac{1}{D_1}}{\displaystyle \frac{1}{D_2}}\right)^2\right\}.`$
The denominator factors are

$`D_1`$ $`=`$ $`\kappa _t^2+\beta (1\beta )Q^2+m_q^2`$
$`D_2`$ $`=`$ $`(𝜿_t𝒌_t)^2+\beta (1\beta )Q^2+m_q^2.`$
We may exploit the symmetry of the integrand in (5) under $`𝜿_t𝜿_t𝒌_t`$ and $`\beta 1\beta `$ to rewrite $`\{\mathrm{}\}`$ as

$$2\left\{\left[\beta ^2+(1\beta )^2\right]\left(\frac{\kappa _t^2}{D_1^2}\frac{\kappa _t^2𝜿_t𝒌_t}{D_1D_2}\right)+\left[m_q^2+\mathrm{\hspace{0.25em}4}Q^2\beta ^2(1\beta )^2\right]\left(\frac{1}{D_1^2}\frac{1}{D_1D_2}\right)\right\}.$$
(20)
The summation in (17) is over massless $`u,d,s`$ quarks and a $`c`$ quark of mass $`m_c=1.4\mathrm{GeV}`$; there is no need to sum $`\overline{u}`$, $`\overline{d}`$, $`\overline{s}`$, $`\overline{c}`$ in addition, because as (5) is written, $`S_u`$, say, is the contribution of a gluon via a $`u`$ quark box of any momentum. The variable $`\beta `$ is the light-cone fraction of the photon momentum carried by the internal quark. The variable $`z`$ is the ratio of Bjorken $`x`$ and the fraction of the proton momentum carried by the gluon. It is specified by the relation

$`{\displaystyle \frac{1}{z}}`$ $`=`$ $`1+{\displaystyle \frac{(𝜿_t(1\beta )𝒌_t)^2+m_q^2}{\beta (1\beta )Q^2}}+{\displaystyle \frac{k_t^2}{Q^2}}`$
$`=`$ $`1+{\displaystyle \frac{\kappa _t^2+m_q^2}{(1\beta )Q^2}}+{\displaystyle \frac{k_t^2+\kappa _t^22𝜿_t𝒌_t+m_q^2}{\beta Q^2}},`$
which is obtained by requiring the outgoing quarks to be on-shell. Following Ref. , we choose the scale $`\mu `$ which controls the unintegrated gluon distribution and the QCD coupling $`\alpha _S`$ to be

$$\mu ^2=k_t^2+\kappa _t^2+m_q^2.$$
(22)
Care is needed in separating this calculation into perturbative and non-perturbative regions. We impose a cut-off $`k_t>k_0`$ for a legitimate perturbative calculation of (5). The smallest cut-off we can choose is the minimum (that is, initial) scale of the function $`h_g(x,k_t^2)`$ from which our two-scale distributions derive. Thus $`k_0`$ is of order 1 GeV. For the contribution from the region of $`k_t<k_0`$ we approximate

$$_0^{k_0^2}\frac{dk_t^2}{k_t^2}f_g(x,k_t^2,\mu ^2)\left[\frac{\mathrm{remainder}}{k_t^2}\right]xg(x,k_0^2)T_g(k_0,\mu )\left[\right]_{k_t=a},$$
(23)
where $`a`$ can be taken to be any value representative of the interval $`(0,k_0)`$. The dependence on the choice of $`a`$ is numerically unimportant.

Now we have to add the direct quark contributions to $`F_2`$, which come from the unintegrated quark distributions $`f_q(x,k_t^2,\mu ^2)`$. If a quark, initially with $`x/z`$ and perturbative transverse momentum $`k_t^{}>k_0`$, splits to a radiated gluon and a ‘final’121212By ‘final’ we mean the struck quark just after the $`P_{qq}`$ splitting. quark with Bjorken $`x`$ and transverse momentum $`\kappa _t`$, then this final quark can couple to the photon and contribute to $`F_2`$ as:

$`F_2^{q(\mathrm{pert})}(x,Q^2)={\displaystyle \underset{q}{}}e_q^2{\displaystyle _{k_0^2}^{Q^2}}{\displaystyle \frac{d\kappa _t^2}{\kappa _t^2}}{\displaystyle \frac{\alpha _S(\kappa _t^2)}{2\pi }}{\displaystyle _{k_0^2}^{\kappa _t^2}}{\displaystyle \frac{dk_t^2}{k_t^2}}{\displaystyle _x^{Q/(Q+k_t)}}𝑑z`$ (24)
$`\left[f_q({\displaystyle \frac{x}{z}},k_t^2,Q^2)+f_{\overline{q}}({\displaystyle \frac{x}{z}},k_t^2,Q^2)\right]P_{qq}(z).`$
where here we have written the antiquark contribution explicitly. As in (14), the upper limit of the $`z`$ integration reflects the angular-ordered constraint of (5) during the quark evolution.

Again we need to account for the non-perturbative domain $`k_t<k_0`$. The initial (integrated) quark distribution $`xq(x,k_0^2)`$ drives our final contribution. Physically the only remaining diagrams that we have not included are those in which a quark (or antiquark) from this initial distribution does not experience real splitting in the perturbative domain, but interacts unchanged with the photon at scale $`Q`$. Hence we write a Sudakov-like factor $`T_q(k_0,Q)`$ to represent the probability of evolution from $`k_0`$ to $`Q`$ without radiation.

$$F_2^{q(\mathrm{non}\mathrm{pert})}(x,Q^2)=\underset{q}{}e_q^2\left(xq(x,k_0^2)+x\overline{q}(x,k_0^2)\right)T_q(k_0,Q).$$
(25)
To avoid double counting, it is important to put a lower limit on $`\kappa _t`$ in both (5) and (24), by enforcing $`\mathrm{\Theta }(\kappa _t^2k_0^2)`$ in the integrations. Without this lower limit on the final transverse momentum $`\kappa _t`$, equations (5) and (24) would partially include low transverse momentum $`\kappa _t`$ quark contributions which are best incorporated in (25), whether they originate from partons with $`k_t>k_0`$ or not.

The structure function $`F_2(x,Q^2)`$ is given by the sum of the gluon-initiated contribution (17), and the quark terms, (24) and (25). In Fig. 5 sample results131313We show in Fig. 5 results obtained using the unintegrated distributions evaluated from (4) with the MRST99 partons , that is, version (b) of $`f_g`$ and $`f_q`$ as discussed in Section 4., shown by the continuous curves, are compared with deep-inelastic structure function data. The gluon and quark components are shown by the dashed and dotted curves respectively. As expected the dominant contribution at small $`x`$ comes from the unintegrated gluon via the quark box and crossed-box contributions, whereas at large $`x`$ the quark terms dominate.

We emphasize that, in the present work, the curves for $`F_2`$ are not the result of a fit to the structure function data. Rather they have been obtained by using single-scale functions, originally fitted more directly to $`F_2`$ data, as plausible input to our ‘last-step’ evolution procedure, which generates two-scale unintegrated distributions $`f_a(x,k_t^2,\mu ^2)`$. Then we use the unintegrated distributions to compute $`F_2`$ via (17), (24) and (25). For the inclusive observable $`F_2`$, we would expect that the insertion of this extra evolution step would not appreciably disturb the description, since it essentially redistributes the distributions in $`k_t`$ space. We see from Fig. 5 that indeed this is the case.

## 6 Relation of $`𝒇`$ to integrated partons

It is important to scrutinise the relationship between the new unintegrated partons $`f_a(x,k_t^2,\mu ^2)`$ and the conventional integrated parton distributions $`a(x,\mu ^2)`$, as obtained in global analyses such as . First we emphasize that we may use either the integrated distributions or the unintegrated distributions to describe both inclusive (like $`F_2`$) and exclusive processes. The framework based on the unintegrated distributions is a bit more complicated. However it accounts for the precise kinematics of the process and an important part of the virtual loop corrections, via the survival factor $`T`$, even at LO. On the other hand, if we work with integrated partons we have to include NLO (and sometimes NNLO) contributions to account for these effects. These differences appear to cause a discrepancy between the integrated and unintegrated approaches. As we explain below, this is to be expected since it arises from simplifications of the LO formalism due to the neglect of terms which are moved into the NLO contribution.

An important equation, sometimes cited as the defining property of unintegrated partons , is

$$a(x,\mu ^2)=^{\mu ^2}\frac{dk_t^2}{k_t^2}f_a(x,k_t^2,\mu ^2),$$
(26)
where $`a`$ represents $`xg`$ or $`xq`$. This is in fact the first equation of Ref. . In the BFKL limit, the $`\mu `$ dependence of $`f`$ vanishes and we have $`f_g(x,k_t^2,\mu ^2)h_g(x,k_t^2)`$ as in (15). In this case, (26) is clearly satisfied. However, in general the situation is complicated by the two separate momentum scales $`k_t`$ and $`\mu `$. The unintegrated partons $`f_g`$ of Ref.  were explicitly constructed to have the property (26), in the sense that the integral of $`f_g`$ over the transverse momentum up to the scale $`\mu `$ would be the same as the integral of the input auxiliary function $`h_g(x,k_t^2)`$ up to the same scale. In contrast, numerical integration over $`k_t`$ of the new unintegrated partons $`f_g`$ and $`f_q`$ presented in this paper (both versions (a) and (b) of Section 4) shows that (26) is only approximately true.141414Note that we cannot compute (26) as it is written, because we can only define the unintegrated function in the regime of perturbative $`k_t>k_0`$. The comparison that is made is between the integral from $`k_0^2`$ to $`\mu ^2`$ and the quantity $`a(x,\mu ^2)a(x,k_0^2)`$. We typically find a discrepancy of order 25% between the right-hand side of (26) and the single-scale distribution that has been used to generate $`f_a`$.

In order to eliminate the discrepancy we may adjust the upper limit $`\mu ^2`$ of the integral in (26) to $`c^2\mu ^2`$. The introduction of $`c`$ is equivalent to a NLO correction for the integrated partons. Typically in the “DGLAP” case (approach (b) of Section 4) we require $`c=0.60.8`$ to reproduce151515This identification is only approximate since the last step of the evolution is not included in the comparison. the original MRST integrated gluon in the domain $`\mu =510\mathrm{GeV}`$ and $`x\stackrel{<}{}0.01`$ (and $`c0.4`$ for approach (a), which embodies BFKL effects). A value $`c<1`$ compensates for the over-large virtual loop correction included in the integrated DGLAP partons on account of the absence of the cut-off $`\mathrm{\Delta }`$ in the integral over $`z^{}`$. At smaller values of $`x`$ the corrections due to the cut-off are mainly cancelled between the virtual and real DGLAP contributions. As $`x`$ increases, the virtual contribution (the second term on the right-hand-side of (1)) increasingly dominates and we have to choose smaller values of $`c`$. Eventually in the domain $`x>0.1`$, $`\mu 10\mathrm{GeV}`$ the main contribution comes from the input and to compensate the absence of the $`z`$ cut-off in the conventional form of the DGLAP equation we must change the input itself.

To summarize, the discrepancy between the integral (26) of the unintegrated parton function and the original integrated distribution is not a cause for concern. Conceptually, there are two different roles for single-scale distributions in the description of data for inclusive observables (where partonic transverse momentum is integrated out). The first role is the traditional one, in the framework of collinear factorization, whereby integrated parton distribution functions are fitted directly to the data. The second role is demonstrated in this paper (following  and ), where we use single-scale functions as *input* to the last-step procedure, see for example (4). We have been forced to introduce a new formalism for calculating $`F_2(x,Q^2)`$, in which unintegrated parton functions are understood to be the fundamental objects; we have emphasized the need to perform a new global fit to data in terms of the new functions $`f_a`$. After this, we do not expect the input single-scale function $`a`$ on the left-hand side of (26) to equal the integral of $`f_a`$ up to $`\mu ^2`$, since $`a`$ itself is not fitted directly to the data, but rather is used as input for the last step of the evolution, which embodies a crucial angular-ordering constraint unique to this last step. Thus our single-scale or ‘auxiliary’ function is not a traditional parton distribution function, but simply an intermediate function.

## 7 Conclusions

Parton distributions, $`f_a(x,k_t^2,\mu ^2)`$ unintegrated over the parton $`k_t`$ are the basic quantities for describing processes initiated by hadrons. An essential ingredient in this description is the existence of the $`k_t`$ factorization theorem . The unintegrated distributions depend on two hard scales — the transverse momentum $`k_t`$ and the factorization scale $`\mu `$. The scale $`\mu `$ drives the angular ordering during the evolution which arises from the coherence of the gluon emissions.

Here we develop a new formalism to determine the unintegrated parton distributions, $`f_a(x,k_t^2,\mu ^2)`$, which embodies both the leading $`\mathrm{ln}Q^2`$ (DGLAP) and $`\mathrm{ln}1/x`$ (BFKL) effects, as well as including a major part of the sub-leading $`\mathrm{ln}1/x`$ contributions. An important observation is that, at leading order, the two-scale functions $`f_a(x,k_t^2,\mu ^2)`$ may be calculated from auxiliary functions $`h_a(x,k_t^2)`$ which satisfy single-scale evolution equations, since the angular ordering restrictions, controlled by $`\mu `$, become important only at the last step of the evolution. The equation for the auxiliary function $`h_g(x,k_t^2)`$ was formulated, and the distributions fitted to the data, in Ref. . These single-scale equations were also devised to include all the leading $`\alpha _S\mathrm{ln}Q^2`$ and $`\alpha _S\mathrm{ln}1/x`$ contributions, and a major part of the sub-leading $`\mathrm{ln}1/x`$ effects.

In other words, the ‘unified’ evolution equations for $`h`$ must be supplemented by a final evolution step in which the $`\mu `$ dependence of the unintegrated parton distributions enters via the angular-ordering constraint. The situation is summarised diagrammatically in Fig. 1. The procedure offers a considerable simplification in the determination of physically realistic unintegrated parton distributions $`f_a(x,k_t^2,\mu ^2)`$ over the full $`x,\mu ^2`$ perturbative domain, including the true kinematics even at leading order. As expected, the gluon and sea quark distributions extend into the $`k_t>\mu `$ region more and more as $`x`$ decreases. We have compared the new unintegrated distributions with those given by previous prescriptions . As compared to Ref. , the new formalism gives a consistent treatment of angular ordering, which leads to the imposition of the $`z`$ integration limit $`z<\mu /(\mu +k_t)`$. As a consequence the distributions $`f_a(x,k_t^2,\mu ^2)`$ decrease faster than those of  for large $`k_t`$, particularly at small $`x`$. An interesting result is that the unintegrated distributions obtained via $`h_a(x,k_t^2)`$ of  are not very different from those obtained via (4) using conventional DGLAP partons — compare the continuous and dotted curves in Fig. 3. It thus appears that the imposition of the angular-ordering constraint is more important than the BFKL effects. This observation has the practical consequence that reasonably accurate predictions for observables can be made using the much simpler, though less complete, prescription of (4).

Finally, we used the new unintegrated distributions to calculate the deep-inelastic structure function $`F_2`$. We also showed the gluon-initiated and quark contributions separately, which, as expected, dominate at small $`x`$ and large $`x`$ respectively. Recall, from , that the rise of the gluon at small $`x`$ is driven by perturbative QCD, which was assumed to have a non-perturbative input which is ‘flat’ in $`x`$. We emphasize that we have not fitted to the deep inelastic data. Nevertheless, Fig. 5 shows that the existing distributions give an adequate description, and therefore they may be used to evaluate other hard processes, such as $`b\overline{b}`$ and large $`q_t`$ prompt photon production in high energy $`p\overline{p}`$ (or $`pp`$) collisions. It is important to use unintegrated distributions for such exclusive reactions.

## Acknowledgements

We thank Jan Kwiecinski and Anna Stasto for many valuable discussions during the course of this research. The work was supported by the UK Particle Physics and Astronomy Research Council (PPARC), and also the Russian Fund for Fundamental Research (98-02-17629). This work was also supported by the EU Framework TMR programme, contract FMRX-CT98- 0194 (DG 12-MIHT).