text:
# Combining WMAP and SDSS Quasar Data on Reionization Constrains Cosmological Parameters and the Star Formation Efficiency

## 1 Introduction

The primieval spectrum of cosmic density fluctuations on large scales is determinable with great precision from cosmic background radiation measurements (Page et al., 2003). In order to determine the total relevant spectrum, information on small scales is also necessary; these become nonlinear at early times, so information concerning reionization in the epochs $`20z6`$ provide vital clues.

Recent observations of high-redshift quasars have provided the first observational signatures of the epoch of reionization. Spectra of quasars at redshift $`z6`$ indicate that the universe was almost fully ionized up to $`z6`$, since even a small neutral fraction in the intergalactic medium (IGM) would have led to complete absorption of a quasar’s continuum radiation. However, the first absorption spectra of quasars at higher redshift indicate that the abundance of neutral hydrogen increases significantly for $`z6`$.

Reionization has generally been assumed to be caused by ionizing photons created in early generations of stars and/or quasars. Given this premise, hydrodynamic simulations as well as semi-analytic calculations seem to indicate that the process of reionization should occur several distinct stages. First, cosmological gas falls into deep enough potential wells (caused by dark matter halos), so that they can cool and collapse to high enough densities to produce stars and/or quasars. The UV photons produced in this process ionizes the local surroundings, first within the halo itself and then outside the halo creating cosmological H II regions. This is generally referred to as the “pre-overlap” stage. As the abundance of these H II regions increases (due to additional “galaxy” formation), they eventually start to “overlap” so that gas in the IGM becomes exposed to multiple sources of ionizing radiation. After this “overlap” stage, the IGM becomes optically thin except for inside self-shielded, high-density clouds (those without ionizing sources).

In this picture, the reionization history of the universe depends on both the growth of density perturbations and the efficiency of star/quasar formation. The former is a complex function of the standard cosmological parameters — the density of the universe $`\mathrm{\Omega }_0`$, the baryon abundance $`\mathrm{\Omega }_b`$, the expansion rate $`H_0`$, and the mass power spectral index $`n`$ and normalization $`\sigma _8`$. The latter can be calculated with some precision based on atomic physics (e.g., cooling rates), but ultimately depends on some unknown parameters relating to the efficiency of turning mass into UV photons that can escape into the IGM.

To this overall picture has been added the recent observations by the Wilkinson Microwave Anisotropy Probe (WMAP) satellite (Bennett et al., 2003). In particular, the high values of the electron optical depth to last scattering $`\tau _{\mathrm{es}}=0.17\pm 0.04`$ (Kogut et al., 2003) and $`\tau _{\mathrm{es}}=0.17\pm 0.07`$ (Spergel et al., 2003) seem to indicate a much earlier epoch or reionization of $`z_{\mathrm{rei}}=17\pm 5`$ (Spergel et al., 2003). How can these measurements be reconciled with the quasar Gunn-Peterson observations? Do they indicate a source of ionizing photons that cannot be accounted for through standard star formation?

In this paper, we first examine how quasar transmission measurements constrain cosmological parameters. In particular, we use the fact that the first generation of UV-generating objects are in fact the tail of the distribution — the rare events that have collapsed to high enough density to produce stars and/or quasars. They thus provide a unique probe of the small scale power spectrum (at $`1`$ Mpc scales) in an analogous way that the X-ray cluster mass function probes $`10`$ Mpc scales. Our constraints are derived from quasars absorption measurements, using theoretical predictions based on the detailed semi-analytic model developed by Chiu and Ostriker (2000), with some improvements derived from recent hydrodynamic and semi-analytical work on reionization.

Second, we address the question of the consistency of the WMAP and quasar measurements, paying particular attention to the degeneracy in the WMAP data between the optical depth $`\tau _{\mathrm{es}}`$ and the spectral index $`n`$. Combining the WMAP constraints with the quasar transmission data, our analysis favors a model with somewhat lower $`\tau _{\mathrm{es}}`$ and $`n`$ than implied by the WMAP data alone, but that are still consistent at the 1-$`\sigma `$ level.

The organization of this paper is as follows. In § 2, we describe the modeling of reionization and Gunn-Peterson absorption. In § 3, we compare observations from both quasars and WMAP to model predictions in order to constrain cosmological parameters. In § 4, we discuss our results and discuss their implications. We summarize and conclude in § 5.

## 2 Modeling Reionization and Gunn-Peterson Absorption

### 2.1 Summary of the Model

The details of the semianalytic model are described in Chiu and Ostriker (2000); the basic principles are summarized here. It is based on a two-phase model of the universe in which a statistical filling factor for ionized gas is self-consistently calculated. It is assumed that the cold, neutral phase has no sources, and evolves passively with the expansion of the universe. The hot, ionized phase contains the ionizing sources and evolves in line with local particle and energy conservation averaged over the phase. The temperature of each phase is calculated using standard physics, with photoheating as the source of heat in the ionized phase and cooling via multiple mechanisms, including H2, atomic lines, and Compton scattering. The time evolution is determined by particle and energy conservation. In each case, we only consider the regions outside of collapsed gas halos, which, of course will rise to the virial temperature. Such collapsed halos are calculated separately, and considered potential sources of ionizing radiation.

The abundance and properties of these potential ionizing sources are calculated on the basis of the Press-Schechter formalism (Press and Schechter, 1974), constrained by the Jeans’ criterion (which utilizes the calculated gas temperatures) and by a cooling criterion (the cooling time must be less than the dynamical time). Cooling in the halos includes the important contributions from H2 cooling. Halos that statisfy the Jeans criterion and that can cool efficiently are sources of ionizing radiation. We calculate the luminosity of each ionizing source using the Schmidt law

$$L(M_b)=(1f_{})M_bc^2ϵ_{\mathrm{esc}}ϵ_{}ϵ_{\mathrm{UV}}t_{\mathrm{dyn}}^1,$$
(1)
where $`f_{}`$ is the fraction of baryons already turned into stars (a small correction), $`M_b`$ is the baryonic mass of a halo, $`ϵ_{\mathrm{esc}}`$ is the escape fraction from the halo, $`ϵ_{}`$ is a resolution factor (related to the fraction of gas that can form stars) determined through calibration (see Chiu and Ostriker (2000) equation (28), and below), $`ϵ_{\mathrm{UV}}`$ is the mass to UV efficiency (where we have absorbed the notation in Chiu and Ostriker (2000) of $`ϵ_{\mathrm{hm}}ϵ_{\mathrm{UV}}`$ into a single efficiency), and $`t_{\mathrm{dyn}}`$ is the dynamical time of the halo. Note that this halo luminosity is somewhat simpler than that used in Chiu and Ostriker (2000), but is consistent with that used in hydrodynamic simulations.

Calibration was done by comparing the redshift of overlap, the ionizing intensity, and the fraction of baryonic mass in stars between the semi-analytic results and the hydrodynamic simulation of Gnedin (2000ab). This simulation was chosen because it is the only published simulation which simulates a statistically “average” universe (as opposed to an individual halo collapse) with sufficient resolution to follow the overlap process, and continues at least redshift 4. The cosmological parameters were $`\mathrm{\Omega }_0=0.3`$, $`\mathrm{\Omega }_\mathrm{\Lambda }=0.7`$, $`h=0.7`$, $`\mathrm{\Omega }_b=0.04`$, $`n=1`$, and $`\sigma _8=0.9`$. The data used for calibration from Gnedin (2000ab) are the redshift of overlap $`z_{}7`$, the ionizing intensity of $`J_{21}0.3`$ at $`z=4`$, and the fraction of baryonic mass in stars $`f_{}0.04`$ at $`z=4`$. The values of $`ϵ_{}`$ and $`ϵ_{\mathrm{UV}}`$ in the semianalytic model were adjusted to best match these three values. The results were $`ϵ_{}=0.03`$ and $`ϵ_{\mathrm{UV}}=1.2\times 10^5`$, which gave $`z_{}=7`$, $`J_{21}=0.6`$, and $`f_{}=0.05`$. Note that Gnedin (2000ab) used $`ϵ_{}=0.05`$ and $`ϵ_{\mathrm{UV}}=4\times 10^5`$, which are remarkably similar to the values in the semianalytic model given the differences in approach. In fact, because we take each halo as a whole, it should not be surprising that $`ϵ_{}`$ is smaller for the semianalytic model than for the hydrodynamic model.

The determination of the escape fraction $`ϵ_{\mathrm{esc}}`$ bears some additional discussion. We use a simple Strömgren sphere approximation (see Appendix) and derive

$$ϵ_{\mathrm{esc}}\left(1\eta \right)^2,$$
(2)
where the quantity $`\eta `$ is given by

$$\eta =\text{Min}[1,\sqrt{\frac{\mathrm{\Delta }_v\overline{\rho _b}\alpha _R}{3ϵ_{}ϵ_{\mathrm{UV}}m_\mathrm{H}^2c^2E_0^1t_{\mathrm{dyn}}^1}}],$$
(3)
where $`\mathrm{\Delta }_v178`$ is the virial overdensity, $`\overline{\rho _b}`$ is the mean baryonic density, $`\alpha _R`$ is the recombination rate, $`m_\mathrm{H}`$ is the hydrogen mass, and $`E_0`$ is 13.6 eV. In this approximation, we assume the halo baryons have a $`r^2`$ density run, the UV output per baryon is constant in the halo, the halo is in ionization equilibrium, and recombinations determine the amount of photons absorbed locally.

The ionizing source luminosity defined above, then, is the luminosity seen by the IGM — i.e., outside of the halos. A further assumption is made in the pre-overlap stage: that at any particular time, the ionized volume around an isolated source is linearly related to its UV luminosity. The constant of proportionality is determined through the global ionization and energy balance, but the size of the H II regions as a function of luminosity is a linear relation. This assumption is consistent with previous work on cosmological H II regions (e.g., Shapiro and Giroux 1987).

The outputs of the model that we use below are the filling factor $`Q`$, the temperature $`T_4`$, and the ionizing intensity $`J_{21}`$.

### 2.2 Density Distribution of Cosmic Gas

One of the most important determinants of how reionization evolves is the degree of gas clumping. The clumping factor, defined by

$$𝒞\rho _b^2/\rho _b^2,$$
(4)
is important not only for determining the ionization balance (and hence the filling factors), but also for determing the *thermal* balance. Note that this clumping factor is only used in the ionized region.

In order to determine the clumping factor, the probability density function (PDF) for the gas density must be known. Miralda-Escudé et al. (2000) found that a good fit to the volume-weighted PDF as seen in hydrodynamical simulations is

$$P_V(\mathrm{\Delta })d\mathrm{\Delta }=A\mathrm{exp}\left[\frac{(\mathrm{\Delta }^{2/3}C_0)^2}{2(\delta _0/3)^2}\right]\mathrm{\Delta }^\beta d\mathrm{\Delta }$$
(5)
where the overdensity $`\mathrm{\Delta }=\rho _b/\overline{\rho _b}`$, $`\delta _0=7.61(1+z)`$ is related to the linear rms gas density fluctuation, and $`\beta `$ is related to the density run at high densities. The parameters $`A`$ and $`C_0`$ are found by requiring the mass and volume to be normalized to unity.

This formula, with the parameterizations of $`\beta `$ and $`\delta _0`$ found in Miralda-Escudé et al. (2000), has been used in many others in calculating characteristics of reionization (e.g., Songaila and Cowie 2002, Fan et al.2002,). However, it is not often noted that the parameterization is based on a single simulation with a particular set of cosmological parameters. In particular, Miralda-Escudé et al. (2000) used the simulation reported in Miralda-Escudé et al. (1996), which was a $`\mathrm{\Lambda }`$-CDM model with $`\mathrm{\Omega }_0=0.4`$, $`\mathrm{\Omega }_\mathrm{\Lambda }=0.6`$, $`\mathrm{\Omega }_bh^2=0.015`$, $`h=0.65`$, $`n=1`$, and $`\sigma _8=0.79`$. Since in this paper we are considering variants in the background cosmology, it is certainly not sufficient to use the Miralda-Escudé et al. (2000) parameterization for the gas PDF. In particular, the amount of power at small scales will depend significantly on all the cosmological parameters (especially $`\sigma _8`$). We therefore analyzed a second simulation (which uses different cosmological parameters) and to generate a second set of parameter fits. The second simulation is one by Cen (2002b) in which $`\mathrm{\Omega }_0=0.3`$, $`\mathrm{\Omega }_\mathrm{\Lambda }=0.7`$, $`\mathrm{\Omega }_bh^2=0.017`$, $`h=0.67`$, $`n=1`$, and $`\sigma _8=0.9`$. The values of $`\beta `$ and $`\delta _0`$ are given for both Miralda-Escudé et al. (2000) and Cen (2002b) are provided in Table 1.

As is clear from Table 1, the value of $`\delta _0`$ and to a lesser extent $`\beta `$ depend on cosmology. This is not surprising given that $`\delta _0`$ *must* depend on the power spectrum. The cosmological dependence on $`\beta `$ is less clear, especially since the best fit value of $`\beta `$ depends significantly on the density run in collapsed halos. Because the Cen (2002b) simulation is at much higher resolution and includes more realistic physics, we simply fit $`\beta `$ to the Cen results as a function of redshift, with a maximum of $`\beta _{\mathrm{max}}=2.5`$, corresponding to an isothermal sphere:

$$\beta \text{Min}[2.5,3.24.73/(1+z)],z4$$
(6)
where we are only considering redshifts $`z4`$.

As for a prescription for finding $`\delta _0`$, we note that the gas PDF also predicts the fraction of mass in collapsed virialized halos:

$$f_{\mathrm{gas}\mathrm{PDF}}(\mathrm{collapsed})_{6\pi ^2}^{\mathrm{}}\mathrm{\Delta }P_V(\mathrm{\Delta })𝑑\mathrm{\Delta }.$$
(7)
The integration point $`6\pi ^2`$ is derived from the fact that for a singular isothermal sphere, the local overdensity density at the virial radius is $`6\pi ^2`$. Therefore, we are approximately taking into account all gas within virialized halos. Similarly, from linear theory, we can use the Press-Schechter formalism (Press and Schechter, 1974) to find the same quantity if we know the correct filtering radius $`R_f`$:

$$f_{\mathrm{PS}}(\mathrm{collapsed})\sqrt{\frac{2}{\pi \sigma _{R_f}^2}}_{\delta _c}^{\mathrm{}}\mathrm{exp}\left[\frac{\delta ^2}{2\sigma _{R_f}}\right]𝑑\delta ,$$
(8)
where $`\delta _c1.69`$ and $`\sigma _{R_f}`$ is the linear RMS mass fluctuation filtered with a tophat of radius $`R_f`$. However, the mass fraction should be equal by the two calculations

$$f_{\mathrm{gas}\mathrm{PDF}}(\mathrm{collapsed})=f_{\mathrm{PS}}(\mathrm{collapsed}).$$
(9)
By analyzing both simulations, we find the following relation leads to satisfactory results:

$$R_fR_J/4,$$
(10)
where $`R_J`$ is the Jeans length defined by

$$R_J=\sqrt{\frac{5\pi k_BT}{12G\overline{\rho }\mu m_Ha^2}}$$
(11)
and $`T`$ is the gas temperature, $`\overline{\rho }`$ is the mean total density, $`\mu `$ is the baryons per particle. With this prescription for finding $`\delta _0`$ and $`\beta `$, we need only use the normalization constraints to find $`A`$ and $`C_0`$ for an arbitrary background cosmology. While by no means perfect, this procedure should account to first order the dependence of the gas PDF on background cosmology, and is certainly more accurate than a treating the evolution of $`\delta _0`$ as independent of cosmology.

Given the gas PDF, one can determine the clumping factor by simple integration over all densities. However, in reality, there is a high density cut-off due to the fact that the very high density regions are not participating in the ionization balance due to self-shielding. However, the determinination of the cutoff in our prescription is not difficult, since we have explicitly made a separation between the “in-halo” and “out-of-halo” calculations. Thus, we calculate the clumping factor by integrating

$$𝒞=\frac{_0^{6\pi ^2}\mathrm{\Delta }^2P_V(\mathrm{\Delta })𝑑\mathrm{\Delta }}{_0^{6\pi ^2}P_V(\mathrm{\Delta })𝑑\mathrm{\Delta }}.$$
(12)
As above, we assume that gas with overdensity $`>6\pi ^2`$ is within collapsed halos.

We should note that this treatment produces smaller clumping factors ($`5`$) than typically used in semi-analytic treatments. This is because we are essentially using a “hybrid” between the reionization treatment of Miralda-Escudé et al. (2000) and the usual clumping factor methods. Conceptually, Miralda-Escudé et al. (2000) considered a density-dependent ionization fraction — they assumed that the universe was ionized up to a critical density, above which was neutral and self-shielded. Here, we are assuming that the collapsed halos are self-shielded, and then treat the rest of the universe using the clumping factor approach. As described in the appendix, at redshifts before full reionization, halos are optically thick if they do not contain sources themselves. By design, this treatment does not include screening and evaporation of mini-halos (e.g., Barkana and Loeb (2002)).

To summarize, the semi-analytic model actually has *three* phases: the “interior of halos,” the hot ionized H II regions, and the cold neutral regions. The “interior of halos” phase includes all the gas in the universe that is in halos with baryonic mass greater than the Jeans mass. We assume that these regions have are virialized with isothermal sphere density runs, and therefore include all gas with overdensity greater than $`6\pi ^2`$. In this phase, the cooling criterion and the Schmidt law are used to calculate the ionizing source function, and the photoionizations and recombinations are treated through the escape fraction calculation (a Strömgren sphere approximation). Those photons that “escape” this phase then enter the “second phase,” which are the ionized H II regions. In this phase, a clumping factor approach is used, with a high density cutoff at $`6\pi ^2`$. The last phase is the cold neutral part of the universe, which uses up ionizing photons only when “converted” into the ionized phase.

### 2.3 Modeling Gunn-Peterson Absorption

The statistics of Gunn-Peterson absorption has been reviewed by numerous others (Miralda-Escudé, Haehnelt, and Rees, 2000; McDonald and Miralda-Escudé 2001; Fan et al.2002; Songaila and Cowie 2002). To summarize, we begin with the stardard optical depth at resonance for a gas with neutral density $`n_{\mathrm{HI}}`$ at redshift $`z`$:

$$\tau =\frac{3\mathrm{\Lambda }_{2\mathrm{p}1\mathrm{s}}\lambda _\alpha ^3n_{\mathrm{HI}}(z)}{8\pi H_0\sqrt{\mathrm{\Omega }_0(1+z)^3+\mathrm{\Omega }_\mathrm{\Lambda }}},$$
(13)
where $`\mathrm{\Lambda }_{2p1s}`$ is the decay rate ($`6.25\times 10^8`$ sec-1), $`\lambda _\alpha `$ is the Lyman-$`\alpha `$ wavelength ($`1.216\times 10^5`$ cm), and the other symbols have their usual meaning in a cosmological context.

Now, we make the standard assumptions of ionization equilibrium and a uniform ionizing intensity in the ionized regions to obtain the neutral fraction. Considering only photoionization and radiative recombination

$$n_{\mathrm{HI}}\mathrm{\Gamma }_{21}J_{21}=X(X+Y/4)/m_p^2R_4T_{4,0}^{0.7}\overline{\rho }_b^2\mathrm{\Delta }^2,$$
(14)
where $`Y=0.24`$ is the helium mass fraction (we assume helium is singly ionized), $`R_44.2\times 10^{13}\mathrm{cm}^3\mathrm{s}^1`$ is the recombination rate at a temperature $`T=10^4\mathrm{K}`$, $`T_{4,0}=T_0/10^4\mathrm{K}`$ is the temperature at mean density in units of $`10^4`$ Kelvin, and $`\mathrm{\Gamma }_{21}`$ is the photoionization rate for $`J_{21}=1`$. We use the method of Hui and Gnedin (1997) to calculate the temperature at mean density after reionization (see also Hui and Haiman 2003). Plugging in values for the various other constants ($`\mathrm{\Gamma }_{21}=4.35\times 10^{12}`$ s-1, $`m_p=1.67\times 10^{24}`$ g, $`\overline{\rho }_b=1.88\times 10^{29}(1+z)^3\mathrm{\Omega }_bh^2`$ g cm-3) gives

$$n_{\mathrm{HI}}=7.63\times 10^{12}\frac{(\mathrm{\Omega }_bh^2)^2T_4^{0.7}\mathrm{\Delta }^2(1+z)^6}{J_{21}}\mathrm{cm}^3$$
(15)
or

$$\tau =\frac{0.316\mathrm{\Omega }_b^2h^3(1+z)^6}{\sqrt{\mathrm{\Omega }_0(1+z)^3+\mathrm{\Omega }_\mathrm{\Lambda }}}\frac{T_4^{0.7}}{J_{21}}\mathrm{\Delta }^2$$
(16)
where $`T_4`$, $`\mathrm{\Delta }`$, and $`J_{21}`$ depend on the redshift $`z`$, but only $`\mathrm{\Delta }`$ also depends on spatial location. It is helpful to reformulate this for $`\mathrm{\Omega }_0(1+z)^3\mathrm{\Omega }_\mathrm{\Lambda }`$, and scaled from $`\mathrm{\Omega }_bh^2=0.02`$ and $`\mathrm{\Omega }_0h^2=0.14`$:

$$\tau =1.536\left(\frac{1+z}{1+5.5}\right)^{4.5}\sqrt{\frac{0.14}{\mathrm{\Omega }_0h^2}}\left(\frac{\mathrm{\Omega }_bh^2}{0.02}\right)^2\frac{T_4^{0.7}}{J_{21}}\mathrm{\Delta }^2.$$
(17)
Defining the optical depth for a uniform medium ($`\mathrm{\Delta }=1`$) $`\tau _u`$ is

$$\tau _u1.536\left(\frac{1+z}{1+5.5}\right)^{4.5}\sqrt{\frac{0.14}{\mathrm{\Omega }_0h^2}}\left(\frac{\mathrm{\Omega }_bh^2}{0.02}\right)^2\frac{T_4^{0.7}}{J_{21}},$$
(18)
we can define the mean transmitted flux ratio at a given redshift

$$𝒯_z=e^{\tau _u\mathrm{\Delta }^2}=Q_0^{\mathrm{}}P_V(\mathrm{\Delta })e^{\tau _u\mathrm{\Delta }^2}𝑑\mathrm{\Delta },$$
(19)
where $`Q`$ is the volume filling factor for the ionized regions (we have now adopted the standard notation — Chiu and Ostriker (2000) used “$`f`$” as the volume filling factor). Note that strictly speaking, there should be a high density cutoff. However, since the exponential has a power of $`\mathrm{\Delta }^2`$, for large values of $`\mathrm{\Delta }`$, the optical depth is very large. The high-density end of the PDF contributes very little to the integral, and the results are thus not sensititive to this cutoff. In fact, at the redshifts of interest, the integral is usually dominated by values of $`\mathrm{\Delta }<1`$, as was noted by Barkana (2001) and others.

Our semi-analytic model provides the values of $`Q`$, $`T_4`$, and $`J_{21}`$, and from the previous section, we have a prescription for determining $`P_V(\mathrm{\Delta })`$. Thus, we have a complete model from which to calculate the expected Gunn-Peterson absorption.

## 3 Comparing Model Predictions and Observations

Having calibrated our model using various hydrodynamic simulations, we are left with one free parameter, $`ϵ_{\mathrm{UV}}`$, in addition to the cosmological parameters. Below, we first consider the case where $`ϵ_{\mathrm{UV}}`$ is constant and derive a constraint on $`\sigma _8\mathrm{\Omega }_0^{0.5}`$, similar to those from rich clusters, and compare the results and predictions to those found by WMAP. We consider joint constraints from WMAP and the SDSS quasars on the spectral index $`n`$ as well as the time-depedence of $`ϵ_{\mathrm{UV}}`$.

### 3.1 Combining SDSS Quasar Data

The quasar data we use combine the compilation at $`z_{\mathrm{abs}}5.5`$ Songaila and Cowie (2002) (Table 2) and based on measurements of six SDSS quasars at $`z>5.7`$ from Becker et al. (2001), Fan et al. (2003), and White et al. (2003) (Table 3). Note that we did not include the J1044-0125 ($`z=5.74`$) data in Table 3 because it is already incorporated into the compilation of Songaila and Cowie (2002).

In order to combine these data, we must take into account the fact that at redshifts $`z5.6`$, the uncertainties in the transmission data are dominated by intrinsic scatter, while at higher redshifts, the uncertainties are dominated by measurement errors. Since our model predicts the *mean* transmission, our likelihood function must properly account for both measurement error as well as the scatter. For a compilation of transmissions, such as that by Songaila and Cowie (2002), the contribution of the data $`D_c=\{(z_j,𝒯_j)\}`$ to the likelihood is simply

$`P(𝒯_z|D_c)`$ $`=`$ $`\mathrm{exp}\left[{\displaystyle \frac{1}{2}}{\displaystyle \underset{j}{}}{\displaystyle \frac{(𝒯_{z_j}𝒯_j)^2}{\sigma _{j,\mathrm{mean}}^2}}\right],`$ (20)
$`\sigma _{j,\mathrm{mean}}^2`$ $`=`$ $`\sigma _{j,\mathrm{scatter}}^2/\sqrt{N}.`$ (21)
The contribution to the likelihood from the compiled data is estimated this way.

For the individual transmission data, especially in the case where both measurement error and scatter are important, the simplest way to do this is to estimate scatter $`\sigma _{\mathrm{scatter}}`$ and add in quadrature with the measurement error $`\sigma _{\mathrm{meas}}`$ to obtain the total uncertainty in the mean for each data point. In particular, for individual transmission data $`D_i=\{(z_j,T_j)\}`$with known (Gaussian) errors $`\sigma _{\mathrm{meas}}`$ and $`\sigma _{\mathrm{scatter}}`$, probability of a mean transmission as a function of redshift $`𝒯_z`$ is

$$P(𝒯_z|D_i)=\mathrm{exp}\left[\frac{1}{2}\underset{j}{}\frac{(𝒯_{z_j}T_j)^2}{\sigma _{j,\mathrm{meas}}^2+\sigma _{j,\mathrm{scatter}}^2}\right].$$
(22)
Note that if the measurement error is negligible, then Equation 22 reduces to Equation (21). To estimate $`\sigma _{\mathrm{scatter}}`$, we use the Songaila and Cowie (2002) compilation for $`z5.6`$ (interpolated), and use the data themselves (binned) to determine $`\sigma _{\mathrm{scatter}}`$ at higher redshifts. These results are also tabluated in Table 3. The total likelihood function, $``$, then is given by the product of the two separate likelihoods.

### 3.2 SDSS Quasar Contraints with a Constant UV-Efficiency Compared with WMAP

For this initial comparison, we fix the other cosmological parameters to their WMAP-only best fit values: $`n=0.99`$, $`\mathrm{\Omega }_bh^2=0.024`$, and $`h=0.72`$. We use the WMAP-only results from Spergel et al. (2003), Table 1, to ensure clarity in the data underlying our analysis. Note, however, that these last two values have additional independent lines of support, as noted in Spergel et al. (2003).

We use a Bayesian method to determine the cosmological constraints. The priors we use are $`\mathrm{\Omega }_0[0.15,0.40]`$ and $`\sigma _8[0.5,1.0]`$. This region essentially bounds the 99% region reported for the combined CMB and cluster analysis of Melchiorri et al. (2002).

In order to constrain $`\mathrm{\Omega }_0`$ and $`\sigma _8`$ separately from $`ϵ_{\mathrm{UV}}`$, we marginalize over the latter parameter to yield the marginalized likelihood distribution:

$$(\mathrm{\Omega }_0,\sigma _8)=(\mathrm{\Omega }_0,\sigma _8,ϵ_{\mathrm{UV}})𝑑ϵ_{\mathrm{UV}}.$$
(23)
The central values and percentile limits are found by integrating over $`(\mathrm{\Omega }_0,\sigma _8)`$.

The result of this analysis is presented in Figure 1, where we plot the likelihood contours in the plane $`\mathrm{\Omega }_0\sigma _8`$. As can be seem by the figure, there is a considerable degeneracy in this plane similar to that derived from rich clusters. Thus, the cosmological constraint from Gunn-Peterson absorption can be summarized by

$$\sigma _8\mathrm{\Omega }_0^{0.5}=0.33\pm 0.01,$$
(24)
where the error term is statistical error only (68%).

Given the complex physics of the reionization model, it is difficult to determine a meaningful systematic error-bar. However, we did investigate the effects of changing $`h`$ and $`n`$. The net result was to shift the right-hand-side of equation (24) so that

$$\sigma _8\mathrm{\Omega }_0^{0.5}(0.33\pm 0.01)\left(\frac{0.72}{h}\right)^{0.9+\frac{n1}{2}}\left(\frac{n}{0.99}\right)^{0.85},$$
(25)
so that lower levels of $`h`$ or $`n`$ increase the right-hand-side. Thus, lowering $`h`$ or $`n`$ slightly would make the Gunn-Peterson constrain more consistent.

These results for a constant UV efficiency are somewhat discordant with WMAP’s marginalized value of $`\sigma _8\mathrm{\Omega }_0^{0.5}=0.48\pm 0.12`$, although only at the 1.25-$`\sigma `$ level. They are consistent with the cluster-determinations of Bahcall et al. (2002) of $`\sigma _8\mathrm{\Omega }_0^{0.6}=0.33\pm 0.03`$ (note different exponent), although again at about the 1-$`\sigma `$ level.

The actual *value* of $`ϵ_{\mathrm{UV}}`$ bears some discussion. Although we have left this parameter as completely free, there certainly are astrophysical constraints on its value. In our 95% region, we find that $`ϵ_{\mathrm{UV}}2\times 10^5`$ (see left panel of Figure 2). For a Scalo mass function with metal-enriched stars (1/20th solar metalicity), the value calculated from population synthesis is $`5\times 10^5`$ (Wiyithe and Loeb, 2002). Given the uncertainties and the relative simplicity of our model, the correspondence is quite remarkable.

The more striking inconsistency is that in all cases, the optical depth to electron scattering,

$$\tau _{\mathrm{es}}_0^{1000}𝑑z\frac{dt}{dz}c\sigma _\mathrm{T}n_e,$$
(26)
where $`\sigma _\mathrm{T}=6.652\times 10^{25}\mathrm{cm}^2`$ is the Thomson cross-section and $`n_e`$ is the electron density, is too low. Our model finds that $`\tau _{\mathrm{es}}0.06`$ (see right panel of Figure 2), whereas Spergel et al. (2003) and Kogut et al. (2003) report WMAP-only best fit values of $`\tau _{\mathrm{es}}=0.17\pm 0.07`$ and $`\tau _{\mathrm{es}}=0.17\pm 0.04`$, respectively. Several articles written in the wake of these results (e.g., Haiman and Holder 2003; Ciardi, Ferrara, and White 2003; Cen 2003) have discussed this inconsistency and indicated that star formation must begun much earlier than previously thought. We consider this issue in the following analysis.

### 3.3 Combined WMAP and SDSS Constraints with a Time-Varying UV-Efficiency

There are good astrophysical reasons to believe that $`ϵ_{\mathrm{UV}}`$ may effectively increase with redshift. For instance, the first generation of stars would have been metal-free, and stellar models predict them to have an significantly higher UV output per baryon (about $`4\times `$ higher was reported by Wiyithe and Loeb (2002)). In addition, the creation of metals and henceforth dust, which would obscure ionizing sources, may also lead to lower *effective* efficiencies at lower redshift than at higher redshift.

Thus we now consider constraints on $`n`$ and the time-dependence of $`ϵ_{\mathrm{UV}}`$ from WMAP and SDSS quasars jointly. To limit the dimensionality of the parameter space, we keep the baryon abundance and the Hubble constant at their WMAP best fit values, and in addition fix $`\mathrm{\Omega }_0h^2=0.14`$ to its WMAP best fit value. For each value of $`n`$ and $`\tau _{\mathrm{es}}`$, WMAP predicts a unique best-fit value of the normalization $`\sigma _8`$ (based on the code provided by Verde et al. 2003 and the accompanying data files from Hinshaw et al. 2003 and Kogut et al. 2003). The normalization given approximately by (for our fixed values of $`h`$, $`\mathrm{\Omega }_0h^2`$, and $`\mathrm{\Omega }_bh^2`$)

$$\sigma _8e^{\tau _{\mathrm{es}}}=0.765+0.6(n1),$$
(27)
where formula (27) is good to better than 0.5%. Note that at the WMAP best-fit values of $`n=0.99`$ and $`\tau _{\mathrm{es}}=0.17`$, this formula gives the WMAP best-fit value of $`\sigma _8=0.9`$.

We assume the following heuristic form for the time-dependence of the UV-efficiency:

$$ϵ_{\mathrm{UV}}=ϵ_{\mathrm{UV},0}\left(1+Ae^{f_{}/f_{,\mathrm{crit}}}\right),$$
(28)
where $`f_{}`$ is the fraction of baryons in stars. Thus, the luminosity of each ionizing source is given by

$$L(M_b)=(1f_{})M_bc^2ϵ_{\mathrm{esc}}ϵ_{}ϵ_{\mathrm{UV},0}t_{\mathrm{dyn}}^1\left(1+Ae^{f_{}/f_{,\mathrm{crit}}}\right).$$
(29)
The new parameters $`A`$ and $`f_{,\mathrm{crit}}`$ define the time-dependence of the efficiency. The essential features are that at very high redshift, the total efficiency is $`(1+A)`$ times greater than at lower redshift, with the transition occuring at $`f_{}f_{,\mathrm{crit}}`$. This form is motivated by the notion that the first stars, being metal-free, should have had a higher effective efficiency, but that as metals were produced by this generation of stars and reinjected back into the IGM, the efficiency would decrease. However, we emphasize that this model was simply selected heuristically, and was not based on any detailed model of star formation feedback on the UV-efficiency. For instance, it should be noted that the “enhancement” factor may be a combination of changes in the intrinsic efficiency $`ϵ_{\mathrm{UV}}`$ and the resolution factor $`ϵ_{}`$ related to the fraction of cooling gas that forms stars. To give a sense of the meaning of $`f_{,\mathrm{crit}}`$, Table 4 presents the transition redshifts $`z_{\mathrm{trans}}`$ where $`f_{}(\mathrm{ln}2)f_{,\mathrm{crit}}`$ so that $`ϵ_{\mathrm{UV}}(z_{\mathrm{trans}})ϵ_{\mathrm{UV},0}(1+A/2)`$.

We fix the value of $`f_{,\mathrm{crit}}=10^4`$ (we discuss the sensitivity to this parameter below). Thus, for each point in the $`n\tau _{\mathrm{es}}`$ plane, our procedure searching parameter space is as follows: (1) Pick value of $`A`$; (2) Adjust $`ϵ_{\mathrm{UV},0}`$ so that $`\tau _{\mathrm{es},\mathrm{model}}`$ (calculated by model) is consistent with $`\tau _{\mathrm{es}}`$; (3) Stop if $`\chi ^2`$ is minimized, otherwise go back to step (1). For each point in the $`n\tau _{\mathrm{es}}`$ plane, we have a value of a best fit $`\stackrel{~}{\chi }^2(n,\tau _{\mathrm{es}})`$. We treat $`\stackrel{~}{\chi }^2`$ as approximately the $`2\mathrm{ln}`$, where $`(n,\tau _{\mathrm{es}}|f_{,\mathrm{crit}})`$ is the likelihood function, and integrate (with uniform priors) to find find the confidence regions (the results are insensitive to the priors).

The effective optical depth for several models is shown in Figure 3. The three illustrative “good-fitting” models are in the 95% region, and the “badly-fitting” model is outside the region. Clearly, high values of $`n`$ and $`\tau _{\mathrm{es}}`$ are ruled out by the data.

Considering all the constraints together, the results are shown in Figure 4. Our analysis strongly favors a narrow range (95% errors):

$`\tau _{\mathrm{es}}`$ $`=`$ $`0.11_{0.03}^{+0.02}`$ (30)
$`n`$ $`=`$ $`0.96_{0.03}^{+0.02}`$ (31)
for $`f_{,\mathrm{crit}}=10^4`$. This range is consistent with the WMAP results, also shown in the Figure. The effect of changing $`f_{,\mathrm{crit}}`$, is only to shift the $`\tau _{\mathrm{es}}`$ constraint up or down. For instance, the constraint for $`f_{,\mathrm{crit}}=3\times 10^5`$ is $`\tau _{\mathrm{es}}=0.12_{0.03}^{+0.02}`$, with no change in the $`n`$ constraint. The results for increasing $`f_{,\mathrm{crit}}`$ are in the opposite direction, decreasing the constraint on $`\tau _{\mathrm{es}}`$ by $`0.02`$ for a factor of 10 increase in $`f_{,\mathrm{crit}}`$.

## 4 Discussion

In order to understand the quasar constraints, let us return to the constant efficiency case. Consider the following “subset” of the quasar data:

Requiring $`z_{\mathrm{overlap}}6`$ and only using the Ly$`\alpha `$ constraint from quasar J1148+5251 ($`z_{\mathrm{em}}=6.43`$);

Using the Ly$`\beta `$ constraint from SDSS 1030+0524 ($`z_{\mathrm{em}}=6.28`$);

Using the Songaila and Cowie (2002) constraints at $`z5`$.

Figure 1 also shows the effects of combining these three constraints. The intersection (Q1+Q2+Q3) subset of data appears to adequately depict the upper limit of the combined dataset. In particular, these three sets of data contain much of the “information” in the full $`\chi ^2`$ treatment. The lower limit of the combined dataset actually also contained within the constraint data in Q3, but changes as a function of $`z_{\mathrm{overlap}}`$. Since Q3 is the union of all these constraints for all values of $`z_{\mathrm{overlap}}`$, the lower limit does not appear.

This illustration shows that given the three degrees for freedom $`\mathrm{\Omega }_0`$, $`\sigma _8`$, and $`ϵ_{\mathrm{UV}}`$, the quasar data constrains to essentially a line in this space. This is because for each $`\mathrm{\Omega }_0`$, the data (Q1+Q2) at $`z6`$ and the data (Q3) at $`z4.5`$ provide strong, essentially independent, constraints on the remaining parameters $`\sigma _8`$ and $`ϵ_{\mathrm{UV}}`$. Of course, the “length” and “width” of the the line (i.e., the exact shape, as illustrated in Figure 1) depend on the likelihood function in detail.

Now let us return to the time-varying efficiency case. For each value of $`f_{,\mathrm{crit}}`$, we have four parameters $`n`$, $`\tau _{\mathrm{es}}`$, $`ϵ_{\mathrm{UV},0}`$, and $`A`$. Heuristically, we now have three constraints on the output: the data (Q1+Q2) at $`z6`$, the data (Q3) at $`z4.5`$, and the self-consistency of $`\tau _{\mathrm{es}}`$ (value “in” = value “out”). Therefore, once again, we expect a “line” in the 4-dimensional parameter space (for fixed $`f_{,\mathrm{crit}}`$). As before, the “length” and “width” of the the line will depend on the likelihood function in detail.

We note also that the implied values for the normalization in our 95% region is $`\sigma _80.83_{0.05}^{+0.03}`$ (based on WMAP best-fit), which for our assumed value of $`\mathrm{\Omega }_0`$ gives $`\sigma _8\mathrm{\Omega }_0^{0.6}=0.38_{0.025}^{+0.015}`$, just overlapping with the cluster measurements $`\sigma _8\mathrm{\Omega }^{0.6}=0.33\pm 0.03`$.

The next question is what are the implications as to the implied values of the UV-efficiency? Do the efficiencies make sense? Table 5 shows a number of parameter combinations from the 95% region. The implied efficiency at high redshift is on the order $`10^4`$; a transition occurs at $`z=1520`$, and the efficiency at low redshift is $`ϵ_{\mathrm{UV},0}=10^{5.55}`$. The span of these efficiencies encompasses those typically calculated through population synthesis of $`10^{54.5}`$. Given uncertainties of factors of a few in the gas collapse fraction (our resolution factor $`ϵ_{}`$), dust absorption (at lower redshift), etc., the values of the efficiencies do not seem unreasonable. They certainly do not approach the upper bound for conversion from nuclear reactions of $`10^3`$. Sokasian et al. (2003) and Wiyithe and Loeb (2003) discuss further the star formation efficiency as relating to Population II and Population III stars.

Our results are consistent with the results of Cen (2003) in that we find that to reach $`\tau _{\mathrm{es}}0.17`$ requires that the spectral index is positively tilted with $`n1.02`$. In this case the effective UV-efficiency was at least $`100\times `$ greater at $`z6`$. However, our calculations indicate that models with such high $`n`$ (and hence high power spectrum normalization $`\sigma _8`$) are inconsistent with quasar transmission measurements at redshift $`z5`$.

Like other authors (e.g., Haiman and Holder 2003, Cen 2003), we find that a value of $`\tau _{\mathrm{es}}=0.17`$ is inconsistent with constraints at $`z6`$ for simple models of reionization, and may require more exotic (though not necessarily implausible) methods for creating ionizing photons, such as mini-quasars or an X-ray background. However, we conclude that WMAP data taken as a whole and quasar observations at $`z6`$ in fact are entirely consistent for reasonable values and time-dependence for the UV-efficiency. Our analysis shows the importance of taking into account the significant degeneracy between $`\tau _{\mathrm{es}}`$ and $`n`$ determined by WMAP.

Having found a consistent set of evolutionary models, let us describe their properties in slightly greater detail. We consider the following two models: a model with constant UV-efficiency with $`h=0.72`$, $`\mathrm{\Omega }_0=0.27`$, $`n=0.99`$, $`\sigma _8=0.64`$, $`ϵ_{\mathrm{UV}}=2.2\times 10^5`$, and $`\tau _{\mathrm{es}}=0.06)`$; and a model with variable efficiency with $`h=0.72`$, $`\mathrm{\Omega }_0=0.27`$, $`n=0.96`$, $`\sigma _8=0.827`$, $`ϵ_{\mathrm{UV}}=7.7\times 10^6(1+17e^{f_{}/10^4})`$, and $`\tau _{\mathrm{es}}=0.11`$. These are both near peak of the likelihood distribution for $`\mathrm{\Omega }_0=0.27`$.

The Gunn-Peterson optical depths are already shown in Figure 3. The best-fit models do not show much difference, as expected since the parameters were fit to these data. Figure 5 shows the reionization properties for the two models. although the last phase of reionization occurs rapidly in both models, the phase in which the average neutral fraction drops from unity to $`0.1`$ takes a significantly longer time in the model with a variable UV-efficiency. This is a “necessary” part of the model in order to achieve a higher electron scattering optical depth. The transition from the higher to the lower UV-efficiency is clearly seen in the stellar baryon fraction $`f_{}`$ at redshift $`z15`$. We note that, even in the case of $`\tau _{\mathrm{es}}=0.11`$, there are not two distinct epochs of ionization, in contrast to the calculations of Cen (2002a). Rather, there is an extended period during which the ionization increases from 1% to 10% before the rapid phase change at $`z6`$.

The thermal properties are show in Figure 6. Again, the rise in the global volume-averaged temperature is much more gradual with the variable UV-efficiency. In addition, in this case, the final temperature is higher, but this is due to the higher power spectrum normalization ($`\sigma _8`$). These temperatures are somewhat lower that the “peak” temperatures derived by Hui and Haiman (2003), but this is probably due to the “sudden” reionization model used in those calculations. A more gradual reionization transition would smooth out the peak, and would appear consistent with our calculations. For instance, Hui and Haiman (2003) consider a “stochastic” reheating process their derived temperatures at $`z4`$ are similar to ours.

Finally, in Figure 7 we show the star formation rate and the electron scattering optical depth to redshift $`z`$. Star formation rates at $`z4`$ from field galaxy measurements have been reported at about $`10^{2.5\pm 0.5}`$ M yr-1 Mpc-1 (Steidel et al., 1999), while a recent lower limit of $`10^{3.1}`$ M yr-1 Mpc-1 has been reported at $`z6`$ (Standway, Bunker, and McMahon, 2003). Our model results imply that there is exists a large population of unobserved sources at $`z6`$. We not that quasar observations alone require that $`\tau _{\mathrm{es}}`$ has reached value of $`0.05`$ at $`z7`$. Finally, the electron scattering optical depth shows that in the constant efficiency case, the full optical depth is reached by $`z10`$, while in the variable efficiency case, the full optical depth is not reached until $`z20`$. This is an observational signature that could be detected by future CMB experiments.

## 5 Summary and Conclusions

We derive constraints on several cosmological parameters based on observations of Gunn-Peterson absorption in high-redshift quasars and WMAP observations (Bennett et al., 2003). We use a semi-analytic model for reionization (Chiu and Ostriker, 2000) that takes into account a number of important physical processes both within collapsing halos (e.g., H2 cooling) and in the intergalactic medium (e.g., H2 cooling, Compton cooling, and photoionization heating). The model is also calibrated to hydrodynamic simulations. We also develop a method for estimating the gas PDF, which is important for properly calculating the mean absorption in the IGM, as a function of cosmological parameters.

We find that the Gunn-Peterson absorption data provide constraints on the power spectrum at small scales in a manner similar to that derived from the cluster mass function. Assuming that the efficiency of producing UV photons per baryon is constant, the constraint takes on the form $`\sigma _8\mathrm{\Omega }_0^{0.5}=0.33\pm 0.01`$ assuming a flat, $`\mathrm{\Lambda }`$-dominated universe with $`h=0.72`$, $`n=0.99`$, and $`\mathrm{\Omega }_bh^2=0.024`$. The best fit for the WMAP data (marginalized over all parameters) is reported as $`\sigma _8\mathrm{\Omega }_0^{0.5}=0.48\pm 0.12`$, which differs by slightly more than $`1\sigma `$. However, the derived value for the optical depth to last-scattering $`\tau _{\mathrm{es}}0.06`$ differs signficantly from the WMAP-determined value of $`\tau _{\mathrm{es}}=0.17\pm 0.04`$.

Since the WMAP constraints on $`\tau _{\mathrm{es}}`$ are somewhat degenerate with the value of the spectral index $`n`$ (Spergel et al., 2003), we then let the primordial spectral index $`n`$ float while fixing the best fit WMAP determination of $`\mathrm{\Omega }_0h^2=0.14`$ and normalizing the power spectrum using WMAP. In addition, we allow the UV-efficiency to have time-dependence. Combining the WMAP constraints with the quasar transmission data, our analysis favors then a model with $`\tau _{\mathrm{es}}=0.11_{0.03}^{+0.02}`$, $`n=0.96_{0.03}^{+0.02}`$, implying $`\sigma _8=0.83_{0.05}^{+0.03}`$ (all at 95% confidence), and an effective UV-efficiency that was at least $`10\times `$ greater at $`z6`$. If future observations confirm this range for the optical depth to electron scattering $`\tau _{\mathrm{es}}`$, then it would appear that no more “exotic” sources of UV-photons are necessary. We are unable to find a model that is consistent with all observational and physical constraints that has an electron scattering optical depth within the 1-$`\sigma `$ range given by the WMAP team. If the EE data finally require a value $`\tau _{\mathrm{es}}0.17`$, then more exotic sources of early ionizing photons will be required that those considered in this paper.

The authors would like to thank Renyue Cen kindly providing his simulation outputs, and Martin Rees, Massimo Ricotti, and David Spergel for useful discussions. XF thanks support from the University of Arizona and a Sloan Research Fellowship

## Appendix A Approximate Model for Escape Fraction

In this appendix, we consider a simple Strömgren sphere approximation for the escape fraction. The basic conditions for star forming halos in our semi-analytic model are as follows:

1. In halos greater than the Jeans mass, the baryonic mass will collapse to a singular isothermal sphere with $`\rho _b=\frac{1}{3}\mathrm{\Delta }_v\overline{\rho }_b(r_v/r)^2`$, where $`\mathrm{\Delta }_v178`$ is the virial overdensity, $`\overline{\rho }_b`$ is the mean baryonic density, and $`r_v`$ is the virial radius. Here the total baryonic mass is $`M_b=\frac{4\pi }{3}\mathrm{\Delta }_v\overline{\rho }_br_v^3`$.
2. In halos where the cooling time $`t_c`$ is less than the dynamical time $`t_d`$, there will be star formation.
3. The rate of UV photon production per unit baryonic mass is $`\theta =c^2ϵ_{}ϵ_{\mathrm{UV}}E_0^1t_d^1`$, where $`ϵ_{\mathrm{UV}}`$ is the mass to UV energy efficiency and $`E_0=13.6`$ eV. Here $`ϵ_{}`$ is the fraction of cooling baryons in the halo that are forming stars (the “resolution” factor), and $`ϵ_{\mathrm{UV}}`$ is the fundamental conversion efficiency from stellar baryons to UV photons.

In regions that are optically thick, so that all photons are absorbed locally, then local ionization equilibrium would imply that $`\theta \rho _b=n_\mathrm{H}^2x^2\alpha `$, where $`n_\mathrm{H}`$ is the total hydrogen density, $`x`$ is the ionization fraction, and $`\alpha `$ is the recombination coefficient. This implies that the ionization fraction is $`x=\sqrt{\theta \rho _b/\alpha }/n_\mathrm{H}r`$. Therefore, as $`r0`$, the ionized fraction also tends to 0.

To check for self consistency, consider the optical depth over a central region. The average density over the sphere $`r`$ is $`\rho _b(r)=\mathrm{\Delta }_v\overline{\rho }_b(r_v/r)^2`$. Optical depth over this region is approximately $`\tau n_\mathrm{H}(1x)\sigma r`$, where $`\sigma `$ is an effective cross section. Since $`n_\mathrm{H}r^2`$, this implies that $`\tau (1x)r^1`$. Therefore, as long as $`(1x)`$ does not vanish, then the optical depth will become very large as $`r0`$. Above, we established that $`x0`$ as $`r0`$ for an optically thick region. We therefore have a self consistent picture for the central part of the halo.

Now what is the condition so that the entire halo is optically thick? For a completely neutral halo, $`\tau \mathrm{\Delta }_v\overline{\rho }_b\sigma r_v/m_{\mathrm{eff}}`$. Plugging in numbers gives $`\tau r_v(1+z)^3/10`$ kpc (in physical units). In comoving units, $`\tau (1+z)^2`$. Therefore at any redshift before full reionization ($`z6`$), all halos of concern are optically thick if fully neutral. They only become optically thin if the neutral fraction is small.

Now consider a model based on a central source approximation. Let $`S(r)`$ be the number of ionizing photons emitted by a central source which pass through a sphere of radius $`r`$. The standard equation for $`S(r)`$ is given by

$$\frac{S(r)}{r}=4\pi r^2n_\mathrm{H}^2x^2\alpha .$$
(A1)
Using the definitions above, we obtain

$$\frac{S(r)}{r}=4\pi \left(\frac{\mathrm{\Delta }_v\overline{\rho }_br_v^2}{3m_{\mathrm{eff}}}\right)^2\frac{x^2\alpha }{r^2}.$$
(A2)
Now consider, the number of photons produced within $`r`$, $`𝒮(r)=\theta M(<r)`$. Differentiating this with respect to $`r`$ gives

$$\frac{𝒮(r)}{r}=\frac{4\pi }{3}\theta \mathrm{\Delta }_v\overline{\rho }_br_v^2.$$
(A3)
Now let us make the ansatz that all these photons can be considered to be radially emitted. Then these two equations can be *combined* to yield

$$\frac{S(r)}{r}=\frac{4\pi }{3}\mathrm{\Delta }_v\overline{\rho }_br_v^2\left(1\frac{\mathrm{\Delta }_v\overline{\rho }_br_v^2}{3\theta m_{\mathrm{eff}}^2}\frac{x^2\alpha }{r^2}\right)$$
(A4)
Because we know that the halo is optically thick if $`x<1`$, let us consider the following prescription: as long as $`x<1`$, all photons are absorbed locally, so $`\frac{S}{r}=0`$. This leads to the expression for the ionized fraction as a function of radius

$$x(r)=\mathrm{Min}[\frac{r}{r_v}\sqrt{\frac{3\theta m_{\mathrm{eff}}^2}{\mathrm{\Delta }_v\overline{\rho }_b\alpha }},1].$$
(A5)
If we define quenching fraction

$$\eta \sqrt{\frac{\mathrm{\Delta }_v\overline{\rho }_b\alpha }{3\theta m_{\mathrm{eff}}^2}},$$
(A6)
then $`\eta 1`$ means that all of the photons are absorbed within the halo. If $`\eta <1`$, then the ionized fraction has reached its maximum value of 1 at the “quenching radius” $`r_q=\eta r_v`$. Note importantly that $`\eta `$ depends only on redshift through the background density $`\overline{\rho }_b`$ and the dynamical time in $`\theta `$, and is essentially *the same for all halos*. (Actually, there is a weak dependence on the halo temperature through the recombination rate $`\alpha `$.)

If $`\eta <1`$, then we use equation (A4) with $`x=1`$ and the boundary condition $`S(r_q)=0`$ to solve for $`S(r)`$ from $`r_q`$ to $`r_v`$. The result is

$$S(r_v)=\frac{4\pi }{3}\theta \mathrm{\Delta }_v\overline{\rho }_br_v^3(1\eta )^2=M_b\theta (1\eta )^2.$$
(A7)
Thus, the factor $`(1\eta )^2`$ is the fraction of UV photons created in the halo that actually escape — the escape fraction $`ϵ_{\mathrm{esc}}`$. At high redshift, when $`\eta 1`$, this factor is 0 — all the photons are consumed within the halo.

To illustrate, consider the “best fit” model with time-variable UV-efficiency described in the main text. Using a temperature of $`10^4`$ Kelvin, $`\mathrm{\Omega }_0h^2=0.14`$, $`\mathrm{\Omega }_bh^2=0.024`$, the quenching fraction is

$$\eta \frac{10^4(1+z)^{0.75}}{\sqrt{ϵ_{}ϵ_{\mathrm{UV}}}}.$$
(A8)
In our “best fit” model, we fixed $`ϵ_{}=0.03`$, while $`ϵ_{\mathrm{UV}}`$ ranged from $`0.8\times 10^5`$ at low redshift to $`1.3\times 10^4`$ at high redshift. At very high redshift, we have

$$\eta \left(\frac{1+z}{57}\right)^{0.75}.$$
(A9)
This implies that at very high redshift, when $`\eta 1`$, all the photons are absorbed locally so that the escape fraction is essentially 0. This is to expected since densities are much higher at high redshift. At redshifts near transition $`15`$, when $`ϵ_{\mathrm{UV}}6\times 10^5`$, we have

$$\eta \left(\frac{1+z}{37}\right)^{0.75},$$
(A10)
implying an escape fraction $`0.2`$, comparable to values used in other semianalytic models. At $`z6`$, we have

$$\eta \left(\frac{1+z}{9}\right)^{0.75},$$
(A11)
implying an escape fraction $`ϵ_{\mathrm{esc}}0.03`$. These values are all within the range used by others in semi-analytic models (e.g., Cen 2003; Wiyithe and Loeb 2002; Haiman and Holder 2003).


format: text

text:
# Entanglement of formation for symmetric Gaussian states

## Abstract

We show that for a fixed amount of entanglement, two–mode squeezed states are those that maximize Einstein–Podolsky–Rosen–like correlations. We use this fact to determine the entanglement of formation for all symmetric Gaussian states corresponding to two modes. This is the first instance in which this measure has been determined for genuine continuous variable systems.

One of the main tasks of Quantum Information Theory is to quantify the entanglement and the quantum correlations that quantum states possess. For that, several entanglement measures have been introduced during the last years MHorodecki . In particular, two of such measures stand out for their well defined physical meaning: the entanglement of distillation and of formation (and the corresponding asymptotic generalization, the entanglement cost) Bennettlongpaper . They quantify the entanglement of a state in terms of the pure state entanglement that can be distilled out of it PRHorodecki  and the one that is needed to prepare it Wootters1 , respectively.

Despite a considerable effort, for the moment we can only evaluate the entanglement of formation (EoF), or the entanglement of distillation for a few sets of mixed states. The reason is that these quantities are defined Bennettlongpaper  in terms of an optimization problem which is extremely difficult to handle analytically. Despite this fact, in a remarkable work Wootters Wootters  managed to derive an analytical expression for the EoF for all two–qubit states. The EoF has been also determined for highly symmetric states (isotropic states Terhal , and Werner states symmetry ). These expressions are important theoretical tools. From a more practical point of view, they can be applied to quantify the entanglement created in current experiments as well as to compare the capability of different experimental set–ups. For low dimensional systems without symmetries, one can still use numerical methods to determine the EoF NumEoF , although they are often not very efficient. For infinite–dimensional systems, however, a numerical approach is not feasible.

Among all quantum states in infinite dimensional systems, Gaussian states play an important role in quantum information. From the experimental point of view, they can be created relatively easily squeezing , and one can use them for quantum cryptography Grangier  and quantum teleportation Kimble ; Bachor . On the theoretical side, separability Geza  and distillability Geza2  criteria for bipartite systems have been fully developed. Moreover, pure Gaussian states are intimately related to Heisenberg’s uncertainty relation since they minimize such a relation for position and momentum operators.

In this work we determine the EoF of all symmetric Gaussian states of two modes. Those states arise naturally in several experimental contexts. For example, when the two output beams of a parametric down converter are sent through optical fibers squeezing , or in atomic ensembles interacting with light expatoms . In order to determine the EoF, we connect the entanglement of pure states, as measured by the von Neumann entropy of the restriction, with the kind of correlations established by Einstein, Podolsky and Rosen (EPR) in their seminal paper EPR . In fact, we show that two–mode squeezed states 2modesq  play a very special role in this relation, since they are the least entangled states for a given correlation of this type. This provides a new characterization of two–mode squeezed states. Finally, we show that the decomposition that leads to the EoF is a decomposition in terms of Gaussian states.

We consider two modes, A and B, with corresponding Hilbert space $`=_A_B`$ and canonical operators $`X_{A,B}`$ and $`P_{A,B}`$. The two–mode squeezed states have the form

$`|\mathrm{\Psi }_s(r):={\displaystyle \frac{1}{\mathrm{cosh}(r)}}{\displaystyle \underset{N=0}{\overset{\mathrm{}}{}}}\mathrm{tanh}^N(r)|N_A|N_B,`$ (1)
where $`r>0`$ is the squeezing parameter, and $`|N`$ denotes the $`N`$-th Fock state, i.e. $`a^{}a|N_A=N|N_A`$, $`b^{}b|N_B=N|N_B`$, where $`a=(X_A+iP_A)/\sqrt{2}`$ and $`b=(X_B+iP_B)/\sqrt{2}`$ are annihilation operators.

In the following, we will denote by $`\psi `$ an arbitrary normalized state in $``$. We define its EPR–uncertainty as follows:

$`\mathrm{\Delta }(\psi ):=\mathrm{min}[1,{\displaystyle \frac{1}{2}}[\mathrm{\Delta }_\psi ^2(X_AX_B)+\mathrm{\Delta }_\psi ^2(P_A+P_B)]],`$ (2)
where, as usual, $`\mathrm{\Delta }_\psi ^2(X):=X\psi |X\psi \psi |X|\psi ^2`$, setting $`\mathrm{\Delta }_\psi ^2(X)=\mathrm{}`$ if $`\psi `$ is not in the domain of $`X`$. Clearly, $`\mathrm{\Delta }(\psi )(0,1]`$. This quantity measures the degree of non–local correlations, and would be zero for the idealized state considered by Einstein, Podolsky, and Rosen EPR . For any state, $`\mathrm{\Delta }(\psi )<1`$ implies the existence of such non–local correlations. Note that this condition is met only if at least one of the uncertainties of $`(X_AX_B)/\sqrt{2}`$ or $`(P_A+P_B)/\sqrt{2}`$ lies below 1 (the standard quantum limit). This implies that the corresponding states must possess a certain squeezing. In fact, the two–mode squeezed states (1) are standard examples of states displaying these correlations since

$`\mathrm{\Delta }[\mathrm{\Psi }_s(r)]=e^{2r}<1.`$ (3)
Any value of $`\mathrm{\Delta }(0,1)`$ is achieved by the two–mode squeezed state with squeezing parameter

$`r_\mathrm{\Delta }:={\displaystyle \frac{1}{2}}\mathrm{ln}(\mathrm{\Delta })`$ (4)
The EPR–uncertainty of a given state $`\psi `$ is certainly related to its entanglement. For pure states this last property is uniquely quantified by the entropy of entanglement, $`E(\psi )`$, which can be determined as follows. Let us write the Schmidt decomposition of $`\psi `$ as notedecomp

$`|\psi ={\displaystyle \underset{N=0}{\overset{\mathrm{}}{}}}c_N|u_N_A|v_N_B,`$ (5)
$`\{u_N\}`$ and $`\{v_N\}`$ are orthonormal bases in $`_{A,B}`$, respectively and $`c=(c_0,c_1,\mathrm{})𝒞`$ where

$$𝒞:=\{cl_{}^2c=1,c_Nc_{N+1}0N\}.$$

Then notelog

$`E(\psi )=𝐞(c):={\displaystyle \underset{N=0}{\overset{\mathrm{}}{}}}c_N^2\mathrm{log}(c_N^2).`$ (6)
Note that this quantity can be infinite for some states. For the two–mode squeezed states (1) we have

$$E[\mathrm{\Psi }_s(r)]=\mathrm{cosh}^2(r)\mathrm{log}[\mathrm{cosh}^2(r)]\mathrm{sinh}^2(r)\mathrm{log}[\mathrm{sinh}^2(r)].$$

With the above definitions we can state the special rôle that two–mode squeezed states play in relation with EPR–correlations and entanglement:

Proposition 1: For all $`\psi `$, $`E(\psi )E[\mathrm{\Psi }_s(r_{\mathrm{\Delta }(\psi )})]`$.

In order to give a clear interpretation of this result, we reformulate it in two equivalent forms: (i) For any given $`\mathrm{\Delta }(0,1)`$, $`E[\mathrm{\Psi }_s(r_\mathrm{\Delta })]=inf_\psi \{E(\psi )\}`$ with $`\psi `$ fulfilling $`\mathrm{\Delta }(\psi )=\mathrm{\Delta }`$; (ii) For any given $`E(0,\mathrm{})`$, $`\mathrm{\Delta }[\mathrm{\Psi }_s(r)]=inf_\psi \{\mathrm{\Delta }(\psi )\}`$ with $`\psi `$ fulfilling $`E(\psi )=E`$ and $`r`$ such that $`E[\mathrm{\Psi }_s(r)]=E`$. Note that the equivalence of this last formulation is ensured by the fact that $`E[\mathrm{\Psi }_s(r)]`$ and $`\mathrm{\Delta }[\mathrm{\Psi }_s(r)]`$ are monotonically increasing and decreasing functions of $`r`$, respectively. The first statement characterizes two–mode squeezed states as the cheapest (regarding entanglement) to achieve a prescribed EPR–uncertainty. The second one characterizes them as those with maximal EPR–correlations (minimal $`\mathrm{\Delta }`$) for any given value of the entanglement.

In order to prove Proposition 1 we introduce two lemmas and the following definition. Given $`c𝒞`$ we define

$`𝜹(c):=1+2{\displaystyle \underset{N=0}{\overset{\mathrm{}}{}}}(c_N^2c_Nc_{N1})N.`$ (7)
We have $`𝜹(c)1`$ and $`𝜹(c)=\mathrm{\Delta }(\psi )`$ whenever $`|u_N=|v_N=|N`$ \[cf. (5)\].

Lemma 1: For all $`\psi `$ with Schmidt decomposition (5), $`\mathrm{\Delta }(\psi )𝜹(c)`$.

Proof: Since $`𝜹(c)1`$ we can restrict ourselves to $`\psi `$ with $`\mathrm{\Delta }(\psi )<1`$. Without loss of generality we can assume that $`\psi |a|\psi =\psi |b|\psi =0`$. Otherwise we can always find $`\psi ^{}`$ fulfilling this condition, with the same Schmidt coefficients as $`\psi `$ and with $`\mathrm{\Delta }(\psi ^{})=\mathrm{\Delta }(\psi )`$ notedispl . We have

$`\mathrm{\Delta }(\psi )`$ $`=`$ $`1+{\displaystyle \underset{N=0}{\overset{\mathrm{}}{}}}c_N^2(u_N|a^{}a|u_N+v_N|b^{}b|v_N)`$
$`{\displaystyle \underset{N,M=0}{\overset{\mathrm{}}{}}}c_Nc_M(u_N|a|u_Mv_N|b|v_M+c.c.).`$
$``$ $`\mathrm{min}[Z(u),Z(v)],`$
where

$`Z(u):=`$ $`1+2{\displaystyle \underset{N=0}{\overset{\mathrm{}}{}}}c_N^2u_N|a^{}a|u_N`$
$`2{\displaystyle \underset{N,M=0}{\overset{\mathrm{}}{}}}c_Nc_M|u_N|a^{}|u_M|^2.`$
Without loss of generality let us assume that $`\mathrm{min}[Z(u),Z(v)]=Z(u)=:Z`$. We can rewrite it as $`Z=_{N=0}^{\mathrm{}}_{M=N+1}^{\mathrm{}}(c_Nc_M)^2X_{N,M}`$, where $`X_{N,M}:=|u_N|a^{}|u_M|^2+|u_M|a^{}|u_N|^2`$. Now, since $`c𝒞`$ we can write $`(c_Nc_M)^2_{R=N}^{M1}(c_Rc_{R+1})^2`$, for $`MN+1`$, so that

$`Z`$ $``$ $`{\displaystyle \underset{R=0}{\overset{\mathrm{}}{}}}(c_Rc_{R+1})^2{\displaystyle \underset{N=0}{\overset{R}{}}}{\displaystyle \underset{M=R+1}{\overset{\mathrm{}}{}}}X_{N,M},`$
$`=`$ $`{\displaystyle \underset{R=0}{\overset{\mathrm{}}{}}}(c_Rc_{R+1})^2(R+1+2Y_R),`$
where

$$Y_R:=\underset{N=0}{\overset{R}{}}\left[u_N|a_N^{}a_N|u_N\underset{0=MN}{\overset{R}{}}|u_N|a^{}|u_M|^2\right],$$

with $`a_N:=au_N|a|u_N`$. Now, using that $`u_Nu_M`$ for $`NM`$ we have $`u_N|a^{}|u_M=u_N|a_N^{}|u_M`$ which, together with $`_{0=MN}^R|u_N|a_N^{}|u_M|^2u_N|a_N^{}a_N|u_N`$, yields that $`Y_R0`$ for all $`R`$ and therefore

$$\mathrm{\Delta }(\psi )Z\underset{R=0}{\overset{\mathrm{}}{}}(c_Rc_{R+1})^2(R+1)𝜹(c).\text{}\text{ }$$

Lemma 1 indicates that for a given set of Schmidt coefficients $`c𝒞`$ EPR–correlations are maximized if the Schmidt vectors are chosen to be Fock states in the right order, i.e. $`|u_N=|v_N=|N`$. Next we will show that for fixed $`\mathrm{\Delta }`$, the choice of Schmidt coefficients minimizing the entropy of entanglement is given by those of a two–mode squeezed state. Since the entropy and the EPR entanglement are explicitly known functionals $`𝐞(c)`$ and $`𝜹(c)`$ on the sequences $`c𝒞`$, this is a classical constrained variational problem.

Lemma 2: For $`\mathrm{\Delta }(0,1)`$, and any sequence $`c𝒞`$ with $`𝜹(c)=\mathrm{\Delta }`$, we have $`𝐞(c)𝐞(c^\mathrm{\Delta })E[\mathrm{\Psi }_s(r_\mathrm{\Delta })]`$, where $`c_N^\mathrm{\Delta }\mathrm{exp}(Nr_\mathrm{\Delta })`$ is the unique geometric sequence in $`𝒞`$ with $`𝜹(c^\mathrm{\Delta })=\mathrm{\Delta }`$.

Sketch of proof: We apply the method of Lagrange multipliers for constrained minima to the infinitely many variables $`c_0,c_1\mathrm{}`$, leaving aside the technicalities of making this rigorous. These involve restricting $`c`$ to finite dimensional spaces, then letting the dimension of the space tend to infinity, and controlling the attained minima in this limit.

With a choice of Langrange multipliers $`\mu `$ and $`\lambda >0`$, designed to simplify the expressions to come, we are thus looking for stationary values $`c𝒞`$ of the functional

$$F(c,\lambda ,\mu ):=𝐞(c)+\frac{\lambda }{2\mathrm{ln}(2)}[𝜹(c)\mathrm{\Delta }]+\frac{(\mu +1)}{\mathrm{ln}(2)}(c1).$$

We obtain

$`2c_N[N\lambda +\mu \mathrm{ln}(c_N^2)]=\lambda [Nc_{N1}+(N+1)c_{N+1}].`$ (8)
where we have defined $`c_1=1`$. One can immediately see that $`c_N>0`$ and thus we can divide (8) by $`c_N`$ and subtract the same expression but for $`N+1`$. Defining $`x_N:=c_{N+1}/c_N=:e^{2r_N}(0,1]`$ for $`N=0,1,\mathrm{}`$ and writing $`\lambda =2r/\mathrm{sinh}^2(r)`$ for some $`r>0`$ we find

$`x_{N+1}=x_NA_NB_N,`$ (9)
where $`N=0,1,\mathrm{}`$ and

$`A_N`$ $`=`$ $`{\displaystyle \frac{4}{N+2}}\left[\mathrm{sinh}^2(r_N){\displaystyle \frac{r_N}{r}}\mathrm{sinh}^2(r)\right],`$ (10a)
$`B_N`$ $`=`$ $`{\displaystyle \frac{N}{N+2}}\left[{\displaystyle \frac{1}{x_N}}{\displaystyle \frac{1}{x_{N1}}}\right].`$ (10b)
If we fix $`r>0`$ and $`x_0`$, we have three possibilities: (i) $`x_0<e^{2r}`$. Then, by induction, $`x_N`$ is decreasing, and will reach some $`x_N<0`$ for finite $`N`$, which is impossible; (ii) $`x_0>e^{2r}`$. Then $`x_N`$ is increasing, and the normalization condition for $`c`$ cannot be fulfilled. Hence we must have the third possibility (iii) $`x_0=e^{2r}`$, which implies that $`x_N=e^{2r}`$ for all $`N`$. Hence $`c_N`$ is a geometric sequence $`\mathrm{exp}(2Nr)`$.

With the help of Lemmas 1 and 2 we are now in the position of proving Proposition 1.

Proof of Proposition 1: Given $`\psi `$, if $`\mathrm{\Delta }(\mathrm{\Psi })=1`$ then it is trivial. Otherwise, using (6) and Lemma 2 we have

$`E(\psi )=𝐞(c)E[\mathrm{\Psi }_s(r_{𝜹(c)})]E[\mathrm{\Psi }_s(r_{\mathrm{\Delta }(\psi )})],`$ (11)
where for the last inequality we have used $`r_{\mathrm{\Delta }(\psi )}r_{𝜹(c)}`$ \[which follows from Lemma 1 and (3)\] and the fact that $`E[\mathrm{\Psi }_s(r)]`$ increases monotonically with $`r`$.

In the following, we will apply Proposition 1 to determine the EoF of symmetric Gaussian states of two modes. For a given density operator, $`\sigma `$, we define its covariance matrix (CM) $`\gamma `$ as usual,

$`\gamma _{ij}:=\mathrm{tr}[(R_iR_j+R_jR_i)\rho ]2\mathrm{t}\mathrm{r}(R_i\rho )\mathrm{tr}(R_j\rho ),`$ (12)
where $`\{R_i,i=1,..,4\}:=\{X_A,P_A,X_B,P_B\}`$. Up to local unitary operations, it can always be written in the standard form Geza

$`\gamma =\left(\begin{array}{cccc}n& 0& k_x& 0\\ 0& n& 0& k_p\\ k_x& 0& m& 0\\ 0& k_p& 0& m\end{array}\right).`$ (13)
We will concentrate here on symmetric states, i.e. those which are invariant under exchange of subindices A and B and therefore fulfilling $`m=n`$. Without loss of generality we can choose $`k_xk_p0`$. In this case, $`\gamma `$ is a CM iff $`n^2k_x^21`$ and describes an entangled state iff $`1>(nk_x)(nk_p)`$ Geza . Next we apply local (unitary) squeezing transformations to yield the state in a more appropriate form without changing its entanglement properties. In the Heisenberg picture, the transformation multiplies (divides) $`X_{A,B}`$ ($`P_{A,B}`$) by $`[(nk_p)/(nk_x)]^{1/4}`$. A simple calculation gives

$`\mathrm{\Delta }(\sigma )=\sqrt{(nk_x)(nk_p)}=:\delta `$ (14)
where $`\mathrm{\Delta }(\sigma )`$ is defined analogously as in (2).

Our goal is to determine the EoF of $`\sigma `$. This is defined as $`E_F(\sigma ):=inf_D(D)`$, where the infimum is taken with respect to all sets of the form $`D=\{p_k,\psi _k\}`$ which give rise to a decomposition of $`\sigma `$, i.e.,

$`\sigma ={\displaystyle \underset{k}{}}p_k|\psi _k\psi _k|,`$ (15)
where the $`\psi _k`$ are normalized and $`p_k0`$. Note that the sum can run over continuous indices. For the set $`D`$ we define

$`(D):={\displaystyle \underset{k}{}}p_kE(\psi _k).`$ (16)
We call the set $`D`$ a decomposition of $`\sigma `$. A particular decomposition $`D_0`$ of $`\sigma `$ is defined through

$$\sigma _^4𝑑\xi W(\xi )|\mathrm{\Psi }_s(r_\delta )\mathrm{\Psi }_s(r_\delta )|W(\xi )^{}e^{\frac{1}{4}\xi ^T(\gamma \gamma _\delta )^1\xi },$$

where $`W(\xi )=e^{i\xi ^TR}`$ is the Weyl displacement operator and $`\gamma _\delta \gamma `$ is the CM of the two–mode squeezed state (1) with squeezing parameter $`r_\delta `$. Since $`W(\chi )`$ are local unitary operators, we have $`(D_0)=E[\mathrm{\Psi }_s(r_\delta )]`$.

We also introduce the auxiliary function $`f:(0,1][0,\mathrm{})`$

$`f(\mathrm{\Delta })=c_+(\mathrm{\Delta })\mathrm{log}[c_+(\mathrm{\Delta })]c_{}(\mathrm{\Delta })\mathrm{log}[c_{}(\mathrm{\Delta })],`$ (17)
where $`c_\pm (\mathrm{\Delta }):=(\mathrm{\Delta }^{1/2}\pm \mathrm{\Delta }^{1/2})^2/4`$. One can readily show that $`f`$ is a convex and decreasing function of $`\mathrm{\Delta }`$ and that

$`E[\mathrm{\Psi }_s(r_\mathrm{\Delta })]=f(\mathrm{\Delta }).`$ (18)
Proposition 2: $`E_F(\sigma )=f[\sqrt{(nk_x)(nk_p)}]`$.

Proof: We just have to prove that for any decomposition $`D`$, $`(D)f(\delta )`$, where $`\delta `$ is given in (14), since the decomposition $`D_0`$ already achieves this value, i.e. $`(D_0)=E[\mathrm{\Psi }_s(r_\delta )]=f(\delta )`$ \[c.f. Eq. (18)\]. For any decomposition we have

$$(D)\underset{k}{}p_kf[\mathrm{\Delta }(\psi _k)]f\left[\underset{k}{}p_k\mathrm{\Delta }(\psi _k)\right]f(\delta ).$$

The first inequality is a consequence of Proposition 1 and (18). The second is due to the convexity of $`f`$. Finally, the last one is a consequence of the fact that $`\delta p_k\mathrm{\Delta }(\psi _k)`$ (which can be easily checked by using the Cauchy–Schwarz inequality) together with the fact that $`f`$ is a decreasing function of its argument.

In summary, we have determined the EoF of symmetric Gaussian states by establishing a connection between EPR–like correlations and the entanglement of a state. The result implies that the measured quantities in some of the recent experiments dealing with atoms expatoms  and photons Glockl  not only qualify entanglement but also quantify it. We expect that the methods introduced here will allow to determine the EoF and other properties of more general Gaussian states. The optimal decomposition $`D_0`$ that gives rise to the EoF is a mixture of Gaussian pure states, which means that those states are the cheapest ones in terms of entanglement to produce symmetric Gaussian states. Thus, it is tempting to conjecture that this is also true for all Gaussian states. Finally, the results presented here provide a new characterization of two–mode squeezed states as those, which achieve a maximal EPR–like correlation for a given value of the present entanglement.

We thank Frank Verstraete for discussions. This work has been supported by the EU IST program (QUPRODIS and RESQ), and the *Kompetenznetzwerk Quanteninformationsverarbeitung der Bayerischen Staatsregierung*. GG acknowledges funding by the *Alexander-von-Humboldt–Stiftung*.